74                 const char* sOutputFileName,
 
   75                 const char* sInputFileName,
 
   76                 bool bAbortAfterInput)
 
   77 : 
DataManager(HP, OF, pS, dInitialTime, sOutputFileName, sInputFileName, bAbortAfterInput)
 
   79         silent_cerr(
"fatal error: you are building SchurDataManager, " 
   80                 "but mbdyn was compiled without MPI. " 
   81                 "Something weird is happening. " 
   82                 "Anyway, please compile with -DUSE_MPI " 
   83                 "to enable parallel solution" << std::endl);
 
  170                 const char* sOutputFileName,
 
  171                 const char* sInputFileName,
 
  172                 bool bAbortAfterInput)
 
  173 :
DataManager(HP, OF, pS, dInitialTime, sOutputFileName, sInputFileName, bAbortAfterInput),
 
  190 wgtflag(WEIGHT_VERTICES),
 
  192 Partitioner(PARTITIONER_DEFAULT),
 
  196 iTotalExpConnections(0)
 
  198         DEBUGCOUT(
"Entering SchurDataManager" << std::endl);
 
  201         DataComm = MBDynComm.Dup();
 
  203         MyRank = DataComm.Get_rank();
 
  211         const char* sKeyWords[] = {
 
  215                 "number" "of" "connections",
 
  225                 "aerodynamic" "element",
 
  286                         pedantic_cerr(
"no explicit connections declared " 
  287                                 "for this input file" << std::endl);
 
  292                 pedantic_cerr(
"no explicit connections declared " 
  293                         "for this input file" << std::endl);
 
  300                 silent_cerr(
"Error: \"begin: parallel;\" expected at line " 
  309                                 silent_cerr(
"Error: Weight flag expected " 
  324                                 } 
else if (HP.
IsKeyWord(
"communication")) {
 
  331                                         silent_cerr(
"invalid weight " 
  343                                         silent_cerr(
"the partition " 
  344                                                 "assignment is not complete; " 
  345                                                 "only " << i << 
" vertices " 
  355                                         silent_cerr(
"illegal value " 
  358                                                 "assignment[" << i << 
"] " 
  360                                                 << 
"; must be between 0 and " 
  372                                 silent_cerr(
"METIS partitioner not available; aborting..." << std::endl);
 
  380                                 silent_cerr(
"CHACO partitioner not available; aborting..." << std::endl);
 
  387                                 silent_cerr(
"unknown partitioner " 
  396                                 silent_cerr(
"Error: \"end: parallel;\" expected " 
  398                                         << 
"; aborting..." << std::endl);
 
  404                         silent_cerr(
"Unknown input at line " 
  409                 case NUMBEROFCONNECTIONS:
 
  420                                         silent_cerr(
"Error: <Connection> " 
  428                                 for (
int k = 0; k < 2; k++) {
 
  460                                                 case AERODYNAMICELEMENT:
 
  497                                                         silent_cerr(
"Error: invalid element type " 
  499                                                                 << 
"; aborting..." << std::endl);
 
  507                                                         silent_cerr(
"Error: label expected " 
  509                                                                 << 
" element type at line " 
  511                                                                 << 
"; aborting..." << std::endl);
 
  518                                                         silent_cerr(
"Error at line " 
  521                                                                 << 
"(" << j << 
") undefined; " 
  522                                                                 "aborting..." << std::endl);
 
  555                                                         silent_cerr(
"Error: invalid node type " 
  557                                                                 << 
"; aborting..." << std::endl);
 
  565                                                         silent_cerr(
"Error: label expected " 
  567                                                                 << 
"; aborting..." << std::endl);
 
  574                                                         silent_cerr(
"Error: at line " 
  577                                                                 << 
"(" << j << 
") undefined; " 
  578                                                                 << 
"aborting..." << std::endl);
 
  585                                                 silent_cerr(
"Unknown input at line " 
  587                                                         << 
"; aborting..." << std::endl);
 
  593                         ASSERT(iNumNodes == iNumElems);
 
  594                         if (iNumNodes + iNumElems !=  2*iTotalExpConnections) {
 
  595                                 silent_cerr(
"Error: total number of nodes and elements" 
  596                                         " in the parallel section at line " 
  598                                         << 
" is not consistent; aborting..." << std::endl);
 
  663                 silent_cerr(
"SchurDataManager::HowManyDofs: " 
  664                         "illegal request (" << 
unsigned(who) << 
")" 
  684                 silent_cerr(
"SchurDataManager::GetDofsList: " 
  685                         "illegal request (" << 
unsigned(who) << 
")" 
  701         DEBUGCOUT(
"Entering SchurDataManager::CreatePartition()" << std::endl);
 
  703         std::vector<std::vector<int> > adj(iTotVertices);       
 
  704         std::vector<int> CommWgts(iTotVertices);                
 
  707         std::vector<int> VertexWgts(iTotVertices);      
 
  720         memset(Vertices.
pXadj, 0, (iTotVertices + 1)*
sizeof(
int));
 
  726         int iMaxConnectionsPerVertex =
 
  729         int GravityPos = 0, AirPropPos = 0;
 
  731         unsigned int* pRotLab = NULL;
 
  736                 memset(pRotPos, 0, iNum*
sizeof(
int));
 
  737                 memset(pRotLab, 0, iNum*
sizeof(
unsigned int));
 
  742         std::vector<const Node *> connectedNodes;
 
  748                 InitList(Vertices.
pAdjncy, iTotVertices*iMaxConnectionsPerVertex, ADJ_UNDEFINED);
 
  751                 for (
unsigned int i = 0; i < 
iTotNodes; i++) {
 
  756                 Node** ppCurrNode = NULL;
 
  763                 for (ElemVecType::const_iterator pTmpEl = 
Elems.begin();
 
  764                         pTmpEl != 
Elems.end();
 
  775                                 pRotLab[iNumRt] = (*pTmpEl)->GetLabel();
 
  779                         (*pTmpEl)->GetConnectedNodes(connectedNodes);
 
  783                         (*pTmpEl)->WorkSpaceDim(&dimA, &dimB);
 
  784                         VertexWgts[iCount] = dimA * dimB;
 
  786                         CommWgts[iCount] = (*pTmpEl)->iGetNumDof()*(*pTmpEl)->iGetNumDof();
 
  788                         for (std::vector<const Node *>::const_iterator i = connectedNodes.begin();
 
  789                                 i != connectedNodes.end();
 
  792                                 Vertices.
pXadj[iCount + 1]++;
 
  796                                 unsigned label = (*i)->GetLabel();
 
  799                                 position = ppCurrNode - ppNodes;
 
  805                                 VertexWgts[iCount] += (*ppCurrNode)->
iGetNumDof();
 
  809                                 if ((iCount*iMaxConnectionsPerVertex) + Vertices.
pXadj[iCount + 1] - 1 < iTotVertices*iMaxConnectionsPerVertex) {
 
  810                                         Vertices.
pAdjncy[(iCount*iMaxConnectionsPerVertex)
 
  811                                                 + Vertices.
pXadj[iCount + 1] - 1] = position;
 
  815                                 Vertices.
pXadj[position + 1]++;
 
  816                                 if ((position*iMaxConnectionsPerVertex) + Vertices.
pXadj[position + 1] - 1 < iTotVertices*iMaxConnectionsPerVertex) {
 
  817                                         Vertices.
pAdjncy[(position*iMaxConnectionsPerVertex)
 
  818                                                 + Vertices.
pXadj[position + 1] - 1] = iCount;
 
  822                                 CommWgts[
position] = (*ppCurrNode)->iGetNumDof();
 
  826                 if (iTotalExpConnections != 0) {
 
  841                                 Vertices.
pXadj[iNdPos + 1]++;
 
  842                                 Vertices.
pXadj[iTotNodes + iElPos + 1]++;
 
  843                                 Vertices.
pAdjncy[(iNdPos*iMaxConnectionsPerVertex)
 
  844                                         + Vertices.
pXadj[iNdPos + 1] - 1] = iElPos + iTotNodes;
 
  845                                 Vertices.
pAdjncy[((iTotNodes + iElPos)*iMaxConnectionsPerVertex)
 
  846                                         + Vertices.
pXadj[iTotNodes + iElPos + 1] - 1] = iNdPos;
 
  852                         if (Vertices.
pXadj[i] > iMaxConnectionsPerVertex) {
 
  853                                 iMax = Vertices.
pXadj[i];
 
  857                 if (iMax <= iMaxConnectionsPerVertex) {
 
  861                 iMaxConnectionsPerVertex = iMax;
 
  871         Pack(Vertices.
pAdjncy, iTotVertices*iMaxConnectionsPerVertex);
 
  880                         std::ofstream ofPartition;
 
  882                         ofPartition.open(
"partition.debug");
 
  883                         ofPartition << 
"# METIS-like Input File" << std::endl
 
  884                                 << 
"# Column 1: computational weight" << std::endl
 
  885                                 << 
"# Column 2: communication weight" << std::endl
 
  886                                 << 
"# Total Vertices: " << iTotVertices << std::endl
 
  887                                 << 
"# Nodes" << std::endl;
 
  888                         for (
unsigned int i = 0; i < 
iTotNodes; i++) {
 
  889                                 ofPartition << 
"# " << i << 
"  Node Type: " 
  890                                         << 
"(" << 
psNodeNames[ppNodes[i]->GetNodeType()] << 
")" 
  891                                         << 
" Label: " << ppNodes[i]->GetLabel() << std::endl
 
  892                                         << VertexWgts[i] << 
" " << CommWgts[i];
 
  893                                 for (
int j = Vertices.
pXadj[i]; j < Vertices.
pXadj[i + 1]; j++) {
 
  894                                         ofPartition << 
" " << Vertices.
pAdjncy[j];
 
  896                                 ofPartition << std::endl;
 
  898                         ofPartition << 
"# Elements" << std::endl;
 
  899                         for (
unsigned int i = 0; i < iTotElem; i++) {
 
  900                                 ofPartition << 
"# " << i + iTotNodes << 
"  Element Type: " 
  902                                         << 
" Label: " << 
Elems[i]->GetLabel() << std::endl
 
  905                                 for (
int j = Vertices.
pXadj[i + iTotNodes];
 
  906                                                 j < Vertices.
pXadj[i + iTotNodes + 1];
 
  908                                         ofPartition << 
" " << Vertices.
pAdjncy[j];
 
  910                                 ofPartition << std::endl;
 
  927                                 vwgt = &VertexWgts[0];
 
  932                                 silent_cout(
"communication weights currently unsupported by CHACO; trying to emulate..." << std::endl);
 
  938                                 silent_cout(
"edges weights currently unsupported by MBDyn" << std::endl);
 
  956                         if (
wgtflag & WEIGHT_VERTICES) {
 
  957                                 vwgt = &VertexWgts[0];
 
  962                                 silent_cout(
"edges weights currently unsupported by MBDyn" << std::endl);
 
  981                         silent_cerr(
"no partition library is available; " 
  982                                 "partition assignments must be provided " 
  983                                 "manually" << std::endl);
 
  998         InterfNodes.
pXadj = NULL;
 
 1002         memset(InterfNodes.
pXadj, 0, DataCommSize*
sizeof(
int));
 
 1008                 memset(InterfNodes.
pAdjncy, 0, DataCommSize*iMaxInterfNodes*8*
sizeof(
int));
 
 1009                 InitList(InterfNodes.
pAdjncy, DataCommSize*iMaxInterfNodes*8, ADJ_UNDEFINED);
 
 1017                 for (
unsigned int i = 0; i < 
iTotNodes; i++) {
 
 1023                                 for (
int j = Vertices.
pXadj[i]; j < Vertices.
pXadj[i + 1]; j++) {
 
 1026                                                 InterfNodes.
pXadj[TmpPrc] += 1;
 
 1027                                                 int k = TmpPrc * iMaxInterfNodes * 2 + InterfNodes.
pXadj[TmpPrc] - 1;
 
 1028                                                 if (k < DataCommSize * iMaxInterfNodes * 2) {
 
 1046                         if (InterfNodes.
pXadj[i] > iMaxInterfNodes) {
 
 1047                                 iMax = InterfNodes.
pXadj[i];
 
 1051                 DataComm.Allreduce(&iMax, &iRMax, 1, MPI::INT, 
MPI::MAX);
 
 1053                 if (iMax <= iMaxInterfNodes) {
 
 1057                 iMaxInterfNodes = iMax;
 
 1063         MPI::Request* pRReq = NULL;
 
 1064         MPI::Request* pSReq = NULL;
 
 1069                 pRReq[i] = MPI::REQUEST_NULL;
 
 1070                 pSReq[i] = MPI::REQUEST_NULL;
 
 1073         const int DIM_TAG = 10;
 
 1080                         pRReq[i] = DataComm.Irecv(InterfNodes.
pAdjncy + iMaxInterfNodes + i*iMaxInterfNodes*2,
 
 1081                                         iMaxInterfNodes, MPI::INT, i, DIM_TAG);
 
 1082                         pSReq[i] = DataComm.Isend(InterfNodes.
pAdjncy + i*iMaxInterfNodes*2,
 
 1083                                         iMaxInterfNodes, MPI::INT, i, DIM_TAG);
 
 1127                                 unsigned int pTmpLab = pIV->
GetLabel();
 
 1129                                 for (
int k = 0; k < iNumRt; k++) {
 
 1130                                         if (pTmpLab == pRotLab[k]) {
 
 1135                                 for (
int j = 0; j < iMyTotRot; j++) {
 
 1136                                         if (pos != pMyRot[j]) {
 
 1137                                                 pMyRot[iMyTotRot] = pos;
 
 1143                                         pMyRot[iMyTotRot] = pos;
 
 1153                 for (
int i = 0; i < iNumRt; i++) {
 
 1154                         if (iMyTotRot == 0) {
 
 1155                                 color = MPI::UNDEFINED;
 
 1160                                         silent_cout(
"InducedVelocity(" << 
Elems[pRotPos[i]]->GetLabel() << 
")" 
 1161                                                 << 
" assigned to process " 
 1169                         pIndVelComm[i] = MBDynComm.Split(color, key);
 
 1171                         IndVelComm[i] = MPI::COMM_WORLD.Split(color, key);
 
 1175                         r->InitializeIndVelComm(pIndVelComm + i);
 
 1178                 for (
int i = 0; i < iMyTotRot; i++) {
 
 1190         for (
unsigned int i = 0; i < iTotVertices - 
iTotNodes; i++) {
 
 1202         MPI::Request::Waitall(DataCommSize, pRReq);
 
 1203         MPI::Request::Waitall(DataCommSize, pSReq);
 
 1206         bool bRecvFlag = 
false, bSentFlag = 
false;
 
 1209                         bRecvFlag = MPI::Request::Testall(DataCommSize, pRReq);
 
 1213                         bSentFlag = MPI::Request::Testall(DataCommSize, pSReq);
 
 1216                 if (bRecvFlag && bSentFlag) {
 
 1225         std::sort(InterfNodes.
pAdjncy,
 
 1226                         InterfNodes.
pAdjncy + iMaxInterfNodes * 2 * DataCommSize);
 
 1229         int* p = std::unique(InterfNodes.
pAdjncy,
 
 1230                         InterfNodes.
pAdjncy + iMaxInterfNodes * 2 * DataCommSize);
 
 1240         unsigned int* llabels = NULL;
 
 1254         int* pPosIntElems = NULL;
 
 1268                                 if (iIVIsMine == 1) {
 
 1270                                         for (std::vector<const Node *>::const_iterator j = connectedNodes.begin();
 
 1271                                                 j != connectedNodes.end();
 
 1274                                                 unsigned int* p = std::find(llabels, llabels + iNumIntNodes, (*j)->GetLabel());
 
 1275                                                 if (p != llabels + iNumIntNodes) {
 
 1295                                 for (std::vector<const Node *>::const_iterator j = connectedNodes.begin();
 
 1296                                         j != connectedNodes.end();
 
 1299                                         unsigned int* p = std::find(llabels, llabels + iNumIntNodes, (*j)->GetLabel());
 
 1300                                         if (p != llabels + iNumIntNodes) {
 
 1338                 if (TmpDofNum != 0) {
 
 1339                         if (i != pPosIntElems[i2Count]) {
 
 1341                                 integer First = (pWithDofs)->iGetFirstIndex();
 
 1345                                 for (
int j = 1; j < TmpDofNum; j++) {
 
 1357         for (
int i = 1; i < iNumIntNodes + 1; i++) {
 
 1367                 if (ppNodes[InterfNodes.
pAdjncy[i]]->iGetNumDof() != 0) {
 
 1368                         integer First = ppNodes[InterfNodes.
pAdjncy[i]]->iGetFirstIndex();
 
 1377                         for (
unsigned int j = 1;
 
 1378                                         j < ppNodes[InterfNodes.
pAdjncy[i]]->iGetNumDof();
 
 1396                 integer First = (pWithDofs)->iGetFirstIndex();
 
 1399                 for (
int j = 1; j < TmpDofNum; j++) {
 
 1407         if ( Vertices.
pXadj != NULL) {
 
 1411         if ( Vertices.
pAdjncy != NULL) {
 
 1415         if ( InterfNodes.
pXadj != NULL) {
 
 1419         if ( InterfNodes.
pAdjncy != NULL) {
 
 1431         if ( pPosIntElems != NULL) {
 
 1435         if ( lTypes != NULL) {
 
 1439         if ( llabels != NULL) {
 
 1449         silent_cout(
"Making partition file...");
 
 1454         time_t tCurrTime(time(NULL));
 
 1456                 << 
"# Partition file for MBDyn. Time: " << ctime(&tCurrTime)  << std::endl
 
 1457                 << 
"# Partition produced ";
 
 1484                 << 
"# Control data to verify the partition " 
 1487                 << 
"Total number of processes: " << DataCommSize << 
";" << std::endl
 
 1488                 << 
"Process #: " << 
MyRank << 
";" << std::endl
 
 1489                 << 
"Total number of Nodes: " << iTotNodes << 
";"  << std::endl
 
 1490                 << 
"Total number of Elements: " << 
Elems.size() << 
";"  << std::endl
 
 1491                 << 
"Total number of Dofs: " << 
iTotDofs << 
";"  << std::endl;
 
 1496                 << 
"Local Dofs number: " << 
iNumLocDofs << std::endl;
 
 1500                 << 
"Local Nodes: " << iNumLocNodes << std::endl;
 
 1514                 << 
"Local Elements: "<< iNumLocElems << std::endl;
 
 1528                 << 
"Interface Dofs number: " << 
iNumIntDofs << std::endl
 
 1530                 << 
"Interface Nodes: " << iNumIntNodes << std::endl;
 
 1544                 << 
"Elements whose internal dofs are interface dofs: " 
 1545                 << iNumIntElems << std::endl;
 
 1559                 << 
"Local Dofs List:" << std::endl;
 
 1563                 if (!(iNumLocDofs%10)) {
 
 1567         if (iNumLocDofs%10) {
 
 1575                 << 
"Local Interface Dofs List:" <<std::endl;
 
 1585         silent_cout(
"... partition file done" << std::endl);
 
 1595         for (; pOld < pList + dim; pOld++) {
 
 1596                 if (*pOld != ADJ_UNDEFINED) {
 
 1607         for (
int i = 0; i < dim; i++) {
 
 1615         for (
int i = 0; i < dim; i++) {
 
 1625         unsigned int* pbegin = 
pLabelsList + (ppFirst - ppNodes);
 
 1626         unsigned int* p = std::find(pbegin, pbegin + dim, label);
 
 1635         DEBUGCOUT(
"Entering SchurDataManager::AssRes()" << std::endl);
 
 1636         ASSERT(pWorkVec != NULL);
 
 1641         } 
catch (ChangedEquationStructure) {
 
 1644                 if (MyElemIter.bGetCurr(pEl) == 
true) {
 
 1645                         silent_cerr(
"Jacobian reassembly requested by " 
 1648                                 "currently unsupported by Schur data manager." 
 1652                         silent_cerr(
"Jacobian reassembly requested " 
 1653                                 "by an element; currently unsupported " 
 1654                                 "by Schur data manager." << std::endl);
 
 1664         DEBUGCOUT(
"Entering SchurDataManager::AssJac()" << std::endl);
 
 1674         DEBUGCOUT(
"Entering SchurDataManager::Update()" << std::endl);
 
 1699         DEBUGCOUT(
"Entering SchurDataManager::DerivativesUpdate()" << std::endl);
 
 1734         DEBUGCOUT(
"Entering SchurDataManager::BeforePredict()" << std::endl);
 
 1759         DEBUGCOUT(
"Entering SchurDataManager::AfterPredict()" << std::endl);
 
 1788         DEBUGCOUT(
"Entering SchurDataManager::AfterConvergence()" << std::endl);
 
 1840         DEBUGCOUT(
"Entering SchurDataManager::Output()" << std::endl);
 
integer * GetDofsList(DofType who) const 
enum SchurDataManager::PartitionLibrary Partitioner
ElemContainerType ElemContainer
#define MBDYN_EXCEPT_ARGS
#define DEBUGCOUTFNAME(fname)
void DerivativesUpdate(void) const 
virtual integer GetInt(integer iDefval=0)
doublereal dGetTime(void) const 
const integer iDefaultMaxInterfNodes
static int position(const MathParser::MathArgs &args)
#define SAFEDELETEARR(pnt)
void BeforePredict(VectorHandler &X, VectorHandler &XP, VectorHandler &XPrev, VectorHandler &XPPrev) const 
int mbdyn_METIS_PartGraph(int iTotVertices, int *pXadj, int *pAdjncy, int *pVertexWgts, int *pCommWgts, int *pEdgeWgts, int DataCommSize, int *pParAmgProcs)
virtual void AssRes(VectorHandler &ResHdl, doublereal dCoef)
Node ** SearchNode(Node **ppFirst, int dim, unsigned int &label)
virtual Elem::Type GetElemType(void) const =0
VecIter< Elem * > MyElemIter
void AssRes(VectorHandler &ResHdl, doublereal dCoef)
void AfterPredict(void) const 
std::ostream & Partition(void) const 
void AssJac(MatrixHandler &JacHdl, doublereal dCoef)
void Output(bool force=false) const 
virtual const InducedVelocity * pGetInducedVelocity(void) const 
DriveCaller * pOutputMeter
unsigned int * pLabelsList
virtual void MakeRestart(void)
virtual bool IsKeyWord(const char *sKeyWord)
doublereal dLastRestartTime
struct DataManager::NodeDataStructure NodeData[Node::LASTNODETYPE]
virtual void AssJac(MatrixHandler &JacHdl, doublereal dCoef)
VectorHandler * pXPrimeCurr
integer iRestartIterations
#define ASSERT(expression)
Elem ** ppFindElem(Elem::Type Typ, unsigned int uElem) const 
virtual void BeforePredict(VectorHandler &, VectorHandler &, VectorHandler &, VectorHandler &) const 
virtual unsigned int iGetNumDof(void) const 
void OutputPartition(void)
virtual doublereal dGet(const doublereal &dVar) const =0
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const 
virtual unsigned int iGetNumDof(void) const =0
void Init(const T *p, unsigned i)
const char * psNodeNames[]
Dof * pGetDofsList(void) const 
virtual Node::Type GetNodeType(void) const =0
#define SAFENEWARRNOFILL(pnt, item, sz)
virtual void Output(OutputHandler &OH) const 
integer HowManyDofs(DofType who) const 
const char * psElemNames[]
SchurDataManager(MBDynParser &HP, unsigned OF, Solver *pS, doublereal dInitialTime, const char *sOutputFileName, const char *sInputFileName, bool bAbortAfterInput)
virtual int GetWord(void)
virtual void AfterConvergence(const VectorHandler &X, const VectorHandler &XP)
void Pack(int *pList, int dim)
void chaco_interface(const int iTotVertices, int *start, int *adjacency, int *vertex_weights, int *comm_weights, int *edge_weights, const int num_processors, int *pParAmgProcs)
#define SAFENEWARR(pnt, item, sz)
virtual void AfterPredict(VectorHandler &X, VectorHandler &XP)
virtual HighParser::ErrOut GetLineData(void) const 
void AfterConvergence(void) const 
unsigned int GetLabel(void) const 
struct DataManager::ElemDataStructure ElemData[Elem::LASTELEMTYPE]
Node * pFindNode(Node::Type Typ, unsigned int uNode) const 
void CreatePartition(void)
const integer iDefaultMaxNodesPerElem
const integer iDefaultMaxConnectionsPerVertex
VariableSubMatrixHandler * pWorkMat
void InitList(int *list, int dim, int value)
virtual void Update(const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)