1 /*  This file is part of MED.
2  *
3  *  COPYRIGHT (C) 1999 - 2019  EDF R&D, CEA/DEN
4  *  MED is free software: you can redistribute it and/or modify
5  *  it under the terms of the GNU Lesser General Public License as published by
6  *  the Free Software Foundation, either version 3 of the License, or
7  *  (at your option) any later version.
8  *
9  *  MED is distributed in the hope that it will be useful,
10  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  *  GNU Lesser General Public License for more details.
13  *
14  *  You should have received a copy of the GNU Lesser General Public License
15  *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
16  */
17 
18 /*
19  *  Use case 12 : read a 2D unstructured mesh with moving grid (generic approach)
20  */
21 
22 #include <med.h>
23 #define MESGERR 1
24 #include <med_utils.h>
25 
26 #include <string.h>
27 
main(int argc,char ** argv)28 int main (int argc, char **argv) {
29   med_idt fid;
30   med_int nmesh;
31   char meshname[MED_NAME_SIZE+1]="";
32   char meshdescription[MED_COMMENT_SIZE+1]="";
33   med_int meshdim;
34   med_int spacedim;
35   med_sorting_type sortingtype;
36   med_int nstep;
37   med_mesh_type meshtype;
38   med_axis_type axistype;
39   char *axisname;
40   char *unitname;
41   char dtunit[MED_SNAME_SIZE+1]="";
42   med_float *coordinates = NULL;
43   med_int ngeo = 0;
44   med_int nnodes = 0;
45   med_int *connectivity = NULL;
46   med_bool coordinatechangement;
47   med_bool geotransformation;
48   med_bool matrixtransformation;
49   med_int matrixsize;
50   med_float matrix[7]={ 0, 0, 0, 0, 0, 0, 0};
51   int i, it, j;
52   med_int profilesize;
53   char profilename[MED_NAME_SIZE+1]="";
54   med_int numdt, numit;
55   med_float dt;
56   med_geometry_type geotype;
57   med_geometry_type *geotypes = MED_GET_CELL_GEOMETRY_TYPE;
58   int ret=-1;
59 
60   /* open MED file with READ ONLY access mode */
61   fid = MEDfileOpen("UsesCase_MEDmesh_9.med",MED_ACC_RDONLY);
62   if (fid < 0) {
63     MESSAGE("ERROR : open file in READ ONLY ACCESS mode ...");
64     goto ERROR;
65   }
66 
67 
68   /* read how many mesh in the file */
69   if ((nmesh = MEDnMesh(fid)) < 0) {
70     MESSAGE("ERROR : read how many mesh ...");
71     goto ERROR;
72   }
73 
74   for (i=0;i<nmesh;i++) {
75 
76     /* read computation space dimension */
77     if ((spacedim = MEDmeshnAxis(fid, i+1)) < 0) {
78       MESSAGE("ERROR : read computation space dimension ...");
79       goto ERROR;
80     }
81 
82     /* memory allocation */
83     if ((axisname  = (char*) malloc(MED_SNAME_SIZE*spacedim+1)) == NULL) {
84       MESSAGE("ERROR : memory allocation ...");
85       goto ERROR;
86     }
87     if ((unitname  = (char*) malloc(MED_SNAME_SIZE*spacedim+1)) == NULL) {
88       MESSAGE("ERROR : memory allocation ...");
89       goto ERROR;
90     }
91 
92     /* read mesh informations : meshname, mesh dimension, mesh type ... */
93     if (MEDmeshInfo(fid, i+1, meshname, &spacedim, &meshdim, &meshtype, meshdescription,
94 		    dtunit, &sortingtype, &nstep,
95 		    &axistype, axisname, unitname) < 0) {
96       MESSAGE("ERROR : mesh info ...");
97       free(axisname);
98       free(unitname);
99       goto ERROR;
100     }
101     free(axisname);
102     free(unitname);
103 
104 
105     /* read how many nodes in the mesh */
106     if ((nnodes = MEDmeshnEntity(fid, meshname, MED_NO_DT, MED_NO_IT, MED_NODE, MED_NONE,
107 				 MED_COORDINATE, MED_NO_CMODE,&coordinatechangement,
108 				 &geotransformation)) < 0) {
109       MESSAGE("ERROR : number of nodes ...");
110       goto ERROR;
111     }
112 
113     /* read mesh nodes coordinates */
114     if ((coordinates = (med_float*) malloc(sizeof(med_float)*nnodes*spacedim)) == NULL) {
115       MESSAGE("ERROR : memory allocation ...");
116       goto ERROR;
117     }
118 
119     if (MEDmeshNodeCoordinateRd(fid, meshname, MED_NO_DT, MED_NO_IT, MED_FULL_INTERLACE,
120 				coordinates) < 0) {
121       MESSAGE("ERROR : nodes coordinates ...");
122       free(coordinates);
123       goto ERROR;
124     }
125 
126     /* read all MED geometry cell types */
127     for (it=1; it<= MED_N_CELL_FIXED_GEO; it++) {
128 
129       geotype = geotypes[it];
130 
131       if ((ngeo = MEDmeshnEntity(fid, meshname, MED_NO_DT, MED_NO_IT, MED_CELL,geotype,
132 				 MED_CONNECTIVITY, MED_NODAL, &coordinatechangement,
133 				 &geotransformation)) < 0) {
134 	MESSAGE("ERROR : number of cell ...");
135 	ISCRUTE(geotype);
136 	goto ERROR;
137       }
138 
139       if (ngeo) {
140 	/* read cells connectivity in the mesh */
141 	if ((connectivity = (med_int *) malloc(sizeof(med_int)*ngeo*(geotype%100))) == NULL) {
142 	  MESSAGE("ERROR : memory allocation ...");
143 	  goto ERROR;
144 	}
145 
146 	if (MEDmeshElementConnectivityRd(fid, meshname, MED_NO_DT, MED_NO_IT, MED_CELL,
147 					 geotype, MED_NODAL, MED_FULL_INTERLACE, connectivity) < 0) {
148 	  MESSAGE("ERROR : cell connectivity ...");
149 	  ISCRUTE(geotype);
150 	  free(connectivity);
151 	  goto ERROR;
152 	}
153 
154 	/* memory deallocation */
155 	free(connectivity);
156 	connectivity = NULL;
157       }
158     }
159 
160 
161     /* read nodes coordinates changements step by step */
162     for (it=1;it<nstep;it++) {
163 
164       if (MEDmeshComputationStepInfo(fid, meshname, it+1,
165 				     &numdt, &numit, &dt) < 0) {
166 	MESSAGE("ERROR : Computing step info ...");
167 	SSCRUTE(meshname);
168 	goto ERROR;
169       }
170 
171       /* test changement : for nodes coordinates */
172       if ((nnodes = MEDmeshnEntityWithProfile(fid, meshname, numdt, numit,
173 					      MED_NODE, MED_NONE,
174 					      MED_COORDINATE, MED_NO_CMODE,
175 					      MED_GLOBAL_STMODE, profilename, &profilesize,
176 					      &coordinatechangement, &geotransformation)) < 0) {
177 	MESSAGE("ERROR : number of nodes ...");
178 	goto ERROR;
179       }
180 
181       /* if only coordinates have changed, then read the new coordinates */
182       if (coordinatechangement && geotransformation) {
183 	if (MEDmeshNodeCoordinateWithProfileRd(fid, meshname, numdt, numit,
184 					       MED_GLOBAL_STMODE,profilename,
185 					       MED_FULL_INTERLACE,MED_ALL_CONSTITUENT,
186 					       coordinates) < 0) {
187 	  MESSAGE("ERROR : nodes coordinates ...");
188 	  free(coordinates);
189 	  goto ERROR;
190 	}
191       }
192 
193       if (coordinatechangement && ! geotransformation) {
194 
195 	matrixsize = MEDmeshnEntity(fid,meshname,numdt,numit,
196 				    MED_NODE,MED_NONE,MED_COORDINATE_TRSF,MED_NODAL,&coordinatechangement,
197 				    &matrixtransformation);
198 
199 	if (matrixsize < 0) {
200 	  MESSAGE("ERROR : matrix transformation ...");
201 	  goto ERROR;
202 	}
203 
204 	if (matrixtransformation) {
205 	  if ( MEDmeshNodeCoordinateTrsfRd(fid, meshname, numdt, numit, matrix) < 0 ) {
206 	    MESSAGE("ERROR : read transformation matrix ...");
207 	    goto ERROR;
208 	  }
209 	}
210       }
211     }
212   }
213   free(coordinates);
214 
215   ret=0;
216  ERROR:
217 
218   /* close MED file */
219   if (MEDfileClose(fid) < 0) {
220     MESSAGE("ERROR : close file");
221     ret=-1;
222   }
223 
224 
225   return ret;
226 }
227 
228