1 /*
2 Copyright (c) 1994 - 2010, Lawrence Livermore National Security, LLC.
3 LLNL-CODE-425250.
4 All rights reserved.
5
6 This file is part of Silo. For details, see silo.llnl.gov.
7
8 Redistribution and use in source and binary forms, with or without
9 modification, are permitted provided that the following conditions
10 are met:
11
12 * Redistributions of source code must retain the above copyright
13 notice, this list of conditions and the disclaimer below.
14 * Redistributions in binary form must reproduce the above copyright
15 notice, this list of conditions and the disclaimer (as noted
16 below) in the documentation and/or other materials provided with
17 the distribution.
18 * Neither the name of the LLNS/LLNL nor the names of its
19 contributors may be used to endorse or promote products derived
20 from this software without specific prior written permission.
21
22 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
23 "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
24 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
25 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL LAWRENCE
26 LIVERMORE NATIONAL SECURITY, LLC, THE U.S. DEPARTMENT OF ENERGY OR
27 CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
28 EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
29 PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
30 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
31 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
32 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
33 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
34
35 This work was produced at Lawrence Livermore National Laboratory under
36 Contract No. DE-AC52-07NA27344 with the DOE. Neither the United
37 States Government nor Lawrence Livermore National Security, LLC nor
38 any of their employees, makes any warranty, express or implied, or
39 assumes any liability or responsibility for the accuracy,
40 completeness, or usefulness of any information, apparatus, product, or
41 process disclosed, or represents that its use would not infringe
42 privately-owned rights. Any reference herein to any specific
43 commercial products, process, or services by trade name, trademark,
44 manufacturer or otherwise does not necessarily constitute or imply its
45 endorsement, recommendation, or favoring by the United States
46 Government or Lawrence Livermore National Security, LLC. The views and
47 opinions of authors expressed herein do not necessarily state or
48 reflect those of the United States Government or Lawrence Livermore
49 National Security, LLC, and shall not be used for advertising or
50 product endorsement purposes.
51 */
52 #include <silo.h>
53 #include <stdio.h>
54 #include <math.h>
55 #include <string.h>
56 #include <std.c>
57
58 #define IND(i,j) i-1][j-1
59
60 #define matrix_assign(matrix,a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34,a41,a42,a43,a44) \
61 { \
62 matrix [IND(1,1)] = a11 ; \
63 matrix [IND(1,2)] = a12 ; \
64 matrix [IND(1,3)] = a13 ; \
65 matrix [IND(1,4)] = a14 ; \
66 matrix [IND(2,1)] = a21 ; \
67 matrix [IND(2,2)] = a22 ; \
68 matrix [IND(2,3)] = a23 ; \
69 matrix [IND(2,4)] = a24 ; \
70 matrix [IND(3,1)] = a31 ; \
71 matrix [IND(3,2)] = a32 ; \
72 matrix [IND(3,3)] = a33 ; \
73 matrix [IND(3,4)] = a34 ; \
74 matrix [IND(4,1)] = a41 ; \
75 matrix [IND(4,2)] = a42 ; \
76 matrix [IND(4,3)] = a43 ; \
77 matrix [IND(4,4)] = a44 ; \
78 }
79
80 #define matrix_mult(matrixa, matrixb, matrixc) \
81 { \
82 for (i = 1; i < 5; i++) { \
83 for (j = 1; j < 5; j++) { \
84 matrixc [IND(i,j)] = matrixa [IND(i,1)] * matrixb [IND(1,j)] + \
85 matrixa [IND(i,2)] * matrixb [IND(2,j)] + \
86 matrixa [IND(i,3)] * matrixb [IND(3,j)] + \
87 matrixa [IND(i,4)] * matrixb [IND(4,j)] ; \
88 } \
89 } \
90 }
91 #ifndef M_PI /* yea, Solaris 5 */
92 #define M_PI 3.14159265358979323846264338327950288 /* pi */
93 #endif
94 #define RAD(deg) M_PI*(deg/180.0)
95
96
97 /*-------------------------------------------------------------------------
98 * Function: main
99 *
100 * Purpose:
101 *
102 * Return: 0
103 *
104 * Programmer:
105 *
106 * Modifications:
107 * Robb Matzke, 1999-04-09
108 * Added argument parsing to control the driver which is used.
109 *
110 * Mark C. Miller, Mon Jan 11 16:27:38 PST 2010
111 * Added missing call to DBFreeFacelist.
112 *
113 * Mark C. Miller, Tue Sep 14 20:50:10 PDT 2010
114 * Fixed uninitialized write for Cvar by assigning a value. Fixed
115 * another missing call to DBFreeFacelist.
116 *-------------------------------------------------------------------------
117 */
118 int
main(int argc,char * argv[])119 main(int argc, char *argv[])
120 {
121
122 DBfile *dbfile = NULL;
123 char *coordnames[3];
124 float *coords[3];
125 int Pnodelist[16];
126 int Cnodelist[8];
127 float x[12], y[12], z[12];
128 int Pshapesize[1];
129 int Pshapecnt[1];
130 int Cshapesize[1];
131 int Cshapecnt[1];
132 DBfacelist *Pfacelist = NULL;
133 DBfacelist *Cfacelist = NULL;
134 int i, j, len;
135 char mesh_command[256];
136 float rot1[4][4], rot2[4][4], final[4][4];
137 float angle;
138 float Cvar[1] = {3.141592678};
139 float Pvar[12];
140 int driver=DB_PDB;
141 char *filename="subhex.silo";
142 int show_all_errors = FALSE;
143
144 /* Parse command-line */
145 for (i=1; i<argc; i++) {
146 if (!strncmp(argv[i], "DB_PDB", 6)) {
147 driver = StringToDriver(argv[i]);
148 filename = "subhex.pdb";
149 } else if (!strncmp(argv[i], "DB_HDF5", 7)) {
150 driver = StringToDriver(argv[i]);
151 filename = "subhex.h5";
152 } else if (!strcmp(argv[i], "show-all-errors")) {
153 show_all_errors = 1;
154 } else if (argv[i][0] != '\0') {
155 fprintf(stderr, "%s: ignored argument `%s'\n", argv[0], argv[i]);
156 }
157 }
158
159 DBShowErrors(show_all_errors?DB_ALL_AND_DRVR:DB_ABORT, NULL);
160 printf("Creating test file \"%s\".\n", filename);
161 dbfile = DBCreate(filename, DB_CLOBBER, DB_LOCAL, "3D ucd hex", driver);
162
163 coordnames[0] = "xcoords";
164 coordnames[1] = "ycoords";
165 coordnames[2] = "zcoords";
166
167 /*
168 * The coordinates of the 12 nodes...
169 */
170
171 x[0] = 0; y[0] = 0; z[0] = 0;
172 x[1] = 1; y[1] = 0; z[1] = 0;
173 x[2] = 1; y[2] = 0; z[2] = 1;
174 x[3] = 0; y[3] = 0; z[3] = 1;
175 x[4] = 0; y[4] = 1; z[4] = 0;
176 x[5] = 1; y[5] = 1; z[5] = 0;
177 x[6] = 1; y[6] = 1; z[6] = 1;
178 x[7] = 0; y[7] = 1; z[7] = 1;
179
180 x[8] = 2; y[8] = 0; z[8] = 0;
181 x[9] = 2; y[9] = 0; z[9] = 1;
182 x[10] = 2; y[10] = 1; z[10] = 0;
183 x[11] = 2; y[11] = 1; z[11] = 1;
184
185 coords[0] = x;
186 coords[1] = y;
187 coords[2] = z;
188
189 angle = 45;
190 angle = M_PI*(angle/180.0);
191 matrix_assign(rot1,
192 1, 0, 0, 0,
193 0, cos(angle), sin(angle), 0,
194 0, -sin(angle), cos(angle), 0,
195 0, 0, 0, 1);
196 matrix_assign(rot2,
197 cos(angle), 0, -sin(angle), 0,
198 0, 1, 0, 0,
199 sin(angle), 0, cos(angle), 0,
200 0, 0, 0, 1);
201 matrix_mult(rot1, rot2, final);
202
203 /*--------------------------------------
204 |
205 | "All" nodes...
206 |
207 +-------------------------------------*/
208
209 for (i = 0; i < 12; i++)
210 {
211 float tx, ty,tz;
212
213 tx = x[i]*final[IND(1,1)] + y[i]*final[IND(1,2)] + z[i]*final[IND(1,3)] + final[IND(1,4)];
214 ty = x[i]*final[IND(2,1)] + y[i]*final[IND(2,2)] + z[i]*final[IND(2,3)] + final[IND(2,4)];
215 tz = x[i]*final[IND(3,1)] + y[i]*final[IND(3,2)] + z[i]*final[IND(3,3)] + final[IND(3,4)];
216
217 x[i] = tx;
218 y[i] = ty;
219 z[i] = tz;
220
221 Pvar[i] = x[i]+y[i]*z[i];
222 }
223
224 /*--------------------------------------
225 |
226 | The parent mesh...
227 |
228 +-------------------------------------*/
229
230 DBPutUcdmesh(dbfile, "parent", 3, (DBCAS_t) coordnames,
231 coords, 12, 2, "Pzonelist", "Pfacelist", DB_FLOAT, NULL);
232
233 DBPutUcdvar1(dbfile, "v", "parent", Pvar, 12, NULL, 0, DB_FLOAT, DB_NODECENT, NULL);
234
235 Pnodelist[ 0] = 0; /* The first hex */
236 Pnodelist[ 1] = 1;
237 Pnodelist[ 2] = 2;
238 Pnodelist[ 3] = 3;
239 Pnodelist[ 4] = 4;
240 Pnodelist[ 5] = 5;
241 Pnodelist[ 6] = 6;
242 Pnodelist[ 7] = 7;
243
244 Pnodelist[ 8] = 1; /* The second hex */
245 Pnodelist[ 9] = 8;
246 Pnodelist[10] = 9;
247 Pnodelist[11] = 2;
248 Pnodelist[12] = 5;
249 Pnodelist[13] = 10;
250 Pnodelist[14] = 11;
251 Pnodelist[15] = 6;
252
253 Pshapecnt[0] = 2; /* There are 2... */
254 Pshapesize[0] = 8; /* ...hexes */
255
256 DBSetDeprecateWarnings(0);
257 DBPutZonelist(dbfile, "Pzonelist", 2, 3, Pnodelist, 12, 0, Pshapesize,
258 Pshapecnt, 1);
259 DBSetDeprecateWarnings(3);
260
261 Pfacelist = DBCalcExternalFacelist(Pnodelist, 12, 0, Pshapesize, Pshapecnt, 1,
262 NULL, 0);
263
264 DBPutFacelist(dbfile, "Pfacelist", Pfacelist->nfaces, Pfacelist->ndims,
265 Pfacelist->nodelist, Pfacelist->lnodelist, Pfacelist->origin,
266 Pfacelist->zoneno, Pfacelist->shapesize, Pfacelist->shapecnt,
267 Pfacelist->nshapes, Pfacelist->types, Pfacelist->typelist,
268 Pfacelist->ntypes);
269
270 DBFreeFacelist(Pfacelist);
271
272 sprintf(mesh_command, "mesh hex; contour v");
273 len = strlen(mesh_command) + 1;
274 DBWrite(dbfile, "_meshtvinfo", mesh_command, &len, 1, DB_CHAR);
275
276 /*--------------------------------------
277 |
278 | The subset...
279 |
280 +-------------------------------------*/
281
282 DBSetDeprecateWarnings(0);
283 DBPutUcdsubmesh(dbfile, "child", "parent", 1, "Czonelist", "Cfacelist", NULL);
284 DBSetDeprecateWarnings(3);
285
286 Cnodelist[ 0] = 1; /* Just one hex, refering to parent nodes */
287 Cnodelist[ 1] = 8;
288 Cnodelist[ 2] = 9;
289 Cnodelist[ 3] = 2;
290 Cnodelist[ 4] = 5;
291 Cnodelist[ 5] = 10;
292 Cnodelist[ 6] = 11;
293 Cnodelist[ 7] = 6;
294
295 Cshapecnt[0] = 1; /* There is 1... */
296 Cshapesize[0] = 8; /* ...hex */
297
298 DBSetDeprecateWarnings(0);
299 DBPutZonelist(dbfile, "Czonelist", 1, 3, Cnodelist, 8, 0, Cshapesize,
300 Cshapecnt, 1);
301 DBSetDeprecateWarnings(3);
302
303 Cfacelist = DBCalcExternalFacelist(Cnodelist, 8, 0, Cshapesize, Cshapecnt, 1,
304 NULL, 0);
305
306 DBPutFacelist(dbfile, "Cfacelist", Cfacelist->nfaces, Cfacelist->ndims,
307 Cfacelist->nodelist, Cfacelist->lnodelist, Cfacelist->origin,
308 Cfacelist->zoneno, Cfacelist->shapesize, Cfacelist->shapecnt,
309 Cfacelist->nshapes, Cfacelist->types, Cfacelist->typelist,
310 Cfacelist->ntypes);
311
312 DBFreeFacelist(Cfacelist);
313
314 DBPutUcdvar1(dbfile, "u", "child", Cvar, 1, NULL, 0, DB_FLOAT, DB_ZONECENT,
315 NULL);
316
317 DBClose(dbfile);
318
319 CleanupDriverStuff();
320 return (0);
321 }
322