00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #define MAX(a,b) ((a) > (b) ? (a) : (b))
00019
00020 #include <med.h>
00021 #define MESGERR 1
00022 #include "med_utils.h"
00023 #include "med_config.h"
00024
00025 #include <stdlib.h>
00026 #include <string.h>
00027 #include <assert.h>
00028
00029 #if TIME_WITH_SYS_TIME
00030 # include <sys/time.h>
00031 # include <time.h>
00032 #else
00033 # if HAVE_SYS_TIME_H
00034 # include <sys/time.h>
00035 # else
00036 # include <time.h>
00037 # endif
00038 #endif
00039
00040 #ifndef HAVE_UNISTD_H
00041 #error "unistd.h reuired."
00042 #endif
00043
00044 #include "getBlocksOfEntitiesPartition.h"
00045 #include "generateDatas.h"
00046 #include "generateFilterArray.h"
00047
00048
00049 #ifdef DEF_LECT_ECR
00050 #define MODE_ACCES MED_ACC_RDWR
00051 #elif DEF_LECT_AJOUT
00052 #define MODE_ACCES MED_ACC_RDEXT
00053 #else
00054 #define MODE_ACCES MED_ACC_CREAT
00055 #endif
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080 int main (int argc, char **argv)
00081
00082
00083 {
00084 med_err _ret=0;
00085 med_idt _fid;
00086
00087 int mpi_size, mpi_rank;
00088 MPI_Comm comm = MPI_COMM_WORLD;
00089 MPI_Info info = MPI_INFO_NULL;
00090
00091 med_int _nentitiesfiltered=0;
00092 med_int *_filterarray=NULL;
00093
00094
00095 MPI_Init(&argc, &argv);
00096 MPI_Comm_size(MPI_COMM_WORLD, &mpi_size);
00097 MPI_Comm_rank(MPI_COMM_WORLD, &mpi_rank);
00098
00099
00100 med_err generateFieldFile( const med_size nentities, const med_size nvaluesperentity, const med_size nconstituentpervalue,
00101 const med_switch_mode constituentmode,GetBlocksOfEntitiesType getBlockOfEntities, const med_int nbblocksperproc,
00102 GenerateDataType generateDatas,
00103 const med_storage_mode storagemode, const med_size profilearraysize, const char * const fieldnameprefix ) {
00104
00105
00106 med_err _ret=-1;
00107 char _filename [255]="";
00108 char _meshname[MED_NAME_SIZE+1]="Empty mesh";
00109 med_int _meshdim=3;
00110 char _meshcomponentname[3*MED_SNAME_SIZE+1] = "x y z ";
00111 char _meshcomponentunit[3*MED_SNAME_SIZE+1] = "cm cm cm ";
00112 char _fieldname [MED_NAME_SIZE+1]="";
00113 char *componentname,*componentunit;
00114 char _profilename[MED_NAME_SIZE+1]=MED_NO_PROFILE;
00115 med_int *_profilearray=0;
00116 int _i=0,_j=0,_k=0, _lastusedrank=0;
00117 med_size _blocksize=0,_lastblocksize=0,_count=0,_stride=0,_start=0,_index=0;
00118 med_float *_arrayvalues;
00119 med_filter filter = MED_FILTER_INIT;
00120 med_size _nusedentities = nentities;
00121 med_size _io_count = nbblocksperproc;
00122 med_idt _fidseq;
00123
00124
00125 char *_MED_MODE_SWITCH_MSG[3]={"MED_FULL_INTERLACE", "MED_NO_INTERLACE","MED_UNDEF_INTERLACE",};
00126 char *_MED_STORAGE_MODE_MSG[3]={"MED_NO_PFLMODE","MED_GLOBAL_PFLMODE", "MED_COMPACT_PFLMODE"};
00127
00128 med_geometry_type _geotype = MED_TRIA6;
00129 med_int _geodim = _geotype/100;
00130 med_int _geonnodes = _geotype%100;
00131 char _ipointname[MED_NAME_SIZE+1];
00132 med_float* _ipointrefcoo = 0;
00133 med_int _ipoint = nvaluesperentity;
00134 med_float* _ipointcoo = 0;
00135 med_float* _ipointwg = 0;
00136
00137 sprintf(_filename,"%s_CPU-%03d_@_%s_%s.med",fieldnameprefix,mpi_size,_MED_MODE_SWITCH_MSG[constituentmode],_MED_STORAGE_MODE_MSG[storagemode]);
00138
00139
00140 if ((_fid = MEDparFileOpen(_filename, MODE_ACCES ,comm, info)) < 0){
00141 MED_ERR_(_ret,MED_ERR_OPEN,MED_ERR_FILE,_filename);
00142 goto ERROR;
00143 }
00144
00145
00146 if (MEDmeshCr( _fid,_meshname,_meshdim,_meshdim, MED_UNSTRUCTURED_MESH,
00147 "Un maillage pour le test parallel","s", MED_SORT_DTIT,
00148 MED_CARTESIAN, _meshcomponentname, _meshcomponentunit) < 0) {
00149 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_MESH,_meshname);
00150 goto ERROR;
00151 };
00152
00153 componentname = (char*) malloc((nconstituentpervalue*MED_SNAME_SIZE+1)*sizeof(char));
00154 componentunit = (char*) malloc((nconstituentpervalue*MED_SNAME_SIZE+1)*sizeof(char));
00155
00156 strcpy(componentname,"");
00157 strcpy(componentunit,"");
00158 strcpy(_fieldname,fieldnameprefix);
00159 if ( MEDfieldCr(_fid,_fieldname,MED_FLOAT64,nconstituentpervalue,componentname,componentunit,"s",_meshname ) < 0) {
00160 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FIELD,_fieldname);
00161 goto ERROR;
00162 };
00163 free(componentname);
00164 free(componentunit);
00165
00166
00167 if ( _ipoint > 1 ) {
00168
00169 MESSAGE("Creating a localization of integration points...");
00170 strcpy(_ipointname,_fieldname);
00171 strcat(_ipointname,"_loc");
00172
00173
00174 _ipointrefcoo = (med_float *) calloc(_geodim*_geonnodes,sizeof(med_float));
00175 _ipointcoo = (med_float *) calloc(_ipoint*_geodim,sizeof(med_float));
00176 _ipointwg = (med_float *) calloc(_ipoint,sizeof(med_float));
00177
00178 if (MEDlocalizationWr(_fid, _ipointname, _geotype, _geotype/100, _ipointrefcoo, constituentmode,
00179 _ipoint, _ipointcoo, _ipointwg, MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT ) < 0) {
00180 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_LOCALIZATION,_ipointname);
00181 ISCRUTE_int(constituentmode);
00182 goto ERROR;
00183 }
00184 free(_ipointrefcoo );
00185 free(_ipointcoo );
00186 free(_ipointwg );
00187
00188 } else {
00189 strcpy(_ipointname,MED_NO_LOCALIZATION);
00190 }
00191
00192 if (profilearraysize) {
00193 MESSAGE("Creating a profile...");
00194
00195 strcpy(_profilename,_fieldname);strcat(_profilename,"_profile");
00196
00197 _profilearray = (med_int*) calloc(profilearraysize,sizeof(med_int));
00198
00199 for (_i=0; _i < profilearraysize; ++_i) _profilearray[_i]=_i;
00200 if ( MEDprofileWr(_fid,_profilename,profilearraysize,_profilearray) < 0) {
00201 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_PROFILE,_profilename);
00202 goto ERROR;
00203 };
00204 _nusedentities = profilearraysize;
00205 } else {
00206 strcpy(_profilename,MED_NO_PROFILE);
00207 }
00208
00209
00210 MESSAGE("Generating partition...");
00211 getBlockOfEntities ( mpi_rank , mpi_size, _nusedentities,
00212 &_start, &_stride, &_io_count, &_blocksize,
00213 &_lastusedrank, &_lastblocksize);
00214
00215 _count=_io_count;
00216 MESSAGE("Generating filter...");
00217 if ( MEDfilterBlockOfEntityCr(_fid, nentities, nvaluesperentity, nconstituentpervalue,
00218 MED_ALL_CONSTITUENT, constituentmode, storagemode, _profilename,
00219 _start,_stride,_count,_blocksize,_lastblocksize, &filter) < 0 ) {
00220 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FILTER,"");
00221 goto ERROR;
00222 }
00223
00224 MESSAGE("Generating datas...");
00225 generateDatas(mpi_rank, _lastusedrank, sizeof(med_float),
00226 storagemode, profilearraysize, _profilearray,
00227 _start, _stride, _count, _blocksize, _lastblocksize,
00228 nentities, nvaluesperentity, nconstituentpervalue,
00229 &_arrayvalues );
00230
00231 MESSAGE("Writing field...");
00232 if ( MEDfieldValueAdvancedWr(_fid,_fieldname,MED_NO_DT,MED_NO_IT,0.0, MED_CELL, _geotype,
00233 _ipointname, &filter, (unsigned char*)_arrayvalues ) < 0) {
00234 MED_ERR_(_ret,MED_ERR_WRITE,MED_ERR_FIELD,_fieldname);
00235 ISCRUTE(mpi_rank);
00236 goto ERROR;
00237 }
00238
00239
00240
00241 H5Fflush(_fid, H5F_SCOPE_GLOBAL );
00242
00243
00244 if (mpi_rank == 0 ) {
00245 MESSAGE("Reading field...");
00246
00247
00248 med_int _nentitiesarrayvalues=0;
00249 med_float *_filteredarrayvalues=NULL;
00250 med_filter filter2=MED_FILTER_INIT;
00251 int _ind=0;
00252 FILE * _asciifile;
00253 char _asciifilename[255]="";
00254
00255
00256 if ((_fidseq = MEDfileOpen(_filename, MED_ACC_RDONLY )) < 0){
00257 MED_ERR_(_ret,MED_ERR_OPEN,MED_ERR_FILE,_filename);
00258 goto ERROR;
00259 }
00260
00261 sprintf(_asciifilename,"%s_CPU-%03d_@_%s_%s.ascii",fieldnameprefix,mpi_size,_MED_MODE_SWITCH_MSG[constituentmode],_MED_STORAGE_MODE_MSG[storagemode]);
00262 _asciifile=fopen(_asciifilename, "w");
00263
00264
00265
00266 if (!_filterarray)
00267 if ( generateFilterArray( nentities, nvaluesperentity, nconstituentpervalue,
00268 profilearraysize, _profilearray,
00269 &_nentitiesfiltered, &_filterarray ) < 0 ) {
00270 goto ERROR;
00271 }
00272
00273 ISCRUTE(_nentitiesfiltered);
00274
00275 for (_i=0; _i < _nentitiesfiltered; ++_i ) {
00276
00277 fprintf(_asciifile,"%d ",_filterarray[_i]) ;
00278 }
00279 fprintf(_asciifile,"\n") ;
00280
00281
00282
00283 if ( profilearraysize ) {
00284 _nentitiesarrayvalues = profilearraysize;
00285 } else {
00286 _nentitiesarrayvalues = nentities;
00287 }
00288
00289
00290
00291
00292
00293 _filteredarrayvalues = (med_float*) malloc(_nentitiesarrayvalues*
00294 nvaluesperentity*
00295 nconstituentpervalue*sizeof(med_float));
00296
00297
00298 for (_i=0;_i<_nentitiesarrayvalues*nvaluesperentity*nconstituentpervalue; ++_i)
00299 _filteredarrayvalues[_i]=-_i;
00300
00301
00302
00303 if ( MEDfilterEntityCr(_fidseq, nentities, nvaluesperentity, nconstituentpervalue,
00304 MED_ALL_CONSTITUENT, constituentmode, storagemode, _profilename,
00305 _nentitiesfiltered,_filterarray, &filter2) < 0 ) {
00306 MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FILTER,"");
00307 goto ERROR;
00308 }
00309
00310 if ( MEDfieldValueAdvancedRd(_fidseq,_fieldname,MED_NO_DT,MED_NO_IT, MED_CELL, _geotype,
00311 &filter2, (unsigned char*)_filteredarrayvalues ) < 0) {
00312 MED_ERR_(_ret,MED_ERR_READ,MED_ERR_FIELD,_fieldname);
00313 ISCRUTE(mpi_rank);
00314 goto ERROR;
00315 }
00316
00317
00318
00319 if ( storagemode == MED_GLOBAL_PFLMODE ) {
00320 switch (constituentmode) {
00321 case MED_FULL_INTERLACE:
00322 for (_i=0; _i < _nentitiesfiltered; ++_i)
00323 for (_j=0; _j < nvaluesperentity; ++_j)
00324 for (_k=0; _k < nconstituentpervalue; ++_k) {
00325 _ind = (_filterarray[_i]-1)*nvaluesperentity*nconstituentpervalue+ _j*nconstituentpervalue+_k;
00326
00327 fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
00328 }
00329 break;
00330 case MED_NO_INTERLACE:
00331 for (_j=0; _j < _nentitiesfiltered; ++_j)
00332 for (_k=0; _k < nvaluesperentity; ++_k)
00333 for (_i=0; _i < nconstituentpervalue; ++_i) {
00334 _ind =_i*nentities*nvaluesperentity+ (_filterarray[_j]-1)*nvaluesperentity +_k;
00335
00336 fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
00337 }
00338 break;
00339 }
00340 } else
00341 switch (constituentmode) {
00342 case MED_FULL_INTERLACE:
00343 for (_i=0; _i < _nentitiesfiltered; ++_i )
00344 for (_j=0; _j < nvaluesperentity; ++_j)
00345 for (_k=0; _k < nconstituentpervalue; ++_k) {
00346 _ind = _i*nvaluesperentity*nconstituentpervalue+_j*nconstituentpervalue+_k;
00347
00348 fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
00349 }
00350 break;
00351 case MED_NO_INTERLACE:
00352 for (_j=0; _j < _nentitiesfiltered; ++_j)
00353 for (_k=0; _k < nvaluesperentity; ++_k)
00354 for (_i=0; _i < nconstituentpervalue; ++_i) {
00355 _ind =_i*_nentitiesfiltered*nvaluesperentity+ _j*nvaluesperentity +_k;
00356
00357
00358 fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
00359 }
00360 break;
00361 }
00362
00363
00364 free(_filteredarrayvalues);
00365
00366 fclose(_asciifile);
00367
00368 if ( MEDfilterClose(&filter2) < 0 ) {
00369 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILTER,"");
00370 goto ERROR;
00371 }
00372
00373 }
00374
00375 if ( MEDfilterClose(&filter) < 0 ) {
00376 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILTER,"");
00377 goto ERROR;
00378 }
00379
00380
00381 _ret=0;
00382 ERROR:
00383 if (_arrayvalues) free(_arrayvalues);
00384 if (profilearraysize) free(_profilearray);
00385
00386 if ( MEDfileClose(_fid) < 0) {
00387 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,""); _ret = -1;
00388 }
00389
00390 if (mpi_rank == 0 ) {
00391 if ( MEDfileClose(_fidseq) < 0) {
00392 MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,""); _ret = -1;
00393 }
00394 }
00395
00396 return _ret;
00397 }
00398
00399
00400
00401
00402
00403 med_size _nbblocksperproc = 0;
00404 int _nentities = 0;
00405 int _nvaluesperentity = 0;
00406 int _nconstituentpervalue = 0;
00407
00408 if (mpi_rank == 0 ) {
00409
00410 struct tm *_tm ;
00411 time_t _tt=time(0);
00412 _tm = localtime(&_tt);
00413
00414 srandom((*_tm).tm_sec * (*_tm).tm_min );
00415 _nbblocksperproc = 1 + (int) (mpi_size * (random() / (RAND_MAX + 1.0)));
00416 _nentities = 1 + (int) (1000.0 * (random() / (RAND_MAX + 1.0)));
00417 _nvaluesperentity = 1 + (int) (11.0 * (random() / (RAND_MAX + 1.0)));
00418 _nconstituentpervalue = 1 + (int) (7.0 * (random() / (RAND_MAX + 1.0)));
00419
00420
00421
00422
00423
00424
00425
00426
00427
00428 }
00429
00430 if ( (sizeof(med_size)%(sizeof(MPI_LONG)))==0 ) {
00431
00432 MPI_Bcast(&_nbblocksperproc , sizeof(med_size)/sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
00433 MPI_Bcast(&_nentities , sizeof(med_size)/sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
00434 MPI_Bcast(&_nvaluesperentity , sizeof(med_size)/sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
00435 MPI_Bcast(&_nconstituentpervalue , sizeof(med_size)/sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
00436 } else {
00437 assert(sizeof(med_size) == (sizeof(MPI_LONG)));
00438 }
00439
00440 char _fieldnameprefix[256] = "";
00441
00442 sprintf(_fieldnameprefix,"NENT-%03d_NVAL-%03d_NCST-%03d_NBL-%03llu",_nentities,_nvaluesperentity,
00443 _nconstituentpervalue,_nbblocksperproc);
00444
00445 GenerateDataType generateDatas = 0;
00446 med_switch_mode _switchmode = MED_UNDEF_INTERLACE;
00447 med_storage_mode _storagemode = MED_UNDEF_PFLMODE;
00448
00449
00450 med_int _profilearraysize=0;
00451
00452
00453 for (_switchmode = MED_FULL_INTERLACE ; _switchmode <= MED_NO_INTERLACE; ++_switchmode) {
00454
00455 if ( _switchmode == MED_FULL_INTERLACE ) generateDatas = generateFullIDatas;
00456 else generateDatas = generateNoIDatas;
00457
00458 for (_storagemode = MED_GLOBAL_PFLMODE ; _storagemode <= MED_COMPACT_PFLMODE; ++_storagemode) {
00459
00460 if ( (_storagemode == MED_GLOBAL_PFLMODE ) && (_profilearraysize) ) _profilearraysize=0;
00461
00462 if ( generateFieldFile( _nentities, _nvaluesperentity, _nconstituentpervalue,
00463 _switchmode, getCyclicBlocksOfEntities, _nbblocksperproc, generateDatas,
00464 _storagemode, _profilearraysize, _fieldnameprefix ) < 0 ) {
00465 MED_ERR_(_ret,MED_ERR_WRITE,MED_ERR_FIELD,_fieldnameprefix);
00466 ISCRUTE(mpi_rank);
00467 goto ERROR;
00468 }
00469
00470 }
00471 }
00472
00473
00474 _ret = 0;
00475 ERROR:
00476
00477 if (_filterarray) free(_filterarray);
00478
00479
00480
00481
00482 MPI_Finalize();
00483
00484
00485 return _ret;
00486 }
00487
00488
00489
00490