1 // Copyright (C) 2010-2013 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 * TestMedParallelWrite.cxx
23 * Created on: 11 mai 2011
31 #include "med_utils.h"
37 int main (int argc, char **argv)
39 const char meshname[MED_NAME_SIZE+1] = "2D unstructured mesh";
40 const med_int spacedim = 2;
41 const med_int meshdim = 2;
42 /* 12345678901234561234567890123456 */
43 const char axisname[2*MED_SNAME_SIZE+1] = "x y ";
44 const char unitname[2*MED_SNAME_SIZE+1] = "cm cm ";
45 med_float coordinates[2222];
46 const med_int nnodes = 1111;
48 med_int* quadconnectivity;
49 const med_int nquad4 = 1000;
53 int mpi_size, mpi_rank;
54 MPI_Comm comm = MPI_COMM_WORLD;
55 MPI_Info info = MPI_INFO_NULL;
57 MPI_Init(&argc, &argv);
58 MPI_Comm_size(MPI_COMM_WORLD, &mpi_size);
59 MPI_Comm_rank(MPI_COMM_WORLD, &mpi_rank);
62 char filename[255]="UsesCase_MEDmesh_parallel.med";
63 /* SSCRUTE(_filename); */
66 printf("mpi_size = %03d\n", mpi_size);
68 /* Ouverture du fichier en mode non-parallel */
69 if ((fid = MEDfileOpen(filename, MED_ACC_CREAT)) < 0){
70 MED_ERR_(_ret, MED_ERR_OPEN,MED_ERR_FILE, filename);
73 /* write a comment in the file */
74 if (MEDfileCommentWr(fid,"A 2D unstructured mesh : 15 nodes, 12 cells") < 0) {
75 MESSAGE("ERROR : write file description ...");
78 /* mesh creation : a 2D unstructured mesh */
79 if (MEDmeshCr(fid, meshname, spacedim, meshdim, MED_UNSTRUCTURED_MESH,
80 "A 2D unstructured mesh","",MED_SORT_DTIT,MED_CARTESIAN, axisname, unitname) < 0) {
81 MESSAGE("ERROR : mesh creation ...");
85 * Building the coordinates of a rectangle of 101 points in the Y-axis,
86 * and 11 in the X-axis
88 for (int j=0; j<11; j++ )
89 for (int i=0; i<101; i++ )
91 coordinates[j*202+i*2] = j+1;
92 coordinates[j*202+i*2+1] = i+1;
95 /* nodes coordinates in a Cartesian axis in full interlace mode
96 (X1,Y1, X2,Y2, X3,Y3, ...) with no iteration and computation step
98 if (MEDmeshNodeCoordinateWr(fid, meshname, MED_NO_DT, MED_NO_IT, 0.0,
99 MED_FULL_INTERLACE, nnodes, coordinates) < 0) {
100 MESSAGE("ERROR : nodes coordinates ...");
103 if ( MEDfileClose( fid ) < 0) {
104 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,"");;
108 } /* End of process ZERO */
114 /* Ouverture du fichier en mode parallel */
115 if ((fid = MEDparFileOpen(filename, MED_ACC_RDWR ,comm, info)) < 0){
116 MED_ERR_(_ret, MED_ERR_OPEN,MED_ERR_FILE, filename);
119 med_int nbofentity = nquad4;
120 med_int nbofvaluesperentity = 1;
121 med_int nbofconstituentpervalue = 4;
122 med_int constituentselect = MED_ALL_CONSTITUENT;
123 med_switch_mode switchmode = MED_FULL_INTERLACE;
124 med_storage_mode storagemode = MED_COMPACT_STMODE;
125 const char *const profilename = MED_NO_PROFILE;
128 * Calculating block sizes
131 int block_size = (100/mpi_size)*10;
132 med_size start = block_size * mpi_rank + 1;
133 med_size stride = block_size;
135 med_size blocksize = block_size;
136 med_size lastblocksize = (100 % mpi_size)*10;
137 if ((mpi_size == mpi_rank+1) && (lastblocksize != 0))
139 blocksize += lastblocksize;
140 stride += lastblocksize;
144 printf("%03d: block_size = %03d\n", mpi_rank, block_size);
145 printf("%03d: start = %03d\n", mpi_rank, start);
146 printf("%03d: stride = %03d\n", mpi_rank, stride);
147 printf("%03d: count = %03d\n", mpi_rank, count);
148 printf("%03d: blocksize = %03d\n", mpi_rank, blocksize);
149 printf("%03d: lastblocksize = %03d\n", mpi_rank, lastblocksize);
150 med_filter filter = MED_FILTER_INIT;
152 if ( MEDfilterBlockOfEntityCr( fid,
155 nbofconstituentpervalue,
167 MESSAGE("ERROR : filter creation ...");
170 // Attention: there is blocksize and block_size and it does not
171 // represent the same quantity, in case we are in the last
172 // block they are different, if not it is the same
173 quadconnectivity = new med_int[blocksize*4];
174 int shift = mpi_rank*block_size;
175 printf("%03d: mpi_rank*block_size = %03d\n", mpi_rank, shift);
176 printf("%03d: block_size = %03d\n", mpi_rank, block_size);
177 int base = shift + shift / 101;
179 for (int i=0; i<blocksize*4; i+=4 )
185 quadconnectivity[i] = base;
186 quadconnectivity[i+1] = base+1;
187 quadconnectivity[i+2] = base+102;
188 quadconnectivity[i+3] = base+101;
191 printf("%03d: number of written quads = %03d\n", mpi_rank, c);
193 if (MEDmeshElementConnectivityAdvancedWr(fid, meshname, MED_NO_DT,
194 MED_NO_IT, 0.0, MED_CELL, MED_QUAD4,
195 MED_NODAL, &filter, quadconnectivity) < 0) {
196 MESSAGE("ERROR : quadrangular cells connectivity ...");
199 if ( MEDfilterClose( &filter ) < 0) {
200 MESSAGE("ERROR : filter closing ...");
203 if ( MEDfileClose( fid ) < 0) {
204 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,""); _ret = -1;
207 /* Barrier before writing family ZERO */
210 if (mpi_rank == 0 ) {
212 if ((fid = MEDfileOpen(filename, MED_ACC_RDWR)) < 0){
213 MED_ERR_(_ret, MED_ERR_OPEN,MED_ERR_FILE, filename);
216 /* create family 0 : by default, all mesh entities family number is 0 */
217 if (MEDfamilyCr(fid, meshname,MED_NO_NAME, 0, 0, MED_NO_GROUP) < 0) {
218 MESSAGE("ERROR : family 0 creation ...");
221 const char familyname_root[MED_NAME_SIZE+1] = "PROCESSOR ";
222 char familyname[MED_NAME_SIZE+1] = " ";
223 for (int i=1; i<mpi_size+1; i++)
225 snprintf(familyname, sizeof familyname, "%s%d", familyname_root, i);
226 if (MEDfamilyCr(fid, meshname,familyname, -i, 0, MED_NO_GROUP) < 0) {
227 MESSAGE("ERROR : family creation ...");
229 printf("%03d: %s\n", mpi_rank, familyname);
232 med_int familynumbers[nquad4];
234 for (int i=0; i<nquad4; i++)
236 if ((i > block_size * l - 1) && (l < mpi_size))
240 familynumbers[i] = -l;
243 if (MEDmeshEntityFamilyNumberWr(fid, meshname, MED_NO_DT, MED_NO_IT,
244 MED_CELL, MED_QUAD4, nquad4, familynumbers) < 0) {
245 MESSAGE("ERROR : nodes family numbers ...");
248 /* Write a Profile */
249 const char profileName[MED_NAME_SIZE+1] = "QUAD4_PROFILE";
250 const med_int profilesize = 9;
251 med_int profilearray[9] = {1, 3, 5, 7, 9, 11, 13, 15, 17};
252 if (MEDprofileWr(fid, profileName, profilesize, profilearray ) < 0) {
253 MESSAGE("ERROR : nodes family numbers ...");
256 /* write localization for integration points */
257 const char localizationName[MED_NAME_SIZE+1] = "QUAD4_INTEGRATION_POINTS_4";
258 const med_float elementcoordinate[6] = {0.0, 0.0, 1.0, 0.0, 0.0,1.0};
259 const med_float iPointCoordinate[8] = {1.0/5, 1.0/5, 3.0/5, 1.0/5, 1.0/5, 3.0/5, 1.0/3, 1.0/3};
260 const med_float weight[4] = {1.0/8, 1.0/8, 1.0/8, 1.0/8};
261 med_int spacedim = 2;
263 if (MEDlocalizationWr(fid, localizationName, MED_QUAD4, spacedim,
264 elementcoordinate, MED_FULL_INTERLACE,
265 nipoint, iPointCoordinate, weight,
266 MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT) < 0) {
267 MESSAGE("ERROR : create family of integration points ...");
270 /* Writing a scalar Field on the Quads right here */
271 const char fieldname[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD";
272 const med_int ncomponent = 1;
273 const char componentname[MED_SNAME_SIZE+1] = "TEMPERATURE";
274 const char componentunit[MED_SNAME_SIZE+1] = "C";
276 if (MEDfieldCr(fid, fieldname, MED_FLOAT64,
277 ncomponent, componentname, componentunit,"",
279 MESSAGE("ERROR : create field");
282 /* write values at cell (QUADS) centers */
283 med_float quad4values[nquad4];
284 for (int i=0; i<nquad4; i++)
285 quad4values[i] = i%100 + 1;
287 med_float quad4values4[nquad4 * 4];
288 long int counter = 0;
289 for (int i=0; i<nquad4; i++)
291 quad4values[i] = i%100 + 1;
292 for (int j=0; j<4; j++)
294 quad4values4[counter] = quad4values[i];
298 if (MEDfieldValueWr(fid, fieldname, MED_NO_DT, MED_NO_IT, 0.0, MED_CELL,
299 MED_QUAD4, MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
300 nquad4, (unsigned char*) quad4values) < 0) {
301 MESSAGE("ERROR : write field values on MED_QUAD4");
304 const char fieldname2[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD_PGAUSS";
305 if (MEDfieldCr(fid, fieldname2, MED_FLOAT64,
306 ncomponent, componentname, componentunit,"",
308 MESSAGE("ERROR : create field");
311 if (MEDfieldValueWithProfileWr(
312 fid, fieldname2, MED_NO_DT, MED_NO_IT, 0.0, MED_CELL, MED_QUAD4,
313 MED_GLOBAL_PFLMODE, MED_NO_PROFILE, localizationName,
314 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
315 nquad4, (unsigned char*) quad4values4) < 0) {
316 MESSAGE("ERROR : write field values on MED_QUAD4");
319 const char fieldname3[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD_ELNO";
320 if (MEDfieldCr(fid, fieldname3, MED_FLOAT64,
321 ncomponent, componentname, componentunit,"",
323 MESSAGE("ERROR : create field");
326 if (MEDfieldValueWithProfileWr(
327 fid, fieldname3, MED_NO_DT, MED_NO_IT, 0.0, MED_CELL, MED_QUAD4,
328 MED_GLOBAL_PFLMODE, MED_NO_PROFILE, MED_GAUSS_ELNO,
329 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
330 nquad4, (unsigned char*) quad4values4) < 0) {
331 MESSAGE("ERROR : write field values on MED_QUAD4");
334 if ( MEDfileClose( fid ) < 0) {
335 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,"");;
338 printf("File UsesCase_MEDmesh_parallel.med has been generated.\n");
339 } /* End of process ZERO */
341 /* MPI_Finalize must be called AFTER MEDclose which may use MPI calls */