Concepts et fonctionnalités avancés sur les champs de résultat

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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/* 
 * Field use case 7 : write a field with computing steps and profiles
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>


int main (int argc, char **argv) {
  med_idt fid;
  const char meshname[MED_NAME_SIZE+1] = "2D unstructured mesh";
  const char fieldname[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD";
  const med_int ncomponent = 1;
  const char componentname[MED_SNAME_SIZE+1] = "TEMPERATURE";
  const char componentunit[MED_SNAME_SIZE+1] = "C";
  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 char profile1name[MED_NAME_SIZE+1] = "MED_TRIA3_PROFILE1";
  const med_int profile1[3] = {1, 4, 8}; 
  const med_int profile1size = 3;
  const char profile2name[MED_NAME_SIZE+1] = "MED_TRIA3_PROFILE2";
  const med_int profile2[5] = {2, 3, 5, 6, 7}; 
  const med_int profile2size = 5;
  const med_int ntria3 = 8;
  const med_int nquad4 = 4;
  int ret=-1;

  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDfield_7.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }

  /* create mesh link */
  if (MEDlinkWr(fid,meshname,"./UsesCase_MEDmesh_1.med") < 0) {
    MESSAGE("ERROR : create mesh link ...");
    goto ERROR;
  }
  
  /* create the profiles in the file */
  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; 
  }
  
  /* field creation : temperature field  : 1 component in celsius degree
   *                  the mesh is the 2D unstructured mesh of UsecaseMEDmesh_1.c
   *                  use case.
   *                  Computation step unit in 'ms'
   */ 
  if (MEDfieldCr(fid, fieldname, MED_FLOAT64, ncomponent, 
                 componentname, componentunit,"ms", meshname) < 0) {
    MESSAGE("ERROR : create field");
    goto ERROR;
  }
  
  /* two computation steps */
  /* write values at cell centers : 8 MED_TRIA3 and 4 MED_QUAD4 */

  /* STEP 1 : dt1 = 5.5, it = 1*/
  /* MED_TRIA3 : with a profile of 3 values in compact memory storage mode */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 1, 1, 5.5, MED_CELL,MED_TRIA3, 
                                 MED_COMPACT_PFLMODE, profile1name, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, ntria3,
                                 (unsigned char*) tria3values_step1_profile1) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3");
    goto ERROR;
  }
   /* MED_QUAD4  : with no profile */ 
  if (MEDfieldValueWithProfileWr(fid, fieldname, 1, 1, 5.5, MED_CELL, MED_QUAD4, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, nquad4,
                                 (unsigned char*) quad4values_step1) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ");
    goto ERROR;
  }

  /* STEP 2 : dt2 = 8.9, it = 1*/
  /* MED_TRIA3 : with a profile of 3 values then a profile of 5 values in global memory storage mode */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2 , 1 , 8.9 , MED_CELL, MED_TRIA3, 
                                 MED_GLOBAL_PFLMODE, profile1name, MED_NO_LOCALIZATION,    
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, ntria3,
                                 (unsigned char*) tria3values_step2_profile1) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3 ...");
    goto ERROR;
  }
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2 , 1 , 8.9 , MED_CELL, MED_TRIA3, 
                                 MED_GLOBAL_PFLMODE, profile2name, MED_NO_LOCALIZATION,    
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, ntria3,
                                 (unsigned char*) tria3values_step2_profile2) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3 ...");
    goto ERROR;
  }

  /* MED_QUAD4 : with no profile */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2, 1, 8.9, MED_CELL, MED_QUAD4, 
                      MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                      MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, nquad4,
                      (unsigned char*) quad4values_step2) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ... ");
    goto ERROR;
  }

  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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_PFLMODE ou MED_COMPACT_PFLMODE). Les routines MEDfieldValueWithProfileRd / mfdrpr, mfdipr permettent de lire les valeurs du champ avec ou sans profil.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

/* 
 * Field use case 9 : read a field with computing steps and profiles (generic approach)
 */


int main (int argc, char **argv) {
  med_idt fid;
  med_idt nfield, i, j;
  char meshname[MED_NAME_SIZE+1]="";
  med_bool localmesh;
  char fieldname[MED_NAME_SIZE+1]="";
  med_field_type fieldtype;
  char *componentname = NULL;
  char *componentunit = NULL;
  char dtunit[MED_SNAME_SIZE+1]="";
  med_float *values = NULL;
  med_int nstep, nvalues;
  med_int ncomponent;
  med_int csit, numit, numdt, it;
  med_float dt;
  med_geometry_type geotype;
  med_geometry_type *geotypes = MED_GET_CELL_GEOMETRY_TYPE;
  med_int nprofile, pit, profilesize;
  char profilename[MED_NAME_SIZE+1]="";
  med_int nintegrationpoint;
  char localizationname[MED_NAME_SIZE+1]="";
  int ret=-1;

  /* open file */
  fid = MEDfileOpen("UsesCase_MEDfield_7.med",MED_ACC_RDONLY);
  if (fid < 0) {
    MESSAGE("ERROR : open file ...");
    return -1;
  }

  /*
   * generic approach : how many fields in the file and identification
   * of each field.
   */
  if ((nfield = MEDnField(fid)) < 0) {
    MESSAGE("ERROR : How many fields in the file ...");
    return -1;
  }

  /* 
   * read values for each field
   */
  for (i=0; i<nfield; i++) {

    if ((ncomponent = MEDfieldnComponent(fid,i+1)) < 0) {
      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) {
      MESSAGE("ERROR : Field info ...");
      free(componentname);
      free(componentunit);
      return -1;
    }

    free(componentname);
    free(componentunit);

    /*
     * Read field values for each computing step 
     */ 
    for (csit=0; csit<nstep; csit++) {

      if (MEDfieldComputingStepInfo(fid, fieldname, csit+1, &numdt, &numit, &dt) < 0) {
        MESSAGE("ERROR : Computing step info ...");
        return -1;
      }

      /* 
       * ... In our case, we suppose that the field values are only defined on cells ... 
       */
      for (it=1; it<=MED_N_CELL_FIXED_GEO; it++) {

        geotype = geotypes[it];
        /*
         * How many profile for each geometry type ? 
         */
        if ((nprofile = MEDfieldnProfile(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                         profilename, localizationname)) < 0) {
          MESSAGE("ERROR : read number of profile ");
          return -1;
        }

        /* 
         * Read values for each profile
         * We suppose that we have no localization name 
         */
        for (pit=0; pit<nprofile; pit++) {
          
          if ((nvalues = MEDfieldnValueWithProfile(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                                   pit+1, MED_COMPACT_PFLMODE, profilename, &profilesize,
                                                   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)) == NULL) {
              MESSAGE("ERROR : memory allocation ...");
              return -1;
            }
            if (MEDfieldValueWithProfileRd(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                           MED_COMPACT_PFLMODE, profilename, 
                                           MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, (unsigned char*)values) < 0) {
              MESSAGE("ERROR : read fields values for cells ..."); 
              free(values);
              return -1; 
            }  

            free(values);
          }
        }
      }
    }
  }

  ret=0;
 ERROR:
  
  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/******************************************************************************
 * - Nom du fichier : test11.c
 *
 * - Description : lecture de champs de resultats MED
 *
 *****************************************************************************/

#include <med.h>
#define MESGERR 1
#include "med_utils.h"
#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_PFLMODE

med_err getFieldsOn(med_idt fid, char * nomcha, med_field_type typcha, med_int ncomp,
                    med_entity_type entite, med_switch_mode stockage, med_int ncstp);

int main (int argc, char **argv)


{
  med_err ret,lret;
  med_idt fid;
  char * fichier = NULL;
  char maa[MED_NAME_SIZE+1]="";
  char desc[MED_COMMENT_SIZE+1]="";
  char pflname[MED_NAME_SIZE+1]="",nomlien[MED_NAME_SIZE+1]="";
  char _meshname [MED_NAME_SIZE+1]="";
  char _dtunit [MED_SNAME_SIZE+1]="";
  char locname[MED_NAME_SIZE+1]="";
  char * lien = NULL;
  char *comp= NULL, *unit= NULL;
  char nomcha  [MED_NAME_SIZE+1]="";
  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;
  med_field_type    typcha;
  med_geometry_type type_geo;
  med_float *refcoo, *gscoo, *wg;
  int i,j;
  med_bool _local;

  char dtunit[MED_SNAME_SIZE+1]="";
  char nomcoo[3*MED_SNAME_SIZE+1]="";
  char unicoo[3*MED_SNAME_SIZE+1]="";
  char geointerpname[MED_NAME_SIZE+1]="";
  char ipointstructmeshname[MED_NAME_SIZE+1]="";
  med_mesh_type type;
  med_sorting_type sort;
  med_int nstep=0;
  med_axis_type rep;
  med_int nsectionmeshcell;
  med_geometry_type sectiongeotype;

  if (argc != 2) {
    MESSAGE("Aucun nom de fichier precise, fichier test10.med utilise ");
    fichier = "test10.med";
  } else {
    fichier = argv[1];
  };


  /* Ouverture du fichier med */
  if ((fid = MEDfileOpen(fichier,MED_ACC_RDONLY)) < 0){
    MESSAGE("Erreur a l'ouverture du fichier : ");SSCRUTE(fichier);
    return -1;
  }

  ret = 0;


 /* Lecture des infos concernant le premier maillage */
  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);
  }


  /* combien de champs dans le fichier */
  if ((ncha = MEDnField(fid)) < 0) {
    MESSAGE("Impossible de lire le nombre de champs : ");ISCRUTE(ncha);
    return ncha;
  }

  printf("Nombre de champs : "IFORMAT" \n",ncha);

  /* lecture de tous les champs  */
  for (i =0;i<ncha;i++) {
    lret = 0;
    printf("\nChamp numero : %d \n",i+1);

    /* Lecture du nombre de composantes */
    if ((ncomp = MEDfieldnComponent(fid,i+1)) < 0) {
      MESSAGE("Erreur a la lecture du nombre de composantes : "); ISCRUTE(ncomp);
      ret = -1; continue;
    }

    /* Lecture du type du champ, des noms des composantes et du nom de l'unite*/
    comp = (char*) malloc(ncomp*MED_SNAME_SIZE+1);
    EXIT_IF(comp == NULL,NULL,NULL);
    unit = (char*) malloc(ncomp*MED_SNAME_SIZE+1);
    EXIT_IF(unit == NULL,NULL,NULL);

    if ( MEDfieldInfo(fid,i+1,nomcha,_meshname,&_local,&typcha,comp,unit,_dtunit,&_ncstp) < 0 ) {
      MESSAGE("Erreur a la demande d'information sur les champs : ");
      ISCRUTE_int(i+1);SSCRUTE(nomcha);ISCRUTE_int(typcha);SSCRUTE(comp);SSCRUTE(unit);
      ISCRUTE(ncomp);
      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 de séquences de calcul |"IFORMAT"|\n",_ncstp);

      /* Le maillage reference est-il porte par un autre fichier */
    if ( !_local ) {

      if ( (lnsize=MEDlinkInfoByName(fid,_meshname) ) < 0 )  {
        MESSAGE("Erreur a la lecture de la taille du lien : ");
        SSCRUTE(_meshname);
        ret = -1;
      } else {

        lien = malloc((lnsize+1)*sizeof(char));
          EXIT_IF(lien == NULL,NULL,NULL);

          if ( MEDlinkRd(fid,_meshname, lien) < 0 )  {
            MESSAGE("Erreur a la lecture du lien : ");
            SSCRUTE(_meshname);SSCRUTE(lien);
            ret = -1;
          } else {
            printf("\tLe maillage |%s| est porte par un fichier distant |%s|\n",_meshname,lien);
          }
          free(lien);
        }
      }
    
    free(comp);
    free(unit);
    
      
    lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_NODE, USER_INTERLACE,_ncstp );
    
    if (lret == 0) lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_CELL, USER_INTERLACE,_ncstp );
    else { MESSAGE("Erreur a la lecture des champs aux noeuds "); ret = -1; continue;}
   
    if (lret == 0) lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_DESCENDING_FACE,USER_INTERLACE,_ncstp);
    else { MESSAGE("Erreur a la lecture des champs aux mailles "); ret = -1; continue;}
   
    if (lret == 0) lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_DESCENDING_EDGE,USER_INTERLACE,_ncstp);
    else {MESSAGE("Erreur a la lecture des champs aux faces "); ret = -1; continue;}
    
    if (lret == 0) lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_NODE_ELEMENT,USER_INTERLACE,_ncstp);
    else {MESSAGE("Erreur a la lecture des champs aux aretes"); ret = -1; continue;}
    
    /*TODO */
