From f4ca8ef1016a9177a5b666b4acda00b03b877c37 Mon Sep 17 00:00:00 2001
From: apo <apo@opencascade.com>
Date: Mon, 3 Oct 2005 06:24:37 +0000
Subject: [PATCH] To introduce quadratic mesh elements

---
 src/CONVERTOR/VISU_Convertor_impl.cxx |  48 +++++--
 src/CONVERTOR/VISU_Convertor_impl.hxx |   7 +
 src/CONVERTOR/VISU_MedConvertor.cxx   | 178 ++++++++++++++++++++++----
 3 files changed, 193 insertions(+), 40 deletions(-)

diff --git a/src/CONVERTOR/VISU_Convertor_impl.cxx b/src/CONVERTOR/VISU_Convertor_impl.cxx
index 196b65bc..94a61b7e 100644
--- a/src/CONVERTOR/VISU_Convertor_impl.cxx
+++ b/src/CONVERTOR/VISU_Convertor_impl.cxx
@@ -841,9 +841,32 @@ namespace VISU
     case VTK_HEXAHEDRON: return 8;
     case VTK_WEDGE: return 6;
     case VTK_PYRAMID: return 5;
+
+    case QUADRATIC_EDGE: return 3;
+    case QUADRATIC_TRIANGLE: return 6;
+    case QUADRATIC_QUAD: return 8;
+    case QUADRATIC_TETRA: return 10;
+    case QUADRATIC_HEXAHEDRON: return 20;
+    case QUADRATIC_WEDGE: return 15;
+    case QUADRATIC_PYRAMID: return 13;
     }
     return -1;
   }
+
+  vtkIdType
+  VTKGeom2VISU(vtkIdType theGeom)
+  { 
+    switch(theGeom){
+    case QUADRATIC_EDGE: return VTK_POLY_LINE;
+    case QUADRATIC_TRIANGLE: return VTK_POLYGON;
+    case QUADRATIC_QUAD: return VTK_POLYGON;
+    case QUADRATIC_TETRA: return VTK_CONVEX_POINT_SET;
+    case QUADRATIC_HEXAHEDRON: return VTK_CONVEX_POINT_SET;
+    case QUADRATIC_WEDGE: return VTK_CONVEX_POINT_SET;
+    case QUADRATIC_PYRAMID: return VTK_CONVEX_POINT_SET;
+    }
+    return theGeom;
+  }
 }
 
 
