7 #include "MEDMEM_Exception.hxx"
8 #include "MEDMEM_Mesh.hxx"
9 #include "MEDMEM_Family.hxx"
10 #include "MEDMEM_Group.hxx"
12 #include "MEDMEM_MedMeshDriver.hxx"
13 #include "MEDMEM_MedFieldDriver.hxx"
14 #include "MEDMEM_Support.hxx"
15 #include "MEDMEM_Field.hxx"
16 #include "MEDMEM_define.hxx"
18 using namespace MEDMEM;
20 double dmax(double x, double y) { return (x>y)?x:y;}
22 double dmin(double x, double y) { return (x>y)?y:x;}
26 void affiche_support(const SUPPORT * mySupport)
28 cout << " - Name : "<<mySupport->getName().c_str()<<endl ;
29 cout << " - Description : "<<mySupport->getDescription().c_str()<<endl ;
30 cout << " - Entity : "<<mySupport->getEntity()<<endl ;
31 cout << " - Entities list : "<<endl ;
32 if (!(mySupport->isOnAllElements())) {
33 int NumberOfTypes = mySupport->getNumberOfTypes() ;
34 cout<<" - NumberOfTypes : "<<NumberOfTypes<<endl;
35 const medGeometryElement * Types = mySupport->getTypes() ;
36 for (int j=0;j<NumberOfTypes;j++) {
37 cout << " * Type "<<Types[j]<<" : " ;
38 int NumberOfElements = mySupport->getNumberOfElements(Types[j]) ;
39 const int * Number = mySupport->getNumber(Types[j]) ;
40 for (int k=0; k<NumberOfElements;k++)
41 cout << Number[k] << " ";
45 cout << " Is on all entities !"<< endl;
49 void affiche_famille(MESH *myMesh,medEntityMesh Entity)
51 int NumberOfFamilies = myMesh->getNumberOfFamilies(Entity) ;
52 cout << "NumberOfFamilies : "<<NumberOfFamilies<<endl;
53 for (int i=1; i<NumberOfFamilies+1;i++) {
54 const FAMILY* myFamily = myMesh->getFamily(Entity,i);
55 affiche_support(myFamily);
56 cout << " - Identifier : "<<myFamily->getIdentifier()<<endl ;
57 int NumberOfAttributes = myFamily->getNumberOfAttributes() ;
58 cout << " - Attributes ("<<NumberOfAttributes<<") :"<<endl;
59 for (int j=1;j<NumberOfAttributes+1;j++)
60 cout << " * "<<myFamily->getAttributeIdentifier(j)<<" : "<<myFamily->getAttributeValue(j)<<", "<<myFamily->getAttributeDescription(j).c_str()<<endl ;
61 int NumberOfGroups = myFamily->getNumberOfGroups() ;
62 cout << " - Groups ("<<NumberOfGroups<<") :"<<endl;
63 for (int j=1;j<NumberOfGroups+1;j++)
64 cout << " * "<<myFamily->getGroupName(j).c_str()<<endl ;
68 void affiche_groupe(MESH *myMesh,medEntityMesh Entity)
70 int NumberOfGroups = myMesh->getNumberOfGroups(Entity) ;
71 cout << "NumberOfGroups : "<<NumberOfGroups<<endl;
72 for (int i=1; i<NumberOfGroups+1;i++) {
73 const GROUP* myGroup = myMesh->getGroup(Entity,i);
74 affiche_support(myGroup);
75 int NumberOfFamillies = myGroup->getNumberOfFamilies() ;
76 cout << " - Families ("<<NumberOfFamillies<<") :"<<endl;
77 for (int j=1;j<NumberOfFamillies+1;j++)
78 cout << " * "<<myGroup->getFamily(j)->getName().c_str()<<endl ;
82 int main (int argc, char ** argv) {
86 if ((argc !=3) && (argc != 4)) {
87 cerr << "Usage : " << argv[0]
88 << " filename meshname [fieldname]" << endl << endl;
92 string filename = argv[1] ;
93 string meshname = argv[2] ;
95 // MESH * myMesh= new MESH(MED_DRIVER,filename,meshname) ;
96 MESH * myMesh= new MESH() ;
97 myMesh->setName(meshname);
98 MED_MESH_RDONLY_DRIVER myMeshDriver(filename,myMesh) ;
99 myMeshDriver.setMeshName(meshname);
100 myMeshDriver.open() ;
101 myMeshDriver.read() ;
102 myMeshDriver.close() ;
104 // int drv = myMesh->addDriver(MED_DRIVER,"sortie.med",meshname);
105 // myMesh->write(drv);
108 int SpaceDimension = myMesh->getSpaceDimension() ;
109 int MeshDimension = myMesh->getMeshDimension() ;
110 int NumberOfNodes = myMesh->getNumberOfNodes() ;
112 cout << "Space Dimension : " << SpaceDimension << endl << endl ;
114 cout << "Mesh Dimension : " << MeshDimension << endl << endl ;
116 const double * Coordinates = myMesh->getCoordinates(MED_FULL_INTERLACE) ;
118 cout << "Show Nodes Coordinates : " << endl ;
120 cout << "Name :" << endl ;
121 const string * CoordinatesNames = myMesh->getCoordinatesNames() ;
122 for(int i=0; i<SpaceDimension ; i++) {
123 cout << " - " << CoordinatesNames[i] << endl ;
125 cout << "Unit :" << endl ;
126 const string * CoordinatesUnits = myMesh->getCoordinatesUnits() ;
127 for(int i=0; i<SpaceDimension ; i++) {
128 cout << " - " << CoordinatesUnits[i] << endl ;
130 for(int i=0; i<NumberOfNodes ; i++) {
131 cout << "Nodes " << i+1 << " : " ;
132 for (int j=0; j<SpaceDimension ; j++)
133 cout << Coordinates[i*SpaceDimension+j] << " " ;
137 int NumberOfTypes = myMesh->getNumberOfTypes(MED_CELL) ;
138 const medGeometryElement * Types = myMesh->getTypes(MED_CELL) ;
140 cout << "Show Connectivity (Nodal) :" << endl ;
141 for (int i=0; i<NumberOfTypes; i++) {
142 cout << "For type " << Types[i] << " : " << endl ;
143 int NumberOfElements = myMesh->getNumberOfElements(MED_CELL,Types[i]);
144 const int * connectivity = myMesh->getConnectivity(MED_FULL_INTERLACE,MED_NODAL,MED_CELL,Types[i]);
145 int NomberOfNodesPerCell = Types[i]%100 ;
146 for (int j=0;j<NumberOfElements;j++){
147 cout << "Element "<< j+1 <<" : " ;
148 for (int k=0;k<NomberOfNodesPerCell;k++)
149 cout << connectivity[j*NomberOfNodesPerCell+k]<<" ";
154 cout << "Show Family :"<<endl ;
155 affiche_famille(myMesh,MED_NODE);
156 affiche_famille(myMesh,MED_CELL);
157 affiche_famille(myMesh,MED_FACE);
158 affiche_famille(myMesh,MED_EDGE);
160 cout << "Show Group :"<<endl ;
161 affiche_groupe(myMesh,MED_NODE);
162 affiche_groupe(myMesh,MED_CELL);
163 affiche_groupe(myMesh,MED_FACE);
164 affiche_groupe(myMesh,MED_EDGE);
166 cout << "Show Reverse Nodal Connectivity :" << endl ;
167 const int * ReverseNodalConnectivity = myMesh->getReverseConnectivity(MED_NODAL) ;
168 const int * ReverseNodalConnectivityIndex = myMesh->getReverseConnectivityIndex(MED_NODAL) ;
169 for (int i=0; i<NumberOfNodes; i++) {
170 cout << "Node "<<i+1<<" : " ;
171 for (int j=ReverseNodalConnectivityIndex[i];j<ReverseNodalConnectivityIndex[i+1];j++)
172 cout << ReverseNodalConnectivity[j-1] << " " ;
176 cout << "Show Connectivity (Descending) :" << endl ;
177 int NumberOfElements ;
178 const int * connectivity ;
179 const int * connectivity_index ;
180 myMesh->calculateConnectivity(MED_FULL_INTERLACE,MED_DESCENDING,MED_CELL);
182 NumberOfElements = myMesh->getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS);
183 connectivity = myMesh->getConnectivity(MED_FULL_INTERLACE,MED_DESCENDING,MED_CELL,MED_ALL_ELEMENTS);
184 connectivity_index = myMesh->getConnectivityIndex(MED_DESCENDING,MED_CELL);
186 catch (MEDEXCEPTION m) {
187 cout << m.what() << endl ;
190 for (int j=0;j<NumberOfElements;j++) {
191 cout << "Element "<<j+1<<" : " ;
192 for (int k=connectivity_index[j];k<connectivity_index[j+1];k++)
193 cout << connectivity[k-1]<<" ";
197 cout << "Show Reverse Descending Connectivity :" << endl ;
198 const int * ReverseDescendingConnectivity = myMesh->getReverseConnectivity(MED_DESCENDING) ;
199 const int * ReverseDescendingConnectivityIndex = myMesh->getReverseConnectivityIndex(MED_DESCENDING) ;
201 int NumberOfConstituents = 0;
203 medEntityMesh constituentEntity ;
205 if (MeshDimension==3) {
206 constituent = "Face" ;
207 constituentEntity = MED_FACE ;
210 if (MeshDimension==2) {
211 constituent = "Edge" ;
212 constituentEntity = MED_EDGE ;
215 if (MeshDimension==1) {
216 INFOS("ERROR : MeshDimension = 1 !");
217 INFOS("We could not see Reverse Descending Connectivity.") ;
219 NumberOfConstituents = myMesh->getNumberOfElements (constituentEntity,MED_ALL_ELEMENTS);
220 for (int i=0; i<NumberOfConstituents; i++) {
221 cout << constituent <<i+1<<" : " ;
222 for (int j=ReverseDescendingConnectivityIndex[i];j<ReverseDescendingConnectivityIndex[i+1];j++)
223 cout << ReverseDescendingConnectivity[j-1] << " " ;
227 cout << "Show "<<constituent<<" Connectivity (Nodal) :" << endl ;
228 const int * face_connectivity = myMesh->getConnectivity(MED_FULL_INTERLACE,MED_NODAL,constituentEntity,MED_ALL_ELEMENTS);
229 const int * face_connectivity_index = myMesh->getConnectivityIndex(MED_NODAL,constituentEntity);
230 for (int i=0; i<NumberOfConstituents; i++) {
231 cout << constituent <<i+1<<" : " ;
232 for (int j=face_connectivity_index[i];j<face_connectivity_index[i+1];j++)
233 cout << face_connectivity[j-1]<<" ";
237 /* test of normal, area, volume, barycenter */
239 SUPPORT * support1 = (SUPPORT*) NULL;
241 //FIELD<double>* normal = new FIELD<double>::FIELD();
242 //FIELD<double>* length = new FIELD<double>::FIELD();
245 string support_name = "Support on all " ;
246 support_name+=constituent;
247 support1 = new SUPPORT(myMesh,support_name,constituentEntity);
248 cout << "Building of the Support on all cells dimensionned (Meshdim-1) of the mesh :"<< endl ;
249 cout << "Face in 3D or Edge in 2D" << endl;
251 cout << "Getting the normal of each face of this support !" << endl ;
253 FIELD<double>* normal = myMesh->getNormal(support1);
255 double normal_square, norm ;
256 double maxnorm=-infty;
257 double minnorm=infty;
259 for (int i = 1; i<=NumberOfConstituents;i++) {
261 cout << "Normal " << i << " " ;
262 for (int j=1; j<=SpaceDimension; j++) {
263 tmp_value = normal->getValueIJ(i,j) ;
264 normal_square += tmp_value*tmp_value ;
265 cout << tmp_value << " " ;
267 norm = sqrt(normal_square);
268 maxnorm = dmax(maxnorm,norm);
269 minnorm = dmin(minnorm,norm);
270 cout << ", Norm = " << norm << endl;
272 cout << "Max Norm " << maxnorm << " Min Norm " << minnorm << endl;
276 if (SpaceDimension == 2)
278 cout << "Getting the length of each edge !" << endl ;
280 FIELD<double>* length = myMesh->getLength(support1);
282 double length_value,maxlength,minlength;
285 for (int i = 1; i<=NumberOfConstituents;i++)
287 length_value = length->getValueIJ(i,1) ;
288 cout << "Length " << i << " " << length_value << endl;
289 maxlength = dmax(maxlength,length_value);
290 minlength = dmin(minlength,length_value);
292 cout << "Max Length " << maxlength << " Min Length " << minlength << endl;
299 cout << "Building of the Support on all space-dimensionned cells of the mesh :"<< endl ;
300 SUPPORT * support = new SUPPORT(myMesh);
302 cout << "Getting the barycenter of each element of this support !" << endl ;
304 //FIELD<double>* barycenter = new FIELD<double>::FIELD();
306 FIELD<double>* barycenter = myMesh->getBarycenter(support);
307 NumberOfElements = myMesh->getNumberOfElements(MED_CELL,MED_ALL_ELEMENTS);
309 for (int i = 1; i<=NumberOfElements;i++)
311 if (SpaceDimension == 3)
312 cout << "Barycenter " << i << " " << barycenter->getValueIJ(i,1) << " " << barycenter->getValueIJ(i,2) << " " << barycenter->getValueIJ(i,3) << endl;
314 if (SpaceDimension == 2)
315 cout << "Barycenter " << i << " " << barycenter->getValueIJ(i,1) << " " << barycenter->getValueIJ(i,2) << endl;
320 //FIELD<double>* volume = new FIELD<double>::FIELD();
321 //FIELD<double>* area = new FIELD<double>::FIELD();
325 if (SpaceDimension == 3)
327 cout << "Getting the Volume of each element of this support which is a 3D one !" << endl;
329 FIELD<double>* volume = myMesh->getVolume(support);
331 double maxvol,minvol,voltot;
335 for (int i = 1; i<=NumberOfElements;i++)
337 cout << "Volume " << i << " " << volume->getValueIJ(i,1) << endl;
338 maxvol = dmax(maxvol,volume->getValueIJ(i,1));
339 minvol = dmin(minvol,volume->getValueIJ(i,1));
340 voltot = voltot + volume->getValueIJ(i,1);
343 cout << "Max Volume " << maxvol << " Min Volume " << minvol << endl;
344 cout << "Support Volume " << voltot << endl;
348 else if (SpaceDimension == 2)
350 cout << "Getting the Area of each element of this support which is a 2D one !" << endl;
352 FIELD<double>* area = myMesh->getArea(support);
354 // cout << "nb of comp "<< area->getNumberOfComponents() << " length " << area->getSupport()->getNumberOfElements(MED_ALL_ELEMENTS) << endl;
356 double maxarea,minarea,areatot;
360 for (int i = 1; i<=NumberOfElements;i++)
362 cout << "Area " << i << " " << area->getValueIJ(i,1) << endl;
363 maxarea = dmax(maxarea,area->getValueIJ(i,1));
364 minarea = dmin(minarea,area->getValueIJ(i,1));
365 areatot = areatot + area->getValueIJ(i,1);
368 cout << "Max Area " << maxarea << " Min Area " << minarea << endl;
369 cout << "Support Area " << areatot << endl;
376 //if (barycenter != NULL) delete barycenter;
377 //if (volume != NULL ) delete volume;
378 //if (area != NULL ) delete area;
381 if (argc < 4) return 0;
385 if (argc != 4) exit(0) ;
386 // else we have a field !
388 string fieldname = argv[3];
390 // SUPPORT * mySupport = new SUPPORT(myMesh,"On_all_node",MED_NODE);
391 SUPPORT * mySupport = new SUPPORT(myMesh,"On_all_cell",MED_CELL);
392 FIELD<double> * myField= new FIELD<double>() ;
394 myField->setName(fieldname);
395 myField->setSupport(mySupport);
396 MED_FIELD_RDONLY_DRIVER<double> myFieldDriver(filename,myField) ;
397 myFieldDriver.setFieldName(fieldname);
398 myFieldDriver.open() ;
401 myFieldDriver.read() ;
404 mySupport = new SUPPORT(myMesh,"On_all_node",MED_NODE);
405 myField->setSupport(mySupport);
407 myFieldDriver.read() ;
409 cout << "Field " << fieldname << " not found !!!" << endl ;
414 myFieldDriver.close() ;
416 cout << "Field "<< myField->getName() << " : " <<myField->getDescription() << endl ;
417 int NumberOfComponents = myField->getNumberOfComponents() ;
418 cout << "- Nombre de composantes : "<< NumberOfComponents << endl ;
419 for (int i=1; i<NumberOfComponents+1; i++) {
420 cout << " - composante "<<i<<" :"<<endl ;
421 cout << " - nom : "<<myField->getComponentName(i)<< endl;
422 cout << " - description : "<<myField->getComponentDescription(i) << endl;
423 cout << " - units : "<<myField->getMEDComponentUnit(i) << endl;
425 cout << "- iteration :" << endl ;
426 cout << " - numero : " << myField->getIterationNumber()<< endl ;
427 cout << " - ordre : " << myField->getOrderNumber()<< endl ;
428 cout << " - temps : " << myField->getTime()<< endl ;
430 cout << "- Valeurs :"<<endl;
431 int NumberOf = mySupport->getNumberOfElements(MED_ALL_ELEMENTS);
432 // for (int i=1; i<NumberOfComponents+1; i++) {
433 // double * value = myField->getValueI(MED_NO_INTERLACE,i) ;
434 // for (int j=0; j<NumberOf; j++)
435 // cout << value[j]<< " ";
438 MEDARRAY<double> * myvalue = myField->getvalue();
439 const double * value ;
440 for (int i=1; i<NumberOf+1; i++) {
441 value = myvalue->getRow(i) ;
442 for (int j=0; j<NumberOfComponents; j++)
443 cout << value[j]<< " ";