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
53 #include <silo.h>
54 #include <stdio.h>
55 #include <math.h>
56 #include <string.h>
57 #include <std.c>
58
59 #define IND(i,j) i-1][j-1
60
61 #define matrix_assign(matrix,a11,a12,a13,a14,a21,a22,a23,a24,a31,a32,a33,a34,a41,a42,a43,a44) \
62 { \
63 matrix [IND(1,1)] = a11 ; \
64 matrix [IND(1,2)] = a12 ; \
65 matrix [IND(1,3)] = a13 ; \
66 matrix [IND(1,4)] = a14 ; \
67 matrix [IND(2,1)] = a21 ; \
68 matrix [IND(2,2)] = a22 ; \
69 matrix [IND(2,3)] = a23 ; \
70 matrix [IND(2,4)] = a24 ; \
71 matrix [IND(3,1)] = a31 ; \
72 matrix [IND(3,2)] = a32 ; \
73 matrix [IND(3,3)] = a33 ; \
74 matrix [IND(3,4)] = a34 ; \
75 matrix [IND(4,1)] = a41 ; \
76 matrix [IND(4,2)] = a42 ; \
77 matrix [IND(4,3)] = a43 ; \
78 matrix [IND(4,4)] = a44 ; \
79 }
80
81 #define matrix_mult(matrixa, matrixb, matrixc) \
82 { \
83 for (i = 1; i < 5; i++) { \
84 for (j = 1; j < 5; j++) { \
85 matrixc [IND(i,j)] = matrixa [IND(i,1)] * matrixb [IND(1,j)] + \
86 matrixa [IND(i,2)] * matrixb [IND(2,j)] + \
87 matrixa [IND(i,3)] * matrixb [IND(3,j)] + \
88 matrixa [IND(i,4)] * matrixb [IND(4,j)] ; \
89 } \
90 } \
91 }
92 #ifndef M_PI /* yea, Solaris 5 */
93 #define M_PI 3.14159265358979323846264338327950288 /* pi */
94 #endif
95 #define RAD(deg) M_PI*(deg/180.0)
96
97 int
main(int argc,char * argv[])98 main(int argc, char *argv[])
99 {
100 DBfile *dbfile = NULL;
101 char *coordnames[3];
102 float *coords[3];
103 int nodelist[4];
104 float x[4], y[4], z[4];
105 int shapesize[1];
106 int shapecnt[1];
107 DBfacelist *facelist = NULL;
108 int matnos[1], matlist[1], dims[1];
109 int i, j, len;
110 char mesh_command[256];
111 float rot1[4][4], rot2[4][4], final[4][4];
112 float angle;
113 float var[4];
114 int driver=DB_PDB;
115 char *filename = "onetet.silo";
116 int show_all_errors = FALSE;
117
118 for (i=1; i<argc; i++) {
119 if (!strncmp(argv[i], "DB_PDB", 6)) {
120 driver = StringToDriver(argv[i]);
121 filename = "onetet.pdb";
122 } else if (!strncmp(argv[i], "DB_HDF5", 7)) {
123 driver = StringToDriver(argv[i]);
124 filename = "onetet.h5";
125 } else if (!strcmp(argv[i], "show-all-errors")) {
126 show_all_errors = 1;
127 } else if (argv[i][0] != '\0') {
128 fprintf(stderr, "%s: ignored argument `%s'\n", argv[0], argv[i]);
129 }
130 }
131
132
133
134 DBShowErrors(show_all_errors?DB_ALL_AND_DRVR:DB_ABORT, NULL);
135 printf("Creating test file \"%s\".\n", filename);
136 dbfile = DBCreate(filename, DB_CLOBBER, DB_LOCAL, "3D ucd tet", driver);
137
138 coordnames[0] = "xcoords";
139 coordnames[1] = "ycoords";
140 coordnames[2] = "zcoords";
141
142 x[0] = 0; y[0] = 0; z[0] = 0;
143 x[1] = 1; y[1] = 0; z[1] = 0;
144 x[2] = 0.5; y[2] = 0; z[2] = 1;
145 x[3] = 0.5; y[3] = 1; z[3] = 1;
146
147 coords[0] = x;
148 coords[1] = y;
149 coords[2] = z;
150
151 angle = 45;
152 angle = M_PI*(angle/180.0);
153 matrix_assign(rot1,
154 1, 0, 0, 0,
155 0, cos(angle), sin(angle), 0,
156 0, -sin(angle), cos(angle), 0,
157 0, 0, 0, 1);
158 matrix_assign(rot2,
159 cos(angle), 0, -sin(angle), 0,
160 0, 1, 0, 0,
161 sin(angle), 0, cos(angle), 0,
162 0, 0, 0, 1);
163 matrix_mult(rot1, rot2, final);
164
165 for (i = 0; i < 4; i++)
166 {
167 float tx, ty,tz;
168
169 tx = x[i]*final[IND(1,1)] + y[i]*final[IND(1,2)] + z[i]*final[IND(1,3)] + final[IND(1,4)];
170 ty = x[i]*final[IND(2,1)] + y[i]*final[IND(2,2)] + z[i]*final[IND(2,3)] + final[IND(2,4)];
171 tz = x[i]*final[IND(3,1)] + y[i]*final[IND(3,2)] + z[i]*final[IND(3,3)] + final[IND(3,4)];
172
173 x[i] = tx;
174 y[i] = ty;
175 z[i] = tz;
176
177 var[i] = x[i]+y[i]*z[i];
178 }
179
180 DBPutUcdmesh(dbfile, "tet", 3, (DBCAS_t) coordnames,
181 coords, 4, 1, "zonelist", "facelist", DB_FLOAT, NULL);
182
183 matnos[0] = 1;
184 matlist[0] = 1;
185 dims[0] = 1;
186
187 DBPutMaterial(dbfile, "mat", "tet", 1, matnos, matlist, dims,
188 1, NULL, NULL, NULL, NULL, 0, DB_FLOAT, NULL);
189
190 DBPutUcdvar1(dbfile, "v", "tet", var, 4, NULL, 0, DB_FLOAT,
191 DB_NODECENT, NULL);
192
193 nodelist[0] = 0;
194 nodelist[1] = 1;
195 nodelist[2] = 2;
196 nodelist[3] = 3;
197
198 shapesize[0] = 4;
199 shapecnt[0] = 1;
200
201 DBSetDeprecateWarnings(0);
202 DBPutZonelist(dbfile, "zonelist", 1, 3, nodelist, 4, 0, shapesize,
203 shapecnt, 1);
204 DBSetDeprecateWarnings(3);
205
206 facelist = DBCalcExternalFacelist(nodelist, 4, 0, shapesize, shapecnt, 1,
207 NULL, 0);
208
209 DBPutFacelist(dbfile, "facelist", facelist->nfaces, facelist->ndims,
210 facelist->nodelist, facelist->lnodelist, facelist->origin,
211 facelist->zoneno, facelist->shapesize, facelist->shapecnt,
212 facelist->nshapes, facelist->types, facelist->typelist,
213 facelist->ntypes);
214
215 sprintf(mesh_command, "mesh tet; contour v");
216 len = strlen(mesh_command) + 1;
217 DBWrite(dbfile, "_meshtvinfo", mesh_command, &len, 1, DB_CHAR);
218
219 DBClose(dbfile);
220
221 DBFreeFacelist(facelist);
222
223 CleanupDriverStuff();
224 return (0);
225 }
226