1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <string.h>
4 #include <math.h>
5 #if defined(_WIN32) && !defined(__NUTC__)
6 # include <io.h>
7 # define unlink _unlink
8 #else
9 # include <unistd.h>
10 #endif
11 #include "cgnslib.h"
12 
13 /* set UNSTRUCTURED_1TO1 to write the unstructured case
14    zone connectivities as 1to1 instead of abutting1to1
15 #define UNSTRUCTURED_1TO1
16 */
17 
18 /* set ABUTTING1TO1_FACES to write the unstructured case
19    zone connectivites as abutting1to1 with faces(elements)
20    instead of points.
21 #define ABUTTING1TO1_FACES
22 */
23 
24 int CellDim = 3, PhyDim = 3;
25 
26 int cgfile, cgbase, cgzone;
27 int CellDim, PhyDim;
28 cgsize_t size[9];
29 
30 #define NUM_ZCONN  5
31 #define ZCONN_NAME "ZoneConnectivity"
32 
33 #define NUM_SIDE 5
34 #define NODE_INDEX(I,J,K) ((I)+NUM_SIDE*(((J)-1)+NUM_SIDE*((K)-1)))
35 #define CELL_INDEX(I,J,K) ((I)+(NUM_SIDE-1)*(((J)-1)+(NUM_SIDE-1)*((K)-1)))
36 
37 int num_coord;
38 float *xcoord, *ycoord, *zcoord;
39 int num_element, num_face;
40 cgsize_t *elements, *faces, *parent;
41 
42 int npts;
43 cgsize_t *pts, *d_pts;
44 
45 char errmsg[256];
46 
47 static float exponents[5] = {0, 1, 0, 0, 0};
48 
49 void init_data();
50 void write_structured(), write_unstructured();
51 void test_zconn(), test_subreg();
52 
error_exit(char * where)53 void error_exit (char *where)
54 {
55     fprintf (stderr, "ERROR:%s:%s\n", where, cg_get_error());
56     exit (1);
57 }
58 
main(int argc,char * argv[])59 int main (int argc, char *argv[])
60 {
61     char outfile[32];
62     int nbases, nzones;
63 
64     strcpy (outfile, "ver31.cgns");
65     if (argc > 1) {
66         int type = 0;
67         int n = 0;
68         if (argv[1][n] == '-') n++;
69         if (argv[1][n] == 'a' || argv[1][n] == 'A' || argv[1][n] == '1') {
70             if (NULL != strchr(argv[1], '2'))
71                 type = CG_FILE_ADF2;
72             else
73                 type = CG_FILE_ADF;
74         }
75         else if (argv[1][n] == 'h' || argv[1][n] == 'H' || argv[1][n] == '2') {
76             type = CG_FILE_HDF5;
77         }
78         else if (argv[1][n] == '3') {
79             type = CG_FILE_ADF2;
80         }
81         else {
82             fprintf(stderr, "unknown option\n");
83             exit (1);
84         }
85         if (cg_set_file_type(type))
86             error_exit("cg_set_file_type");
87     }
88     init_data();
89     unlink(outfile);
90     if (cg_open(outfile, CG_MODE_WRITE, &cgfile)) error_exit("cg_open");
91     write_structured();
92     write_unstructured();
93     if (cg_close(cgfile)) error_exit("cg_close");
94 
95     if (cg_open(outfile, CG_MODE_READ, &cgfile)) error_exit("cg_open");
96     if (cg_nbases(cgfile, &nbases)) error_exit("cg_nbases");
97     for (cgbase = 1; cgbase <= nbases; cgbase++) {
98         if (cg_nzones(cgfile, cgbase, &nzones)) error_exit("cg_nzones");
99         for (cgzone = 1; cgzone <= nzones; cgzone++) {
100             test_zconn();
101             test_subreg();
102         }
103     }
104     if (cg_close(cgfile)) error_exit("cg_close");
105 
106     return 0;
107 }
108 
init_data()109 void init_data()
110 {
111     int n, i, j, k, nn, nf, np;
112 
113     /* compute coordinates - make it twice as big for use with cylindrical */
114 
115     num_coord = NUM_SIDE * NUM_SIDE * NUM_SIDE;
116     xcoord = (float *) malloc (6 * num_coord * sizeof(float));
117     if (NULL == xcoord) {
118         fprintf(stderr, "malloc failed for coordinates\n");
119         exit(1);
120     }
121     ycoord = xcoord + 2 * num_coord;
122     zcoord = ycoord + 2 * num_coord;
123     for (n = 0, k = 0; k < NUM_SIDE; k++) {
124         for (j = 0; j < NUM_SIDE; j++) {
125             for (i = 0; i < NUM_SIDE; i++, n++) {
126                 xcoord[n] = (float)i;
127                 ycoord[n] = (float)j;
128                 zcoord[n] = (float)k;
129             }
130         }
131     }
132 
133     /* compute elements */
134 
135     num_element = (NUM_SIDE - 1) * (NUM_SIDE - 1) * (NUM_SIDE - 1);
136     elements = (cgsize_t *) malloc (8 * num_element * sizeof(cgsize_t));
137     if (NULL == elements) {
138         fprintf(stderr, "malloc failed for elements");
139         exit(1);
140     }
141     for (n = 0, k = 1; k < NUM_SIDE; k++) {
142         for (j = 1; j < NUM_SIDE; j++) {
143             for (i = 1; i < NUM_SIDE; i++) {
144                 nn = NODE_INDEX(i, j, k);
145                 elements[n++] = nn;
146                 elements[n++] = nn + 1;
147                 elements[n++] = nn + 1 + NUM_SIDE;
148                 elements[n++] = nn + NUM_SIDE;
149                 nn += NUM_SIDE * NUM_SIDE;
150                 elements[n++] = nn;
151                 elements[n++] = nn + 1;
152                 elements[n++] = nn + 1 + NUM_SIDE;
153                 elements[n++] = nn + NUM_SIDE;
154             }
155         }
156     }
157 
158     /* compute outside face elements */
159 
160     num_face = 6 * (NUM_SIDE - 1) * (NUM_SIDE - 1);
161     faces = (cgsize_t *) malloc (4 * num_face * sizeof(cgsize_t));
162     parent = (cgsize_t *) malloc (4 * num_face * sizeof(cgsize_t));
163     if (NULL == faces || NULL == parent) {
164         fprintf(stderr, "malloc failed for elements");
165         exit(1);
166     }
167     for (n = 0; n < 4*num_face; n++)
168         parent[n] = 0;
169     nf = np = 0;
170     n = 2 * num_face;
171     i = 1;
172     for (k = 1; k < NUM_SIDE; k++) {
173         for (j = 1; j < NUM_SIDE; j++) {
174             nn = NODE_INDEX(i, j, k);
175             faces[nf++]  = nn;
176             faces[nf++]  = nn + NUM_SIDE * NUM_SIDE;
177             faces[nf++]  = nn + NUM_SIDE * (NUM_SIDE + 1);
178             faces[nf++]  = nn + NUM_SIDE;
179             parent[np]   = CELL_INDEX(i, j, k);
180             parent[np+n] = 5;
181             np++;
182         }
183     }
184     i = NUM_SIDE;
185     for (k = 1; k < NUM_SIDE; k++) {
186         for (j = 1; j < NUM_SIDE; j++) {
187             nn = NODE_INDEX(i, j, k);
188             faces[nf++]  = nn;
189             faces[nf++]  = nn + NUM_SIDE;
190             faces[nf++]  = nn + NUM_SIDE * (NUM_SIDE + 1);
191             faces[nf++]  = nn + NUM_SIDE * NUM_SIDE;
192             parent[np]   = CELL_INDEX(i-1, j, k);
193             parent[np+n] = 3;
194             np++;
195         }
196     }
197     j = 1;
198     for (k = 1; k < NUM_SIDE; k++) {
199         for (i = 1; i < NUM_SIDE; i++) {
200             nn = NODE_INDEX(i, j, k);
201             faces[nf++]  = nn;
202             faces[nf++]  = nn + 1;
203             faces[nf++]  = nn + 1 + NUM_SIDE * NUM_SIDE;
204             faces[nf++]  = nn + NUM_SIDE * NUM_SIDE;
205             parent[np]   = CELL_INDEX(i, j, k);
206             parent[np+n] = 2;
207             np++;
208         }
209     }
210     j = NUM_SIDE;
211     for (k = 1; k < NUM_SIDE; k++) {
212         for (i = 1; i < NUM_SIDE; i++) {
213             nn = NODE_INDEX(i, j, k);
214             faces[nf++]  = nn;
215             faces[nf++]  = nn + NUM_SIDE * NUM_SIDE;
216             faces[nf++]  = nn + 1 + NUM_SIDE * NUM_SIDE;
217             faces[nf++]  = nn + 1;
218             parent[np]   = CELL_INDEX(i, j-1, k);
219             parent[np+n] = 4;
220             np++;
221         }
222     }
223     k = 1;
224     for (j = 1; j < NUM_SIDE; j++) {
225         for (i = 1; i < NUM_SIDE; i++) {
226             nn = NODE_INDEX(i, j, k);
227             faces[nf++]  = nn;
228             faces[nf++]  = nn + NUM_SIDE;
229             faces[nf++]  = nn + NUM_SIDE + 1;
230             faces[nf++]  = nn + 1;
231             parent[np]   = CELL_INDEX(i, j, k);
232             parent[np+n] = 1;
233             np++;
234         }
235     }
236     k = NUM_SIDE;
237     for (j = 1; j < NUM_SIDE; j++) {
238         for (i = 1; i < NUM_SIDE; i++) {
239             nn = NODE_INDEX(i, j, k);
240             faces[nf++]  = nn;
241             faces[nf++]  = nn + 1;
242             faces[nf++]  = nn + NUM_SIDE + 1;
243             faces[nf++]  = nn + NUM_SIDE;
244             parent[np]   = CELL_INDEX(i, j, k-1);
245             parent[np+n] = 6;
246             np++;
247         }
248     }
249 
250     /* connectivity points - make it big enough to hold 4 surfaces */
251 
252     npts = NUM_SIDE * NUM_SIDE;
253     pts = (cgsize_t *) malloc (12 * npts * sizeof(cgsize_t));
254     if (NULL == pts) {
255         fprintf(stderr, "malloc failed for connectivity points");
256         exit(1);
257     }
258     d_pts = pts + 6 * npts;
259 }
260 
write_coords(int nz)261 void write_coords(int nz)
262 {
263     int k, nn, n, nij, koff, cgcoord;
264 
265     koff = nz == 1 ? 1 - NUM_SIDE : 0;
266     nij = NUM_SIDE * NUM_SIDE;
267     for (n = 0, k = 0; k < NUM_SIDE; k++) {
268         for (nn = 0; nn < nij; nn++)
269             zcoord[n++] = (float)(k + koff);
270     }
271 
272     if (cg_coord_write(cgfile, cgbase, nz, CGNS_ENUMV(RealSingle),
273             "CoordinateX", xcoord, &cgcoord) ||
274         cg_coord_write(cgfile, cgbase, nz, CGNS_ENUMV(RealSingle),
275             "CoordinateY", ycoord, &cgcoord) ||
276         cg_coord_write(cgfile, cgbase, nz, CGNS_ENUMV(RealSingle),
277             "CoordinateZ", zcoord, &cgcoord)) {
278         sprintf (errmsg, "zone %d coordinates", nz);
279         error_exit(errmsg);
280     }
281     for (n = 1; n <= 3; n++) {
282         if (cg_goto(cgfile, cgbase, "Zone_t", nz, "GridCoordinates_t", 1,
283                 "DataArray_t", n, NULL) ||
284             cg_exponents_write(CGNS_ENUMV(RealSingle), exponents))
285             error_exit("coordinate exponents");
286     }
287 }
288 
write_elements(int nz)289 void write_elements(int nz)
290 {
291     int cgsect;
292 
293     if (cg_section_write(cgfile, cgbase, nz, "Elements", CGNS_ENUMV(HEXA_8),
294             1, num_element, 0, elements, &cgsect) ||
295         cg_section_write(cgfile, cgbase, nz, "Faces", CGNS_ENUMV(QUAD_4),
296             num_element+1, num_element+num_face, 0, faces, &cgsect) ||
297         cg_parent_data_write(cgfile, cgbase, nz, cgsect, parent)) {
298         sprintf (errmsg, "zone %d elements", nz);
299         error_exit(errmsg);
300     }
301 }
302 
303 /*------------------------------------------------------------------------
304  * structured grid
305  *------------------------------------------------------------------------*/
306 
write_structured()307 void write_structured()
308 {
309     int i, j, n, nc, cgfam, cgz, cgbc, cgconn, cghole, cgsr;
310     cgsize_t range[6], d_range[6], dims[2];
311     int transform[3], rind[6], iter[NUM_ZCONN];
312     char name[33], cname[33], hname[33], sname[33], zcname[33];
313     char cpointers[32*NUM_ZCONN+1], spointers[32*NUM_ZCONN+1];
314 
315     if (cg_base_write(cgfile, "Structured", CellDim, PhyDim, &cgbase) ||
316         cg_goto(cgfile, cgbase, "end") ||
317         cg_descriptor_write("Descriptor", "Multi-block Structured Grid") ||
318         cg_dataclass_write(CGNS_ENUMV(Dimensional)) ||
319         cg_units_write(CGNS_ENUMV(Kilogram), CGNS_ENUMV(Meter),
320             CGNS_ENUMV(Second), CGNS_ENUMV(Kelvin), CGNS_ENUMV(Radian)))
321         error_exit("structured base");
322     if (cg_simulation_type_write(cgfile, cgbase, CGNS_ENUMV(TimeAccurate)))
323         error_exit("simulation type");
324 
325     /* write a FamilyBCDataset node and children */
326 
327     if (cg_family_write(cgfile, cgbase, "Family", &cgfam) ||
328         cg_fambc_write(cgfile, cgbase, cgfam, "Wall",
329             CGNS_ENUMV(BCWall), &cgbc) ||
330         cg_goto(cgfile, cgbase, "Family_t", 1, "FamilyBC_t", 1, NULL) ||
331         cg_bcdataset_write("Dataset", CGNS_ENUMV(BCWall),
332             CGNS_ENUMV(Dirichlet)) ||
333         cg_gorel(cgfile, "Dataset", 0, NULL) ||
334 /*        cg_goto(cgfile, cgbase, "Family_t", 1, "FamilyBC_t", 1,
335             "FamilyBCDataSet_t", 1, NULL) ||*/
336         cg_descriptor_write("descriptor", "this is a descriptor") ||
337         cg_dataclass_write(CGNS_ENUMV(NormalizedByDimensional)) ||
338         cg_units_write(CGNS_ENUMV(Kilogram), CGNS_ENUMV(Meter),
339             CGNS_ENUMV(Second), CGNS_ENUMV(Kelvin), CGNS_ENUMV(Radian)) ||
340         cg_state_write("reference state") ||
341         cg_user_data_write("Userdata") ||
342         cg_gorel(cgfile, "UserDefinedData_t", 1, NULL) ||
343         cg_famname_write("Family"))
344         error_exit("family bc");
345 
346     /* write zones */
347 
348     for (n = 0; n < 3; n++) {
349         size[n]   = NUM_SIDE;
350         size[n+3] = NUM_SIDE - 1;
351         size[n+6] = 0;
352     }
353     for (n = 1; n <= 2; n++) {
354         sprintf(name, "Zone%d", n);
355         if (cg_zone_write(cgfile, cgbase, name, size,
356                 CGNS_ENUMV(Structured), &cgzone)) {
357             sprintf (errmsg, "structured zone %d", n);
358             error_exit(errmsg);
359         }
360         write_coords(n);
361     }
362 
363     /* write inlet BC (zone 1) as point range */
364 
365     for (n = 0; n < 3; n++) {
366         range[n] = 1;
367         range[n+3] = NUM_SIDE;
368     }
369     range[5] = 1;
370     if (cg_boco_write(cgfile, cgbase, 1, "Inlet", CGNS_ENUMV(BCInflow),
371             CGNS_ENUMV(PointRange), 2, range, &cgbc))
372         error_exit("inlet boco");
373 
374     /* write outlet BC (zone 2) as point list */
375 
376     for (n = 0, j = 1; j <= NUM_SIDE; j++) {
377         for (i = 1; i <= NUM_SIDE; i++) {
378             pts[n++] = i;
379             pts[n++] = j;
380             pts[n++] = NUM_SIDE;
381         }
382     }
383     if (cg_boco_write(cgfile, cgbase, 2, "Outlet", CGNS_ENUMV(BCOutflow),
384             CGNS_ENUMV(PointList), npts, pts, &cgbc))
385         error_exit("outlet boco");
386 
387     /* write connectivities in multiple ZoneGridConnectivity_t nodes */
388 
389     for (nc = 1; nc <= NUM_ZCONN; nc++) {
390         /* create ZoneGridConnectivity_t node */
391         sprintf(zcname, "%s%d", ZCONN_NAME, nc);
392         if (cg_zconn_write(cgfile, cgbase, 1, zcname, &cgz) ||
393             cg_zconn_write(cgfile, cgbase, 2, zcname, &cgz))
394             error_exit(name);
395         sprintf(cname, "conn%d", nc);
396         sprintf(hname, "hole%d", nc);
397 
398         /* write zone 1 to zone 2 connectivity as 1to1 */
399 
400         for (n = 0; n < 3; n++) {
401             range[n] = d_range[n] = 1;
402             range[n+3] = d_range[n+3] = NUM_SIDE;
403             transform[n] = n + 1;
404         }
405         range[2] = NUM_SIDE;
406         d_range[5] = 1;
407         if (cg_1to1_write(cgfile, cgbase, 1, cname, "Zone2",
408                 range, d_range, transform, &cgconn)) {
409             sprintf(errmsg, "Zone1/%s/%s", zcname, cname);
410             error_exit(errmsg);
411         }
412 
413         /* write zone 2 to zone 1 connectivity as Abutting1to1 */
414 
415         for (n = 0; n < 3; n++) {
416             range[n] = 1;
417             range[n+3] = NUM_SIDE;
418         }
419         range[5] = 1;
420         for (n = 0, j = 1; j <= NUM_SIDE; j++) {
421             for (i = 1; i <= NUM_SIDE; i++) {
422                 d_pts[n++] = i;
423                 d_pts[n++] = j;
424                 d_pts[n++] = 1;
425             }
426         }
427         if (cg_conn_write(cgfile, cgbase, 2, cname,
428                 CGNS_ENUMV(Vertex), CGNS_ENUMV(Abutting1to1),
429                 CGNS_ENUMV(PointRange), 2, range, "Zone1",
430                 CGNS_ENUMV(Structured), CGNS_ENUMV(PointListDonor),
431                 CGNS_ENUMV(Integer), npts, d_pts, &cgconn)) {
432             sprintf(errmsg, "Zone2/%s/%s", zcname, cname);
433             error_exit(errmsg);
434         }
435 
436         /* write hole in zone 1 as PointList */
437 
438         if (cg_hole_write(cgfile, cgbase, 1, hname, CGNS_ENUMV(Vertex),
439                 CGNS_ENUMV(PointList), 1, npts, d_pts, &cghole)) {
440             sprintf(errmsg, "Zone1/%s/%s", zcname, hname);
441             error_exit(errmsg);
442         }
443 
444         /* write hole in zone 2 as PointRange */
445 
446         if (cg_hole_write(cgfile, cgbase, 2, hname, CGNS_ENUMV(Vertex),
447                 CGNS_ENUMV(PointRange), 1, 2, range, &cghole)) {
448             sprintf(errmsg, "Zone2/%s/%s", zcname, hname);
449             error_exit(errmsg);
450         }
451     }
452 
453     /* write SubRegion 1 as BCRegionName */
454 
455     if (cg_subreg_bcname_write(cgfile, cgbase, 1, "SubRegion1",
456             2, "Inlet", &cgsr))
457         error_exit("cg_subreg_bcname_write(Inlet)");
458     if (cg_goto(cgfile, cgbase, "Zone_t", 1, "ZoneSubRegion_t", 1, NULL))
459         error_exit("cg_goto ZoneSubRegion_t 1");
460     if (cg_dataclass_write(CGNS_ENUMV(Dimensional)))
461         error_exit("cg_dataclass_write Dimensional");
462     if (cg_units_write(CGNS_ENUMV(Kilogram), CGNS_ENUMV(Meter),
463             CGNS_ENUMV(Second), CGNS_ENUMV(Kelvin), CGNS_ENUMV(Degree)))
464         error_exit("cg_units_write");
465     if (cg_descriptor_write("descriptor", "this is zone 1, subregion 1"))
466         error_exit("cg_descriptor_write");
467 
468     if (cg_subreg_bcname_write(cgfile, cgbase, 2, "SubRegion1",
469             2, "Outlet", &cgsr))
470         error_exit("cg_subreg_bcname_write(Outlet)");
471     if (cg_goto(cgfile, cgbase, "Zone_t", 2, "SubRegion1", 0, NULL))
472         error_exit("cg_goto SubRegion1");
473     if (cg_famname_write("Family"))
474         error_exit("cg_famname_write SubRegion1");
475     if (cg_user_data_write("Userdata"))
476         error_exit("cg_user_data_write");
477     if (cg_gorel(cgfile, "UserDefinedData_t", 1, NULL))
478         error_exit("cg_gorel UserDefinedData");
479     if (cg_famname_write("Family"))
480         error_exit("cg_famname_write SubRegion1/UserDefinedData");
481 
482     /* write SubRegion 2 as ConnectivityRegionName */
483 
484     sprintf(cname, "%s2/conn2", ZCONN_NAME);
485     strcpy(sname, "SubRegion2");
486     if (cg_subreg_gcname_write(cgfile, cgbase, 1, sname, 2, cname, &cgsr) ||
487         cg_subreg_gcname_write(cgfile, cgbase, 2, sname, 2, cname, &cgsr)) {
488         sprintf(errmsg, "cg_subreg_gcname_write(%s)", sname);
489         error_exit(errmsg);
490     }
491 
492     /* write subregion as PointRange in Zone1 and PointList in Zone2 */
493 
494     for (n = 0; n < 3; n++) {
495         range[n] = 1;
496         range[n+3] = NUM_SIDE;
497     }
498     range[5] = 1;
499     for (n = 0, j = 1; j < NUM_SIDE; j++) {
500         for (i = 1; i < NUM_SIDE; i++) {
501             pts[n++] = i;
502             pts[n++] = j;
503             pts[n++] = NUM_SIDE;
504         }
505     }
506     for (n = 0; n < 6; n++)
507         rind[n] = 1;
508     dims[0] = NUM_SIDE * NUM_SIDE;
509     dims[1] = (NUM_SIDE - 1) * (NUM_SIDE - 1);
510 
511     for (nc = 3; nc <= NUM_ZCONN; nc++) {
512         sprintf(sname, "SubRegion%d", nc);
513         if (cg_subreg_ptset_write(cgfile, cgbase, 1, sname, 2,
514                 CGNS_ENUMV(Vertex), CGNS_ENUMV(PointRange),
515                 2, range, &cgsr)) {
516             sprintf(errmsg, "cg_subreg_ptset_write(%s PointRange)", sname);
517             error_exit(errmsg);
518         }
519         if (cg_goto(cgfile, cgbase, "Zone_t", 1, sname, 0, NULL)) {
520             sprintf(errmsg, "cg_goto %s", sname);
521             error_exit(errmsg);
522         }
523         if (cg_rind_write(rind)) {
524             sprintf(errmsg, "cg_rind_write %s", sname);
525             error_exit(errmsg);
526         }
527         if (cg_array_write("DataArray", CGNS_ENUMV(RealSingle), 1,
528                 dims, xcoord)) {
529             sprintf(errmsg, "cg_array_write %s", sname);
530             error_exit(errmsg);
531         }
532         if (cg_gorel(cgfile, "DataArray", 0, NULL) ||
533             cg_exponents_write(CGNS_ENUMV(RealSingle), exponents)) {
534             sprintf(errmsg, "cg_exponents_write %s", sname);
535             error_exit(errmsg);
536         }
537 
538         if (cg_subreg_ptset_write(cgfile, cgbase, 2, sname, 2,
539                 CGNS_ENUMV(KFaceCenter), CGNS_ENUMV(PointList),
540                 dims[1], pts, &cgsr)) {
541             sprintf(errmsg, "cg_subreg_ptset_write(%s PointList)", sname);
542             error_exit(errmsg);
543         }
544         if (cg_goto(cgfile, cgbase, "Zone_t", 2, sname, 0, NULL)) {
545             sprintf(errmsg, "cg_goto %s", sname);
546             error_exit(errmsg);
547         }
548         if (cg_array_write("DataArray", CGNS_ENUMV(RealSingle), 1,
549                 &dims[1], xcoord)) {
550             sprintf(errmsg, "cg_array_write %s", sname);
551             error_exit(errmsg);
552         }
553         if (cg_gorel(cgfile, "DataArray", 0, NULL) ||
554             cg_exponents_write(CGNS_ENUMV(RealSingle), exponents)) {
555             sprintf(errmsg, "cg_exponents_write %s", sname);
556             error_exit(errmsg);
557         }
558     }
559 
560     /* create BaseIterativeData_t node */
561 
562     for (n = 0; n < NUM_ZCONN; n++)
563         iter[n] = n + 1;
564     dims[0] = NUM_ZCONN;
565 
566     if (cg_biter_write(cgfile, cgbase, "BaseIterativeData", NUM_ZCONN))
567         error_exit("cg_biter_write");
568     if (cg_goto(cgfile, cgbase, "BaseIterativeData", 0, NULL))
569         error_exit("cg_goto BaseIterativeData");
570     if (cg_array_write("IterationValues", CGNS_ENUMV(Integer),
571         1, dims, iter)) error_exit("cg_array_write IterationValues");
572 
573     /* create the ZoneIterativeData_t node */
574 
575     dims[0] = 32;
576     dims[1] = NUM_ZCONN;
577     n = 32 * NUM_ZCONN;
578     for (i = 0; i < n; i++) {
579         cpointers[i] = ' ';
580         spointers[i] = ' ';
581     }
582     cpointers[n] = 0;
583     spointers[n] = 0;
584 
585     for (n = 0, nc = 1; nc <= NUM_ZCONN; nc++, n+=32) {
586         sprintf(zcname, "%s%d", ZCONN_NAME, nc);
587         i = strlen(zcname);
588         strcpy(&cpointers[n], zcname);
589         cpointers[n+i] = ' ';
590         sprintf(sname, "SubRegion%d", nc);
591         i = strlen(sname);
592         strcpy(&spointers[n], sname);
593         spointers[n+i] = ' ';
594     }
595 
596     for (n = 1; n <= 2; n++) {
597         sprintf(name, "Zone%dInterativeData", n);
598         if (cg_ziter_write(cgfile, cgbase, n, name)) {
599             sprintf(errmsg, "%s:cg_ziter_write", name);
600             error_exit(errmsg);
601         }
602         if (cg_goto(cgfile, cgbase, "Zone_t", n, name, 0, NULL)) {
603             sprintf(errmsg, "%s:cg_goto", name);
604             error_exit(errmsg);
605         }
606         if (cg_array_write("ZoneGridConnectivityPointers",
607                 CGNS_ENUMV(Character), 2, dims, cpointers)) {
608             sprintf(errmsg, "%s:cg_array_write ZoneGridConnectivityPointers",
609                 name);
610             error_exit(errmsg);
611         }
612         if (cg_array_write("ZoneSubRegionPointers",
613                 CGNS_ENUMV(Character), 2, dims, spointers)) {
614             sprintf(errmsg, "%s:cg_array_write ZoneSubRegionPointers",
615                 name);
616             error_exit(errmsg);
617         }
618     }
619 }
620 
621 /*------------------------------------------------------------------------
622  * unstructured grid
623  *------------------------------------------------------------------------*/
624 
write_unstructured()625 void write_unstructured()
626 {
627     int n, nc, cgconn, cgz, cghole;
628     int iter[NUM_ZCONN];
629     cgsize_t range[2], dims[2];
630     char name[33], zcname[33], pointers[32*NUM_ZCONN+1];
631 #ifdef UNSTRUCTURED_1TO1
632     int d_range[2], transform;
633 #else
634 # ifdef ABUTTING1TO1_FACES
635     GridLocation_t location;
636     int nelem;
637 # else
638     int i, j;
639 # endif
640 #endif
641 
642     if (cg_base_write(cgfile, "Unstructured", CellDim, PhyDim, &cgbase) ||
643         cg_goto(cgfile, cgbase, "end") ||
644         cg_descriptor_write("Descriptor", "Multi-block Unstructured Grid") ||
645         cg_dataclass_write(CGNS_ENUMV(NormalizedByDimensional)) ||
646         cg_units_write(CGNS_ENUMV(Kilogram), CGNS_ENUMV(Meter),
647             CGNS_ENUMV(Second), CGNS_ENUMV(Kelvin), CGNS_ENUMV(Radian)))
648         error_exit("unstructured base");
649 
650     /* write zones */
651 
652     for (n = 0; n < 9; n++)
653         size[n] = 0;
654     size[0] = num_coord;
655     size[1] = num_element;
656     for (n = 1; n <= 2; n++) {
657         sprintf(name, "Zone%d", n);
658         if (cg_zone_write(cgfile, cgbase, name, size,
659                 CGNS_ENUMV(Unstructured), &cgzone)) {
660             sprintf (errmsg, "unstructured zone %d", n);
661             error_exit(errmsg);
662         }
663         write_coords(n);
664         write_elements(n);
665     }
666 
667 #ifdef ABUTTING1TO1_FACES
668     int nelem = (NUM_SIDE - 1) * (NUM_SIDE - 1);
669 #endif
670 
671     /* write connectivities in multiple ZoneGridConnectivity_t nodes */
672     for (nc = 1; nc <= NUM_ZCONN; nc++) {
673         /* create ZoneGridConnectivity_t node */
674         sprintf(zcname, "%s%d", ZCONN_NAME, nc);
675         if (cg_zconn_write(cgfile, cgbase, 1, zcname, &cgz) ||
676             cg_zconn_write(cgfile, cgbase, 2, zcname, &cgz))
677             error_exit(name);
678 
679         sprintf(name, "conn%d", nc);
680 
681 #ifdef UNSTRUCTURED_1TO1
682 
683         /* write connectivities as 1to1 */
684 
685         range[0] = NODE_INDEX(1, 1, NUM_SIDE);
686         range[1] = NODE_INDEX(NUM_SIDE, NUM_SIDE, NUM_SIDE);
687         d_range[0] = NODE_INDEX(1, 1, 1);
688         d_range[1] = NODE_INDEX(NUM_SIDE, NUM_SIDE, 1);
689         transform = 1;
690         if (cg_1to1_write(cgfile, cgbase, 1, name, "Zone2",
691                 range, d_range, &transform, &cgconn)) {
692             sprintf(errmsg, "Zone1/%s/%s", name);
693             error_exit(errmsg);
694         }
695 
696         if (cg_1to1_write(cgfile, cgbase, 2, name, "Zone1",
697                 d_range, range, &transform, &cgconn)) {
698             sprintf(errmsg, "Zone2/%s/%s", zcname, name);
699             error_exit(errmsg);
700         }
701 
702 #else
703 # ifdef ABUTTING1TO1_FACES
704 
705         nelem = (NUM_SIDE - 1) * (NUM_SIDE - 1);
706 
707         /* zone 1 to zone 2 connectivity as Abutting1to1 with element range */
708 
709         range[0] = num_element + num_face - nelem + 1;
710         range[1] = num_element + num_face;
711         for (n = 0; n < nelem; n++)
712             d_pts[n] = range[0] + n - nelem;
713         location = FaceCenter;
714 
715         if (cg_conn_write(cgfile, cgbase, 1, name,
716                 location, CGNS_ENUMV(Abutting1to1),
717                 CGNS_ENUMV(PointRange), 2, range, "Zone2",
718                 CGNS_ENUMV(Unstructured), CGNS_ENUMV(PointListDonor),
719                 CGNS_ENUMV(Integer), nelem, d_pts, &cgconn)) {
720             sprintf(errmsg, "Zone1/%s/%s", zcname, name);
721             error_exit(errmsg);
722         }
723 
724         /* zone 2 to zone 1 connectivity as Abutting1to1 with element list */
725 
726         for (n = 0; n < nelem; n++) {
727             pts[n]   = num_element + num_face - 2 * nelem + 1 + n;
728             d_pts[n] = pts[n] + nelem;
729         }
730         if (cg_conn_write(cgfile, cgbase, 2, name,
731                 location, CGNS_ENUMV(Abutting1to1),
732                 CGNS_ENUMV(PointList), nelem, pts, "Zone1",
733                 CGNS_ENUMV(Unstructured), CGNS_ENUMV(PointListDonor),
734                 CGNS_ENUMV(Integer), nelem, d_pts, &cgconn)) {
735             sprintf(errmsg, "Zone2/%s/%s", zcname, name);
736             error_exit(errmsg);
737         }
738 
739 # else
740 
741         /* zone 1 to zone 2 connectivity as Abutting1to1 with point range */
742 
743         range[0] = NODE_INDEX(1, 1, NUM_SIDE);
744         range[1] = NODE_INDEX(NUM_SIDE, NUM_SIDE, NUM_SIDE);
745         for (n = 0, j = 1; j <= NUM_SIDE; j++) {
746             for (i = 1; i <= NUM_SIDE; i++)
747                 d_pts[n++] = NODE_INDEX(i, j, 1);
748         }
749         if (cg_conn_write(cgfile, cgbase, 1, name,
750                 CGNS_ENUMV(Vertex), CGNS_ENUMV(Abutting1to1),
751                 CGNS_ENUMV(PointRange), 2, range, "Zone2",
752                 CGNS_ENUMV(Unstructured), CGNS_ENUMV(PointListDonor),
753                 CGNS_ENUMV(Integer), npts, d_pts, &cgconn)) {
754             sprintf(errmsg, "Zone1/%s/%s", zcname, name);
755             error_exit(errmsg);
756         }
757 
758         /* zone 2 to zone 1 connectivity as Abutting1to1 with point list */
759 
760         for (n = 0, j = 1; j <= NUM_SIDE; j++) {
761             for (i = 1; i <= NUM_SIDE; i++) {
762                 pts[n] = d_pts[n];
763                 d_pts[n++] = NODE_INDEX(i, j, NUM_SIDE);
764             }
765         }
766         if (cg_conn_write(cgfile, cgbase, 2, name,
767                 CGNS_ENUMV(Vertex), CGNS_ENUMV(Abutting1to1),
768                 CGNS_ENUMV(PointList), npts, pts, "Zone1",
769                 CGNS_ENUMV(Unstructured), CGNS_ENUMV(PointListDonor),
770                 CGNS_ENUMV(Integer), npts, d_pts, &cgconn)) {
771             sprintf(errmsg, "Zone2/%s/%s", zcname, name);
772             error_exit(errmsg);
773         }
774 
775         sprintf(name, "hole%d", nc);
776 
777         /* write hole in zone 1 as PointRange */
778 
779         if (cg_hole_write(cgfile, cgbase, 1, name, CGNS_ENUMV(Vertex),
780                 CGNS_ENUMV(PointRange), 1, 2, range, &cghole)) {
781             sprintf(errmsg, "Zone1/%s/%s", zcname, name);
782             error_exit(errmsg);
783         }
784 
785         /* write hole in zone 2 as PointList */
786 
787         if (cg_hole_write(cgfile, cgbase, 2, name, CGNS_ENUMV(Vertex),
788                 CGNS_ENUMV(PointList), 1, npts, pts, &cghole)) {
789             sprintf(errmsg, "Zone2/%s/%s", zcname, name);
790             error_exit(errmsg);
791         }
792 
793 # endif
794 #endif
795     }
796 
797     /* create BaseIterativeData_t node */
798 
799 
800     for (n = 0; n < NUM_ZCONN; n++)
801         iter[n] = n + 1;
802     dims[0] = NUM_ZCONN;
803 
804     if (cg_biter_write(cgfile, cgbase, "BaseIterativeData", NUM_ZCONN))
805         error_exit("cg_biter_write");
806     if (cg_goto(cgfile, cgbase, "BaseIterativeData", 0, NULL))
807         error_exit("cg_goto BaseIterativeData");
808     if (cg_array_write("IterationValues", CGNS_ENUMV(Integer),
809         1, dims, iter)) error_exit("cg_array_write IterationValues");
810 
811     /* create the ZoneIterativeData_t node */
812 
813     dims[0] = 32;
814     dims[1] = NUM_ZCONN;
815     n = 32 * NUM_ZCONN;
816     for (i = 0; i < n; i++)
817         pointers[i] = ' ';
818     pointers[n] = 0;
819 
820     for (n = 0, nc = 1; nc <= NUM_ZCONN; nc++, n+=32) {
821         sprintf(zcname, "%s%d", ZCONN_NAME, nc);
822         i = strlen(zcname);
823         strcpy(&pointers[n], zcname);
824         pointers[n+i] = ' ';
825     }
826 
827     for (n = 1; n <= 2; n++) {
828         sprintf(name, "Zone%dInterativeData", n);
829         if (cg_ziter_write(cgfile, cgbase, n, name)) {
830             sprintf(errmsg, "%s:cg_ziter_write", name);
831             error_exit(errmsg);
832         }
833         if (cg_goto(cgfile, cgbase, "Zone_t", n, name, 0, NULL)) {
834             sprintf(errmsg, "%s:cg_goto", name);
835             error_exit(errmsg);
836         }
837         if (cg_array_write("ZoneGridConnectivityPointers",
838                 CGNS_ENUMV(Character), 2, dims, pointers)) {
839             sprintf(errmsg, "%s:cg_array_write", name);
840             error_exit(errmsg);
841         }
842     }
843 }
844 
845 /*------------------------------------------------------------------------
846  * test Zone connectivity data
847  *------------------------------------------------------------------------*/
848 
set_zconn(const char * zcname)849 void set_zconn(const char *zcname)
850 {
851     int n;
852     char name[33];
853 
854     for (n = 1; n <= NUM_ZCONN; n++) {
855         if (cg_zconn_read(cgfile, cgbase, cgzone, n, name)) {
856             sprintf(errmsg, "%s:cg_zconn_read", zcname);
857             error_exit(errmsg);
858         }
859         if (0 == strcmp(name, zcname)) {
860             if (cg_zconn_set(cgfile, cgbase, cgzone, n)) {
861                 sprintf(errmsg, "%s:cg_zconn_set", zcname);
862                 error_exit(errmsg);
863             }
864             return;
865         }
866     }
867     sprintf(errmsg, "zconn %s not found", zcname);
868     error_exit(errmsg);
869 }
870 
test_zconn()871 void test_zconn()
872 {
873     int n, nz, nc;
874     int nconn, n1to1, nholes;
875     char name[33], zcname[33], expected[33], pointers[32*NUM_ZCONN+1];
876     CGNS_ENUMT(GridLocation_t) location;
877     CGNS_ENUMT(PointSetType_t) ptype;
878     cgsize_t np;
879 
880     if (cg_nzconns(cgfile, cgbase, cgzone, &nz) || nz != NUM_ZCONN)
881         error_exit("cg_nzconns");
882     if (cg_ziter_read(cgfile, cgbase, cgzone, name))
883         error_exit("cg_ziter_read");
884     if (cg_goto(cgfile, cgbase, "Zone_t", cgzone, name, 0, NULL))
885         error_exit("cg_goto ziter");
886     if (cg_array_read(1, pointers)) error_exit("cg_array_read");
887 
888     for (nz = 0; nz < NUM_ZCONN; nz++) {
889         strncpy(zcname, &pointers[32*nz], 32);
890         zcname[32] = 0;
891         for (n = strlen(zcname)-1; n >= 0 && zcname[n] == ' '; n--)
892             ;
893         zcname[++n] = 0;
894         set_zconn(zcname);
895         n = strlen(ZCONN_NAME);
896         nc = atoi(&zcname[n]);
897 
898         if (cg_nconns(cgfile, cgbase, cgzone, &nconn)) {
899             sprintf(errmsg, "%s:cg_nconns", zcname);
900             error_exit(errmsg);
901         }
902         if (cg_n1to1(cgfile, cgbase, cgzone, &n1to1)) {
903             sprintf(errmsg, "%s:cg_n1to1", zcname);
904             error_exit(errmsg);
905         }
906         if ((nconn + n1to1) != 1) {
907             sprintf(errmsg, "%s:nconns=%d", zcname, (nconn + n1to1));
908             error_exit(errmsg);
909         }
910         if (n1to1) {
911             cgsize_t range[6], d_range[6];
912             int transform[3];
913             if (cg_1to1_read(cgfile, cgbase, cgzone, 1, name, expected,
914                     range, d_range, transform)) {
915                 sprintf(errmsg, "%s:cg_1to1_read", zcname);
916                 error_exit(errmsg);
917             }
918         }
919         else {
920             CGNS_ENUMT(GridConnectivityType_t) type;
921             CGNS_ENUMT(PointSetType_t) d_ptype;
922             CGNS_ENUMT(ZoneType_t) d_ztype;
923             CGNS_ENUMT(DataType_t) d_dtype;
924             cgsize_t d_npnts;
925            if (cg_conn_info(cgfile, cgbase, cgzone, 1, name, &location,
926                     &type, &ptype, &np, expected, &d_ztype, &d_ptype,
927                     &d_dtype, &d_npnts)) {
928                 sprintf(errmsg, "%s:cg_conn_info", zcname);
929                 error_exit(errmsg);
930             }
931         }
932         sprintf(expected, "conn%d", nc);
933         if (strcmp(expected, name)) {
934             sprintf(errmsg, "%s:conn name %s != %s", zcname, name, expected);
935             error_exit(errmsg);
936         }
937 
938         if (cg_nholes(cgfile, cgbase, cgzone, &nholes)) {
939             sprintf(errmsg, "%s:cg_zconn_set", zcname);
940             error_exit(errmsg);
941         }
942         if (nholes != 1) {
943             sprintf(errmsg, "%s:nholes=%d", zcname, nholes);
944             error_exit(errmsg);
945         }
946         if (cg_hole_info(cgfile, cgbase, cgzone, 1, name, &location,
947                 &ptype, &n, &np)) {
948             sprintf(errmsg, "%s:cg_hole_info", zcname);
949             error_exit(errmsg);
950         }
951         sprintf(expected, "hole%d", nc);
952         if (strcmp(expected, name)) {
953             sprintf(errmsg, "%s:hole name %s != %s", zcname, name, expected);
954             error_exit(errmsg);
955         }
956     }
957 }
958 
959 /*------------------------------------------------------------------------
960  * test sub region data
961  *------------------------------------------------------------------------*/
962 
test_subreg()963 void test_subreg()
964 {
965     int ns, dim, rind[6];
966     int bclen, gclen, ndescr;
967     char srname[33], name[33], text[65], *descr;
968     CGNS_ENUMT(GridLocation_t) location;
969     CGNS_ENUMT(PointSetType_t) ptype;
970     CGNS_ENUMT(DataClass_t) dclass;
971     CGNS_ENUMT(MassUnits_t) mass;
972     CGNS_ENUMT(LengthUnits_t) length;
973     CGNS_ENUMT(TimeUnits_t) time;
974     CGNS_ENUMT(TemperatureUnits_t) temp;
975     CGNS_ENUMT(AngleUnits_t) ang;
976     cgsize_t np;
977 
978     if (cg_nsubregs(cgfile, cgbase, cgzone, &ns))
979         error_exit("cg_nsubregs");
980     if (ns == 0 && cgbase == 2) return;
981     if (ns != NUM_ZCONN) error_exit("invalid nsubregs");
982 
983     for (ns = 1; ns < NUM_ZCONN; ns++) {
984         if (cg_subreg_info(cgfile, cgbase, cgzone, ns, srname, &dim,
985                 &location, &ptype, &np, &bclen, &gclen))
986             error_exit("cg_subreg_info");
987         if (dim != 2) error_exit("dim not 2");
988 
989         if (bclen) {
990             if (gclen || ptype != CGNS_ENUMV(PointSetTypeNull) || np)
991                 error_exit("bclen");
992             if (strcmp(srname, "SubRegion1")) error_exit("SubRegion1");
993             /* these should be errors */
994             if (!cg_subreg_gcname_read(cgfile, cgbase, cgzone, ns, text) ||
995                 !cg_subreg_ptset_read(cgfile, cgbase, cgzone, ns, pts))
996                 error_exit("cg_subreg_bcname_read:multiple");
997             /* valid */
998             if (cg_subreg_bcname_read(cgfile, cgbase, cgzone, ns, text))
999                 error_exit("gc_subreg_bcname_read");
1000 
1001             if (cgzone == 1) {
1002                 if (strcmp(text, "Inlet")) error_exit("Inlet");
1003                 if (cg_goto(cgfile, cgbase, "Zone_t", 1, "ZoneSubRegion_t", ns, NULL))
1004                     error_exit("cg_goto ZoneSubRegion_t");
1005                 if (cg_dataclass_read(&dclass) ||
1006                     cg_units_read(&mass, &length, &time, &temp, &ang) ||
1007                     cg_ndescriptors(&ndescr))
1008                     error_exit("dataclass,units or ndescriptors");
1009                 if (ndescr != 1) error_exit("ndescriptors not 1");
1010                 if (cg_descriptor_read(1, name, &descr))
1011                     error_exit("cg_descriptor_read");
1012                 if (strcmp(name, "descriptor") ||
1013                     strcmp(descr, "this is zone 1, subregion 1"))
1014                     error_exit("descriptor");
1015                 cg_free(descr);
1016             }
1017             else {
1018                 if (strcmp(text, "Outlet")) error_exit("Outlet");
1019                 if (cg_goto(cgfile, cgbase, "Zone_t", 2, srname, 0, NULL))
1020                     error_exit("cg_goto SubRegion1");
1021                 if (cg_famname_read(name) || cg_user_data_read(1, text))
1022                     error_exit("famname or user data");
1023                 if (strcmp(name, "Family")) error_exit("Family");
1024                 if (strcmp(text, "Userdata")) error_exit("Userdata");
1025                 if (cg_gorel(cgfile, "UserDefinedData_t", 1, NULL))
1026                     error_exit("cg_gorel SubRegion1/UserDefinedData");
1027                 if (cg_famname_read(name) || strcmp(name, "Family"))
1028                     error_exit("Userdata/Family");
1029             }
1030         }
1031 
1032         else if (gclen) {
1033             if (bclen || ptype != CGNS_ENUMV(PointSetTypeNull) || np)
1034                 error_exit("gclen");
1035             if (strcmp(srname, "SubRegion2")) error_exit("SubRegion2");
1036             /* these should be errors */
1037             if (!cg_subreg_bcname_read(cgfile, cgbase, cgzone, ns, text) ||
1038                 !cg_subreg_ptset_read(cgfile, cgbase, cgzone, ns, pts))
1039                 error_exit("cg_subreg_bcname_read:multiple");
1040             /* valid */
1041             if (cg_subreg_gcname_read(cgfile, cgbase, cgzone, ns, text))
1042                 error_exit("cg_subreg_gcname_read");
1043             if (strcmp(text, "ZoneConnectivity2/conn2"))
1044                 error_exit("ZoneConnectivity2/conn2");
1045         }
1046 
1047         else {
1048             if (bclen || gclen) error_exit("ptset");
1049             if (strcmp(srname, "SubRegion3") &&
1050                 strcmp(srname, "SubRegion4") &&
1051                 strcmp(srname, "SubRegion5")) error_exit("SubRegion[345]");
1052             /* these should be errors */
1053             if (!cg_subreg_bcname_read(cgfile, cgbase, cgzone, ns, text) ||
1054                 !cg_subreg_gcname_read(cgfile, cgbase, cgzone, ns, text))
1055                 error_exit("cg_subreg_ptset_read:multiple");
1056 
1057             if (cgzone == 1) {
1058                 if (location != CGNS_ENUMV(Vertex) ||
1059                     ptype != CGNS_ENUMV(PointRange) || np != 2)
1060                     error_exit("not Vertex,PointRange or 2");
1061                 if (cg_goto(cgfile, cgbase, "Zone_t", 1, srname, 0, NULL))
1062                     error_exit("cg_goto SubRegion[345]");
1063                 if (cg_rind_read(rind)) error_exit("cg_rind_read");
1064             }
1065             else {
1066                 if (location != CGNS_ENUMV(KFaceCenter) ||
1067                     ptype != CGNS_ENUMV(PointList) ||
1068                     np != (NUM_SIDE-1)*(NUM_SIDE-1))
1069                     error_exit("not EdgeCenter,PointList or npts");
1070             }
1071             if (cg_subreg_ptset_read(cgfile, cgbase, cgzone, ns, pts))
1072                 error_exit("cg_subreg_ptset_read");
1073         }
1074     }
1075 }
1076 
1077