@@ -1245,15 +1268,16 @@ namespace
   {
     vtkIdType aNbCells = theSubProfile->myNbCells;
     vtkIdType aCellsSize = theSubProfile->myCellsSize;
-    vtkIdType aGeom = theSubProfile->myGeom;
-    vtkIdType aNbNodes = VTKGeom2NbNodes(aGeom);
-    
+    vtkIdType aVGeom = theSubProfile->myGeom;
+    vtkIdType aNbNodes = VTKGeom2NbNodes(aVGeom);
+    vtkIdType aGeom = VTKGeom2VISU(aVGeom);
+
     const TSubMeshID& aSubMeshID = theSubProfile->mySubMeshID;
 
     const VISU::TGeom2SubMesh& aGeom2SubMesh = theMeshOnEntity->myGeom2SubMesh;
-    VISU::TGeom2SubMesh::const_iterator anIter = aGeom2SubMesh.find(aGeom);
+    VISU::TGeom2SubMesh::const_iterator anIter = aGeom2SubMesh.find(aVGeom);
     if(anIter == aGeom2SubMesh.end())
-      EXCEPTION(runtime_error,"GetCells >> There is no elements for the GEOM("<<aGeom<<")");
+      EXCEPTION(runtime_error,"GetCells >> There is no elements for the GEOM("<<aVGeom<<")");
     
     const VISU::TSubMeshImpl& aSubMesh = anIter->second;
     const TCell2Connect& aCell2Connect = aSubMesh.myCell2Connect;
@@ -1308,7 +1332,7 @@ namespace
 		      const PProfileImpl& theProfile,
 		      const PSubProfileImpl& theSubProfile)
   {
-    INITMSG(MYDEBUG,"GetMeshOnSubProfile - aGeom = "<<theSubProfile->myGeom<<endl);
+    INITMSG(MYDEBUG,"GetMeshOnSubProfile - aVGeom = "<<theSubProfile->myGeom<<endl);
     
     const TVTKSource& aSource = theSubProfile->GetSource();
     if(theSubProfile->myIsVTKDone)
@@ -1422,12 +1446,12 @@ namespace
     TGeom2GaussSubMesh::const_iterator anIter = aGeom2GaussSubMesh.begin();
     TVector<float> aDataValues(aNbComp);
     for(int aTupleId = 0; anIter != aGeom2GaussSubMesh.end(); anIter++){
-      vtkIdType aGeom = anIter->first;
+      vtkIdType aVGeom = anIter->first;
       PGaussSubMeshImpl aGaussSubMesh = anIter->second;
       if(!aGaussSubMesh->myIsDone)
 	continue;
       
-      TGeom2Value::const_iterator anIter2 = aGeom2Value.find(aGeom);
+      TGeom2Value::const_iterator anIter2 = aGeom2Value.find(aVGeom);
       if(anIter2 == aGeom2Value.end()){
 	EXCEPTION(runtime_error,
 		  "GetTimeStampOnGaussMesh >> Can't find values for corresponding Gauss Points SubMesh");
@@ -1439,7 +1463,7 @@ namespace
       if(aNbGauss < 1)
 	continue;
 
-      INITMSG(MYDEBUG,"- aGeom = "<<aGeom<<
+      INITMSG(MYDEBUG,"- aVGeom = "<<aVGeom<<
 	      "; aNbElem = "<<aNbElem<<
 	      "; aNbComp = "<<aNbComp<<
 	      "; aNbGauss = "<<aNbGauss<<
@@ -1532,7 +1556,7 @@ namespace
 		  const PGaussSubMeshImpl& theGaussSubMesh)
   {
     PGaussImpl aGauss = theGaussSubMesh->myGauss;
-    INITMSG(MYDEBUG,"GetGaussSubMesh - aGeom = "<<aGauss->myGeom<<endl);
+    INITMSG(MYDEBUG,"GetGaussSubMesh - aVGeom = "<<aGauss->myGeom<<endl);
     
     if(!theGaussSubMesh->myIsDone)
       return;
@@ -1657,11 +1681,11 @@ VISU_Convertor_impl
       aSubMeshArr.resize(aGeom2SubMesh.size());
 
       for(vtkIdType anID = 0, aCellID = 0; anIter != aGeom2SubMesh.end(); anIter++, anID++){
-	const vtkIdType& aGeom = anIter->first;
+	const vtkIdType& aVGeom = anIter->first;
 	PSubMeshImpl aSubMesh = anIter->second;
 	const TVTKSource& aSource = aSubMesh->GetSource();
 	aSource->SetPoints(GetPoints(aMesh));
-	GetCellsOnSubMesh(aSource,aMeshOnEntity,aSubMesh,aGeom);
+	GetCellsOnSubMesh(aSource,aMeshOnEntity,aSubMesh,aVGeom);
 	anAppendFilter->AddInput(aSource.GetPointer());
 
 	aSubMesh->myStartID = aCellID;
diff --git a/src/CONVERTOR/VISU_Convertor_impl.hxx b/src/CONVERTOR/VISU_Convertor_impl.hxx
index c7d231a0..eae3b60b 100644
--- a/src/CONVERTOR/VISU_Convertor_impl.hxx
+++ b/src/CONVERTOR/VISU_Convertor_impl.hxx
@@ -24,6 +24,13 @@ class VTKViewer_AppendFilter;
 
 namespace VISU
 {
+  const vtkIdType QUADRATIC_EDGE = 30;
+  const vtkIdType QUADRATIC_TRIANGLE = 31;
+  const vtkIdType QUADRATIC_QUAD = 32;
+  const vtkIdType QUADRATIC_TETRA = 33;
+  const vtkIdType QUADRATIC_HEXAHEDRON = 34;
+  const vtkIdType QUADRATIC_WEDGE = 35;
+  const vtkIdType QUADRATIC_PYRAMID = 36;
 
   using MED::TCSlice;
   using MED::TSlice;
diff --git a/src/CONVERTOR/VISU_MedConvertor.cxx b/src/CONVERTOR/VISU_MedConvertor.cxx
index c0530b21..54d2037e 100644
--- a/src/CONVERTOR/VISU_MedConvertor.cxx
+++ b/src/CONVERTOR/VISU_MedConvertor.cxx
@@ -70,7 +70,7 @@ static int MY_GROUP_DEBUG = 0;
 namespace
 {
   //---------------------------------------------------------------
-  int
+  vtkIdType
   MEDGeom2NbNodes(MED::EGeometrieElement theMEDGeomType)
   { 
     return theMEDGeomType % 100;
@@ -78,25 +78,27 @@ namespace
 
 
   //---------------------------------------------------------------
-  int
+  vtkIdType
   MEDGeomToVTK(MED::EGeometrieElement theMEDGeomType)
   { 
     switch(theMEDGeomType){
     case MED::ePOINT1: return VTK_VERTEX;
     case MED::eSEG2: return VTK_LINE;
-    case MED::eSEG3: return VTK_LINE;
     case MED::eTRIA3: return VTK_TRIANGLE;
-    case MED::eTRIA6: return VTK_TRIANGLE;
     case MED::eQUAD4: return VTK_QUAD;
-    case MED::eQUAD8: return VTK_QUAD;
     case MED::eTETRA4: return VTK_TETRA;
-    case MED::eTETRA10: return VTK_TETRA;
     case MED::eHEXA8: return VTK_HEXAHEDRON;
-    case MED::eHEXA20: return VTK_HEXAHEDRON;
     case MED::ePENTA6: return VTK_WEDGE;
-    case MED::ePENTA15: return VTK_WEDGE;
     case MED::ePYRA5: return VTK_PYRAMID;
-    case MED::ePYRA13: return VTK_PYRAMID;
+
+    case MED::eSEG3: return VISU::QUADRATIC_EDGE;
+    case MED::eTRIA6: return VISU::QUADRATIC_TRIANGLE;
+    case MED::eQUAD8: return VISU::QUADRATIC_QUAD;
+    case MED::eTETRA10: return VISU::QUADRATIC_TETRA;
+    case MED::eHEXA20: return VISU::QUADRATIC_HEXAHEDRON;
+    case MED::ePENTA15: return VISU::QUADRATIC_WEDGE;
+    case MED::ePYRA13: return VISU::QUADRATIC_PYRAMID;
+
     case MED::ePOLYGONE: return VTK_POLYGON;
     case MED::ePOLYEDRE: return VTK_CONVEX_POINT_SET;
     }
@@ -117,6 +119,15 @@ namespace
     case VTK_HEXAHEDRON: return MED::eHEXA8;
     case VTK_WEDGE: return MED::ePENTA6;
     case VTK_PYRAMID: return MED::ePYRA5;
+
+    case VISU::QUADRATIC_EDGE: return MED::eSEG3;
+    case VISU::QUADRATIC_TRIANGLE: return MED::eTRIA6;
+    case VISU::QUADRATIC_QUAD: return MED::eQUAD8;
+    case VISU::QUADRATIC_TETRA: return MED::eTETRA10;
+    case VISU::QUADRATIC_HEXAHEDRON: return MED::eHEXA20;
+    case VISU::QUADRATIC_WEDGE: return MED::ePENTA15;
+    case VISU::QUADRATIC_PYRAMID: return MED::ePYRA13;
+
     case VTK_POLYGON: return MED::ePOLYGONE;
     case VTK_CONVEX_POINT_SET: return MED::ePOLYEDRE;
     }
@@ -154,12 +165,15 @@ namespace
 
   //---------------------------------------------------------------
   PMEDSubProfile
-  CrSubProfile(MED::EGeometrieElement theMGeom,
+  CrSubProfile(const MED::PWrapper& theMEDWrapper,
+	       const MED::PMeshInfo& theMeshInfo,
+	       MED::EEntiteMaillage theMEntity,
+	       MED::EGeometrieElement theMGeom,
 	       const MED::TGeom2Size& theGeom2Size,
 	       const MED::TGeom2Profile& theGeom2Profile)
   {
     vtkIdType aVGeom = MEDGeomToVTK(theMGeom);
-    vtkIdType aVNbNodes = VTKGeom2NbNodes(aVGeom);
+    vtkIdType aVNbNodes = VTKGeom2NbNodes(theMGeom);
 
     PMEDSubProfile aSubProfile(new TMEDSubProfile());
     aSubProfile->myGeom = aVGeom;
@@ -170,18 +184,81 @@ namespace
     if(aTimeStampIter == theGeom2Size.end())
       aSubProfile->myStatus = eRemoveAll;
     else{
-      aSubProfile->myNbCells = aTimeStampIter->second;
       MED::TGeom2Profile::const_iterator aProfileIter = theGeom2Profile.find(theMGeom);
       if(aProfileIter != theGeom2Profile.end()){
 	MED::PProfileInfo aProfileInfo = aProfileIter->second;
 
-	aSubProfile->myNbCells = aProfileInfo->myElemNum.size();
 	aSubProfile->myName = aProfileInfo->GetName();
 	aSubProfile->myStatus = eAddPart;
+
+	const MED::TElemNum& anElemNum = aProfileInfo->myElemNum;
+	TInt aNbElem = anElemNum.size();
+	aSubProfile->myNbCells = aNbElem;
+	switch(theMGeom){
+	case MED::ePOLYGONE: {
+	  MED::PPolygoneInfo aPolygoneInfo = 
+	    theMEDWrapper->GetPPolygoneInfo(theMeshInfo,theMEntity,theMGeom);
+	  for(TInt anElemId = 0; anElemId < aNbElem ; anElemId++){
+	    TInt aNbConn = aPolygoneInfo->GetNbConn(anElemNum[anElemId]);
+	    aSubProfile->myCellsSize += aNbConn;
+	  }
+	  break;
+	}
+	case MED::ePOLYEDRE: {
+	  MED::PPolyedreInfo aPolyedreInfo = 
+	    theMEDWrapper->GetPPolyedreInfo(theMeshInfo,theMEntity,theMGeom);
+	  for(TInt anElemId = 0; anElemId < aNbElem ; anElemId++){
+	    MED::TCConnSliceArr aConnSliceArr = 
+	      aPolyedreInfo->GetConnSliceArr(anElemNum[anElemId]);
+	    TInt aNbFaces = aConnSliceArr.size();
+	    TInt aCellSize = 0;
+	    for(TInt iFace = 0; iFace < aNbFaces; iFace++){
+	      MED::TCConnSlice aConnSlice = aConnSliceArr[iFace];
+	      TInt aNbConn = aConnSlice.size();
+	      aCellSize += aNbConn;
+	    }
+	    aSubProfile->myCellsSize += aCellSize;
+	  }
+	  break;
+	}
+	default: {
+	  aSubProfile->myCellsSize = aSubProfile->myNbCells*aVNbNodes;
+	}}
+      }else{
+	TInt aNbElem = aTimeStampIter->second;
+	aSubProfile->myNbCells = aNbElem;
+	switch(theMGeom){
+	case MED::ePOLYGONE: {
+	  MED::PPolygoneInfo aPolygoneInfo = 
+	    theMEDWrapper->GetPPolygoneInfo(theMeshInfo,theMEntity,theMGeom);
+	  for(TInt anElemId = 0; anElemId < aNbElem ; anElemId++){
+	    TInt aNbConn = aPolygoneInfo->GetNbConn(anElemId);
+	    aSubProfile->myCellsSize += aNbConn;
+	  }
+	  break;
+	}
+	case MED::ePOLYEDRE: {
+	  MED::PPolyedreInfo aPolyedreInfo = 
+	    theMEDWrapper->GetPPolyedreInfo(theMeshInfo,theMEntity,theMGeom);
+	  for(TInt anElemId = 0; anElemId < aNbElem ; anElemId++){
+	    MED::TCConnSliceArr aConnSliceArr = 
+	      aPolyedreInfo->GetConnSliceArr(anElemId);
+	    TInt aNbFaces = aConnSliceArr.size();
+	    TInt aCellSize = 0;
+	    for(TInt iFace = 0; iFace < aNbFaces; iFace++){
+	      MED::TCConnSlice aConnSlice = aConnSliceArr[iFace];
+	      TInt aNbConn = aConnSlice.size();
+	      aCellSize += aNbConn;
+	    }
+	    aSubProfile->myCellsSize += aCellSize;
+	  }
+	  break;
+	}
+	default: {
+	  aSubProfile->myCellsSize = aSubProfile->myNbCells*aVNbNodes;
+	}}
       }
     }
-    aSubProfile->myCellsSize = aSubProfile->myNbCells*aVNbNodes;
-    
     INITMSG(MYDEBUG,
 	    "- aVGeom = "<<aVGeom<<
 	    "; aName = '"<<aSubProfile->myName<<"'"<<
@@ -196,8 +273,11 @@ namespace
 
   //---------------------------------------------------------------
   TProfileKey
-  GetProfileKey(const MED::TTimeStampVal& theTimeStampVal,
+  GetProfileKey(const MED::PWrapper& theMEDWrapper,
+		const MED::PMeshInfo& theMeshInfo,
+		const MED::TTimeStampVal& theTimeStampVal,
 		const VISU::TMEDMeshOnEntity& theMeshOnEntity,
+		MED::EEntiteMaillage theMEntity,
 		const MED::TGeom2Size& theGeom2Size)
   {
     INITMSG(MYDEBUG,"GetProfileKey"<<endl);
@@ -209,7 +289,12 @@ namespace
     MED::TGeom2Size::const_iterator anIter = aGeom2Size.begin();
     for(; anIter != aGeom2Size.end(); anIter++){
       MED::EGeometrieElement aMGeom = anIter->first;
-      PSubProfile aSubProfile = CrSubProfile(aMGeom,theGeom2Size,aGeom2Profile);
+      PSubProfile aSubProfile = CrSubProfile(theMEDWrapper,
+					     theMeshInfo,
+					     theMEntity,
+					     aMGeom,
+					     theGeom2Size,
+					     aGeom2Profile);
       aProfileKey.insert(aSubProfile);
     }
     
@@ -219,8 +304,11 @@ namespace
 
   //---------------------------------------------------------------
   void
-  InitProfile(MED::TTimeStampVal& theTimeStampVal,
+  InitProfile(const MED::PWrapper& theMEDWrapper,
+	      const MED::PMeshInfo& theMeshInfo,
+	      MED::TTimeStampVal& theTimeStampVal,
 	      VISU::TMEDMeshOnEntity& theMeshOnEntity,
+	      MED::EEntiteMaillage theMEntity,
 	      const MED::TGeom2Size& theGeom2Size,
 	      VISU::TMEDValForTime& theValForTime)
   {
@@ -228,8 +316,11 @@ namespace
     
     TProfileMap& aProfileMap = theMeshOnEntity.myProfileMap;
     
-    TProfileKey aProfileKey = GetProfileKey(theTimeStampVal,
+    TProfileKey aProfileKey = GetProfileKey(theMEDWrapper,
+					    theMeshInfo,
+					    theTimeStampVal,
 					    theMeshOnEntity,
+					    theMEntity,
 					    theGeom2Size);
     
     TProfileMap::const_iterator anIter = aProfileMap.find(aProfileKey);
@@ -383,15 +474,29 @@ namespace
 
   //---------------------------------------------------------------
   void
-  InitGaussProfile(MED::TTimeStampVal& theTimeStampVal,
+  InitGaussProfile(const MED::PWrapper& theMEDWrapper,
+		   const MED::PMeshInfo& theMeshInfo,
+		   MED::TTimeStampVal& theTimeStampVal,
 		   VISU::TMEDMeshOnEntity& theMeshOnEntity,
+		   MED::EEntiteMaillage theMEntity,
 		   const MED::TGeom2Size& theGeom2Size,
 		   VISU::TMEDValForTime& theValForTime)
   {
     INITMSG(MYDEBUG,"InitGaussProfile"<<endl);
+
     // The order of the function calls is important
-    InitProfile(theTimeStampVal,theMeshOnEntity,theGeom2Size,theValForTime);
-    InitGaussMesh(theTimeStampVal,theMeshOnEntity,theGeom2Size,theValForTime);
+    InitProfile(theMEDWrapper,
+		theMeshInfo,
+		theTimeStampVal,
+		theMeshOnEntity,
+		theMEntity,
+		theGeom2Size,
+		theValForTime);
+
+    InitGaussMesh(theTimeStampVal,
+		  theMeshOnEntity,
+		  theGeom2Size,
+		  theValForTime);
   }
 
 
@@ -495,7 +600,7 @@ namespace
 	  }
 	  default: {
 	    vtkIdType aVGeom = MEDGeomToVTK(aMGeom);
-	    int aVNbNodes = VTKGeom2NbNodes(aVGeom);
+	    vtkIdType aVNbNodes = VTKGeom2NbNodes(aVGeom);
 	    MED::PCellInfo aCellInfo = theMEDWrapper->GetPCellInfo(aMeshInfo,aMEntity,aMGeom);
 	    TInt aNbElem = aCellInfo->GetNbElem();
 	    aMeshOnEntity->myNbCells += aNbElem;
@@ -809,8 +914,11 @@ namespace
 	aValForTime->myTime = VISU::TTime(aDt,anUnitDt);
 	INITMSG(MYDEBUG,"aDt = '"<<aDt<<", "<<anUnitDt<<"'\n");
 	
-	InitGaussProfile(aTimeStampVal,
+	InitGaussProfile(theMEDWrapper,
+			 aMeshInfo,
+			 aTimeStampVal,
 			 aMeshOnEntity,
+			 aMEntity,
 			 aGeom2Size,
 			 aValForTime);
       }
@@ -1412,7 +1520,7 @@ VISU_MedConvertor
 
     for(; aGeom2SizeIter != aGeom2Size.end(); aGeom2SizeIter++){
       const MED::EGeometrieElement& aMGeom = aGeom2SizeIter->first;
-      int aVGeom = MEDGeomToVTK(aMGeom);
+      vtkIdType aVGeom = MEDGeomToVTK(aMGeom);
       INITMSG(MYDEBUG,"aMGeom = "<<aMGeom<<"\n");
       switch(aMGeom){
       case MED::ePOLYGONE: {
@@ -1480,7 +1588,7 @@ VISU_MedConvertor
 	break;
       }
       default: {
-	TInt aVNbNodes = VTKGeom2NbNodes(aVGeom);
+	vtkIdType aVNbNodes = VTKGeom2NbNodes(aVGeom);
 	
 	MED::PCellInfo aCellInfo = theMed->GetPCellInfo(aMeshInfo,aMEntity,aMGeom);
 	TInt aNbElem = aCellInfo->GetNbElem();
@@ -1514,15 +1622,29 @@ VISU_MedConvertor
 	    }
 	    
 	    switch(aMGeom){
-	    case MED::eTETRA4:
 	    case MED::eTETRA10:
+	      anArray[4] = aConnect[4];
+	      anArray[5] = aConnect[5];
+	      anArray[6] = aConnect[7];  
+	      anArray[7] = aConnect[6];  
+	      anArray[8] = aConnect[8];  
+	      anArray[9] = aConnect[9];  
+	    case MED::eTETRA4:
 	      anArray[0] = aConnect[0];
 	      anArray[1] = aConnect[1];
 	      anArray[2] = aConnect[3];  
 	      anArray[3] = aConnect[2];  
 	      break;
-	    case MED::ePYRA5:
 	    case MED::ePYRA13:
+	      anArray[5] = aConnect[5];
+	      anArray[6] = aConnect[8];  
+	      anArray[7] = aConnect[7];
+	      anArray[8] = aConnect[6];  
+	      anArray[9] = aConnect[9];
+	      anArray[10] = aConnect[10];
+	      anArray[11] = aConnect[11];
+	      anArray[12] = aConnect[12];
+	    case MED::ePYRA5:
 	      anArray[0] = aConnect[0];
 	      anArray[1] = aConnect[3];  
 	      anArray[2] = aConnect[2];
-- 
2.39.2