FEI Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
snl_fei_tester.cpp
Go to the documentation of this file.
1/*--------------------------------------------------------------------*/
2/* Copyright 2005 Sandia Corporation. */
3/* Under the terms of Contract DE-AC04-94AL85000, there is a */
4/* non-exclusive license for use of this work by or on behalf */
5/* of the U.S. Government. Export of this program may require */
6/* a license from the United States Government. */
7/*--------------------------------------------------------------------*/
8
9#include <fei_macros.hpp>
10
12
14#include <fei_ArrayUtils.hpp>
16
17#include <fei_base.hpp>
18
19#ifdef HAVE_FEI_FETI
20#include <FETI_DP_FiniteElementData.h>
21#endif
22
25
26#undef fei_file
27#define fei_file "snl_fei_tester.cpp"
28#include <fei_ErrMacros.hpp>
29
30//----------------------------------------------------------------------------
32 MPI_Comm comm, int localProc, int numProcs)
33 : comm_(comm),
34 factory_(),
35 vecSpace_(),
36 matrixGraph_(),
37 A_(),
38 x_(),
39 b_(),
40 linSys_(NULL),
41 linSysCore_(NULL),
42 feData_(NULL),
43 data_(data_reader),
44 idTypes_(),
45 numPatterns_(0),
46 localProc_(localProc),
47 numProcs_(numProcs)
48{
49}
50
51//----------------------------------------------------------------------------
53{
54 delete linSysCore_;
55 delete feData_;
56}
57
58//----------------------------------------------------------------------------
60{
61 if (factory_.get() == NULL) {
62 try {
64 }
65 catch (std::runtime_error& exc) {
66 fei::console_out() << exc.what()<<FEI_ENDL;
67 return(1);
68 }
69 if (factory_.get() == NULL) {
70 ERReturn(-1);
71 }
72 }
73
74 std::vector<std::string> stdstrings;
76 stdstrings);
77 fei::ParameterSet paramset;
78 fei::utils::parse_strings(stdstrings, " ", paramset);
79
80 if (!path_.empty()) {
81 paramset.add(fei::Param("debugOutput", path_.c_str()));
82 }
83
84 factory_->parameters(paramset);
85
87
88 vecSpace_->setParameters(paramset);
89
91
94
95 matrixGraph_->setParameters(paramset);
96
98
100
101 int i;
102 for(i=0; i<data_->numSharedNodeSets_; ++i) {
103 CommNodeSet& nodeSet = data_->sharedNodeSets_[i];
104
107 nodeSet.nodeIDs_,
108 nodeSet.procsPerNode_,
109 nodeSet.procs_) );
110 }
111
113
114 return(0);
115}
116
117//----------------------------------------------------------------------------
119{
120 FEI_OSTRINGSTREAM osstr;
121 osstr << "A_" << A_->typeName() << ".np"<<numProcs_;
122 std::string str = osstr.str();
123 A_->writeToFile(str.c_str());
124}
125
126//----------------------------------------------------------------------------
127void snl_fei_tester::setParameter(const char* param)
128{
129 std::vector<std::string> stdstrings;
130 fei::utils::char_ptrs_to_strings(1, &param, stdstrings);
131 fei::ParameterSet paramset;
132 fei::utils::parse_strings(stdstrings, " ", paramset);
133 factory_->parameters(paramset);
134 vecSpace_->setParameters(paramset);
135 matrixGraph_->setParameters(paramset);
136
137 linSys_->parameters(1, &param);
138 A_->parameters(paramset);
139}
140
141//----------------------------------------------------------------------------
143{
145
149
151
153
154 std::vector<std::string> stdstrings;
156 fei::ParameterSet paramset;
157 fei::utils::parse_strings(stdstrings, " ", paramset);
158 CHK_ERR( A_->parameters(paramset) );
159
161 linSys_->setRHS(b_);
163
164 CHK_ERR( A_->putScalar(0.0) );
165 CHK_ERR( b_->putScalar(0.0) );
166
168
170
172
173 int i;
174 for(i=0; i<data_->numBCNodeSets_; ++i) {
175 BCNodeSet& bcSet = data_->bcNodeSets_[i];
176 int fieldSize = data_->getFieldSize(bcSet.fieldID_);
177 if (fieldSize < 1) {
178 continue;
179 }
180
182 bcSet.nodeIDs_,
184 bcSet.fieldID_,
185 bcSet.offsetsIntoField_,
186 bcSet.prescribed_values_) );
187 }
188
190
191 return(0);
192}
193
194//----------------------------------------------------------------------------
196{
198
199 std::vector<std::string> stdstrings;
201 fei::ParameterSet paramset;
202 fei::utils::parse_strings(stdstrings," ",paramset);
203
204 int status, itersTaken = 0;
205 CHK_ERR( solver->solve(linSys_.get(),
206 NULL, //preconditioningMatrix
207 paramset, itersTaken, status) );
208
210
211 return(0);
212}
213
214//----------------------------------------------------------------------------
216{
218 numProcs_, localProc_, 1));
219
221 numProcs_, localProc_, 1));
222
224 numProcs_, localProc_, 1));
225
227 data_->checkFileName_.c_str(), "node", 1);
228
230 data_->checkFileName_.c_str(), "elem", 1);
231
233 data_->checkFileName_.c_str(), "mult", 1);
234 int globalErr = err;
235#ifndef FEI_SER
236 if (MPI_SUCCESS != MPI_Allreduce(&err, &globalErr, 1, MPI_INT, MPI_SUM,
237 comm_)) return(-1);
238#endif
239 if (globalErr != 0) return(-1);
240 return(0);
241}
242
243//----------------------------------------------------------------------------
245{
247
248 //nodeIDType == 0
249 idTypes_.push_back(0);
250
251 //constraintIDType == 1
252 idTypes_.push_back(1);
253
254 //elemDofIDType == 2
255 idTypes_.push_back(2);
256
258
259 nodeTypeOffset_ = 0;
261 elemTypeOffset_ = 2;
262}
263
264//----------------------------------------------------------------------------
266{
267 for(int i=0; i<data_->numElemBlocks_; ++i) {
268 ElemBlock& eb = data_->elemBlocks_[i];
269
270 int patternID;
271 definePattern(eb, patternID);
272
274 eb.numElements_,
275 patternID) );
276
277 for(int j=0; j<eb.numElements_; ++j) {
278 std::vector<int> conn(eb.numNodesPerElement_);
279 for(int ii=0; ii<eb.numNodesPerElement_; ++ii) {
280 conn[ii] = eb.elemConn_[j][ii];
281 }
282
284 eb.elemIDs_[j],
285 &conn[0]) );
286 }
287 }
288
289 return(0);
290}
291
292//----------------------------------------------------------------------------
294{
295 int i;
296 for(i=0; i<data_->numElemBlocks_; ++i) {
297 ElemBlock& eb = data_->elemBlocks_[i];
298
299 if (eb.numElements_ < 1) {
300 continue;
301 }
302
303 int numIndices = matrixGraph_->getConnectivityNumIndices(eb.blockID_);
304
305 std::vector<int> indices(numIndices);
306
307 for(int j=0; j<eb.numElements_; ++j) {
308 int checkNum;
310 eb.elemIDs_[j],
311 numIndices,
312 &indices[0],
313 checkNum) );
314 if (numIndices != checkNum) {
315 ERReturn(-1);
316 }
317
318 CHK_ERR( A_->sumIn(eb.blockID_, eb.elemIDs_[j],
319 eb.elemStiff_[j]) );
320
321 CHK_ERR( b_->sumIn(numIndices, &indices[0],
322 eb.elemLoad_[j], 0) );
323 }
324 }
325
326 return(0);
327}
328
329//----------------------------------------------------------------------------
331{
332 std::vector<int> idTypes;
333 int constraintID = localProc_*100000;
334 int i;
335 for(i=0; i<data_->numCRMultSets_; ++i) {
336 CRSet& crSet = data_->crMultSets_[i];
337
338 for(int j=0; j<1; ++j) {
339 idTypes.assign(crSet.numNodes_, idTypes_[nodeTypeOffset_]);
340
341 crSet.crID_ = constraintID++;
342 int constraintIDType = idTypes_[constraintTypeOffset_];
344 constraintIDType,
345 crSet.numNodes_,
346 &idTypes[0],
347 crSet.nodeIDs_[j],
348 crSet.fieldIDs_) );
349 }
350 }
351
352 for(i=0; i<data_->numCRPenSets_; ++i) {
353 CRSet& crSet = data_->crPenSets_[i];
354
355 for(int j=0; j<1; ++j) {
356 idTypes.assign(crSet.numNodes_, idTypes_[nodeTypeOffset_]);
357
358 crSet.crID_ = constraintID++;
359 int constraintIDType = idTypes_[constraintTypeOffset_];
361 constraintIDType,
362 crSet.numNodes_,
363 &idTypes[0],
364 crSet.nodeIDs_[j],
365 crSet.fieldIDs_) );
366 }
367 }
368
369 std::map<int,int> fieldDB;
370 for(i=0; i<data_->numFields_; ++i) {
371 fieldDB.insert(std::pair<int,int>(data_->fieldIDs_[i], data_->fieldSizes_[i]));
372 }
373
374 std::vector<int> nodeIDs;
375 std::vector<int> fieldIDs;
376 std::vector<double> weights;
377
378 for(i=0; i<data_->numSlaveVars_; i++) {
379 int ii;
380 CRSet& crSet = data_->slaveVars_[i];
381
382 nodeIDs.resize(crSet.numNodes_+1);
383 nodeIDs[0] = crSet.slaveNodeID_;
384 fieldIDs.resize(0);
385 fieldIDs.push_back(crSet.slaveFieldID_);
386
387 for(ii=0; ii<crSet.numNodes_; ++ii) {
388 nodeIDs[ii+1] = crSet.nodeIDs_[0][ii];
389 fieldIDs.push_back(crSet.fieldIDs_[ii]);
390 }
391
392 idTypes.assign(crSet.numNodes_+1, idTypes_[nodeTypeOffset_]);
393
394 int fieldSize = fieldDB[crSet.slaveFieldID_];
395 weights.resize(0);
396 for(ii=0; ii<fieldSize; ++ii) weights.push_back(0.0);
397 weights[crSet.slaveOffset_] = -1.0;
398 int offset = 0;
399 for(ii=0; ii<crSet.numNodes_; ++ii) {
400 fieldSize = fieldDB[crSet.fieldIDs_[ii]];
401 for(int jj=0; jj<fieldSize; ++jj) {
402 weights.push_back(crSet.weights_[offset++]);
403 }
404 }
405
407 &idTypes[0],
408 &nodeIDs[0],
409 &fieldIDs[0],
410 0,
411 crSet.slaveOffset_,
412 &weights[0],
413 crSet.values_[0]));
414 }
415
416 return(0);
417}
418
419//----------------------------------------------------------------------------
421{
422 int i;
423 for(i=0; i<data_->numCRMultSets_; ++i) {
424 CRSet& crSet = data_->crMultSets_[i];
425
426 for(int j=0; j<1; ++j) {
428 crSet.weights_,
429 crSet.values_[j]) );
430 }
431 }
432
433 for(i=0; i<data_->numCRPenSets_; ++i) {
434 CRSet& crSet = data_->crPenSets_[i];
435
436 for(int j=0; j<1; ++j) {
438 crSet.weights_,
439 crSet.penValues_[j],
440 crSet.values_[j]) );
441 }
442 }
443
444 return(0);
445}
446
447//----------------------------------------------------------------------------
449{
450 int i, j, numIDTypes = 1;
451 numIDTypes += eb.numElemDOF_>0 ? 1 : 0;
452
453 //find out how many nodal fields there are, total.
454 std::vector<int> nodalFieldIDs;
455 std::vector<int> flatFieldIDsArray;
456 for(i=0; i<eb.numNodesPerElement_; ++i) {
457 for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
458 fei::sortedListInsert(eb.nodalFieldIDs_[i][j], nodalFieldIDs);
459 flatFieldIDsArray.push_back(eb.nodalFieldIDs_[i][j]);
460 }
461 }
462
463 patternID = numPatterns_++;
464
465 if (numIDTypes == 1 && nodalFieldIDs.size() == 1) {
466 //This is a very simple pattern
469 nodalFieldIDs[0]);
470 }
471 else if (numIDTypes == 1) {
472 std::vector<int> numFieldsPerID(eb.numNodesPerElement_);
473
477 &flatFieldIDsArray[0]);
478 }
479 else {
480 std::vector<int> idTypes(eb.numNodesPerElement_+1, idTypes_[nodeTypeOffset_]);
481 idTypes[idTypes.size()-1] = idTypes_[elemTypeOffset_];
482 std::vector<int> numFieldsPerID(idTypes.size());
483 std::vector<int> fieldIDs;
484 for(i=0; i<eb.numNodesPerElement_; ++i) {
485 numFieldsPerID[i] = eb.numFieldsPerNode_[i];
486 for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
487 fieldIDs.push_back(eb.nodalFieldIDs_[i][j]);
488 }
489 }
490 numFieldsPerID[idTypes.size()-1] = eb.numElemDOF_;
491 for(i=0; i<eb.numElemDOF_; ++i) {
492 fieldIDs.push_back(eb.elemDOFFieldIDs_[i]);
493 }
494
495 patternID = matrixGraph_->definePattern(idTypes.size(),
496 &idTypes[0],
497 &numFieldsPerID[0],
498 &fieldIDs[0]);
499 }
500}
501
502//----------------------------------------------------------------------------
504 const char* solnFileName, int numProcs,
505 int localProc, int solveCounter)
506{
507 (void)solveCounter;
508
510
511 int* nodeList = new int[numLocalNodes];
512
513 int checkNum = 0;
515 numLocalNodes, nodeList, checkNum);
516 if (err != 0) {
517 ERReturn(-1);
518 }
519
520 FEI_OSTRINGSTREAM fileName;
521 fileName<< solnFileName<<".node."<<solveCounter<<"."<<numProcs<<"."<<localProc;
522 std::string str = fileName.str();
523 FEI_OFSTREAM outfile(str.c_str());
524
525 if (!outfile || outfile.bad()) {
526 fei::console_out() << "ERROR opening solution output file " << fileName.str() << FEI_ENDL;
527 return(-1);
528 }
529
530 outfile.setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
531
532 std::vector<double> solnData;
533 std::vector<int> fieldList;
534
535 int totalSize = 0;
536
537 for(int i=0; i<numLocalNodes; i++) {
538 int idType = idTypes_[nodeTypeOffset_];
539 int ID = nodeList[i];
540
541 int numDOF = vecSpace_->getNumDegreesOfFreedom(idType, ID);
542 int numFields = vecSpace_->getNumFields(idType, ID);
543 solnData.resize(numDOF);
544 vecSpace_->getFields(idType, ID, fieldList);
545
546 outfile << ID << " " << numDOF << FEI_ENDL;
547 for(int j=0; j<numFields; ++j) {
548 int fieldSize = vecSpace_->getFieldSize(fieldList[j]);
549 totalSize += fieldSize;
550
551 CHK_ERR( vec->copyOutFieldData(fieldList[j], idType,
552 1, &ID, &solnData[0]) );
553
554 for(int k=0; k<fieldSize; ++k) {
555 outfile << solnData[k] << " ";
556 }
557 }
558 outfile << FEI_ENDL;
559 }
560
561 FEI_COUT << "save-node-soln: wrote " << totalSize << " entries for " << numLocalNodes << " nodes to " << str << FEI_ENDL;
562
563 delete [] nodeList;
564
565 outfile.close();
566 return(0);
567}
568
569//----------------------------------------------------------------------------
571 const char* solnFileName,
572 int numProcs,
573 int localProc, int solveCounter)
574{
575 (void)solveCounter;
576
578
579 int* elemList = new int[numLocalElems];
580
581 int checkNum = 0;
583 numLocalElems, elemList, checkNum);
584 if (err != 0) {
585 ERReturn(-1);
586 }
587
588 FEI_OSTRINGSTREAM fileName;
589 fileName<< solnFileName<<".elem."<<solveCounter<<"."<<numProcs<<"."<<localProc;
590 std::string str = fileName.str();
591 FEI_OFSTREAM outfile(str.c_str());
592
593 if (!outfile || outfile.bad()) {
594 fei::console_out() << "ERROR opening solution output file " << fileName.str() << FEI_ENDL;
595 return(-1);
596 }
597
598 std::vector<double> solnData;
599 std::vector<int> fieldList;
600
601 for(int i=0; i<numLocalElems; i++) {
602 int idType = idTypes_[elemTypeOffset_];
603 int ID = elemList[i];
604
605 int numDOF = vecSpace_->getNumDegreesOfFreedom(idType, ID);
606 int numFields = vecSpace_->getNumFields(idType, ID);
607 solnData.resize(numDOF);
608 vecSpace_->getFields(idType, ID, fieldList);
609
610 outfile << ID << " " << numDOF << FEI_ENDL;
611 for(int j=0; j<numFields; ++j) {
612 int fieldSize = vecSpace_->getFieldSize(fieldList[j]);
613
614 CHK_ERR( vec->copyOutFieldData(fieldList[j], idType,
615 1, &ID, &solnData[0]) );
616
617 for(int k=0; k<fieldSize; ++k) {
618 outfile << solnData[k] << " ";
619 }
620 }
621 outfile << FEI_ENDL;
622 }
623
624 delete [] elemList;
625
626 outfile.close();
627 return(0);
628}
629
630//----------------------------------------------------------------------------
632 const char* solnFileName,
633 int numProcs, int localProc,
634 int solveCounter)
635{
636 (void)solveCounter;
637
639
640 int* globalNumCRs = new int[numProcs];
641#ifndef FEI_SER
642 if (MPI_Allgather(&numLocalCRs, 1, MPI_INT, globalNumCRs, 1, MPI_INT,
643 comm_) != MPI_SUCCESS) {
644 ERReturn(-1);
645 }
646#endif
647
648 int localCRStart = 0;
649#ifndef FEI_SER
650 for(int p=0; p<localProc; p++) localCRStart += globalNumCRs[p];
651#endif
652
653 delete [] globalNumCRs;
654
655 std::vector<int> crList(numLocalCRs);
656
657 int checkNum = 0;
659 idTypes_[constraintTypeOffset_], numLocalCRs,
660 numLocalCRs ? &crList[0] : 0, checkNum);
661 if (err != 0) {
662 ERReturn(-1);
663 }
664
665 FEI_OSTRINGSTREAM fileName;
666 fileName<< solnFileName<<".mult."<<solveCounter<<"."<<numProcs<<"."<<localProc;
667 std::string str = fileName.str();
668 FEI_OFSTREAM outfile(str.c_str());
669
670 if (!outfile || outfile.bad()) {
671 fei::console_out() << "ERROR opening solution output file " << fileName.str() << FEI_ENDL;
672 return(-1);
673 }
674
675 std::vector<double> solnData;
676 std::vector<int> fieldList;
677
678 for(int i=0; i<numLocalCRs; i++) {
679 int idType = idTypes_[constraintTypeOffset_];
680 int ID = crList[i];
681
682 solnData.resize(1);
683
684 outfile << localCRStart++ << " " << 1 << FEI_ENDL;
685 for(int j=0; j<1; ++j) {
686 int globalIndex = -1;
687 CHK_ERR( vecSpace_->getGlobalIndex(idType, ID, globalIndex) );
688
689 CHK_ERR( vec->copyOut(1, &globalIndex, &solnData[0]) );
690
691 for(int k=0; k<1; ++k) {
692 outfile << solnData[k] << " ";
693 }
694 }
695 outfile << FEI_ENDL;
696 }
697
698 outfile.close();
699 return(0);
700}
int numNodes_
Definition BCNodeSet.hpp:21
int * offsetsIntoField_
Definition BCNodeSet.hpp:24
GlobalID * nodeIDs_
Definition BCNodeSet.hpp:22
double * prescribed_values_
Definition BCNodeSet.hpp:25
int * fieldIDs_
Definition CRSet.hpp:57
int slaveOffset_
Definition CRSet.hpp:53
int crID_
Definition CRSet.hpp:37
double * weights_
Definition CRSet.hpp:59
int slaveFieldID_
Definition CRSet.hpp:48
GlobalID slaveNodeID_
Definition CRSet.hpp:45
GlobalID ** nodeIDs_
Definition CRSet.hpp:55
double * values_
Definition CRSet.hpp:60
double * penValues_
Definition CRSet.hpp:61
int numNodes_
Definition CRSet.hpp:40
GlobalID * nodeIDs_
int * procsPerNode_
char ** paramStrings_
BCNodeSet * bcNodeSets_
CRSet * crMultSets_
CRSet * slaveVars_
int * fieldIDs_
std::string solnFileName_
ElemBlock * elemBlocks_
CommNodeSet * sharedNodeSets_
int numCRPenSets_
int * fieldSizes_
int numSharedNodeSets_
int numElemBlocks_
int getFieldSize(int fieldID)
int numCRMultSets_
std::string solverLibraryName_
int numBCNodeSets_
std::string checkFileName_
int numSlaveVars_
CRSet * crPenSets_
int * numFieldsPerNode_
Definition ElemBlock.hpp:20
GlobalID blockID_
Definition ElemBlock.hpp:17
double *** elemStiff_
Definition ElemBlock.hpp:26
int numNodesPerElement_
Definition ElemBlock.hpp:19
GlobalID ** elemConn_
Definition ElemBlock.hpp:23
GlobalID * elemIDs_
Definition ElemBlock.hpp:22
double ** elemLoad_
Definition ElemBlock.hpp:27
int numElements_
Definition ElemBlock.hpp:18
int numElemDOF_
Definition ElemBlock.hpp:28
int ** nodalFieldIDs_
Definition ElemBlock.hpp:21
int * elemDOFFieldIDs_
Definition ElemBlock.hpp:29
virtual void parameters(const fei::ParameterSet &paramset)
virtual fei::SharedPtr< fei::LinearSystem > createLinearSystem(fei::SharedPtr< fei::MatrixGraph > &matrixGraph)
virtual void setRHS(fei::SharedPtr< fei::Vector > &rhs)
virtual int loadEssentialBCs(int numIDs, const int *IDs, int idType, int fieldID, int offsetIntoField, const double *prescribedValues)
virtual void setSolutionVector(fei::SharedPtr< fei::Vector > &soln)
virtual int loadPenaltyConstraint(int constraintID, const double *weights, double penaltyValue, double rhsValue)=0
virtual int loadLagrangeConstraint(int constraintID, const double *weights, double rhsValue)=0
virtual int parameters(int numParams, const char *const *paramStrings)=0
virtual void setMatrix(fei::SharedPtr< fei::Matrix > &matrix)
virtual int loadComplete(bool applyBCs=true, bool globalAssemble=true)=0
virtual fei::SharedPtr< fei::MatrixGraph > createMatrixGraph(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > columnSpace, const char *name)=0
virtual int definePattern(int numIDs, int idType)=0
virtual int initPenaltyConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)=0
virtual int initComplete()=0
virtual void setIndicesMode(int mode)=0
virtual int createSlaveMatrices()=0
virtual int initConnectivity(int blockID, int connectivityID, const int *connectedIdentifiers)=0
virtual int getConnectivityNumIndices(int blockID) const =0
virtual int initSlaveConstraint(int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs, int offsetOfSlave, int offsetIntoSlaveField, const double *weights, double rhsValue)=0
virtual int initLagrangeConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)=0
virtual int getConnectivityIndices(int blockID, int connectivityID, int indicesAllocLen, int *indices, int &numIndices)=0
virtual void setParameters(const fei::ParameterSet &params)=0
virtual int initConnectivityBlock(int blockID, int numConnectivityLists, int patternID, bool diagonal=false)=0
virtual fei::SharedPtr< fei::Matrix > createMatrix(fei::SharedPtr< fei::MatrixGraph > matrixGraph)=0
virtual int writeToFile(const char *filename, bool matrixMarketFormat=true)=0
virtual const char * typeName()=0
virtual int parameters(const fei::ParameterSet &paramset)=0
virtual int sumIn(int numRows, const int *rows, int numCols, const int *cols, const double *const *values, int format=0)=0
virtual int putScalar(double scalar)=0
void add(const Param &param, bool maintain_unique_keys=true)
virtual fei::SharedPtr< fei::Solver > createSolver(const char *name=0)=0
virtual fei::SharedPtr< VectorSpace > createVectorSpace(MPI_Comm, const char *name)
unsigned getFieldSize(int fieldID)
int getNumOwnedAndSharedIDs(int idType)
int initSharedIDs(int numShared, int idType, const int *sharedIDs, const int *numSharingProcsPerID, const int *sharingProcs)
void setParameters(const fei::ParameterSet &paramset)
void getFields(std::vector< int > &fieldIDs)
int getOwnedAndSharedIDs(int idtype, int lenList, int *IDs, int &numOwnedAndSharedIDs)
int getNumDegreesOfFreedom(int idType, int ID)
int getGlobalIndex(int idType, int ID, int fieldID, int fieldOffset, int whichComponentOfField, int &globalIndex)
void defineFields(int numFields, const int *fieldIDs, const int *fieldSizes, const int *fieldTypes=NULL)
void defineIDTypes(int numIDTypes, const int *idTypes)
virtual fei::SharedPtr< fei::Vector > createVector(fei::SharedPtr< fei::VectorSpace > vecSpace, int numVectors=1)=0
virtual int copyOut(int numValues, const int *indices, double *values, int vectorIndex=0) const =0
virtual int putScalar(double scalar)=0
virtual int scatterToOverlap()=0
virtual int sumIn(int numValues, const int *indices, const double *values, int vectorIndex=0)=0
virtual int copyOutFieldData(int fieldID, int idType, int numIDs, const int *IDs, double *data, int vectorIndex=0)=0
std::string path_
Definition feitester.hpp:60
int save_block_node_soln(DataReader &data, fei::Vector *vec, const char *solnFileName, int numProcs, int localProc, int solveCounter)
fei::SharedPtr< fei::VectorSpace > vecSpace_
void setParameter(const char *param)
fei::SharedPtr< DataReader > data_
fei::SharedPtr< fei::Factory > factory_
fei::SharedPtr< fei::LinearSystem > linSys_
std::vector< int > idTypes_
LinearSystemCore * linSysCore_
fei::SharedPtr< fei::Vector > x_
void definePattern(ElemBlock &eb, int &patternID)
int save_multiplier_soln(DataReader &data, fei::Vector *vec, const char *solnFileName, int numProcs, int localProc, int solveCounter)
snl_fei_tester(fei::SharedPtr< DataReader > data_reader, MPI_Comm comm, int localProc, int numProcs)
fei::SharedPtr< fei::Vector > b_
int save_block_elem_soln(DataReader &data, fei::Vector *vec, const char *solnFileName, int numProcs, int localProc, int solveCounter)
FiniteElementData * feData_
void defineFieldsAndIDTypes()
fei::SharedPtr< fei::MatrixGraph > matrixGraph_
fei::SharedPtr< fei::Matrix > A_
#define ERReturn(a)
#define CHK_ERR(a)
#define FEI_OFSTREAM
#define FEI_ENDL
#define FEI_COUT
#define IOS_FLOATFIELD
#define IOS_SCIENTIFIC
#define MPI_SUCCESS
Definition fei_mpi.h:62
#define MPI_Comm
Definition fei_mpi.h:56
#define FEI_OSTRINGSTREAM
int checkSolution(int localProc, int numProcs, const char *solnFileName, const char *checkFileName, const char *extension, int solveCounter)
Definition SolnCheck.cpp:63
void parse_strings(std::vector< std::string > &stdstrings, const char *separator_string, fei::ParameterSet &paramset)
void char_ptrs_to_strings(int numStrings, const char *const *charstrings, std::vector< std::string > &stdstrings)
fei::SharedPtr< fei::Factory > create_fei_Factory(MPI_Comm comm, const char *libraryName)
std::ostream & console_out()
int sortedListInsert(const T &item, std::vector< T > &list)