1 /* ----------------------------------------------------------------  */
2 /* readFoam                                        7.07.2005 Wittig  */
3 /* ----------------------------------------------------------------  */
4 
5 
6 #include <cgx.h>
7 #include <dirent.h>
8 
9 #define kappa    1.4
10 
11 #define TEST     0
12 
13 #define MAX_ENTITIES      100
14 #define CHECK_INPUT       0
15 
16 extern Sets      *setx;
17 extern Summen    *anzx;
18 extern Alias     *alias;
19 extern SumGeo    anzGeo[1];
20 extern SumAsci   sumAsci[1];
21 
22 
23 typedef struct {
24   int anz_f;
25   int *face;
26   int *ori;
27 }Elemf;
28 
29 typedef struct {
30   int p;
31   int u;
32   int t;
33   int rho;
34 }Ds;
35 
36 
get_values_in_string(char * handle,double * fval)37 int get_values_in_string( char *handle, double *fval)
38 {
39   int i=0, pos=0;
40   static int parenthesis=0;
41   char byte;
42   char string[MAX_LINE_LENGTH];
43 
44 
45   /* count open and closed parenthesis */
46   do
47   {
48     byte = handle[pos++];
49     if(byte=='(') parenthesis++;
50     else if(byte==')')
51     {
52       parenthesis--;
53 
54       /* return if the first open '(' is closed */
55       if(parenthesis==0) break;
56 
57       /* if one "(" is open get the values and return */
58       if(parenthesis==1)
59       {
60         string[i++]=0;
61         sscanf(string, "%lf %lf %lf", &fval[0], &fval[1], &fval[2]);
62         return(pos);
63       }
64     }
65     else if (parenthesis==2)
66     {
67       if(i<MAX_LINE_LENGTH-2) string[i++]=byte;
68     }
69 
70   }while(byte!=0);
71   fval[0]=fval[1]=fval[2]=0.;
72   return(-parenthesis);
73 }
74 
defineElemNodes(Summen * anzx,Faces * face,Elements * elem,int i,int * eface,int anz_f)75 int defineElemNodes(Summen *anzx, Faces *face, Elements *elem, int i, int *eface, int anz_f)
76 {
77   int f, l, n, j, k, ni[8];
78   int fbuf[6], slavnod[2];
79 
80   if (anz_f==6) /* HEXA8  */
81   {
82     /* mark the faces as used */
83     for(f=0; f<anz_f; f++) face[eface[f]].group++;
84 
85     elem[i].type = 1;
86     anzx->etype[elem[i].type]++;
87 
88     /* bottom */
89     /* owner */
90     if(face[eface[0]].group==1)
91     {
92       elem[i].nod[0] = face[eface[0]].nod[0];
93       elem[i].nod[4] = face[eface[0]].nod[1];
94       elem[i].nod[7] = face[eface[0]].nod[2];
95       elem[i].nod[3] = face[eface[0]].nod[3];
96     }
97     /* neighbour */
98     else if(face[eface[0]].group==2)
99     {
100       elem[i].nod[0] = face[eface[0]].nod[0];
101       elem[i].nod[4] = face[eface[0]].nod[3];
102       elem[i].nod[7] = face[eface[0]].nod[2];
103       elem[i].nod[3] = face[eface[0]].nod[1];
104     }
105     else { printf("ERROR in group:%d\n",face[eface[0]].group); exit(1); }
106 
107     /* search face on node 0 and 1 from face0 */
108     for(f=2; f<anz_f; f++)
109     {
110       for(l=0; l<4; l++)
111       {
112         if(face[eface[f]].nod[l]==face[eface[0]].nod[0])
113         {
114           slavnod[0]=l;
115           goto found_face;
116         }
117       }
118     }
119     printf("ERROR: found no slavface for elem:%d\n",i ); exit(1);
120   found_face:;
121     for(l=0; l<4; l++)
122     {
123       for(n=1; n<4; n++)
124       {
125         if(face[eface[f]].nod[l]==face[eface[0]].nod[n])
126         {
127           if((l==0)&&(slavnod[0]==1)) slavnod[1]=2;
128         else if((l==1)&&(slavnod[0]==2)) slavnod[1]=3;
129         else if((l==2)&&(slavnod[0]==3)) slavnod[1]=0;
130         else if((l==3)&&(slavnod[0]==0)) slavnod[1]=1;
131 
132         else if((l==2)&&(slavnod[0]==1)) slavnod[1]=0;
133         else if((l==3)&&(slavnod[0]==2)) slavnod[1]=1;
134         else if((l==0)&&(slavnod[0]==3)) slavnod[1]=2;
135         else if((l==1)&&(slavnod[0]==0)) slavnod[1]=3;
136         else { printf("ERROR: slavnod nod found l:%d slavnod[0]:%d\n",l,slavnod[0]); exit(1); }
137           goto found_nodn;
138         }
139       }
140     }
141     printf("ERROR: found no slavnod[1] for elem:%d\n",i ); exit(1);
142   found_nodn:;
143 
144     /* determine the node from face1 at the slavnod[1] */
145     for(l=0; l<4; l++)
146     {
147       if(face[eface[f]].nod[slavnod[1]]==face[eface[1]].nod[l])
148       {
149         /* top */
150         if(face[eface[1]].group==1)
151         {
152           if(l==0)
153           {
154             elem[i].nod[1] = face[eface[1]].nod[0];
155             elem[i].nod[2] = face[eface[1]].nod[1];
156             elem[i].nod[6] = face[eface[1]].nod[2];
157             elem[i].nod[5] = face[eface[1]].nod[3];
158           }
159           else if(l==1)
160           {
161             elem[i].nod[1] = face[eface[1]].nod[1];
162             elem[i].nod[2] = face[eface[1]].nod[2];
163             elem[i].nod[6] = face[eface[1]].nod[3];
164             elem[i].nod[5] = face[eface[1]].nod[0];
165           }
166           else if(l==2)
167           {
168             elem[i].nod[1] = face[eface[1]].nod[2];
169             elem[i].nod[2] = face[eface[1]].nod[3];
170             elem[i].nod[6] = face[eface[1]].nod[0];
171             elem[i].nod[5] = face[eface[1]].nod[1];
172           }
173           else if(l==3)
174           {
175             elem[i].nod[1] = face[eface[1]].nod[3];
176             elem[i].nod[2] = face[eface[1]].nod[0];
177             elem[i].nod[6] = face[eface[1]].nod[1];
178             elem[i].nod[5] = face[eface[1]].nod[2];
179           }
180           else { printf("ERROR: l nod found\n"); exit(1); }
181         }
182         else if(face[eface[1]].group==2)
183         {
184           if(l==0)
185           {
186             elem[i].nod[1] = face[eface[1]].nod[0];
187             elem[i].nod[2] = face[eface[1]].nod[3];
188             elem[i].nod[6] = face[eface[1]].nod[2];
189             elem[i].nod[5] = face[eface[1]].nod[1];
190           }
191           else if(l==1)
192           {
193             elem[i].nod[1] = face[eface[1]].nod[1];
194             elem[i].nod[2] = face[eface[1]].nod[0];
195             elem[i].nod[6] = face[eface[1]].nod[3];
196             elem[i].nod[5] = face[eface[1]].nod[2];
197           }
198           else if(l==2)
199           {
200             elem[i].nod[1] = face[eface[1]].nod[2];
201             elem[i].nod[2] = face[eface[1]].nod[1];
202             elem[i].nod[6] = face[eface[1]].nod[0];
203             elem[i].nod[5] = face[eface[1]].nod[3];
204           }
205           else if(l==3)
206           {
207             elem[i].nod[1] = face[eface[1]].nod[3];
208             elem[i].nod[2] = face[eface[1]].nod[2];
209             elem[i].nod[6] = face[eface[1]].nod[1];
210             elem[i].nod[5] = face[eface[1]].nod[0];
211           }
212           else { printf("ERROR: l nod found\n"); exit(1); }
213         }
214         else { printf("ERROR in group:%d\n",face[eface[1]].group); exit(1); }
215         break;
216       }
217     }
218 
219     for(f=0; f<anz_f; f++)
220     {
221       face[eface[f]].elem_nr=elem[i].nr;
222 
223       n=0;
224       for(j=0; j<8; j++)
225         for(k=0; k<4; k++)
226           if(elem[i].nod[j]==face[eface[f]].nod[k])
227             ni[n++]=j+1;
228       if (n==4)
229       {
230         /* check which sides are involved */
231         for (j=0; j<6; j++) fbuf[j]=0;
232         for (j=0; j<4; j++)
233         {
234           if ((ni[j]==1)||(ni[j]==2)||(ni[j]==3)||(ni[j]==4)) fbuf[0]++;
235           if ((ni[j]==5)||(ni[j]==8)||(ni[j]==7)||(ni[j]==6)) fbuf[1]++;
236             if ((ni[j]==1)||(ni[j]==5)||(ni[j]==6)||(ni[j]==2)) fbuf[2]++;
237           if ((ni[j]==2)||(ni[j]==6)||(ni[j]==7)||(ni[j]==3)) fbuf[3]++;
238           if ((ni[j]==3)||(ni[j]==7)||(ni[j]==8)||(ni[j]==4)) fbuf[4]++;
239           if ((ni[j]==4)||(ni[j]==8)||(ni[j]==5)||(ni[j]==1)) fbuf[5]++;
240         }
241       }
242       for (j=0; j<6; j++)
243       {
244         /* printf(" j:%d fbuf[j]:%d\n", j, fbuf[j]); */
245         if(fbuf[j]==4)
246         {
247           face[eface[f]].nr=j;
248           if( anzx->emax>anzx->f)
249           {
250             if ( (face = (Faces *)realloc((Faces *)face, (anzx->emax+2) * sizeof(Faces))) == NULL )
251             {  printf("\n\n ERROR: malloc failed, face\n\n") ; exit(1); }
252 	  }
253           face[elem[i].nr].indx[j]=eface[f];
254 	}
255       }
256     }
257   }
258   else if (anz_f==5) /* PE6    */
259   {
260     elem[i].type = 2; printf("elem-type:pe6 not installed\n"); return(0);
261   }
262   else if (anz_f==4) /* TET4   */
263   {
264     elem[i].type = 3; printf("elem-type:tet4 not installed\n"); return(0);
265   }
266   else
267   {
268     printf("\n\n ERROR: number of faces:%d of element:%d do not match available types of elements\n", anz_f, i);
269     return(0);
270   }
271 
272 #if TEST
273           printf (" e:%d t:%d f:",elem[i].nr, elem[i].type);
274           for(j=0; j<anz_f; j++)  printf (" %d",eface[j]);
275           printf ("\nnode:");
276           for(j=0; j<8; j++)  printf (" %d",elem[i].nod[j]);
277           printf ("\n");
278 #endif
279    return(1);
280 }
281 
282 
283 
284 
285 
readFoam(char * datin,Summen * apre,Sets ** sptr,Nodes ** nptr,Elements ** eptr,Datasets ** lptr)286 int readFoam(char *datin, Summen *apre, Sets **sptr, Nodes **nptr, Elements **eptr, Datasets **lptr  )
287 {
288   int i,j,k=0,l,n,f,f1=0,f2, anz_f, length, pos;
289   int eface[8], sum_dir, ncomps;
290   FILE *handle;
291   DIR *dirp;
292   struct dirent *dp;
293   Rsort *rsort=NULL;
294   Ds ds;
295 
296   double cp, Rg;
297   double ps, us, ts, pt, tt, rho, ma;
298 
299 
300   char file[MAX_LINE_LENGTH], rec_str[MAX_LINE_LENGTH];
301   char dat[2][MAX_LINE_LENGTH];
302   char resultname[MAX_LINE_LENGTH];
303   char buffer[MAX_LINE_LENGTH];
304   char **resultdir=NULL;
305   char **resultdir2=NULL;
306 
307   Nodes     *node = NULL;
308   Faces     *face = NULL;
309   Elements  *elem = NULL;
310   Datasets  *lcase= NULL;
311   Elemf     *elemf= NULL;
312 
313   int anz_dat=0, ipuf=0, *vpn=NULL;
314   double  fval[10];
315 
316   int startFace, nFaces, setNr=-1, formFlag=0;
317 
318   anzx=apre;
319   setx=*sptr;
320 
321   anzx->n=-1;
322   anzx->f=-1;
323   anzx->e=-1;
324   anzx->l=-1;
325 
326   /* Open the files and check to see that it was opened correctly */
327 
328   /* mesh */
329   strcpy(anzx->model, "foam");
330   printf (" MODEL NAME:  %s\n", anzx->model);
331 
332   anzx->n=anzx->f=anzx->e=0;
333   anzx->nmax=-MAX_INTEGER;  anzx->nmin=MAX_INTEGER;
334   anzx->emax=-MAX_INTEGER;  anzx->emin=MAX_INTEGER;
335 
336   printf ("\n Read the Mesh  \n");
337 
338   /* nodes */
339   /* located in "datin"/constant/polymesh/points */
340   sprintf( file, "%s/constant/polyMesh/points", datin);
341   handle = fopen (file, "r");
342   if ( handle== NULL )
343   {
344     printf ("\nThe input file \"%s\" could not be opened.\n", file);
345   }
346   else printf ("\n%s opened\n",file);
347   while(1)
348   {
349     length = frecord( handle, rec_str);
350     if (rec_str[length] == (char)EOF) break;
351 
352     /* read data if a block starts */
353     if(rec_str[0]=='(')
354     {
355       if ( (node = (Nodes *)malloc( (anzx->n+1) * sizeof(Nodes))) == NULL )
356       {  printf("\n\n ERROR: malloc failed, node\n\n") ; exit(1); }
357       for(i=0; i<anzx->n; i++)
358       {
359         node[i].nr=i+1;
360         node[node[i].nr].indx=i;
361         length = frecord( handle, rec_str);
362         sscanf(&rec_str[1], "%lf %lf %lf", &node[node[i].nr].nx, &node[node[i].nr].ny, &node[node[i].nr].nz);
363 #if TEST
364         printf (" n=%d x=%lf y=%lf z=%lf \n", node[i].nr, node[node[i].nr].nx, node[node[i].nr].ny, node[node[i].nr].nz);
365 #endif
366       }
367       break;
368     }
369     else if(length>0) anzx->n=atoi(rec_str);
370   }
371   fclose(handle);
372   anzx->nmax=anzx->n;
373   anzx->nmin=1;
374 
375   /* faces */
376   /* located in "datin"/constant/polymesh/faces */
377   sprintf( file, "%s/constant/polyMesh/faces", datin);
378   handle = fopen (file, "r");
379   if ( handle== NULL )
380   {
381     printf ("\nThe input file \"%s\" could not be opened.\n", file);
382   }
383   else printf ("\n%s opened\n",file);
384   while(1)
385   {
386     length = frecord( handle, rec_str);
387     if (rec_str[length] == (char)EOF) break;
388 
389     /* read data if a block starts */
390     if(rec_str[0]=='(')
391     {
392       if ( (face = (Faces *)malloc( (anzx->f+2) * sizeof(Faces))) == NULL )
393       {  printf("\n\n ERROR: malloc failed, face\n\n") ; exit(1); }
394       for(i=0; i<anzx->f; i++)
395       {
396         face[i].group=0;
397         face[i].mat=-1;
398         face[i].side=NULL;
399         length = frecord( handle, rec_str);
400         sscanf(rec_str, "%d",&n);
401         sscanf(&rec_str[2], "%d %d %d %d %d %d %d %d",&face[i].nod[0],&face[i].nod[1],&face[i].nod[2],&face[i].nod[3],&face[i].nod[4],&face[i].nod[5],&face[i].nod[6],&face[i].nod[7]);
402         for(j=0; j<n; j++)  face[i].nod[j]++; /* nodenr is one higher! */
403 
404         if(n==3)      face[i].type = 7;  /* TRI3  */
405         else if(n==6) face[i].type = 8;  /* TRI6  */
406         else if(n==4) face[i].type = 9;  /* QUAD4 */
407         else if(n==8) face[i].type = 10; /* QUAD8 */
408         else
409         {
410           printf("\n\n ERROR: number of nodes:%d of face:%d do not match available types of elements\n",n,i);
411           exit(1);
412         }
413 
414 #if TEST
415         printf (" f:%d t:%d",i, face[i].type);
416         for(j=0; j<n; j++)  printf (" %d",face[i].nod[j]);
417         printf ("\n");
418 #endif
419       }
420       break;
421     }
422     else if(length>0) anzx->f=atoi(rec_str);
423   }
424   fclose(handle);
425 
426   /* elements */
427   /* located in "datin"/constant/polymesh/cells */
428   sprintf( file, "%s/constant/polyMesh/cells", datin);
429   handle = fopen (file, "r");
430   if ( handle== NULL )
431   {
432     //printf ("\nThe input file \"%s\" could not be opened.\n", file);
433 
434     /* no cells found, try to read owner and neighbour */
435     sprintf( file, "%s/constant/polyMesh/owner", datin);
436     handle = fopen (file, "r");
437     if ( handle== NULL )
438     {
439       printf("ERROR: found no owner file\n"); exit(1);
440     }
441 
442     /* create a list of faces per element */
443     printf ("\n%s opened\n",file);
444     while(1)
445     {
446       length = frecord( handle, rec_str);
447       if (rec_str[length] == (char)EOF) break;
448 
449       if(rec_str[0]=='(') break;
450     }
451     /* read data if a block starts */
452     i=0;
453     while(1)
454     {
455       length = frecord( handle, rec_str);
456       if (rec_str[length] == (char)EOF) break;
457       if(rec_str[0]==')') break;
458       sscanf(rec_str, "%d",&face[i].elem_nr);
459       if(face[i].elem_nr>anzx->e) anzx->e=face[i].elem_nr;
460       i++;
461     }
462     fclose(handle);
463     anzx->e++;
464 
465     sprintf( file, "%s/constant/polyMesh/neighbour", datin);
466     handle = fopen (file, "r");
467     if ( handle== NULL )
468     {
469       printf("ERROR: found no owner file\n"); exit(1);
470     }
471 
472     /* create a list of faces per element */
473     printf ("\n%s opened\n",file);
474     while(1)
475     {
476       length = frecord( handle, rec_str);
477       if (rec_str[length] == (char)EOF) break;
478 
479       if(rec_str[0]=='(') break;
480     }
481     /* read data if a block starts */
482     i=0;
483     while(1)
484     {
485       length = frecord( handle, rec_str);
486       if (rec_str[length] == (char)EOF) break;
487       if(rec_str[0]==')') break;
488       sscanf(rec_str, "%d",&face[i++].mat);
489     }
490     fclose(handle);
491 
492     /* define the cells */
493     if ( (elemf = (Elemf *)malloc( (anzx->e+1) * sizeof(Elemf))) == NULL )
494     {  printf("\n\n ERROR: malloc failed, elemf\n\n") ; exit(1); }
495     for(i=0; i<anzx->e; i++) elemf[i].face=NULL;
496     for(i=0; i<anzx->e; i++) elemf[i].anz_f=0;
497     for(i=0; i<anzx->f; i++)
498     {
499       if ( (elemf[face[i].elem_nr].face = (int *)realloc( (int *)elemf[face[i].elem_nr].face, (elemf[face[i].elem_nr].anz_f+1) * sizeof(int))) == NULL )
500       {  printf("\n\n ERROR: malloc failed, elemf\n\n") ; exit(1); }
501 
502       elemf[face[i].elem_nr].face[elemf[face[i].elem_nr].anz_f]=i;
503       elemf[face[i].elem_nr].anz_f++;
504       if(face[i].mat>-1)
505       {
506         if ( (elemf[face[i].mat].face = (int *)realloc( (int *)elemf[face[i].mat].face, (elemf[face[i].mat].anz_f+1) * sizeof(int))) == NULL )
507         {  printf("\n\n ERROR: malloc failed, elemf\n\n") ; exit(1); }
508         elemf[face[i].mat].face[elemf[face[i].mat].anz_f]=i;
509         elemf[face[i].mat].anz_f++;
510       }
511     }
512 
513     /* define the elements */
514     if ( (elem = (Elements *)malloc( (anzx->e+1) * sizeof(Elements))) == NULL )
515     {  printf("\n\n ERROR: malloc failed, elem\n\n") ; exit(1); }
516     for(i=0; i<anzx->e; i++)
517     {
518       /* search two disjunct faces */
519       eface[0]=elemf[i].face[0];
520       f2=2;
521       if     (face[eface[0]].type == 7) f1=3;  /* TRI3  */
522       else if(face[eface[0]].type == 8) f1=6;  /* TRI6  */
523       else if(face[eface[0]].type == 9) f1=4;  /* QUAD4 */
524       else if(face[eface[0]].type == 10) f1=8; /* QUAD8 */
525       for(f=1; f<elemf[i].anz_f; f++)
526       {
527 	if     (face[elemf[i].face[f]].type == 7) k=3;  /* TRI3  */
528         else if(face[elemf[i].face[f]].type == 8) k=6;  /* TRI6  */
529         else if(face[elemf[i].face[f]].type == 9) k=4;  /* QUAD4 */
530         else if(face[elemf[i].face[f]].type == 10) k=8; /* QUAD8 */
531         for(n=0; n<k; n++)
532 	{
533 	  for(j=0; j<f1; j++)
534           {
535 	    if(face[elemf[i].face[f]].nod[n]==face[eface[0]].nod[j])
536             {
537               eface[f2++]=elemf[i].face[f];
538               goto isconnected;
539             }
540           }
541         }
542         eface[1]=elemf[i].face[f];
543       isconnected:;
544       }
545 #if TEST
546       printf("face order (%d): ", elemf[i].anz_f);
547       for(f=0; f<elemf[i].anz_f; f++) printf(" %d", eface[f]);
548       printf("\n");
549 #endif
550 
551       elem[i].nr = i+1;
552       elem[i].group= 0;
553       elem[i].mat = 1;
554 
555       defineElemNodes(anzx, face, elem, i, eface, elemf[i].anz_f);
556     }
557   }
558   else
559   {
560     printf ("\n%s opened\n",file);
561     while(1)
562     {
563       length = frecord( handle, rec_str);
564       if (rec_str[length] == (char)EOF) break;
565 
566       /* read data if a block starts */
567       if(rec_str[0]=='(')
568       {
569         if ( (elem = (Elements *)malloc( (anzx->e+1) * sizeof(Elements))) == NULL )
570         {  printf("\n\n ERROR: malloc failed, elem\n\n") ; exit(1); }
571         for(i=0; i<anzx->e; i++)
572         {
573     	  elem[i].nr = i+1;
574           elem[i].group= 0;
575           elem[i].mat = 1;
576           length = frecord( handle, rec_str);
577           sscanf(rec_str, "%d",&anz_f);
578           sscanf(&rec_str[2], "%d %d %d %d %d %d",&eface[0],&eface[1],&eface[2],&eface[3],&eface[4],&eface[5]);
579           defineElemNodes(anzx, face, elem, i, eface, anz_f);
580         }
581         break;
582       }
583       else if(length>0) anzx->e=atoi(rec_str);
584     }
585     fclose(handle);
586   }
587   anzx->emax=anzx->e;
588   anzx->emin=1;
589 
590 
591   /* boundaries */
592   /* located in "datin"/constant/polymesh/boundary */
593   sprintf( file, "%s/constant/polyMesh/boundary", datin);
594   handle = fopen (file, "r");
595   if ( handle== NULL )
596   {
597     printf ("\nThe input file \"%s\" could not be opened.\n", file);
598   }
599   else printf ("\n%s opened\n",file);
600   while(1)
601   {
602     length = frecord( handle, rec_str);
603     if (rec_str[length] == (char)EOF) break;
604 
605     /* read data if a block starts */
606     if(rec_str[0]=='(')
607     {
608       while(1)
609       {
610         length = frecord( handle, rec_str);
611         if (rec_str[length] == (char)EOF) break;
612         length = strfind(rec_str, ")");
613         if ( length!=-1 ) break;
614 
615         length = strfind(rec_str, "{");
616         if ( length!=-1 )
617         {
618           printf("\n   setname:%s\n",  buffer);
619           setNr=getSetNrx(buffer);
620           if(setNr<0) setNr=pre_setax( buffer, "is", 0);
621         }
622         else sscanf(rec_str, "%s", buffer);
623 
624         length = strfind(rec_str, "nFaces");
625         if ( length!=-1 )
626         {
627           sscanf(rec_str,"%*s %d", &nFaces);
628           printf("   nFaces:%d\n", nFaces );
629         }
630 
631         length = strfind(rec_str, "startFace");
632         if ( length!=-1 )
633         {
634           sscanf(rec_str,"%*s %d", &startFace);
635           printf("   startFace:%d\n", startFace );
636         }
637 
638         /* store the faces in the set if the defining block was closed */
639         length = strfind(rec_str, "}");
640         if ( length!=-1 )
641         {
642           for(j=0; j<nFaces; j++)
643           {
644             if(face[j+startFace].nr>-1)
645             {
646               setax(setNr,"f", j+startFace);
647               i=setax(setNr,"j",0);
648               if(i>-1)
649               {
650                 setx[setNr].elf[i].e=face[j+startFace].elem_nr;
651                 setx[setNr].elf[i].f=face[j+startFace].nr;
652               }
653 	    }
654           }
655         }
656       }
657     }
658   }
659   fclose(handle);
660 
661   /* scan over all project directories and search for results */
662   sum_dir=0;
663   dirp = opendir(datin);
664   printf("seach for results in %s\n", datin );
665   if (dirp != NULL) while ((dp = readdir(dirp)) != NULL)
666   {
667 
668     /* if the dir-name starts with a number, its a result dir */
669     if((dp->d_name[0]>=48)&&(dp->d_name[0]<=57))
670     {
671       if ( (resultdir2 = (char **)realloc( resultdir2, (sum_dir+2) * sizeof(char *))) == NULL )
672         printf("\n\n ERROR: realloc failure\n\n" );
673       if ( (resultdir2[sum_dir] = (char *)malloc( (MAX_LINE_LENGTH) * sizeof(char))) == NULL )
674         printf("\n\n ERROR: realloc failure\n\n" );
675       sprintf( resultdir2[sum_dir], "%s/%s", datin, dp->d_name );
676 
677       if ( (rsort = (Rsort *)realloc( rsort, (sum_dir+1) * sizeof(Rsort))) == NULL )
678         printf("ERROR: realloc failed: Rsort\n\n" );
679       rsort[sum_dir].r=atof(dp->d_name);
680       rsort[sum_dir].i=sum_dir;
681       sum_dir++;
682     }
683   }
684   closedir(dirp);
685 
686   /* sort the resultdirs according to their name-value (its the time-step value) */
687   qsort( rsort, sum_dir, sizeof(Rsort), (void *)compareRsort );
688   if ( (resultdir = (char **)realloc( resultdir, (sum_dir+1) * sizeof(char *))) == NULL )
689     printf("\n\n ERROR: realloc failure\n\n" );
690   for( i=0; i<sum_dir; i++)
691   {
692     if ( (resultdir[i] = (char *)malloc( (MAX_LINE_LENGTH) * sizeof(char))) == NULL )
693       printf("\n\n ERROR: realloc failure\n\n" );
694     strcpy( resultdir[i], resultdir2[rsort[i].i]);
695     printf(" found results in dir:%s\n",resultdir[i]);
696   }
697 
698 
699   if ( (vpn = (int *)malloc( (anzx->nmax+1) * sizeof(int))) == NULL )
700     printf("\n\n ERROR: malloc failure\n\n" );
701   for(i=0; i<=anzx->nmax; i++) vpn[i]=0;
702 
703   for(l=0; l<sum_dir; l++)
704   {
705     /* search for files in dir */
706     dirp = opendir(resultdir[l]);
707     if (dirp != NULL)
708     {
709      ds.p=ds.t=ds.u=ds.rho=-1;
710 
711      /* scan all files */
712      while ((dp = readdir(dirp)) != NULL)
713      {
714       if ((strlen(dp->d_name) >0 )&&(dp->d_name[strlen(dp->d_name)-1]!='~')&&(dp->d_name[0]!='.'))
715       {
716         /* open the dir/file */
717         sprintf(file, "%s/%s",resultdir[l],dp->d_name);
718         handle = fopen (file, "r");
719         if ( handle== NULL )
720         {
721           printf ("\nThe input file \"%s\" could not be opened.\n", file);
722           goto noFile;
723         }
724         else printf("file:%s opened\n", file);
725 
726         /* determine the type of data and create the entity names */
727         ncomps=0;
728         while(1)
729         {
730           length = frecord( handle, rec_str);
731           // printf ("record: %s", rec_str);
732           if (rec_str[length] == (char)EOF) break;
733 
734           /* read data when the header starts */
735           else if(rec_str[0]=='{')
736           {
737             while(1)
738             {
739               length = frecord( handle, rec_str);
740               /* printf ("record: %s", rec_str); */
741               if (rec_str[length] == (char)EOF) { printf("ERROR\n"); exit(1); }
742               if (rec_str[0] == '}') goto foundEntities;
743 
744               sscanf(rec_str, "%s %s", dat[0], dat[1]);
745               // printf ("dat:%s| dat:%s|\n",dat[0], dat[1]);
746 
747               if(compareStrings(dat[0], "class")>0)
748               {
749                 if( compare(dat[1], "volScalarField", 14)==14) ncomps=1;
750                 else if( compare(dat[1], "volVectorField", 14)==14) ncomps=3;
751                 else { printf("Warning; class:%s not known\n",dat[1]); ncomps=0; }
752               }
753               if(compareStrings(dat[0], "object")>0)
754               { strcpy(resultname,dat[1]); resultname[strlen(dat[1])]=0; }
755             }
756           }
757         } foundEntities:;
758         printf ("  resultname:%s entities:%d\n", resultname, ncomps);
759 
760         if(ncomps<1) goto foundNoEntities;
761 
762         /* define the lcases */
763         anzx->l++;
764 
765         /* store Datasets */
766         if ( (lcase = (Datasets *)realloc(lcase, (anzx->l+1) * sizeof(Datasets))) == NULL )
767           printf("\n\n ERROR: malloc failed, lcase\n\n") ;
768 
769         i=strlen(resultdir[l]); while((i)&&(resultdir[l][i]!='/')) i--;
770         lcase[anzx->l].value=atof(&resultdir[l][i+1]);
771         if(compareStrings(&resultdir[l][i+1],"0")==1) sprintf( lcase[anzx->l].dataset_name,"%s", "initial");
772         else sprintf( lcase[anzx->l].dataset_name,"%s    ", "result");
773 
774         lcase[anzx->l].irtype = 1;
775         lcase[anzx->l].npheader = 0;
776         strcpy(lcase[anzx->l].analysis_name,"");
777         strcpy(lcase[anzx->l].name, dp->d_name );
778         strcpy(lcase[anzx->l].dataset_text,"");
779         lcase[anzx->l].step_number=l+1;
780         lcase[anzx->l].analysis_type=1;
781 
782         printf ("lcase.name[%d]= %s value:%f\n", anzx->l, lcase[anzx->l].name, lcase[anzx->l].value);
783         if(lcase[anzx->l].name[0]=='p') ds.p=anzx->l;
784         if(lcase[anzx->l].name[0]=='U') ds.u=anzx->l;
785         if(lcase[anzx->l].name[0]=='T') ds.t=anzx->l;
786         if(lcase[anzx->l].name[0]=='r') ds.rho=anzx->l;
787         lcase[anzx->l].ncomps  = ncomps;
788 
789 
790         if ( (lcase[anzx->l].nmax = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
791           printf("\n\n ERROR: malloc failure\n\n" );
792         if ( (lcase[anzx->l].nmin = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
793           printf("\n\n ERROR: malloc failure\n\n" );
794         if ( (lcase[anzx->l].max = (float *)malloc( lcase[anzx->l].ncomps * sizeof(float))) == NULL )
795           printf("\n\n ERROR: malloc failure\n\n" );
796         if ( (lcase[anzx->l].min = (float *)malloc( lcase[anzx->l].ncomps * sizeof(float))) == NULL )
797           printf("\n\n ERROR: malloc failure\n\n" );
798         if ( (lcase[anzx->l].dat = (float **)malloc( lcase[anzx->l].ncomps * sizeof(float *))) == NULL )
799           printf("\n\n ERROR: malloc failure\n\n" );
800         if ( (lcase[anzx->l].icname = (char **)malloc( lcase[anzx->l].ncomps * sizeof(char *))) == NULL )
801           printf("\n\n ERROR: malloc failure\n\n" );
802         if ( (lcase[anzx->l].compName = (char **)malloc( lcase[anzx->l].ncomps * sizeof(char *))) == NULL )
803           printf("\n\n ERROR: malloc failure\n\n" );
804         for(i=0; i<lcase[anzx->l].ncomps; i++)
805         {
806           if ( (lcase[anzx->l].compName[i] = (char *)malloc( MAX_LINE_LENGTH * sizeof(char))) == NULL )
807             printf("\n\n ERROR: malloc failed\n\n" );
808           if ( (lcase[anzx->l].dat[i] = (float *)malloc( (anzx->nmax+1) * sizeof(float))) == NULL )
809             printf("\n\n ERROR: malloc failed\n\n" );
810           if ( (lcase[anzx->l].icname[i] = (char *)malloc( MAX_LINE_LENGTH * sizeof(char))) == NULL )
811              printf("\n\n ERROR: malloc failed\n\n" );
812           lcase[anzx->l].max[i]=-MAX_FLOAT;
813           lcase[anzx->l].min[i]=MAX_FLOAT;
814         }
815         if ( (lcase[anzx->l].menu = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
816           printf("\n\n ERROR: malloc failure\n\n" );
817         if ( (lcase[anzx->l].ictype = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
818           printf("\n\n ERROR: malloc failure\n\n" );
819         if ( (lcase[anzx->l].icind1 = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
820           printf("\n\n ERROR: malloc failure\n\n" );
821         if ( (lcase[anzx->l].icind2 = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
822           printf("\n\n ERROR: malloc failure\n\n" );
823         if ( (lcase[anzx->l].iexist = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
824           printf("\n\n ERROR: malloc failure\n\n" );
825 
826         for(i=0; i<lcase[anzx->l].ncomps; i++)
827         {
828           sprintf(lcase[anzx->l].compName[i],"%s%d ",resultname, i+1);
829           lcase[anzx->l].menu[i] = 1;
830           if(ncomps==1)
831 	  {
832             lcase[anzx->l].ictype[i] = 1;
833             lcase[anzx->l].icind1[i] = 0;
834 	  }
835           else if(ncomps==3)
836           {
837             lcase[anzx->l].ictype[i]=2;
838             lcase[anzx->l].icind1[i]=i+1;
839           }
840           lcase[anzx->l].icind2[i] = 0;
841           lcase[anzx->l].iexist[i] = 0;
842           for(j=0; j<=anzx->nmax; j++) lcase[anzx->l].dat[i][j]=0.;
843         }
844 
845         while(1)
846         {
847           length = frecord( handle, rec_str);
848           if (rec_str[length] == (char)EOF) break;
849 
850           /* read "internalField" data */
851           if ( strfind(rec_str, "internalField") !=-1 )
852           {
853             sscanf(rec_str, "%*s %s", buffer);
854             if(compareStrings(buffer,"uniform")>0)
855             {
856               formFlag=1;
857               if     (ncomps==1) sscanf(rec_str, "%*s %*s %lf", &fval[0]);
858               else if(ncomps==3)
859               {
860                 length = strfind(rec_str, "(");
861                  sscanf(&rec_str[length+1], "%lf %lf %lf", &fval[0], &fval[1], &fval[2]);
862 	      }
863               else { printf("ERROR1; class not known\n"); exit(1); }
864 		//printf("%s vals:%lf %lf %lf\n", &rec_str[length], fval[0], fval[1], fval[2]);
865 	    }
866             else
867             {
868               formFlag=0;
869               length = frecord( handle, rec_str);
870               anz_dat=atoi(rec_str);
871               if(anzx->e!=anz_dat)
872               { printf("ERROR: sum-data:%d do not match sum-elem:%d\n",anz_dat,anzx->e); exit(1); }
873               length = frecord( handle, rec_str);
874 	    }
875 
876             for(i=0; i<anzx->e; i++)
877             {
878 	      if(formFlag==0)
879 	      {
880                 length = frecord( handle, rec_str);
881                 if     (ncomps==1) sscanf(rec_str, "%lf", &fval[0]);
882                 else if(ncomps==3)
883                 {
884                   sscanf(&rec_str[1], "%lf %lf %lf", &fval[0], &fval[1], &fval[2]);
885 	        }
886                 else { printf("ERROR2; class not known (ncomps:%d)\n",ncomps); exit(1); }
887 	      }
888               if (elem[i].type == 1) ipuf = 8;  /* HEXA8 */
889               if (elem[i].type == 2) ipuf = 6;  /* PENTA6 */
890               if (elem[i].type == 3) ipuf = 4;  /* TET4 */
891               for(n=0; n<lcase[anzx->l].ncomps; n++)
892               {
893                 for(j=0; j<ipuf; j++)
894                 {
895                   lcase[anzx->l].dat[n][elem[i].nod[j]]+=fval[n];
896                   if((vpn[0]==0)&&(n==0))
897                   {
898                     vpn[elem[i].nod[j]]++; /* count the hits */
899 	          }
900                 }
901               }
902             }
903 	    vpn[0]=1;
904           }
905 
906           /* read "boundaryField" data */
907           if ( strfind(rec_str, "boundaryField") !=-1 )
908           {
909             length = frecord( handle, rec_str);
910             while(1)
911             {
912               length = frecord( handle, rec_str);
913               if (rec_str[length] == (char)EOF) break;
914 
915               length = strfind(rec_str, "{");
916               if ( length!=-1 )
917               {
918                 printf("\n   setname:%s\n",  buffer);
919 
920                 setNr=getSetNrx(buffer);
921                 if(setNr<0) { printf("ERROR; Set:%s in boundaryField not known\n",buffer); exit(1); }
922 
923 		while(1)
924                 {
925                   length = frecord( handle, rec_str);
926                   if (rec_str[length] == (char)EOF) break;
927                   if ( strfind(rec_str, "}")!=-1) break;
928 
929 
930                   if ( strfind(rec_str, "fixedValue")!=-1)
931                   {
932 		    do
933                     {
934                       length = frecord( handle, rec_str); if (rec_str[length] == (char)EOF) break;
935                     }while(strfind(rec_str, "value")==-1);
936 
937                     /* reset all affected nodes to 0. */
938                     for (j=0; j<setx[setNr].anz_f; j++)
939                     {
940                       f=setx[setNr].face[j];
941                       if (face[f].type == 7) ipuf = 3;  /* TRI3 */
942                       else if (face[f].type== 8) ipuf = 6;  /* TRI6  */
943                       else if (face[f].type == 9) ipuf = 4;  /* QU4 */
944                       else if (face[f].type == 10) ipuf = 8;  /* QU4 */
945                       else
946                       {
947                         printf("\n\n ERROR: face.type:%d not known\n",face[f].type);
948                         exit(1);
949                       }
950                       for(n=0; n<lcase[anzx->l].ncomps; n++)
951                       {
952                         for(k=0; k<ipuf; k++)
953                         {
954     	                  lcase[anzx->l].dat[n][face[f].nod[k]]=0.;
955                           vpn[face[f].nod[k]]=0;
956                         }
957                       }
958 
959     		    }
960 
961                     if ( strfind(rec_str, "nonuniform")!=-1)
962                     {
963                       i=0;
964 		      f=-1;
965                       do
966 		      {
967   		        pos=0;
968                         do
969 		        {
970 		          length=get_values_in_string( &rec_str[pos], &fval[0]);
971                           if (length<1) break;
972                           pos+=length;
973 			  f=face[setx[setNr].elf[i].e].indx[setx[setNr].elf[i].f]; i++;
974 
975                           if (face[f].type == 7) ipuf = 3;  /* TRI3 */
976                           else if (face[f].type== 8) ipuf = 6;  /* TRI6  */
977                           else if (face[f].type == 9) ipuf = 4;  /* QU4 */
978                           else if (face[f].type == 10) ipuf = 8;  /* QU4 */
979                           else
980                           {
981                             printf("\n\n ERROR: face.type:%d not known\n",face[f].type);
982                             exit(1);
983                           }
984                           for(n=0; n<lcase[anzx->l].ncomps; n++)
985                           {
986                             for(j=0; j<ipuf; j++)
987                             {
988                               lcase[anzx->l].dat[n][face[f].nod[j]]+=fval[n];
989                               if(n==0) vpn[face[f].nod[j]]++; /* count the hits */
990                             }
991                           }
992 #if TEST
993 			  if(ncomps==1) printf("face:%d e:%d val:%lf \n", f, face[f].elem_nr, fval[0]);
994 			  if(ncomps==3) printf("face:%d e:%d vals:%lf %lf %lf\n", f, face[f].elem_nr, fval[0], fval[1], fval[2]);
995 #endif
996 			}while(1);
997 
998                         if((length<0)||(f==-1)) { length = frecord( handle, rec_str); if (rec_str[length] == (char)EOF) break; }
999                         else break;
1000 		      }while(1); /* as long as a parenthesis is open */
1001 
1002 		    }
1003                     else if ( strfind(rec_str, "uniform")!=-1)
1004 		    {
1005                       if     (ncomps==1) sscanf(rec_str, "%*s %*s %lf", &fval[0]);
1006                       else if(ncomps==3)
1007                       {
1008                         length = strfind(rec_str, "(");
1009                         sscanf(&rec_str[length+1], "%lf %lf %lf", &fval[0], &fval[1], &fval[2]);
1010   	              }
1011                       else { printf("ERROR3; class not known\n"); exit(1); }
1012 #if TEST
1013   		      if(ncomps==1)  printf("%s vals:%lf\n", &rec_str[length], fval[0]);
1014   		      if(ncomps==3)  printf("%s vals:%lf %lf %lf\n", &rec_str[length], fval[0], fval[1], fval[2]);
1015 #endif
1016                       for(i=0; i<setx[setNr].anz_elf; i++)
1017                       {
1018                         f=face[setx[setNr].elf[i].e].indx[setx[setNr].elf[i].f];
1019 
1020                         if (face[f].type == 7) ipuf = 3;  /* TRI3 */
1021                         else if (face[f].type== 8) ipuf = 6;  /* TRI6  */
1022                         else if (face[f].type == 9) ipuf = 4;  /* QU4 */
1023                         else if (face[f].type == 10) ipuf = 8;  /* QU4 */
1024                         else
1025                         {
1026                           printf("\n\n ERROR: face.type:%d not known\n",face[f].type);
1027                           exit(1);
1028                         }
1029                         for(n=0; n<lcase[anzx->l].ncomps; n++)
1030                         {
1031                           for(j=0; j<ipuf; j++)
1032                           {
1033 			    lcase[anzx->l].dat[n][face[f].nod[j]]+=fval[n];
1034                             vpn[face[f].nod[j]]++; /* count the hits */
1035                           }
1036                         }
1037     		      }
1038 		    }
1039 		  }
1040 
1041 
1042                 }
1043               }
1044               else sscanf(rec_str, "%s", buffer);
1045 
1046             }
1047 
1048           }
1049         }
1050 
1051         /* divide the nodal values by the amount of contributor elements */
1052         for(i=0; i<anzx->n; i++)
1053         {
1054           for(n=0; n<lcase[anzx->l].ncomps; n++)
1055           {
1056             lcase[anzx->l].dat[n][node[i].nr]/=vpn[node[i].nr];
1057 
1058             if(lcase[anzx->l].dat[n][node[i].nr] > lcase[anzx->l].max[n])
1059             {
1060               lcase[anzx->l].max[n]=lcase[anzx->l].dat[n][node[i].nr];
1061               lcase[anzx->l].nmax[n]=node[i].nr;
1062             }
1063             if(lcase[anzx->l].dat[n][node[i].nr] < lcase[anzx->l].min[n])
1064             {
1065               lcase[anzx->l].min[n]=lcase[anzx->l].dat[n][node[i].nr];
1066               lcase[anzx->l].nmin[n]=node[i].nr;
1067             }
1068           }
1069         }
1070 
1071         foundNoEntities:;
1072         fclose(handle);
1073         noFile:;
1074       }
1075      }
1076      closedir(dirp);
1077 
1078      /* if data were found calc derived data like Ma etc. (based on pre-defined kappa ) */
1079      printf("p:%d u:%d t:%d rho:%d\n",ds.p, ds.u, ds.t, ds.rho) ;
1080      if((ds.p!=-1)&&(ds.u!=-1)&&(ds.t!=-1)&&(ds.rho!=-1))
1081      {
1082        anzx->l++;
1083 
1084        /* store Datasets */
1085        if ( (lcase = (Datasets *)realloc(lcase, (anzx->l+1) * sizeof(Datasets))) == NULL )
1086          printf("\n\n ERROR: malloc failed, lcase\n\n") ;
1087 
1088        lcase[anzx->l].value=lcase[anzx->l-1].value;
1089        sprintf( lcase[anzx->l].dataset_name,"%s", lcase[anzx->l-1].dataset_name);
1090 
1091        lcase[anzx->l].irtype = 1;
1092        lcase[anzx->l].npheader = 0;
1093        strcpy(lcase[anzx->l].analysis_name,"");
1094        strcpy(lcase[anzx->l].name, "derived");
1095        strcpy(lcase[anzx->l].dataset_text,"");
1096        lcase[anzx->l].step_number=lcase[anzx->l-1].step_number;
1097        lcase[anzx->l].analysis_type=1;
1098 
1099        printf ("lcase.name[%d]= %s value:%f\n", anzx->l, lcase[anzx->l].name, lcase[anzx->l].value);
1100        lcase[anzx->l].ncomps  = 3;
1101 
1102        if ( (lcase[anzx->l].nmax = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1103          printf("\n\n ERROR: malloc failure\n\n" );
1104        if ( (lcase[anzx->l].nmin = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1105          printf("\n\n ERROR: malloc failure\n\n" );
1106        if ( (lcase[anzx->l].max = (float *)malloc( lcase[anzx->l].ncomps * sizeof(float))) == NULL )
1107          printf("\n\n ERROR: malloc failure\n\n" );
1108        if ( (lcase[anzx->l].min = (float *)malloc( lcase[anzx->l].ncomps * sizeof(float))) == NULL )
1109          printf("\n\n ERROR: malloc failure\n\n" );
1110        if ( (lcase[anzx->l].dat = (float **)malloc( lcase[anzx->l].ncomps * sizeof(float *))) == NULL )
1111          printf("\n\n ERROR: malloc failure\n\n" );
1112        if ( (lcase[anzx->l].icname = (char **)malloc( lcase[anzx->l].ncomps * sizeof(char *))) == NULL )
1113          printf("\n\n ERROR: malloc failure\n\n" );
1114        if ( (lcase[anzx->l].compName = (char **)malloc( lcase[anzx->l].ncomps * sizeof(char *))) == NULL )
1115          printf("\n\n ERROR: malloc failure\n\n" );
1116 
1117        for(i=0; i<lcase[anzx->l].ncomps; i++)
1118        {
1119          if ( (lcase[anzx->l].compName[i] = (char *)malloc( MAX_LINE_LENGTH * sizeof(char))) == NULL )
1120            printf("\n\n ERROR: malloc failed\n\n" );
1121          if ( (lcase[anzx->l].dat[i]      = (float *)malloc( (anzx->nmax+1) * sizeof(float))) == NULL )
1122            printf("\n\n ERROR: malloc failed\n\n" );
1123          if ( (lcase[anzx->l].icname[i]   = (char *)malloc( MAX_LINE_LENGTH * sizeof(char))) == NULL )
1124            printf("\n\n ERROR: malloc failed\n\n" );
1125          lcase[anzx->l].max[i]=-MAX_FLOAT;
1126          lcase[anzx->l].min[i]=MAX_FLOAT;
1127        }
1128        if ( (lcase[anzx->l].menu = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1129          printf("\n\n ERROR: malloc failure\n\n" );
1130        if ( (lcase[anzx->l].ictype = (int *)malloc(lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1131          printf("\n\n ERROR: malloc failure\n\n" );
1132        if ( (lcase[anzx->l].icind1 = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1133          printf("\n\n ERROR: malloc failure\n\n" );
1134        if ( (lcase[anzx->l].icind2 = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1135          printf("\n\n ERROR: malloc failure\n\n" );
1136        if ( (lcase[anzx->l].iexist = (int *)malloc( lcase[anzx->l].ncomps * sizeof(int))) == NULL )
1137          printf("\n\n ERROR: malloc failure\n\n" );
1138 
1139        for(i=0; i<lcase[anzx->l].ncomps; i++)
1140        {
1141          lcase[anzx->l].menu[i] = 1;
1142          lcase[anzx->l].ictype[i] = 1;
1143          lcase[anzx->l].icind1[i] = 0;
1144          lcase[anzx->l].icind2[i] = 0;
1145          lcase[anzx->l].iexist[i] = 0;
1146        }
1147 
1148        sprintf(lcase[anzx->l].compName[0],"TT ");
1149        sprintf(lcase[anzx->l].compName[1],"PT ");
1150        sprintf(lcase[anzx->l].compName[2],"Ma ");
1151 
1152        for(i=0; i<anzx->n; i++)
1153        {
1154          ps=lcase[ds.p].dat[0][node[i].nr];
1155          us=sqrt(lcase[ds.u].dat[0][node[i].nr]*lcase[ds.u].dat[0][node[i].nr]+lcase[ds.u].dat[1][node[i].nr]*lcase[ds.u].dat[1][node[i].nr]+lcase[ds.u].dat[2][node[i].nr]*lcase[ds.u].dat[2][node[i].nr]);
1156          ts=lcase[ds.t].dat[0][node[i].nr];
1157          rho=lcase[ds.rho].dat[0][node[i].nr];
1158          Rg=ps/rho/ts;
1159          cp=Rg/(1-1/kappa);
1160 	 //printf("R:%f cp:%f\n", Rg, cp);
1161          tt=us*us*.5/cp + ts;
1162 	 pt = ps/pow( (ts/tt), (kappa/(kappa-1)) );
1163          ma=us/sqrt(kappa*ps/rho);
1164          lcase[anzx->l].dat[0][node[i].nr]=tt;
1165          lcase[anzx->l].dat[1][node[i].nr]=pt;
1166          lcase[anzx->l].dat[2][node[i].nr]=ma;
1167 
1168          for(j=0; j<lcase[anzx->l].ncomps; j++)
1169 	 {
1170            if(lcase[anzx->l].dat[j][node[i].nr] > lcase[anzx->l].max[j])
1171            {
1172              lcase[anzx->l].max[j]=lcase[anzx->l].dat[j][node[i].nr];
1173              lcase[anzx->l].nmax[j]=node[i].nr;
1174            }
1175            if(lcase[anzx->l].dat[j][node[i].nr] < lcase[anzx->l].min[j])
1176            {
1177              lcase[anzx->l].min[j]=lcase[anzx->l].dat[j][node[i].nr];
1178              lcase[anzx->l].nmin[j]=node[i].nr;
1179            }
1180          }
1181        }
1182      }
1183 
1184     }
1185   }
1186   free(vpn);
1187   anzx->l++;
1188 
1189 
1190   /* free the temporary faces */
1191   free(face);
1192   anzx->f=0;
1193 
1194   /* set all .loaded to 1 to indicate that the data are available */
1195   for (i=0; i<anzx->l; i++) { lcase[i].fileptr=NULL; lcase[i].loaded=1; }
1196 
1197   /* delete all faces in sets */
1198   for(i=0; i<anzx->sets; i++) setx[i].anz_f=0;
1199 
1200   *sptr=setx; *nptr =  node; *eptr = elem; *lptr = lcase;
1201   return (1);
1202 }
1203