src/MEDWrapper/Base/MED_Algorithm.cxx

   1 // Copyright (C) 2007-2015  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License, or (at your option) any later version.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22 #include "MED_Algorithm.hxx"
  23 #include "MED_Wrapper.hxx"
  24
  25 #include "MED_Utilities.hxx"
  26
  27 #ifdef _DEBUG_
  28 static int MYDEBUG = 0;
  29 static int MYVALUEDEBUG = 0;
  30 #else
  31 // static int MYDEBUG = 0;
  32 // static int MYVALUEDEBUG = 0;
  33 #endif
  34
  35 namespace MED
  36 {
  37   //---------------------------------------------------------------
  38   TEntity2TGeom2ElemInfo
  39   GetEntity2TGeom2ElemInfo(const PWrapper& theWrapper,
  40                            const PMeshInfo& theMeshInfo,
  41                            const MED::TEntityInfo& theEntityInfo)
  42   {
  43     MSG(MYDEBUG,"GetElemsByEntity(...)");
  44     TEntity2TGeom2ElemInfo anEntity2TGeom2ElemInfo;
  45     MED::TEntityInfo::const_iterator anIter = theEntityInfo.begin();
  46     PElemInfo anElemInfo;
  47     TErr anErr;
  48     for(; anIter != theEntityInfo.end(); anIter++){
  49       const EEntiteMaillage& anEntity = anIter->first;
  50       const TGeom2Size& aGeom2Size = anIter->second;
  51       TGeom2ElemInfo& aGeom2ElemInfo = anEntity2TGeom2ElemInfo[anEntity];
  52
  53       if(anEntity == eNOEUD){
  54         aGeom2ElemInfo[ePOINT1] = theWrapper->GetPElemInfo(theMeshInfo);
  55         continue;
  56       }
  57
  58       TGeom2Size::const_iterator anIter2 = aGeom2Size.begin();
  59       for(; anIter2 != aGeom2Size.end(); anIter2++){
  60         const EGeometrieElement& aGeom = anIter2->first;
  61         aGeom2ElemInfo[aGeom] = theWrapper->GetPElemInfo(theMeshInfo,anEntity,aGeom,MED::eNOD,&anErr);
  62       }
  63     }
  64     ADDMSG(MYDEBUG,"\n");
  65     return anEntity2TGeom2ElemInfo;
  66   }
  67
  68
  69   //---------------------------------------------------------------
  70   TFamilyInfoSet
  71   GetFamilyInfoSet(const PWrapper& theWrapper,
  72                    const PMeshInfo& theMeshInfo)
  73   {
  74     MSG(MYDEBUG,"GetFamilies(...)");
  75     TErr anErr;
  76     TFamilyInfoSet aFamilyInfoSet;
  77     TInt aNbFam = theWrapper->GetNbFamilies(*theMeshInfo);
  78     INITMSG(MYDEBUG,"GetNbFamilies() = "<<aNbFam<<"\n");
  79     for(TInt iFam = 1; iFam <= aNbFam; iFam++){
  80       PFamilyInfo aFamilyInfo = theWrapper->GetPFamilyInfo(theMeshInfo,iFam,&anErr);
  81       if(anErr >= 0)
  82         aFamilyInfoSet.insert(aFamilyInfo);
  83     }
  84     ADDMSG(MYDEBUG,"\n");
  85     return aFamilyInfoSet;
  86   }
  87
  88
  89   //---------------------------------------------------------------
  90   TGroupInfo
  91   GetGroupInfo(const TFamilyInfoSet& theFamilyInfoSet)
  92   {
  93     MSG(MYDEBUG,"GetFamiliesByGroup(...)");
  94     TGroupInfo aGroup;
  95     TFamilyInfoSet::const_iterator anIter = theFamilyInfoSet.begin();
  96     for(; anIter != theFamilyInfoSet.end(); anIter++){
  97       const PFamilyInfo& aFamilyInfo = *anIter;
  98       TInt aNbGroup = aFamilyInfo->GetNbGroup();
  99       for(TInt iGroup = 0; iGroup < aNbGroup; iGroup++){
 100         aGroup[aFamilyInfo->GetGroupName(iGroup)].insert(aFamilyInfo);
 101       }
 102     }
 103
 104 #ifdef _DEBUG_
 105     if(MYDEBUG){
 106       TGroupInfo::const_iterator anIter = aGroup.begin();
 107       for(; anIter != aGroup.end(); anIter++){
 108         const std::string& aName = anIter->first;
 109         INITMSG(MYDEBUG,"aGroupName = '"<<aName<<"'\n");
 110         const TFamilyInfoSet& aFamilyInfoSet = anIter->second;
 111         TFamilyInfoSet::const_iterator anFamIter = aFamilyInfoSet.begin();
 112         for(; anFamIter != aFamilyInfoSet.end(); anFamIter++){
 113           const PFamilyInfo& aFamilyInfo = *anFamIter;
 114           INITMSG(MYDEBUG,"aFamilyName = '"<<aFamilyInfo->GetName()<<"'\n");
 115         }
 116       }
 117       ADDMSG(MYDEBUG,"\n");
 118     }
 119 #endif
 120
 121     return aGroup;
 122   }
 123
 124
 125   //---------------------------------------------------------------
 126   TFieldInfo2TimeStampInfoSet
 127   GetFieldInfo2TimeStampInfoSet(const PWrapper& theWrapper,
 128                                 const PMeshInfo& theMeshInfo,
 129                                 const MED::TEntityInfo& theEntityInfo)
 130   {
 131     MSG(MYDEBUG,"GetFieldsByEntity(...)");
 132     TFieldInfo2TimeStampInfoSet aFieldInfo2TimeStampInfoSet;
 133     TInt aNbFields = theWrapper->GetNbFields();
 134     INITMSG(MYDEBUG,"GetNbFields() = "<<aNbFields<<"\n");
 135     for(TInt iField = 1; iField <= aNbFields; iField++){
 136       PFieldInfo aFieldInfo = theWrapper->GetPFieldInfo(theMeshInfo,iField);
 137       INITMSG(MYDEBUG,"aFieldName = '"<<aFieldInfo->GetName()<<
 138               "'; aNbComp = "<<aFieldInfo->GetNbComp()<<"; ");
 139       TGeom2Size aGeom2Size;
 140       EEntiteMaillage anEntity = EEntiteMaillage(-1);
 141       TInt aNbTimeStamps = theWrapper->GetNbTimeStamps(aFieldInfo,theEntityInfo,anEntity,aGeom2Size);
 142       ADDMSG(MYDEBUG,"anEntity = "<<anEntity<<"; GetNbTimeStamps = "<<aNbTimeStamps<<"\n");
 143       for(TInt iTimeStamp = 1; iTimeStamp <= aNbTimeStamps; iTimeStamp++){
 144         PTimeStampInfo aTimeStamp =
 145           theWrapper->GetPTimeStampInfo(aFieldInfo,anEntity,aGeom2Size,iTimeStamp);
 146         aFieldInfo2TimeStampInfoSet[aFieldInfo].insert(aTimeStamp);
 147         INITMSG(MYDEBUG,
 148                 "aDt = "<<aTimeStamp->GetDt()<<
 149                 ", Unit = \'"<<aTimeStamp->GetUnitDt()<<"\n");
 150       }
 151     }
 152     ADDMSG(MYDEBUG,"\n");
 153     return aFieldInfo2TimeStampInfoSet;
 154   }
 155
 156
 157   //---------------------------------------------------------------
 158   TEntite2TFieldInfo2TimeStampInfoSet
 159   GetEntite2TFieldInfo2TimeStampInfoSet(const TFieldInfo2TimeStampInfoSet& theFieldInfo2TimeStampInfoSet)
 160   {
 161     TEntite2TFieldInfo2TimeStampInfoSet anEntite2TFieldInfo2TimeStampInfoSet;
 162     TFieldInfo2TimeStampInfoSet::const_iterator anIter = theFieldInfo2TimeStampInfoSet.begin();
 163     for(; anIter != theFieldInfo2TimeStampInfoSet.end(); anIter++){
 164       const TTimeStampInfoSet& aTimeStampInfoSet = anIter->second;
 165       //const PFieldInfo& aFieldInfo = anIter->first;
 166       if(aTimeStampInfoSet.empty())
 167         continue;
 168       const PTimeStampInfo& aTimeStampInfo = *aTimeStampInfoSet.begin();
 169       anEntite2TFieldInfo2TimeStampInfoSet[ConvertEntity(aTimeStampInfo->GetEntity())].insert(*anIter);
 170     }
 171     return anEntite2TFieldInfo2TimeStampInfoSet;
 172   }
 173
 174
 175   //---------------------------------------------------------------
 176   bool
 177   operator<(const TFamilyTSize& theLeft, const TFamilyTSize& theRight)
 178   {
 179     const MED::PFamilyInfo& aLeftInfo = boost::get<0>(theLeft);
 180     const MED::PFamilyInfo& aRightInfo = boost::get<0>(theRight);
 181     return aLeftInfo->GetId() < aRightInfo->GetId();
 182   }
 183
 184
 185   //---------------------------------------------------------------
 186   TEntity2FamilySet
 187   GetEntity2FamilySet(const PWrapper& theWrapper,
 188                       const TEntity2TGeom2ElemInfo& theEntity2TGeom2ElemInfo,
 189                       const TFamilyInfoSet& theFamilyInfoSet)
 190   {
 191     MSG(MYDEBUG,"GetFamiliesByEntity(...)");
 192     TEntity2FamilySet anEntity2FamilySet;
 193
 194     typedef std::map<TInt,PFamilyInfo> TId2Family;
 195     TId2Family anId2Family;
 196     TFamilyInfoSet::const_iterator anIter = theFamilyInfoSet.begin();
 197     for(; anIter != theFamilyInfoSet.end(); anIter++){
 198       const PFamilyInfo& aFamilyInfo = *anIter;
 199       anId2Family.insert(TId2Family::value_type(aFamilyInfo->GetId(),aFamilyInfo));
 200     }
 201
 202     if(!anId2Family.empty()){
 203       typedef std::map<TInt,TInt> TFamilyID2Size;
 204       typedef std::map<EEntiteMaillage,TFamilyID2Size> TEntity2FamilyID;
 205       TEntity2FamilyID anEntity2FamilyID;
 206
 207       if(!theEntity2TGeom2ElemInfo.empty()){
 208         TEntity2TGeom2ElemInfo::const_iterator anIter = theEntity2TGeom2ElemInfo.begin();
 209         for(; anIter != theEntity2TGeom2ElemInfo.end(); anIter++){
 210           const EEntiteMaillage& anEntity = anIter->first;
 211           TFamilyID2Size& aFamilyID2Size = anEntity2FamilyID[anEntity];
 212           const TGeom2ElemInfo& aGeom2ElemInfo = anIter->second;
 213           TGeom2ElemInfo::const_iterator aGeom2ElemInfoIter = aGeom2ElemInfo.begin();
 214           for(; aGeom2ElemInfoIter != aGeom2ElemInfo.end(); aGeom2ElemInfoIter++){
 215             const PElemInfo& aElemInfo = aGeom2ElemInfoIter->second;
 216             if(TInt aNbElem = aElemInfo->GetNbElem()){
 217               for(TInt i = 0; i < aNbElem; i++){
 218                 aFamilyID2Size[aElemInfo->GetFamNum(i)] += 1;
 219               }
 220             }
 221           }
 222         }
 223       }
 224
 225       if(!anEntity2FamilyID.empty()){
 226         TEntity2FamilyID::const_iterator anIter = anEntity2FamilyID.begin();
 227         for(; anIter != anEntity2FamilyID.end(); anIter++){
 228           const EEntiteMaillage& anEntity = anIter->first;
 229           INITMSG(MYDEBUG,"anEntity = "<<anEntity<<":\n");
 230           const TFamilyID2Size& aFamilyID2Size = anIter->second;
 231           TFamilyID2Size::const_iterator anIter2 = aFamilyID2Size.begin();
 232           for(; anIter2 != aFamilyID2Size.end(); anIter2++){
 233             TInt anId = anIter2->first;
 234             TInt aSize = anIter2->second;
 235             TId2Family::const_iterator anIter3 = anId2Family.find(anId);
 236             if(anIter3 != anId2Family.end()){
 237               const PFamilyInfo& aFamilyInfo = anIter3->second;
 238               anEntity2FamilySet[anEntity].insert(TFamilyTSize(aFamilyInfo,aSize));
 239               INITMSG(MYDEBUG,
 240                       "aFamilyName = '"<<aFamilyInfo->GetName()<<
 241                       "' anId = "<<aFamilyInfo->GetId()<<"\n");
 242             }
 243           }
 244         }
 245       }
 246     }
 247     ADDMSG(MYDEBUG,"\n");
 248     return anEntity2FamilySet;
 249   }
 250
 251
 252   //---------------------------------------------------------------
 253   TKey2Gauss
 254   GetKey2Gauss(const PWrapper& theWrapper,
 255                TErr* theErr,
 256                EModeSwitch theMode)
 257   {
 258     INITMSG(MYDEBUG,"GetKey2Gauss - theMode = "<<theMode<<std::endl);
 259     TKey2Gauss aKey2Gauss;
 260     TInt aNbGauss = theWrapper->GetNbGauss(theErr);
 261     for(TInt anId = 1; anId <= aNbGauss; anId++){
 262       TGaussInfo::TInfo aPreInfo = theWrapper->GetGaussPreInfo(anId);
 263       PGaussInfo anInfo = theWrapper->CrGaussInfo(aPreInfo,theMode);
 264       theWrapper->GetGaussInfo(anId,anInfo,theErr);
 265       TGaussInfo::TKey aKey = boost::get<0>(aPreInfo);
 266       aKey2Gauss[aKey] = anInfo;
 267
 268 #ifdef _DEBUG_
 269       const EGeometrieElement& aGeom = boost::get<0>(aKey);
 270       const std::string& aName = boost::get<1>(aKey);
 271       INITMSG(MYDEBUG,
 272               "- aGeom = "<<aGeom<<
 273               "; aName = '"<<aName<<"'"<<
 274               std::endl);
 275 #endif
 276
 277     }
 278     return aKey2Gauss;
 279   }
 280
 281
 282   //---------------------------------------------------------------
 283   PProfileInfo
 284   GetProfileInfo(const PWrapper& theWrapper,
 285                  const std::string& theProfileName,
 286                  TErr* theErr,
 287                  EModeProfil theMode)
 288   {
 289     PProfileInfo anInfo;
 290     TInt aNbProfiles = theWrapper->GetNbProfiles(theErr);
 291     for(TInt anId = 1; anId <= aNbProfiles; anId++){
 292       TProfileInfo::TInfo aPreInfo = theWrapper->GetProfilePreInfo(anId);
 293       const std::string& aName = boost::get<0>(aPreInfo);
 294       if(aName == theProfileName)
 295         return theWrapper->GetPProfileInfo(anId,theMode,theErr);
 296     }
 297     return anInfo;
 298   }
 299
 300
 301   //---------------------------------------------------------------
 302   TMKey2Profile
 303   GetMKey2Profile(const PWrapper& theWrapper,
 304                   TErr* theErr,
 305                   EModeProfil theMode)
 306   {
 307     INITMSG(MYDEBUG,"GetMKey2Profile - theMode = "<<theMode<<std::endl);
 308     TKey2Profile aKey2Profile;
 309     TInt aNbProfiles = theWrapper->GetNbProfiles(theErr);
 310     for(TInt anId = 1; anId <= aNbProfiles; anId++){
 311       TProfileInfo::TInfo aPreInfo = theWrapper->GetProfilePreInfo(anId);
 312       PProfileInfo anInfo = theWrapper->GetPProfileInfo(anId,theMode,theErr);
 313       const std::string& aName = boost::get<0>(aPreInfo);
 314       aKey2Profile[aName] = anInfo;
 315
 316 #ifdef _DEBUG_
 317       INITMSG(MYDEBUG,
 318               "- aName = '"<<aName<<"'"<<
 319               " : "<<
 320               std::endl);
 321       TInt aNbElem = anInfo->GetSize();
 322       for(TInt iElem = 0; iElem < aNbElem; iElem++){
 323         ADDMSG(MYVALUEDEBUG,anInfo->GetElemNum(iElem)<<", ");
 324       }
 325       ADDMSG(MYVALUEDEBUG, std::endl);
 326 #endif
 327
 328     }
 329     return TMKey2Profile(theMode,aKey2Profile);
 330   }
 331
 332   //---------------------------------------------------------------
 333   EEntiteMaillage
 334   GetEntityByFamilyId(PGrilleInfo& theInfo,TInt theId){
 335     TElemNum::iterator aNodeFamIter = (theInfo->myFamNumNode).begin();
 336     for(;aNodeFamIter != (theInfo->myFamNumNode).end(); aNodeFamIter++){
 337       if(theId == *aNodeFamIter)
 338         return eNOEUD;
 339     }
 340     TElemNum::iterator aCellFamIter = (theInfo->myFamNum).begin();
 341     for(;aCellFamIter != (theInfo->myFamNum).end(); aCellFamIter++){
 342       if(theId == *aCellFamIter)
 343         return eMAILLE;
 344     }
 345     EXCEPTION(std::runtime_error, "GetEntityByFamilyId - fails");
 346     return EEntiteMaillage(-1);
 347   }
 348
 349   TFamilyID2NbCells
 350   GetFamilyID2NbCells(PGrilleInfo& theInfo){
 351     TFamilyID2NbCells aFamily2NbCells;
 352     TInt aNbNodes = theInfo->myFamNumNode.size();
 353     TInt aNbCells = theInfo->myFamNum.size();
 354     for(TInt i=0; i<aNbNodes; i++) aFamily2NbCells[theInfo->GetFamNumNode(i)] = 0;
 355     for(TInt i=0; i<aNbCells; i++) aFamily2NbCells[theInfo->GetFamNum(i)] = 0;
 356     for(TInt i=0; i<aNbNodes; i++) aFamily2NbCells[theInfo->GetFamNumNode(i)] += 1;
 357     for(TInt i=0; i<aNbCells; i++) aFamily2NbCells[theInfo->GetFamNum(i)] += 1;
 358     return aFamily2NbCells;
 359   }
 360
 361   EEntiteMaillage ConvertEntity(const EEntiteMaillage& aEntity){
 362     switch( aEntity ){
 363
 364     case eNOEUD_ELEMENT:
 365     case eMAILLE: return eMAILLE; //eNOEUD_ELEMENT it is eMAILLE
 366
 367     case eFACE:
 368     case eARETE:
 369     case eNOEUD: return aEntity; break;
 370     default: return EEntiteMaillage(-1);
 371
 372     }
 373   }
 374 }