Créer et lire des champs avec gestion des profils
Lire des champs créés avec MED 2.x et définis sur plusieurs maillages
Créer et lire un champ aux points d'intégration des éléments
Créer et lire un champ aux noeuds des éléments
Créer et lire un champ lié à un maillage contenant des éléments de structure
Ecrire et lire un maillage ou un champ en filtrant les données
Ecriture et lecture en parallèle dans un seul fichier
Créer et lire des champs avec gestion des profils
Au même titre que pour un maillage évolutif, il est possible d'utiliser les profils pour sélectionner les valeurs des champs de résultats. Les champs ne peuvent être definis que sur une partie du maillage. Le profil d'un champ indique sur quelles entités du maillage se trouvent les valeurs.
A l'écriture des champs, l'utilisation des profils est possible pour l'écriture des valeurs du champ avec les routines MEDfieldValueWithProfileWr / mfdrpw, mfdipw.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const med_float tria3values_step1_profile1[3] = {1000., 4000., 8000.};
const med_float tria3values_step2_profile1[8] = {1500., 0., 0., 4500., 0., 0., 0., 8500.};
const med_float tria3values_step2_profile2[8] = { 0., 2500., 3500., 0., 5500., 6500., 7500., 0.};
const med_float quad4values_step1[4] = {10000., 20000., 30000., 40000.};
const med_float quad4values_step2[4] = {15000., 25000., 35000., 45000.};
const med_int profile1[3] = {1, 4, 8};
const med_int profile2[5] = {2, 3, 5, 6, 7};
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
if (
MEDlinkWr(fid,meshname,
"./UsesCase_MEDmesh_1.med") < 0) {
MESSAGE(
"ERROR : create mesh link ...");
goto ERROR;
}
if (
MEDprofileWr(fid, profile1name, profile1size, profile1 ) < 0) {
MESSAGE(
"ERROR : create profile ...");
goto ERROR;
}
if (
MEDprofileWr(fid, profile2name, profile2size, profile2 ) < 0) {
MESSAGE(
"ERROR : create profile ...");
goto ERROR;
}
componentname, componentunit,"ms", meshname) < 0) {
goto ERROR;
}
(unsigned char*) tria3values_step1_profile1) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3");
goto ERROR;
}
(unsigned char*) quad4values_step1) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ");
goto ERROR;
}
(unsigned char*) tria3values_step2_profile1) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3 ...");
goto ERROR;
}
(unsigned char*) tria3values_step2_profile2) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3 ...");
goto ERROR;
}
(unsigned char*) quad4values_step2) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ... ");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
A la lecture, la routine MEDfieldnProfile / mfdnpf permet de lire le nombre de profil dans un champ pour une séquence de calcul donnée. La routine MEDfieldnValueWithProfile / mfdnvp permet de lire le nombre de valeurs à lire en tenant compte du profil et mode de stockage des données en mémoire (MED_GLOBAL_STMODE ou MED_COMPACT_STMODE). Les routines MEDfieldValueWithProfileRd / mfdrpr, mfdipr permettent de lire les valeurs du champ avec ou sans profil.
#define MESGERR 1
#include <string.h>
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char *componentname = NULL;
char *componentunit = NULL;
med_int nprofile, pit, profilesize;
int ret=-1;
if (fid < 0) {
return -1;
}
MESSAGE(
"ERROR : How many fields in the file ...");
return -1;
}
for (i=0; i<nfield; i++) {
MESSAGE(
"ERROR : number of field component ...");
return -1;
}
if ((componentname = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if ((componentunit = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if (
MEDfieldInfo(fid, i+1, fieldname, meshname, &localmesh, &fieldtype,
componentname, componentunit, dtunit, &nstep) < 0) {
free(componentname);
free(componentunit);
return -1;
}
free(componentname);
free(componentunit);
for (csit=0; csit<nstep; csit++) {
MESSAGE(
"ERROR : Computing step info ...");
return -1;
}
geotype = geotypes[it];
profilename, localizationname)) < 0) {
MESSAGE(
"ERROR : read number of profile ");
return -1;
}
for (pit=0; pit<nprofile; pit++) {
localizationname, &nintegrationpoint)) < 0) {
MESSAGE(
"ERROR : read number of values with a profile ...");
return -1;
}
if (nvalues) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
MESSAGE(
"ERROR : read fields values for cells ...");
free(values);
return -1;
}
free(values);
}
}
}
}
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Lire des champs créés avec MED 2.x et définis sur plusieurs maillages
Dans MED 2.x, les champs pouvaient être définis sur plusieurs maillages. Cette possibilité disparaît avec MED 3.0 et l'apparition des maillages évolutifs. MED 3.0 fournit néanmoins des routines permettant la lecture de fichier MED 2.x pouvant contenir des champs définis sur plusieurs maillages. Ces routines spécifiques sont à utiliser en lieu et place de celles utilisées dans le cas d'utilisation précédent (ces routines peuvent évidemment également lire des champs dans des fichiers MED 3.0), il s'agit de MEDfield23ComputingStepMeshInfo / mfdoci, MEDfield23nValueWithProfile / mfdonv, MEDfield23ValueWithProfileRd / mfdorr, mfdoir.
Le cas d'utilisation suivant propose une approche générique pour lire les champs d'un fichier MED, basés sur un ou plusieurs maillages.
#define MESGERR 1
#include <string.h>
#ifdef DEF_LECT_ECR
#define MODE_ACCES MED_ACC_RDWR
#elif DEF_LECT_AJOUT
#define MODE_ACCES MED_ACC_RDEXT
#else
#define MODE_ACCES MED_ACC_CREAT
#endif
#ifndef USER_INTERLACE
#define USER_INTERLACE MED_FULL_INTERLACE
#endif
#define USER_MODE MED_COMPACT_STMODE
int main (
int argc,
char **argv)
{
char * fichier = NULL;
char * lien = NULL;
char *comp= NULL, *unit= NULL;
med_int mdim=0,sdim=0,ncomp,ncha,npro,nln,pflsize,*pflval,nval;
med_int _ncstp=0,ngauss=0,nloc=0,locsdim=0,lnsize=0;
int t1,t2,t3;
int i,j;
if (argc != 2) {
MESSAGE(
"Aucun nom de fichier precise, fichier test10.med utilise ");
fichier = "test10.med";
} else {
fichier = argv[1];
};
return -1;
}
ret = 0;
if (
MEDmeshInfo( fid, 1, maa, &sdim, &mdim, &type, desc, dtunit, &sort,
&nstep, &rep, nomcoo,unicoo) < 0 ) {
MESSAGE(
"Erreur a la lecture des informations sur le maillage : ");
SSCRUTE(maa);
return -1;
} else {
printf(
"Maillage de nom : |%s| , de dimension : "IFORMAT" , et de type %d\n",maa,mdim,type);
printf(
"\t -Dimension de l'espace : "IFORMAT"\n",sdim);
printf("\t -Description du maillage : %s\n",desc);
printf("\t -Noms des axes : |%s|\n",nomcoo);
printf("\t -Unités des axes : |%s|\n",unicoo);
printf("\t -Type de repère : %d\n",rep);
printf(
"\t -Nombre d'étapes de calcul : "IFORMAT"\n",nstep);
printf("\t -Unité des dates : |%s|\n",dtunit);
}
return ncha;
}
printf(
"Nombre de champs : "IFORMAT" \n",ncha);
for (i =0;i<ncha;i++) {
lret = 0;
printf("\nChamp numero : %d \n",i+1);
MESSAGE(
"Erreur a la lecture du nombre de composantes : ");
ISCRUTE(ncomp);
ret = -1; continue;
}
if (
MEDfieldInfo(fid,i+1,nomcha,_meshname,&_local,&typcha,comp,unit,_dtunit,&_ncstp) < 0 ) {
MESSAGE(
"Erreur a la demande d'information sur les champs : ");
ret = -1; continue;
}
printf("Nom du champ : |%s| de type %d\n",nomcha,typcha);
printf("Nom des composantes : |%s|\n",comp);
printf("Unites des composantes : |%s| \n",unit);
printf("Unites des dates : |%s| \n",_dtunit);
printf("Le maillage associé est |%s|\n",_meshname);
printf(
"Nombre d'étapes de calcul |"IFORMAT"|\n",_ncstp);
if ( !_local ) {
MESSAGE(
"Erreur a la lecture de la taille du lien : ");
ret = -1;
} else {
lien = malloc((lnsize+1)*sizeof(char));
MESSAGE(
"Erreur a la lecture du lien : ");
ret = -1;
} else {
printf("\tLe maillage |%s| est porte par un fichier distant |%s|\n",_meshname,lien);
}
free(lien);
}
}
free(comp);
free(unit);
else {
MESSAGE(
"Erreur a la lecture des champs aux noeuds "); ret = -1;
continue;}
else {
MESSAGE(
"Erreur a la lecture des champs aux mailles "); ret = -1;
continue;}
else {
MESSAGE(
"Erreur a la lecture des champs aux faces "); ret = -1;
continue;}
else {
MESSAGE(
"Erreur a la lecture des champs aux aretes"); ret = -1;
continue;}
if (lret != 0) {
MESSAGE(
"Erreur a la lecture des champs aux noeuds des mailles "); ret = -1;};
}
printf(
"\nNombre de profils stockes : "IFORMAT"\n\n",npro);
for (i=1 ; i <= npro ; i++ ) {
ret = -1;continue;
}
printf(
"\t- Profil n°%i de nom |%s| et de taille "IFORMAT"\n",i,pflname,nval);
MESSAGE(
"Erreur a la lecture des valeurs du profil : ");
ret = -1;
} else {
printf("\t");
for (j=0;j<nval;j++) printf(
" "IFORMAT" ",*(pflval+j));
printf("\n\n");
}
free(pflval);
}
printf(
"\nNombre de liens stockes : "IFORMAT"\n\n",nln);
for (i=1 ; i <= nln ; i++ ) {
ret = -1;continue;
}
printf(
"\t- Lien n°%i de nom |%s| et de taille "IFORMAT"\n",i,nomlien,nval);
lien = (char * ) malloc((nval+1)*sizeof(char));
MESSAGE(
"Erreur a la lecture du lien : ");
ret = -1;
} else {
lien[nval] = '\0';
printf("\t\t|%s|\n\n",lien);
}
free(lien);
}
printf(
"\nNombre de localisations stockees : "IFORMAT"\n\n",nloc);
for (i=1 ; i <= nloc ; i++ ) {
geointerpname, ipointstructmeshname,&nsectionmeshcell,
§iongeotype) < 0) {
ret = -1;continue;
}
printf(
"\t- Loc. n°%i de nom |%s| de dimension "IFORMAT" avec "IFORMAT" pts de GAUSS \n",i,locname,locsdim,ngauss);
t1 = (type_geo%100)*(type_geo/100);
t2 = ngauss*(type_geo/100);
t3 = ngauss;
MESSAGE(
"Erreur a la lecture des valeurs de la localisation : ");
ret = -1;
} else {
printf("\t Coordonnees de l'element de reference de type %i :\n\t\t",type_geo);
for (j=0;j<t1;j++) printf(" %f ",*(refcoo+j));
printf("\n");
printf("\t Localisation des points de GAUSS : \n\t\t");
for (j=0;j<t2;j++) printf(" %f ",*(gscoo+j));
printf("\n");
printf("\t Poids associes aux points de GAUSS :\n\t\t");
for (j=0;j<t3;j++) printf(" %f ",*(wg+j));
printf("\n\n");
}
free(refcoo);
free(gscoo);
free(wg);
}
return ret;
}
int i,j,k,l,m,n,nb_geo=0;
med_int nbpdtnor=0,pflsize,*pflval,ngauss=0,ngroup,*vale=NULL,nval;
char * lien = NULL;
const char * const * AFF;
switch (entite) {
break;
break;
break;
break;
}
for (k=1;k<=nb_geo;k++) {
nbpdtnor = ncstp;
if (nbpdtnor < 1 ) continue;
for (j=0;j<nbpdtnor;j++) {
&nmesh, meshname,&localmesh, &meshnumdt, &meshnumit ) <0) {
MESSAGE(
"Erreur a la demande d'information sur (pdt,nor) : ");
ret = -1; continue;
}
for (i=0;i< nmesh;++i) {
if ( (_nprofile =
MEDfield23nProfile(fid,nomcha,numdt,numo,entite,type_geo[k],i+1,meshname,
pflname,locname ) ) < 0 ) {
MESSAGE(
"Erreur a la demande du nombre de profils referencés par le champ : ");
ret = -1; continue;
};
for (l=0;l<_nprofile;l++) {
locname, &ngauss) ) < 0 ) {
MESSAGE(
"Erreur a la lecture du nombre de valeurs du champ : ");
ret = -1; continue;
};
printf(
"\n +Pas de Temps n."IFORMAT" (%f) [%s], n. d'ordre "IFORMAT", avec "IFORMAT" valeur(s) par entité.\n",numdt,dt,dt_unit,numo,ngauss);
printf(
"\t- Il y a "IFORMAT" entités qui portent des valeurs en mode %i. Chaque entite %s\
de type geometrique %s associes au profile |%s| a "IFORMAT" valeurs associées \n",
nval,
USER_MODE,AFF_ENT[(
int)entite],AFF[k],pflname,ngauss);
printf("\t- Le maillage associé est |%s|\n",meshname);
(unsigned char*) valr) < 0 ) {
MESSAGE(
"Erreur a la lecture des valeurs du champ : ");
ret = -1;
}
} else {
(unsigned char*) vale) < 0 ) {
MESSAGE(
"Erreur a la lecture des valeurs du champ : ");
ret = -1;
};
}
if ( strlen(locname) )
printf("\t- Modèle de localisation des points de Gauss de nom |%s|\n",locname);
ngroup = (type_geo[k] % 100);
else
ngroup = ngauss;
switch (stockage) {
printf("\t- Valeurs :\n\t");
for (m=0;m<(nval*ngauss)/ngroup;m++) {
printf("|");
for (n=0;n<ngroup*ncomp;n++)
printf(" %f ",*(valr+(m*ngroup*ncomp)+n));
else
printf(
" "IFORMAT" ",*(vale+(m*ngroup*ncomp)+n));
}
break;
printf("\t- Valeurs :\n\t");
for (m=0;m<ncomp;m++) {
printf("|");
for (n=0;n<(nval*ngauss);n++)
printf(" %f ",*(valr+(m*nval)+n));
else
printf(
" "IFORMAT" ",*(vale+(m*nval)+n));
}
break;
}
printf("|\n");
if ( valr ) {free(valr);valr = NULL;}}
else
if (vale) { free(vale);vale = NULL; }
printf("\t- Profil : MED_NO_PROFILE\n");
else {
MESSAGE(
"Erreur a la lecture du nombre de valeurs du profil : ");
ret = -1; continue;
}
printf(
"\t- Profil : |%s| de taille "IFORMAT"\n",pflname,pflsize);
MESSAGE(
"Erreur a la lecture des valeurs du profil : ");
ret = -1;
}
printf("\t");
for (m=0;m<pflsize;m++) printf(
" "IFORMAT" ",*(pflval+m));
printf("\n");
free(pflval);
}
}
}
}
}
return ret;
}
Créer et lire un champ aux points d'intégration des éléments
MED fournit la possibilité d'exprimer les champs de résultat sur les points de Gauss (ou points d'intégration) des éléments d'un maillage. Dans ce cadre, il est possible de localiser ces points sur des éléments de référence en des lieux différents selon la modélisation numérique choisie. Pour chaque type de modélisation, il est possible de spécifier cette localisation sur des éléments de référence. On distingue différentes familles de points de Gauss en fonction du nombre de points d'intégration. Chaque point d'intégration est localisé au sein d'un élément de référence par ses coordonnées et se voit associer un poids.
La localisation des points de Gauss pour un élément de référence nécessite donc de connaître le type géométrique de l'élément, les coordonnées des noeuds de l'élément, les coordonnées et le poids de chaque point de Gauss. L'expression des coordonnées d'un élément de référence peut se faire dans un repère de coordonnées dont la dimension est supérieure à celle de l'élément de référence. La référence à un élément de référence se fait à l'appel des routines d'écriture et lecture des valeurs des champs.
Si les points de Gauss se confondent avec les noeuds de l'élément, il est inutile de créer une localisation factice avec des poids qui ne signifient rien et des coordonnées des points de Gauss identiques à celles des noeuds. Dans ce cas de figure, il faut utiliser le mot clé réservé MED_GAUSS_ELNO à l'écriture des valeurs d'un champ pour indiquer le type de localisation.
Il est possible pour une même type géométrique d'élément d'associer plusieurs éléments de référence. Il suffit pour cela d'associer chaque localisation avec un profil lors de l'appel d'écriture et lecture des valeurs des champs.
MED permet également la prise en compte des fonctions de forme et des familles d'interpolation. Une famille d'interpolation est l'ensemble d'interpolations disponibles pour un champ donné. Chacune des interpolations est décrite par des fonctions de forme polynomiales. Le nombre de variable des fonctions de forme est égale à la dimension de l'espace de la maille de référence utilisée pour le construire. Il est donc possible d'associer à un champ plusieurs interpolations définies sur des mailles de référence des différents types géométriques du maillage de calcul sur lequel repose le champ de résultat.
L'écriture des fonctions de forme et des familles d'interpolation est optionnelle pour échanger des champs aux points d'intégration.
A l'écriture, la création d'un élément de référence se fait avec la routine MEDlocalizationWr / mlclow. Un élément de référence peut être associé à une famille d'interpolation. La création d'une fonction d'interpolation se fait avec la routine MEDinterpBaseFunctionWr / mipcre. L'écriture d'une fonction de forme se fait avec la routine MEDinterpBaseFunctionWr / mipbfw. La routine MEDfieldInterpWr / mfdinw permet d'associer une famille d'interpolation à un champ résultat.
Le cas d'utilisation suivant montre un exemple d'écriture des valeurs des champs aux points de Gauss avec définition des éléments de référence et référence à une famille d'interpolation.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const med_float tria3values_step1_profile1[9] = {1000.,1010.,1020.,
4000.,4010.,4020.,
8000.,8010.,8020. };
const med_float tria3values_step2_profile1[24] = {1500.,1510.,1520.,
0., 0., 0.,
0., 0., 0.,
4500.,4510,4520.,
0., 0., 0.,
0., 0., 0.,
0., 0., 0.,
8500., 8510, 8520 };
const med_float tria3values_step2_profile2[32] = { 0., 0., 0., 0.,
2500.,2510.,2520,2530.,
3500.,3510.,3520.,3530.,
0., 0., 0., 0.,
5500.,5510.,5520.,5530.,
6500.,6510.,6520.,6530.,
7500.,7510.,7520.,7530.,
0., 0., 0., 0. };
const med_float quad4values_step1[4] = {10000., 20000., 30000., 40000.};
const med_float quad4values_step2[4] = {15000., 25000., 35000., 45000.};
const med_int profile1[3] = {1, 4, 8};
const med_int profile2[5] = {2, 3, 5, 6, 7};
const char localization1name[
MED_NAME_SIZE+1] =
"TRIA3_INTEGRATION_POINTS_3";
const med_float weight1[3] = {1.0/6, 1.0/6, 1.0/6};
const med_float elementcoordinate[6] = {0.0, 0.0, 1.0, 0.0, 0.0,1.0};
const med_float ipoint1coordinate[6] = {1.0/6, 1.0/6, 2.0/3, 1.0/6, 1.0/6, 2.0/6};
const char localization2name[
MED_NAME_SIZE+1] =
"TRIA3_INTEGRATION_POINTS_4";
const med_float weight2[6] = {25.0/(24*4), 25.0/(24*4), 25.0/(24*4), -27.0/(24*4)};
const med_float ipoint2coordinate[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};
const char interpname[
MED_NAME_SIZE+1] =
"MED_TRIA3 interpolation family";
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
if (
MEDlinkWr(fid,meshname,
"./UsesCase_MEDmesh_1.med") < 0) {
MESSAGE(
"ERROR : create mesh link ...");
goto ERROR;
}
if (
MEDprofileWr(fid, profile1name, profile1size, profile1 ) < 0) {
MESSAGE(
"ERROR : create profile ...");
goto ERROR;
}
if (
MEDprofileWr(fid, profile2name, profile2size, profile2 ) < 0) {
MESSAGE(
"ERROR : create profile ...");
goto ERROR;
}
spacedim = 2;
nipoint = 3;
nipoint, ipoint1coordinate, weight1,
MESSAGE(
"ERROR : create famlily of integration points ...");
goto ERROR;
}
spacedim = 2;
nipoint = 4;
nipoint, ipoint2coordinate, weight2,
MESSAGE(
"ERROR : create famlily of integration points ...");
goto ERROR;
}
ncomponent, componentname, componentunit,
"ms", meshname) < 0) {
goto ERROR;
}
MESSAGE(
"ERROR : write field interpolation family name ...");
goto ERROR;
}
ntria3, (unsigned char*) tria3values_step1_profile1) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3");
goto ERROR;
}
nquad4, (unsigned char*) quad4values_step1) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ");
goto ERROR;
}
ntria3, (unsigned char*) tria3values_step2_profile1) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3 ...");
goto ERROR;
}
ntria3, (unsigned char*) tria3values_step2_profile2) < 0) {
MESSAGE(
"ERROR : write field values on MED_TRIA3 ...");
goto ERROR;
}
nquad4, (unsigned char*) quad4values_step2) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ... ");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
A la lecture, la routine MEDlocalizationInfoByName / mlclni permet de lire les informations relatives à un élément de référence dont on connaît le nom. Une autre possibilité est de lire le nombre d'élément de référence avec la routine MEDnLocalization / mlcnlc et d'itérer afin de récupérer à chaque itération les informations relatives à l'élément de référence avec la routine MEDlocalizationInfo / mlclci et lire l'élément de référence avec la routine MEDlocalizationRd / mlclor.
Pour les fonctions d'interpolation, la routine MEDinterpInfoByName / mipiin informe des caractéristiques de la fonction d'interpolation dont on connaît le nom. Il est également possible de lire le nombre de famille d'interpolation avec la routine MEDnInterp / mipnip et d'itérer sur ces familles. Ces deux fonctions renvoie le nombre de fonctions de forme. Il reste donc à itérer sur chacune d'entre elles et d'appeler la routine MEDinterpBaseFunctionRd / mipbfr pour lire chaque polynôme.
Dans un champ de résultat, il est possible de lire le nombre de famille d'interpolation associé aux champ avec la routine MEDfieldnInterp / mfdnin . En itérant sur toutes ces familles, on peut lire le nom de chacune d'entre elles avec la routine MEDfieldInterpInfo / mfdini.
Le cas d'utilisation suivant propose une approche générique pour lire des champs aux points de Gauss définies sur les mailles d'un maillage.
#define MESGERR 1
#include <string.h>
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char *componentname = NULL;
char *componentunit = NULL;
med_int nprofile, pit, profilesize;
int k;
int ret=-1;
if (fid < 0) {
return -1;
}
MESSAGE(
"ERROR : How many fields in the file ...");
return -1;
}
for (i=0; i<nfield; i++) {
MESSAGE(
"ERROR : number of field component ...");
return -1;
}
if ((componentname = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if ((componentunit = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if (
MEDfieldInfo(fid, i+1, fieldname, meshname, &localmesh, &fieldtype,
componentname, componentunit, dtunit, &nstep) < 0) {
free(componentname);
free(componentunit);
return -1;
}
free(componentname);
free(componentunit);
MESSAGE(
"ERROR : Read how many interpolation functions for the field ...");
return -1;
}
for (it=0; it<ninterp; it++) {
MESSAGE(
"ERROR : read interpolation family name ...");
return -1;
}
}
for (csit=0; csit<nstep; csit++) {
MESSAGE(
"ERROR : Computing step info ...");
return -1;
}
geotype = geotypes[it];
profilename, localizationname)) < 0) {
MESSAGE(
"ERROR : read number of profile ");
return -1;
}
for (pit=0; pit<nprofile; pit++) {
localizationname, &nintegrationpoint)) < 0) {
MESSAGE(
"ERROR : read number of values with a profile ...");
return -1;
}
if (nvalues) {
if ((values = (
med_float *) malloc(
sizeof(
med_float)*nvalues*ncomponent*nintegrationpoint)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
(unsigned char*) values) < 0) {
MESSAGE(
"ERROR : read fields values for cells ...");
free(values);
return -1;
}
free(values);
}
}
}
}
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Le cas d'utilisation suivant montre un exemple de création d'une famille d'interpolation.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char interpname[
MED_NAME_SIZE+1] =
"MED_TRIA3 interpolation family";
const med_int const power1_1[] = {0,0,1,0,0,1};
const med_float const coefficient1_1[] = {1,-1,-1};
const med_int const power1_2[] = {1,0};
const med_float const coefficient1_2[] = {1};
const med_int const power1_3[] = {0,1};
const med_float const coefficient1_3[] = {1};
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : interpolation family creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : first base function creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : second base function creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : third base function creation ...");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Le cas d'utilisation suivant montre un exemple de lecture de famille d'interpolation dont on connaît le nom.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const char interpname[
MED_NAME_SIZE+1] =
"MED_TRIA3 interpolation family";
int basisfuncit =0;
int powerit =0;
int coefficientit =0;
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
&nvariable,&maxdegree,&nmaxcoefficient) < 0) {
MESSAGE(
"ERROR : interpolation function information ...");
goto ERROR;
}
for ( basisfuncit=1; basisfuncit<= nbasisfunc; ++basisfuncit) {
MESSAGE(
"ERROR : read number of coefficient in the base function ...");
goto ERROR;
}
MESSAGE(
"ERROR : read base function ...");
free(coefficient); free(power);
goto ERROR;
}
free(coefficient);
free(power);
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Le cas d'utilisation suivant montre un exemple de lecture de famille d'interpolation par une approche itérative.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
int it =0;
int basisfuncit =0;
int powerit =0;
int coefficientit =0;
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : read number of interpolation ...");
goto ERROR;
}
for (it=1; it<= ninterp; it++) {
if (
MEDinterpInfo(fid,it,interpname,&geotype,&cellnodes,&nbasisfunc,
&nvariable,&maxdegree,&nmaxcoefficient) < 0) {
MESSAGE(
"ERROR : interpolation function information ...");
goto ERROR;
}
for ( basisfuncit=1; basisfuncit<= nbasisfunc; ++basisfuncit) {
MESSAGE(
"ERROR : read number of coefficient in the base function ...");
goto ERROR;
}
MESSAGE(
"ERROR : read base function ...");
free(coefficient); free(power);
goto ERROR;
}
free(coefficient);
free(power);
}
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Créer et lire un champ aux noeuds des éléments
MED permet d'écrire et lire des champs aux noeuds des éléments d'un maillage. Pour cela, il suffit d'indiquer le mot clé MED_NODE_ELEMENT comme type d'entité lors de l'appel des routines d'écriture / lecture des valeurs des champs comme le montrent les deux cas d'utilisation suivants.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const med_float quad4values_step1[4*4] = { 10000., 20000., 30000., 40000.,
50000., 60000., 70000., 80000.,
90000., 100000., 110000., 120000.,
130000., 140000., 150000., 160000. };
const med_float quad4values_step2[4*4] = { 100., 200., 300., 400.,
500., 600., 700., 800.,
900., 1000., 1100., 1200.,
1300., 1400., 1500., 1600. };
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
if (
MEDlinkWr(fid,meshname,
"./UsesCase_MEDmesh_1.med") < 0) {
MESSAGE(
"ERROR : create mesh link ...");
goto ERROR;
}
ncomponent, componentname, componentunit,
"ms", meshname) < 0) {
goto ERROR;
}
nquad4, (unsigned char*) quad4values_step1) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ");
goto ERROR;
}
nquad4, (unsigned char*) quad4values_step2) < 0) {
MESSAGE(
"ERROR : write field values on MED_QUAD4 ... ");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
#define MESGERR 1
#include <string.h>
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char *componentname = NULL;
char *componentunit = NULL;
med_int nprofile, pit, profilesize;
int k;
int ret=-1;
if (fid < 0) {
goto ERROR;
}
MESSAGE(
"ERROR : How many fields in the file ...");
goto ERROR;
}
for (i=0; i<nfield; i++) {
MESSAGE(
"ERROR : number of field component ...");
goto ERROR;
}
if ((componentname = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
if ((componentunit = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
if (
MEDfieldInfo(fid, i+1, fieldname, meshname, &localmesh, &fieldtype,
componentname, componentunit, dtunit, &nstep) < 0) {
free(componentname);
free(componentunit);
goto ERROR;
}
free(componentname);
free(componentunit);
for (csit=0; csit<nstep; csit++) {
MESSAGE(
"ERROR : Computing step info ...");
goto ERROR;
}
geotype = geotypes[it];
profilename, localizationname)) < 0) {
MESSAGE(
"ERROR : read number of profile ");
goto ERROR;
}
for (pit=0; pit<nprofile; pit++) {
localizationname, &nintegrationpoint)) < 0) {
MESSAGE(
"ERROR : read number of values with a profile ...");
goto ERROR;
}
if (nvalues) {
if ((values = (
med_float *) malloc(
sizeof(
med_float)*nvalues*ncomponent*nintegrationpoint)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
(unsigned char*) values) < 0) {
MESSAGE(
"ERROR : read fields values for node elements ...");
free(values);
goto ERROR;
}
free(values);
}
}
}
}
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
Créer et lire un champ lié à un maillage contenant des éléments de structure
Au même titre que pour types d'éléments pré-définis dans le modèle MED, MED permet l'écriture et la lecture de champs sur les éléments de structure. Pour cela, il suffit d'indiquer le mot clé MED_STRUCT_ELEMENT comme type d'entité lors de l'appel des routines d'écriture / lecture des valeurs des champs comme le montrent les deux cas d'utilisation suivants.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const med_float part_speed1[3*3] = { 1.1, 2.2, 3.3,
4.4, 5.5, 6.6,
7.7, 8.8, 9.9 };
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
if (
MEDlinkWr(fid,meshname,
"./UsesCase_MEDstructElement_1.med") < 0) {
MESSAGE(
"ERROR : create mesh link ...");
goto ERROR;
}
MESSAGE(
"ERROR : file mounting ...");
goto ERROR;
}
ncomponent, componentname, componentunit,
"ms", meshname) < 0) {
goto ERROR;
}
npart, (unsigned char*) part_speed1) < 0) {
MESSAGE(
"ERROR : write field values on MED_PARTICLE ");
goto ERROR;
}
MESSAGE(
"ERROR : file unmounting ...");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
#define MESGERR 1
#include <string.h>
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char *componentname = NULL;
char *componentunit = NULL;
med_int nprofile, pit, profilesize;
int k;
int ret=-1;
if (fid < 0) {
goto ERROR;
}
MESSAGE(
"ERROR : file mounting ...");
goto ERROR;
}
MESSAGE(
"ERROR : file mounting ...");
goto ERROR;
}
MESSAGE(
"ERROR : read number of struct element models ...");
goto ERROR;
}
for (it=0; it<nmodels; it++) {
supportmeshname, &entitype, &nnode, &ncell,
&geocelltype, &nconstatt, &anyprofile, nvaratt+it) < 0) {
MESSAGE(
"ERROR : struct element models information ...");
free(nvaratt);
goto ERROR;
}
}
free(nvaratt);
MESSAGE(
"ERROR : How many fields in the file ...");
goto ERROR;
}
for (i=0; i<nfield; i++) {
MESSAGE(
"ERROR : number of field component ...");
goto ERROR;
}
if ((componentname = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
if ((componentunit = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
if (
MEDfieldInfo(fid, i+1, fieldname, meshname, &localmesh, &fieldtype,
componentname, componentunit, dtunit, &nstep) < 0) {
free(componentname);
free(componentunit);
goto ERROR;
}
free(componentname);
free(componentunit);
for (csit=0; csit<nstep; csit++) {
MESSAGE(
"ERROR : Computing step info ...");
goto ERROR;
}
for (it=0; it<nmodels; it++) {
geotype = *(geotypes+it);
profilename, localizationname)) < 0) {
MESSAGE(
"ERROR : read number of profile ");
goto ERROR;
}
for (pit=0; pit<nprofile; pit++) {
localizationname, &nintegrationpoint)) < 0) {
MESSAGE(
"ERROR : read number of values with a profile ...");
goto ERROR;
}
if (nvalues) {
if ((values = (
med_float *) malloc(
sizeof(
med_float)*nvalues*ncomponent*nintegrationpoint)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
goto ERROR;
}
(unsigned char*) values) < 0) {
MESSAGE(
"ERROR : read fields values for cells ...");
goto ERROR;
free(values);
}
free(values);
}
}
}
}
}
ret=0;
ERROR :
free(geotypes);
MESSAGE(
"ERROR : file unmounting ...");
ret=-1;
}
MESSAGE(
"ERROR : file unmounting ...");
ret= -1;
}
ret= -1;
}
return ret;
}
Il est également possible de définir des champs aux points d'intégration des éléments de structure. Dans ce cas là, l'élément de référence est déjà décrit à la définition de l'élément de structure. Les coordonnées des points d'intégration y sont relatives. Il est possible d'indiquer l'utilisation d'un maillage support définissant une section du modèle d'élément de structure. Ce maillage support est alors utilisé comme section de l'élément de structure à chaque point d'intégration. Auquel cas un champ utilisant cette localisation définira autant de valeur par élément qu'il y a de maille dans le maillage section de chaque point d'intégration comme montrent les deux cas d'utilisation suivants.
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
const med_float tempvalue[3*1*4] = { 1.1, 2.2, 3.3, 4.4,
5.5, 6.6, 7.7, 8.8,
9.9, 10.1,11.11, 12.12};
const char localization[
MED_NAME_SIZE+1] =
"BEAM_INTEGRATION_POINTS";
const char localization2[
MED_NAME_SIZE+1] =
"BEAM_INTEGRATION_TRANSF";
const med_float elementcoordinate[3*3] = { 0.0,0.0,0.0,
0.0,0.0,0.0,
0.0,0.0,0.0,};
const med_float ipointcoordinate[3*3] = { 0.0,0.0,2.5,
0.0,0.0,3.5,
0.0,0.0,4.5};
const med_float weight[4] = {1.0/4, 1.0/4, 1.0/4, 1.0/4};
const med_int const power1_1[] = {0,0,1,0,0,1};
const med_float const coefficient1_1[] = {1,-1,-1};
const med_int const power1_2[] = {1,0};
const med_float const coefficient1_2[] = {1};
const med_int const power1_3[] = {0,1};
const med_float const coefficient1_3[] = {1};
int ret=-1;
if (fid < 0) {
MESSAGE(
"ERROR : file creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : file mounting ...");
goto ERROR;
}
MESSAGE(
"ERROR : file mounting ...");
goto ERROR;
}
if (
MEDlinkWr(fid,meshname,
"./UsesCase_MEDstructElement_1.med") < 0) {
MESSAGE(
"ERROR : create mesh link ...");
goto ERROR;
}
MESSAGE(
"ERROR : interpolation family creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : first base function creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : second base function creation ...");
goto ERROR;
}
MESSAGE(
"ERROR : third base function creation ...");
goto ERROR;
}
nipoint, ipointcoordinate, weight,
MESSAGE(
"ERROR : create famlily of integration points ...");
goto ERROR;
}
nipoint, ipointcoordinate, weight,
interpname, beamsectionname) < 0) {
MESSAGE(
"ERROR : create famlily of integration points ...");
goto ERROR;
}
ncomponent, componentname, componentunit,
"ms", meshname) < 0) {
goto ERROR;
}
nbeam, (unsigned char*) tempvalue) < 0) {
MESSAGE(
"ERROR : write field values on MED_BEAM ");
goto ERROR;
}
nbeam, (unsigned char*) tempvalue) < 0) {
MESSAGE(
"ERROR : write field values on MED_BEAM ");
goto ERROR;
}
MESSAGE(
"ERROR : file unmounting ...");
goto ERROR;
}
MESSAGE(
"ERROR : file unmounting ...");
goto ERROR;
}
ret=0;
ERROR:
ret=-1;
}
return ret;
}
#define MESGERR 1
#include <string.h>
#define MESGERR 1
#include <string.h>
int main (
int argc,
char **argv) {
char *componentname = NULL;
char *componentunit = NULL;
med_int nprofile, pit, profilesize;
int k;
int ret=-1;
if (fid < 0) {
return -1;
}
MESSAGE(
"ERROR : file mounting ...");
return -1;
}
MESSAGE(
"ERROR : file mounting ...");
return -1;
}
MESSAGE(
"ERROR : file mounting ...");
return -1;
}
MESSAGE(
"ERROR : How many fields in the file ...");
return -1;
}
for (i=0; i<nfield; i++) {
MESSAGE(
"ERROR : number of field component ...");
return -1;
}
if ((componentname = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if ((componentunit = (
char *) malloc(ncomponent*
MED_SNAME_SIZE+1)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
if (
MEDfieldInfo(fid, i+1, fieldname, meshname, &localmesh, &fieldtype,
componentname, componentunit, dtunit, &nstep) < 0) {
free(componentname);
free(componentunit);
return -1;
}
free(componentname);
free(componentunit);
&geotransformation)) < 0) {
MESSAGE(
"ERROR : number of nodes ...");
return -1;
}
for (it=0; it<nmodels; it++) {
supportmeshname, &entitype, &nnode, &ncell,
&geocelltype, &nconstatt, &anyprofile, nvaratt+it) < 0) {
MESSAGE(
"ERROR : struct element models information ...");
return -1;
}
}
for (csit=0; csit<nstep; csit++) {
MESSAGE(
"ERROR : Computing step info ...");
return -1;
}
for (it=0; it<nmodels; it++) {
geotype = *(geotypes+it);
profilename, localizationname)) < 0) {
MESSAGE(
"ERROR : read number of profile ");
return -1;
}
for (pit=0; pit<nprofile; pit++) {
localizationname, &nintegrationpoint)) < 0) {
MESSAGE(
"ERROR : read number of values with a profile ...");
return -1;
}
if (nvalues) {
if ((values = (
med_float *) malloc(
sizeof(
med_float)*nvalues*ncomponent*nintegrationpoint)) == NULL) {
MESSAGE(
"ERROR : memory allocation ...");
return -1;
}
(unsigned char*) values) < 0) {
MESSAGE(
"ERROR : read fields values for cells ...");
free(values);
return -1;
}
free(values);
}
}
}
}
ret=0;
ERROR:
if (nvaratt)
free(nvaratt);
if (geotypes)
free(geotypes);
}
MESSAGE(
"ERROR : file unmounting ...");
ret=-1;
}
MESSAGE(
"ERROR : file unmounting ...");
ret=-1;
}
MESSAGE(
"ERROR : file unmounting ...");
ret=-1;
}
ret=-1;
}
return ret;
}
Ecrire et lire un maillage ou un champ en filtrant les données
Dans MED, un filtre de données peut être utilisé pour contrôler un sous-ensemble d'entités concernées par un appel à l'API.
En mode séquentiel, la routine MEDfilterEntityCr / mfrcre permet de créer un filtre élémentaire en sélectionnant les entités pour lesquelles on veut lire/écrire des valeurs. Cette sélection permet une lecture/écriture avancée vers/depuis les emplacements mémoire sélectionnés. Elle s'utilise uniquement en mode séquentiel (un seul processus).
Une fois créé un filtre peut être passé en arguments aux routines avancées de l'API. Ces routines ne manipuleront que les données filtrées. Ces routines avancées sont celles qui permettent d'écrire et lire les coordonnées des noeuds, la connectivité des éléments d'un maillage non structuré, les valeurs des champs de résultats : MEDfieldValueAdvancedWr / mfdraw, mfdiaw, MEDfieldValueAdvancedRd / mfdrar, mfdiar, MEDmeshElementConnectivityAdvancedWr / mmhyaw, MEDmeshElementConnectivityAdvancedRd / mmhyar, MEDmeshNodeCoordinateAdvancedWr / mmhcaw, MEDmeshNodeCoordinateAdvancedRd / mmhcar.
Il est possible d'allouer un tableau de filtre avec la routine MEDfilterAllocate / mfrall, puis de le dé-allouer avec la routine MEDfilterDeAllocate / mfrdea.
Ecriture et lecture en parallèle dans un seul fichier
En mode parallèle, le filtre est un moyen de définir le domaine concerné par chacun des processeurs. Les écritures et lectures peuvent donc se faire en parallèle entre plusieurs processeurs dans un seul et même fichier MED.
La routine MEDfilterBlockOfEntityCr / mfrblc permet de créer un filtre en sélectionnant les entités par blocs continus de taille constante pour lesquelles on veut lire/écrire des valeurs. Cette sélection permet une lecture/écriture avancée vers/depuis les emplacements mémoire sélectionnés. Elle s'utilise aussi bien en mode séquentiel qu'en mode parallèle (un ou plusieurs processus).
Dans ce cadre d'utilisation, le fichier doit être ouvert avec la routine MEDparFileOpen / mfipfo.
L'exemple suivant montre un cas d'utilisation en parallèle en écriture et lecture de champs de résultats.
#define MAX(a,b) ((a) > (b) ? (a) : (b))
#define MESGERR 1
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#ifdef PPRO_NT_CALL
# include <windows.h>
# include <Lmcons.h>
# include <sys/timeb.h>
# include <time.h>
#else
# if HAVE_SYS_TIME_H
# include <sys/time.h>
# endif
#include <time.h>
#ifndef HAVE_UNISTD_H
#error "unistd.h required."
#endif
#endif
#ifdef DEF_LECT_ECR
#define MODE_ACCES MED_ACC_RDWR
#elif DEF_LECT_AJOUT
#define MODE_ACCES MED_ACC_RDEXT
#else
#define MODE_ACCES MED_ACC_CREAT
#endif
typedef struct {
MPI_Info info;
MPI_Comm comm;
int mpi_size;
int mpi_rank;
char _filename [255]="";
char *componentname,*componentunit;
int _i=0,_j=0,_k=0, _lastusedrank=0;
med_size _blocksize=0,_lastblocksize=0,_count=0,_stride=0,_start=0,_index=0;
MPI_Info info = cominfo->
info;
MPI_Comm comm = cominfo->
comm;
char *_MED_MODE_SWITCH_MSG[3]={"MED_FULL_INTERLACE", "MED_NO_INTERLACE","MED_UNDEF_INTERLACE",};
char *_MED_STORAGE_MODE_MSG[3]={"MED_NO_STMODE","MED_GLOBAL_STMODE", "MED_COMPACT_STMODE"};
med_int _ipoint = nvaluesperentity;
sprintf(_filename,"%s_CPU-%03d_@_%s_%s.med",fieldnameprefix,mpi_size,_MED_MODE_SWITCH_MSG[constituentmode],_MED_STORAGE_MODE_MSG[storagemode]);
goto ERROR;
}
goto ERROR;
};
componentname = (
char*) malloc((nconstituentpervalue*
MED_SNAME_SIZE+1)*
sizeof(char));
componentunit = (
char*) malloc((nconstituentpervalue*
MED_SNAME_SIZE+1)*
sizeof(char));
strcpy(componentname,"");
strcpy(componentunit,"");
strcpy(_fieldname,fieldnameprefix);
if (
MEDfieldCr(_fid,_fieldname,
MED_FLOAT64,nconstituentpervalue,componentname,componentunit,
"s",_meshname ) < 0) {
goto ERROR;
};
free(componentname);
free(componentunit);
if ( _ipoint > 1 ) {
MESSAGE(
"Creating a localization of integration points...");
strcpy(_ipointname,_fieldname);
strcat(_ipointname,"_loc");
if (
MEDlocalizationWr(_fid, _ipointname, _geotype, _geotype/100, _ipointrefcoo, constituentmode,
goto ERROR;
}
free(_ipointrefcoo );
free(_ipointcoo );
free(_ipointwg );
} else {
}
if (profilearraysize) {
strcpy(_profilename,_fieldname);strcat(_profilename,"_profile");
for (_i=0; _i < profilearraysize; ++_i) _profilearray[_i]=_i;
if (
MEDprofileWr(_fid,_profilename,profilearraysize,_profilearray) < 0) {
goto ERROR;
};
_nusedentities = profilearraysize;
} else {
}
MESSAGE(
"Generating partition...");
getBlockOfEntities ( mpi_rank , mpi_size, _nusedentities,
&_start, &_stride, &_io_count, &_blocksize,
&_lastusedrank, &_lastblocksize);
_count=_io_count;
_start,_stride,_count,_blocksize,_lastblocksize, &filter) < 0 ) {
goto ERROR;
}
generateDatas(mpi_rank, _lastusedrank,
sizeof(
med_float),
storagemode, profilearraysize, _profilearray,
_start, _stride, _count, _blocksize, _lastblocksize,
nentities, nvaluesperentity, nconstituentpervalue,
&_arrayvalues );
_ipointname, &filter, (unsigned char*)_arrayvalues ) < 0) {
goto ERROR;
}
H5Fflush(_fid, H5F_SCOPE_GLOBAL );
if (mpi_rank == 0 ) {
int _ind=0;
FILE * _asciifile;
char _asciifilename[255]="";
goto ERROR;
}
sprintf(_asciifilename,"%s_CPU-%03d_@_%s_%s.ascii",fieldnameprefix,mpi_size,_MED_MODE_SWITCH_MSG[constituentmode],_MED_STORAGE_MODE_MSG[storagemode]);
_asciifile=fopen(_asciifilename, "w");
profilearraysize, _profilearray,
goto ERROR;
}
}
fprintf(_asciifile,"\n") ;
if ( profilearraysize ) {
_nentitiesarrayvalues = profilearraysize;
} else {
_nentitiesarrayvalues = nentities;
}
_filteredarrayvalues = (
med_float*) malloc(_nentitiesarrayvalues*
nvaluesperentity*
for (_i=0;_i<_nentitiesarrayvalues*nvaluesperentity*nconstituentpervalue; ++_i)
_filteredarrayvalues[_i]=-_i;
goto ERROR;
}
&filter2, (unsigned char*)_filteredarrayvalues ) < 0) {
goto ERROR;
}
switch (constituentmode) {
for (_j=0; _j < nvaluesperentity; ++_j)
for (_k=0; _k < nconstituentpervalue; ++_k) {
_ind = (cominfo->
filterarray[_i]-1)*nvaluesperentity*nconstituentpervalue+ _j*nconstituentpervalue+_k;
fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
}
break;
for (_k=0; _k < nvaluesperentity; ++_k)
for (_i=0; _i < nconstituentpervalue; ++_i) {
_ind =_i*nentities*nvaluesperentity+ (cominfo->
filterarray[_j]-1)*nvaluesperentity +_k;
fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
}
break;
}
} else
switch (constituentmode) {
for (_j=0; _j < nvaluesperentity; ++_j)
for (_k=0; _k < nconstituentpervalue; ++_k) {
_ind = _i*nvaluesperentity*nconstituentpervalue+_j*nconstituentpervalue+_k;
fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
}
break;
for (_k=0; _k < nvaluesperentity; ++_k)
for (_i=0; _i < nconstituentpervalue; ++_i) {
fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
}
break;
}
free(_filteredarrayvalues);
fclose(_asciifile);
goto ERROR;
}
}
goto ERROR;
}
_ret=0;
ERROR:
if (_arrayvalues) free(_arrayvalues);
if (profilearraysize) free(_profilearray);
}
if (mpi_rank == 0 ) {
}
}
return _ret;
}
int main (
int argc,
char **argv)
{
_cominfo.
comm = MPI_COMM_WORLD;
_cominfo.
info = MPI_INFO_NULL;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &(_cominfo.
mpi_size));
MPI_Comm_rank(MPI_COMM_WORLD, &(_cominfo.
mpi_rank));
int _nentities = 0;
int _nvaluesperentity = 0;
int _nconstituentpervalue = 0;
struct tm *_tm ;
time_t _tt=time(0);
_tm = localtime(&_tt);
srandom((*_tm).tm_sec * (*_tm).tm_min );
_nbblocksperproc = 1 + (int) (_cominfo.
mpi_size * (random() / (RAND_MAX + 1.0)));
_nentities = 1 + (int) (1000.0 * (random() / (RAND_MAX + 1.0)));
_nvaluesperentity = 1 + (int) (11.0 * (random() / (RAND_MAX + 1.0)));
_nconstituentpervalue = 1 + (int) (7.0 * (random() / (RAND_MAX + 1.0)));
}
if ( (
sizeof(
med_size)%(
sizeof(MPI_LONG)))==0 ) {
MPI_Bcast(&_nbblocksperproc ,
sizeof(
med_size)/
sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
MPI_Bcast(&_nentities ,
sizeof(
med_size)/
sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
MPI_Bcast(&_nvaluesperentity ,
sizeof(
med_size)/
sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
MPI_Bcast(&_nconstituentpervalue ,
sizeof(
med_size)/
sizeof(MPI_LONG), MPI_LONG, 0, MPI_COMM_WORLD);
} else {
assert(
sizeof(
med_size) == (
sizeof(MPI_LONG)));
}
char _fieldnameprefix[256] = "";
sprintf(_fieldnameprefix,"NENT-%03d_NVAL-%03d_NCST-%03d_NBL-%03llu",_nentities,_nvaluesperentity,
_nconstituentpervalue,_nbblocksperproc);
if ( (_storagemode ==
MED_GLOBAL_STMODE ) && (_profilearraysize) ) _profilearraysize=0;
_storagemode, _profilearraysize, _fieldnameprefix, & _cominfo) < 0 ) {
goto ERROR;
}
}
}
_ret = 0;
ERROR:
MPI_Finalize();
return _ret;
}