FEI Version of the Day
Loading...
Searching...
No Matches
fei_Trilinos_Helpers.hpp
1/*
2// @HEADER
3// ************************************************************************
4// FEI: Finite Element Interface to Linear Solvers
5// Copyright (2005) Sandia Corporation.
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
8// U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Alan Williams (william@sandia.gov)
38//
39// ************************************************************************
40// @HEADER
41*/
42
43
44#ifndef _fei_Trilinos_Helpers_hpp_
45#define _fei_Trilinos_Helpers_hpp_
46
47#include "fei_trilinos_macros.hpp"
48#include "fei_fwd.hpp"
49
50#include <fei_Include_Trilinos.hpp>
51
52#include <fei_mpi.h>
53#include <fei_SharedPtr.hpp>
54
55#include <fei_LinearProblemManager.hpp>
56#include <fei_VectorSpace.hpp>
57#include <fei_Reducer.hpp>
58#include <fei_MatrixGraph.hpp>
59
60namespace Trilinos_Helpers {
61
62#ifdef HAVE_FEI_EPETRA
63
69 Epetra_Map create_Epetra_Map(MPI_Comm comm,
70 const std::vector<int>& local_eqns);
71
78 create_Epetra_BlockMap(const fei::SharedPtr<fei::VectorSpace>& vecspace);
79
81 create_Epetra_CrsGraph(const fei::SharedPtr<fei::MatrixGraph>& matgraph,
82 bool blockEntries,
83 bool orderRowsWithLocalColsFirst=false);
84
86 create_from_Epetra_Matrix(fei::SharedPtr<fei::MatrixGraph> matrixGraph,
87 bool blockEntryMatrix,
89 bool orderRowsWithLocalColsFirst=false);
90
92 create_from_LPM_EpetraBasic(fei::SharedPtr<fei::MatrixGraph> matrixGraph,
93 bool blockEntryMatrix,
96 lpm_epetrabasic);
97#endif
98
102 void copy_parameterset(const fei::ParameterSet& paramset,
103 Teuchos::ParameterList& paramlist);
104
108 void copy_parameterlist(const Teuchos::ParameterList& paramlist,
109 fei::ParameterSet& paramset);
110
111#ifdef HAVE_FEI_EPETRA
116 get_Epetra_MultiVector(fei::Vector* feivec, bool soln_vec);
117
121 Epetra_VbrMatrix* get_Epetra_VbrMatrix(fei::Matrix* feimat);
122
126 Epetra_CrsMatrix* get_Epetra_CrsMatrix(fei::Matrix* feimat);
127
132 void get_Epetra_pointers(fei::SharedPtr<fei::Matrix> feiA,
135 Epetra_CrsMatrix*& crsA,
136 Epetra_Operator*& opA,
139
141 int zero_Epetra_VbrMatrix(Epetra_VbrMatrix* mat);
142
146 inline
147 double* getBeginPointer(const fei::SharedPtr<fei::Matrix>& feiA)
148 {
149 Epetra_CrsMatrix* A = get_Epetra_CrsMatrix(feiA.get());
150 return (*A)[0];
151 }
152
157 inline
158 int getOffsetG(const fei::SharedPtr<fei::Matrix>& feiA,
159 int global_row_index, int global_col_index)
160 {
161 Epetra_CrsMatrix* A = get_Epetra_CrsMatrix(feiA.get());
162 const Epetra_Map& erowmap = A->RowMap();
163 const Epetra_Map& ecolmap = A->ColMap();
164 int local_row = erowmap.LID(global_row_index);
165 int local_col = ecolmap.LID(global_col_index);
166
167 int* rowOffsets;
168 int* colIndices;
169 double* coefs;
170 A->ExtractCrsDataPointers(rowOffsets, colIndices, coefs);
171
172 int* row_ptr = &colIndices[rowOffsets[local_row]];
173 int* end_row = &colIndices[rowOffsets[local_row+1]];
174
175 int col_offset = 0;
176 for(; row_ptr != end_row; ++row_ptr) {
177 if (*row_ptr == local_col) break;
178 ++col_offset;
179 }
180
181 return rowOffsets[local_row] + col_offset;
182 }
183
188 inline
189 int getOffsetL(const fei::SharedPtr<fei::Matrix>& feiA,
190 int local_row_index, int local_col_index)
191 {
192 Epetra_CrsMatrix* A = get_Epetra_CrsMatrix(feiA.get());
193
194 int* rowOffsets;
195 int* colIndices;
196 double* coefs;
197 A->ExtractCrsDataPointers(rowOffsets, colIndices, coefs);
198
199 int* row_ptr = &colIndices[rowOffsets[local_row_index]];
200 int* end_row = &colIndices[rowOffsets[local_row_index+1]];
201
202 int col_offset = 0;
203 for(; row_ptr != end_row; ++row_ptr) {
204 if (*row_ptr == local_col_index) break;
205 ++col_offset;
206 }
207
208 return rowOffsets[local_row_index] + col_offset;
209 }
210
211#endif // HAVE_FEI_EPETRA
212
213}//namespace Trilinos_Helpers
214
215#endif // _Trilinos_Helpers_hpp_
216
int LID(int GID) const
const Epetra_Map & RowMap() const
int ExtractCrsDataPointers(int *&IndexOffset, int *&Indices, double *&Values_in) const
const Epetra_Map & ColMap() const
T * get() const