9#include <fei_macros.hpp>
12#include <test_utils/fei_test_utils.hpp>
13#include <test_utils/LibraryFactory.hpp>
14#include <test_utils/test_Matrix.hpp>
15#include <test_utils/test_VectorSpace.hpp>
16#include <test_utils/test_MatrixGraph.hpp>
17#include <fei_MatrixGraph_Impl2.hpp>
18#include <fei_Factory.hpp>
20#include <snl_fei_Factory.hpp>
21#include <fei_Vector_Impl.hpp>
22#include <fei_Matrix_Impl.hpp>
24#ifdef HAVE_FEI_AZTECOO
25#include <fei_Aztec_LinSysCore.hpp>
27#include <fei_Factory_Trilinos.hpp>
30#include <FETI_DP_FiniteElementData.h>
34#define fei_file "test_Matrix.cpp"
35#include <fei_ErrMacros.hpp>
37int test_matrix_unit1()
39 std::vector<double> data1D;
40 std::vector<const double*>data2D;
44 double** values =
new double*[numRows];
46 for(i=0; i<numRows; ++i) {
47 values[i] =
new double[numCols];
48 for(j=0; j<numCols; ++j) {
53 fei::Matrix_core::copyTransposeToWorkArrays(numRows, numCols,
54 values, data1D, data2D);
56 for(i=0; i<numRows; ++i) {
57 for(j=0; j<numCols; ++j) {
58 if (std::abs(values[i][j] - (data2D[j])[i]) > 1.e-49) {
70void test_Matrix_unit2(MPI_Comm comm,
int numProcs,
int localProc)
76 FEI_COUT <<
"testing fei::Matrix_Impl...";
81 int rowfield = 0, rowfieldsize = 1;
83 rowspace->defineFields(1, &rowfield, &rowfieldsize);
84 rowspace->defineIDTypes(1, &idType);
88 int patternID1 = mgraph->definePattern(2, idType, rowfield);
90 fei::Pattern* rowpattern = mgraph->getPattern(patternID1);
92 mgraph->initConnectivityBlock(0, 1, patternID1);
94 std::vector<int> ids(2);
95 ids[0] = 0; ids[1] = 1;
97 int err = mgraph->initConnectivity(0, 0, &ids[0]);
99 FEI_OSTRINGSTREAM osstr;
100 osstr <<
"test_Matrix_unit2, initConnectivity returned err="<<err;
101 throw std::runtime_error(osstr.str());
104 err = mgraph->initComplete();
106 FEI_OSTRINGSTREAM osstr;
107 osstr <<
"test_Matrix_unit2, initComplete returned err="<<err;
108 throw std::runtime_error(osstr.str());
111 bool factory_created =
false;
115 factory_created =
true;
119 if (!factory_created) {
120 FEI_COUT <<
"failed to create Trilinos factory."<<FEI_ENDL;
128 std::vector<double> coefs(numrowindices*numrowindices, 1.0);
129 std::vector<double*> coefs_2D(numrowindices);
130 for(
int i=0; i<numrowindices; ++i) {
131 coefs_2D[i] = &(coefs[i*numrowindices]);
134 err = feimat->sumIn(0, 0, &coefs_2D[0]);
136 FEI_OSTRINGSTREAM osstr;
137 osstr <<
"test_Matrix_unit2, feimat->sumIn returned err="<<err;
138 throw std::runtime_error(osstr.str());
141 err = feimat->globalAssemble();
143 FEI_OSTRINGSTREAM osstr;
144 osstr <<
"test_Matrix_unit2, feimat->globalAssemble returned err="<<err;
145 throw std::runtime_error(osstr.str());
148 err = feimat->writeToFile(
"feimat2.mtx",
false);
150 FEI_OSTRINGSTREAM osstr;
151 osstr <<
"test_Matrix_unit2, feimat->writeToFile returned err="<<err;
152 throw std::runtime_error(osstr.str());
155 fei::FillableMat feimat_ss;
156 err = fei_test_utils::copy_feiMatrix_to_FillableMat(*feimat, feimat_ss);
158 FEI_OSTRINGSTREAM osstr;
159 osstr <<
"test_Matrix_unit2, copy_feiMatrix_to_FillableMat returned err="<<err;
160 throw std::runtime_error(osstr.str());
163 fei_test_utils::writeMatrix(
"feimat_ss2.mtx", feimat_ss);
165 FEI_COUT <<
"ok"<<FEI_ENDL;
168void test_Matrix_unit4(MPI_Comm comm,
int numProcs,
int localProc)
174 FEI_COUT <<
"testing fei::Matrix_Impl with FEI_BLOCK_DIAGONAL_ROW...";
179 int rowfield = 0, rowfieldsize = 2;
181 rowspace->defineFields(1, &rowfield, &rowfieldsize);
182 rowspace->defineIDTypes(1, &idType);
186 int patternID1 = mgraph->definePattern(2, idType, rowfield);
188 fei::Pattern* rowpattern = mgraph->getPattern(patternID1);
190 bool diagonal =
true;
191 mgraph->initConnectivityBlock(0, 1, patternID1, diagonal);
193 std::vector<int> ids(2);
194 ids[0] = 0; ids[1] = 1;
196 int err = mgraph->initConnectivity(0, 0, &ids[0]);
198 FEI_OSTRINGSTREAM osstr;
199 osstr <<
"test_Matrix_unit4, initConnectivity returned err="<<err;
200 throw std::runtime_error(osstr.str());
203 err = mgraph->initComplete();
205 FEI_OSTRINGSTREAM osstr;
206 osstr <<
"test_Matrix_unit4, initComplete returned err="<<err;
207 throw std::runtime_error(osstr.str());
215 FEI_COUT <<
"Trilinos not available."<<FEI_ENDL;
223 paramset.
add(blktrue);
224 factory->parameters(paramset);
228 paramset.
add(blkfalse);
229 factory->parameters(paramset);
235 std::vector<double> coefs(numrowindices*rowfieldsize*rowfieldsize, 1.0);
236 std::vector<double*> coefs_2D(numrowindices*rowfieldsize);
238 for(
int i=0; i<numrowindices*rowfieldsize; ++i) {
239 coefs_2D[i] = &(coefs[offset]);
240 offset += rowfieldsize;
243 err = feimat->sumIn(0, 0, &coefs_2D[0], FEI_BLOCK_DIAGONAL_ROW);
245 FEI_OSTRINGSTREAM osstr;
246 osstr <<
"test_Matrix_unit4, feimat->sumIn returned err="<<err;
247 throw std::runtime_error(osstr.str());
250 err = feiblkmat->sumIn(0, 0, &coefs_2D[0], FEI_BLOCK_DIAGONAL_ROW);
252 FEI_OSTRINGSTREAM osstr;
253 osstr <<
"test_Matrix_unit4, feiblkmat->sumIn returned err="<<err;
254 throw std::runtime_error(osstr.str());
257 err = feimat->globalAssemble();
259 FEI_OSTRINGSTREAM osstr;
260 osstr <<
"test_Matrix_unit4, feimat->globalAssemble returned err="<<err;
261 throw std::runtime_error(osstr.str());
264 err = feiblkmat->globalAssemble();
266 FEI_OSTRINGSTREAM osstr;
267 osstr <<
"test_Matrix_unit4, feimat->globalAssemble returned err="<<err;
268 throw std::runtime_error(osstr.str());
271 feimat->writeToFile(
"feimat_blkdiag.mtx");
272 feiblkmat->writeToFile(
"feiblkmat_blkdiag.mtx");
274 FEI_COUT <<
"ok"<<FEI_ENDL;
277test_Matrix::test_Matrix(MPI_Comm comm)
282test_Matrix::~test_Matrix()
286int test_Matrix::runtests()
291 if (localProc_==0) FEI_COUT <<
"Testing Factory_Trilinos fei::Matrix..." << FEI_ENDL;
298 if (localProc_==0) FEI_COUT << FEI_ENDL;
300#ifdef HAVE_FEI_AZTECOO
307 if (localProc_==0) FEI_COUT <<
"Testing 'old' factory_aztec fei::Matrix..." << FEI_ENDL;
322 testData test_data(localProc_, numProcs_);
325 test_VectorSpace::create_VectorSpace(comm_, &test_data, localProc_, numProcs_,
326 false,
false, (
const char*)0, factory);
327 int err = vspace->initComplete();
329 FEI_COUT <<
"ERROR, failed to create valid fei::VectorSpace." << FEI_ENDL;
330 throw std::runtime_error(
"test_Vector::vector_test1: ERROR, failed to create valid fei::VectorSpace.");
334 factory->createMatrixGraph(vspace, vspace, NULL);
336 std::vector<int>& fieldIDs = test_data.fieldIDs;
337 std::vector<int>& idTypes = test_data.idTypes;
338 std::vector<int>& ids = test_data.ids;
340 int numIDs = ids.size();
341 int fieldID = fieldIDs[0];
342 int idType = idTypes[0];
344 int patternID = mgraph->definePattern(numIDs, idType, fieldID);
346 mgraph->initConnectivityBlock(0, 1, patternID);
348 mgraph->initConnectivity(0, 0, &ids[0]);
350 mgraph->initComplete();
359 FEI_COUT <<
" matrix_test1: testing fei::Matrix with type '"
360 << mat->typeName() <<
"':"<<FEI_ENDL;
367 FEI_COUT <<
" testing get{Global/Local}NumRows,getRowLength...";
369 int mglobalrows = mat->getGlobalNumRows();
370 int vglobaleqns = rspace->getGlobalNumIndices();
372 if (mglobalrows != vglobaleqns) {
373 throw std::runtime_error(
"mat reports different num rows than vector-space eqns");
376 std::vector<int> global_offsets;
377 rspace->getGlobalIndexOffsets(global_offsets);
379 int my_num_rows = mat->getLocalNumRows();
380 if (my_num_rows != global_offsets[localProc_+1]-global_offsets[localProc_]) {
381 throw std::runtime_error(
"num-local-rows mis-match between mat and vector-space");
384 int i, my_first_row = global_offsets[localProc_];
385 std::vector<int> row_lengths(my_num_rows);
387 for(i=0; i<my_num_rows; ++i) {
388 int errcode = mat->getRowLength(i+my_first_row, row_lengths[i]);
390 throw std::runtime_error(
"nonzero errcode from mat->getRowLength");
394 if (localProc_==0) FEI_COUT <<
"ok" << FEI_ENDL;
397int test_Matrix::serialtest1()
400 std::vector<int>& fieldIDs = testdata->fieldIDs;
401 std::vector<int>& fieldSizes = testdata->fieldSizes;
402 std::vector<int>& idTypes = testdata->idTypes;
403 std::vector<int>& ids = testdata->ids;
407 vspc->defineFields(fieldIDs.size(),
411 vspc->defineIDTypes(idTypes.size(), &idTypes[0]);
416 int numIDs = ids.size();
417 int fieldID = fieldIDs[0];
418 int idType = idTypes[0];
420 int patternID = matgraph->definePattern(numIDs, idType, fieldID);
422 CHK_ERR( matgraph->initConnectivityBlock(0, 1, patternID) );
424 CHK_ERR( matgraph->initConnectivity(0, 0, &ids[0]) );
426 CHK_ERR( matgraph->initComplete() );
429 int localsize = matgraph->getRowSpace()->getNumIndices_Owned();
434 std::vector<int> indices(numIDs);
435 CHK_ERR( matgraph->getConnectivityIndices(0, 0, numIDs, &indices[0], numIDs) );
437 std::vector<double> data1(numIDs*numIDs);
438 std::vector<double*> data2d(numIDs);
441 for(i=0; i<numIDs; ++i) {
442 data2d[i] = &(data1[i*numIDs]);
445 for(i=0; i<numIDs*numIDs; ++i) {
449 CHK_ERR( matrix->sumIn(numIDs, &indices[0], numIDs, &indices[0],
452 CHK_ERR( matrix->sumIn(0, 0, &data2d[0], 0) );
454 CHK_ERR( matrixT->sumIn(numIDs, &indices[0],
455 numIDs, &indices[0], &data2d[0], 3) );
457 CHK_ERR( matrixT->sumIn(0, 0, &data2d[0], 3) );
459 if (*ssmat != *ssmatT) {
466int test_Matrix::serialtest2()
469 std::vector<int>& idTypes = testdata->idTypes;
470 std::vector<int>& ids = testdata->ids;
474 vspc->defineIDTypes(idTypes.size(), &idTypes[0]);
478 int numIDs = ids.size();
479 int idType = idTypes[0];
481 int patternID = matgraph->definePattern(numIDs, idType);
483 CHK_ERR( matgraph->initConnectivityBlock(0, 1, patternID) );
485 CHK_ERR( matgraph->initConnectivity(0, 0, &ids[0]) );
487 CHK_ERR( matgraph->initComplete() );
490 int localsize = matgraph->getRowSpace()->getNumIndices_Owned();
495 std::vector<int> indices(numIDs);
496 CHK_ERR( matgraph->getConnectivityIndices(0, 0, numIDs, &indices[0], numIDs) );
498 std::vector<double> data1(numIDs*numIDs);
499 std::vector<double*> data2d(numIDs);
502 for(i=0; i<numIDs; ++i) {
503 data2d[i] = &(data1[i*numIDs]);
506 for(i=0; i<numIDs*numIDs; ++i) {
510 CHK_ERR( matrix->
sumIn(numIDs, &indices[0],
511 numIDs, &indices[0], &data2d[0], 0) );
513 CHK_ERR( matrixT->
sumIn(numIDs, &indices[0],
514 numIDs, &indices[0], &data2d[0], 3) );
516 if (*ssmat != *ssmatT) {
527int test_Matrix::serialtest3()
530 std::vector<int>& fieldIDs = testdata->fieldIDs;
531 std::vector<int>& fieldSizes = testdata->fieldSizes;
532 std::vector<int>& idTypes = testdata->idTypes;
533 std::vector<int>& ids = testdata->ids;
537 vspc->defineFields(fieldIDs.size(), &fieldIDs[0], &fieldSizes[0]);
539 vspc->defineIDTypes(idTypes.size(), &idTypes[0]);
544 int numIDs = ids.size();
545 int fieldID = fieldIDs[0];
546 int idType = idTypes[0];
548 int patternID = matgraph->definePattern(numIDs, idType, fieldID);
550 CHK_ERR( matgraph->initConnectivityBlock(0, 1, patternID) );
552 CHK_ERR( matgraph->initConnectivity(0, 0, &ids[0]) );
556 int offsetOfSlave = 1;
557 int offsetIntoSlaveField = 0;
558 std::vector<double> weights(2);
561 double rhsValue = 0.0;
562 std::vector<int> cr_idtypes(2, idTypes[0]);
563 std::vector<int> cr_fieldIDs(2, fieldIDs[0]);
565 CHK_ERR( matgraph->initSlaveConstraint(2,
570 offsetIntoSlaveField,
574 CHK_ERR( matgraph->initComplete() );
577 int localsize = matgraph->getRowSpace()->getNumIndices_Owned();
581 if (matrix == NULL) {
585 std::vector<int> indices(numIDs);
586 CHK_ERR( matgraph->getConnectivityIndices(0, 0, numIDs,
587 &indices[0], numIDs) );
589 std::vector<double> data1(numIDs*numIDs);
590 std::vector<double*> data2d(numIDs);
593 for(i=0; i<numIDs; ++i) {
594 data2d[i] = &(data1[i*numIDs]);
597 for(i=0; i<numIDs*numIDs; ++i) {
601 CHK_ERR( matrix->
sumIn(numIDs, &indices[0],
602 numIDs, &indices[0], &data2d[0], 0) );
604 CHK_ERR( matrix->
sumIn(0, 0, &data2d[0], 0) );
612int test_Matrix::test1()
617int test_Matrix::test2()
622int test_Matrix::test3()
626 std::vector<int>& idTypes = testdata->idTypes;
627 std::vector<int>& ids = testdata->ids;
631 std::string paramstr(
"debugOutput .");
632 char* param =
const_cast<char*
>(paramstr.c_str());
634 CHK_ERR( fedata->parameters(1, ¶m) );
639 test_VectorSpace::create_VectorSpace(comm_,
640 testdata, localProc_, numProcs_,
641 false,
false,
"U_FEMat", factory);
644 test_MatrixGraph::create_MatrixGraph(testdata, localProc_, numProcs_,
645 false,
false,
"U_FEMat", vectorSpacePtr,
648 CHK_ERR( matrixGraphPtr->initComplete() );
661 int numIndices = matrixGraphPtr->getConnectivityNumIndices(blockID);
663 std::vector<int> indicesArray(numIndices);
664 int* indicesPtr = &indicesArray[0];
666 int checkNumIndices = 0;
667 CHK_ERR( matrixGraphPtr->getConnectivityIndices(blockID, 0,
668 numIndices, indicesPtr,
671 std::vector<double> data(ids.size(), 1.0);
672 double* dptr = &data[0];
673 std::vector<double*> coefPtrs(ids.size(), dptr);
675 CHK_ERR( mat_fed->sumIn(blockID, 0, &coefPtrs[0]) );
677 CHK_ERR( vec_fed->sumIn(blockID, 0, &data[0]) );
679 CHK_ERR( mat_fed->gatherFromOverlap() );
681 CHK_ERR( fedata->loadComplete() );
693int test_Matrix::test4()
virtual int sumIn(int numRows, const int *rows, int numCols, const int *cols, const double *const *values, int format=0)=0
void add(const Param ¶m, bool maintain_unique_keys=true)
int getNumIndices() const
fei::SharedPtr< fei::Factory > create_fei_Factory(MPI_Comm comm, const char *libraryName)