FEI Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
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_sstream.hpp>
10#include <fei_fstream.hpp>
11
13
15
18#include <snl_fei_Utils.hpp>
19
20#include <fei_FEI_Impl.hpp>
21
23
24#ifdef HAVE_FEI_FETI
25#include <FETI_DP_FiniteElementData.h>
26#endif
27
30
31#undef fei_file
32#define fei_file "FEI_tester.cpp"
33
34#include <fei_ErrMacros.hpp>
35
36//----------------------------------------------------------------------------
38 MPI_Comm comm, int localProc, int numProcs, bool useNewFEI)
39 : comm_(comm),
40 fei_(),
41 wrapper_(),
42 data_(data_reader),
43 localProc_(localProc),
44 numProcs_(numProcs),
45 numMatrices(1),
46 matrixIDs(NULL),
47 numRHSs(1),
48 rhsIDs(NULL),
49 useNewFEI_(useNewFEI)
50{
51}
52
53//----------------------------------------------------------------------------
55{
56 delete [] matrixIDs;
57 delete [] rhsIDs;
58}
59
60//----------------------------------------------------------------------------
62{
63 if (data_.get() == NULL) {
64 ERReturn(-1);
65 }
66
68
69 const char* feiVersionString;
70 CHK_ERR( fei_->version(feiVersionString) );
71
72 FEI_COUT << "FEI version: " << feiVersionString << FEI_ENDL;
73
75
76 //Now we check the solveType. A regular Ax=b solve corresponds to
77 //solveType==FEI_SINGLE_SYSTEM. The aggregate stuff (solveType==
78 //FEI_AGGREGATE_SUM) is for power users who
79 //want to assemble more than one matrix system and solve a linear
80 //combination of them, an aggregate system.
81
84 }
85
87
88 int numBlkActNodes;
89 for(int i=0; i<data_->numElemBlocks_; ++i) {
90 ElemBlock& eblk = data_->elemBlocks_[i];
91 int elemBlockID = eblk.blockID_;
92 CHK_ERR( fei_->getNumBlockActNodes(elemBlockID, numBlkActNodes) );
93 }
94
95 //************** Initialization Phase is now complete *****************
96
97 return(0);
98}
99
100//----------------------------------------------------------------------------
102{
104 CHK_ERR(fei_->resetMatrix()); //try out all of these reset functions,
105 CHK_ERR(fei_->resetSystem()); //just for the sake of coverage...
106
107 // ************ Load Phase (delegated to a function) ***************
108
109 //obviously only one of these load-phases should be
110 //performed.
111
114 }
117 }
118
120
122
123 //**** let's try out the 'put' functions. ******************
125
126 return(0);
127}
128
129//----------------------------------------------------------------------------
131{
132 //If solverLibraryName is TEST_LSC, then we're not running a real solver so
133 //just return.
134 std::string sname(data_->solverLibraryName_);
135 if (sname == "TEST_LSC") {
136 return( 0 );
137 }
138
139 int status;
140 int err = fei_->solve(status);
141
142 //FEI_COUT << "fei->solve, err = " << err << ", status = " << status << FEI_ENDL;
143
144 if (err != 0 || status != 0) {
145 FEI_COUT << "!!!! solve returned: err: "<<err<<", status: "<<status<<FEI_ENDL;
146 return(err);
147 }
148
149 if (localProc_ == 0) {
150 int itersTaken;
151 CHK_ERR( fei_->iterations(itersTaken));
152 //FEI_COUT << "iterations: " << itersTaken << FEI_ENDL;
153 }
154
156
157 return(0);
158}
159
160//----------------------------------------------------------------------------
162{
163}
164
165//----------------------------------------------------------------------------
167{
168}
169
170//----------------------------------------------------------------------------
172{
173 //If solverLibraryName is TEST_LSC, then we're not running a real solver so
174 //just check the matrix.
175 std::string sname(data_->solverLibraryName_);
176 if (sname == "TEST_LSC") {
177 return( lsc_matrix_check() );
178 }
179
181 numProcs_, localProc_, 1));
182
184 numProcs_, localProc_, 1));
185
187 numProcs_, localProc_, 1));
188
190 data_->checkFileName_.c_str(), "node", 1);
191
193 data_->checkFileName_.c_str(), "elem", 1);
194
196 data_->checkFileName_.c_str(), "mult", 1);
197 int globalErr = err;
198#ifndef FEI_SER
199 if (MPI_SUCCESS != MPI_Allreduce(&err, &globalErr, 1, MPI_INT, MPI_SUM,
200 comm_)) return(-1);
201#endif
202 if (globalErr != 0) {
203 ERReturn(-1);
204 }
205
206 return(0);
207}
208
209//----------------------------------------------------------------------------
210int FEI_tester::createFEIinstance(const char* solverName)
211{
212 try {
214 }
215 catch(std::runtime_error& exc) {
216 fei::console_out() << exc.what()<<FEI_ENDL;
217 return(1);
218 }
219
220 if (wrapper_.get() == NULL) ERReturn(-1);
221
222 if (useNewFEI_) {
224 }
225 else {
227 }
228
229 if (fei_.get() == NULL) ERReturn(-1);
230
231 return(0);
232}
233
234//----------------------------------------------------------------------------
236{
237 snl_fei::getIntParamValue("numMatrices",
241
242 matrixIDs = new int[numMatrices];
243 numRHSs = 1;
244 rhsIDs = new int[1];
245 rhsIDs[0] = 0;
246
247 for(int i=0; i<numMatrices; i++) {
248 matrixIDs[i] = i;
249 }
250
252 return(0);
253}
254
255//----------------------------------------------------------------------------
257{
260 FEI_COUT << "FEI_tester: bad solveType: " << data_->solveType_ << FEI_ENDL;
261 return(-1);
262 }
263
265
267
268 int i;
269 for(i=0; i<data_->numElemBlocks_; i++) {
270 ElemBlock& block = data_->elemBlocks_[i];
271
273 block.numElements_,
275 block.numFieldsPerNode_,
276 block.nodalFieldIDs_,
277 block.numElemDOF_,
278 block.elemDOFFieldIDs_,
279 block.interleaveStrategy_) );
280
281 for(int el=0; el<block.numElements_; el++) {
283 block.elemIDs_[el],
284 block.elemConn_[el]));
285 }
286 }
287
288 for(i=0; i<data_->numSharedNodeSets_; i++) {
289 CommNodeSet& shNodeSet = data_->sharedNodeSets_[i];
290
291 CHK_ERR(fei_->initSharedNodes(shNodeSet.numNodes_, shNodeSet.nodeIDs_,
292 shNodeSet.procsPerNode_, shNodeSet.procs_));
293 }
294
295 //********* Initialize any slave variables *****************************
296
297 for(i=0; i<data_->numSlaveVars_; i++) {
298 CRSet& crSet = data_->slaveVars_[i];
299
301 crSet.slaveFieldID_,
302 crSet.slaveOffset_,
303 crSet.numNodes_,
304 crSet.nodeIDs_[0],
305 crSet.fieldIDs_,
306 crSet.weights_,
307 crSet.values_[0]) );
308 }
309
310 //*********** Initialize Constraint Relation Equations *****************
311
312 for(i=0; i<data_->numCRMultSets_; i++) {
313 CRSet& crSet = data_->crMultSets_[i];
314
315 for(int j=0; j<1; j++) {
317 crSet.nodeIDs_[j],
318 crSet.fieldIDs_,
319 crSet.crID_));
320 }
321 }
322
323 for(i=0; i<data_->numCRPenSets_; i++) {
324 CRSet& crSet = data_->crPenSets_[i];
325
326 for(int j=0; j<1; j++) {
328 crSet.nodeIDs_[j],
329 crSet.fieldIDs_,
330 crSet.crID_));
331 }
332 }
333
335 return(0);
336}
337
338//----------------------------------------------------------------------------
340{
341 int i;
342
343 //*************** Boundary Condition Node Sets *************************
344
345 for(i=0; i<data_->numBCNodeSets_; i++) {
346 BCNodeSet& bcSet = data_->bcNodeSets_[i];
347
349 bcSet.nodeIDs_,
350 bcSet.fieldID_,
351 bcSet.offsetsIntoField_,
352 bcSet.prescribed_values_));
353 }
354
355 for(i=0; i<data_->numElemBlocks_; i++) {
356 ElemBlock& block = data_->elemBlocks_[i];
357
358 for(int el=0; el<block.numElements_; el++) {
359
361 block.elemIDs_[el],
362 block.elemConn_[el],
363 block.elemStiff_[el],
364 block.elemFormat_));
365 }
366 }
367
368 for(i=0; i<data_->numElemBlocks_; i++) {
369 ElemBlock& block = data_->elemBlocks_[i];
370
371 for(int el=0; el<block.numElements_; el++) {
372
374 block.elemIDs_[el],
375 block.elemConn_[el],
376 block.elemLoad_[el]));
377 }
378 }
379
380 //******** Load Constraint Relation Equations ***********************
381
382 for(i=0; i<data_->numCRMultSets_; i++) {
383 CRSet& crSet = data_->crMultSets_[i];
384
385 for(int j=0; j<1; j++) {
387 crSet.numNodes_,
388 crSet.nodeIDs_[j],
389 crSet.fieldIDs_,
390 crSet.weights_,
391 crSet.values_[j]))
392 }
393 }
394
395 for(i=0; i<data_->numCRPenSets_; i++) {
396 CRSet& crSet = data_->crPenSets_[i];
397
398 for(int j=0; j<1; j++) {
400 crSet.numNodes_,
401 crSet.nodeIDs_[j],
402 crSet.fieldIDs_,
403 crSet.weights_,
404 crSet.values_[j],
405 crSet.penValues_[j]))
406 }
407 }
408 return(0);
409}
410
411//----------------------------------------------------------------------------
413{
414 int i;
415
416 for(i=0; i<numMatrices; i++) {
418
419 for(int j=0; j<data_->numElemBlocks_; j++) {
420 ElemBlock& block = data_->elemBlocks_[j];
421
422 for(int el=0; el<block.numElements_; el++) {
423
425 block.elemIDs_[el],
426 block.elemConn_[el],
427 block.elemStiff_[el],
428 block.elemFormat_))
429 }
430 }
431 }
432
433 for(i=0; i<numRHSs; i++) {
435
436 for(int j=0; j<data_->numElemBlocks_; j++) {
437 ElemBlock& block = data_->elemBlocks_[j];
438
439 for(int el=0; el<block.numElements_; el++) {
441 block.elemIDs_[el],
442 block.elemConn_[el],
443 block.elemLoad_[el]))
444 }
445 }
446 }
447
448 //*************** Boundary Condition Node Sets *************************
449
450 for(i=0; i<data_->numBCNodeSets_; i++) {
451 BCNodeSet& bcSet = data_->bcNodeSets_[i];
452
454 bcSet.nodeIDs_,
455 bcSet.fieldID_,
456 bcSet.offsetsIntoField_,
457 bcSet.prescribed_values_))
458 }
459
460 double* matScalars = new double[numMatrices];
461 for(i=0; i<numMatrices; i++) {
462 matScalars[i] = 1.0;
463 }
464
465 int rhsScaleID = rhsIDs[0];
466 double rhsScalar = 1.0;
467
469 CHK_ERR(fei_->setRHSScalars(1, &rhsScaleID, &rhsScalar))
470
471 delete [] matScalars;
472 return(0);
473}
474
475//----------------------------------------------------------------------------
477{
478 std::string sname(data_->solverLibraryName_);
479 if (sname == "TEST_LSC") {
480 return( 0 );
481 }
482 //FEI_COUT << "numFields: " << data_->numFields_ << FEI_ENDL;
483 double* norms = new double[data_->numFields_];
484 int *fields = new int[data_->numFields_];
485 for(int i=0; i<data_->numFields_; ++i) {
486 fields[i] = data_->fieldIDs_[i];
487 }
488
489 CHK_ERR( fei_->residualNorm(1, data_->numFields_, fields, norms) );
490 //int n;
491 //for(n=0; n<data_->numFields_; n++) {
492 // FEI_COUT << " field["<<n<<"]: " << fields[n]
493 // << ", 1-norm: " << norms[n] << FEI_ENDL;
494 //}
495
496 delete [] fields;
497 delete [] norms;
498 return(0);
499}
500
501//----------------------------------------------------------------------------
503{
504 //let's exercise the putNodalFieldData function.
505
506 int numNodes;
507 CHK_ERR( fei_->getNumLocalNodes(numNodes) );
508 std::vector<int> nodeIDs(numNodes);
509 int* nodeIDsPtr = &nodeIDs[0];
510 int checkNumNodes;
511 CHK_ERR( fei_->getLocalNodeIDList(checkNumNodes, nodeIDsPtr, numNodes) );
512
513 for(int i=0; i<data_->numFields_; ++i) {
514 int fieldID = data_->fieldIDs_[i];
515 int fieldSize = data_->fieldSizes_[i];
516 std::vector<double> data(numNodes*fieldSize, 0.0001);
517
518 CHK_ERR( fei_->putNodalFieldData(fieldID, numNodes, nodeIDsPtr,
519 &data[0]) );
520 }
521
522 return(0);
523}
524
525//----------------------------------------------------------------------------
527 const char* solnFileName, int numProcs,
528 int localProc, int solveCounter)
529{
530 (void)solveCounter;
531 int i;
532
533 int maxNumEqnsPerNode = 0;
534 for(i=0; i<data.numFields_; ++i) {
535 maxNumEqnsPerNode += data.fieldSizes_[i];
536 }
537
538 std::vector<double> soln(maxNumEqnsPerNode);
539
540 int numNodes;
541 CHK_ERR( fei.getNumLocalNodes(numNodes) );
542
543 std::vector<GlobalID> nodes(numNodes);
544 int* nodesPtr = &nodes[0];
545
546 int checkNumNodes;
547 CHK_ERR( fei.getLocalNodeIDList( checkNumNodes, nodesPtr, numNodes) );
548
549 if (checkNumNodes != numNodes) {
550 ERReturn(-1);
551 }
552
553 FEI_OSTRINGSTREAM fileName;
554 fileName << solnFileName<<".node."<<solveCounter<<"."<<numProcs<<"."<<localProc;
555 FEI_OFSTREAM outfile(fileName.str().c_str());
556
557 if (!outfile || outfile.bad()) {
558 FEI_COUT << "ERROR opening solution output file " << fileName.str() << FEI_ENDL;
559 return(1);
560 }
561
562 outfile.setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
563
564 std::vector<int> offsets(2);
565
566 for(i=0; i<numNodes; ++i) {
567 CHK_ERR( fei.getNodalSolution(1, &(nodesPtr[i]),
568 &offsets[0], &soln[0]) );
569
570 int numDOF = offsets[1];
571
572 outfile << nodesPtr[i] << " " << numDOF << FEI_ENDL;
573 for(int j=0; j<numDOF; j++) {
574 outfile << soln[j] << " ";
575 }
576 outfile << FEI_ENDL;
577 }
578
579 return(0);
580}
581
582//----------------------------------------------------------------------------
584 const char* solnFileName,
585 int numProcs, int localProc,
586 int solveCounter)
587{
588 int returnValue = 0;
589 FEI_OSTRINGSTREAM fileName;
590 fileName << solnFileName<<".elem."<<solveCounter<<"."<<numProcs<<"."<<localProc;
591 FEI_OFSTREAM outfile(fileName.str().c_str());
592
593 if (!outfile || outfile.bad()) {
594 FEI_COUT << "ERROR opening elem-solution output file " << fileName.str() << FEI_ENDL;
595 return(1);
596 }
597
598 for(int i=0; i<data.numElemBlocks_; i++) {
599 if (returnValue != 0) break;
600
601 ElemBlock& eb = data.elemBlocks_[i];
602
603 GlobalID blockID = eb.blockID_;
604 int numElems;
605 CHK_ERR( fei.getNumBlockElements(blockID, numElems))
606
607 int dofPerElem;
608 CHK_ERR( fei.getNumBlockElemDOF(blockID, dofPerElem))
609 int totalNumElemDOF = numElems*dofPerElem;
610
611 if (totalNumElemDOF < 1) {
612 continue;
613 }
614
615 GlobalID* elemIDs = new GlobalID[numElems];
616 if (elemIDs==NULL) return(-1);
617
618 int err = fei.getBlockElemIDList(blockID, numElems, elemIDs);
619 if (err) returnValue = 1;
620
621 int* offsets = new int[numElems+1];
622 if (offsets == NULL) return(-1);
623
624 if (totalNumElemDOF > 0) {
625 double* solnValues = new double[totalNumElemDOF];
626 if (solnValues == NULL) return(-1);
627
628 err = fei.getBlockElemSolution(blockID, numElems, elemIDs,
629 dofPerElem, solnValues);
630 if (err) returnValue = 1;
631
632 if (!err) {
633 for(int j=0; j<numElems; j++) {
634
635 outfile << (int)elemIDs[j] << " " << dofPerElem << FEI_ENDL << " ";
636 for(int k=0; k<dofPerElem; k++) {
637 outfile << solnValues[j*dofPerElem + k] << " ";
638 }
639 outfile << FEI_ENDL;
640 }
641 }
642
643 delete [] solnValues;
644 }
645
646 delete [] elemIDs;
647 delete [] offsets;
648 }
649
650 return(0);
651}
652
653//----------------------------------------------------------------------------
655 const char* solnFileName,
656 int numProcs, int localProc, int solveCounter)
657{
658 int numCRs = 0;
659
660 CHK_ERR( fei.getNumCRMultipliers(numCRs) );
661
662 int* globalNumCRs = new int[numProcs];
663#ifndef FEI_SER
664 if (MPI_Allgather(&numCRs, 1, MPI_INT, globalNumCRs, 1, MPI_INT,
665 comm_) != MPI_SUCCESS) {
666 ERReturn(-1);
667 }
668#endif
669
670 int localCRStart = 0;
671 for(int p=0; p<localProc; p++) localCRStart += globalNumCRs[p];
672
673 delete [] globalNumCRs;
674
675 FEI_OSTRINGSTREAM fileName;
676 fileName << solnFileName<<".mult."<<solveCounter<<"."<<numProcs<<"."<<localProc;
677 FEI_OFSTREAM outfile(fileName.str().c_str());
678
679 if (!outfile || outfile.bad()) {
680 FEI_COUT << "ERROR opening mult-solution output file " << fileName.str() << FEI_ENDL;
681 return(-1);
682 }
683
684 int* CRIDs = numCRs > 0 ? new int[numCRs] : NULL;
685 double* results = numCRs > 0 ? new double[numCRs] : NULL;
686
687 if (numCRs > 0 && (CRIDs==NULL || results==NULL)) {
688 ERReturn(-1);
689 }
690
691 if (numCRs < 1) {
692 return(0);
693 }
694
695 CHK_ERR( fei.getCRMultIDList(numCRs, CRIDs) );
696
697 std::string sname(data_->solverLibraryName_);
698 if (sname == "FETI") {
699 for(int ii=0; ii<numCRs; ++ii) results[ii] = -999.99;
700 }
701 else {
702 CHK_ERR( fei.getCRMultipliers(numCRs, CRIDs, results));
703 }
704
705 for(int i=0; i<numCRs; i++) {
706 outfile << localCRStart++ << " " << 1 << FEI_ENDL;
707
708 outfile << " " << results[i] << FEI_ENDL;
709 }
710
711 delete [] CRIDs;
712 delete [] results;
713
714 return(0);
715}
716
717//----------------------------------------------------------------------------
719{
720 if (localProc_ == 0) {
721 char* current_dir = NULL;
722 CHK_ERR( fei_test_utils::dirname(data_->solnFileName_.c_str(), current_dir));
723
724 FEI_OSTRINGSTREAM solnMtxName;
725 solnMtxName<< current_dir<<"/A_TLSC.mtx";
726 fei::FillableMat solnMtx, checkMtx;
727 CHK_ERR( SolnCheck::readMatrix(solnMtxName.str().c_str(), numProcs_, solnMtx) );
729 int err = SolnCheck::compareMatrices(solnMtx, checkMtx);
730 delete [] current_dir;
731 if (err == 0) {
732 FEI_COUT << "Utst_fei_lsc: TEST PASSED" << FEI_ENDL;
733 }
734 else {
735 FEI_COUT << "Utst_fei_lsc: TEST FAILED" << FEI_ENDL;
736 return(-1);
737 }
738 }
739
740 return(0);
741}
int numNodes_
Definition: BCNodeSet.hpp:21
int * offsetsIntoField_
Definition: BCNodeSet.hpp:24
GlobalID * nodeIDs_
Definition: BCNodeSet.hpp:22
int fieldID_
Definition: BCNodeSet.hpp:23
double * prescribed_values_
Definition: BCNodeSet.hpp:25
Definition: CRSet.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
int ** procs_
Definition: CommNodeSet.hpp:19
GlobalID * nodeIDs_
Definition: CommNodeSet.hpp:18
int * procsPerNode_
Definition: CommNodeSet.hpp:20
char ** paramStrings_
Definition: DataReader.hpp:39
BCNodeSet * bcNodeSets_
Definition: DataReader.hpp:60
CRSet * crMultSets_
Definition: DataReader.hpp:51
CRSet * slaveVars_
Definition: DataReader.hpp:54
int solveType_
Definition: DataReader.hpp:28
int * fieldIDs_
Definition: DataReader.hpp:35
std::string solnFileName_
Definition: DataReader.hpp:31
ElemBlock * elemBlocks_
Definition: DataReader.hpp:42
CommNodeSet * sharedNodeSets_
Definition: DataReader.hpp:63
int numCRPenSets_
Definition: DataReader.hpp:56
int * fieldSizes_
Definition: DataReader.hpp:36
int numSharedNodeSets_
Definition: DataReader.hpp:62
int numElemBlocks_
Definition: DataReader.hpp:41
int numCRMultSets_
Definition: DataReader.hpp:50
std::string solverLibraryName_
Definition: DataReader.hpp:30
int numFields_
Definition: DataReader.hpp:34
int numBCNodeSets_
Definition: DataReader.hpp:59
int numParams_
Definition: DataReader.hpp:38
std::string checkFileName_
Definition: DataReader.hpp:32
int numSlaveVars_
Definition: DataReader.hpp:53
CRSet * crPenSets_
Definition: DataReader.hpp:57
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 elemFormat_
Definition: ElemBlock.hpp:25
int interleaveStrategy_
Definition: ElemBlock.hpp:30
int ** nodalFieldIDs_
Definition: ElemBlock.hpp:21
int * elemDOFFieldIDs_
Definition: ElemBlock.hpp:29
int localProc_
Definition: FEI_tester.hpp:84
int createFEIinstance(const char *solverName)
Definition: FEI_tester.cpp:210
int testCheckResult()
Definition: FEI_tester.cpp:171
int initializationPhase()
Definition: FEI_tester.cpp:256
int numMatrices
Definition: FEI_tester.hpp:86
int testInitialization()
Definition: FEI_tester.cpp:61
int lsc_matrix_check()
Definition: FEI_tester.cpp:718
int save_multiplier_soln(DataReader &data, FEI &fei, const char *solnFileName, int numProcs, int localProc, int solveCounter)
Definition: FEI_tester.cpp:654
int save_block_node_soln(DataReader &data, FEI &fei, const char *solnFileName, int numProcs, int localProc, int solveCounter)
Definition: FEI_tester.cpp:526
int save_block_elem_soln(DataReader &data, FEI &fei, const char *solnFileName, int numProcs, int localProc, int solveCounter)
Definition: FEI_tester.cpp:583
void dumpMatrixFiles()
Definition: FEI_tester.cpp:161
fei::SharedPtr< FEI > fei_
Definition: FEI_tester.hpp:78
bool useNewFEI_
Definition: FEI_tester.hpp:90
int aggregateLoadPhase()
Definition: FEI_tester.cpp:412
int exerciseResidualNorm()
Definition: FEI_tester.cpp:476
int setIDlists()
Definition: FEI_tester.cpp:235
int testLoading()
Definition: FEI_tester.cpp:101
fei::SharedPtr< DataReader > data_
Definition: FEI_tester.hpp:82
int normalLoadPhase()
Definition: FEI_tester.cpp:339
int testSolve()
Definition: FEI_tester.cpp:130
int exercisePutFunctions()
Definition: FEI_tester.cpp:502
fei::SharedPtr< LibraryWrapper > wrapper_
Definition: FEI_tester.hpp:80
FEI_tester(fei::SharedPtr< DataReader > data_reader, MPI_Comm comm, int localProc, int numProcs, bool useNewFEI=false)
Definition: FEI_tester.cpp:37
MPI_Comm comm_
Definition: FEI_tester.hpp:76
void setParameter(const char *param)
Definition: FEI_tester.cpp:166
int * matrixIDs
Definition: FEI_tester.hpp:87
int * rhsIDs
Definition: FEI_tester.hpp:89
Definition: FEI.hpp:144
virtual int initComplete()=0
virtual int getLocalNodeIDList(int &numNodes, GlobalID *nodeIDs, int lenNodeIDs)=0
virtual int residualNorm(int whichNorm, int numFields, int *fieldIDs, double *norms)=0
virtual int setCurrentMatrix(int matrixID)=0
virtual int initElemBlock(GlobalID elemBlockID, int numElements, int numNodesPerElement, const int *numFieldsPerNode, const int *const *nodalFieldIDs, int numElemDofFieldsPerElement, const int *elemDOFFieldIDs, int interleaveStrategy)=0
virtual int sumInElemMatrix(GlobalID elemBlockID, GlobalID elemID, const GlobalID *elemConn, const double *const *elemStiffness, int elemFormat)=0
virtual int iterations(int &itersTaken) const =0
virtual int loadCRPen(int CRPenID, int numCRNodes, const GlobalID *CRNodeIDs, const int *CRFieldIDs, const double *CRWeights, double CRValue, double penValue)=0
virtual int initCRPen(int numCRNodes, const GlobalID *CRNodeIDs, const int *CRFieldIDs, int &CRID)=0
virtual int putNodalFieldData(int fieldID, int numNodes, const GlobalID *nodeIDs, const double *data)=0
virtual int setMatScalars(int numScalars, const int *IDs, const double *scalars)=0
virtual int initCRMult(int numCRNodes, const GlobalID *CRNodeIDs, const int *CRFieldIDs, int &CRID)=0
virtual int loadCRMult(int CRMultID, int numCRNodes, const GlobalID *CRNodeIDs, const int *CRFieldIDs, const double *CRWeights, double CRValue)=0
virtual int initFields(int numFields, const int *fieldSizes, const int *fieldIDs, const int *fieldTypes=NULL)=0
virtual int loadComplete(bool applyBCs=true, bool globalAssemble=true)=0
virtual int getNumLocalNodes(int &numNodes)=0
virtual int sumInElemRHS(GlobalID elemBlockID, GlobalID elemID, const GlobalID *elemConn, const double *elemLoad)=0
virtual int setIDLists(int numMatrices, const int *matrixIDs, int numRHSs, const int *rhsIDs)=0
virtual int setCurrentRHS(int rhsID)=0
virtual int resetMatrix(double s=0.0)=0
virtual int version(const char *&versionString)=0
virtual int setRHSScalars(int numScalars, const int *IDs, const double *scalars)=0
virtual int parameters(int numParams, const char *const *paramStrings)=0
virtual int resetRHSVector(double s=0.0)=0
virtual int initSharedNodes(int numSharedNodes, const GlobalID *sharedNodeIDs, const int *numProcsPerNode, const int *const *sharingProcIDs)=0
virtual int initSlaveVariable(GlobalID slaveNodeID, int slaveFieldID, int offsetIntoSlaveField, int numMasterNodes, const GlobalID *masterNodeIDs, const int *masterFieldIDs, const double *weights, double rhsValue)=0
virtual int resetSystem(double s=0.0)=0
virtual int getNumBlockActNodes(GlobalID elemBlockID, int &numNodes) const =0
virtual int setSolveType(int solveType)=0
virtual int solve(int &status)=0
virtual int loadNodeBCs(int numNodes, const GlobalID *nodeIDs, int fieldID, const int *offsetsIntoField, const double *prescribedValues)=0
virtual int initElem(GlobalID elemBlockID, GlobalID elemID, const GlobalID *elemConn)=0
void reset(T *p=0)
T * get() const
#define ERReturn(a)
#define CHK_ERR(a)
#define FEI_AGGREGATE_SUM
Definition: fei_defs.h:67
int GlobalID
Definition: fei_defs.h:60
#define FEI_SINGLE_SYSTEM
Definition: fei_defs.h:65
#define FEI_OFSTREAM
Definition: fei_fstream.hpp:14
#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
Definition: fei_sstream.hpp:32
#define FEI_Implementation
Definition: fei_version.h:66
int compareMatrices(fei::FillableMat &mat1, fei::FillableMat &mat2)
Definition: SolnCheck.cpp:57
int checkSolution(int localProc, int numProcs, const char *solnFileName, const char *checkFileName, const char *extension, int solveCounter)
Definition: SolnCheck.cpp:63
int readMatrix(const char *baseName, int np, fei::FillableMat &matrix)
Definition: SolnCheck.cpp:51
int dirname(const char *name, char *&dir)
std::ostream & console_out()
fei::SharedPtr< LibraryWrapper > create_LibraryWrapper(MPI_Comm comm, const char *libraryName)
int getIntParamValue(const char *key, int numParams, const char *const *params, int &paramValue)