/*     if (lret == 0) lret = getFieldsOn(fid, nomcha, typcha, ncomp, MED_STRUCT_ELEMENT,USER_INTERLACE,_ncstp); */
/*     else {MESSAGE("Erreur a la lecture des champs aux éléments de structure"); ret = -1; continue;} */
    
    if  (lret != 0) {MESSAGE("Erreur a la lecture des champs aux noeuds des mailles "); ret = -1;};
  }


  /* Interrogation des profils */
  npro = MEDnProfile(fid);

  printf("\nNombre de profils stockes : "IFORMAT"\n\n",npro);
  for (i=1 ; i <= npro ; i++ ) {
    if ( MEDprofileInfo(fid, i, pflname, &nval) < 0)  {
      MESSAGE("Erreur a la demande d'information sur le profil n° : "); ISCRUTE_int(i);
      ret = -1;continue;
    }
    printf("\t- Profil n°%i de nom |%s| et de taille "IFORMAT"\n",i,pflname,nval);
    pflval = (med_int*) malloc(sizeof(med_int)*nval);
    if ( MEDprofileRd(fid, pflname, pflval) < 0) {
      MESSAGE("Erreur a la lecture des valeurs du profil : ");
      SSCRUTE(pflname);
      ret = -1;
    } else {
      printf("\t");
      for (j=0;j<nval;j++) printf(" "IFORMAT" ",*(pflval+j));
      printf("\n\n");
    }
    free(pflval);
  }

  /* Interrogation des liens */
  nln = MEDnLink(fid);

  printf("\nNombre de liens stockes : "IFORMAT"\n\n",nln);
  for (i=1 ; i <= nln ; i++ ) {
    if ( MEDlinkInfo(fid, i, nomlien, &nval) < 0)  {
      MESSAGE("Erreur a la demande d'information sur le lien n° : "); ISCRUTE_int(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));
    EXIT_IF(lien == NULL,NULL,NULL);

    if ( MEDlinkRd(fid, nomlien, lien ) < 0 )  {
      MESSAGE("Erreur a la lecture du lien : ");
      SSCRUTE(nomlien);SSCRUTE(lien);
      ret = -1;
    } else {
      lien[nval] = '\0';
      printf("\t\t|%s|\n\n",lien);
    }
    free(lien);
  }

  /* Interrogation des localisations des points de GAUSS */
  nloc = MEDnLocalization(fid);

  printf("\nNombre de localisations stockees : "IFORMAT"\n\n",nloc);
  for (i=1 ; i <= nloc ; i++ ) {
    if ( MEDlocalizationInfo(fid, i, locname, &type_geo, &locsdim,&ngauss,
                             geointerpname, ipointstructmeshname,&nsectionmeshcell,
                             &sectiongeotype) < 0)  {
      MESSAGE("Erreur a la demande d'information sur la localisation n° : "); ISCRUTE_int(i);
      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;
    refcoo = (med_float *) malloc(sizeof(med_float)*t1 );
    gscoo  = (med_float *) malloc(sizeof(med_float)*t2 );
    wg     = (med_float *) malloc(sizeof(med_float)*t3 );

    if ( MEDlocalizationRd(fid, locname, USER_INTERLACE, refcoo, gscoo, wg  ) < 0) {
      MESSAGE("Erreur a la lecture des valeurs de la localisation : ");
      SSCRUTE(locname);
      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);
  }


  /* fermeture du fichier */
  if ( MEDfileClose(fid) < 0) return -1;

  return ret;
}

