src/MedClient/src/create_mesh_c3h8.c

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 /*
  24   creation d'une geometrie 3d : un cube [0,1]^3
  25   maillé uniformement en hexahedres reguliers;
  26   avec n (=argv[1]) noeuds dans chaque direction.
  27   2 champs:
  28   - DbleVectNode champ vectoriel reel sur les noeuds
  29   - DbleVectCell champ vectoriel reel sur les cellules
  30
  31   En sortie, il y aura production d'un fichier MED
  32   cube_hexa8_n.med qui contiendra un seul maillage et 2 champs
  33   avec une seule famille la FAMILLE_0
  34 */
  35
  36 #include <med.h>
  37 #include <string.h>
  38 #include <stdlib.h>
  39
  40 #include <time.h>
  41
  42 int main (int argc, char **argv)
  43 {
  44   med_err ret;
  45   med_idt fid;
  46   char maa[MED_TAILLE_NOM+1] = "cube_hexa8";
  47   med_int mdim = 3;
  48
  49   int nnoe_dir;
  50   int nelt_dir;
  51   med_int nnoe;
  52   int i, j, k, ijk;
  53
  54   med_float * coo;
  55   med_int * numnoe;
  56   med_int * nufano;
  57
  58   med_float hxsize;
  59   med_float hysize;
  60   med_float hzsize;
  61
  62   med_float * DbleVectNode;
  63   med_float * DbleVectCell;
  64
  65   time_t t1;
  66
  67   /*
  68     Le maillage
  69    */
  70
  71   char nomcoo[3*MED_TAILLE_PNOM+1] = "x       y       z       ";
  72   char unicoo[3*MED_TAILLE_PNOM+1] = "cm      cm      cm      ";
  73   /*  char nomnoe[19*MED_TAILLE_PNOM+1] = "nom1    nom2    nom3    nom4";*/
  74   char *nomnoe ;
  75
  76   med_int nhexa8;
  77   med_int * hexa8;
  78   med_int * numhexa8;
  79   med_int * nufahexa8;
  80
  81   char * nomhexa8;
  82   int indexN1, indexN2, indexN3, indexN4, indexN5, indexN6, indexN7, indexN8;
  83
  84   char nomfam[MED_TAILLE_NOM+1];
  85   med_int numfam;
  86   char attdes[MED_TAILLE_DESC+1];
  87   med_int natt;
  88   med_int attide;
  89   med_int attval;
  90   med_int ngro;
  91   char gro[MED_TAILLE_LNOM+1];
  92   int nfame = 0;
  93   int nfamn = 0;
  94
  95   char MedFile[100] = "cube_hexa8_";
  96   char buff[100];
  97
  98   /*
  99     Les champs
 100   */
 101
 102   char champDbleVectNode[MED_TAILLE_NOM+1] = "DbleVectNode";
 103   char compDbleVectNode[MED_TAILLE_PNOM*3+1] = "comp1   comp2   comp3   " ;
 104   char unitDbleVectNode[MED_TAILLE_PNOM*3+1] = "unit1   unit2   unit3   " ;
 105
 106   char champDbleVectCell[MED_TAILLE_NOM+1] = "DbleVectCell";
 107   char compDbleVectCell[MED_TAILLE_PNOM*3+1] = "comp1   comp2   comp3   " ;
 108   char unitDbleVectCell[MED_TAILLE_PNOM*3+1] = "unit1   unit2   unit3   " ;
 109
 110   if (argc != 2)
 111     {
 112       printf("Usage: %s <n> \n",argv[0]);
 113       printf("       where\n");
 114       printf("       - n is the number of nodes in each direction.\n");
 115       printf("\n");
 116       printf("This program will produce a MED file cube_hexa8_n.med\n");
 117       exit(0);
 118     }
 119
 120   nnoe_dir = atoi(argv[1]);
 121   nelt_dir = nnoe_dir-1;
 122   nnoe = nnoe_dir*nnoe_dir*nnoe_dir;
 123
 124   coo = malloc(mdim*nnoe*sizeof(med_float));
 125   numnoe = malloc(nnoe*sizeof(med_int));
 126   nufano = malloc(nnoe*sizeof(med_int));
 127   nomnoe = malloc((MED_TAILLE_PNOM*nnoe+1)*sizeof(char));
 128
 129   hxsize = 1./((med_float) (nnoe_dir - 1));
 130   hysize = hxsize;
 131   hzsize = hxsize;
 132
 133   nhexa8 = nelt_dir*nelt_dir*nelt_dir;
 134   hexa8 = malloc(8*nhexa8*sizeof(med_int));
 135   numhexa8 = malloc(nhexa8*sizeof(med_int));
 136   nufahexa8 = malloc(nhexa8*sizeof(med_int));
 137   nomhexa8 = malloc((MED_TAILLE_PNOM*nhexa8+1)*sizeof(char));
 138
 139   DbleVectNode = malloc(mdim*nnoe*sizeof(med_float));
 140   DbleVectCell = malloc(mdim*nhexa8*sizeof(med_float));
 141
 142   /*
 143     les noeuds:
 144   */
 145
 146   for(k=0;k<nnoe_dir;k++)
 147     {
 148       for(j=0;j<nnoe_dir;j++)
 149         {
 150           for (i=0;i<nnoe_dir;i++)
 151             {
 152               int ijk = k*nnoe_dir*nnoe_dir+j*nnoe_dir+i;
 153
 154               numnoe[ijk] = ijk+1;
 155               nufano[ijk] = 0;
 156
 157               coo[mdim*ijk] = ((med_float) i)*hxsize;
 158               coo[mdim*ijk+1] = ((med_float) j)*hysize;
 159               coo[mdim*ijk+2] = ((med_float) k)*hzsize;
 160
 161               /*
 162               printf("Coordonnées %d   X = %lf  Y = %lf  Z = %lf\n",(ijk+1),coo[mdim*ijk],coo[mdim*ijk+1],coo[mdim*ijk+2]);
 163               */
 164             }
 165         }
 166     }
 167
 168   /*
 169     les elements:
 170   */
 171
 172   for(k=0;k<nelt_dir;k++)
 173     {
 174       for(j=0;j<nelt_dir;j++)
 175         {
 176           for (i=0;i<nelt_dir;i++)
 177             {
 178               int ijk = k*nelt_dir*nelt_dir+j*nelt_dir+i;
 179
 180               numhexa8[ijk] = ijk+1;
 181               nufahexa8[ijk] = 0;
 182
 183               indexN5 = k*nnoe_dir*nnoe_dir+j*nnoe_dir+i+1;
 184               indexN8 = indexN5+1;
 185               indexN1 = indexN5+nnoe_dir;
 186               indexN4 = indexN8+nnoe_dir;
 187
 188               indexN6 = indexN5+nnoe_dir*nnoe_dir;
 189               indexN7 = indexN8+nnoe_dir*nnoe_dir;
 190               indexN2 = indexN1+nnoe_dir*nnoe_dir;
 191               indexN3 = indexN4+nnoe_dir*nnoe_dir;
 192
 193               hexa8[8*ijk] = indexN1;
 194               hexa8[8*ijk+1] = indexN2;
 195               hexa8[8*ijk+2] = indexN3;
 196               hexa8[8*ijk+3] = indexN4;
 197               hexa8[8*ijk+4] = indexN5;
 198               hexa8[8*ijk+5] = indexN6;
 199               hexa8[8*ijk+6] = indexN7;
 200               hexa8[8*ijk+7] = indexN8;
 201
 202               /*
 203               printf("Connectivitée %d  i1 = %d  i2 = %d  i3 = %d  i4 = %d  i5 = %d  i6 = %d  i7 = %d  i8 = %d\n",(ijk+1),hexa8[8*ijk],hexa8[8*ijk+1],hexa8[8*ijk+2],hexa8[8*ijk+3],hexa8[8*ijk+4],hexa8[8*ijk+5],hexa8[8*ijk+6],hexa8[8*ijk+7]);
 204               */
 205             }
 206         }
 207     }
 208
 209   /*
 210     Les champs
 211   */
 212
 213   (void) time(&t1);
 214
 215    srand((int) t1); /* use time in seconds to set seed */
 216
 217    for(i=0;i<nnoe;i++)
 218      {
 219        DbleVectNode[mdim*i] = (med_float)
 220          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 221
 222        DbleVectNode[mdim*i+1] = (med_float)
 223          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 224
 225        DbleVectNode[mdim*i+2] = (med_float)
 226          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 227
 228        /*
 229          printf("i %d DbleVectNode %lf %lf\n",i,DbleVectNode[mdim*i],
 230          DbleVectNode[mdim*i+1],DbleVectNode[mdim*i+2]);
 231        */
 232      }
 233
 234    for(i=0;i<nhexa8;i++)
 235      {
 236        DbleVectCell[mdim*i] = (med_float)
 237          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 238
 239        DbleVectCell[mdim*i+1] = (med_float)
 240          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 241
 242        DbleVectCell[mdim*i+2] = (med_float)
 243          (1+(int) (100.0*rand()/(RAND_MAX+1.0)));
 244
 245        /*
 246          printf("i %d DbleVectCell %lf %lf\n",i,DbleVectCell[mdim*i],
 247          DbleVectCell[mdim*i+1],DbleVectCell[mdim*i+2]);
 248        */
 249      }
 250
 251   /***************************************************************************/
 252
 253   sprintf(buff,"%d",nnoe_dir);
 254   strcat(MedFile,buff);
 255   strcat(MedFile,".med");
 256
 257   fid = MEDouvrir(MedFile,RDWR);
 258
 259   if (fid < 0)
 260     ret = -1;
 261   else
 262     ret = 0;
 263   printf("%d\n",ret);
 264
 265   /***************************************************************************/
 266
 267   if (ret == 0)
 268     ret = MEDmaaCr(fid,maa,mdim);
 269   printf("%d\n",ret);
 270
 271   if (ret == 0)
 272     ret = MEDunvCr(fid,maa);
 273   printf("%d\n",ret);
 274
 275   /***************************************************************************/
 276
 277   if (ret == 0)
 278     ret = MEDnoeudsEcr(fid,maa,mdim,coo,MED_FULL_INTERLACE,MED_CART,
 279                        nomcoo,unicoo,nomnoe,MED_FAUX,numnoe,MED_VRAI,
 280                        nufano,nnoe,WRONLY);
 281   printf("%d\n",ret);
 282
 283   /*
 284     ecriture des mailles MED_HEXA8 :
 285     - connectivite
 286     - noms (optionnel)
 287     - numeros (optionnel)
 288     - numeros des familles
 289   */
 290
 291   if (ret == 0)
 292     ret = MEDelementsEcr(fid,maa,mdim,hexa8,MED_FULL_INTERLACE,
 293                          nomhexa8,MED_FAUX,numhexa8,MED_VRAI,nufahexa8,nhexa8,
 294                          MED_MAILLE,MED_HEXA8,MED_NOD,WRONLY);
 295   printf("%d \n",ret);
 296
 297   /***************************************************************************/
 298   /* ecriture des familles */
 299   /* Conventions :
 300      - toujours creer une famille de numero 0 ne comportant aucun attribut
 301        ni groupe (famille de reference pour les noeuds ou les elements
 302        qui ne sont rattaches a aucun groupe ni attribut)
 303      - les numeros de familles de noeuds sont > 0
 304      - les numeros de familles des elements sont < 0
 305      - rien d'imposer sur les noms de familles
 306    */
 307
 308   /* la famille 0 */
 309   if (ret == 0)
 310     {
 311       strcpy(nomfam,"FAMILLE_0");
 312       numfam = 0;
 313       ret = MEDfamCr(fid,maa,nomfam,numfam,&attide,&attval,attdes,0,
 314                      gro,0);
 315     }
 316   printf("%d \n",ret);
 317
 318   /***************************************************************************/
 319   /*
 320     Les Champs
 321   */
 322
 323   if (ret == 0)
 324     {
 325       ret = MEDchampCr(fid,champDbleVectNode,MED_REEL64,compDbleVectNode,
 326                        unitDbleVectNode,mdim);
 327
 328       printf("MEDchampCr DbleVectNode : %d \n",ret);
 329
 330       if (ret == 0)
 331         {
 332           ret = MEDchampEcr(fid, maa, champDbleVectNode,
 333                             (unsigned char *)DbleVectNode,
 334                             MED_NO_INTERLACE, nnoe,
 335                             MED_NOPG, MED_ALL, MED_NOPFL, WRONLY, MED_NOEUD,
 336                             0, MED_NOPDT,"        ", 0., MED_NONOR);
 337
 338           printf("MEDchampEcr DbleVectNode : %d \n",ret);
 339         }
 340     }
 341
 342   if (ret == 0)
 343     {
 344       ret = MEDchampCr(fid,champDbleVectCell,MED_REEL64,compDbleVectCell,
 345                        unitDbleVectCell,mdim);
 346
 347       printf("MEDchampCr DbleVectCell : %d \n",ret);
 348
 349       if (ret == 0)
 350         {
 351           ret = MEDchampEcr(fid, maa, champDbleVectCell,
 352                             (unsigned char *)DbleVectCell,
 353                             MED_NO_INTERLACE, nhexa8,
 354                             MED_NOPG, MED_ALL, MED_NOPFL, WRONLY, MED_MAILLE,
 355                             MED_HEXA8, MED_NOPDT,"        ", 0., MED_NONOR);
 356
 357           printf("MEDchampEcr DbleVectCell : %d \n",ret);
 358         }
 359     }
 360
 361   /***************************************************************************/
 362
 363   ret = MEDfermer(fid);
 364   printf("%d\n",ret);
 365
 366   free(coo);
 367   free(numnoe);
 368   free(nufano);
 369   free(nomnoe);
 370   free(hexa8);
 371   free(numhexa8);
 372   free(nufahexa8);
 373   free(nomhexa8);
 374   free(DbleVectNode);
 375   free(DbleVectCell);
 376
 377   return 0;
 378 }
 379
 380