FEI Version of the Day
Loading...
Searching...
No Matches
fei_MatrixGraph_Impl2.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#include <fei_sstream.hpp>
10
11#include <limits>
12#include <cmath>
13
14#include <fei_MatrixGraph_Impl2.hpp>
15
16#include <fei_utils.hpp>
17#include "fei_TemplateUtils.hpp"
18
19#include <fei_Pattern.hpp>
20#include <fei_LogManager.hpp>
21#include <fei_TemplateUtils.hpp>
22#include <fei_impl_utils.hpp>
23#include <snl_fei_Utils.hpp>
24#include <fei_FieldMask.hpp>
25#include <fei_Record.hpp>
26#include <snl_fei_RecordCollection.hpp>
27#include <fei_VectorSpace.hpp>
28#include <fei_ParameterSet.hpp>
29#include <fei_ostream_ops.hpp>
30#include <fei_Reducer.hpp>
31#include <fei_GraphReducer.hpp>
32#include <fei_ConnectivityBlock.hpp>
33#include <snl_fei_BlkSizeMsgHandler.hpp>
34#include <fei_Graph_Impl.hpp>
35#include <snl_fei_Constraint.hpp>
36
37#include <fei_EqnBuffer.hpp>
38#include <fei_EqnCommMgr.hpp>
39#include <SNL_FEI_Structure.hpp>
40
41#undef fei_file
42#define fei_file "fei_MatrixGraph.cpp"
43
44#include <fei_ErrMacros.hpp>
45
46namespace snl_fei {
47static unsigned getFieldSize(int fieldID,
48 fei::VectorSpace* space1,
49 fei::VectorSpace* space2)
50{
51 unsigned fieldsize = 0;
52 bool foundfield = false;
53 if (space1 != NULL) {
54 try {
55 fieldsize = space1->getFieldSize(fieldID);
56 foundfield = true;
57 }
58 catch (...) {
59 foundfield = false;
60 }
61 }
62
63 if (!foundfield && space2 != NULL) {
64 try {
65 fieldsize = space2->getFieldSize(fieldID);
66 }
67 catch (std::runtime_error& exc) {
68 std::string msg("snl_fei::getFieldSize: ");
69 msg += exc.what();
70 throw std::runtime_error(msg);
71 }
72 }
73
74 return(fieldsize);
75}
76
77}//namespace snl_fei
78
79//------------------------------------------------------------------------------
83 const char* name)
84{
86 columnSpace,
87 name));
88
89 return(mgptr);
90}
91
92//=====Constructor==============================================================
95 const char* name)
96 : comm_(rowSpace->getCommunicator()),
97 rowSpace_(rowSpace),
98 colSpace_(colSpace),
99 haveRowSpace_(false),
100 haveColSpace_(false),
101 symmetric_(false),
102 remotelyOwnedGraphRows_(NULL),
103 simpleProblem_(false),
104 blockEntryGraph_(false),
105 patterns_(),
106 connectivityBlocks_(),
107 arbitraryBlockCounter_(1),
108 sparseBlocks_(),
109 lagrangeConstraints_(),
110 penaltyConstraints_(),
111 slaveConstraints_(),
112 ptEqualBlk_(false),
113 newSlaveData_(false),
114 localNumSlaves_(0),
115 globalNumSlaves_(0),
116 D_(NULL),
117 g_(),
118 g_nonzero_(false),
119 reducer_(),
120 name_(),
121 dbgprefix_("MatGrph: "),
122 tmpIntArray1_(),
123 tmpIntArray2_(),
124 includeAllSlaveConstraints_(false)
125{
126 localProc_ = fei::localProc(comm_);
127 numProcs_ = fei::numProcs(comm_);
128
129 if (rowSpace_.get() == NULL) {
130 voidERReturn;
131 }
132 else haveRowSpace_ = true;
133
134 if (colSpace_.get() != NULL) haveColSpace_ = true;
135 else colSpace_ = rowSpace_;
136
137 setName(name);
138}
139
140//------------------------------------------------------------------------------
142{
143 fei::destroyValues(patterns_);
144 patterns_.clear();
145
146 fei::destroyValues(connectivityBlocks_);
147 connectivityBlocks_.clear();
148
149 int i, len = sparseBlocks_.size();
150 for(i=0; i<len; ++i) {
151 delete sparseBlocks_[i];
152 }
153
154 fei::destroyValues(lagrangeConstraints_);
155 lagrangeConstraints_.clear();
156
157 fei::destroyValues(penaltyConstraints_);
158 penaltyConstraints_.clear();
159
160 fei::destroyValues(slaveConstraints_);
161 slaveConstraints_.clear();
162}
163
164//------------------------------------------------------------------------------
166{
167 const fei::Param* param = 0;
168 fei::Param::ParamType ptype = fei::Param::BAD_TYPE;
169
170 param = params.get("FEI_OUTPUT_LEVEL");
171 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
172 if (ptype == fei::Param::STRING) {
174 log_manager.setOutputLevel(param->getStringValue().c_str());
175 setOutputLevel(fei::utils::string_to_output_level(param->getStringValue()));
176 }
177
178 param = params.get("FEI_LOG_EQN");
179 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
180 if (ptype == fei::Param::INT) {
181 addLogEqn(param->getIntValue());
182 }
183
184 param = params.get("FEI_LOG_ID");
185 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
186 if (ptype == fei::Param::INT) {
187 addLogID(param->getIntValue());
188 }
189
190 param = params.get("name");
191 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
192 if (ptype == fei::Param::STRING) {
193 setName(param->getStringValue().c_str());
194 }
195
196 param = params.get("BLOCK_GRAPH");
197 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
198 if (ptype != fei::Param::BAD_TYPE) {
199 blockEntryGraph_ = true;
200 }
201
202 param = params.get("BLOCK_MATRIX");
203 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
204 if (ptype != fei::Param::BAD_TYPE) {
205 if (ptype == fei::Param::BOOL) {
206 blockEntryGraph_ = param->getBoolValue();
207 }
208 else {
209 blockEntryGraph_ = true;
210 }
211 }
212
213 param = params.get("INCLUDE_ALL_SLAVE_CONSTRAINTS");
214 ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
215 if (ptype != fei::Param::BAD_TYPE) {
216 includeAllSlaveConstraints_ = true;
217 }
218
219}
220
221//----------------------------------------------------------------------------
223{
224 rowSpace_ = rowSpace;
225 haveRowSpace_ = true;
226 if (!haveColSpace_) symmetric_ = true;
227}
228
229//----------------------------------------------------------------------------
231{
232 colSpace_ = colSpace;
233 haveColSpace_ = true;
234 if (colSpace_ == rowSpace_) symmetric_ = true;
235 else symmetric_ = false;
236}
237
238//----------------------------------------------------------------------------
239int fei::MatrixGraph_Impl2::addPattern(fei::Pattern* pattern)
240{
241 std::map<int,fei::Pattern*>::iterator
242 p_iter = patterns_.begin(), p_iter_end = patterns_.end();
243
244 bool matches_existing_pattern = false;
245 int pattern_id = -1;
246 while(p_iter != p_iter_end && !matches_existing_pattern) {
247 if (*pattern == *(p_iter->second)) {
248 matches_existing_pattern = true;
249 pattern_id = p_iter->first;
250 }
251 else ++p_iter;
252 }
253
254 if (!matches_existing_pattern) {
255 pattern_id = patterns_.size();
256 patterns_.insert(std::make_pair(pattern_id, pattern));
257 }
258 else {
259 delete pattern;
260 }
261
262 return pattern_id;
263}
264
265//----------------------------------------------------------------------------
267 int idType)
268{
269 snl_fei::RecordCollection* rec_coll = NULL;
270 rowSpace_->getRecordCollection(idType, rec_coll);
271
272 fei::Pattern* pattern = new fei::Pattern(numIDs, idType, rec_coll);
273 return addPattern(pattern);
274}
275
276//----------------------------------------------------------------------------
278 int idType,
279 int fieldID)
280{
281 unsigned fieldsize = 0;
282 try {
283 fieldsize = snl_fei::getFieldSize(fieldID, rowSpace_.get(), colSpace_.get());
284 }
285 catch (std::runtime_error& exc) {
286 FEI_OSTRINGSTREAM osstr;
287 osstr << "fei::MatrixGraph_Impl2::definePattern caught error: "<<exc.what();
288 throw std::runtime_error(osstr.str());
289 }
290
291 snl_fei::RecordCollection* rec_coll = NULL;
292 rowSpace_->getRecordCollection(idType, rec_coll);
293
294 fei::Pattern* pattern =
295 new fei::Pattern(numIDs, idType, rec_coll, fieldID, fieldsize);
296 return addPattern(pattern);
297}
298
299//----------------------------------------------------------------------------
301 int idType,
302 const int* numFieldsPerID,
303 const int* fieldIDs)
304{
305 std::vector<int> fieldSizes;
306 try {
307 int offset = 0;
308 for(int i=0; i<numIDs; ++i) {
309 for(int j=0; j<numFieldsPerID[i]; ++j) {
310 fieldSizes.push_back(snl_fei::getFieldSize(fieldIDs[offset++],
311 rowSpace_.get(), colSpace_.get()));
312 }
313 }
314 }
315 catch (std::runtime_error& exc) {
316 FEI_OSTRINGSTREAM osstr;
317 osstr << "fei::MatrixGraph_Impl2::definePattern caught error: "<<exc.what();
318 throw std::runtime_error(osstr.str());
319 }
320
321 snl_fei::RecordCollection* rec_coll = NULL;
322 rowSpace_->getRecordCollection(idType, rec_coll);
323
324 fei::Pattern* pattern = new fei::Pattern(numIDs, idType, rec_coll,
325 numFieldsPerID, fieldIDs, &(fieldSizes[0]));
326 return addPattern(pattern);
327}
328
329//----------------------------------------------------------------------------
331 const int* idTypes,
332 const int* numFieldsPerID,
333 const int* fieldIDs)
334{
335 std::vector<int> fieldSizes;
336 try {
337 int offset = 0;
338 for(int i=0; i<numIDs; ++i) {
339 for(int j=0; j<numFieldsPerID[i]; ++j) {
340 fieldSizes.push_back(snl_fei::getFieldSize(fieldIDs[offset++],
341 rowSpace_.get(), colSpace_.get()));
342 }
343 }
344 }
345 catch (std::runtime_error& exc) {
346 FEI_OSTRINGSTREAM osstr;
347 osstr << "fei::MatrixGraph_Impl2::definePattern caught error: "<<exc.what();
348 throw std::runtime_error(osstr.str());
349 }
350
351 std::vector<snl_fei::RecordCollection*> recordCollections(numIDs);
352 for(int i=0; i<numIDs; ++i) {
353 rowSpace_->getRecordCollection(idTypes[i], recordCollections[i]);
354 }
355
356 fei::Pattern* pattern = new fei::Pattern(numIDs, idTypes, &recordCollections[0],
357 numFieldsPerID, fieldIDs, &(fieldSizes[0]));
358 return addPattern(pattern);
359}
360
361//------------------------------------------------------------------------------
363{
364 std::map<int,fei::Pattern*>::iterator
365 p_iter = patterns_.find(patternID);
366
367 if (p_iter == patterns_.end()) {
368 return NULL;
369 }
370
371 fei::Pattern* pattern = (*p_iter).second;
372 return pattern;
373}
374
375//------------------------------------------------------------------------------
377 int numConnectivityLists,
378 int patternID,
379 bool diagonal)
380{
381 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
382 (*output_stream_) <<dbgprefix_<< "initConnectivityBlock symm, blkID="
383 << blockID<<", num="<<numConnectivityLists<<", patternID="
384 <<patternID<<FEI_ENDL;
385 }
386
387 if (blockID < 0) {
388 fei::console_out() << "fei::MatrixGraph_Impl2::initConnectivityBlock: blockID ("
389 << blockID << ") must be non-negative." << FEI_ENDL;
390 ERReturn(-1);
391 }
392
393 std::map<int,fei::Pattern*>::iterator
394 p_iter = patterns_.find(patternID);
395
396 if (p_iter == patterns_.end()) {
397 ERReturn(-1);
398 }
399
400 fei::Pattern* pattern = (*p_iter).second;
401
402 if (getConnectivityBlock(blockID) != NULL) ERReturn(-1);
403
404 fei::ConnectivityBlock* cblock =
405 new fei::ConnectivityBlock(blockID, pattern, numConnectivityLists);
406
407 connectivityBlocks_.insert(std::pair<int,fei::ConnectivityBlock*>(blockID, cblock));
408 if (diagonal) cblock->setIsDiagonal(diagonal);
409
410 return(0);
411}
412
413//------------------------------------------------------------------------------
415 int patternID,
416 bool diagonal)
417{
418 int blockID = connectivityBlocks_.size();
419 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
420 (*output_stream_) <<dbgprefix_<< "initConnectivityBlock symm, blkID="
421 << blockID<<", num="<<numConnectivityLists<<", patternID="
422 <<patternID<<FEI_ENDL;
423 }
424
425 std::map<int,fei::Pattern*>::iterator
426 p_iter = patterns_.find(patternID);
427
428 if (p_iter == patterns_.end()) {
429 FEI_OSTRINGSTREAM osstr;
430 osstr << "fei::MatrixGraph_Impl2::initConnectivityBlock: ERROR, patternID "
431 << patternID << " not found.";
432 throw std::runtime_error(osstr.str());
433 }
434
435 fei::Pattern* pattern = (*p_iter).second;
436
437 fei::ConnectivityBlock* cblock =
438 new fei::ConnectivityBlock(blockID, pattern, numConnectivityLists);
439
440 connectivityBlocks_.insert(std::pair<int,fei::ConnectivityBlock*>(blockID, cblock));
441 if (diagonal) cblock->setIsDiagonal(diagonal);
442
443 return(blockID);
444}
445
446//------------------------------------------------------------------------------
448 int numConnectivityLists,
449 int rowPatternID,
450 int colPatternID)
451{
452 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
453 (*output_stream_) <<dbgprefix_<< "initConnectivityBlock nonsymm, blkID="
454 << blockID<<", num="<<numConnectivityLists<<", row-patternID="
455 <<rowPatternID<<", col-patternID="<<colPatternID<<FEI_ENDL;
456 }
457
458 if (blockID < 0) {
459 FEI_OSTRINGSTREAM osstr;
460 osstr << "fei::MatrixGraph_Impl2::initConnectivityBlock: blockID ("
461 << blockID << ") must be non-negative." << FEI_ENDL;
462 throw std::runtime_error(osstr.str());
463 }
464
465 std::map<int,fei::Pattern*>::iterator
466 rp_iter = patterns_.find(rowPatternID);
467 if (rp_iter == patterns_.end()) ERReturn(-1);
468
469 fei::Pattern* rowpattern = (*rp_iter).second;
470
471 std::map<int,fei::Pattern*>::iterator
472 cp_iter = patterns_.find(colPatternID);
473 if (cp_iter == patterns_.end()) ERReturn(-1);
474
475 fei::Pattern* colpattern = (*cp_iter).second;
476
477 if (getConnectivityBlock(blockID) != NULL) ERReturn(-1);
478
479 fei::ConnectivityBlock* cblock =
480 new fei::ConnectivityBlock(blockID, rowpattern,
481 colpattern, numConnectivityLists);
482
483 connectivityBlocks_.insert(std::pair<int,fei::ConnectivityBlock*>(blockID, cblock));
484
485 return(0);
486}
487
488//------------------------------------------------------------------------------
490 int connectivityID,
491 const int* connectedIdentifiers)
492{
493 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
494 (*output_stream_) <<dbgprefix_<< "initConnectivity blkID="
495 <<blockID<<", connID="<<connectivityID<<FEI_ENDL;
496 }
497
498 //for each connected identifier, get a pointer to its record from the
499 //solution-space.
500 //store those pointers in the appropriate connectivity-table, mapped to
501 //the user-provided connectivityID.
502
503 fei::ConnectivityBlock* connblk = getConnectivityBlock(blockID);
504 if (connblk == NULL) {
505 FEI_OSTRINGSTREAM osstr;
506 osstr << "fei::MatrixGraph_Impl2::initConnectivity ERROR, blockID " << blockID
507 << " doesn't correspond to an existing ConnectivityBlock.";
508 throw std::runtime_error(osstr.str());
509 }
510
511 fei::Pattern* pattern = connblk->getRowPattern();
512 int numIDs = pattern->getNumIDs();
513
514 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
515 for(int ii=0; ii<numIDs; ++ii) {
516 if (isLogID(connectedIdentifiers[ii])) {
517 FEI_OSTREAM& os = *output_stream_;
518 os << dbgprefix_<<"initConnectivity blockID " << blockID << ", connID "
519 << connectivityID << ": ";
520 for(int jj=0; jj<numIDs; ++jj) {
521 os << connectedIdentifiers[jj] << " ";
522 }
523 os << FEI_ENDL;
524 }
525 }
526 }
527
528 if (rowSpace_.get() == NULL) ERReturn(-1);
529
530 std::map<int,int>& connectivityIDs = connblk->getConnectivityIDs();
531
532 int idOffset = -1;
533 std::map<int,int>::iterator
534 iter = connectivityIDs.find(connectivityID);
535 if (iter == connectivityIDs.end()) {
536 idOffset = connectivityIDs.size();
537 connectivityIDs.insert(iter, std::make_pair(connectivityID,idOffset));
538 }
539 else {
540 idOffset = iter->second;
541 }
542
543 idOffset *= pattern->getNumIDs();
544
545 int* rlist = &(connblk->getRowConnectivities()[idOffset]);
546
547 CHK_ERR( getConnectivityRecords(pattern, rowSpace_.get(),
548 connectedIdentifiers, rlist) );
549
550 for(int i=0; i<numIDs; ++i) {
551 if (pattern->getNumFieldsPerID()[i] > 0) {
552 pattern->getRecordCollections()[i]->getRecordWithLocalID(rlist[i])->isInLocalSubdomain_ = true;
553 }
554 }
555
556 return(0);
557}
558
559//------------------------------------------------------------------------------
561 int numRows,
562 const int* rowIDs,
563 const int* rowOffsets,
564 const int* packedColumnIDs)
565{
567 rowIDs, rowOffsets);
568
569 int max_row_len = 0;
570 for(int i=0; i<numRows; ++i) {
571 int row_len = rowOffsets[i+1]-rowOffsets[i];
572 if (row_len > max_row_len) max_row_len = row_len;
573 }
574
575 int patternID = definePattern(max_row_len, idType);
576 fei::Pattern* pattern = getPattern(patternID);
577 block->setRowPattern(pattern);
578
579 sparseBlocks_.push_back(block);
580
581 int* row_records = &(block->getRowConnectivities()[0]);
582
583 CHK_ERR( getConnectivityRecords(rowSpace_.get(),
584 idType, numRows, rowIDs, row_records) );
585
586 fei::VectorSpace* colSpace = rowSpace_.get();
587 if (colSpace_.get() != NULL) {
588 colSpace = colSpace_.get();
589 }
590
591 int* col_records = &(block->getColConnectivities()[0]);
592
593 CHK_ERR( getConnectivityRecords(colSpace, idType, rowOffsets[numRows],
594 packedColumnIDs, col_records));
595
596 return(0);
597}
598
599//------------------------------------------------------------------------------
601 int numRows,
602 const int* rowIDs,
603 const int* rowLengths,
604 const int*const* columnIDs)
605{
607 rowIDs, rowLengths, true);
608
609 int max_row_len = 0;
610 for(int i=0; i<numRows; ++i) {
611 int row_len = rowLengths[i];
612 if (row_len > max_row_len) max_row_len = row_len;
613 }
614
615 int patternID = definePattern(max_row_len, idType);
616 fei::Pattern* pattern = getPattern(patternID);
617 block->setRowPattern(pattern);
618
619 sparseBlocks_.push_back(block);
620
621 int* row_records = &(block->getRowConnectivities()[0]);
622
623 CHK_ERR( getConnectivityRecords(rowSpace_.get(),
624 idType, numRows, rowIDs, row_records) );
625
626 fei::VectorSpace* colSpace = rowSpace_.get();
627 if (colSpace_.get() != NULL) {
628 colSpace = colSpace_.get();
629 }
630
631 int* col_records = &(block->getColConnectivities()[0]);
632
633 int offset = 0;
634 for(int i=0; i<numRows; ++i) {
635 CHK_ERR( getConnectivityRecords(colSpace, idType, rowLengths[i],
636 columnIDs[i], &(col_records[offset])));
637 offset += rowLengths[i];
638 }
639
640 return(0);
641}
642
643//------------------------------------------------------------------------------
644int fei::MatrixGraph_Impl2::getConnectivityRecords(fei::VectorSpace* vecSpace,
645 int idType,
646 int numIDs,
647 const int* IDs,
648 int* records)
649{
650 snl_fei::RecordCollection* collection = NULL;
651 CHK_ERR( vecSpace->getRecordCollection(idType, collection) );
652
653 for(int i=0; i<numIDs; ++i) {
654 records[i] = collection->getLocalID(IDs[i]);
655 if (records[i] == -1) {
656 CHK_ERR( vecSpace->addDOFs(idType, 1, &IDs[i], &records[i]) );
657 }
658 }
659
660 return(0);
661}
662
663//------------------------------------------------------------------------------
664int fei::MatrixGraph_Impl2::getConnectivityRecords(fei::VectorSpace* vecSpace,
665 int idType,
666 int fieldID,
667 int numIDs,
668 const int* IDs,
669 int* records)
670{
671 snl_fei::RecordCollection* collection = NULL;
672 CHK_ERR( vecSpace->getRecordCollection(idType, collection) );
673
674 for(int i=0; i<numIDs; ++i) {
675 records[i] = collection->getLocalID(IDs[i]);
676 if (records[i] == -1) {
677 CHK_ERR( vecSpace->addDOFs(fieldID, idType, 1, &IDs[i], &records[i]));
678 }
679 }
680
681 return(0);
682}
683
684//------------------------------------------------------------------------------
685int fei::MatrixGraph_Impl2::getConnectivityRecords(fei::Pattern* pattern,
686 fei::VectorSpace* vecSpace,
687 const int* connectedIdentifiers,
688 int* recordList)
689{
690 fei::Pattern::PatternType pType = pattern->getPatternType();
691 int i, numIDs = pattern->getNumIDs();
692 const int* numFieldsPerID = pattern->getNumFieldsPerID();
693 const int* fieldIDs = pattern->getFieldIDs();
694
695 if (pType == fei::Pattern::GENERAL) {
696 const int* idTypes = pattern->getIDTypes();
697
698 int fieldOffset = 0;
699 for(i=0; i<numIDs; ++i) {
700 int id = connectedIdentifiers[i];
701
702 for(int nf=0; nf<numFieldsPerID[i]; ++nf) {
703 CHK_ERR( vecSpace->addDOFs(fieldIDs[fieldOffset++],
704 idTypes[i], 1, &id,
705 &recordList[i]));
706 }
707 }
708 }
709 else if (pType == fei::Pattern::SINGLE_IDTYPE ||
710 pType == fei::Pattern::SIMPLE ||
711 pType == fei::Pattern::NO_FIELD) {
712
713 int idType = pattern->getIDTypes()[0];
714
715 switch(pType) {
716 case fei::Pattern::SINGLE_IDTYPE:
717 {
718 int fieldOffset = 0;
719 for(i=0; i<numIDs; ++i) {
720 int id = connectedIdentifiers[i];
721 for(int nf=0; nf<numFieldsPerID[i]; ++nf) {
722 CHK_ERR( vecSpace->addDOFs(fieldIDs[fieldOffset++],
723 idType, 1, &id,
724 &(recordList[i])));
725 }
726 }
727 }
728 break;
729 case fei::Pattern::SIMPLE:
730 {
731 CHK_ERR( vecSpace->addDOFs(fieldIDs[0], idType,
732 numIDs, connectedIdentifiers,
733 recordList) );
734 }
735 break;
736 case fei::Pattern::NO_FIELD:
737 CHK_ERR( vecSpace->addDOFs(idType, numIDs,
738 connectedIdentifiers, recordList));
739 break;
740 case fei::Pattern::GENERAL:
741 //Include a stub for this case to avoid compiler warning...
742 std::abort(); break;
743 }
744 }
745 else {
746 ERReturn(-1);
747 }
748
749 return(0);
750}
751
752//------------------------------------------------------------------------------
754 int connectivityID,
755 const int* rowConnectedIdentifiers,
756 const int* colConnectedIdentifiers)
757{
758 //for each connected identifier, get a pointer to its record from the
759 //solution-space.
760 //store those pointers in the appropriate connectivity-table, mapped to
761 //the user-provided connectivityID.
762
763 fei::ConnectivityBlock* connblk = getConnectivityBlock(blockID);
764 if (connblk == NULL) {
765 FEI_OSTRINGSTREAM osstr;
766 osstr << "fei::MatrixGraph_Impl2::initConnectivity ERROR, blockID " << blockID
767 << " doesn't correspond to an existing ConnectivityBlock.";
768 throw std::runtime_error(osstr.str());
769 }
770
771 fei::Pattern* pattern = connblk->getRowPattern();
772 int numIDs = pattern->getNumIDs();
773 fei::Pattern* colPattern = connblk->getColPattern();
774 int numColIDs = colPattern->getNumIDs();
775 if (rowSpace_.get() == NULL) {
776 ERReturn(-1);
777 }
778 if (colSpace_.get() == NULL) {
779 ERReturn(-1);
780 }
781
782 std::map<int,int>& connectivityIDs = connblk->getConnectivityIDs();
783
784 int i, idOffset = -1;
785 std::map<int,int>::iterator
786 iter = connectivityIDs.find(connectivityID);
787 if (iter == connectivityIDs.end()) {
788 idOffset = connectivityIDs.size();
789 connectivityIDs.insert(iter, std::make_pair(connectivityID,idOffset));
790 }
791 else {
792 idOffset = iter->second;
793 }
794 int* row_rlist =
795 &(connblk->getRowConnectivities()[idOffset*numIDs]);
796 int* col_rlist =
797 &(connblk->getColConnectivities()[idOffset*numColIDs]);
798
799 CHK_ERR( getConnectivityRecords(pattern, rowSpace_.get(),
800 rowConnectedIdentifiers, row_rlist) );
801
802 for(i=0; i<numIDs; ++i)
803 pattern->getRecordCollections()[i]->getRecordWithLocalID(row_rlist[i])->isInLocalSubdomain_ = true;
804 CHK_ERR( getConnectivityRecords(colPattern, colSpace_.get(),
805 colConnectedIdentifiers, col_rlist) );
806
807 for(i=0; i<numColIDs; ++i)
808 colPattern->getRecordCollections()[i]->getRecordWithLocalID(col_rlist[i])->isInLocalSubdomain_ = true;
809
810 return(0);
811}
812
813//------------------------------------------------------------------------------
815 const int* connectedIdentifiers)
816{
817 std::map<int,fei::Pattern*>::iterator
818 p_iter = patterns_.find(patternID);
819 if (p_iter == patterns_.end()) ERReturn(-1);
820
821 fei::Pattern* pattern = p_iter->second;
822
823 int blockID = -arbitraryBlockCounter_++;
824 fei::ConnectivityBlock* cblock = new fei::ConnectivityBlock(blockID, pattern, 1);
825
826 connectivityBlocks_.insert(std::pair<int,fei::ConnectivityBlock*>(blockID, cblock));
827
828 CHK_ERR( initConnectivity(blockID, 0, connectedIdentifiers) );
829 return(0);
830}
831
832//------------------------------------------------------------------------------
834 const int* rowConnectedIdentifiers,
835 int colPatternID,
836 const int* colConnectedIdentifiers)
837{
838 std::map<int,fei::Pattern*>::iterator
839 rp_iter = patterns_.find(rowPatternID);
840 if (rp_iter == patterns_.end()) ERReturn(-1);
841
842 fei::Pattern* rowpattern = (*rp_iter).second;
843
844 std::map<int,fei::Pattern*>::iterator
845 cp_iter = patterns_.find(colPatternID);
846 if (cp_iter == patterns_.end()) ERReturn(-1);
847
848 fei::Pattern* colpattern = (*cp_iter).second;
849
850 int blockID = -arbitraryBlockCounter_++;
851 fei::ConnectivityBlock* cblock = new fei::ConnectivityBlock(blockID, rowpattern,
852 colpattern, 1);
853
854 connectivityBlocks_.insert(std::pair<int,fei::ConnectivityBlock*>(blockID, cblock));
855
856 CHK_ERR( initConnectivity(blockID, 0,
857 rowConnectedIdentifiers,
858 colConnectedIdentifiers) );
859 return(0);
860}
861
862//------------------------------------------------------------------------------
864 int& numIndices)
865{
866 std::map<int,fei::Pattern*>::iterator
867 p_iter = patterns_.find(patternID);
868 if (p_iter == patterns_.end()) ERReturn(-1);
869 fei::Pattern* pattern = (*p_iter).second;
870
871 numIndices = pattern->getNumIndices();
872
873 return(0);
874}
875
876//------------------------------------------------------------------------------
878 const int* IDs,
879 std::vector<int>& indices)
880{
881 std::map<int,fei::Pattern*>::iterator
882 p_iter = patterns_.find(patternID);
883 if (p_iter == patterns_.end()) ERReturn(-1);
884
885 fei::Pattern* pattern = (*p_iter).second;
886
887 indices.resize(pattern->getNumIndices());
888
889 int numIDs = pattern->getNumIDs();
890 const int* idTypes = pattern->getIDTypes();
891 const int* numFieldsPerID = pattern->getNumFieldsPerID();
892 const int* fieldIDs = pattern->getFieldIDs();
893
894 int offset = 0, fieldOffset = 0;
895 for(int i=0; i<numIDs; ++i) {
896 for(int j=0; j<numFieldsPerID[i]; ++j) {
897 CHK_ERR( rowSpace_->getGlobalIndices(1, &(IDs[i]), idTypes[i],
898 fieldIDs[fieldOffset],
899 &(indices[offset])) );
900 offset += snl_fei::getFieldSize(fieldIDs[fieldOffset++],
901 rowSpace_.get(), colSpace_.get());
902 }
903 }
904
905 return(0);
906}
907
908//------------------------------------------------------------------------------
910 int fieldID,
911 int numRows,
912 const int* rowIDs,
913 const int* rowOffsets,
914 const int* packedColumnIDs)
915{
916 fei::ConnectivityBlock* block = new fei::ConnectivityBlock(fieldID, numRows,
917 rowIDs, rowOffsets);
918
919 sparseBlocks_.push_back(block);
920
921 int* row_records = &(block->getRowConnectivities()[0]);
922
923 CHK_ERR( getConnectivityRecords(rowSpace_.get(), idType, fieldID,
924 numRows, rowIDs, row_records) );
925
926 fei::VectorSpace* colSpace = rowSpace_.get();
927 if (colSpace_.get() != NULL) {
928 colSpace = colSpace_.get();
929 }
930
931 int* col_records = &(block->getColConnectivities()[0]);
932
933 CHK_ERR( getConnectivityRecords(colSpace, idType, fieldID, rowOffsets[numRows],
934 packedColumnIDs, col_records));
935
936 return(0);
937}
938
939//------------------------------------------------------------------------------
941 int constraintIDType,
942 int numIDs,
943 const int* idTypes,
944 const int* IDs,
945 const int* fieldIDs)
946{
947 if (rowSpace_.get() == NULL) ERReturn(-1);
948
949 ConstraintType* constraint = getLagrangeConstraint(constraintID);
950
951 CHK_ERR( rowSpace_->addDOFs(constraintIDType, 1, &constraintID) );
952
953 if (haveColSpace_) {
954 if (colSpace_.get() == NULL) {
955 ERReturn(-1);
956 }
957 CHK_ERR( colSpace_->addDOFs(constraintIDType,
958 1, &constraintID) );
959 }
960
961 ConstraintType* newconstraint =
962 new ConstraintType(constraintID, constraintIDType,
963 false, //isSlave
964 false, //isPenalty
965 numIDs, idTypes, IDs, fieldIDs,
966 0, //slaveOffset
967 0, //offsetIntoSlaveField
968 NULL, //weights
969 0.0, //rhsValue
970 rowSpace_.get());
971
972 if (constraint != NULL) {
973 if (*constraint != *newconstraint) {
974 delete constraint;
975 }
976 else {
977 delete newconstraint;
978 return(0);
979 }
980 }
981
982 lagrangeConstraints_[constraintID] = newconstraint;
983
984 return(0);
985}
986
987//------------------------------------------------------------------------------
989 int constraintIDType,
990 int numIDs,
991 const int* idTypes,
992 const int* IDs,
993 const int* fieldIDs)
994{
995 if (rowSpace_.get() == NULL) ERReturn(-1);
996
997 ConstraintType* constraint = getPenaltyConstraint(constraintID);
998
999 ConstraintType* newconstraint =
1000 new ConstraintType(constraintID, constraintIDType,
1001 false, //isSlave
1002 true, //isPenalty
1003 numIDs, idTypes, IDs, fieldIDs,
1004 0, //slaveOffset
1005 0, //offsetIntoSlaveField
1006 NULL, //weights
1007 0.0, //rhsValue
1008 rowSpace_.get());
1009
1010 if (constraint != NULL) {
1011 if (*constraint != *newconstraint) {
1012 delete constraint;
1013 }
1014 else {
1015 delete newconstraint;
1016 return(0);
1017 }
1018 }
1019
1020 penaltyConstraints_[constraintID] = newconstraint;
1021
1022 return(0);
1023}
1024
1025//------------------------------------------------------------------------------
1027{
1028 snl_fei::RecordCollection* collection = NULL;
1029 rowSpace_->getRecordCollection(idType, collection);
1030 if (collection == NULL) {
1031 throw std::runtime_error("fei::MatrixGraph_Impl2::hasSlaveDof: ERROR, unknown idType");
1032 }
1033
1034 fei::Record<int>* rec = collection->getRecordWithID(ID);
1035
1036 if (rec == NULL) {
1037 FEI_OSTRINGSTREAM osstr;
1038 osstr << "fei::MatrixGraph_Impl2::hasSlaveDof: ERROR, specified ID ("
1039 << ID << ") not found.";
1040 throw std::runtime_error(osstr.str());
1041 }
1042
1043 return( rec->hasSlaveDof() );
1044}
1045
1046//------------------------------------------------------------------------------
1048 const int* idTypes,
1049 const int* IDs,
1050 const int* fieldIDs,
1051 int offsetOfSlave,
1052 int offsetIntoSlaveField,
1053 const double* weights,
1054 double rhsValue)
1055{
1056 if (rowSpace_.get() == NULL) ERReturn(-1);
1057
1058 FEI_OSTRINGSTREAM idosstr;
1059 idosstr << IDs[offsetOfSlave] << fieldIDs[offsetOfSlave] << offsetIntoSlaveField;
1060 FEI_ISTRINGSTREAM idisstr(idosstr.str());
1061 int crID = IDs[offsetOfSlave];
1062 idisstr >> crID;
1063
1064 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
1065 FEI_OSTREAM& dbgos = *output_stream_;
1066
1067 dbgos << dbgprefix_<<"initSlaveConstraint crID=" << crID << ", slaveID="
1068 << IDs[offsetOfSlave];
1069
1070 if (output_level_ >= fei::FULL_LOGS) {
1071 dbgos << " { ";
1072 for(int ii=0; ii<numIDs; ++ii) {
1073 dbgos << "("<<IDs[ii] << ","<<fieldIDs[ii]<<") ";
1074 }
1075 dbgos << "}"<<FEI_ENDL;
1076 }
1077 else dbgos << FEI_ENDL;
1078 }
1079
1080 ConstraintType* constraint = getSlaveConstraint(crID);
1081
1082 ConstraintType* newconstraint =
1083 new ConstraintType(crID, 0,
1084 true, //isSlave
1085 false, //isPenalty
1086 numIDs, idTypes, IDs, fieldIDs,
1087 offsetOfSlave,
1088 offsetIntoSlaveField,
1089 weights,
1090 rhsValue,
1091 rowSpace_.get());
1092
1093 if (constraint != NULL) {
1094 if (*constraint != *newconstraint) {
1095 if (!constraint->structurallySame(*newconstraint)) {
1096 FEI_OSTRINGSTREAM osstr;
1097 osstr << "fei::MatrixGraph_Impl2::initSlaveConstraint: slave ID "<<IDs[offsetOfSlave]
1098 << " is already constrained, with different connectivity. Changing the"
1099 << " the structure of an existing constraint is not allowed.";
1100 throw std::runtime_error(osstr.str());
1101 }
1102 newSlaveData_ = true;
1103 delete constraint;
1104 }
1105 else {
1106 delete newconstraint;
1107 return(0);
1108 }
1109 }
1110 else {
1111 newSlaveData_ = true;
1112 }
1113
1114 slaveConstraints_[crID] = newconstraint;
1115
1116 return(0);
1117}
1118
1119//------------------------------------------------------------------------------
1120bool
1121fei::MatrixGraph_Impl2::newSlaveData()
1122{
1123 bool globalBool = false;
1124 fei::Allreduce(comm_, newSlaveData_, globalBool);
1125 newSlaveData_ = globalBool;
1126
1127 return(newSlaveData_);
1128}
1129
1130//------------------------------------------------------------------------------
1133{
1134 std::map<int,ConstraintType*>::iterator
1135 c_iter = lagrangeConstraints_.find(constraintID);
1136 if (c_iter == lagrangeConstraints_.end()) return(NULL);
1137
1138 return( (*c_iter).second );
1139}
1140
1141//------------------------------------------------------------------------------
1144{
1145 std::map<int,ConstraintType*>::iterator
1146 c_iter = slaveConstraints_.find(constraintID);
1147 if (c_iter == slaveConstraints_.end()) return(NULL);
1148
1149 return( (*c_iter).second );
1150}
1151
1152//------------------------------------------------------------------------------
1155{
1156 std::map<int,ConstraintType*>::iterator
1157 c_iter = penaltyConstraints_.find(constraintID);
1158 if (c_iter == penaltyConstraints_.end()) return(NULL);
1159
1160 return( (*c_iter).second );
1161}
1162
1163//------------------------------------------------------------------------------
1165{
1166 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
1167 (*output_stream_) <<dbgprefix_<< "initComplete" << FEI_ENDL;
1168 }
1169 if (rowSpace_.get() == NULL) {
1170 ERReturn(-1);
1171 }
1172 else {
1173 CHK_ERR( rowSpace_->initComplete() );
1174 }
1175
1176 if (haveColSpace_ && colSpace_.get() == NULL) ERReturn(-1);
1177 if (haveColSpace_) {
1178 CHK_ERR( colSpace_->initComplete() );
1179 }
1180
1181 if (rowSpace_->fieldMasks_.size() == 1 &&
1182 rowSpace_->fieldMasks_[0]->getNumFields() == 1) {
1183 simpleProblem_ = true;
1184 }
1185
1186 std::vector<int>& eqnNums = rowSpace_->getEqnNumbers();
1187 vspcEqnPtr_ = eqnNums.size() > 0 ? &eqnNums[0] : NULL;
1188
1189 //If there are slave constraints (on any processor), we need to create
1190 //a slave dependency matrix.
1191 localNumSlaves_ = slaveConstraints_.size();
1192 globalNumSlaves_ = 0;
1193 CHK_ERR( fei::GlobalSum(comm_, localNumSlaves_, globalNumSlaves_) );
1194
1195 if (globalNumSlaves_ > 0) {
1196 //we need to allocate the slave dependency and reduction matrices too.
1197 CHK_ERR( createSlaveMatrices() );
1198 }
1199
1200 return(0);
1201}
1202
1203//----------------------------------------------------------------------------
1204int fei::MatrixGraph_Impl2::createAlgebraicGraph(bool blockEntryGraph,
1205 fei::Graph* graph,
1206 bool gatherFromOverlap)
1207{
1208 //Now run the connectivity blocks, adding vertices (locations of nonzeros) to
1209 //the graph.
1210
1211 bool wasAlreadyBlockEntry = blockEntryGraph_;
1212
1213 blockEntryGraph_ = blockEntryGraph;
1214
1215 for(unsigned i=0; i<sparseBlocks_.size(); ++i) {
1216 fei::ConnectivityBlock& cblock = *(sparseBlocks_[i]);
1217 CHK_ERR( addBlockToGraph_sparse(graph, &cblock) );
1218 }
1219
1220 std::map<int,fei::ConnectivityBlock*>::const_iterator
1221 cb_iter = connectivityBlocks_.begin(),
1222 cb_end = connectivityBlocks_.end();
1223
1224 while(cb_iter != cb_end) {
1225 fei::ConnectivityBlock& cblock = *((*cb_iter).second);
1226
1227 bool symmetricBlock = cblock.isSymmetric();
1228
1229 if (symmetricBlock) {
1230 fei::Pattern* pattern = cblock.getRowPattern();
1231 if (pattern == NULL) {
1232 ERReturn(-1);
1233 }
1234
1235 fei::Pattern::PatternType pType = pattern->getPatternType();
1236 if (pType == fei::Pattern::GENERAL || pType == fei::Pattern::SINGLE_IDTYPE) {
1237 CHK_ERR( addBlockToGraph_multiField_symmetric(graph, &cblock) );
1238 }
1239 else if (pType == fei::Pattern::SIMPLE) {
1240 CHK_ERR( addBlockToGraph_singleField_symmetric(graph, &cblock) );
1241 }
1242 else if (pType == fei::Pattern::NO_FIELD) {
1243 CHK_ERR( addBlockToGraph_noField_symmetric(graph, &cblock) );
1244 }
1245 }
1246 else {
1247 fei::Pattern* pattern = cblock.getRowPattern();
1248 fei::Pattern* colpattern = cblock.getColPattern();
1249 if (pattern == NULL || colpattern == NULL) {
1250 ERReturn(-1);
1251 }
1252
1253 fei::Pattern::PatternType pType = pattern->getPatternType();
1254 fei::Pattern::PatternType colPType = colpattern->getPatternType();
1255 if (pType == fei::Pattern::SIMPLE && colPType == fei::Pattern::SIMPLE) {
1256 CHK_ERR( addBlockToGraph_singleField_nonsymmetric(graph, &cblock) );
1257 }
1258 else {
1259 CHK_ERR( addBlockToGraph_multiField_nonsymmetric(graph, &cblock) );
1260 }
1261 }
1262 ++cb_iter;
1263 }
1264
1265 CHK_ERR( addLagrangeConstraintsToGraph(graph) );
1266
1267 CHK_ERR( addPenaltyConstraintsToGraph(graph) );
1268
1269 if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
1270 FEI_OSTREAM& os = *output_stream_;
1271 os << dbgprefix_<<"# before graph->gatherFromOverlap()" << FEI_ENDL;
1272 CHK_ERR( graph->writeLocalGraph(os) );
1273 CHK_ERR( graph->writeRemoteGraph(os) );
1274 }
1275
1276 if (gatherFromOverlap) {
1277 CHK_ERR( graph->gatherFromOverlap() );
1278 }
1279
1280 if (blockEntryGraph) {
1281 CHK_ERR( exchangeBlkEqnSizes(graph) );
1282 }
1283
1284 if (output_level_ >= fei::FULL_LOGS && fei::numProcs(comm_)>1 &&
1285 output_stream_ != NULL && gatherFromOverlap) {
1286 FEI_OSTREAM& os = *output_stream_;
1287 os << dbgprefix_<<" after graph->gatherFromOverlap()" << FEI_ENDL;
1288 CHK_ERR( graph->writeLocalGraph(*output_stream_) );
1289 }
1290
1291 if (!wasAlreadyBlockEntry) {
1292 setIndicesMode(POINT_ENTRY_GRAPH);
1293 }
1294
1295 return(0);
1296}
1297
1298//----------------------------------------------------------------------------
1301 bool localRowGraph_includeSharedRows)
1302{
1304
1305 std::vector<int> globalOffsets;
1306
1307 if (blockEntryGraph) {
1308 if (reducer_.get() != NULL) {
1309 throw std::runtime_error("fei::MatrixGraph_Impl2::createGraph ERROR, can't specify both block-entry assembly and slave-constraint reduction.");
1310 }
1311
1312 rowSpace_->getGlobalBlkIndexOffsets(globalOffsets);
1313 }
1314 else {
1315 rowSpace_->getGlobalIndexOffsets(globalOffsets);
1316 }
1317
1318 if ((int)globalOffsets.size() < localProc_+2) return localRows;
1319
1320 int firstOffset = globalOffsets[localProc_];
1321 int lastOffset = globalOffsets[localProc_+1] - 1;
1322
1323 if (reducer_.get() != NULL) {
1324 std::vector<int>& reduced_eqns = reducer_->getLocalReducedEqns();
1325 if (!reduced_eqns.empty()) {
1326 firstOffset = reduced_eqns[0];
1327 lastOffset = reduced_eqns[reduced_eqns.size()-1];
1328 }
1329 }
1330
1331 fei::SharedPtr<fei::Graph> inner_graph(new fei::Graph_Impl(comm_, firstOffset, lastOffset) );
1333
1334 if (reducer_.get() == NULL) {
1335 graph = inner_graph;
1336 }
1337 else {
1338 graph.reset( new fei::GraphReducer(reducer_, inner_graph) );
1339 }
1340
1341 bool gatherFromOverlap = !localRowGraph_includeSharedRows;
1342 int err = createAlgebraicGraph(blockEntryGraph, graph.get(),
1343 gatherFromOverlap);
1344 if (err != 0) {
1345 return(localRows);
1346 }
1347
1348 localRows = fei::createSparseRowGraph(*(inner_graph->getLocalGraph()));
1349
1350 remotelyOwnedGraphRows_ = fei::createSparseRowGraph(inner_graph->getRemoteGraph());
1351
1352 remotelyOwnedGraphRows_->blockEntries = blockEntryGraph;
1353
1354 if (localRowGraph_includeSharedRows &&
1355 remotelyOwnedGraphRows_->rowNumbers.size() > 0) {
1356
1357 fei::SharedPtr<fei::SparseRowGraph> localPlusShared =
1358 snl_fei::mergeSparseRowGraphs(localRows.get(), remotelyOwnedGraphRows_.get());
1359 localRows = localPlusShared;
1360 }
1361
1362 return(localRows);
1363}
1364
1365//----------------------------------------------------------------------------
1366int fei::MatrixGraph_Impl2::exchangeBlkEqnSizes(fei::Graph* graph)
1367{
1368 //If point-equals-block (meaning block-equations are of size 1) or
1369 //if blockEntryGraph_ is false (meaning we haven't constructed a block-entry-
1370 //graph) then there is nothing to do in this function.
1371 if ( rowSpace_->getPointBlockMap()->ptEqualBlk() ) {
1372 return(0);
1373 }
1374
1375 snl_fei::BlkSizeMsgHandler blkHandler(rowSpace_.get(), graph, comm_);
1376
1377 CHK_ERR( blkHandler.do_the_exchange() );
1378
1379 return(0);
1380}
1381
1382//----------------------------------------------------------------------------
1384{
1385 //In this function we extract the dependency matrix of equation numbers
1386 //from the slave-constraints locally on each processor, then gather those
1387 //slave-equations onto each processor (so that each processor holds ALL
1388 //slave-equations).
1389 //Then we remove any couplings that may exist among the slave-equations and
1390 //finally create a FillableMat object to hold the final dependency matrix D_.
1391 //
1392
1393 if (!newSlaveData()) {
1394 return(0);
1395 }
1396
1397 std::vector<int>& eqnNums = rowSpace_->getEqnNumbers();
1398 vspcEqnPtr_ = eqnNums.size() > 0 ? &eqnNums[0] : NULL;
1399
1400 fei::SharedPtr<fei::FillableMat> local_D(new fei::FillableMat);
1402
1403 std::vector<int> masterEqns;
1404 std::vector<double> masterCoefs;
1405
1406 std::map<int,ConstraintType*>::const_iterator
1407 cr_iter = slaveConstraints_.begin(),
1408 cr_end = slaveConstraints_.end();
1409
1410 for(; cr_iter != cr_end; ++cr_iter) {
1411 ConstraintType* cr = (*cr_iter).second;
1412
1413 fei::Record<int>* slaveRecord = cr->getSlave();
1414 int slaveID = slaveRecord->getID();
1415 int slaveIDType = cr->getIDType();
1416 int slaveFieldID = cr->getSlaveFieldID();
1417 int offsetIntoSlaveField = cr->getOffsetIntoSlaveField();
1418 int slaveEqn = -1;
1419 CHK_ERR( rowSpace_->getGlobalIndex(slaveIDType, slaveID,
1420 slaveFieldID, 0, offsetIntoSlaveField,
1421 slaveEqn) );
1422
1423 std::vector<int>& masterRecords_vec = cr->getMasters();
1424 int* masterRecords = &masterRecords_vec[0];
1425 std::vector<snl_fei::RecordCollection*>& masterRecColls = cr->getMasterRecordCollections();
1426 std::vector<int>& masterIDTypes = cr->getMasterIDTypes();
1427 std::vector<int>& masterFieldIDs = cr->getMasterFieldIDs();
1428 std::vector<double>& masterWeights = cr->getMasterWeights();
1429 double* masterWtPtr = &masterWeights[0];
1430
1431 masterEqns.resize(masterWeights.size());
1432 masterCoefs.resize(masterWeights.size());
1433
1434 int* masterEqnsPtr = &(masterEqns[0]);
1435 double* masterCoefsPtr = &(masterCoefs[0]);
1436
1437 int offset = 0;
1438 for(size_t j=0; j<masterIDTypes.size(); ++j) {
1439 fei::Record<int>* masterRecord = masterRecColls[j]->getRecordWithLocalID(masterRecords[j]);
1440 int* eqnNumbers = vspcEqnPtr_+masterRecord->getOffsetIntoEqnNumbers();
1441 fei::FieldMask* mask = masterRecord->getFieldMask();
1442 int eqnOffset = 0;
1443 if (!simpleProblem_) {
1444 int err = mask->getFieldEqnOffset(masterFieldIDs[j], eqnOffset);
1445 if (err != 0) {
1446 throw std::runtime_error("FEI ERROR, failed to get eqn-offset for constraint master-field.");
1447 }
1448 }
1449
1450 unsigned fieldSize = rowSpace_->getFieldSize(masterFieldIDs[j]);
1451
1452 int eqn = eqnNumbers[eqnOffset];
1453 for(unsigned k=0; k<fieldSize; ++k) {
1454 masterEqnsPtr[offset++] = eqn+k;
1455 }
1456 }
1457
1458 double fei_eps = 1.e-49;
1459
1460 offset = 0;
1461 for(size_t jj=0; jj<masterEqns.size(); ++jj) {
1462 if ( includeAllSlaveConstraints_ ||
1463 (std::abs(masterWtPtr[jj]) > fei_eps) ) {
1464 masterCoefsPtr[offset] = masterWtPtr[jj];
1465 masterEqnsPtr[offset++] = masterEqnsPtr[jj];
1466 }
1467 }
1468
1469 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
1470 bool log_eqn = isLogEqn(slaveEqn);
1471 if (!log_eqn) {
1472 for(unsigned ii=0; ii<masterEqns.size(); ++ii) {
1473 if (isLogEqn(masterEqnsPtr[ii])) {
1474 log_eqn = true;
1475 break;
1476 }
1477 }
1478 }
1479
1480 if (log_eqn) {
1481 FEI_OSTREAM& os = *output_stream_;
1482 os << "createSlaveMatrices: " << slaveEqn << " = ";
1483 for(unsigned ii=0; ii<masterEqns.size(); ++ii) {
1484 if (ii!=0) os << "+ ";
1485 os << masterCoefsPtr[ii]<<"*"<<masterEqnsPtr[ii]<<" ";
1486 }
1487 os << FEI_ENDL;
1488 }
1489 }
1490
1491 local_D->putRow(slaveEqn, masterEqnsPtr, masterCoefsPtr, offset);
1492
1493 if ( includeAllSlaveConstraints_ ||
1494 (std::abs(cr->getRHSValue()) > fei_eps) ) {
1495 fei::put_entry(*local_g, slaveEqn, cr->getRHSValue());
1496 }
1497 }
1498
1499 if (D_.get() == NULL) {
1500 D_.reset(new fei::FillableMat);
1501 }
1502
1503 if (g_.get() == NULL) {
1504 g_.reset(new fei::CSVec);
1505 }
1506
1507#ifndef FEI_SER
1508 if (numProcs_ > 1) {
1509 fei::impl_utils::global_union(comm_, *local_D, *D_);
1510 fei::impl_utils::global_union(comm_, *local_g, *g_);
1511 }
1512 else {
1513 *D_ = *local_D;
1514 *g_ = *local_g;
1515 }
1516#else
1517 *D_ = *local_D;
1518 *g_ = *local_g;
1519#endif
1520
1521 if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
1522 (*output_stream_) << "# D_ (pre-removeCouplings):"<<FEI_ENDL;
1523 (*output_stream_) << *D_;
1524 }
1525
1526 int levelsOfCoupling = fei::impl_utils::remove_couplings(*D_);
1527 (void)levelsOfCoupling;
1528
1529 if (reducer_.get() == NULL) {
1530 reducer_.reset(new fei::Reducer(D_, g_, comm_));
1531
1532 std::vector<int> indices;
1533 rowSpace_->getIndices_Owned(indices);
1534
1535 reducer_->setLocalUnreducedEqns(indices);
1536 }
1537 else {
1538 reducer_->initialize();
1539 }
1540
1541 if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
1542 (*output_stream_) << "# D_:"<<FEI_ENDL;
1543 (*output_stream_) << *D_;
1544 }
1545
1546 double fei_eps = 1.e-49;
1547
1548 g_nonzero_ = false;
1549 for(size_t j=0; j<g_->size(); ++j) {
1550 double coef = g_->coefs()[j];
1551 if ( includeAllSlaveConstraints_ ||
1552 (std::abs(coef) > fei_eps) ) {
1553 g_nonzero_ = true;
1554 }
1555 }
1556
1557 newSlaveData_ = false;
1558
1559 return(0);
1560}
1561
1562//----------------------------------------------------------------------------
1565{
1566 return( reducer_ );
1567}
1568
1569//----------------------------------------------------------------------------
1572{
1573 return( remotelyOwnedGraphRows_ );
1574}
1575
1576//----------------------------------------------------------------------------
1577void fei::MatrixGraph_Impl2::getConstrainedIndices(std::vector<int>& crindices) const
1578{
1579 if (constrained_indices_.empty()) {
1580 crindices.clear();
1581 return;
1582 }
1583
1584 std::set<int>::const_iterator
1585 s_iter = constrained_indices_.begin(),
1586 s_end = constrained_indices_.end();
1587
1588 crindices.resize(constrained_indices_.size());
1589
1590 int offset = 0;
1591 for(; s_iter != s_end; ++s_iter) {
1592 crindices[offset++] = *s_iter;
1593 }
1594}
1595
1596//----------------------------------------------------------------------------
1597int fei::MatrixGraph_Impl2::addLagrangeConstraintsToGraph(fei::Graph* graph)
1598{
1599 std::vector<int> indices;
1600 std::map<int,ConstraintType*>::const_iterator
1601 cr_iter = lagrangeConstraints_.begin(),
1602 cr_end = lagrangeConstraints_.end();
1603
1604 constrained_indices_.clear();
1605
1606 while(cr_iter != cr_end) {
1607 ConstraintType* cr = (*cr_iter).second;
1608 int crID = cr->getConstraintID();
1609
1610 CHK_ERR( getConstraintConnectivityIndices(cr, indices) );
1611
1612 int numIndices = indices.size();
1613 int* indicesPtr = &(indices[0]);
1614
1615 for(int i=0; i<numIndices; ++i) {
1616 constrained_indices_.insert(indicesPtr[i]);
1617 }
1618
1619 int crEqnRow = -1, tmp;
1620 if (blockEntryGraph_) {
1621 CHK_ERR( rowSpace_->getGlobalBlkIndex(cr->getIDType(),
1622 crID, crEqnRow) );
1623 cr->setBlkEqnNumber(crEqnRow);
1624 CHK_ERR( rowSpace_->getGlobalIndex(cr->getIDType(),
1625 crID, tmp) );
1626 cr->setEqnNumber(tmp);
1627 }
1628 else {
1629 CHK_ERR( rowSpace_->getGlobalIndex(cr->getIDType(),
1630 crID, crEqnRow) );
1631 cr->setEqnNumber(crEqnRow);
1632 CHK_ERR( rowSpace_->getGlobalBlkIndex(cr->getIDType(),
1633 crID, tmp) );
1634 cr->setBlkEqnNumber(tmp);
1635 }
1636
1637 //now add the row contribution
1638 CHK_ERR( graph->addIndices(crEqnRow, numIndices, &(indices[0])) );
1639
1640 //Let's add a diagonal entry to the graph for this constraint-equation,
1641 //just in case we need to fiddle with this equation later (e.g. discard the
1642 //constraint equation and replace it with a dirichlet boundary condition...).
1643 CHK_ERR( graph->addIndices(crEqnRow, 1, &crEqnRow) );
1644
1645 //and finally, add the column contribution (which is simply the transpose
1646 //of the row contribution)
1647 for(int k=0; k<numIndices; ++k) {
1648 CHK_ERR( graph->addIndices(indicesPtr[k], 1, &crEqnRow) );
1649 }
1650
1651 ++cr_iter;
1652 }
1653
1654 return(0);
1655}
1656
1657//----------------------------------------------------------------------------
1660 std::vector<int>& globalIndices)
1661{
1662 std::vector<int>& fieldSizes = tmpIntArray1_;
1663 std::vector<int>& ones = tmpIntArray2_;
1664
1665 std::vector<int>& constrainedRecords = cr->getMasters();
1666 std::vector<int>& constrainedFieldIDs = cr->getMasterFieldIDs();
1667 std::vector<snl_fei::RecordCollection*>& recordCollections = cr->getMasterRecordCollections();
1668
1669 int len = constrainedRecords.size();
1670 fieldSizes.resize(len);
1671
1672 ones.assign(len, 1);
1673
1674 int numIndices = 0;
1675 if (blockEntryGraph_) {
1676 numIndices = len;
1677 }
1678 else {
1679 for(int j=0; j<len; ++j) {
1680 unsigned fieldSize = rowSpace_->getFieldSize(constrainedFieldIDs[j]);
1681 fieldSizes[j] = fieldSize;
1682 numIndices += fieldSize;
1683 }
1684 }
1685
1686 globalIndices.resize(numIndices);
1687
1688 int checkNum;
1689 CHK_ERR( getConnectivityIndices_multiField(&recordCollections[0],
1690 &constrainedRecords[0],
1691 len, &ones[0],
1692 &constrainedFieldIDs[0],
1693 &fieldSizes[0],
1694 numIndices, &globalIndices[0],
1695 checkNum) );
1696 if (numIndices != checkNum) {
1697 ERReturn(-1);
1698 }
1699
1700 return(0);
1701}
1702
1703//----------------------------------------------------------------------------
1704int fei::MatrixGraph_Impl2::addPenaltyConstraintsToGraph(fei::Graph* graph)
1705{
1706 std::vector<int> indices;
1707 std::map<int,ConstraintType*>::const_iterator
1708 cr_iter = penaltyConstraints_.begin(),
1709 cr_end = penaltyConstraints_.end();
1710
1711 while(cr_iter != cr_end) {
1712 ConstraintType* cr = (*cr_iter).second;
1713
1714 CHK_ERR( getConstraintConnectivityIndices(cr, indices) );
1715
1716 int numIndices = indices.size();
1717
1718 //now add the symmetric contributions to the graph
1719 CHK_ERR( graph->addSymmetricIndices(numIndices, &(indices[0])) );
1720
1721 ++cr_iter;
1722 }
1723
1724 return(0);
1725}
1726
1727//----------------------------------------------------------------------------
1729 bool& equivalent) const
1730{
1731 equivalent = false;
1732 int localCode = 1;
1733 int globalCode = 0;
1734
1735 int myNumBlocks = getNumConnectivityBlocks();
1736 int numBlocks = matrixGraph.getNumConnectivityBlocks();
1737
1738 if (myNumBlocks != numBlocks) {
1739 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1740 equivalent = globalCode > 0 ? false : true;
1741 return(0);
1742 }
1743
1744 if (numBlocks > 0) {
1745 std::vector<int> myBlockIDs;
1746 std::vector<int> blockIDs;
1747
1748 CHK_ERR( getConnectivityBlockIDs(myBlockIDs) );
1749 CHK_ERR( matrixGraph.getConnectivityBlockIDs(blockIDs) );
1750
1751 for(int i=0; i<numBlocks; ++i) {
1752 if (myBlockIDs[i] != blockIDs[i]) {
1753 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1754 equivalent = globalCode > 0 ? false : true;
1755 return(0);
1756 }
1757
1758 const fei::ConnectivityBlock* mycblock = getConnectivityBlock(myBlockIDs[i]);
1759 const fei::ConnectivityBlock* cblock = matrixGraph.getConnectivityBlock(blockIDs[i]);
1760
1761 int myNumLists = mycblock->getConnectivityIDs().size();
1762 int numLists = cblock->getConnectivityIDs().size();
1763
1764 if (myNumLists != numLists ||
1765 mycblock->isSymmetric() != cblock->isSymmetric()) {
1766 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1767 equivalent = globalCode > 0 ? false : true;
1768 return(0);
1769 }
1770
1771 int myNumIDsPerList = mycblock->getRowPattern()->getNumIDs();
1772 int numIDsPerList = cblock->getRowPattern()->getNumIDs();
1773
1774 if (myNumIDsPerList != numIDsPerList) {
1775 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1776 equivalent = globalCode > 0 ? false : true;
1777 return(0);
1778 }
1779 }
1780 }
1781
1782 int numMyLagrangeConstraints = getLocalNumLagrangeConstraints();
1783 int numMySlaveConstraints = getGlobalNumSlaveConstraints();
1784 int numLagrangeConstraints = matrixGraph.getLocalNumLagrangeConstraints();
1785 int numSlaveConstraints = matrixGraph.getGlobalNumSlaveConstraints();
1786
1787 if (numMyLagrangeConstraints != numLagrangeConstraints ||
1788 numMySlaveConstraints != numSlaveConstraints) {
1789 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1790 equivalent = globalCode > 0 ? false : true;
1791 return(0);
1792 }
1793
1794 localCode = 0;
1795 CHK_ERR( fei::GlobalMax(comm_, localCode, globalCode) );
1796 equivalent = globalCode > 0 ? false : true;
1797 return(0);
1798}
1799
1800//----------------------------------------------------------------------------
1802{
1803 return(connectivityBlocks_.size());
1804}
1805
1806//----------------------------------------------------------------------------
1807int fei::MatrixGraph_Impl2::getConnectivityBlockIDs(std::vector<int>& blockIDs) const
1808{
1809 blockIDs.resize(connectivityBlocks_.size());
1810
1811 std::map<int,fei::ConnectivityBlock*>::const_iterator
1812 cdb_iter = connectivityBlocks_.begin();
1813
1814 for(size_t i=0; i<blockIDs.size(); ++i, ++cdb_iter) {
1815 int blockID = (*cdb_iter).first;
1816 blockIDs[i] = blockID;
1817 }
1818
1819 return(0);
1820}
1821
1822//----------------------------------------------------------------------------
1824{
1825 const fei::ConnectivityBlock* cblock = getConnectivityBlock(blockID);
1826 if (cblock == NULL) return(-1);
1827
1828 const fei::Pattern* pattern = cblock->getRowPattern();
1829 return(pattern->getNumIDs());
1830}
1831
1832//----------------------------------------------------------------------------
1834{
1835 const fei::ConnectivityBlock* cblock = getConnectivityBlock(blockID);
1836 if (cblock == NULL) return(-1);
1837
1838 const fei::Pattern* pattern = cblock->getRowPattern();
1839 return( blockEntryGraph_ ?
1840 pattern->getNumIDs() : pattern->getNumIndices());
1841}
1842
1843//----------------------------------------------------------------------------
1845 int& numRowIndices,
1846 int& numColIndices)
1847{
1848 fei::ConnectivityBlock* cblock = getConnectivityBlock(blockID);
1849 if (cblock == NULL) return(-1);
1850
1851 fei::Pattern* pattern = cblock->getRowPattern();
1852 numRowIndices = pattern->getNumIndices();
1853 fei::Pattern* colpattern = cblock->isSymmetric() ?
1854 pattern : cblock->getColPattern();
1855 numColIndices = colpattern->getNumIndices();
1856
1857 return(0);
1858}
1859
1860//----------------------------------------------------------------------------
1862{
1863 std::map<int,fei::ConnectivityBlock*>::const_iterator
1864 c_iter = connectivityBlocks_.find(blockID);
1865 if (c_iter == connectivityBlocks_.end()) return(0);
1866
1867 return( (*c_iter).second );
1868}
1869
1870//----------------------------------------------------------------------------
1872{
1873 std::map<int,fei::ConnectivityBlock*>::iterator
1874 c_iter = connectivityBlocks_.find(blockID);
1875 if (c_iter == connectivityBlocks_.end()) return(0);
1876
1877 return( (*c_iter).second );
1878}
1879
1880//----------------------------------------------------------------------------
1882 int connectivityID,
1883 int indicesAllocLen,
1884 int* indices,
1885 int& numIndices)
1886{
1887 fei::ConnectivityBlock* cblock = getConnectivityBlock(blockID);
1888 if (cblock == NULL) return(-1);
1889
1890 std::vector<int>& eqnNums = rowSpace_->getEqnNumbers();
1891 vspcEqnPtr_ = eqnNums.size() > 0 ? &eqnNums[0] : NULL;
1892
1893 fei::Pattern* pattern = cblock->getRowPattern();
1894 numIndices = pattern->getNumIndices();
1895
1896 int len = numIndices > indicesAllocLen ? indicesAllocLen : numIndices;
1897
1898 int* records = cblock->getRowConnectivity(connectivityID);
1899 if (records == NULL) {
1900 ERReturn(-1);
1901 }
1902
1903 fei::Pattern::PatternType pType = pattern->getPatternType();
1904
1905 if (pType == fei::Pattern::GENERAL || pType == fei::Pattern::SINGLE_IDTYPE) {
1906 const int* numFieldsPerID = pattern->getNumFieldsPerID();
1907 const int* fieldIDs = pattern->getFieldIDs();
1908 int totalNumFields = pattern->getTotalNumFields();
1909
1910 std::vector<int> fieldSizes(totalNumFields);
1911
1912 for(int j=0; j<totalNumFields; ++j) {
1913 fieldSizes[j] = snl_fei::getFieldSize(fieldIDs[j], rowSpace_.get(),
1914 colSpace_.get());
1915 }
1916
1917 CHK_ERR( getConnectivityIndices_multiField(pattern->getRecordCollections(),
1918 records, pattern->getNumIDs(),
1919 numFieldsPerID,
1920 fieldIDs, &fieldSizes[0],
1921 len, indices, numIndices) );
1922 }
1923 else if (pType == fei::Pattern::SIMPLE) {
1924 const int* fieldIDs = pattern->getFieldIDs();
1925
1926 int fieldID = fieldIDs[0];
1927 unsigned fieldSize = snl_fei::getFieldSize(fieldID,
1928 rowSpace_.get(),
1929 colSpace_.get());
1930
1931 CHK_ERR( getConnectivityIndices_singleField(pattern->getRecordCollections(),
1932 records, pattern->getNumIDs(),
1933 fieldID, fieldSize,
1934 len, indices, numIndices) );
1935 }
1936 else if (pType == fei::Pattern::NO_FIELD) {
1937 CHK_ERR( getConnectivityIndices_noField(pattern->getRecordCollections(),
1938 records, pattern->getNumIDs(),
1939 len, indices, numIndices) );
1940 }
1941
1942 return(0);
1943}
1944
1945//----------------------------------------------------------------------------
1947 int connectivityID,
1948 int rowIndicesAllocLen,
1949 int* rowIndices,
1950 int& numRowIndices,
1951 int colIndicesAllocLen,
1952 int* colIndices,
1953 int& numColIndices)
1954{
1955 fei::ConnectivityBlock* cblock = getConnectivityBlock(blockID);
1956 if (cblock == NULL) return(-1);
1957
1958 std::vector<int>& eqnNums = rowSpace_->getEqnNumbers();
1959 vspcEqnPtr_ = eqnNums.size() > 0 ? &eqnNums[0] : NULL;
1960
1961 fei::Pattern* pattern = cblock->getRowPattern();
1962 numRowIndices = pattern->getNumIndices();
1963
1964 int len = numRowIndices > rowIndicesAllocLen ?
1965 rowIndicesAllocLen : numRowIndices;
1966
1967 int* records = cblock->getRowConnectivity(connectivityID);
1968 if (records == NULL) {
1969 ERReturn(-1);
1970 }
1971
1972 fei::Pattern::PatternType pType = pattern->getPatternType();
1973
1974 if (pType == fei::Pattern::GENERAL || pType == fei::Pattern::SINGLE_IDTYPE) {
1975 const int* numFieldsPerID = pattern->getNumFieldsPerID();
1976 const int* fieldIDs = pattern->getFieldIDs();
1977 int totalNumFields = pattern->getTotalNumFields();
1978
1979 std::vector<int> fieldSizes(totalNumFields);
1980
1981 for(int j=0; j<totalNumFields; ++j) {
1982 fieldSizes[j] = snl_fei::getFieldSize(fieldIDs[j], rowSpace_.get(),
1983 colSpace_.get());
1984 }
1985
1986 CHK_ERR( getConnectivityIndices_multiField(pattern->getRecordCollections(),
1987 records, pattern->getNumIDs(),
1988 numFieldsPerID,
1989 fieldIDs, &fieldSizes[0],
1990 len, rowIndices, numRowIndices) );
1991 }
1992 else if (pType == fei::Pattern::SIMPLE) {
1993 const int* fieldIDs = pattern->getFieldIDs();
1994
1995 int fieldID = fieldIDs[0];
1996 unsigned fieldSize = snl_fei::getFieldSize(fieldID, rowSpace_.get(),
1997 colSpace_.get());
1998
1999 CHK_ERR( getConnectivityIndices_singleField(pattern->getRecordCollections(),
2000 records, pattern->getNumIDs(),
2001 fieldID, fieldSize,
2002 len, rowIndices, numRowIndices) );
2003 }
2004 else if (pType == fei::Pattern::NO_FIELD) {
2005 CHK_ERR( getConnectivityIndices_noField(pattern->getRecordCollections(),
2006 records, pattern->getNumIDs(),
2007 len, rowIndices, numRowIndices) );
2008 }
2009
2010 fei::Pattern* colpattern = cblock->isSymmetric() ?
2011 pattern : cblock->getColPattern();
2012 pType = colpattern->getPatternType();
2013 numColIndices = colpattern->getNumIndices();
2014 len = numColIndices > colIndicesAllocLen ?
2015 colIndicesAllocLen : numColIndices;
2016
2017 if (!cblock->isSymmetric()) {
2018 records = cblock->getColConnectivity(connectivityID);
2019 }
2020 if (records == NULL) {
2021 return(-1);
2022 }
2023
2024 if (pType == fei::Pattern::GENERAL || pType == fei::Pattern::SINGLE_IDTYPE) {
2025 const int* numFieldsPerID = colpattern->getNumFieldsPerID();
2026 const int* fieldIDs = colpattern->getFieldIDs();
2027 int totalNumFields = colpattern->getTotalNumFields();
2028
2029 std::vector<int> fieldSizes(totalNumFields);
2030
2031 for(int j=0; j<totalNumFields; ++j) {
2032 fieldSizes[j] = snl_fei::getFieldSize(fieldIDs[j], rowSpace_.get(),
2033 colSpace_.get());
2034 }
2035
2036 CHK_ERR( getConnectivityIndices_multiField(colpattern->getRecordCollections(),
2037 records, colpattern->getNumIDs(),
2038 numFieldsPerID,
2039 fieldIDs, &fieldSizes[0],
2040 len, colIndices, numColIndices) );
2041 }
2042 else if (pType == fei::Pattern::SIMPLE) {
2043 const int* fieldIDs = colpattern->getFieldIDs();
2044
2045 int fieldID = fieldIDs[0];
2046 unsigned fieldSize = snl_fei::getFieldSize(fieldID, rowSpace_.get(),
2047 colSpace_.get());
2048
2049 CHK_ERR( getConnectivityIndices_singleField(colpattern->getRecordCollections(),
2050 records, colpattern->getNumIDs(),
2051 fieldID, fieldSize,
2052 len, colIndices, numColIndices) );
2053 }
2054 else if (pType == fei::Pattern::NO_FIELD) {
2055 CHK_ERR( getConnectivityIndices_noField(colpattern->getRecordCollections(),
2056 records, colpattern->getNumIDs(),
2057 len, colIndices, numColIndices) );
2058 }
2059
2060 return(0);
2061}
2062
2063//----------------------------------------------------------------------------
2065{
2066 return( lagrangeConstraints_.size() );
2067}
2068
2069//----------------------------------------------------------------------------
2070int fei::MatrixGraph_Impl2::addBlockToGraph_multiField_symmetric(fei::Graph* graph,
2071 fei::ConnectivityBlock* cblock)
2072{
2073 fei::Pattern* pattern = cblock->getRowPattern();
2074 int j;
2075 int numIDs = pattern->getNumIDs();
2076 int numIndices = blockEntryGraph_ ?
2077 pattern->getNumIDs() : pattern->getNumIndices();
2078 int checkNumIndices = numIndices;
2079 std::vector<int> indices(numIndices);
2080 int* indicesPtr = &indices[0];
2081
2082 const int* numFieldsPerID = pattern->getNumFieldsPerID();
2083 const int* fieldIDs = pattern->getFieldIDs();
2084 int totalNumFields = pattern->getTotalNumFields();
2085
2086 std::vector<int> fieldSizes(totalNumFields);
2087
2088 for(j=0; j<totalNumFields; ++j) {
2089 fieldSizes[j] = snl_fei::getFieldSize(fieldIDs[j], rowSpace_.get(),
2090 colSpace_.get());
2091 }
2092
2093 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2094 std::vector<int>& values = cblock->getRowConnectivities();
2095
2096 std::map<int,int>::iterator
2097 c_iter = connIDs.begin(),
2098 c_end = connIDs.end();
2099
2100 for(; c_iter != c_end; ++c_iter) {
2101 int offset = c_iter->second;
2102 int* records = &values[offset*numIDs];
2103
2104 CHK_ERR( getConnectivityIndices_multiField(pattern->getRecordCollections(),
2105 records, numIDs,
2106 numFieldsPerID,
2107 fieldIDs,
2108 &fieldSizes[0],
2109 numIndices,
2110 indicesPtr,
2111 checkNumIndices) );
2112
2113 if (checkNumIndices != numIndices) {
2114 ERReturn(-1);
2115 }
2116
2117 if (output_level_ > fei::BRIEF_LOGS && output_stream_ != NULL) {
2118 FEI_OSTREAM& os = *output_stream_;
2119
2120 const snl_fei::RecordCollection*const* recordColls = pattern->getRecordCollections();
2121 unsigned thisoffset = 0;
2122 for(int ii=0; ii<numIDs; ++ii) {
2123 const fei::Record<int>* record = recordColls[ii]->getRecordWithLocalID(records[ii]);
2124 int ID = record->getID();
2125 os << dbgprefix_<<"scatterIndices: ID=" <<ID<<": ";
2126 int num = pattern->getNumIndicesPerID()[ii];
2127 for(int jj=0; jj<num; ++jj) {
2128 os << indicesPtr[thisoffset++] << " ";
2129 }
2130 os << FEI_ENDL;
2131 }
2132
2133 for(int ii=0; ii<numIndices; ++ii) {
2134 if (isLogEqn(indicesPtr[ii])) {
2135 os << "adding Symm inds: ";
2136 for(int jj=0; jj<numIndices; ++jj) {
2137 os << indicesPtr[jj] << " ";
2138 }
2139 os << FEI_ENDL;
2140 break;
2141 }
2142 }
2143 }
2144
2145 //now we have the indices array, so we're ready to push it into
2146 //the graph container
2147 if (numIndices == numIDs || !cblock->isDiagonal()) {
2148 CHK_ERR( graph->addSymmetricIndices(numIndices, indicesPtr,
2149 cblock->isDiagonal()) );
2150 }
2151 else {
2152 int ioffset = 0;
2153 int* fieldSizesPtr = &fieldSizes[0];
2154 for(int i=0; i<numIDs; ++i) {
2155 CHK_ERR( graph->addSymmetricIndices(fieldSizesPtr[i], &(indicesPtr[ioffset])));
2156 ioffset += fieldSizesPtr[i];
2157 }
2158 }
2159 }
2160
2161 return(0);
2162}
2163
2164//----------------------------------------------------------------------------
2165int fei::MatrixGraph_Impl2::addBlockToGraph_multiField_nonsymmetric(fei::Graph* graph,
2166 fei::ConnectivityBlock* cblock)
2167{
2168 fei::Pattern* pattern = cblock->getRowPattern();
2169 fei::Pattern* colpattern = cblock->getColPattern();
2170 int j;
2171 int numIDs = pattern->getNumIDs();
2172 int numIndices = blockEntryGraph_ ? numIDs : pattern->getNumIndices();
2173 int checkNumIndices = numIndices;
2174 std::vector<int> indices(numIndices);
2175 int* indicesPtr = &indices[0];
2176
2177 int numColIDs = colpattern->getNumIDs();
2178 int numColIndices = blockEntryGraph_ ? numColIDs : colpattern->getNumIndices();
2179 int checkNumColIndices = numColIndices;
2180 std::vector<int> colindices(numColIndices);
2181 int* colindicesPtr = &colindices[0];
2182
2183 const int* numFieldsPerID = pattern->getNumFieldsPerID();
2184 const int* fieldIDs = pattern->getFieldIDs();
2185 int totalNumFields = pattern->getTotalNumFields();
2186
2187 const int* numFieldsPerColID = colpattern->getNumFieldsPerID();
2188 const int* colfieldIDs = colpattern->getFieldIDs();
2189 int totalNumColFields = colpattern->getTotalNumFields();
2190
2191 std::vector<int> fieldSizes(totalNumFields);
2192 std::vector<int> colfieldSizes(totalNumColFields);
2193
2194 for(j=0; j<totalNumFields; ++j) {
2195 fieldSizes[j] = snl_fei::getFieldSize(fieldIDs[j], rowSpace_.get(),
2196 colSpace_.get());
2197 }
2198 for(j=0; j<totalNumColFields; ++j) {
2199 colfieldSizes[j] = snl_fei::getFieldSize(colfieldIDs[j], rowSpace_.get(),
2200 colSpace_.get());
2201 }
2202
2203 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2204 std::vector<int>& rowrecords = cblock->getRowConnectivities();
2205 std::vector<int>& colrecords = cblock->getColConnectivities();
2206
2207 std::map<int,int>::iterator
2208 c_iter = connIDs.begin(),
2209 c_end = connIDs.end();
2210
2211 for(; c_iter != c_end; ++c_iter) {
2212 int offset = c_iter->second;
2213 int* records = &rowrecords[offset*numIDs];
2214
2215 int* colRecords = &colrecords[offset*numColIDs];
2216
2217 CHK_ERR( getConnectivityIndices_multiField(pattern->getRecordCollections(),
2218 records, numIDs,
2219 numFieldsPerID,
2220 fieldIDs,
2221 &fieldSizes[0],
2222 numIndices,
2223 indicesPtr,
2224 checkNumIndices) );
2225
2226 if (checkNumIndices != numIndices) {
2227 ERReturn(-1);
2228 }
2229
2230
2231 CHK_ERR( getConnectivityIndices_multiField(colpattern->getRecordCollections(),
2232 colRecords, numColIDs,
2233 numFieldsPerColID,
2234 colfieldIDs,
2235 &colfieldSizes[0],
2236 numColIndices,
2237 colindicesPtr,
2238 checkNumColIndices) );
2239
2240 if (checkNumColIndices != numColIndices) {
2241 ERReturn(-1);
2242 }
2243
2244 //now we have the indices arrays, so we're ready to push them into
2245 //the graph container
2246 for(int r=0; r<numIndices; ++r) {
2247 CHK_ERR( graph->addIndices(indicesPtr[r], numColIndices, colindicesPtr));
2248 }
2249 }
2250
2251 return(0);
2252}
2253
2254//----------------------------------------------------------------------------
2255int fei::MatrixGraph_Impl2::getConnectivityIndices_multiField(const snl_fei::RecordCollection*const* recordCollections,
2256 int* records,
2257 int numRecords,
2258 const int* numFieldsPerID,
2259 const int* fieldIDs,
2260 const int* fieldSizes,
2261 int indicesAllocLen,
2262 int* indices,
2263 int& numIndices)
2264{
2265 numIndices = 0;
2266 int fld_offset = 0;
2267
2268 for(int i=0; i<numRecords; ++i) {
2269 const fei::Record<int>* record = recordCollections[i]->getRecordWithLocalID(records[i]);
2270 if (record==NULL) continue;
2271
2272 if (blockEntryGraph_) {
2273 indices[numIndices++] = record->getNumber();
2274 continue;
2275 }
2276
2277 const fei::FieldMask* fieldMask = record->getFieldMask();
2278 int* eqnNumbers = vspcEqnPtr_ + record->getOffsetIntoEqnNumbers();
2279
2280 for(int nf=0; nf<numFieldsPerID[i]; ++nf) {
2281 int eqnOffset = 0;
2282 if (!simpleProblem_) {
2283 int err = fieldMask->getFieldEqnOffset(fieldIDs[fld_offset], eqnOffset);
2284 if (err != 0) {
2285 for(int fs=0; fs<fieldSizes[fld_offset]; ++fs) {
2286 indices[numIndices++] = -1;
2287 }
2288 continue;
2289 }
2290 }
2291
2292 for(int fs=0; fs<fieldSizes[fld_offset]; ++fs) {
2293 indices[numIndices++] = eqnNumbers[eqnOffset+fs];
2294 }
2295
2296 ++fld_offset;
2297 }
2298 }
2299
2300 return(0);
2301}
2302
2303//----------------------------------------------------------------------------
2304int fei::MatrixGraph_Impl2::addBlockToGraph_singleField_symmetric(fei::Graph* graph,
2305 fei::ConnectivityBlock* cblock)
2306{
2307 fei::Pattern* pattern = cblock->getRowPattern();
2308 int numIDs = pattern->getNumIDs();
2309 int numIndices = blockEntryGraph_ ? numIDs : pattern->getNumIndices();
2310 int checkNumIndices = numIndices;
2311 std::vector<int> indices(numIndices);
2312 int* indicesPtr = &indices[0];
2313
2314 const int* fieldIDs = pattern->getFieldIDs();
2315
2316 int fieldID = fieldIDs[0];
2317 unsigned fieldSize = snl_fei::getFieldSize(fieldID, rowSpace_.get(),
2318 colSpace_.get());
2319
2320 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2321 std::vector<int>& rowrecords = cblock->getRowConnectivities();
2322
2323 std::map<int,int>::iterator
2324 c_iter = connIDs.begin(),
2325 c_end = connIDs.end();
2326
2327 for(; c_iter != c_end; ++c_iter) {
2328 int offset = c_iter->second;
2329 int* records = &rowrecords[offset*numIDs];
2330
2331 CHK_ERR( getConnectivityIndices_singleField(pattern->getRecordCollections(),
2332 records, numIDs,
2333 fieldID, fieldSize,
2334 checkNumIndices,
2335 indicesPtr,
2336 numIndices) );
2337
2338 if (checkNumIndices != numIndices) {
2339 ERReturn(-1);
2340 }
2341
2342 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
2343 for(int ii=0; ii<numIndices; ++ii) {
2344 if (isLogEqn(indicesPtr[ii])) {
2345 FEI_OSTREAM& os = *output_stream_;
2346 os << "adding Symm inds: ";
2347 for(int jj=0; jj<numIndices; ++jj) {
2348 os << indicesPtr[jj] << " ";
2349 }
2350 os << FEI_ENDL;
2351 break;
2352 }
2353 }
2354 }
2355
2356 //now we have the indices array, so we're ready to push it into
2357 //the graph container
2358 if (numIndices == numIDs || !cblock->isDiagonal()) {
2359 CHK_ERR( graph->addSymmetricIndices(numIndices, indicesPtr,
2360 cblock->isDiagonal()));
2361 }
2362 else {
2363 int ioffset = 0;
2364 for(int i=0; i<numIDs; ++i) {
2365 CHK_ERR( graph->addSymmetricIndices(fieldSize, &(indicesPtr[ioffset])));
2366 ioffset += fieldSize;
2367 }
2368 }
2369 }
2370
2371 return(0);
2372}
2373
2374//----------------------------------------------------------------------------
2375int fei::MatrixGraph_Impl2::addBlockToGraph_singleField_nonsymmetric(fei::Graph* graph,
2376 fei::ConnectivityBlock* cblock)
2377{
2378 fei::Pattern* pattern = cblock->getRowPattern();
2379 fei::Pattern* colpattern = cblock->getColPattern();
2380 int numIDs = pattern->getNumIDs();
2381 int numIndices = blockEntryGraph_ ? numIDs : pattern->getNumIndices();
2382 int checkNumIndices = numIndices;
2383 std::vector<int> indices(numIndices);
2384 int* indicesPtr = &indices[0];
2385
2386 int numColIDs = colpattern->getNumIDs();
2387 int numColIndices = blockEntryGraph_ ?
2388 numColIDs : colpattern->getNumIndices();
2389 int checkNumColIndices = numColIndices;
2390 std::vector<int> colindices(numColIndices);
2391 int* colindicesPtr = &colindices[0];
2392
2393 const int* fieldIDs = pattern->getFieldIDs();
2394
2395 int rowFieldID = fieldIDs[0];
2396 int rowFieldSize = snl_fei::getFieldSize(rowFieldID, rowSpace_.get(),
2397 colSpace_.get());
2398
2399 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2400 std::vector<int>& rowrecords = cblock->getRowConnectivities();
2401 std::vector<int>& colrecords = cblock->getColConnectivities();
2402
2403 int colFieldID = colpattern->getFieldIDs()[0];
2404 int colFieldSize = snl_fei::getFieldSize(colFieldID, rowSpace_.get(),
2405 colSpace_.get());
2406
2407 std::map<int,int>::iterator
2408 c_iter = connIDs.begin(),
2409 c_end = connIDs.end();
2410
2411 for(; c_iter != c_end; ++c_iter) {
2412 int offset = c_iter->second;
2413 int* records = &rowrecords[offset*numIDs];
2414
2415 int* colRecords = &colrecords[offset*numColIDs];
2416
2417 if (blockEntryGraph_) {
2418 rowSpace_->getGlobalBlkIndicesL(numIDs, pattern->getRecordCollections(),
2419 records, checkNumIndices,
2420 indicesPtr, numIndices);
2421
2422 colSpace_->getGlobalBlkIndicesL(numColIDs, colpattern->getRecordCollections(),
2423 colRecords, checkNumColIndices,
2424 colindicesPtr, numColIndices);
2425 }
2426 else {
2427 rowSpace_->getGlobalIndicesL(numIDs, pattern->getRecordCollections(),
2428 records, rowFieldID,
2429 rowFieldSize, checkNumIndices,
2430 indicesPtr, numIndices);
2431
2432 colSpace_->getGlobalIndicesL(numColIDs, colpattern->getRecordCollections(),
2433 colRecords, colFieldID,
2434 colFieldSize, checkNumColIndices,
2435 colindicesPtr, numColIndices);
2436 }
2437
2438 if (checkNumIndices != numIndices || checkNumColIndices != numColIndices) {
2439 ERReturn(-1);
2440 }
2441
2442 for(int r=0; r<numIndices; ++r) {
2443 CHK_ERR( graph->addIndices(indicesPtr[r], numColIndices, colindicesPtr));
2444 }
2445 }
2446
2447 return(0);
2448}
2449
2450//----------------------------------------------------------------------------
2451int fei::MatrixGraph_Impl2::getConnectivityIndices_singleField(const snl_fei::RecordCollection*const* recordCollections,
2452 int* records,
2453 int numRecords,
2454 int fieldID,
2455 int fieldSize,
2456 int indicesAllocLen,
2457 int* indices,
2458 int& numIndices)
2459{
2460 numIndices = 0;
2461
2462 for(int i=0; i<numRecords; ++i) {
2463 if (numIndices == indicesAllocLen) break;
2464
2465 const fei::Record<int>* record = recordCollections[i]->getRecordWithLocalID(records[i]);
2466
2467 if (blockEntryGraph_) {
2468 indices[numIndices++] = record->getNumber();
2469 continue;
2470 }
2471
2472 int* eqnNumbers = vspcEqnPtr_+record->getOffsetIntoEqnNumbers();
2473
2474 int eqnOffset = 0;
2475 if (!simpleProblem_) {
2476 const fei::FieldMask* fieldMask = record->getFieldMask();
2477 int err = fieldMask->getFieldEqnOffset(fieldID, eqnOffset);
2478 if (err != 0) {
2479 indices[numIndices++] = -1;
2480 if (fieldSize > 1) {
2481 for(int fs=1; fs<fieldSize; ++fs) {
2482 indices[numIndices++] = -1;
2483 }
2484 }
2485 continue;
2486 }
2487 }
2488
2489 indices[numIndices++] = eqnNumbers[eqnOffset];
2490 if (fieldSize > 1) {
2491 for(int fs=1; fs<fieldSize; ++fs) {
2492 indices[numIndices++] = eqnOffset >= 0 ? eqnNumbers[eqnOffset+fs] : -1;
2493 }
2494 }
2495 }
2496
2497 return(0);
2498}
2499
2500//----------------------------------------------------------------------------
2501int fei::MatrixGraph_Impl2::getConnectivityIndices_noField(const snl_fei::RecordCollection*const* recordCollections,
2502 int* records,
2503 int numRecords,
2504 int indicesAllocLen,
2505 int* indices,
2506 int& numIndices)
2507{
2508 numIndices = 0;
2509
2510 for(int i=0; i<numRecords; ++i) {
2511
2512 const fei::Record<int>* record = recordCollections[i]->getRecordWithLocalID(records[i]);
2513 int* eqnNumbers = vspcEqnPtr_+record->getOffsetIntoEqnNumbers();
2514
2515 if (blockEntryGraph_) {
2516 indices[numIndices++] = record->getNumber();
2517 }
2518 else {
2519 indices[numIndices++] = eqnNumbers[0];
2520 }
2521 }
2522
2523 return(0);
2524}
2525
2526//----------------------------------------------------------------------------
2527int fei::MatrixGraph_Impl2::addBlockToGraph_noField_symmetric(fei::Graph* graph,
2528 fei::ConnectivityBlock* cblock)
2529{
2530 fei::Pattern* pattern = cblock->getRowPattern();
2531 int numIDs = pattern->getNumIDs();
2532 int numIndices = pattern->getNumIndices();
2533 std::vector<int> indices(numIndices);
2534 int* indicesPtr = &indices[0];
2535
2536 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2537 std::vector<int>& rowrecords = cblock->getRowConnectivities();
2538
2539 std::map<int,int>::iterator
2540 c_iter = connIDs.begin(),
2541 c_end = connIDs.end();
2542
2543 for(; c_iter != c_end; ++c_iter) {
2544 int offset = c_iter->second;
2545 int* records = &rowrecords[offset*numIDs];
2546
2547 int checkNumIndices;
2548 CHK_ERR( getConnectivityIndices_noField(pattern->getRecordCollections(), records, numIDs, numIndices,
2549 indicesPtr, checkNumIndices) );
2550
2551 if (checkNumIndices != numIndices) {
2552 ERReturn(-1);
2553 }
2554
2555 if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
2556 for(int ii=0; ii<numIndices; ++ii) {
2557 if (isLogEqn(indicesPtr[ii])) {
2558 FEI_OSTREAM& os = *output_stream_;
2559 os << "adding Symm inds: ";
2560 for(int jj=0; jj<numIndices; ++jj) {
2561 os << indicesPtr[jj] << " ";
2562 }
2563 os << FEI_ENDL;
2564 break;
2565 }
2566 }
2567 }
2568
2569 //now we have the indices array, so we're ready to push it into
2570 //the graph container
2571 CHK_ERR( graph->addSymmetricIndices(numIndices, indicesPtr) );
2572 }
2573
2574 return(0);
2575}
2576
2577//----------------------------------------------------------------------------
2578int fei::MatrixGraph_Impl2::addBlockToGraph_sparse(fei::Graph* graph,
2579 fei::ConnectivityBlock* cblock)
2580{
2581 std::vector<int> row_indices;
2582 std::vector<int> indices;
2583
2584 fei::Pattern* pattern = cblock->getRowPattern();
2585 const snl_fei::RecordCollection*const* recordCollections = pattern->getRecordCollections();
2586
2587 std::map<int,int>& connIDs = cblock->getConnectivityIDs();
2588 std::vector<int>& connOffsets = cblock->getConnectivityOffsets();
2589 int* rowrecords = &(cblock->getRowConnectivities()[0]);
2590 int* colrecords = &(cblock->getColConnectivities()[0]);
2591 bool haveField = cblock->haveFieldID();
2592 int fieldID = cblock->fieldID();
2593 int fieldSize = 1;
2594
2595 if (haveField) {
2596 fieldSize = snl_fei::getFieldSize(fieldID, rowSpace_.get(),
2597 colSpace_.get());
2598 }
2599
2600 std::map<int,int>::iterator
2601 c_iter = connIDs.begin(),
2602 c_end = connIDs.end();
2603
2604 for(; c_iter != c_end; ++c_iter) {
2605 int offset = c_iter->second;
2606 int rowlen = connOffsets[offset+1] - offset;
2607
2608 int* records = &(rowrecords[offset]);
2609
2610 int checkNumIndices;
2611 row_indices.resize(fieldSize);
2612
2613 if (haveField) {
2614 CHK_ERR( getConnectivityIndices_singleField(recordCollections, records, 1,
2615 fieldID, fieldSize,
2616 fieldSize,
2617 &row_indices[0],
2618 checkNumIndices) );
2619 }
2620 else {
2621 CHK_ERR( getConnectivityIndices_noField(recordCollections, records, 1, fieldSize,
2622 &row_indices[0],
2623 checkNumIndices) );
2624 }
2625
2626 indices.resize(fieldSize*rowlen);
2627 int* indicesPtr = &indices[0];
2628 int* crecords = &(colrecords[offset]);
2629
2630 if (haveField) {
2631 CHK_ERR( getConnectivityIndices_singleField(recordCollections, crecords, rowlen,
2632 fieldID, fieldSize,
2633 fieldSize*rowlen,
2634 indicesPtr,
2635 checkNumIndices) );
2636 }
2637 else {
2638 CHK_ERR( getConnectivityIndices_noField(recordCollections, crecords, rowlen, rowlen,
2639 indicesPtr, checkNumIndices) );
2640 }
2641 if (checkNumIndices != rowlen) {
2642 ERReturn(-1);
2643 }
2644
2645 //now we have the indices, so we're ready to push them into
2646 //the graph container
2647 int* row_ind_ptr = &row_indices[0];
2648 for(int i=0; i<fieldSize; ++i) {
2649 CHK_ERR( graph->addIndices(row_ind_ptr[i], fieldSize*rowlen,
2650 indicesPtr) );
2651 }
2652 }
2653
2654 return(0);
2655}
2656
2657//----------------------------------------------------------------------------
2658void fei::MatrixGraph_Impl2::setName(const char* name)
2659{
2660 if (name == NULL) return;
2661
2662 if (name_ == name) return;
2663
2664 name_ = name;
2665 dbgprefix_ = "MatGraph_"+name_+": ";
2666}
2667
2668//----------------------------------------------------------------------------
2670{
2671 if (mode == BLOCK_ENTRY_GRAPH) {
2672 blockEntryGraph_ = true;
2673 return;
2674 }
2675
2676 if (mode == POINT_ENTRY_GRAPH) {
2677 blockEntryGraph_ = false;
2678 return;
2679 }
2680
2681 voidERReturn;
2682}
2683
2684//----------------------------------------------------------------------------
2690
const fei::Pattern * getColPattern() const
std::vector< int > & getConnectivityOffsets()
const int * getColConnectivity(int ID) const
std::vector< int > & getColConnectivities()
const int * getRowConnectivity(int ID) const
std::vector< int > & getRowConnectivities()
const std::map< int, int > & getConnectivityIDs() const
const fei::Pattern * getRowPattern() const
int getFieldEqnOffset(int fieldID, int &offset) const
virtual int addIndices(int row, int len, const int *indices)=0
virtual int gatherFromOverlap()=0
virtual int writeLocalGraph(FEI_OSTREAM &os, bool debug=false, bool prefixLinesWithPoundSign=true)=0
virtual int addSymmetricIndices(int numIndices, int *indices, bool diagonal=false)=0
virtual int writeRemoteGraph(FEI_OSTREAM &os)=0
void setOutputLevel(OutputLevel olevel)
static LogManager & getLogManager()
virtual fei::SharedPtr< fei::MatrixGraph > createMatrixGraph(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > columnSpace, const char *name)
void setParameters(const fei::ParameterSet &params)
MatrixGraph_Impl2(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > colSpace, const char *name=NULL)
const fei::ConnectivityBlock * getConnectivityBlock(int blockID) const
ConstraintType * getLagrangeConstraint(int constraintID)
int getConnectivityIndices(int blockID, int connectivityID, int indicesAllocLen, int *indices, int &numIndices)
int getConnectivityNumIndices(int blockID) const
ConstraintType * getSlaveConstraint(int constraintID)
int getConnectivityBlockIDs(std::vector< int > &blockIDs) const
int definePattern(int numIDs, int idType)
void setColumnSpace(fei::SharedPtr< fei::VectorSpace > columnSpace)
ConstraintType * getPenaltyConstraint(int constraintID)
int initPenaltyConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)
int initConnectivityBlock(int blockID, int numConnectivityLists, int patternID, bool diagonal=false)
void setRowSpace(fei::SharedPtr< fei::VectorSpace > rowSpace)
bool hasSlaveDof(int ID, int idType)
int getPatternIndices(int patternID, const int *IDs, std::vector< int > &indices)
void getConstrainedIndices(std::vector< int > &crindices) const
int initSlaveConstraint(int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs, int offsetOfSlave, int offsetIntoSlaveField, const double *weights, double rhsValue)
int initLagrangeConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)
fei::Pattern * getPattern(int patternID)
int compareStructure(const fei::MatrixGraph &matrixGraph, bool &equivalent) const
fei::SharedPtr< fei::Reducer > getReducer()
fei::SharedPtr< FillableMat > getSlaveDependencyMatrix()
int getConstraintConnectivityIndices(ConstraintType *cr, std::vector< int > &globalIndices)
int initConnectivity(int blockID, int connectivityID, const int *connectedIdentifiers)
fei::SharedPtr< fei::SparseRowGraph > getRemotelyOwnedGraphRows()
fei::SharedPtr< fei::SparseRowGraph > createGraph(bool blockEntryGraph, bool localRowGraph_includeSharedRows=false)
int getNumIDsPerConnectivityList(int blockID) const
int getPatternNumIndices(int patternID, int &numIndices)
virtual const fei::ConnectivityBlock * getConnectivityBlock(int blockID) const =0
virtual int getLocalNumLagrangeConstraints() const =0
virtual int getGlobalNumSlaveConstraints() const =0
virtual int getConnectivityBlockIDs(std::vector< int > &blockIDs) const =0
virtual int getNumConnectivityBlocks() const =0
const std::string & getStringValue() const
ParamType getType() const
Definition fei_Param.hpp:98
bool getBoolValue() const
int getIntValue() const
const Param * get(const char *name) const
const int * getNumFieldsPerID() const
int getNumIDs() const
snl_fei::RecordCollection *const * getRecordCollections() const
PatternType getPatternType() const
const int * getNumIndicesPerID() const
const int * getIDTypes() const
int getNumIndices() const
int getTotalNumFields() const
const int * getFieldIDs() const
fei::FieldMask * getFieldMask()
GlobalIDType getID() const
int getOffsetIntoEqnNumbers() const
GlobalIDType getNumber() const
void reset(T *p=0)
unsigned getFieldSize(int fieldID)
int addDOFs(int fieldID, int idType, int numIDs, const int *IDs)
int getRecordCollection(int idType, snl_fei::RecordCollection *&records)
std::vector< int > & getMasterIDTypes()
std::vector< int > & getMasterFieldIDs()
bool structurallySame(const Constraint< RecordType > &rhs)
std::vector< snl_fei::RecordCollection * > & getMasterRecordCollections()
void setBlkEqnNumber(int blkEqn)
std::vector< int > & getMasters()
std::vector< double > & getMasterWeights()
fei::Record< int > * getRecordWithID(int ID)
fei::OutputLevel string_to_output_level(const std::string &str)
Definition fei_utils.cpp:58
int localProc(MPI_Comm comm)
void destroyValues(MAP_TYPE &map_obj)
fei::SharedPtr< fei::SparseRowGraph > createSparseRowGraph(const std::vector< snl_fei::RaggedTable< MAP_TYPE, SET_TYPE > * > &tables)
std::ostream & console_out()
int Allreduce(MPI_Comm comm, bool localBool, bool &globalBool)
int GlobalMax(MPI_Comm comm, std::vector< T > &local, std::vector< T > &global)
int numProcs(MPI_Comm comm)
int GlobalSum(MPI_Comm comm, std::vector< T > &local, std::vector< T > &global)
snl_fei::Constraint< fei::Record< int > * > ConstraintType