med_err getFieldsOn(med_idt fid, char * nomcha, med_field_type typcha, med_int ncomp,
                    med_entity_type entite, med_switch_mode stockage, med_int ncstp) {

  int i,j,k,l,m,n,nb_geo=0;
  med_int nbpdtnor=0,pflsize,*pflval,ngauss=0,ngroup,*vale=NULL,nval;
  med_int numdt=0,numo=0,_nprofile;
  med_int meshnumdt=0,meshnumit=0;
  med_float *valr=NULL,dt=0.0;
  med_err ret=0;
  char pflname [MED_NAME_SIZE+1]="";
  char locname [MED_NAME_SIZE+1]="";
  char meshname [MED_NAME_SIZE+1]="";
  char * lien = NULL;
  char dt_unit [MED_SNAME_SIZE+1]="unknown";
  med_bool localmesh;
  med_int nmesh=0;

  med_geometry_type * type_geo;

  const char * const * AFF;
  const char * const * AFF_ENT=MED_GET_ENTITY_TYPENAME+1;
  switch (entite) {
  case MED_NODE :
    type_geo = MED_GET_NODE_GEOMETRY_TYPE;
    nb_geo   = MED_N_NODE_FIXED_GEO;
    AFF      = MED_GET_NODE_GEOMETRY_TYPENAME;
    break;
  case  MED_CELL :
  case  MED_NODE_ELEMENT :
    type_geo = MED_GET_CELL_GEOMETRY_TYPE;
    nb_geo   = MED_N_CELL_FIXED_GEO;
    AFF      = MED_GET_CELL_GEOMETRY_TYPENAME;
    break;
  case  MED_DESCENDING_FACE :
    type_geo = MED_GET_FACE_GEOMETRY_TYPE;
    nb_geo   = MED_N_FACE_FIXED_GEO;
    AFF      = MED_GET_FACE_GEOMETRY_TYPENAME;
    break;
  case  MED_DESCENDING_EDGE :
    type_geo = MED_GET_EDGE_GEOMETRY_TYPE;
    nb_geo   = MED_N_EDGE_FIXED_GEO;
    AFF      = MED_GET_EDGE_GEOMETRY_TYPENAME;
    break;
  }

  for (k=1;k<=nb_geo;k++) {

    /* Combien de (PDT,NOR) a lire */
    nbpdtnor = ncstp;
    if (nbpdtnor < 1 ) continue;

    for (j=0;j<nbpdtnor;j++) {

      if ( MEDfield23ComputingStepMeshInfo(fid,nomcha,j+1, &numdt, &numo, &dt,
                                           &nmesh, meshname,&localmesh, &meshnumdt, &meshnumit ) <0) {
        MESSAGE("Erreur a la demande d'information sur (pdt,nor) : ");
        ISCRUTE(numdt); ISCRUTE(numo);ISCRUTE(nmesh);SSCRUTE(meshname);ISCRUTE_int(localmesh);
        ISCRUTE(meshnumdt);ISCRUTE(meshnumit);
        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 : ");
          SSCRUTE(nomcha); ISCRUTE(numdt); ISCRUTE(numo);SSCRUTE(meshname);
          ISCRUTE_int(entite);ISCRUTE_int(type_geo[k]);SSCRUTE(pflname);SSCRUTE(locname);
          SSCRUTE(AFF_ENT[(int)entite]);SSCRUTE(AFF[k]);
          ret = -1; continue;
        };

        for (l=0;l<_nprofile;l++) {


          if ( (nval = MEDfield23nValueWithProfile(fid, nomcha, numdt, numo,  entite, type_geo[k],meshname,
                                                   l+1,  USER_MODE, pflname,&pflsize,
                                                   locname, &ngauss) ) < 0 ) {
            MESSAGE("Erreur a la lecture du nombre de valeurs du champ : ");
            SSCRUTE(nomcha);ISCRUTE(numdt);ISCRUTE(numo);SSCRUTE(meshname);
            ISCRUTE_int(entite);ISCRUTE_int(type_geo[k]);
            ISCRUTE_int(USER_MODE);
            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);

          /*Lecture des valeurs du champ */
          if (typcha == MED_FLOAT64) {

            valr = (med_float*) calloc(ncomp*nval*ngauss,sizeof(med_float));
            EXIT_IF(valr == NULL,NULL,NULL);

            if (MEDfield23ValueWithProfileRd(fid, nomcha, numdt,numo, entite,type_geo[k],meshname,
                                             USER_MODE, pflname, stockage,MED_ALL_CONSTITUENT,
                                           (unsigned char*) valr) < 0 ) {
              MESSAGE("Erreur a la lecture des valeurs du champ : ");
              SSCRUTE(nomcha);ISCRUTE_int(entite);ISCRUTE_int(type_geo[k]);
              ISCRUTE(numdt);ISCRUTE(numo);
              ret = -1;
            }
          } else {

            vale = (med_int*) calloc(ncomp*nval*ngauss,sizeof(med_int));
            EXIT_IF(vale == NULL,NULL,NULL);

            if (MEDfield23ValueWithProfileRd(fid, nomcha, numdt,numo, entite,type_geo[k],meshname,
                                           USER_MODE, pflname, stockage,MED_ALL_CONSTITUENT,
                                           (unsigned char*) vale) < 0 ) {
              MESSAGE("Erreur a la lecture des valeurs du champ : ");
              SSCRUTE(nomcha);ISCRUTE_int(entite);ISCRUTE_int(type_geo[k]);
              ISCRUTE(numdt);ISCRUTE(numo);
              ret = -1;
            };
          }

          if ( strlen(locname) )
            printf("\t- Modèle de localisation des points de Gauss de nom |%s|\n",locname);

          if (entite == MED_NODE_ELEMENT)
            ngroup = (type_geo[k] % 100);
          else
            ngroup = ngauss;

          switch (stockage) {

          case MED_FULL_INTERLACE :
            printf("\t- Valeurs :\n\t");
            for (m=0;m<(nval*ngauss)/ngroup;m++) {
              printf("|");
              for (n=0;n<ngroup*ncomp;n++)
                if (typcha == MED_FLOAT64)
                  printf(" %f ",*(valr+(m*ngroup*ncomp)+n));
                else
                  printf(" "IFORMAT" ",*(vale+(m*ngroup*ncomp)+n));

            }
            break;

            /*Affichage en fonction du profil à traiter*/
          case MED_NO_INTERLACE :
            printf("\t- Valeurs :\n\t");
            for (m=0;m<ncomp;m++) {
              printf("|");
              for (n=0;n<(nval*ngauss);n++)
                if (typcha == MED_FLOAT64)
                  printf(" %f ",*(valr+(m*nval)+n));
                else
                  printf(" "IFORMAT" ",*(vale+(m*nval)+n));
            }
            break;
          }

          printf("|\n");
          if (typcha == MED_FLOAT64) {
            if ( valr ) {free(valr);valr = NULL;}}
          else
            if (vale) { free(vale);vale = NULL; }

          /*Lecture du profil associe */
          if (strcmp(pflname,MED_NO_PROFILE) == 0 )
            printf("\t- Profil : MED_NO_PROFILE\n");
          else {
            if ( (pflsize = MEDprofileSizeByName(fid,pflname)) <0 )  {
              MESSAGE("Erreur a la lecture du nombre de valeurs du profil : ");
              SSCRUTE(pflname);
              ret = -1; continue;
            }

            printf("\t- Profil : |%s| de taille "IFORMAT"\n",pflname,pflsize);

            pflval = (med_int*) malloc(sizeof(med_int)*pflsize);
            EXIT_IF(pflval == NULL,NULL,NULL);
            if ( MEDprofileRd(fid,pflname,pflval) <0) {
              MESSAGE("Erreur a la lecture des valeurs du profil : ");
              SSCRUTE(pflname);
              ret = -1;
            }
            printf("\t");
            for (m=0;m<pflsize;m++) printf(" "IFORMAT" ",*(pflval+m));
            printf("\n");
            free(pflval);
          }
        }
      }
    }
  } /* fin for sur les mailles*/

  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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/* 
 * Field use case 10 : write a field in a MED file with computing steps,
 *                     profiles, integration points and interpolation families
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>


int main (int argc, char **argv) {
  med_idt fid;
  const char meshname[MED_NAME_SIZE+1] = "2D unstructured mesh";
  const char fieldname[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD";
  const med_int ncomponent = 1;
  const char componentname[MED_SNAME_SIZE+1] = "TEMPERATURE";
  const char componentunit[MED_SNAME_SIZE+1] = "C";
  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_int ntria3 = 8;
  const med_int nquad4 = 4;
  const med_float quad4values_step1[4] = {10000., 20000., 30000., 40000.};
  const med_float quad4values_step2[4] = {15000., 25000., 35000., 45000.};
  const char profile1name[MED_NAME_SIZE+1] = "MED_TRIA3_PROFILE1";
  const med_int profile1[3] = {1, 4, 8}; 
  const med_int profile1size = 3;
  const char profile2name[MED_NAME_SIZE+1] = "MED_TRIA3_PROFILE2";
  const med_int profile2[5] = {2, 3, 5, 6, 7}; 
  const med_int profile2size = 5;
  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};
  med_int nipoint, spacedim;
  const char interpname[MED_NAME_SIZE+1] = "MED_TRIA3 interpolation family";
  int ret=-1;

  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDfield_10.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }

  /* create mesh link */
  if (MEDlinkWr(fid,meshname,"./UsesCase_MEDmesh_1.med") < 0) {
    MESSAGE("ERROR : create mesh link ...");
    goto ERROR;
  }
  
  /* create the profiles in the file */
  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; 
  }

  /* create the localization elements for integration points */
  spacedim = 2;
  nipoint = 3;
  if (MEDlocalizationWr(fid, localization1name, MED_TRIA3, spacedim, 
                        elementcoordinate, MED_FULL_INTERLACE, 
                        nipoint, ipoint1coordinate, weight1,
                        MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT) < 0) {
    MESSAGE("ERROR : create famlily of integration points ...");
    goto ERROR; 
  }

  spacedim = 2;
  nipoint = 4;
  if (MEDlocalizationWr(fid, localization2name, MED_TRIA3, spacedim, 
                        elementcoordinate, MED_FULL_INTERLACE, 
                        nipoint, ipoint2coordinate, weight2, 
                        MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT) < 0) {
    MESSAGE("ERROR : create famlily of integration points ...");
    goto ERROR; 
  }

  /* 
   * Temperature field  creation : 
   * - 1 component 
   * - component unit : celsius degree
   * - mesh is the 2D unstructured mesh of UsecaseMEDmesh_1.c use case.
   * - computation step unit in 'ms'
   */ 
  if (MEDfieldCr(fid, fieldname, MED_FLOAT64, 
                 ncomponent, componentname, componentunit,
                 "ms", meshname) < 0) {
    MESSAGE("ERROR : create field");
    goto ERROR;
  }

  /* write interpolation family name for MED_TRIA3 cell type */
  /* The interpolation family "interpname" is created in the UsesCase_MEDinterp_1 
     use case */
  if (MEDfieldInterpWr(fid,fieldname,interpname) <0) {
    MESSAGE("ERROR : write field interpolation family name ...");
    goto ERROR;
  }  
  
  /* two computation steps */
  /* write values at cell centers : 8 MED_TRIA3 and 4 MED_QUAD4 */

  /* STEP 1 : dt1 = 5.5, it = 1*/
  /* MED_TRIA3 : with a profile of 3 values in compact memory storage mode 
     and a family of 3 integration points */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 1, 1, 5.5, MED_CELL,MED_TRIA3, 
                                 MED_COMPACT_PFLMODE, profile1name, localization1name,    
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                 ntria3, (unsigned char*) tria3values_step1_profile1) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3");
    goto ERROR;
  }
   /* MED_QUAD4  : with no profile */ 
  if (MEDfieldValueWithProfileWr(fid, fieldname, 1, 1, 5.5, MED_CELL, MED_QUAD4, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                 nquad4, (unsigned char*) quad4values_step1) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ");
    goto ERROR;
  }

  /* STEP 2 : dt2 = 8.9, it = 1*/
  /* MED_TRIA3 : with a profile of 3 values then a profile of 5 values in global memory storage mode 
   * For each profile, a family of 3 and then 4 integration points */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2 , 1 , 8.9 , MED_CELL, MED_TRIA3, 
                                 MED_GLOBAL_PFLMODE, profile1name, localization1name,    
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,  
                                 ntria3, (unsigned char*) tria3values_step2_profile1) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3 ...");
    goto ERROR;
  }
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2 , 1 , 8.9 , MED_CELL, MED_TRIA3, 
                                 MED_GLOBAL_PFLMODE, profile2name, localization2name,    
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,  
                                 ntria3, (unsigned char*) tria3values_step2_profile2) < 0) {
    MESSAGE("ERROR : write field values on MED_TRIA3 ...");
    goto ERROR;
  }

  /* MED_QUAD4 : with no profile */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2, 1, 8.9, MED_CELL, MED_QUAD4, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,  
                                 nquad4, (unsigned char*) quad4values_step2) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ... ");
    goto ERROR;
  }

  ret=0;
 ERROR:
  
  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

/* 
 * Field use case 12 : read a field (generic approach) in a MED file with computing steps,
 *                     profiles, integration points and interpolation families   
 */


int main (int argc, char **argv) {
  med_idt fid;
  med_idt nfield, i, j;
  char meshname[MED_NAME_SIZE+1]="";
  med_bool localmesh;
  char fieldname[MED_NAME_SIZE+1]="";
  med_field_type fieldtype;
  char *componentname = NULL;
  char *componentunit = NULL;
  char dtunit[MED_SNAME_SIZE+1]="";
  med_float *values = NULL;
  med_int nstep, nvalues;
  med_int ncomponent;
  med_int csit, numit, numdt, it;
  med_float dt;
  med_geometry_type geotype;
  med_geometry_type *geotypes = MED_GET_CELL_GEOMETRY_TYPE;
  med_int nprofile, pit, profilesize;
  char profilename[MED_NAME_SIZE+1]="";
  med_int nintegrationpoint;
  char localizationname[MED_NAME_SIZE+1]="";
  int k;
  med_int ninterp;
  char interpname[MED_NAME_SIZE+1];
  int ret=-1;


  /* open file */
  fid = MEDfileOpen("UsesCase_MEDfield_10.med",MED_ACC_RDONLY);
  if (fid < 0) {
    MESSAGE("ERROR : open file ...");
    return -1;
  }

  /*
   * generic approach : how many fields in the file and identification
   * of each field.
   */
  if ((nfield = MEDnField(fid)) < 0) {
    MESSAGE("ERROR : How many fields in the file ...");
    return -1;
  }

  /* 
   * read values for each field
   */
  for (i=0; i<nfield; i++) {

    if ((ncomponent = MEDfieldnComponent(fid,i+1)) < 0) {
      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) {
      MESSAGE("ERROR : Field info ...");
      free(componentname);
      free(componentunit);
      return -1;
    }

    free(componentname);
    free(componentunit);

    /* 
     * Read how many interpolation family name in the field ?  
     */
    if ((ninterp = MEDfieldnInterp(fid, fieldname)) < 0) {
      MESSAGE("ERROR : Read how many interpolation functions for the field ...");
      return -1;
    }
    /* - Read each interlolation family name 
     * - The way to read an interploation family and it's basis functions 
     *  is described in UsesCase_MEDinterp_2 and UsesCase_MEDinterp_3 uses case
     */
    for (it=0; it<ninterp; it++) {
      if (MEDfieldInterpInfo(fid,fieldname,it+1,interpname) < 0) {
        MESSAGE("ERROR : read interpolation family name ...");
        return -1;
      }
    }

    /*
     * Read field values for each computing step 
     */ 
    for (csit=0; csit<nstep; csit++) {

      if (MEDfieldComputingStepInfo(fid, fieldname, csit+1, &numdt, &numit, &dt) < 0) {
        MESSAGE("ERROR : Computing step info ...");
        return -1;
      }

      /* 
       * ... In our case, we suppose that the field values are only defined on cells ... 
       */
      for (it=1; it<=MED_N_CELL_FIXED_GEO; it++) {

        geotype = geotypes[it];

        /*
         * How many profile for each geometry type ? 
         */
        if ((nprofile = MEDfieldnProfile(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                         profilename, localizationname)) < 0) {
          MESSAGE("ERROR : read number of profile ");
          return -1;
        }

        /* 
         * Read values for each profile
         */
        for (pit=0; pit<nprofile; pit++) {
          
          if ((nvalues = MEDfieldnValueWithProfile(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                                   pit+1, MED_COMPACT_PFLMODE, profilename, &profilesize,
                                                   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;
            }
            if (MEDfieldValueWithProfileRd(fid, fieldname, numdt, numit, MED_CELL, geotype,
                                           MED_COMPACT_PFLMODE, profilename,
                                           MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                           (unsigned char*) values) < 0) {
              MESSAGE("ERROR : read fields values for cells ..."); 
              free(values);
              return -1; 
            }  
            free(values);
          }
        }
      }
    }
  }
  
  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    ret=-1; 
  } 
  
  return ret;
}

