FEI Version of the Day
Loading...
Searching...
No Matches
test_VectorSpace.cpp
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
10#include <fei_macros.hpp>
11
12#include <test_utils/test_VectorSpace.hpp>
13
14#include <snl_fei_Factory.hpp>
15#include <fei_ParameterSet.hpp>
16#include <fei_LibraryWrapper.hpp>
17#include <assert.h>
18
19#undef fei_file
20#define fei_file "test_VectorSpace.cpp"
21#include <fei_ErrMacros.hpp>
22
23test_VectorSpace::test_VectorSpace(MPI_Comm comm)
24 : tester(comm)
25{
26}
27
28test_VectorSpace::~test_VectorSpace()
29{
30}
31
32int test_VectorSpace::runtests()
33{
34 CHK_ERR( test0() );
35 CHK_ERR( test1() );
36 CHK_ERR( test2() );
37 CHK_ERR( test3() );
38 CHK_ERR( test4() );
39 return(0);
40}
41
42int test_VectorSpace::test0()
43{
44 int i, ID = 0;
45 int idType = 0;
46 std::vector<int> fieldIDs(numProcs_);
47 int fieldSize = 2;
48 std::vector<int> fieldSizes(numProcs_, fieldSize);
49
50 for(i=0; i<numProcs_; ++i) fieldIDs[i] = i;
51
52 fei::VectorSpace vspace(comm_);
53
54 vspace.defineIDTypes(1, &idType);
55 vspace.defineFields(fieldIDs.size(), &fieldIDs[0], &fieldSizes[0]);
56
57 CHK_ERR( vspace.addDOFs(idType, 1, &ID) );
58
59 for(i=0; i<numProcs_; ++i) {
60 if (i == localProc_) continue;
61
62 int numSharingProcsPerID = 1;
63 int sharingProc = i;
64
65 CHK_ERR( vspace.initSharedIDs(1, idType, &ID,
66 &numSharingProcsPerID, &sharingProc) );
67 }
68
69 CHK_ERR( vspace.addDOFs(fieldIDs[localProc_], idType,
70 1, &ID) );
71
72 CHK_ERR( vspace.initComplete() );
73
74 std::vector<int> globalIndexOffsets;
75
76 vspace.getGlobalIndexOffsets(globalIndexOffsets);
77
78 return(0);
79}
80
81int test_VectorSpace::test1()
82{
83 //A general test of fei::VectorSpace, including defining fields and
84 //idTypes, initializing solution entries, shared ids if numProcs_ > 1, then
85 //testing getGlobalIndexOffsets and getGlobalIndex several times.
86
87 testData* testdata = new testData(localProc_, numProcs_);
88 std::vector<int>& fieldIDs = testdata->fieldIDs;
89 std::vector<int>& idTypes = testdata->idTypes;
90 std::vector<int>& ids = testdata->ids;
91
92 int numDOFsPerID = 4; //1 scalar field and 1 3-vector field
93
95 fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, wrapper));
96
98 test_VectorSpace::create_VectorSpace(comm_,
99 testdata, localProc_, numProcs_,
100 true, //defineBothFields
101 true, //initSolnBothFields
102 "U", //name
103 factory);
104
105 CHK_ERR( vectorSpacePtr->initComplete() );
106
107 fei::VectorSpace& vectorSpace = *vectorSpacePtr;
108
109 std::vector<int> globalOffsets;
110
111 vectorSpace.getGlobalIndexOffsets(globalOffsets);
112
113 if (localProc_ > 0) {
114 if (globalOffsets[localProc_] != (localProc_+1)*2*numDOFsPerID) {
115 ERReturn(-1);
116 }
117 }
118
119 if (globalOffsets[localProc_+1] != (localProc_+2)*2*numDOFsPerID) {
120 ERReturn(-1);
121 }
122
123 std::vector<int> globalBlkOffsets;
124 vectorSpace.getGlobalBlkIndexOffsets(globalBlkOffsets);
125 if (localProc_ > 0) {
126 if (globalBlkOffsets[localProc_] != (localProc_+1)*2) {
127 ERReturn(-1);
128 }
129 }
130
131 if (globalBlkOffsets[localProc_+1] != (localProc_+2)*2) {
132 ERReturn(-1);
133 }
134
135 int len = ids.size();
136 int globalIndex = 0;
137 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[0], fieldIDs[0],
138 0, 0, globalIndex) );
139
140 int correctIndex = globalOffsets[localProc_];
141 int correctBlkIndex = globalBlkOffsets[localProc_];
142 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
143 if (localProc_ != 0) correctBlkIndex -= 2;
144
145 if (globalIndex != correctIndex) {
146 ERReturn(-1);
147 }
148
149 int globalBlkIndex = 0;
150 CHK_ERR( vectorSpace.getGlobalBlkIndex(idTypes[0], ids[0], globalBlkIndex) );
151 if (globalBlkIndex != correctBlkIndex) {
152 fei::console_out() << "localProc_: " << localProc_ << ", globalBlkIndex: "
153 << globalBlkIndex << ", correctBlkIndex: " << correctBlkIndex << FEI_ENDL;
154 ERReturn(-1);
155 }
156
157 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[1], fieldIDs[0],
158 0, 0, globalIndex) );
159
160 correctIndex = globalOffsets[localProc_] + 4;
161 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
162
163 if (globalIndex != correctIndex) {
164 ERReturn(-1);
165 }
166
167 CHK_ERR( vectorSpace.getGlobalBlkIndex(idTypes[0], ids[1], globalBlkIndex) );
168 correctBlkIndex = globalBlkOffsets[localProc_]+1;
169 if (localProc_ != 0) correctBlkIndex -= 2;
170
171 if (globalBlkIndex != correctBlkIndex) {
172 fei::console_out() << "2localProc_: " << localProc_ << ", globalBlkIndex: "
173 << globalBlkIndex << ", correctBlkIndex: " << correctBlkIndex << FEI_ENDL;
174 ERReturn(-1);
175 }
176
177 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[len-1], fieldIDs[0],
178 0, 0, globalIndex) );
179 correctIndex = globalOffsets[localProc_] + 12;
180 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
181
182 if (globalIndex != correctIndex) {
183 ERReturn(-1);
184 }
185
186 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[0], fieldIDs[1],
187 0, 0, globalIndex) );
188 correctIndex = globalOffsets[localProc_] + 1;
189 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
190
191 if (globalIndex != correctIndex) {
192 ERReturn(-1);
193 }
194
195 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[1], fieldIDs[1],
196 0, 0, globalIndex) );
197 correctIndex = globalOffsets[localProc_] + 5;
198 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
199
200 if (globalIndex != correctIndex) {
201 ERReturn(-1);
202 }
203
204 CHK_ERR( vectorSpace.getGlobalIndex(idTypes[0], ids[len-1], fieldIDs[1],
205 0, 0, globalIndex) );
206 correctIndex = globalOffsets[localProc_] + 13;
207 if (localProc_ != 0) correctIndex -= 2*numDOFsPerID;
208
209 if (globalIndex != correctIndex) {
210 ERReturn(-1);
211 }
212
213 std::vector<int> globalIndices(ids.size()*numDOFsPerID);
214
215 CHK_ERR( vectorSpace.getGlobalIndices(ids.size(),
216 &(ids[0]),
217 idTypes[0], fieldIDs[0],
218 &globalIndices[0] ));
219
220 std::vector<int> idFieldIDs(ids.size(), fieldIDs[1]);
221 std::vector<int> idIDTypes(ids.size(), idTypes[0]);
222
223 CHK_ERR( vectorSpace.getGlobalIndices(ids.size(),
224 &ids[0],
225 idTypes[0], fieldIDs[0],
226 &globalIndices[0] ));
227
228 CHK_ERR( vectorSpace.getGlobalBlkIndices(ids.size(),
229 &ids[0],
230 idTypes[0],
231 &globalIndices[0] ));
232
233 CHK_ERR( vectorSpace.getGlobalIndices(ids.size(),
234 &ids[0],
235 &idIDTypes[0],
236 &idFieldIDs[0],
237 &globalIndices[0]) );
238
239 unsigned numFields = vectorSpace.getNumFields();
240 if (numFields != fieldIDs.size()) {
241 ERReturn(-1);
242 }
243
244 std::vector<int> testgetfields;
245 vectorSpace.getFields(testgetfields);
246 if (testgetfields.size() != numFields) {
247 ERReturn(-1);
248 }
249
250 if (fieldIDs != testgetfields) {
251 ERReturn(-1);
252 }
253
254 delete testdata;
255
256 return(0);
257}
258
259int test_VectorSpace::test2()
260{
261 //only run this test if asserts are enabled (i.e., if NDEBUG is not defined)
262 bool run_this_test = false;
263 assert( (run_this_test = true) == true);
264 if (!run_this_test) return(0);
265
266 //Initializes a fei::VectorSpace object using the same data as used
267 //above, but then adds a globally inconsistent shared-id
268 //before calling initComplete. The goal is to make sure that the shared-id
269 //safety-check in fei::VectorSpace will catch the incorrect shared-id
270 //data and return an error-code.
271
272 testData* testdata = new testData(localProc_, numProcs_);
273 std::vector<int>& idTypes = testdata->idTypes;
274 std::vector<int>& ids = testdata->ids;
275
277 fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, wrapper));
278
279 fei::SharedPtr<fei::VectorSpace> vectorSpacePtr =
280 test_VectorSpace::create_VectorSpace(comm_,
281 testdata, localProc_,
282 numProcs_,
283 true, //defineBothFields
284 true, //initSolnBothFields
285 "U2", factory);
286
287 if (localProc_ < numProcs_-1) {
288 int numSharingProcsPerID = 1;
289 int sharingProc = localProc_+1;
290
291 CHK_ERR( vectorSpacePtr->initSharedIDs(1, idTypes[0], &ids[0],
292 &numSharingProcsPerID,
293 &sharingProc) );
294 }
295
296 fei::ParameterSet paramset;
297 paramset.add(fei::Param("FEI_CHECK_SHARED_IDS", true));
298
299 vectorSpacePtr->setParameters(paramset);
300
301 //we just provided globally inconsistent shared-id info, so we expect the
302 //initComplete method to return(-1) after the shared-id safety-check fails,
303 //if we're running on more than one processor.
304 int err = vectorSpacePtr->initComplete();
305 if (numProcs_ > 1) if (err == 0) return(-1);
306 delete testdata;
307
308 return(0);
309}
310
311int test_VectorSpace::test3()
312{
313 testData* testdata = new testData(localProc_, numProcs_);
314
316 fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, wrapper));
317
318 fei::SharedPtr<fei::VectorSpace> vectorSpacePtr =
319 test_VectorSpace::create_VectorSpace(comm_,
320 testdata, localProc_, numProcs_,
321 true, //defineBothFields
322 true, //initSolnBothFields
323 "U3", factory);
324
326 factory->createVectorSpace(comm_, "U3copy");
327
328 fei::Param param("debugOutput", ".");
329 fei::ParameterSet paramset;
330 paramset.add(param);
331 copy->setParameters(paramset);
332
333
334 CHK_ERR( copy->addVectorSpace(vectorSpacePtr.get()) );
335
336 CHK_ERR( vectorSpacePtr->initComplete() );
337
338 CHK_ERR( copy->initComplete() );
339
340 std::vector<int> globalOffsets;
341 std::vector<int> globalOffsetsCopy;
342
343 vectorSpacePtr->getGlobalIndexOffsets(globalOffsets);
344
345 copy->getGlobalIndexOffsets(globalOffsetsCopy);
346
347 for(size_t i=0; i<globalOffsets.size(); ++i) {
348 if (globalOffsets[i] != globalOffsetsCopy[i]) {
349 ERReturn(-1);
350 }
351 }
352
353 CHK_ERR( copy->initComplete() );
354
355 copy->getGlobalIndexOffsets(globalOffsetsCopy);
356
357 if (globalOffsetsCopy[numProcs_] != globalOffsets[numProcs_]) {
358 ERReturn(-1);
359 }
360
361 delete testdata;
362
363 return(0);
364}
365
366int test_VectorSpace::test4()
367{
368 return(0);
369}
370
372test_VectorSpace::create_VectorSpace(MPI_Comm comm)
373{
374 int localProc = 0, numProcs = 1;
375#ifndef FEI_SER
376 MPI_Comm_rank(comm, &localProc);
377 MPI_Comm_size(comm, &numProcs);
378#endif
379
380 testData test_data(localProc, numProcs);
382
384 test_VectorSpace::create_VectorSpace(comm, &test_data, localProc, numProcs,
385 false, false, (const char*)0, factory);
386 int err = vspace->initComplete();
387 if (err != 0) {
388 FEI_COUT << "ERROR, failed to create valid fei::VectorSpace." << FEI_ENDL;
389 throw std::runtime_error("test_Vector::vector_test1: ERROR, failed to create valid fei::VectorSpace.");
390 }
391
392 return(vspace);
393}
394
396test_VectorSpace::create_VectorSpace(MPI_Comm comm,
397 testData* testdata,
398 int localProc,
399 int numProcs,
400 bool defineBothFields,
401 bool initSolnBothFields,
402 const char* name,
404 bool turnOnDebugOutput)
405{
406 //
407 //This function creates a VectorSpace object, then initializes it as follows:
408 //
409 //defineFields testdata->fieldIDs, testdata->fieldSizes
410 //defineIDTypes testdata->idTypes
411 //addDOFs testdata->ids with associated testdata->fieldIDs[0]
412 //addDOFs testdata->ids with testdata->fieldIDs[1] if bothFields
413 //initSharedIDs testdata->[shared-id-data] with testdata->idTypes[0]
414 //
416
417 if (factory.get() == NULL) {
418 vsptr.reset(new fei::VectorSpace(comm, name));
419 }
420 else {
421 vsptr = factory->createVectorSpace(comm, name);
422 }
423
424 fei::VectorSpace& vectorSpace = *vsptr;
425
426 std::vector<int>& fieldIDs = testdata->fieldIDs;
427 std::vector<int>& fieldSizes = testdata->fieldSizes;
428 std::vector<int>& idTypes = testdata->idTypes;
429 std::vector<int>& ids = testdata->ids;
430 std::vector<int>& sharedIDs = testdata->sharedIDs;
431 std::vector<int>& numSharingProcsPerID = testdata->numSharingProcsPerID;
432 std::vector<int>& sharingProcs = testdata->sharingProcs;
433
434 fei::ParameterSet paramset;
435 fei::Param param1("name", name);
436 paramset.add(param1);
437 if (turnOnDebugOutput) {
438 fei::Param param2("debugOutput", ".");
439 paramset.add(param2);
440 }
441
442 vectorSpace.setParameters(paramset);
443
444 int numFields = defineBothFields ? 2 : 1;
445 vectorSpace.defineFields(numFields,
446 &fieldIDs[0],
447 &fieldSizes[0]);
448
449 vectorSpace.defineIDTypes(idTypes.size(),
450 &idTypes[0]);
451
452 vectorSpace.addDOFs(fieldIDs[0],
453 idTypes[0],
454 ids.size(),
455 &ids[0]);
456
457 if (initSolnBothFields) {
458 vectorSpace.addDOFs(fieldIDs[1],
459 idTypes[0],
460 ids.size(),
461 &ids[0]);
462 }
463
464 vectorSpace.initSharedIDs(sharedIDs.size(),
465 idTypes[0],
466 sharedIDs.size() ? &sharedIDs[0] : 0,
467 numSharingProcsPerID.size() ? &numSharingProcsPerID[0] : 0,
468 sharingProcs.size() ? &sharingProcs[0] : 0
469 );
470
471 return(vsptr);
472}
void add(const Param &param, bool maintain_unique_keys=true)
void reset(T *p=0)
T * get() const
void getGlobalIndexOffsets(std::vector< int > &globalOffsets) const
int getGlobalIndices(int numIDs, const int *IDs, int idType, int fieldID, int *globalIndices)
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)
void getGlobalBlkIndexOffsets(std::vector< int > &globalBlkOffsets) const
int getGlobalBlkIndex(int idType, int ID, int &globalBlkIndex)
int addDOFs(int fieldID, int idType, int numIDs, const int *IDs)
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)
int getGlobalBlkIndices(int numIDs, const int *IDs, int idType, int *globalBlkIndices)
int localProc(MPI_Comm comm)
std::ostream & console_out()
int numProcs(MPI_Comm comm)