FEI  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
FEData.hpp
1 #ifndef _FEData_h_
2 #define _FEData_h_
3 
4 /*--------------------------------------------------------------------*/
5 /* Copyright 2005 Sandia Corporation. */
6 /* Under the terms of Contract DE-AC04-94AL85000, there is a */
7 /* non-exclusive license for use of this work by or on behalf */
8 /* of the U.S. Government. Export of this program may require */
9 /* a license from the United States Government. */
10 /*--------------------------------------------------------------------*/
11 
12 #include <fei_iostream.hpp>
13 #include <fei_fstream.hpp>
14 
15 #include <fei_mpi.h>
16 
17 #define dbgOut() if (debugOutputLevel_ > 0) *dbgOStreamPtr_
18 
19 #include <fei_defs.h>
20 #include <fei_FiniteElementData.hpp>
21 
22 #include <snl_fei_Utils.hpp>
23 
28 class FEData : public virtual FiniteElementData {
29  public:
31  FEData(MPI_Comm comm)
32  :
33  comm_(comm), numProcs_(1), localProc_(0),
34  debugOutputLevel_(0),
35  dbgPath_(NULL),
36  dbgOStreamPtr_(NULL),
37  dbgFileOpened_(false)
38  {
39 #ifndef FEI_SER
40  if (MPI_Comm_rank(comm_, &localProc_) != MPI_SUCCESS) MPI_Abort(comm_,-1);
41  if (MPI_Comm_size(comm_, &numProcs_) != MPI_SUCCESS) MPI_Abort(comm_,-1);
42 #endif
43  setDebugLog(0, ".");
44  }
45 
46  virtual ~FEData()
47  {
48  if (dbgFileOpened_ == true) { dbgFStreamPtr_->close(); }
49 
50  delete [] dbgPath_;
51  delete dbgOStreamPtr_;
52  }
53 
59  int parameters(int numParams, char** params);
60 
68  int setLookup(Lookup& lookup)
69  {
70  dbgOut() << "setLookup" << FEI_ENDL;
71  return(0);
72  }
73 
74 
78  int describeStructure(int numElemBlocks,
79  const int* numElemsPerBlock,
80  const int* numNodesPerElem,
81  const int* elemMatrixSizePerBlock,
82  int totalNumNodes,
83  int numSharedNodes,
84  int numMultCRs)
85  {
86  dbgOut() << "describeStructure" << FEI_ENDL
87  << " numElemBlocks: " << numElemBlocks << FEI_ENDL;
88  for(int i=0; i<numElemBlocks; i++) {
89  dbgOut() << " elem-block " << i << ": " << FEI_ENDL
90  << " number of elements: " << numElemsPerBlock[i] << FEI_ENDL
91  << " nodes per element: " << numNodesPerElem[i] << FEI_ENDL;
92  dbgOut() << " elemMatrixSizePerBlock: "
93  << elemMatrixSizePerBlock[i] << FEI_ENDL;
94  }
95  return(0);
96  }
97 
100  int setConnectivity(int elemBlockID,
101  int elemID,
102  int numNodes,
103  const int* nodeNumbers,
104  const int* numDofPerNode,
105  const int* dof_ids)
106  {
107  dbgOut() << "setConnectivity" << FEI_ENDL
108  << " elemBlockID: " << elemBlockID << ", elemID: " << elemID
109  << ", numNodes: " << numNodes << FEI_ENDL
110  << " nodeNumbers: ";
111  for(int i=0; i<numNodes; i++) {
112  dbgOut() << nodeNumbers[i] << " ";
113  }
114  dbgOut() << FEI_ENDL;
115 
116  dbgOut() << " numDOFPerNode: ";
117  for(int j=0; j<numNodes; j++) {
118  dbgOut() << numDofPerNode[j] << " ";
119  }
120  dbgOut() << FEI_ENDL;
121 
122  return(0);
123  }
124 
125 
126  int setElemMatrix(int elemBlockID,
127  int elemID,
128  int numNodes,
129  const int* nodeNumbers,
130  const int* dofPerNode,
131  const int* dof_ids,
132  const double *const * coefs)
133  {
134  dbgOut() << "setElemMatrix" << FEI_ENDL
135  << " elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
136  << " numNodes: " << numNodes << FEI_ENDL;
137  int i;
138  dbgOut() << " nodeNumbers: ";
139  for(i=0; i<numNodes; i++) {
140  dbgOut() << nodeNumbers[i] << " ";
141  }
142  dbgOut() << FEI_ENDL << " dofPerNode: ";
143  int numRows = 0;
144  for(i=0; i<numNodes; i++) {
145  dbgOut() << dofPerNode[i] << " ";
146  numRows += dofPerNode[i];
147  }
148  dbgOut() << FEI_ENDL << " coefs:" << FEI_ENDL;
149  for(i=0; i<numRows; i++) {
150  dbgOut() << " ";
151  for(int j=0; j<numRows; j++) {
152  dbgOut() << coefs[i][j] << " ";
153  }
154  dbgOut() << FEI_ENDL;
155  }
156 
157  return(0);
158  }
159 
160 
161  int setElemVector(int elemBlockID,
162  int elemID,
163  int numNodes,
164  const int* nodeNumbers,
165  const int* dofPerNode,
166  const int* dof_ids,
167  const double* coefs)
168  {
169  dbgOut() << "setElemVector" << FEI_ENDL
170  << " elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
171  << " numNodes: " << numNodes << FEI_ENDL;
172  int i;
173  dbgOut() << " nodeNumbers: ";
174  for(i=0; i<numNodes; i++) {
175  dbgOut() << nodeNumbers[i] << " ";
176  }
177  dbgOut() << FEI_ENDL << " dofPerNode: ";
178  int numRows = 0;
179  for(i=0; i<numNodes; i++) {
180  dbgOut() << dofPerNode[i] << " ";
181  numRows += dofPerNode[i];
182  }
183  dbgOut() << FEI_ENDL << " coefs:" << FEI_ENDL << " ";
184  for(i=0; i<numRows; i++) {
185  dbgOut() << coefs[i] << " ";
186  }
187  dbgOut() << FEI_ENDL;
188  return(0);
189  }
190 
191  int setDirichletBCs(int numBCs,
192  const int* nodeNumbers,
193  const int* dofOffsets,
194  const double* values)
195  {
196  dbgOut() << "setDirichletBCs" << FEI_ENDL
197  << " numBCs: " << numBCs << FEI_ENDL;
198  for(int i=0; i<numBCs; i++) {
199  dbgOut() << " nodeNumber: " << nodeNumbers[i] << ", "
200  << "dof-offset: " << dofOffsets[i] << ", value: " << values[i]<<FEI_ENDL;
201  }
202 
203  return(0);
204  }
205 
206  int sumIntoMatrix(int numRowNodes,
207  const int* rowNodeNumbers,
208  const int* rowDofOffsets,
209  const int* numColNodesPerRow,
210  const int* colNodeNumbers,
211  const int* colDofOffsets,
212  const double* coefs)
213  {
214  dbgOut() << "sumIntoMatrix, numRowNodes: " << numRowNodes << FEI_ENDL;
215  int offset = 0;
216  for(int i=0; i<numRowNodes; i++) {
217  dbgOut() << " rowNodeNumber " << rowNodeNumbers[i]
218  << ", rowDofOffset " << rowDofOffsets[i] << FEI_ENDL;
219  for(int j=0; j<numColNodesPerRow[i]; j++) {
220  dbgOut() << " colNodeNumber " << colNodeNumbers[offset]
221  << ", colDofOffset " << colDofOffsets[offset]
222  << ", value: " << coefs[offset]<<FEI_ENDL;
223  offset++;
224  }
225  }
226 
227  return(0);
228  }
229 
230  int sumIntoRHSVector(int numNodes,
231  const int* nodeNumbers,
232  const int* dofOffsets,
233  const double* coefs)
234  {
235  dbgOut() << "sumIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
236  for(int i=0; i<numNodes; i++) {
237  dbgOut() << " nodeNumber " << nodeNumbers[i]
238  << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
239  }
240 
241  return(0);
242  }
243 
244  int putIntoRHSVector(int numNodes,
245  const int* nodeNumbers,
246  const int* dofOffsets,
247  const double* coefs)
248  {
249  dbgOut() << "putIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
250  for(int i=0; i<numNodes; i++) {
251  dbgOut() << " nodeNumber " << nodeNumbers[i]
252  << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
253  }
254 
255  return(0);
256  }
257 
258  int loadComplete()
259  {
260  dbgOut() << "loadComplete" << FEI_ENDL;
261 
262  return(0);
263  }
264 
274  int launchSolver(int& solveStatus, int& iterations)
275  {
276  dbgOut() << "launchSolver" << FEI_ENDL;
277 
278  solveStatus = 0;
279  iterations = 0;
280 
281  return(0);
282  }
283 
284  int reset()
285  {
286  dbgOut() << "reset" << FEI_ENDL;
287  return(0);
288  }
289 
290  int resetRHSVector()
291  {
292  dbgOut() << "resetRHSVector" << FEI_ENDL;
293  return(0);
294  }
295 
296  int resetMatrix()
297  {
298  dbgOut() << "resetMatrix" << FEI_ENDL;
299  return(0);
300  }
301 
302  int deleteConstraints()
303  {
304  dbgOut() << "deleteConstraints" << FEI_ENDL;
305  return(0);
306  }
307 
308  int getSolnEntry(int nodeNumber,
309  int dofOffset,
310  double& value)
311  {
312  dbgOut() << "getSolnEntry, nodeNumber: " << nodeNumber
313  << ", dofOffset: " << dofOffset << FEI_ENDL;
314 
315  value = -999.99;
316 
317  return(0);
318  }
319 
320  int getMultiplierSoln(int CRID, double& lagrangeMultiplier)
321  {
322  lagrangeMultiplier = -999.99;
323  return(0);
324  }
325 
337  int putNodalFieldData(int fieldID,
338  int fieldSize,
339  int numNodes,
340  const int* nodeNumbers,
341  const double* coefs)
342  {
343  dbgOut() << "putNodalFieldData, fieldID: " << fieldID << ", fieldSize: "
344  << fieldSize << FEI_ENDL;
345  int offset = 0;
346  for(int i=0; i<numNodes; i++) {
347  dbgOut() << " nodeNumber " << nodeNumbers[i] << ", coefs: ";
348  for(int j=0; j<fieldSize; j++) {
349  dbgOut() << coefs[offset++] << " ";
350  }
351  dbgOut() << FEI_ENDL;
352  }
353 
354  return(0);
355  }
356 
357  int setMultiplierCR(int CRID,
358  int numNodes,
359  const int* nodeNumbers,
360  const int* dofOffsets,
361  const double* coefWeights,
362  double rhsValue)
363  {
364  dbgOut() << "setMultiplierCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
365  for(int i=0; i<numNodes; i++) {
366  dbgOut() << " nodeNumber " << nodeNumbers[i] << ", dof-offset "
367  << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
368  }
369 
370  dbgOut() << " rhsValue: " << rhsValue << FEI_ENDL;
371 
372  return(0);
373  }
374 
375  int setPenaltyCR(int CRID,
376  int numNodes,
377  const int* nodeNumbers,
378  const int* dofOffsets,
379  const double* coefWeights,
380  double penaltyValue,
381  double rhsValue)
382  {
383  dbgOut() << "setPenaltyCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
384  for(int i=0; i<numNodes; i++) {
385  dbgOut() << " nodeNumber " << nodeNumbers[i] << ", dof-offset "
386  << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
387  }
388 
389  dbgOut() << " penaltyValue: " << penaltyValue << FEI_ENDL;
390  dbgOut() << " rhsValue: " << rhsValue << FEI_ENDL;
391 
392  return(0);
393  }
394 
395  private:
396  int setDebugLog(int debugOutputLevel, const char* path);
397 
398  MPI_Comm comm_;
399  int numProcs_, localProc_;
400 
401  int debugOutputLevel_;
402  char* dbgPath_;
403  FEI_OSTREAM* dbgOStreamPtr_;
404  bool dbgFileOpened_;
405  FEI_OFSTREAM* dbgFStreamPtr_;
406 };
407 
408 #endif // _FEData_h_
FEData(MPI_Comm comm)
int getMultiplierSoln(int CRID, double &lagrangeMultiplier)
int setElemMatrix(int elemBlockID, int elemID, int numNodes, const int *nodeNumbers, const int *dofPerNode, const int *dof_ids, const double *const *coefs)
int launchSolver(int &solveStatus, int &iterations)
int describeStructure(int numElemBlocks, const int *numElemsPerBlock, const int *numNodesPerElem, const int *elemMatrixSizePerBlock, int totalNumNodes, int numSharedNodes, int numMultCRs)
int setMultiplierCR(int CRID, int numNodes, const int *nodeNumbers, const int *dofOffsets, const double *coefWeights, double rhsValue)
int loadComplete()
int parameters(int numParams, char **params)
Definition: FEData.cpp:20
int reset()
int setConnectivity(int elemBlockID, int elemID, int numNodes, const int *nodeNumbers, const int *numDofPerNode, const int *dof_ids)
int setPenaltyCR(int CRID, int numNodes, const int *nodeNumbers, const int *dofOffsets, const double *coefWeights, double penaltyValue, double rhsValue)
int getSolnEntry(int nodeNumber, int dofOffset, double &value)
int setDirichletBCs(int numBCs, const int *nodeNumbers, const int *dofOffsets, const double *values)
int sumIntoMatrix(int numRowNodes, const int *rowNodeNumbers, const int *rowDofOffsets, const int *numColNodesPerRow, const int *colNodeNumbers, const int *colDofOffsets, const double *coefs)
int setElemVector(int elemBlockID, int elemID, int numNodes, const int *nodeNumbers, const int *dofPerNode, const int *dof_ids, const double *coefs)
int deleteConstraints()
int setLookup(Lookup &lookup)
int putNodalFieldData(int fieldID, int fieldSize, int numNodes, const int *nodeNumbers, const double *coefs)