Le cas d'utilisation suivant montre un exemple de création d'une famille d'interpolation.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * Interp use case 1 : write an interpolation family
 * In this example, the interpolation family can be used with
 * the TEMPERATURE field of UsesCase_MEDfield_10 use case
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

int main (int argc, char **argv) {
  med_idt fid;
  char interpname[MED_NAME_SIZE+1] = "MED_TRIA3 interpolation family";
  const med_int nvariable=2;
  const med_int maxdegree=1;
  const med_int nmaxcoefficient=3;
  const med_int         ncoefficient1_1 = 3;
  const med_int   const power1_1[]         = {0,0,1,0,0,1};
  const med_float const coefficient1_1[]   = {1,-1,-1};

  const med_int         ncoefficient1_2 = 1;
  const med_int   const power1_2[]         = {1,0};
  const med_float const coefficient1_2[]   = {1};

  const med_int         ncoefficient1_3 = 1;
  const med_int   const power1_3[]         = {0,1};
  const med_float const coefficient1_3[]   = {1};
  int ret=-1;

  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDinterp_1.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }
 
  /* Family interpolation creation :
     - reference element = MED_TRIA3 
     - integration points of UsesCase_MEDfield_10 use case 
     are used to build the interpolation                 
     - basis functions are P1(X)= 1-X1-X2, P2(X)= X1, P3(X)= X2 
   */
  if (MEDinterpCr(fid, interpname, MED_TRIA3, MED_FALSE, nvariable, maxdegree, nmaxcoefficient) < 0) {
    MESSAGE("ERROR : interpolation family creation ...");
    goto ERROR;
  }

  /* Basis functions creation */
  if (MEDinterpBaseFunctionWr(fid,interpname,1,ncoefficient1_1,power1_1,coefficient1_1) < 0) {
    MESSAGE("ERROR : first base function creation ...");
    goto ERROR;
  }
  
  if (MEDinterpBaseFunctionWr(fid,interpname,2,ncoefficient1_2,power1_2,coefficient1_2) < 0) {
    MESSAGE("ERROR : second base function creation ...");
    goto ERROR;
  }

  if (MEDinterpBaseFunctionWr(fid,interpname,3,ncoefficient1_3,power1_3,coefficient1_3) < 0) {
    MESSAGE("ERROR : third base function creation ...");
    goto ERROR;
  }

  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    ret=-1; 
  } 
  
  return ret;
}

Le cas d'utilisation suivant montre un exemple de lecture de famille d'interpolation dont on connaît le nom.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * Interp use case 2 : read an interpolation family with a direct access by name
 * In this example, the interpolation family can be used with
 * the TEMPERATURE field of UsesCase_MEDfield_10 use case
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

int main (int argc, char **argv) {
  med_idt fid;
  const char interpname[MED_NAME_SIZE+1] = "MED_TRIA3 interpolation family";
  med_geometry_type geotype                    =MED_NONE;
  med_bool          cellnodes                  =MED_FALSE;
  med_int           nbasisfunc              =0;
  med_int           nvariable               =0;
  med_int           maxdegree                  =0;
  med_int           nmaxcoefficient            =0;
  int               basisfuncit                =0;
  int               powerit                    =0;
  med_int           ncoefficient            =0;
  med_int*          power                      =NULL;
  med_float*        coefficient                =NULL;
  int               coefficientit              =0;
  int ret=-1;


  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDinterp_1.med",MED_ACC_RDONLY);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }
 
  /* with direct access by the family name */
  if (MEDinterpInfoByName(fid,interpname,&geotype,&cellnodes,&nbasisfunc,
                          &nvariable,&maxdegree,&nmaxcoefficient) < 0) {
    MESSAGE("ERROR : interpolation function information ...");
    goto ERROR;
  }

  /* read each basis function */
  for ( basisfuncit=1; basisfuncit<= nbasisfunc; ++basisfuncit) {

    if ((ncoefficient = MEDinterpBaseFunctionCoefSize(fid,interpname,basisfuncit) ) <0 ) {
      MESSAGE("ERROR : read number of coefficient in the base function ...");
      goto ERROR;
    }

    coefficient = (med_float*) calloc(sizeof(med_float),ncoefficient);
    power       = (med_int*)   calloc(sizeof(med_int),nvariable*ncoefficient);
    
    if (MEDinterpBaseFunctionRd(fid,interpname,basisfuncit,&ncoefficient,power,coefficient) < 0) {
      MESSAGE("ERROR : read base function ...");
      free(coefficient); free(power);
      goto ERROR;
    }

    free(coefficient);
    free(power);
  }

  ret=0;
 ERROR:
  
  
  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    ret=-1; 
  } 
  
  return ret;
}

Le cas d'utilisation suivant montre un exemple de lecture de famille d'interpolation par une approche itérative.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * Interp use case 1 : read an interpolation family with a generic access by an iterator.
 * In this example, the interpolation family can be used with
 * the TEMPERATURE field of UsesCase_MEDfield_10 use case
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

int main (int argc, char **argv) {
  med_idt fid;
  char              interpname[MED_NAME_SIZE+1]="";
  med_int           ninterp                    =0;
  int               it                         =0;
  med_geometry_type geotype                    =MED_NONE;
  med_bool          cellnodes                  =MED_FALSE;
  med_int           nbasisfunc              =0;
  med_int           nvariable               =0;
  med_int           maxdegree                  =0;
  med_int           nmaxcoefficient            =0;
  int               basisfuncit                =0;
  int               powerit                    =0;
  med_int           ncoefficient            =0;
  med_int*          power                      =NULL;
  med_float*        coefficient                =NULL;
  int               coefficientit              =0;
  int ret=-1;


  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDinterp_1.med",MED_ACC_RDONLY);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }
 
  /* how many interpolation family in the file ? */
  if ((ninterp = MEDnInterp(fid)) < 0) {
    MESSAGE("ERROR : read number of interpolation ...");
    goto ERROR;
  }

  /* read each interpolation family */
  /* with an access by an iterator */
  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;
    }

    /* read each basis function */
    for ( basisfuncit=1; basisfuncit<= nbasisfunc; ++basisfuncit) {
      if ((ncoefficient = MEDinterpBaseFunctionCoefSize(fid,interpname,basisfuncit) ) <0 ) {
        MESSAGE("ERROR : read number of coefficient in the base function ...");
        goto ERROR;
      }
      
      coefficient = (med_float*) calloc(sizeof(med_float),ncoefficient);
      power       = (med_int*)   calloc(sizeof(med_int),nvariable*ncoefficient);
      
      if (MEDinterpBaseFunctionRd(fid,interpname,basisfuncit,&ncoefficient,power,coefficient) < 0) {
        MESSAGE("ERROR : read base function ...");
        free(coefficient); free(power);
        goto ERROR;
      }

      free(coefficient);
      free(power);
    }

  }

  ret=0;
 ERROR:
  
  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/* 
 * Field use case 13 : write a field on nodes elements in a MED file
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>


int main (int argc, char **argv) {
  med_idt fid;
  const char meshname[MED_NAME_SIZE+1] = "2D unstructured mesh";
  const char fieldname[MED_NAME_SIZE+1] = "TEMPERATURE_FIELD";
  const med_int ncomponent = 1;
  const char componentname[MED_SNAME_SIZE+1] = "TEMPERATURE";
  const char componentunit[MED_SNAME_SIZE+1] = "C";
  const med_int nquad4 = 4;
  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;
  
  /* file creation */
  fid = MEDfileOpen("UsesCase_MEDfield_13.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }

  /* create mesh link */
  if (MEDlinkWr(fid,meshname,"./UsesCase_MEDmesh_1.med") < 0) {
    MESSAGE("ERROR : create mesh link ...");
    goto ERROR;
  }
  
  /* 
   * Temperature field  creation : 
   * - 1 component 
   * - component unit : celsius degree
   * - mesh is the 2D unstructured mesh of UsecaseMEDmesh_1.c use case.
   * - computation step unit in 'ms'
   */ 
  if (MEDfieldCr(fid, fieldname, MED_FLOAT64, 
                 ncomponent, componentname, componentunit,
                 "ms", meshname) < 0) {
    MESSAGE("ERROR : create field");
    goto ERROR;
  }
  
  /* two computation steps */
  /* write values at nodes elements : 4 MED_QUAD4 */

  /* STEP 1 : dt1 = 5.5, it = 1*/
  /* MED_QUAD4  : with no profile */ 
  if (MEDfieldValueWithProfileWr(fid, fieldname, 1, 1, 5.5, MED_NODE_ELEMENT, MED_QUAD4, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                 nquad4, (unsigned char*) quad4values_step1) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ");
    goto ERROR;
  }

  /* STEP 2 : dt2 = 8.9, it = 1*/
  /* MED_QUAD4 : with no profile */
  if (MEDfieldValueWithProfileWr(fid, fieldname, 2, 1, 8.9, MED_NODE_ELEMENT, MED_QUAD4, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,  
                                 nquad4, (unsigned char*) quad4values_step2) < 0) {
    MESSAGE("ERROR : write field values on MED_QUAD4 ... ");
    goto ERROR;
  }

  ret=0;
 ERROR:
  
  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    ret=-1; 
  } 
  
  return ret;
}

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

/* 
 * Field use case 14 : read a field on nodes elements (generic approach) in a MED file
 */


int main (int argc, char **argv) {
  med_idt fid;
  med_idt nfield, i, j;
  char meshname[MED_NAME_SIZE+1]="";
  med_bool localmesh;
  char fieldname[MED_NAME_SIZE+1]="";
  med_field_type fieldtype;
  char *componentname = NULL;
  char *componentunit = NULL;
  char dtunit[MED_SNAME_SIZE+1]="";
  med_float *values = NULL;
  med_int nstep, nvalues;
  med_int ncomponent;
  med_int csit, numit, numdt, it;
  med_float dt;
  med_geometry_type geotype;
  med_geometry_type *geotypes = MED_GET_CELL_GEOMETRY_TYPE;
  med_int nprofile, pit, profilesize;
  char profilename[MED_NAME_SIZE+1]="";
  med_int nintegrationpoint;
  char localizationname[MED_NAME_SIZE+1]="";
  int k;
  int ret=-1;

  /* open file */
  fid = MEDfileOpen("UsesCase_MEDfield_13.med",MED_ACC_RDONLY);
  if (fid < 0) {
    MESSAGE("ERROR : open file ...");
    goto ERROR;
  }

  /*
   * generic approach : how many fields in the file and identification
   * of each field.
   */
  if ((nfield = MEDnField(fid)) < 0) {
    MESSAGE("ERROR : How many fields in the file ...");
    goto ERROR;
  }

  /* 
   * read values for each field
   */
  for (i=0; i<nfield; i++) {

    if ((ncomponent = MEDfieldnComponent(fid,i+1)) < 0) {
      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) {
      MESSAGE("ERROR : Field info ...");
      free(componentname);
      free(componentunit);
      goto ERROR;
    }

    free(componentname);
    free(componentunit);

    /*
     * Read field values for each computing step 
     */ 
    for (csit=0; csit<nstep; csit++) {

      if (MEDfieldComputingStepInfo(fid, fieldname, csit+1, &numdt, &numit, &dt) < 0) {
        MESSAGE("ERROR : Computing step info ...");
        goto ERROR;
      }

      /* 
       * ... In our case, we suppose that the field values are only defined on nodes element ... 
       */
      for (it=1; it<= MED_N_CELL_FIXED_GEO; it++) {

        geotype = geotypes[it];
        /*
         * How many profile for each geometry type ? 
         */
        if ((nprofile = MEDfieldnProfile(fid, fieldname, numdt, numit, MED_NODE_ELEMENT, geotype,
                                         profilename, localizationname)) < 0) {
          MESSAGE("ERROR : read number of profile ");
          goto ERROR;
        }

        /* 
         * Read values for each profile
         */
        for (pit=0; pit<nprofile; pit++) {
          
          if ((nvalues = MEDfieldnValueWithProfile(fid, fieldname, numdt, numit, MED_NODE_ELEMENT, geotype,
                                                   pit+1, MED_COMPACT_PFLMODE, profilename, &profilesize,
                                                   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;
            }
            if (MEDfieldValueWithProfileRd(fid, fieldname, numdt, numit, MED_NODE_ELEMENT, geotype,
                                           MED_COMPACT_PFLMODE, profilename,
                                           MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                           (unsigned char*) values) < 0) {
              MESSAGE("ERROR : read fields values for cells ..."); 
              free(values);
              goto ERROR; 
            }  
            free(values);
          }
        }
      }
    }
  }

  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * Field use case 15 : write a field in a MED file with 
 * values defined on struct elements
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>


int main (int argc, char **argv) {
  med_idt fid=0,mfid=0;
  const char meshname[MED_NAME_SIZE+1] = "COMPUT_MESH";
  const char fieldname[MED_NAME_SIZE+1] = "SPEED";
  const med_int ncomponent = 3;
  /*                                              123456789012345612345678901234561234567890123456 */   
  const char componentname[3*MED_SNAME_SIZE+1] = "Vx              Vy              Vz";
  const char componentunit[3*MED_SNAME_SIZE+1] = "m/s             m/s             m/s";
  med_geometry_type geotype=MED_NONE;
  const med_int npart = 3;
  char structelementname[MED_NAME_SIZE+1]="";
  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;


  /* File creation to write the field */
  fid = MEDfileOpen("UsesCase_MEDfield_15.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }

  /* Create mesh link */
  if (MEDlinkWr(fid,meshname,"./UsesCase_MEDstructElement_1.med") < 0) {
    MESSAGE("ERROR : create mesh link ...");
    goto ERROR;
  }

  /*
   * Read struct element geometric type
   */
  if (( mfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_ELSTRUCT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    goto ERROR;
  }

  strcpy(structelementname,MED_PARTICLE_NAME);
  geotype = MEDstructElementGeotype(fid,structelementname);


  /*
   * Speed field  creation for particles :
   * - 3 component
   * - component unit : m/s
   * - mesh is the 3D computation mesh of UsesCase_MEDstructElement_1 use case.
   * - computation step unit in 'ms'
   */
  if (MEDfieldCr(fid, fieldname, MED_FLOAT64,
                 ncomponent, componentname, componentunit,
                 "ms", meshname) < 0) {
    MESSAGE("ERROR : create field");
    goto ERROR;
  }

  if (MEDfieldValueWithProfileWr(fid, fieldname, MED_NO_DT, MED_NO_IT, MED_UNDEF_DT, MED_STRUCT_ELEMENT, geotype, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, MED_NO_LOCALIZATION,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
                                 npart, (unsigned char*) part_speed1) < 0) {
    MESSAGE("ERROR : write field values on MED_PARTICLE ");
    goto ERROR;
  }

  if ( MEDfileObjectsUnmount(fid, mfid, MED_ELSTRUCT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    goto ERROR;
  }

  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");
    ret=-1;
  }
  
  return ret;
}

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

/* 
 * Field use case 16 : read a field (generic approach) in a MED file with 
 * - values defined on struct elements
 */


int main (int argc, char **argv) {
  med_idt fid,mfid,sfid;
  med_idt nfield, i, j;
  char meshname[MED_NAME_SIZE+1]="";
  med_bool localmesh;
  char fieldname[MED_NAME_SIZE+1]="";
  med_field_type fieldtype;
  char *componentname = NULL;
  char *componentunit = NULL;
  char dtunit[MED_SNAME_SIZE+1]="";
  med_float *values = NULL;
  med_int nstep, nvalues;
  med_int ncomponent;
  med_int csit, numit, numdt, it;
  med_float dt;
  med_geometry_type *geotypes, geotype;
  med_int nmodels;
  med_int nprofile, pit, profilesize;
  char profilename[MED_NAME_SIZE+1]="";
  med_int nintegrationpoint;
  char localizationname[MED_NAME_SIZE+1]="";
  int k;
  char elementname[MED_NAME_SIZE+1];
  char supportmeshname[MED_NAME_SIZE+1];
  med_entity_type entitype;
  med_int elementdim;
  med_int nnode,ncell;
  med_geometry_type geocelltype;
  med_bool anyprofile=0;
  med_int nconstatt, *nvaratt;
  int ret=-1;

  /* open file */
  fid = MEDfileOpen("UsesCase_MEDfield_15.med",MED_ACC_RDWR);
  if (fid < 0) {
    MESSAGE("ERROR : open file ...");
    goto ERROR;
  }

  if (( mfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_MESH_SUPPORT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    goto ERROR;
  }

  if (( sfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_ELSTRUCT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    goto ERROR;
  }

  /* How many struct element models ? */
  if ((nmodels =  MEDnStructElement(mfid)) < 0) {
    MESSAGE("ERROR : read number of struct element models ...");
    goto ERROR;
  }
  geotypes = (med_geometry_type *) malloc(sizeof(med_geometry_type)*nmodels);
  nvaratt = (med_int *) malloc(sizeof(med_int)*nmodels);

  for (it=0; it<nmodels; it++) {
    if (MEDstructElementInfo(mfid, it+1, elementname, geotypes+it, &elementdim,
                             supportmeshname, &entitype, &nnode, &ncell,
                             &geocelltype, &nconstatt, &anyprofile, nvaratt+it) < 0) {
      MESSAGE("ERROR : struct element models information ...");
      free(nvaratt);
      goto ERROR;
    }
  }
  free(nvaratt);


  /*
   * generic approach : how many fields in the file and identification
   * of each field.
   */
  if ((nfield = MEDnField(fid)) < 0) {
    MESSAGE("ERROR : How many fields in the file ...");
    goto ERROR;
  }

  /*
   * read values for each field
   */
  for (i=0; i<nfield; i++) {

    if ((ncomponent = MEDfieldnComponent(fid,i+1)) < 0) {
      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) {
      MESSAGE("ERROR : Field info ...");
      free(componentname);
      free(componentunit);
      goto ERROR;
    }

    free(componentname);
    free(componentunit);

    /*
     * Read field values for each computing step 
     */ 
    for (csit=0; csit<nstep; csit++) {

      if (MEDfieldComputingStepInfo(fid, fieldname, csit+1, &numdt, &numit, &dt) < 0) {
        MESSAGE("ERROR : Computing step info ...");
        goto ERROR;
      }

      /* 
       * ... In our case, we suppose that the field values are only defined on struct element ... 
       */
      for (it=0; it<nmodels; it++) {

        geotype = *(geotypes+it);
     
        /*
         * How many profile for each geometry type ? 
         */
        if ((nprofile = MEDfieldnProfile(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                         profilename, localizationname)) < 0) {
          MESSAGE("ERROR : read number of profile ");
          goto ERROR;
        }
        
        /* 
         * Read values for each profile
         */
        for (pit=0; pit<nprofile; pit++) {
          
          if ((nvalues = MEDfieldnValueWithProfile(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                                   pit+1, MED_COMPACT_PFLMODE, profilename, &profilesize,
                                                   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;
            }
            if (MEDfieldValueWithProfileRd(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                           MED_COMPACT_PFLMODE, profilename,
                                           MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                           (unsigned char*) values) < 0) {
              MESSAGE("ERROR : read fields values for cells ..."); 
              goto ERROR; 
              free(values);
            }  
            free(values);
          }
        }
      }
    }
  }
  
  ret=0;
 ERROR :

  free(geotypes);

  if ( MEDfileObjectsUnmount(fid, sfid, MED_ELSTRUCT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    ret=-1;
  }

  if ( MEDfileObjectsUnmount(fid, mfid, MED_MESH_SUPPORT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    ret= -1;
  }

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * Field use case 17 : write a field in a MED file with 
 *                     values defined on integration points
 *                     on struct elements 
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>


int main (int argc, char **argv) {
  med_idt fid=0,mfid=0,sfid=0;
  const med_int spacedim = 3;
  const char meshname[MED_NAME_SIZE+1] = "COMPUT_MESH";
  const char fieldname[MED_NAME_SIZE+1] = "TEMPERATURE";
  const med_int ncomponent = 1;
  /*                                      1234567890123456123456*/   
  const char componentname[MED_SNAME_SIZE+1] = "TEMP            ";
  const char componentunit[MED_SNAME_SIZE+1] = "CELSIUS         ";
  med_geometry_type geotype=MED_NONE;
  const med_int nbeam = 1;
  char structelementname[MED_NAME_SIZE+1];
  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 char beamsectionname[MED_NAME_SIZE+1]="BEAM_SECTION_MESH";
  const med_int nipoint = 3;
  char interpname[MED_NAME_SIZE+1] = "geometrical transformation";
  const med_int nvariable=2;
  const med_int maxdegree=1;
  const med_int nmaxcoefficient=3;
  const med_int         ncoefficient1_1 = 3;
  const med_int   const power1_1[]         = {0,0,1,0,0,1};
  const med_float const coefficient1_1[]   = {1,-1,-1};
  const med_int         ncoefficient1_2 = 1;
  const med_int   const power1_2[]         = {1,0};
  const med_float const coefficient1_2[]   = {1};
  const med_int         ncoefficient1_3 = 1;
  const med_int   const power1_3[]         = {0,1};
  const med_float const coefficient1_3[]   = {1};
  int ret=-1;


  /* Open file to write the field */
  fid = MEDfileOpen("UsesCase_MEDfield_17.med",MED_ACC_CREAT);
  if (fid < 0) {
    MESSAGE("ERROR : file creation ...");
    goto ERROR;
  }

  if (( mfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_MESH_SUPPORT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    goto ERROR;
  }

  if (( sfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_ELSTRUCT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    goto ERROR;
  }

  /* Create mesh link */
  if (MEDlinkWr(fid,meshname,"./UsesCase_MEDstructElement_1.med") < 0) {
    MESSAGE("ERROR : create mesh link ...");
    goto ERROR;
  }

  /*
   * Read struct element geometric type
   */
  strcpy(structelementname,MED_BEAM_NAME);
  geotype = MEDstructElementGeotype(fid,structelementname);


  /* create a geometrical transformation fonction */
  if (MEDinterpCr(fid, interpname, geotype, MED_FALSE, nvariable, maxdegree, nmaxcoefficient) < 0) {
    MESSAGE("ERROR : interpolation family creation ...");
    goto ERROR;
  }
  /* Basis functions creation */
  if (MEDinterpBaseFunctionWr(fid,interpname,1,ncoefficient1_1,power1_1,coefficient1_1) < 0) {
    MESSAGE("ERROR : first base function creation ...");
    goto ERROR;
  }
  
  if (MEDinterpBaseFunctionWr(fid,interpname,2,ncoefficient1_2,power1_2,coefficient1_2) < 0) {
    MESSAGE("ERROR : second base function creation ...");
    goto ERROR;
  }

  if (MEDinterpBaseFunctionWr(fid,interpname,3,ncoefficient1_3,power1_3,coefficient1_3) < 0) {
    MESSAGE("ERROR : third base function creation ...");
    goto ERROR;
  }


  /* create the families of integration points 
     for the struct element */
  if (MEDlocalizationWr(fid, localization, geotype, spacedim, 
                        elementcoordinate, MED_FULL_INTERLACE, 
                        nipoint, ipointcoordinate, weight, 
                        MED_NO_INTERPOLATION, beamsectionname) < 0) {
    MESSAGE("ERROR : create famlily of integration points ...");
    goto ERROR; 
  }

  if (MEDlocalizationWr(fid, localization2, geotype, spacedim, 
                        elementcoordinate, MED_FULL_INTERLACE, 
                        nipoint, ipointcoordinate, weight, 
                        interpname, beamsectionname) < 0) {
    MESSAGE("ERROR : create famlily of integration points ...");
    goto ERROR; 
  }

  /*
   * Temperature  field creation for beam struct element :
   * - 1 component, 3 integration points, 4 cells in the support mesh for each 
   *   the section of each integration point 
   * - mesh is the 3D computation mesh of UsesCase_MEDstructElement_1 use case.
   */
  if (MEDfieldCr(fid, fieldname, MED_FLOAT64,
                 ncomponent, componentname, componentunit,
                 "ms", meshname) < 0) {
    MESSAGE("ERROR : create field");
    goto ERROR;
  }

  if (MEDfieldValueWithProfileWr(fid, fieldname, MED_NO_DT, MED_NO_IT, MED_UNDEF_DT, 
                                 MED_STRUCT_ELEMENT, geotype, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, localization,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
                                 nbeam, (unsigned char*) tempvalue) < 0) {
    MESSAGE("ERROR : write field values on MED_BEAM ");
    goto ERROR;
  }

  if (MEDfieldValueWithProfileWr(fid, fieldname, MED_NO_DT, 1, MED_UNDEF_DT, 
                                 MED_STRUCT_ELEMENT, geotype, 
                                 MED_COMPACT_PFLMODE, MED_NO_PROFILE, localization2,
                                 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
                                 nbeam, (unsigned char*) tempvalue) < 0) {
    MESSAGE("ERROR : write field values on MED_BEAM ");
    goto ERROR;
  }

  /* unmount file objects */


  if ( MEDfileObjectsUnmount(fid, mfid, MED_MESH_SUPPORT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    goto ERROR;
  }

  if ( MEDfileObjectsUnmount(fid, sfid, MED_ELSTRUCT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    goto ERROR;
  }

  ret=0;
 ERROR:

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");
    ret=-1;
  }
  
  return ret;
}

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

#include <med.h>
#define MESGERR 1
#include <med_utils.h>

#include <string.h>

/* 
 * Field use case 18 : read a field (generic approach) in a MED file with
 *                     values defined on integration points on struct elements
 */


int main (int argc, char **argv) {
  med_idt fid,mfid,sfid,cmfid;
  med_idt nfield, i, j;
  char meshname[MED_NAME_SIZE+1]="";
  med_bool localmesh;
  char fieldname[MED_NAME_SIZE+1]="";
  med_field_type fieldtype;
  char *componentname = NULL;
  char *componentunit = NULL;
  char dtunit[MED_SNAME_SIZE+1]="";
  med_float *values = NULL;
  med_int nstep, nvalues;
  med_int ncomponent;
  med_int csit, numit, numdt, it;
  med_float dt;
  med_geometry_type *geotypes=NULL, geotype;
  med_int nmodels;
  med_int nprofile, pit, profilesize;
  char profilename[MED_NAME_SIZE+1]="";
  med_int nintegrationpoint;
  char localizationname[MED_NAME_SIZE+1]="";
  int k;
  char elementname[MED_NAME_SIZE+1]="";
  char supportmeshname[MED_NAME_SIZE+1]="";
  med_entity_type entitype;
  med_int elementdim;
  med_int nnode,ncell;
  med_geometry_type geocelltype;
  med_bool anyprofile=0;
  med_int nconstatt, *nvaratt=NULL;
  med_bool coordinatechangement;
  med_bool geotransformation;
  int ret=-1;


  /* open file */
  fid = MEDfileOpen("UsesCase_MEDfield_17.med",MED_ACC_RDWR);
  if (fid < 0) {
    MESSAGE("ERROR : open file ...");
    return -1;
  }

  if (( mfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_MESH_SUPPORT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    return -1;
  }

  if (( sfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_ELSTRUCT)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    return -1;
  }

  if (( cmfid=MEDfileObjectsMount(fid,  "UsesCase_MEDstructElement_1.med",MED_MESH)) < 0 ) {
    MESSAGE("ERROR : file mounting ...");
    return -1;
  }


  /*
   * generic approach : how many fields in the file and identification
   * of each field.
   */
  if ((nfield = MEDnField(fid)) < 0) {
    MESSAGE("ERROR : How many fields in the file ...");
    return -1;
  }

  /*
   * read values for each field
   */
  for (i=0; i<nfield; i++) {

    if ((ncomponent = MEDfieldnComponent(fid,i+1)) < 0) {
      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) {
      MESSAGE("ERROR : Field info ...");
      free(componentname);
      free(componentunit);
      return -1;
    }

    free(componentname);
    free(componentunit);

    /* read how many struct element models in the mesh ? */
    if ((nmodels = MEDmeshnEntity(fid, meshname, MED_NO_DT, MED_NO_IT, MED_STRUCT_ELEMENT, MED_GEO_ALL,
                                  MED_UNDEF_DATATYPE, MED_NO_CMODE,&coordinatechangement,
                                  &geotransformation)) < 0) {
      MESSAGE("ERROR : number of nodes ...");
      return -1;
    }
    geotypes = (med_geometry_type *) malloc(sizeof(med_geometry_type)*nmodels);
    nvaratt = (med_int *) malloc(sizeof(med_int)*nmodels);
    
    for (it=0; it<nmodels; it++) {
      if (MEDstructElementInfo(mfid, it+1, elementname, geotypes+it, &elementdim,
                               supportmeshname, &entitype, &nnode, &ncell,
                               &geocelltype, &nconstatt, &anyprofile, nvaratt+it) < 0) {
        MESSAGE("ERROR : struct element models information ...");
        return -1;
      }
    }


    /*
     * Read field values for each computing step 
     */ 
    for (csit=0; csit<nstep; csit++) {

      if (MEDfieldComputingStepInfo(fid, fieldname, csit+1, &numdt, &numit, &dt) < 0) {
        MESSAGE("ERROR : Computing step info ...");
        return -1;
      }

      /* 
       * ... In our case, we suppose that the field values are only defined on struct element ... 
       */
      for (it=0; it<nmodels; it++) {
        geotype = *(geotypes+it);
        
        /*
         * How many profile for each geometry type ? 
         */
        if ((nprofile = MEDfieldnProfile(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                         profilename, localizationname)) < 0) {
          MESSAGE("ERROR : read number of profile ");
          return -1;
        }

        /* 
         * Read values for each profile
         */
        for (pit=0; pit<nprofile; pit++) {
          
          if ((nvalues = MEDfieldnValueWithProfile(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                                   pit+1, MED_COMPACT_PFLMODE, profilename, &profilesize,
                                                   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;
            }
            if (MEDfieldValueWithProfileRd(fid, fieldname, numdt, numit, MED_STRUCT_ELEMENT, geotype,
                                           MED_COMPACT_PFLMODE, profilename,
                                           MED_FULL_INTERLACE, MED_ALL_CONSTITUENT, 
                                           (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);

  }



  if ( MEDfileObjectsUnmount(fid, mfid, MED_MESH_SUPPORT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    ret=-1;
  }

  if ( MEDfileObjectsUnmount(fid, sfid, MED_ELSTRUCT) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    ret=-1;
  }


  if ( MEDfileObjectsUnmount(fid, cmfid, MED_MESH) < 0 ) {
    MESSAGE("ERROR : file unmounting ...");
    ret=-1;
  }

  /* close file */
  if (MEDfileClose(fid) < 0) {
    MESSAGE("ERROR : close file ...");             
    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.

/*  This file is part of MED.
 *
 *  COPYRIGHT (C) 1999 - 2015  EDF R&D, CEA/DEN
 *  MED is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU Lesser General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  MED is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Lesser General Public License for more details.
 *
 *  You should have received a copy of the GNU Lesser General Public License
 *  along with MED.  If not, see <http://www.gnu.org/licenses/>.
 */

#define MAX(a,b) ((a) > (b) ? (a) : (b))

#include <med.h>
#define MESGERR 1
#include "med_utils.h"
#include "med_config.h"

#include <stdlib.h>
#include <string.h>
#include <assert.h>

#if TIME_WITH_SYS_TIME
# include <sys/time.h>
# include <time.h>
#else
# if HAVE_SYS_TIME_H
#  include <sys/time.h>
# else
#  include <time.h>
# endif
#endif

#ifndef HAVE_UNISTD_H
#error "unistd.h reuired."
#endif

#include "getBlocksOfEntitiesPartition.h"
#include "generateDatas.h"
#include "generateFilterArray.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

/*Valeur de l'enum dans med.h.in*/
/* #define MED_NO_INTERLACE 2 */
/* #define MED_FULL_INTERLACE 1 */


/* #ifndef USER_INTERLACE */
/* #define USER_INTERLACE MED_FULL_INTERLACE */
/* #warning "Defining MED_FULL_INTERLACE mode..." */
/* #endif */

/* #if USER_INTERLACE == MED_NO_INTERLACE */
/* #define generateDatas generateNoIDatas */
/* #warning "Using generateNoIDatas..." */
/* #elif USER_INTERLACE == MED_FULL_INTERLACE */
/* #define generateDatas generateFullIDatas */
/* #warning "Using generateFullIDatas..." */
/* #else */
/* #error "The USER_INTERLACE macro value match neither MED_NO_INTERLACE nor MED_FULL_INTERLACE" */
/* #endif */

/* #define USER_MODE MED_COMPACT */


int main (int argc, char **argv)


{
  med_err _ret=0;
  med_idt _fid;

  int mpi_size, mpi_rank;
  MPI_Comm comm = MPI_COMM_WORLD;
  MPI_Info info = MPI_INFO_NULL;

  med_int    _nentitiesfiltered=0;
  med_int    *_filterarray=NULL;


  MPI_Init(&argc, &argv);
  MPI_Comm_size(MPI_COMM_WORLD, &mpi_size);
  MPI_Comm_rank(MPI_COMM_WORLD, &mpi_rank);


  med_err generateFieldFile( const med_size nentities, const med_size nvaluesperentity, const med_size nconstituentpervalue,
                             const med_switch_mode constituentmode,GetBlocksOfEntitiesType getBlockOfEntities, const med_int nbblocksperproc,
                             GenerateDataType generateDatas,
                             const med_storage_mode storagemode, const med_size profilearraysize,  const char * const fieldnameprefix ) {

/*     static int   _fileno=0; */
    med_err      _ret=-1;
    char         _filename   [255]="";
    char         _meshname[MED_NAME_SIZE+1]="Empty mesh";
    med_int      _meshdim=3;
    char         _meshcomponentname[3*MED_SNAME_SIZE+1] = "x               y               z               ";
    char         _meshcomponentunit[3*MED_SNAME_SIZE+1] = "cm              cm              cm              ";
    char         _fieldname  [MED_NAME_SIZE+1]="";
    char         *componentname,*componentunit;
    char         _profilename[MED_NAME_SIZE+1]=MED_NO_PROFILE;
    med_int       *_profilearray=0;
    int          _i=0,_j=0,_k=0, _lastusedrank=0;
    med_size     _blocksize=0,_lastblocksize=0,_count=0,_stride=0,_start=0,_index=0;
    med_float    *_arrayvalues;
    med_filter   filter = MED_FILTER_INIT;
    med_size     _nusedentities        = nentities;
    med_size     _io_count                = nbblocksperproc;
    med_idt      _fidseq;

    /*TODO : EXTERNALISER CES DEFINITIONS ET GENERALISER LE PRINCIPE */
    char         *_MED_MODE_SWITCH_MSG[3]={"MED_FULL_INTERLACE", "MED_NO_INTERLACE","MED_UNDEF_INTERLACE",};
    char         *_MED_STORAGE_MODE_MSG[3]={"MED_NO_PFLMODE","MED_GLOBAL_PFLMODE", "MED_COMPACT_PFLMODE"};

    med_geometry_type     _geotype       = MED_TRIA6;
    med_int               _geodim        = _geotype/100;
    med_int               _geonnodes     = _geotype%100;
    char       _ipointname[MED_NAME_SIZE+1];
    med_float* _ipointrefcoo = 0;
    med_int    _ipoint       = nvaluesperentity;
    med_float* _ipointcoo    = 0;
    med_float* _ipointwg     = 0;

    sprintf(_filename,"%s_CPU-%03d_@_%s_%s.med",fieldnameprefix,mpi_size,_MED_MODE_SWITCH_MSG[constituentmode],_MED_STORAGE_MODE_MSG[storagemode]);
/*     SSCRUTE(_filename); */
    /* Ouverture du fichier en mode parallel */
    if ((_fid = MEDparFileOpen(_filename, MODE_ACCES ,comm, info)) < 0){
      MED_ERR_(_ret,MED_ERR_OPEN,MED_ERR_FILE,_filename);
      goto ERROR;
    }

/*     SSCRUTE(_meshname); */
    if (MEDmeshCr( _fid,_meshname,_meshdim,_meshdim, MED_UNSTRUCTURED_MESH,
                   "Un maillage pour le test parallel","s", MED_SORT_DTIT,
                   MED_CARTESIAN, _meshcomponentname, _meshcomponentunit) < 0) {
      MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_MESH,_meshname);
      goto ERROR;
    };

    componentname = (char*) malloc((nconstituentpervalue*MED_SNAME_SIZE+1)*sizeof(char));
    componentunit = (char*) malloc((nconstituentpervalue*MED_SNAME_SIZE+1)*sizeof(char));
    /*TODO : Compléter le nom */
    strcpy(componentname,"");
    strcpy(componentunit,"");
    strcpy(_fieldname,fieldnameprefix);
    if ( MEDfieldCr(_fid,_fieldname,MED_FLOAT64,nconstituentpervalue,componentname,componentunit,"s",_meshname ) < 0) {
      MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FIELD,_fieldname);
      goto ERROR;
    };
    free(componentname);
    free(componentunit);


    if ( _ipoint > 1 ) {

      MESSAGE("Creating a localization of integration points...");
      strcpy(_ipointname,_fieldname);
      strcat(_ipointname,"_loc");

      /*Attention ancienne spec*/
      _ipointrefcoo = (med_float *) calloc(_geodim*_geonnodes,sizeof(med_float));
      _ipointcoo    = (med_float *) calloc(_ipoint*_geodim,sizeof(med_float));
      _ipointwg     = (med_float *) calloc(_ipoint,sizeof(med_float));

      if (MEDlocalizationWr(_fid, _ipointname, _geotype, _geotype/100, _ipointrefcoo, constituentmode,
                            _ipoint, _ipointcoo, _ipointwg, MED_NO_INTERPOLATION, MED_NO_MESH_SUPPORT ) < 0) {
        MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_LOCALIZATION,_ipointname);
        ISCRUTE_int(constituentmode);
        goto ERROR;
      }
      free(_ipointrefcoo );
      free(_ipointcoo    );
      free(_ipointwg     );

    } else {
      strcpy(_ipointname,MED_NO_LOCALIZATION);
    }

    if (profilearraysize) {
      MESSAGE("Creating a profile...");

      strcpy(_profilename,_fieldname);strcat(_profilename,"_profile");

      _profilearray = (med_int*) calloc(profilearraysize,sizeof(med_int));

      for (_i=0; _i < profilearraysize; ++_i) _profilearray[_i]=_i;
      if ( MEDprofileWr(_fid,_profilename,profilearraysize,_profilearray) < 0) {
        MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_PROFILE,_profilename);
        goto ERROR;
      };
      _nusedentities = profilearraysize;
    } else {
      strcpy(_profilename,MED_NO_PROFILE);
    }


    MESSAGE("Generating partition...");
    getBlockOfEntities ( mpi_rank , mpi_size, _nusedentities,
                         &_start, &_stride, &_io_count, &_blocksize,
                         &_lastusedrank, &_lastblocksize);

    _count=_io_count;
    MESSAGE("Generating filter...");
    if ( MEDfilterBlockOfEntityCr(_fid, nentities, nvaluesperentity, nconstituentpervalue,
                                  MED_ALL_CONSTITUENT, constituentmode, storagemode, _profilename,
                                  _start,_stride,_count,_blocksize,_lastblocksize,  &filter) < 0 ) {
        MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FILTER,"");
        goto ERROR;
    }

    MESSAGE("Generating datas...");
    generateDatas(mpi_rank, _lastusedrank, sizeof(med_float),
                  storagemode, profilearraysize, _profilearray,
                  _start, _stride, _count, _blocksize, _lastblocksize,
                  nentities, nvaluesperentity, nconstituentpervalue,
                  &_arrayvalues );

    MESSAGE("Writing field...");
    if ( MEDfieldValueAdvancedWr(_fid,_fieldname,MED_NO_DT,MED_NO_IT,0.0, MED_CELL, _geotype,
                                 _ipointname, &filter, (unsigned char*)_arrayvalues ) < 0) {
      MED_ERR_(_ret,MED_ERR_WRITE,MED_ERR_FIELD,_fieldname);
      ISCRUTE(mpi_rank);
      goto ERROR;
    }

    /* Test de lecture du même fichier avec filtre simple par un seul processeur */
    /* TODO : Créer MEDflush */
    H5Fflush(_fid, H5F_SCOPE_GLOBAL );

    /*Le flush suffit pas besoin de synchroniser les processus : MPI_Barrier(MPI_COMM_WORLD); */
    if (mpi_rank == 0 ) {
      MESSAGE("Reading field...");


      med_int    _nentitiesarrayvalues=0;
      med_float  *_filteredarrayvalues=NULL;
      med_filter filter2=MED_FILTER_INIT;
      int        _ind=0;
      FILE *     _asciifile;
      char       _asciifilename[255]="";


      if ((_fidseq = MEDfileOpen(_filename, MED_ACC_RDONLY )) < 0){
        MED_ERR_(_ret,MED_ERR_OPEN,MED_ERR_FILE,_filename);
        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");

      /*Génère un filtre de selection simple s'il n'a pas déjà été généré lors d'un précédent appel */
      /*TODO : Déplacer cette appel dans le main après avoir externaliser la génération du profile */
      if (!_filterarray)
        if ( generateFilterArray(  nentities,  nvaluesperentity, nconstituentpervalue,
                                   profilearraysize, _profilearray,
                                   &_nentitiesfiltered, &_filterarray ) < 0 ) {
          goto ERROR;
        }

      ISCRUTE(_nentitiesfiltered);
      /*Stocke le filtre utilisé dans le fichier .ascii*/
      for (_i=0; _i < _nentitiesfiltered; ++_i ) {
/*      ISCRUTE(_filterarray[_i]); */
        fprintf(_asciifile,"%d ",_filterarray[_i]) ;
      }
      fprintf(_asciifile,"\n") ;


      /*Pas de profile possible (profilearraysize == 0) en MED_GLOBAL_PFLMODE sur un fichier géré en parallel */
      if ( profilearraysize ) {
        _nentitiesarrayvalues = profilearraysize;
      } else {
        _nentitiesarrayvalues = nentities;
      }

      /*Attention allocation mémoire potentiellement grosse car réalisée uniquement par le processus 0
       qui rassemble les données.*/
      /* C'est une taille maxi qui ne prend pas en compte le COMPACT+filter */
      /* TODO : Ajuster la taille au storage_mode*/
      _filteredarrayvalues = (med_float*) malloc(_nentitiesarrayvalues*
                                                 nvaluesperentity*
                                                 nconstituentpervalue*sizeof(med_float));

      /* Permet de vérifier une erreur d'indiçage après la lecture */
      for (_i=0;_i<_nentitiesarrayvalues*nvaluesperentity*nconstituentpervalue; ++_i)
        _filteredarrayvalues[_i]=-_i;


      /*Création d'un filtre de sélection simple, pour une lecture séquentielle par le processys 0*/
      if ( MEDfilterEntityCr(_fidseq, nentities, nvaluesperentity, nconstituentpervalue,
                             MED_ALL_CONSTITUENT, constituentmode, storagemode, _profilename,
                             _nentitiesfiltered,_filterarray, &filter2) < 0 ) {
        MED_ERR_(_ret,MED_ERR_CREATE,MED_ERR_FILTER,"");
        goto ERROR;
      }

      if ( MEDfieldValueAdvancedRd(_fidseq,_fieldname,MED_NO_DT,MED_NO_IT, MED_CELL, _geotype,
                                   &filter2, (unsigned char*)_filteredarrayvalues ) < 0) {
        MED_ERR_(_ret,MED_ERR_READ,MED_ERR_FIELD,_fieldname);
        ISCRUTE(mpi_rank);
        goto ERROR;
      }

      /*AFFICHAGE TOUJOURS EN FULL INTERLACE QUELQUES SOIENT LES COMBINAISONS*/
      /*TODO : Externaliser l'affichage*/
      if ( storagemode == MED_GLOBAL_PFLMODE ) {
        switch (constituentmode) {
        case MED_FULL_INTERLACE:
          for (_i=0; _i < _nentitiesfiltered; ++_i)
            for (_j=0; _j < nvaluesperentity; ++_j)
              for (_k=0; _k < nconstituentpervalue; ++_k) {
                _ind = (_filterarray[_i]-1)*nvaluesperentity*nconstituentpervalue+ _j*nconstituentpervalue+_k;
/*              fprintf(stdout,"%s%3d%s = %f\n","_filteredarrayvaluesFULLGLB[",_ind,"]",_filteredarrayvalues[_ind]) ; */
                fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
              }
          break;
        case MED_NO_INTERLACE:
          for (_j=0; _j < _nentitiesfiltered; ++_j)
            for (_k=0; _k < nvaluesperentity; ++_k)
              for (_i=0; _i < nconstituentpervalue; ++_i) {
                _ind =_i*nentities*nvaluesperentity+ (_filterarray[_j]-1)*nvaluesperentity +_k;
/*              fprintf(stdout,"%s%3d%s = %f\n","_filteredarrayvaluesNOGLB[",_ind,"]",_filteredarrayvalues[_ind]); */
                fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
              }
          break;
        }
      }  else
        switch (constituentmode) {
        case MED_FULL_INTERLACE:
          for (_i=0; _i < _nentitiesfiltered; ++_i )
            for (_j=0; _j < nvaluesperentity; ++_j)
              for (_k=0; _k < nconstituentpervalue; ++_k) {
                _ind = _i*nvaluesperentity*nconstituentpervalue+_j*nconstituentpervalue+_k;
/*              fprintf(stdout,"%s%3d%s = %f\n","_filteredarrayvaluesFULLCP[",_ind,"]",_filteredarrayvalues[_ind]) ; */
                fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]) ;
          }
          break;
        case MED_NO_INTERLACE:
          for (_j=0; _j < _nentitiesfiltered; ++_j)
            for (_k=0; _k < nvaluesperentity; ++_k)
              for (_i=0; _i < nconstituentpervalue; ++_i) {
                _ind =_i*_nentitiesfiltered*nvaluesperentity+ _j*nvaluesperentity +_k;
                /* _ind =_i*_nentitiesarrayvalues*nvaluesperentity+ (_filterarray[_j]-1)*nvaluesperentity +_k; */
/*              fprintf(stdout,"%s%3d%s = %f\n","_filteredarrayvaluesNOCP[",_ind,"]",_filteredarrayvalues[_ind]); */
                fprintf(_asciifile,"%f\n",_filteredarrayvalues[_ind]);
              }
          break;
        }


      free(_filteredarrayvalues);

      fclose(_asciifile);

      if ( MEDfilterClose(&filter2) < 0 ) {
        MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILTER,"");
        goto ERROR;
      }

    } /*fin if (mpi_rank == 0) */

  if ( MEDfilterClose(&filter) < 0 ) {
    MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILTER,"");
    goto ERROR;
  }


    _ret=0;
  ERROR:
    if (_arrayvalues)     free(_arrayvalues);
    if (profilearraysize) free(_profilearray);

    if (  MEDfileClose(_fid) < 0) {
      MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,""); _ret = -1;
    }

    if (mpi_rank == 0 ) {
      if (  MEDfileClose(_fidseq) < 0) {
        MED_ERR_(_ret,MED_ERR_CLOSE,MED_ERR_FILE,""); _ret = -1;
      }
    }

    return _ret;
  }




  /* MAIN */
  med_size            _nbblocksperproc    = 0;
  int           _nentities             = 0;
  int           _nvaluesperentity      = 0;
  int           _nconstituentpervalue  = 0;

  if (mpi_rank == 0 ) {

    struct tm *_tm ;
    time_t _tt=time(0);
    _tm = localtime(&_tt);

    srandom((*_tm).tm_sec * (*_tm).tm_min );
    _nbblocksperproc         = 1 + (int) (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)));
/*     _nbblocksperproc         = 1 + (int) (mpi_size * (random() / (RAND_MAX + 1.0))); */
/*     _nentities            = 1 + (int) (5.0 * (random() / (RAND_MAX + 1.0))); */
/*     _nvaluesperentity     = 1 + (int) (3.0 * (random() / (RAND_MAX + 1.0))); */
/*     _nconstituentpervalue = 1 + (int) (3.0 * (random() / (RAND_MAX + 1.0))); */
/*     _nbblocksperproc         = 2; */
/*     _nentities            = 4; */
/*     _nvaluesperentity     = 1; */
/*     _nconstituentpervalue = 1; */

  }

  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);

  GenerateDataType generateDatas = 0;
  med_switch_mode  _switchmode  = MED_UNDEF_INTERLACE;
  med_storage_mode _storagemode = MED_UNDEF_PFLMODE;
  /*Pour que les 4 fichiers générés soient identiques, on désactive l'utilisation des profils
    qui n'est pas utilisable en mode GLOBAL et // */
  med_int          _profilearraysize=0;
  /*   med_int       _profilearraysize=_nentities/2; */

  for (_switchmode = MED_FULL_INTERLACE ; _switchmode <= MED_NO_INTERLACE; ++_switchmode) {

    if ( _switchmode == MED_FULL_INTERLACE ) generateDatas = generateFullIDatas;
    else generateDatas = generateNoIDatas;

    for (_storagemode = MED_GLOBAL_PFLMODE ; _storagemode <= MED_COMPACT_PFLMODE; ++_storagemode) {

      if ( (_storagemode == MED_GLOBAL_PFLMODE ) && (_profilearraysize) ) _profilearraysize=0;

      if ( generateFieldFile( _nentities,  _nvaluesperentity, _nconstituentpervalue,
                              _switchmode, getCyclicBlocksOfEntities, _nbblocksperproc, generateDatas,
                              _storagemode, _profilearraysize,  _fieldnameprefix ) < 0 ) {
        MED_ERR_(_ret,MED_ERR_WRITE,MED_ERR_FIELD,_fieldnameprefix);
        ISCRUTE(mpi_rank);
        goto ERROR;
      }

    }
  }


  _ret = 0;
 ERROR:

  if (_filterarray) free(_filterarray);

  /*pour arch. BLueGene : Sync entre GPFS et LSF : sleep(360) */

  /* MPI_Finalize must be called AFTER MEDclose which may use MPI calls */
  MPI_Finalize();

  /* Catcher l'erreur en retour mpirun  et .sh*/
  return _ret;
}







Généré le Thu Oct 8 14:27:13 2015 pour MED fichier par  doxygen 1.6.1