Branch data Line data Source code
1 : : //- Class: Cholla
2 : : //- Description: C Interface for the Cholla module
3 : : //- Facet-based geometry definition and evaluation
4 : : //- Owner: Steven J. Owen
5 : : //- Checked by:
6 : : //- Version:
7 : :
8 : : #include "time.h"
9 : : #include "Cholla.h"
10 : : #include "ChollaEngine.hpp"
11 : : #include "CubitPoint.hpp"
12 : : #include "CubitPointData.hpp"
13 : : #include "CubitFacetEdge.hpp"
14 : : #include "CubitFacetEdgeData.hpp"
15 : : #include "CubitFacet.hpp"
16 : : #include "CubitFacetData.hpp"
17 : : #include "CubitQuadFacet.hpp"
18 : : #include "CubitQuadFacetData.hpp"
19 : : #include "CastTo.hpp"
20 : : #include "TDFacetBoundaryPoint.hpp"
21 : : #include "FacetEvalTool.hpp"
22 : : #include "FacetDataUtil.hpp"
23 : :
24 : : #include <iostream>
25 : : #include <fstream>
26 : : using std::ifstream;
27 : :
28 : : #define END_OF_BLOCKS "END_OF_BLOCKS"
29 : :
30 : : // Internal Static function prototypes
31 : :
32 : : static void time_stamp( FILE *fp );
33 : : static CubitStatus build_facets(int numTri, int numQuad, int numEdge,
34 : : int numVert, int* triEdge,
35 : : int *quadEdge,
36 : : int* edgeVert, double* vert,
37 : : CubitPoint **point_array,
38 : : CubitFacetEdge **edge_array,
39 : : CubitFacet **facet_array,
40 : : CubitQuadFacet **qfacet_array,
41 : : DLIList <FacetEntity *> &facet_list);
42 : : static CubitStatus extractEdgeTangsFromFile( int num_file_edges,
43 : : int *edge_nodes_file,
44 : : double *edge_node_tangs_file,
45 : : int numTri,
46 : : int *edgeVert,
47 : : double *edgeVertTang );
48 : : static CubitStatus extractTriNormalsFromFile( int num_file_tris,
49 : : int *tri_nodes_file,
50 : : double *tri_node_norms_file,
51 : : int numTri,
52 : : int *edgeVert,
53 : : int *triEdge,
54 : : double *triVertNorms );
55 : : static CubitStatus extractQuadNormalsFromFile( int num_file_quads,
56 : : int *quad_nodes_file,
57 : : double *quad_node_norms_file,
58 : : int numQuad,
59 : : int *edgeVert,
60 : : int *quadEdge,
61 : : double *quadVertNorms );
62 : : static void constructTriVerts( int triEdge[3],
63 : : int *edgeVert,
64 : : int this_tri[3] );
65 : : static void constructQuadVerts( int quadEdge[3],
66 : : int *edgeVert,
67 : : int this_quad[4] );
68 : : static void checkMemoryAllocations( int num_nodes_per_elem,
69 : : int additional_num_elems,
70 : : int *num_elems,
71 : : int **nodes,
72 : : double **node_normals );
73 : :
74 : : static int classify_local_convexity_at_edge(CubitFacetEdge *edge_ptr);
75 : :
76 : : //===========================================================================
77 : : // Function: constructBezier
78 : : // Purpose: Return control points for a quartic Bezier for each face, using
79 : : // a given feature angle tolerance to distinguish C0 and C1 edges.
80 : : // Separate arrays of tris and quad faces should be supplied
81 : : // Date: 10/28/2002
82 : : // Author: sjowen
83 : : //===========================================================================
84 : 0 : void constructBezier( double angle, int numTri, int numQuad, int numEdge,
85 : : int numVert, int* triEdge, int* quadEdge,
86 : : int* edgeVert, double* vert,
87 : : double* edgeCtrlPts, double* triCtrlPts,
88 : : double* quadCtrlPts )
89 : : {
90 : :
91 : : // create arrays of facet entities from the arrays
92 : :
93 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
94 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
95 : 0 : CubitFacet **facet_array = NULL;
96 [ # # ]: 0 : if (numTri > 0)
97 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numTri];
98 : 0 : CubitQuadFacet **qfacet_array = NULL;
99 [ # # ]: 0 : if (numQuad > 0)
100 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numQuad];
101 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
102 : 0 : CubitStatus stat = CUBIT_SUCCESS;
103 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
104 : : edgeVert, vert, point_array, edge_array, facet_array,
105 [ # # ]: 0 : qfacet_array, facet_list);
106 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
107 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
108 : 0 : return;
109 : : }
110 : :
111 : : // create lists of points and edges
112 : :
113 : : int ii;
114 [ # # ][ # # ]: 0 : DLIList<CubitPoint *> point_list;
[ # # ][ # # ]
115 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
116 [ # # ]: 0 : point_list.append(point_array[ii]);
117 : :
118 [ # # ][ # # ]: 0 : DLIList<CubitFacetEdge *> edge_list;
119 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
120 [ # # ]: 0 : edge_list.append(edge_array[ii]);
121 : :
122 : : // prepare the ChollaEngine
123 : :
124 [ # # ][ # # ]: 0 : DLIList<FacetEntity *> fpoint_list, fedge_list;
[ # # ][ # # ]
125 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( point_list, fpoint_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
126 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( edge_list, fedge_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
127 [ # # ][ # # ]: 0 : ChollaEngine *ceng = new ChollaEngine( facet_list, fedge_list, fpoint_list );
128 : :
129 : : CubitBoolean use_feature_angle;
130 [ # # ][ # # ]: 0 : if (angle < 0.00001 || angle > 179.9999)
131 : 0 : use_feature_angle = CUBIT_FALSE;
132 : : else
133 : 0 : use_feature_angle = CUBIT_TRUE;
134 : 0 : int interp_order = 4; // assume we are creating beziers
135 : 0 : CubitBoolean smooth_non_manifold = CUBIT_TRUE;
136 : 0 : CubitBoolean split_surfaces = CUBIT_FALSE;
137 : :
138 : : // create the geometry
139 : :
140 : : //CubitStatus rv =
141 : : ceng->create_geometry(use_feature_angle,angle,
142 : : interp_order,smooth_non_manifold,
143 [ # # ]: 0 : split_surfaces);
144 : :
145 : : // get the control points on edges
146 : :
147 : : int jj;
148 : : int index;
149 : : CubitVector *control_points;
150 : : CubitFacetEdge *edge_ptr;
151 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
152 : : {
153 : 0 : edge_ptr = edge_array[ii];
154 [ # # ]: 0 : control_points = edge_ptr->control_points();
155 [ # # ]: 0 : assert(control_points != NULL);
156 [ # # ]: 0 : for (jj=0; jj<NUM_EDGE_CPTS; jj++)
157 : : {
158 : 0 : index = 3 * (ii * NUM_EDGE_CPTS + jj);
159 [ # # ]: 0 : edgeCtrlPts[index] = control_points[jj].x();
160 [ # # ]: 0 : edgeCtrlPts[index+1] = control_points[jj].y();
161 [ # # ]: 0 : edgeCtrlPts[index+2] = control_points[jj].z();
162 : : }
163 : : }
164 : :
165 : : // get the control points on triangles
166 : :
167 : : CubitFacet *facet_ptr;
168 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
169 : : {
170 : 0 : facet_ptr = facet_array[ii];
171 [ # # ]: 0 : control_points = facet_ptr->control_points();
172 [ # # ]: 0 : assert(control_points != NULL);
173 [ # # ]: 0 : for (jj=0; jj<NUM_TRI_CPTS; jj++)
174 : : {
175 : 0 : index = 3 * (ii * NUM_TRI_CPTS + jj);
176 [ # # ]: 0 : triCtrlPts[index] = control_points[jj].x();
177 [ # # ]: 0 : triCtrlPts[index+1] = control_points[jj].y();
178 [ # # ]: 0 : triCtrlPts[index+2] = control_points[jj].z();
179 : : }
180 : : }
181 : :
182 : : // get the control points on quads
183 : :
184 : : //CubitVector qctrl_points[NUM_QUAD_CPTS];
185 : : CubitQuadFacet *qfacet_ptr;
186 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
187 : : {
188 : 0 : qfacet_ptr = qfacet_array[ii];
189 : 0 : index = 3 * (ii * NUM_QUAD_CPTS);
190 [ # # ]: 0 : qfacet_ptr->get_control_points( &(quadCtrlPts[index]) );
191 : : }
192 : :
193 : : // delete the ChollaEngine
194 : :
195 [ # # ]: 0 : ceng->delete_eval_tools();
196 [ # # ]: 0 : ceng->delete_me();
197 [ # # ][ # # ]: 0 : delete ceng;
198 : :
199 : : // delete temp arrays (facets are deleted with the eval tool)
200 : :
201 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
202 : : {
203 [ # # ]: 0 : qfacet_array[ii]->remove_tri_facets();
204 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
205 : : }
206 [ # # ]: 0 : delete [] point_array;
207 [ # # ]: 0 : delete [] edge_array;
208 [ # # ]: 0 : if (facet_array)
209 [ # # ]: 0 : delete [] facet_array;
210 [ # # ]: 0 : if (qfacet_array)
211 [ # # ]: 0 : delete [] qfacet_array;
212 : :
213 : : }
214 : :
215 : : //===========================================================================
216 : : // Function: constructTriNormals
217 : : // Purpose: Return normals and tangents for a set of triangles using
218 : : // a given feature angle tolerance to distinguish C0 and C1 edges.
219 : : // Date: 10/10/2003
220 : : // Author: sjowen
221 : : //===========================================================================
222 : 0 : void constructTriNormals( double angle, int numTri, int numEdge,
223 : : int numVert, int* triEdge,
224 : : int* edgeVert, double* vert,
225 : : double* triVertNorms, double* edgeVertTang,
226 : : int* edgeTypeFlags)
227 : : {
228 : :
229 : : // create arrays of facet entities from the arrays
230 : :
231 : 0 : int* quadEdge = NULL;
232 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
233 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
234 : 0 : CubitFacet **facet_array = NULL;
235 [ # # ]: 0 : if (numTri > 0)
236 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numTri];
237 : 0 : int numQuad = 0;
238 : 0 : CubitQuadFacet **qfacet_array = NULL;
239 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
240 : 0 : CubitStatus stat = CUBIT_SUCCESS;
241 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
242 : : edgeVert, vert, point_array, edge_array, facet_array,
243 [ # # ]: 0 : qfacet_array, facet_list);
244 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
245 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
246 : 0 : return;
247 : : }
248 : :
249 : : // create lists of points and edges
250 : :
251 : : int ii;
252 [ # # ][ # # ]: 0 : DLIList<CubitPoint *> point_list;
[ # # ][ # # ]
253 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
254 [ # # ]: 0 : point_list.append(point_array[ii]);
255 : :
256 [ # # ][ # # ]: 0 : DLIList<CubitFacetEdge *> edge_list;
257 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++){
258 [ # # ]: 0 : edge_list.append(edge_array[ii]);
259 : : }
260 : :
261 : :
262 : : // prepare the ChollaEngine
263 : :
264 [ # # ][ # # ]: 0 : DLIList<FacetEntity *> fpoint_list, fedge_list;
[ # # ][ # # ]
265 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( point_list, fpoint_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
266 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( edge_list, fedge_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
267 [ # # ][ # # ]: 0 : ChollaEngine *ceng = new ChollaEngine( facet_list, fedge_list, fpoint_list );
268 : :
269 : : CubitBoolean use_feature_angle;
270 [ # # ][ # # ]: 0 : if (angle < 0.00001 || angle > 179.9999)
271 : 0 : use_feature_angle = CUBIT_FALSE;
272 : : else
273 : 0 : use_feature_angle = CUBIT_TRUE;
274 : 0 : int interp_order = 4; // assume we are creating beziers
275 : 0 : CubitBoolean smooth_non_manifold = CUBIT_TRUE;
276 : 0 : CubitBoolean split_surfaces = CUBIT_FALSE;
277 : :
278 : : // create the geometry
279 : :
280 : : //CubitStatus rv =
281 : : ceng->create_geometry(use_feature_angle,angle,
282 : : interp_order,smooth_non_manifold,
283 [ # # ]: 0 : split_surfaces);
284 : :
285 : : // get the control points on edges and infer the tangents
286 : :
287 : 0 : int index = 0;
288 : : CubitVector *control_points;
289 : : CubitFacetEdge *edge_ptr;
290 [ # # ][ # # ]: 0 : CubitVector v0, v1, t0, t1;
[ # # ][ # # ]
291 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
292 : : {
293 : 0 : edge_ptr = edge_array[ii];
294 [ # # ][ # # ]: 0 : v0 = edge_ptr->point(0)->coordinates();
[ # # ]
295 [ # # ][ # # ]: 0 : v1 = edge_ptr->point(1)->coordinates();
[ # # ]
296 [ # # ]: 0 : control_points = edge_ptr->control_points();
297 [ # # ]: 0 : assert(control_points != NULL);
298 : 0 : index = 6 * ii;
299 [ # # ][ # # ]: 0 : t0 = control_points[0] - v0;
300 [ # # ][ # # ]: 0 : t1 = v1 - control_points[2];
301 [ # # ][ # # ]: 0 : if(edge_ptr->is_flipped()){
302 [ # # ]: 0 : t0 *= (-1.0);
303 [ # # ]: 0 : t1 *= (-1.0);
304 : : }
305 [ # # ]: 0 : t0.normalize();
306 [ # # ]: 0 : t1.normalize();
307 [ # # ]: 0 : edgeVertTang[index ] = t0.x();
308 [ # # ]: 0 : edgeVertTang[index+1] = t0.y();
309 [ # # ]: 0 : edgeVertTang[index+2] = t0.z();
310 [ # # ]: 0 : edgeVertTang[index+3] = t1.x();
311 [ # # ]: 0 : edgeVertTang[index+4] = t1.y();
312 [ # # ]: 0 : edgeVertTang[index+5] = t1.z();
313 [ # # ]: 0 : edgeTypeFlags[ii]=classify_local_convexity_at_edge(edge_ptr);
314 : : }
315 : :
316 : : // get the normal on triangles
317 : :
318 [ # # ][ # # ]: 0 : CubitVector n0, n1, n2;
[ # # ]
319 : : CubitFacet *facet_ptr;
320 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
321 : : {
322 : 0 : facet_ptr = facet_array[ii];
323 : 0 : index = 9 * ii;
324 [ # # ][ # # ]: 0 : n0 = facet_ptr->point(0)->normal(facet_ptr);
[ # # ]
325 [ # # ][ # # ]: 0 : n1 = facet_ptr->point(1)->normal(facet_ptr);
[ # # ]
326 [ # # ][ # # ]: 0 : n2 = facet_ptr->point(2)->normal(facet_ptr);
[ # # ]
327 [ # # ]: 0 : triVertNorms[index ] = n0.x();
328 [ # # ]: 0 : triVertNorms[index + 1] = n0.y();
329 [ # # ]: 0 : triVertNorms[index + 2] = n0.z();
330 [ # # ]: 0 : triVertNorms[index + 3] = n1.x();
331 [ # # ]: 0 : triVertNorms[index + 4] = n1.y();
332 [ # # ]: 0 : triVertNorms[index + 5] = n1.z();
333 [ # # ]: 0 : triVertNorms[index + 6] = n2.x();
334 [ # # ]: 0 : triVertNorms[index + 7] = n2.y();
335 [ # # ]: 0 : triVertNorms[index + 8] = n2.z();
336 : : }
337 : :
338 : : // delete the ChollaEngine
339 : :
340 [ # # ]: 0 : ceng->delete_eval_tools();
341 [ # # ]: 0 : ceng->delete_me();
342 [ # # ][ # # ]: 0 : delete ceng;
343 : :
344 : : // delete temp arrays (facets are deleted with the eval tool)
345 : :
346 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
347 : : {
348 [ # # ]: 0 : qfacet_array[ii]->remove_tri_facets();
349 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
350 : : }
351 [ # # ]: 0 : delete [] point_array;
352 [ # # ]: 0 : delete [] edge_array;
353 [ # # ]: 0 : if (facet_array)
354 [ # # ]: 0 : delete [] facet_array;
355 [ # # ]: 0 : if (qfacet_array)
356 [ # # ]: 0 : delete [] qfacet_array;
357 : : }
358 : :
359 : : //===========================================================================
360 : : // Function: constructQuadNormals
361 : : // Purpose: Return normals and tangents for a set of quads using
362 : : // a given feature angle tolerance to distinguish C0 and C1 edges.
363 : : // Date: 10/10/2003
364 : : // Author: sjowen
365 : : //===========================================================================
366 : 0 : void constructQuadNormals( double angle, int numQuad, int numEdge,
367 : : int numVert, int* quadEdge,
368 : : int* edgeVert, double* vert,
369 : : double* quadVertNorms, double* edgeVertTang,
370 : : int* edgeTypeFlags )
371 : : {
372 : :
373 : : // create arrays of facet entities from the arrays
374 : :
375 : 0 : int *triEdge = NULL;
376 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
377 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
378 : 0 : CubitFacet **facet_array = NULL;
379 : 0 : CubitQuadFacet **qfacet_array = NULL;
380 : 0 : int numTri = 0;
381 [ # # ]: 0 : if (numQuad > 0)
382 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numQuad];
383 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
384 : 0 : CubitStatus stat = CUBIT_SUCCESS;
385 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
386 : : edgeVert, vert, point_array, edge_array, facet_array,
387 [ # # ]: 0 : qfacet_array, facet_list);
388 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
389 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
390 : 0 : return;
391 : : }
392 : : // create lists of points and edges
393 : :
394 : : int ii;
395 : :
396 [ # # ][ # # ]: 0 : DLIList<CubitPoint *> point_list;
[ # # ][ # # ]
397 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
398 [ # # ]: 0 : point_list.append(point_array[ii]);
399 : :
400 [ # # ][ # # ]: 0 : DLIList<CubitFacetEdge *> edge_list;
401 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
402 [ # # ]: 0 : edge_list.append(edge_array[ii]);
403 : :
404 : : // prepare the ChollaEngine
405 : :
406 [ # # ][ # # ]: 0 : DLIList<FacetEntity *> fpoint_list, fedge_list;
[ # # ][ # # ]
407 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( point_list, fpoint_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
408 [ # # ][ # # ]: 0 : CAST_LIST_TO_PARENT( edge_list, fedge_list );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
409 [ # # ][ # # ]: 0 : ChollaEngine *ceng = new ChollaEngine( facet_list, fedge_list, fpoint_list );
410 : :
411 : : CubitBoolean use_feature_angle;
412 [ # # ][ # # ]: 0 : if (angle < 0.00001 || angle > 179.9999)
413 : 0 : use_feature_angle = CUBIT_FALSE;
414 : : else
415 : 0 : use_feature_angle = CUBIT_TRUE;
416 : 0 : int interp_order = 4; // assume we are creating beziers
417 : 0 : CubitBoolean smooth_non_manifold = CUBIT_TRUE;
418 : 0 : CubitBoolean split_surfaces = CUBIT_FALSE;
419 : :
420 : : // create the geometry
421 : :
422 : : ceng->create_geometry(use_feature_angle,angle,
423 : : interp_order,smooth_non_manifold,
424 [ # # ]: 0 : split_surfaces);
425 : :
426 : : // get the control points on edges and infer the tangents
427 : :
428 : 0 : int mydebug = 0;
429 : 0 : FILE *fp = NULL;
430 [ # # ]: 0 : if (mydebug)
431 : : {
432 [ # # ]: 0 : fp = fopen("edges.txt", "w");
433 : : }
434 : 0 : int index = 0;
435 : : CubitVector *control_points;
436 : : CubitFacetEdge *edge_ptr;
437 [ # # ][ # # ]: 0 : CubitVector v0, v1, t0, t1;
[ # # ][ # # ]
438 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
439 : : {
440 : 0 : edge_ptr = edge_array[ii];
441 [ # # ][ # # ]: 0 : v0 = edge_ptr->point(0)->coordinates();
[ # # ]
442 [ # # ][ # # ]: 0 : v1 = edge_ptr->point(1)->coordinates();
[ # # ]
443 [ # # ]: 0 : control_points = edge_ptr->control_points();
444 [ # # ]: 0 : assert(control_points != NULL);
445 [ # # ]: 0 : if (mydebug)
446 : : {
447 : : int kk;
448 [ # # ]: 0 : for(kk=0; kk<3; kk++);
449 [ # # ][ # # ]: 0 : fprintf(fp, "(%d) %.6f %.6f %.6f\n", edge_ptr->id(), v0.x(), v0.y(), v0.z());
[ # # ][ # # ]
[ # # ]
450 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", control_points[0].x(), control_points[0].y(), control_points[0].z());
[ # # ][ # # ]
451 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", control_points[1].x(), control_points[1].y(), control_points[1].z());
[ # # ][ # # ]
452 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", control_points[2].x(), control_points[2].y(), control_points[2].z());
[ # # ][ # # ]
453 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", v1.x(), v1.y(), v1.z());
[ # # ][ # # ]
454 : : }
455 : 0 : index = 6 * ii;
456 [ # # ][ # # ]: 0 : t0 = control_points[0] - v0;
457 [ # # ][ # # ]: 0 : t1 = v1 - control_points[2];
458 [ # # ][ # # ]: 0 : if(edge_ptr->is_flipped()){
459 [ # # ]: 0 : t0 *= (-1.0);
460 [ # # ]: 0 : t1 *= (-1.0);
461 : : }
462 [ # # ]: 0 : t0.normalize();
463 [ # # ]: 0 : t1.normalize();
464 [ # # ]: 0 : edgeVertTang[index ] = t0.x();
465 [ # # ]: 0 : edgeVertTang[index+1] = t0.y();
466 [ # # ]: 0 : edgeVertTang[index+2] = t0.z();
467 [ # # ]: 0 : edgeVertTang[index+3] = t1.x();
468 [ # # ]: 0 : edgeVertTang[index+4] = t1.y();
469 [ # # ]: 0 : edgeVertTang[index+5] = t1.z();
470 [ # # ]: 0 : edgeTypeFlags[ii]=classify_local_convexity_at_edge(edge_ptr);
471 : : }
472 [ # # ]: 0 : if (mydebug)
473 : : {
474 [ # # ]: 0 : fclose(fp);
475 : : }
476 : :
477 : : // get the normal on quad
478 : :
479 [ # # ][ # # ]: 0 : CubitVector n0, n1, n2, n3;
[ # # ][ # # ]
480 : : CubitQuadFacet *qfacet_ptr;
481 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
482 : : {
483 : 0 : qfacet_ptr = qfacet_array[ii];
484 : 0 : index = 12 * ii;
485 [ # # ][ # # ]: 0 : n0 = qfacet_ptr->point(0)->normal(qfacet_ptr);
[ # # ]
486 [ # # ][ # # ]: 0 : n1 = qfacet_ptr->point(1)->normal(qfacet_ptr);
[ # # ]
487 [ # # ][ # # ]: 0 : n2 = qfacet_ptr->point(2)->normal(qfacet_ptr);
[ # # ]
488 [ # # ][ # # ]: 0 : n3 = qfacet_ptr->point(3)->normal(qfacet_ptr);
[ # # ]
489 [ # # ]: 0 : quadVertNorms[index ] = n0.x();
490 [ # # ]: 0 : quadVertNorms[index + 1 ] = n0.y();
491 [ # # ]: 0 : quadVertNorms[index + 2 ] = n0.z();
492 [ # # ]: 0 : quadVertNorms[index + 3 ] = n1.x();
493 [ # # ]: 0 : quadVertNorms[index + 4 ] = n1.y();
494 [ # # ]: 0 : quadVertNorms[index + 5 ] = n1.z();
495 [ # # ]: 0 : quadVertNorms[index + 6 ] = n2.x();
496 [ # # ]: 0 : quadVertNorms[index + 7 ] = n2.y();
497 [ # # ]: 0 : quadVertNorms[index + 8 ] = n2.z();
498 [ # # ]: 0 : quadVertNorms[index + 9 ] = n3.x();
499 [ # # ]: 0 : quadVertNorms[index + 10] = n3.y();
500 [ # # ]: 0 : quadVertNorms[index + 11] = n3.z();
501 : : }
502 : :
503 : : // delete the ChollaEngine
504 : :
505 [ # # ]: 0 : ceng->delete_eval_tools();
506 [ # # ]: 0 : ceng->delete_me();
507 [ # # ][ # # ]: 0 : delete ceng;
508 : :
509 : : // delete temp arrays
510 : :
511 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
512 : : {
513 [ # # ]: 0 : qfacet_array[ii]->remove_tri_facets();
514 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
515 : : }
516 [ # # ]: 0 : delete [] point_array;
517 [ # # ]: 0 : delete [] edge_array;
518 [ # # ]: 0 : if (facet_array)
519 [ # # ]: 0 : delete [] facet_array;
520 [ # # ]: 0 : if (qfacet_array)
521 [ # # ]: 0 : delete [] qfacet_array;
522 : : }
523 : :
524 : : //===========================================================================
525 : : // Function: build_facets
526 : : // Purpose: static function that creates facet entities
527 : : // from initial the arrays
528 : : // Date: 10/28/2002
529 : : // Author: sjowen
530 : : //===========================================================================
531 : 0 : static CubitStatus build_facets(int numTri, int numQuad, int numEdge,
532 : : int numVert, int* triEdge, int *quadEdge,
533 : : int* edgeVert, double* vert,
534 : : CubitPoint **point_array,
535 : : CubitFacetEdge **edge_array,
536 : : CubitFacet **facet_array,
537 : : CubitQuadFacet **qfacet_array,
538 : : DLIList <FacetEntity *> &facet_list)
539 : : {
540 : : // create the points
541 : :
542 : : double x, y, z;
543 : : int ii;
544 : :
545 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
546 : : {
547 : 0 : x = vert[ii*3];
548 : 0 : y = vert[ii*3 + 1];
549 : 0 : z = vert[ii*3 + 2];
550 [ # # ][ # # ]: 0 : point_array[ii] = (CubitPoint *) new CubitPointData( x, y, z );
551 [ # # ]: 0 : point_array[ii]->set_id(ii);
552 : : }
553 : :
554 : : // create the edges
555 : :
556 : : int ip, jp;
557 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
558 : : {
559 : 0 : ip = edgeVert[ii*2];
560 : 0 : jp = edgeVert[ii*2 + 1];
561 [ # # ][ # # ]: 0 : assert(ip < numVert && jp < numVert && ip >= 0 && jp >= 0);
[ # # ][ # # ]
562 : 0 : edge_array[ii] = (CubitFacetEdge *) new CubitFacetEdgeData( point_array[ip],
563 [ # # ][ # # ]: 0 : point_array[jp] );
564 [ # # ]: 0 : edge_array[ii]->set_id(ii);
565 : : }
566 : :
567 : : // create the tri facets
568 : :
569 : : int begin_face;
570 : : int jj, iedge;
571 : : CubitFacetEdge *edges[4];
572 : 0 : CubitFacet *facet_ptr = NULL;
573 : :
574 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
575 : : {
576 : 0 : begin_face = 3 * ii;
577 [ # # ]: 0 : for(jj=0; jj<3; jj++)
578 : : {
579 : 0 : iedge = triEdge[begin_face+jj];
580 : 0 : edges[jj] = edge_array[iedge];
581 : : }
582 : : facet_ptr = (CubitFacet *)
583 [ # # ][ # # ]: 0 : new CubitFacetData(edges[0], edges[1], edges[2]);
584 [ # # ]: 0 : facet_list.append(facet_ptr);
585 : 0 : facet_array[ii] = facet_ptr;
586 : : }
587 : :
588 : : // create the quad facets
589 : :
590 : 0 : CubitQuadFacet *qfacet_ptr = NULL;
591 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
592 : : {
593 : 0 : begin_face = 4 * ii;
594 [ # # ]: 0 : for(jj=0; jj<4; jj++)
595 : : {
596 : 0 : iedge = quadEdge[begin_face+jj];
597 : 0 : edges[jj] = edge_array[iedge];
598 : : }
599 : : qfacet_ptr = (CubitQuadFacet *)
600 [ # # ][ # # ]: 0 : new CubitQuadFacetData(edges[0], edges[1], edges[2], edges[3]);
601 [ # # ][ # # ]: 0 : facet_list.append(qfacet_ptr->get_tri_facet(0));
602 [ # # ][ # # ]: 0 : facet_list.append(qfacet_ptr->get_tri_facet(1));
603 : 0 : qfacet_array[ii] = qfacet_ptr;
604 : : }
605 : :
606 : 0 : int mydebug = 0;
607 [ # # ]: 0 : if (mydebug)
608 : : {
609 [ # # ]: 0 : DLIList<CubitFacet*> myfacet_list;
610 [ # # ][ # # ]: 0 : CAST_LIST(facet_list, myfacet_list, CubitFacet);
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
611 [ # # ][ # # ]: 0 : FacetDataUtil::write_facets("C:\\CUBIT\\cholla_apps\\examples\\sierra\\test.facets",myfacet_list);
612 : : }
613 : :
614 : 0 : return CUBIT_SUCCESS;
615 : :
616 : : }
617 : :
618 : : //===========================================================================
619 : : // Function: evalBezierFace
620 : : // Purpose: Evaluate a set of quartic Bezier patches at a specified
621 : : // parametric location given its control points. Evaluates
622 : : // location, normal and/or derivative. Expects list of either
623 : : // quads or tris (not both)
624 : : //
625 : : // Date: 10/28/2002
626 : : // Author: sjowen
627 : : //===========================================================================
628 : 0 : void evalBezierFace( int numFace, int numEdge, int numVert, int numVertPerFace,
629 : : int* faceEdge, int* edgeVert,
630 : : double* vert, double* faceCtrlPts,
631 : : double* edgeCtrlPts,
632 : : int numLocs, double* paramLocation,
633 : : double* location, double* normal, double* deriv )
634 : : {
635 : : // create arrays of facet entities from the arrays
636 : :
637 : 0 : int numTri = 0;
638 : 0 : int numQuad = 0;
639 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
640 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
641 : 0 : CubitFacet **facet_array = NULL;
642 [ # # ]: 0 : int *tri_edge = (numVertPerFace == 3) ? faceEdge : NULL;
643 [ # # ]: 0 : int *quad_edge = (numVertPerFace == 4) ? faceEdge : NULL;
644 [ # # ]: 0 : if (numVertPerFace == 3)
645 : : {
646 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numFace];
647 : 0 : numTri = numFace;
648 : : }
649 : 0 : CubitQuadFacet **qfacet_array = NULL;
650 [ # # ]: 0 : if (numVertPerFace == 4)
651 : : {
652 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numFace];
653 : 0 : numQuad = numFace;
654 : : }
655 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
656 : 0 : CubitStatus stat = CUBIT_SUCCESS;
657 : : stat = build_facets(numTri, numQuad, numEdge, numVert, tri_edge, quad_edge,
658 : : edgeVert, vert, point_array, edge_array, facet_array,
659 [ # # ]: 0 : qfacet_array, facet_list);
660 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
661 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
662 : 0 : return;
663 : : }
664 : :
665 : : // set the control points on edges
666 : :
667 : : int index, ii, jj;
668 : : CubitFacetEdge *edge_ptr;
669 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
670 : : {
671 : 0 : edge_ptr = edge_array[ii];
672 : 0 : index = 3 * (ii * NUM_EDGE_CPTS);
673 [ # # ]: 0 : edge_ptr->set_control_points( &(edgeCtrlPts[index]) );
674 : : }
675 : :
676 : : // set the control points on triangles
677 : :
678 : : CubitFacet *facet_ptr;
679 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
680 : : {
681 : 0 : facet_ptr = facet_array[ii];
682 : 0 : index = 3 * (ii * NUM_TRI_CPTS);
683 [ # # ]: 0 : facet_ptr->set_control_points( &(faceCtrlPts[index]) );
684 : : }
685 : :
686 : : // set the control points on quads
687 : :
688 : : CubitQuadFacet *qfacet_ptr;
689 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
690 : : {
691 : 0 : qfacet_ptr = qfacet_array[ii];
692 : 0 : index = 3 * (ii * NUM_QUAD_CPTS);
693 [ # # ]: 0 : qfacet_ptr->set_control_points( &(faceCtrlPts[index]) );
694 : : }
695 : :
696 : : // Do the evaluations
697 : :
698 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
699 : : {
700 [ # # ]: 0 : for (jj=0; jj<numLocs; jj++)
701 : : {
702 : 0 : facet_ptr = facet_array[ii];
703 : 0 : CubitVector areacoord(paramLocation[jj*3],
704 : 0 : paramLocation[jj*3+1],
705 [ # # ]: 0 : paramLocation[jj*3+2]);
706 [ # # ]: 0 : CubitVector eval_point;
707 : 0 : CubitVector *eval_point_ptr = &eval_point;
708 [ # # ]: 0 : CubitVector eval_normal;
709 : 0 : CubitVector *eval_normal_ptr = NULL;
710 [ # # ]: 0 : if (normal != NULL)
711 : : {
712 : 0 : eval_normal_ptr = &eval_normal;
713 : : }
714 [ # # ]: 0 : facet_ptr->evaluate( areacoord, eval_point_ptr, eval_normal_ptr );
715 : 0 : index = 3 * ii * jj;
716 [ # # ]: 0 : location[index] = eval_point.x();
717 [ # # ]: 0 : location[index+1] = eval_point.y();
718 [ # # ]: 0 : location[index+2] = eval_point.z();
719 [ # # ]: 0 : if (normal != NULL)
720 : : {
721 [ # # ]: 0 : normal[index] = eval_normal.x();
722 [ # # ]: 0 : normal[index+1] = eval_normal.y();
723 [ # # ]: 0 : normal[index+2] = eval_normal.z();
724 : : }
725 : : }
726 : : }
727 : :
728 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
729 : : {
730 [ # # ]: 0 : for (jj=0; jj<numLocs; jj++)
731 : : {
732 : 0 : qfacet_ptr = qfacet_array[ii];
733 [ # # ]: 0 : CubitVector eval_point;
734 : 0 : CubitVector *eval_point_ptr = &eval_point;
735 [ # # ]: 0 : CubitVector eval_normal;
736 : 0 : CubitVector *eval_normal_ptr = NULL;
737 [ # # ]: 0 : if (normal != NULL)
738 : 0 : eval_normal_ptr = &eval_normal;
739 [ # # ][ # # ]: 0 : CubitVector eval_du, eval_dv;
740 : 0 : CubitVector *eval_du_ptr = NULL;
741 : 0 : CubitVector *eval_dv_ptr = NULL;
742 [ # # ]: 0 : if (deriv != NULL)
743 : : {
744 : 0 : eval_du_ptr = &eval_du;
745 : 0 : eval_dv_ptr = &eval_dv;
746 : : }
747 : :
748 : 0 : qfacet_ptr->evaluate( paramLocation[2*jj], paramLocation[2*jj+1],
749 : : eval_point_ptr, eval_normal_ptr,
750 [ # # ]: 0 : eval_du_ptr, eval_dv_ptr);
751 : 0 : index = 3 * ii * jj;
752 [ # # ]: 0 : location[index] = eval_point.x();
753 [ # # ]: 0 : location[index+1] = eval_point.y();
754 [ # # ]: 0 : location[index+2] = eval_point.z();
755 [ # # ]: 0 : if (normal != NULL)
756 : : {
757 [ # # ]: 0 : normal[index] = eval_normal.x();
758 [ # # ]: 0 : normal[index+1] = eval_normal.y();
759 [ # # ]: 0 : normal[index+2] = eval_normal.z();
760 : : }
761 [ # # ]: 0 : if (deriv != NULL)
762 : : {
763 : 0 : index = 6 * ii * jj;
764 [ # # ]: 0 : deriv[index] = eval_du.x();
765 [ # # ]: 0 : deriv[index+1] = eval_du.y();
766 [ # # ]: 0 : deriv[index+2] = eval_du.z();
767 [ # # ]: 0 : deriv[index+3] = eval_dv.x();
768 [ # # ]: 0 : deriv[index+4] = eval_dv.y();
769 [ # # ]: 0 : deriv[index+5] = eval_dv.z();
770 : : }
771 : : }
772 : : }
773 : :
774 : : // delete the temp arrays
775 : :
776 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
777 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
778 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
779 [ # # ][ # # ]: 0 : delete facet_array[ii];
780 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
781 [ # # ][ # # ]: 0 : delete edge_array[ii];
782 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
783 [ # # ][ # # ]: 0 : delete point_array[ii];
784 [ # # ]: 0 : delete [] point_array;
785 [ # # ]: 0 : delete [] edge_array;
786 [ # # ]: 0 : if (facet_array)
787 [ # # ]: 0 : delete [] facet_array;
788 [ # # ]: 0 : if (qfacet_array)
789 [ # # ][ # # ]: 0 : delete [] qfacet_array;
[ # # ]
790 : : }
791 : :
792 : : //===========================================================================
793 : : // Function: evalBezierFaceFromNorms
794 : : // Purpose: This is the same as the evalBezierFace except it differes by
795 : : // the input arguments. This function takes a list of normals
796 : : // and tangents at the face vertices and computes bezier control
797 : : // points internally. Normals and tangents should have been
798 : : // computed previously in constructTriNormals or constructQuadNormals.
799 : : // This function is not as computationally efficient as evalBezierFace
800 : : // since it requires Bezier control points to be computed
801 : : // as part of the call - however, it is more memory efficient,
802 : : // requiring fewer variables to be stored with the calling application.
803 : : // The following argument list describes only those arguments that
804 : : // differ from evalBezierFace above.
805 : : //
806 : : // Date: 10/15/2003
807 : : // Author: sjowen
808 : : //===========================================================================
809 : 0 : void evalBezierFaceFromNorms( int numFace, int numEdge, int numVert,
810 : : int numVertPerFace, int* faceEdge, int* edgeVert,
811 : : double* vert, double* vertNorms, double* edgeVertTang,
812 : : int numLocs, double* paramLocation,
813 : : double* location, double* normal, double* deriv )
814 : : {
815 : : // create arrays of facet entities from the arrays
816 : :
817 : 0 : int numTri = 0;
818 : 0 : int numQuad = 0;
819 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
820 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
821 : 0 : CubitFacet **facet_array = NULL;
822 [ # # ]: 0 : int *triEdge = (numVertPerFace == 3) ? faceEdge : NULL;
823 [ # # ]: 0 : int *quadEdge = (numVertPerFace == 4) ? faceEdge : NULL;
824 [ # # ]: 0 : if (numVertPerFace == 3)
825 : : {
826 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numFace];
827 : 0 : numTri = numFace;
828 : : }
829 : 0 : CubitQuadFacet **qfacet_array = NULL;
830 [ # # ]: 0 : if (numVertPerFace == 4)
831 : : {
832 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numFace];
833 : 0 : numQuad = numFace;
834 : : }
835 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
836 : 0 : CubitStatus stat = CUBIT_SUCCESS;
837 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
838 : : edgeVert, vert, point_array, edge_array, facet_array,
839 [ # # ]: 0 : qfacet_array, facet_list);
840 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
841 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
842 : 0 : return;
843 : : }
844 : :
845 : : // set the normals on triangle vertices
846 : :
847 : : int index, ii, jj;
848 : : CubitPoint *pt;
849 : : CubitFacet *facet_ptr;
850 [ # # ]: 0 : CubitVector pt_normal;
851 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
852 : : {
853 : 0 : facet_ptr = facet_array[ii];
854 [ # # ]: 0 : for(jj=0; jj<3; jj++)
855 : : {
856 : 0 : index = ii * 9 + (3 * jj);
857 [ # # ]: 0 : pt_normal.x( vertNorms[ index ] );
858 [ # # ]: 0 : pt_normal.y( vertNorms[ index + 1 ] );
859 [ # # ]: 0 : pt_normal.z( vertNorms[ index + 2 ] );
860 [ # # ]: 0 : pt = facet_ptr->point(jj);
861 [ # # ]: 0 : TDFacetBoundaryPoint::add_facet_boundary_point( pt, facet_ptr, pt_normal );
862 : : }
863 : : }
864 : :
865 : : // set the normals on quad vertices
866 : :
867 : : CubitQuadFacet *qfacet_ptr;
868 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
869 : : {
870 : 0 : qfacet_ptr = qfacet_array[ii];
871 [ # # ]: 0 : for(jj=0; jj<4; jj++)
872 : : {
873 : 0 : index = ii * 12 + (3 * jj);
874 [ # # ]: 0 : pt_normal.x( vertNorms[ index ] );
875 [ # # ]: 0 : pt_normal.y( vertNorms[ index + 1 ] );
876 [ # # ]: 0 : pt_normal.z( vertNorms[ index + 2 ] );
877 [ # # ]: 0 : pt = qfacet_ptr->point(jj);
878 [ # # ]: 0 : TDFacetBoundaryPoint::add_facet_boundary_point( pt, qfacet_ptr, pt_normal );
879 : : }
880 : : }
881 : :
882 : : // set the control points on edges
883 : :
884 : 0 : int mydebug = 0;
885 : 0 : FILE *fp = NULL;
886 [ # # ]: 0 : if (mydebug)
887 [ # # ]: 0 : fp = fopen("edges.txt", "w");
888 : : CubitFacetEdge *edge_ptr;
889 [ # # ][ # # ]: 0 : CubitVector N0, N1, T0, T1, P0, P1;
[ # # ][ # # ]
[ # # ][ # # ]
890 [ # # ][ # # ]: 0 : CubitVector Pi[3];
891 : : CubitPoint *pt0, *pt1;
892 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
893 : : {
894 : 0 : index = 6 * ii;
895 : 0 : edge_ptr = edge_array[ii];
896 [ # # ]: 0 : pt0 = edge_ptr->point(0);
897 [ # # ]: 0 : pt1 = edge_ptr->point(1);
898 [ # # ][ # # ]: 0 : P0 = pt0->coordinates();
899 [ # # ][ # # ]: 0 : P1 = pt1->coordinates();
900 [ # # ][ # # ]: 0 : N0 = pt0->normal( edge_ptr );
901 [ # # ][ # # ]: 0 : N1 = pt1->normal( edge_ptr );
902 [ # # ]: 0 : T0.x( edgeVertTang[ index ] );
903 [ # # ]: 0 : T0.y( edgeVertTang[ index + 1] );
904 [ # # ]: 0 : T0.z( edgeVertTang[ index + 2] );
905 [ # # ]: 0 : T1.x( edgeVertTang[ index + 3] );
906 [ # # ]: 0 : T1.y( edgeVertTang[ index + 4] );
907 [ # # ]: 0 : T1.z( edgeVertTang[ index + 5] );
908 [ # # ]: 0 : FacetEvalTool::init_edge_control_points(P0, P1, N0, N1, T0, T1, Pi);
909 [ # # ]: 0 : if (mydebug)
910 : : {
911 : : int kk;
912 [ # # ]: 0 : for(kk=0; kk<3; kk++);
913 [ # # ][ # # ]: 0 : fprintf(fp, "(%d) %.6f %.6f %.6f\n", edge_ptr->id(), P0.x(), P0.y(), P0.z());
[ # # ][ # # ]
[ # # ]
914 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[0].x(), Pi[0].y(), Pi[0].z());
[ # # ][ # # ]
915 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[1].x(), Pi[1].y(), Pi[1].z());
[ # # ][ # # ]
916 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[2].x(), Pi[2].y(), Pi[2].z());
[ # # ][ # # ]
917 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", P1.x(), P1.y(), P1.z());
[ # # ][ # # ]
918 : : }
919 [ # # ]: 0 : edge_ptr->control_points( Pi, 4 );
920 : : }
921 [ # # ]: 0 : if (mydebug)
922 [ # # ]: 0 : fclose(fp);
923 : :
924 : : // set the control points on triangles
925 : :
926 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
927 : : {
928 : 0 : facet_ptr = facet_array[ii];
929 [ # # ]: 0 : facet_ptr->init_patch();
930 : : }
931 : :
932 : : // set the control points on quads
933 : :
934 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
935 : : {
936 : 0 : qfacet_ptr = qfacet_array[ii];
937 [ # # ]: 0 : qfacet_ptr->init_patch();
938 : : }
939 : :
940 : : // Do the evaluations
941 : :
942 : 0 : index = 0;
943 : 0 : int nindex = 0;
944 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
945 : : {
946 : 0 : facet_ptr = facet_array[ii];
947 [ # # ]: 0 : for (jj=0; jj<numLocs; jj++)
948 : : {
949 : 0 : CubitVector areacoord(paramLocation[jj*3],
950 : 0 : paramLocation[jj*3+1],
951 [ # # ]: 0 : paramLocation[jj*3+2]);
952 [ # # ]: 0 : CubitVector eval_point;
953 : 0 : CubitVector *eval_point_ptr = &eval_point;
954 [ # # ]: 0 : CubitVector eval_normal;
955 : 0 : CubitVector *eval_normal_ptr = NULL;
956 [ # # ]: 0 : if (normal != NULL)
957 : : {
958 : 0 : eval_normal_ptr = &eval_normal;
959 : : }
960 [ # # ]: 0 : facet_ptr->evaluate( areacoord, eval_point_ptr, eval_normal_ptr );
961 [ # # ]: 0 : location[index++] = eval_point.x();
962 [ # # ]: 0 : location[index++] = eval_point.y();
963 [ # # ]: 0 : location[index++] = eval_point.z();
964 [ # # ]: 0 : if (normal != NULL)
965 : : {
966 [ # # ]: 0 : normal[nindex++] = eval_normal.x();
967 [ # # ]: 0 : normal[nindex++] = eval_normal.y();
968 [ # # ]: 0 : normal[nindex++] = eval_normal.z();
969 : : }
970 : : }
971 : : }
972 : :
973 : 0 : int dindex = 0;
974 : 0 : nindex = index = 0;
975 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
976 : : {
977 [ # # ]: 0 : for (jj=0; jj<numLocs; jj++)
978 : : {
979 : 0 : qfacet_ptr = qfacet_array[ii];
980 [ # # ]: 0 : CubitVector eval_point;
981 : 0 : CubitVector *eval_point_ptr = &eval_point;
982 [ # # ]: 0 : CubitVector eval_normal;
983 : 0 : CubitVector *eval_normal_ptr = NULL;
984 [ # # ]: 0 : if (normal != NULL)
985 : 0 : eval_normal_ptr = &eval_normal;
986 [ # # ][ # # ]: 0 : CubitVector eval_du, eval_dv;
987 : 0 : CubitVector *eval_du_ptr = NULL;
988 : 0 : CubitVector *eval_dv_ptr = NULL;
989 [ # # ]: 0 : if (deriv != NULL)
990 : : {
991 : 0 : eval_du_ptr = &eval_du;
992 : 0 : eval_dv_ptr = &eval_dv;
993 : : }
994 : :
995 : 0 : qfacet_ptr->evaluate( paramLocation[2*jj], paramLocation[2*jj+1],
996 : : eval_point_ptr, eval_normal_ptr,
997 [ # # ]: 0 : eval_du_ptr, eval_dv_ptr);
998 [ # # ]: 0 : location[index++] = eval_point.x();
999 [ # # ]: 0 : location[index++] = eval_point.y();
1000 [ # # ]: 0 : location[index++] = eval_point.z();
1001 [ # # ]: 0 : if (normal != NULL)
1002 : : {
1003 [ # # ]: 0 : normal[nindex++] = eval_normal.x();
1004 [ # # ]: 0 : normal[nindex++] = eval_normal.y();
1005 [ # # ]: 0 : normal[nindex++] = eval_normal.z();
1006 : : }
1007 [ # # ]: 0 : if (deriv != NULL)
1008 : : {
1009 : 0 : index = 6 * ii * jj;
1010 [ # # ]: 0 : deriv[dindex++] = eval_du.x();
1011 [ # # ]: 0 : deriv[dindex++] = eval_du.y();
1012 [ # # ]: 0 : deriv[dindex++] = eval_du.z();
1013 [ # # ]: 0 : deriv[dindex++] = eval_dv.x();
1014 [ # # ]: 0 : deriv[dindex++] = eval_dv.y();
1015 [ # # ]: 0 : deriv[dindex++] = eval_dv.z();
1016 : : }
1017 : : }
1018 : : }
1019 : :
1020 : : // delete the temp arrays
1021 : :
1022 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1023 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
1024 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1025 [ # # ][ # # ]: 0 : delete facet_array[ii];
1026 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1027 [ # # ][ # # ]: 0 : delete edge_array[ii];
1028 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
1029 [ # # ][ # # ]: 0 : delete point_array[ii];
1030 [ # # ]: 0 : delete [] point_array;
1031 [ # # ]: 0 : delete [] edge_array;
1032 [ # # ]: 0 : if (facet_array)
1033 [ # # ]: 0 : delete [] facet_array;
1034 [ # # ]: 0 : if (qfacet_array)
1035 [ # # ][ # # ]: 0 : delete [] qfacet_array;
[ # # ]
1036 : : }
1037 : :
1038 : : //===========================================================================
1039 : : // Function: projToBezierFace
1040 : : // Purpose: Project a set of x-y-z locations to Bezier patches. Finds the
1041 : : // closest point on one of the patches. Evaluates location,
1042 : : // normal and/or derivative. Expects list of either quads or
1043 : : // tris (not both)
1044 : : //
1045 : : // Date: 12/08/2003
1046 : : // Author: sjowen
1047 : : //===========================================================================
1048 : 0 : void projToBezierFace( int numFace, int numEdge, int numVert, int numVertPerFace,
1049 : : int* faceEdge, int* edgeVert,
1050 : : double* vert, double* faceCtrlPts,
1051 : : double* edgeCtrlPts,
1052 : : int numLocs, double* xyz,
1053 : : int specify_tol, double converge_tol,
1054 : : double* xyzOnFace, double* normal, double* deriv )
1055 : : {
1056 : : // create arrays of facet entities from the arrays
1057 : :
1058 : 0 : int numTri = 0;
1059 : 0 : int numQuad = 0;
1060 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
1061 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
1062 : 0 : CubitFacet **facet_array = NULL;
1063 [ # # ]: 0 : int *tri_edge = (numVertPerFace == 3) ? faceEdge : NULL;
1064 [ # # ]: 0 : int *quad_edge = (numVertPerFace == 4) ? faceEdge : NULL;
1065 [ # # ]: 0 : if (numVertPerFace == 3)
1066 : : {
1067 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numFace];
1068 : 0 : numTri = numFace;
1069 : : }
1070 : 0 : CubitQuadFacet **qfacet_array = NULL;
1071 [ # # ]: 0 : if (numVertPerFace == 4)
1072 : : {
1073 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numFace];
1074 : 0 : numQuad = numFace;
1075 : : }
1076 [ # # ]: 0 : DLIList<FacetEntity *> facetentity_list;
1077 : 0 : CubitStatus stat = CUBIT_SUCCESS;
1078 : : stat = build_facets(numTri, numQuad, numEdge, numVert, quad_edge, tri_edge,
1079 : : edgeVert, vert, point_array, edge_array, facet_array,
1080 [ # # ]: 0 : qfacet_array, facetentity_list);
1081 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
1082 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
1083 : 0 : return;
1084 : : }
1085 [ # # ][ # # ]: 0 : DLIList<CubitFacet *> facet_list;
[ # # ][ # # ]
1086 [ # # ][ # # ]: 0 : CAST_LIST(facetentity_list, facet_list, CubitFacet);
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
1087 : :
1088 : : // set the control points on edges
1089 : :
1090 : : double edgelen;
1091 : 0 : double compare_tol = 0.0;
1092 : : int index, ii;
1093 : : CubitFacetEdge *edge_ptr;
1094 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1095 : : {
1096 : 0 : edge_ptr = edge_array[ii];
1097 : 0 : index = 3 * (ii * NUM_EDGE_CPTS);
1098 [ # # ]: 0 : edge_ptr->set_control_points( &(edgeCtrlPts[index]) );
1099 [ # # ]: 0 : if ( !specify_tol )
1100 : : {
1101 [ # # ]: 0 : edgelen = edge_ptr->length();
1102 : 0 : compare_tol += edgelen;
1103 : : }
1104 : : }
1105 [ # # ]: 0 : if ( specify_tol )
1106 : : {
1107 : 0 : compare_tol = converge_tol;
1108 : : }
1109 : : else
1110 : : {
1111 : 0 : compare_tol /= numEdge;
1112 : 0 : compare_tol *= 1.0e-3;
1113 : : }
1114 : :
1115 : : // set the control points on triangles
1116 : :
1117 : : CubitFacet *facet_ptr;
1118 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1119 : : {
1120 : 0 : facet_ptr = facet_array[ii];
1121 : 0 : index = 3 * (ii * NUM_TRI_CPTS);
1122 [ # # ]: 0 : facet_ptr->set_control_points( &(faceCtrlPts[index]) );
1123 : : }
1124 : :
1125 : : // set the control points on quads
1126 : :
1127 : : CubitQuadFacet *qfacet_ptr;
1128 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1129 : : {
1130 : 0 : qfacet_ptr = qfacet_array[ii];
1131 : 0 : index = 3 * (ii * NUM_QUAD_CPTS);
1132 [ # # ]: 0 : qfacet_ptr->set_control_points( &(faceCtrlPts[index]) );
1133 : : }
1134 : :
1135 : : // Do the projections
1136 : :
1137 : 0 : CubitFacet *last_facet = NULL;
1138 : 0 : int interp_order = 4;
1139 : 0 : CubitBoolean trim = CUBIT_TRUE;
1140 : : CubitBoolean outside;
1141 [ # # ]: 0 : for (ii=0; ii<numLocs; ii++)
1142 : : {
1143 : 0 : CubitVector cur_location(xyz[ii*3],
1144 : 0 : xyz[ii*3+1],
1145 [ # # ]: 0 : xyz[ii*3+2]);
1146 [ # # ]: 0 : CubitVector eval_point;
1147 : 0 : CubitVector *eval_point_ptr = &eval_point;
1148 [ # # ]: 0 : CubitVector eval_normal;
1149 : 0 : CubitVector *eval_normal_ptr = NULL;
1150 [ # # ]: 0 : if (normal != NULL)
1151 : : {
1152 : 0 : eval_normal_ptr = &eval_normal;
1153 : : }
1154 : : FacetEvalTool::project_to_facets(facet_list, last_facet,
1155 : : interp_order, compare_tol,
1156 : : cur_location, trim,
1157 : : &outside, eval_point_ptr,
1158 [ # # ]: 0 : eval_normal_ptr);
1159 : 0 : index = 3 * ii;
1160 [ # # ]: 0 : xyzOnFace[index] = eval_point.x();
1161 [ # # ]: 0 : xyzOnFace[index+1] = eval_point.y();
1162 [ # # ]: 0 : xyzOnFace[index+2] = eval_point.z();
1163 [ # # ]: 0 : if (normal != NULL)
1164 : : {
1165 [ # # ]: 0 : normal[index] = eval_normal.x();
1166 [ # # ]: 0 : normal[index+1] = eval_normal.y();
1167 [ # # ]: 0 : normal[index+2] = eval_normal.z();
1168 : : }
1169 : : }
1170 : :
1171 : : // no derivatives yet!!
1172 : :
1173 : :
1174 : : // delete the temp arrays
1175 : :
1176 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1177 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
1178 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1179 [ # # ][ # # ]: 0 : delete facet_array[ii];
1180 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1181 [ # # ][ # # ]: 0 : delete edge_array[ii];
1182 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
1183 [ # # ][ # # ]: 0 : delete point_array[ii];
1184 [ # # ]: 0 : delete [] point_array;
1185 [ # # ]: 0 : delete [] edge_array;
1186 [ # # ]: 0 : if (facet_array)
1187 [ # # ]: 0 : delete [] facet_array;
1188 [ # # ]: 0 : if (qfacet_array)
1189 [ # # ]: 0 : delete [] qfacet_array;
1190 : : }
1191 : :
1192 : : //===========================================================================
1193 : : // Function: projToBezierFaceFromNorms
1194 : : // Purpose: This is the same as the projToBezierFace except it differes by
1195 : : // the input arguments. This function takes a list of normals and tangents
1196 : : // at the face vertices and computes bezier control points internally.
1197 : : // Normals and tangents should have been computed previously in
1198 : : // constructTriNormals or constructQuadNormals. This function is not as
1199 : : // computationally efficient as evalBezierFace since it requires Bezier
1200 : : // control points to be computed as part of the call - however, it is more
1201 : : // memory efficient, requiring fewer variables to be stored with the calling
1202 : : // application. The following argument list describes only those arguments
1203 : : // that differ from projToBezierFace above.
1204 : : //
1205 : : // Date: 12/08/2003
1206 : : // Author: sjowen
1207 : : //===========================================================================
1208 : 0 : void projToBezierFaceFromNorms( int numFace, int numEdge, int numVert,
1209 : : int numVertPerFace, int* faceEdge, int* edgeVert,
1210 : : double* vert, double* vertNorms, double* edgeVertTang,
1211 : : int numLocs, double* xyz,
1212 : : int specify_tol, double converge_tol,
1213 : : double* xyzOnFace, double* normal, double* deriv )
1214 : : {
1215 : : // create arrays of facet entities from the arrays
1216 : :
1217 : 0 : int numTri = 0;
1218 : 0 : int numQuad = 0;
1219 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
1220 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
1221 : 0 : CubitFacet **facet_array = NULL;
1222 [ # # ]: 0 : int *triEdge = (numVertPerFace == 3) ? faceEdge : NULL;
1223 [ # # ]: 0 : int *quadEdge = (numVertPerFace == 4) ? faceEdge : NULL;
1224 [ # # ]: 0 : if (numVertPerFace == 3)
1225 : : {
1226 [ # # ][ # # ]: 0 : facet_array = new CubitFacet * [numFace];
1227 : 0 : numTri = numFace;
1228 : : }
1229 : 0 : CubitQuadFacet **qfacet_array = NULL;
1230 [ # # ]: 0 : if (numVertPerFace == 4)
1231 : : {
1232 [ # # ][ # # ]: 0 : qfacet_array = new CubitQuadFacet * [numFace];
1233 : 0 : numQuad = numFace;
1234 : : }
1235 [ # # ]: 0 : DLIList<FacetEntity *> facetentity_list;
1236 : 0 : CubitStatus stat = CUBIT_SUCCESS;
1237 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
1238 : : edgeVert, vert, point_array, edge_array, facet_array,
1239 [ # # ]: 0 : qfacet_array, facetentity_list);
1240 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
1241 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
1242 : 0 : return;
1243 : : }
1244 [ # # ][ # # ]: 0 : DLIList<CubitFacet *> facet_list;
[ # # ][ # # ]
1245 [ # # ][ # # ]: 0 : CAST_LIST(facetentity_list, facet_list, CubitFacet);
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
1246 : :
1247 : : // set the normals on triangle vertices
1248 : :
1249 : : int index, ii, jj;
1250 : : CubitPoint *pt;
1251 : : CubitFacet *facet_ptr;
1252 [ # # ]: 0 : CubitVector pt_normal;
1253 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1254 : : {
1255 : 0 : facet_ptr = facet_array[ii];
1256 [ # # ]: 0 : for(jj=0; jj<3; jj++)
1257 : : {
1258 : 0 : index = ii * 9 + (3 * jj);
1259 [ # # ]: 0 : pt_normal.x( vertNorms[ index ] );
1260 [ # # ]: 0 : pt_normal.y( vertNorms[ index + 1 ] );
1261 [ # # ]: 0 : pt_normal.z( vertNorms[ index + 2 ] );
1262 [ # # ]: 0 : pt = facet_ptr->point(jj);
1263 [ # # ]: 0 : TDFacetBoundaryPoint::add_facet_boundary_point( pt, facet_ptr, pt_normal );
1264 : : }
1265 : : }
1266 : :
1267 : : // set the normals on quad vertices
1268 : :
1269 : : CubitQuadFacet *qfacet_ptr;
1270 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1271 : : {
1272 : 0 : qfacet_ptr = qfacet_array[ii];
1273 [ # # ]: 0 : for(jj=0; jj<4; jj++)
1274 : : {
1275 : 0 : index = ii * 12 + (3 * jj);
1276 [ # # ]: 0 : pt_normal.x( vertNorms[ index ] );
1277 [ # # ]: 0 : pt_normal.y( vertNorms[ index + 1 ] );
1278 [ # # ]: 0 : pt_normal.z( vertNorms[ index + 2 ] );
1279 [ # # ]: 0 : pt = qfacet_ptr->point(jj);
1280 [ # # ]: 0 : TDFacetBoundaryPoint::add_facet_boundary_point( pt, qfacet_ptr, pt_normal );
1281 : : }
1282 : : }
1283 : :
1284 : : // set the control points on edges
1285 : :
1286 : : double edgelen;
1287 : 0 : double compare_tol = 0.0;
1288 : : CubitFacetEdge *edge_ptr;
1289 [ # # ][ # # ]: 0 : CubitVector N0, N1, T0, T1, P0, P1;
[ # # ][ # # ]
[ # # ][ # # ]
1290 [ # # ][ # # ]: 0 : CubitVector Pi[3];
1291 : : CubitPoint *pt0, *pt1;
1292 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1293 : : {
1294 : 0 : index = 6 * ii;
1295 : 0 : edge_ptr = edge_array[ii];
1296 [ # # ]: 0 : pt0 = edge_ptr->point(0);
1297 [ # # ]: 0 : pt1 = edge_ptr->point(1);
1298 [ # # ][ # # ]: 0 : P0 = pt0->coordinates();
1299 [ # # ][ # # ]: 0 : P1 = pt1->coordinates();
1300 [ # # ][ # # ]: 0 : N0 = pt0->normal( edge_ptr );
1301 [ # # ][ # # ]: 0 : N1 = pt1->normal( edge_ptr );
1302 [ # # ]: 0 : T0.x( edgeVertTang[ index ] );
1303 [ # # ]: 0 : T0.y( edgeVertTang[ index + 1] );
1304 [ # # ]: 0 : T0.z( edgeVertTang[ index + 2] );
1305 [ # # ]: 0 : T1.x( edgeVertTang[ index + 3] );
1306 [ # # ]: 0 : T1.y( edgeVertTang[ index + 4] );
1307 [ # # ]: 0 : T1.z( edgeVertTang[ index + 5] );
1308 [ # # ]: 0 : FacetEvalTool::init_edge_control_points(P0, P1, N0, N1, T0, T1, Pi);
1309 [ # # ]: 0 : edge_ptr->control_points( Pi, 4 );
1310 [ # # ]: 0 : if ( !specify_tol )
1311 : : {
1312 [ # # ]: 0 : edgelen = edge_ptr->length();
1313 : 0 : compare_tol += edgelen;
1314 : : }
1315 : : }
1316 [ # # ]: 0 : if ( specify_tol )
1317 : : {
1318 : 0 : compare_tol = converge_tol;
1319 : : }
1320 : : else
1321 : : {
1322 : 0 : compare_tol /= numEdge;
1323 : 0 : compare_tol *= 1.0e-3;
1324 : : }
1325 : :
1326 : : // set the control points on triangles
1327 : :
1328 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1329 : : {
1330 : 0 : facet_ptr = facet_array[ii];
1331 [ # # ]: 0 : facet_ptr->init_patch();
1332 : : }
1333 : :
1334 : : // set the control points on quads
1335 : :
1336 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1337 : : {
1338 : 0 : qfacet_ptr = qfacet_array[ii];
1339 [ # # ]: 0 : qfacet_ptr->init_patch();
1340 : : }
1341 : :
1342 : : // Do the projections
1343 : :
1344 : 0 : CubitFacet *last_facet = NULL;
1345 : 0 : int interp_order = 4;
1346 : 0 : CubitBoolean trim = CUBIT_TRUE;
1347 : : CubitBoolean outside;
1348 : :
1349 [ # # ]: 0 : for (ii=0; ii<numLocs; ii++)
1350 : : {
1351 : 0 : CubitVector cur_location(xyz[ii*3],
1352 : 0 : xyz[ii*3+1],
1353 [ # # ]: 0 : xyz[ii*3+2]);
1354 [ # # ]: 0 : CubitVector eval_point;
1355 : 0 : CubitVector *eval_point_ptr = &eval_point;
1356 [ # # ]: 0 : CubitVector eval_normal;
1357 : 0 : CubitVector *eval_normal_ptr = NULL;
1358 [ # # ]: 0 : if (normal != NULL)
1359 : : {
1360 : 0 : eval_normal_ptr = &eval_normal;
1361 : : }
1362 : : FacetEvalTool::project_to_facets(facet_list, last_facet,
1363 : : interp_order, compare_tol,
1364 : : cur_location, trim,
1365 : : &outside, eval_point_ptr,
1366 [ # # ]: 0 : eval_normal_ptr);
1367 : : //double dist = cur_location.distance_between(eval_point);
1368 : 0 : index = 3 * ii;
1369 [ # # ]: 0 : xyzOnFace[index] = eval_point.x();
1370 [ # # ]: 0 : xyzOnFace[index+1] = eval_point.y();
1371 [ # # ]: 0 : xyzOnFace[index+2] = eval_point.z();
1372 [ # # ]: 0 : if (normal != NULL)
1373 : : {
1374 [ # # ]: 0 : normal[index] = eval_normal.x();
1375 [ # # ]: 0 : normal[index+1] = eval_normal.y();
1376 [ # # ]: 0 : normal[index+2] = eval_normal.z();
1377 : : }
1378 : : }
1379 : :
1380 : : // no derivatives yet!!
1381 : :
1382 : : // delete the temp arrays
1383 : :
1384 [ # # ]: 0 : for (ii=0; ii<numQuad; ii++)
1385 [ # # ][ # # ]: 0 : delete qfacet_array[ii];
1386 [ # # ]: 0 : for (ii=0; ii<numTri; ii++)
1387 [ # # ][ # # ]: 0 : delete facet_array[ii];
1388 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1389 [ # # ][ # # ]: 0 : delete edge_array[ii];
1390 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
1391 [ # # ][ # # ]: 0 : delete point_array[ii];
1392 [ # # ]: 0 : delete [] point_array;
1393 [ # # ]: 0 : delete [] edge_array;
1394 [ # # ]: 0 : if (facet_array)
1395 [ # # ]: 0 : delete [] facet_array;
1396 [ # # ]: 0 : if (qfacet_array)
1397 [ # # ]: 0 : delete [] qfacet_array;
1398 : : }
1399 : :
1400 : : //---------------------------------------------------------------------------
1401 : : // Functions for debugging
1402 : : //---------------------------------------------------------------------------
1403 : :
1404 : : //===========================================================================
1405 : : // Function: dumpMesh
1406 : : // Purpose: dump the face mesh to a file that can be read by CUBIT
1407 : : // Use the same definition of parameters as used with
1408 : : // resolveFaceVectors
1409 : : // Date: 10/28/2002
1410 : : // Author: sjowen
1411 : : //===========================================================================
1412 : 0 : void dumpMesh(const char *fileName, int includeResults, double angle,
1413 : : int numTri, int numQuad, int numEdge,
1414 : : int numVert, int* triEdge, int* quadEdge,
1415 : : int* edgeVert, double* vert,
1416 : : double* edgeCtrlPts, double* triCtrlPts,
1417 : : double* quadCtrlPts)
1418 : : {
1419 : 0 : FILE *fp = fopen(fileName, "w");
1420 [ # # ]: 0 : assert(fp != NULL); // couldn't open file for writing
1421 : :
1422 : : // write the header info
1423 : :
1424 : 0 : fprintf(fp, "Cholla version 1.0\n");
1425 : 0 : time_stamp( fp );
1426 : 0 : fprintf(fp, "%d %d %d %d %d\n", numTri, numQuad, numEdge, numVert, includeResults);
1427 : 0 : fprintf(fp, "%f\n", angle);
1428 : :
1429 : : // faceEdgeBegin
1430 : :
1431 : : int ii;
1432 : :
1433 : : // faceEdge
1434 : :
1435 : 0 : fprintf(fp, "triEdge\n");
1436 [ # # ]: 0 : for(ii=0; ii<numTri; ii++)
1437 : : {
1438 : 0 : fprintf(fp, "%d %d %d\n", triEdge[ii*3], triEdge[ii*3+1], triEdge[ii*3+2]);
1439 : : }
1440 : :
1441 : 0 : fprintf(fp, "quadEdge\n");
1442 [ # # ]: 0 : for(ii=0; ii<numQuad; ii++)
1443 : : {
1444 : 0 : fprintf(fp, "%d %d %d %d\n", quadEdge[ii*4], quadEdge[ii*4+1],
1445 : 0 : quadEdge[ii*4+2], quadEdge[ii*4+3]);
1446 : : }
1447 : :
1448 : : // edgeVert
1449 : :
1450 : 0 : fprintf(fp, "edgeVert\n");
1451 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1452 : : {
1453 : 0 : fprintf(fp, "%d %d\n", edgeVert[ii*2], edgeVert[ii*2+1]);
1454 : : }
1455 : :
1456 : : // vert
1457 : :
1458 : 0 : fprintf(fp, "vert\n");
1459 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
1460 : : {
1461 : 0 : fprintf(fp, "%f %f %f\n", vert[ii*3], vert[ii*3+1], vert[ii*3+2]);
1462 : : }
1463 : :
1464 [ # # ]: 0 : if (includeResults != 0)
1465 : : {
1466 : : // edge control points
1467 : :
1468 : : int jj;
1469 : 0 : fprintf(fp, "edgeCtrlPts\n");
1470 [ # # ]: 0 : for(ii=0; ii<numEdge; ii++)
1471 : : {
1472 [ # # ]: 0 : for (jj=0; jj<NUM_EDGE_CPTS; jj++)
1473 : : {
1474 : 0 : fprintf(fp, "%f %f %f\n", edgeCtrlPts[ii*jj*3],
1475 : 0 : edgeCtrlPts[ii*jj*3+1],
1476 : 0 : edgeCtrlPts[ii*jj*3+2]);
1477 : : }
1478 : : }
1479 : :
1480 : : // triangle control points
1481 : :
1482 : 0 : fprintf(fp, "triCtrlPts\n");
1483 [ # # ]: 0 : for(ii=0; ii<numTri; ii++)
1484 : : {
1485 [ # # ]: 0 : for (jj=0; jj<NUM_TRI_CPTS; jj++)
1486 : : {
1487 : 0 : fprintf(fp, "%f %f %f\n", triCtrlPts[ii*jj*3],
1488 : 0 : triCtrlPts[ii*jj*3+1],
1489 : 0 : triCtrlPts[ii*jj*3+2]);
1490 : : }
1491 : : }
1492 : :
1493 : : // quad control points
1494 : :
1495 : 0 : fprintf(fp, "quadCtrlPts\n");
1496 [ # # ]: 0 : for(ii=0; ii<numQuad; ii++)
1497 : : {
1498 [ # # ]: 0 : for (jj=0; jj<NUM_QUAD_CPTS; jj++)
1499 : : {
1500 : 0 : fprintf(fp, "%f %f %f\n", quadCtrlPts[ii*jj*3],
1501 : 0 : quadCtrlPts[ii*jj*3+1],
1502 : 0 : quadCtrlPts[ii*jj*3+2]);
1503 : : }
1504 : : }
1505 : : }
1506 : 0 : fclose(fp);
1507 : 0 : }
1508 : :
1509 : : //===========================================================================
1510 : : // Function: dumpResults
1511 : : // Purpose: dump only the results to a file
1512 : : // Date: 10/28/2002
1513 : : // Author: sjowen
1514 : : //===========================================================================
1515 : 0 : void dumpResults(const char *fileName, int numEdge, int numTri, int numQuad,
1516 : : double* edgeCtrlPts, double* triCtrlPts,
1517 : : double* quadCtrlPts )
1518 : : {
1519 : 0 : FILE *fp = fopen( fileName, "w" );
1520 [ # # ]: 0 : assert(fp != NULL);
1521 : :
1522 : : // write the header info
1523 : :
1524 : 0 : fprintf(fp, "Cholla version 1.0 Results\n");
1525 : 0 : time_stamp( fp );
1526 : 0 : fprintf(fp, "%d %d %d\n", numEdge, numTri, numQuad );
1527 : :
1528 : : int ii,jj;
1529 : 0 : fprintf(fp, "edgeCtrlPts\n");
1530 [ # # ]: 0 : for(ii=0; ii<numEdge; ii++)
1531 : : {
1532 [ # # ]: 0 : for (jj=0; jj<NUM_EDGE_CPTS; jj++)
1533 : : {
1534 : 0 : fprintf(fp, "%f %f %f\n", edgeCtrlPts[ii*jj*3],
1535 : 0 : edgeCtrlPts[ii*jj*3+1],
1536 : 0 : edgeCtrlPts[ii*jj*3+2]);
1537 : : }
1538 : : }
1539 : :
1540 : : // triangle control points
1541 : :
1542 : 0 : fprintf(fp, "triCtrlPts\n");
1543 [ # # ]: 0 : for(ii=0; ii<numTri; ii++)
1544 : : {
1545 [ # # ]: 0 : for (jj=0; jj<NUM_TRI_CPTS; jj++)
1546 : : {
1547 : 0 : fprintf(fp, "%f %f %f\n", triCtrlPts[ii*jj*3],
1548 : 0 : triCtrlPts[ii*jj*3+1],
1549 : 0 : triCtrlPts[ii*jj*3+2]);
1550 : : }
1551 : : }
1552 : :
1553 : : // quad control points
1554 : :
1555 : 0 : fprintf(fp, "quadCtrlPts\n");
1556 [ # # ]: 0 : for(ii=0; ii<numQuad; ii++)
1557 : : {
1558 [ # # ]: 0 : for (jj=0; jj<NUM_QUAD_CPTS; jj++)
1559 : : {
1560 : 0 : fprintf(fp, "%f %f %f\n", quadCtrlPts[ii*jj*3],
1561 : 0 : quadCtrlPts[ii*jj*3+1],
1562 : 0 : quadCtrlPts[ii*jj*3+2]);
1563 : : }
1564 : : }
1565 : :
1566 : 0 : fclose(fp);
1567 : 0 : }
1568 : :
1569 : : //===========================================================================
1570 : : // Function: importMesh
1571 : : // Purpose: import the face mesh from the specified file
1572 : : // allocate the arrays and pass them back
1573 : : // Date: 11/28/2002
1574 : : // Author: sjowen
1575 : : //===========================================================================
1576 : 0 : void importMesh(const char *fileName, int *includeResults,
1577 : : double *angle, int *numTri, int *numQuad, int *numEdge,
1578 : : int *numVert, int** triEdge, int** quadEdge,
1579 : : int** edgeVert, double** vert,
1580 : : double** edgeCtrlPts, double** triCtrlPts,
1581 : : double** quadCtrlPts)
1582 : : {
1583 [ # # ]: 0 : FILE *fp = fopen(fileName, "r");
1584 [ # # ]: 0 : assert(fp != NULL); // couldn't open file for reading
1585 : :
1586 : : // read the header info
1587 : :
1588 : : char version[128];
1589 [ # # ]: 0 : char* s = fgets(version, 128, fp);
1590 [ # # ]: 0 : if (NULL == s) {
1591 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1592 [ # # ]: 0 : fclose(fp);
1593 : 0 : return;
1594 : : }
1595 [ # # ]: 0 : assert(strcmp(version, "Cholla version 1.0\n") == 0); // don't recognize file
1596 : : char filetime[128];
1597 [ # # ]: 0 : s = fgets(filetime, 128, fp);
1598 [ # # ]: 0 : if (NULL == s) {
1599 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1600 [ # # ]: 0 : fclose(fp);
1601 : 0 : return;
1602 : : }
1603 : :
1604 : : char line[256];
1605 : : int num_tri;
1606 : : int num_quad;
1607 : : int num_edge;
1608 : : int num_vert;
1609 [ # # ]: 0 : s = fgets(line,256,fp);
1610 [ # # ]: 0 : if (NULL == s) {
1611 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1612 [ # # ]: 0 : fclose(fp);
1613 : 0 : return;
1614 : : }
1615 : 0 : sscanf(line, "%d %d %d %d %d\n", &num_tri, &num_quad, &num_edge, &num_vert, includeResults);
1616 : 0 : *numTri = num_tri;
1617 : 0 : *numQuad = num_quad;
1618 : 0 : *numEdge = num_edge;
1619 : 0 : *numVert = num_vert;
1620 [ # # ]: 0 : int nread = fscanf(fp, "%lf\n", angle);
1621 [ # # ]: 0 : if (1 != nread) {
1622 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1623 [ # # ]: 0 : fclose(fp);
1624 : 0 : return;
1625 : : }
1626 : :
1627 : : // triEdge
1628 : :
1629 [ # # ][ # # ]: 0 : *triEdge = new int [num_tri*3];
1630 : 0 : int *tri_edge = *triEdge;
1631 : : char array_name[128];
1632 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1633 [ # # ]: 0 : if (1 != nread) {
1634 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1635 [ # # ]: 0 : fclose(fp);
1636 : 0 : return;
1637 : : }
1638 [ # # ]: 0 : assert(strcmp(array_name, "triEdge") == 0); // check start of array
1639 : : int ii;
1640 [ # # ]: 0 : for(ii=0; ii<num_tri; ii++)
1641 : : {
1642 [ # # ]: 0 : nread = fscanf(fp, "%d %d %d\n", &tri_edge[ii*3], &tri_edge[ii*3+1], &tri_edge[ii*3+2]);
1643 [ # # ]: 0 : if (3 != nread) {
1644 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1645 [ # # ]: 0 : fclose(fp);
1646 : 0 : return;
1647 : : }
1648 : : }
1649 : :
1650 : : // quadEdge
1651 : :
1652 [ # # ][ # # ]: 0 : *quadEdge = new int [num_quad*4];
1653 : 0 : int *quad_edge = *quadEdge;
1654 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1655 [ # # ]: 0 : if (1 != nread) {
1656 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1657 [ # # ]: 0 : fclose(fp);
1658 : 0 : return;
1659 : : }
1660 [ # # ]: 0 : assert(strcmp(array_name, "quadEdge") == 0); // check start of array
1661 [ # # ]: 0 : for(ii=0; ii<num_quad; ii++)
1662 : : {
1663 : 0 : nread = fscanf(fp, "%d %d %d %d\n", &quad_edge[ii*4], &quad_edge[ii*4+1],
1664 [ # # ]: 0 : &quad_edge[ii*4+2], &quad_edge[ii*4+3]);
1665 [ # # ]: 0 : if (4 != nread) {
1666 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1667 [ # # ]: 0 : fclose(fp);
1668 : 0 : return;
1669 : : }
1670 : : }
1671 : :
1672 : : // edgeVert
1673 : :
1674 [ # # ][ # # ]: 0 : *edgeVert = new int [2*num_edge];
1675 : 0 : int *edge_vert = *edgeVert;
1676 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1677 [ # # ]: 0 : if (1 != nread) {
1678 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1679 [ # # ]: 0 : fclose(fp);
1680 : 0 : return;
1681 : : }
1682 [ # # ]: 0 : assert(strcmp(array_name, "edgeVert") == 0); // check start of array
1683 [ # # ]: 0 : for (ii=0; ii<num_edge; ii++)
1684 : : {
1685 [ # # ]: 0 : nread = fscanf(fp, "%d %d\n", &edge_vert[ii*2], &edge_vert[ii*2+1]);
1686 [ # # ]: 0 : if (2 != nread) {
1687 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1688 [ # # ]: 0 : fclose(fp);
1689 : 0 : return;
1690 : : }
1691 : : }
1692 : :
1693 : : // vert
1694 : :
1695 [ # # ][ # # ]: 0 : *vert = new double [3*num_vert];
1696 : 0 : double *my_vert = *vert;
1697 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1698 [ # # ]: 0 : if (1 != nread) {
1699 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1700 [ # # ]: 0 : fclose(fp);
1701 : 0 : return;
1702 : : }
1703 [ # # ]: 0 : assert(strcmp(array_name, "vert") == 0); // check start of array
1704 [ # # ]: 0 : for (ii=0; ii<num_vert; ii++)
1705 : : {
1706 [ # # ]: 0 : nread = fscanf(fp, "%lf %lf %lf\n", &my_vert[ii*3], &my_vert[ii*3+1], &my_vert[ii*3+2]);
1707 [ # # ]: 0 : if (3 != nread) {
1708 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1709 [ # # ]: 0 : fclose(fp);
1710 : 0 : return;
1711 : : }
1712 : : }
1713 : :
1714 [ # # ][ # # ]: 0 : *edgeCtrlPts = new double [3*num_edge*NUM_EDGE_CPTS];
1715 [ # # ]: 0 : if (numTri != 0)
1716 [ # # ][ # # ]: 0 : *triCtrlPts = new double [3*num_tri*NUM_TRI_CPTS];
1717 [ # # ]: 0 : if (numQuad != 0)
1718 [ # # ][ # # ]: 0 : *quadCtrlPts = new double [3*num_quad*NUM_QUAD_CPTS];
1719 : :
1720 : 0 : double *edge_ctrl_pts = *edgeCtrlPts;
1721 : 0 : double *tri_ctrl_pts = *triCtrlPts;
1722 : 0 : double *quad_ctrl_pts = *quadCtrlPts;
1723 [ # # ]: 0 : if (*includeResults != 0)
1724 : : {
1725 : : // edge control points
1726 : :
1727 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1728 [ # # ]: 0 : if (1 != nread) {
1729 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1730 [ # # ]: 0 : fclose(fp);
1731 : 0 : return;
1732 : : }
1733 [ # # ]: 0 : assert(strcmp(array_name, "edgeCtrlPts") == 0); // check start of array
1734 [ # # ]: 0 : for(ii=0; ii<num_edge*NUM_EDGE_CPTS; ii++)
1735 : : {
1736 : : nread = fscanf(fp, "%lf %lf %lf\n", &edge_ctrl_pts[ii*3],
1737 : 0 : &edge_ctrl_pts[ii*3+1],
1738 [ # # ]: 0 : &edge_ctrl_pts[ii*3+2]);
1739 [ # # ]: 0 : if (3 != nread) {
1740 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1741 [ # # ]: 0 : fclose(fp);
1742 : 0 : return;
1743 : : }
1744 : : }
1745 : :
1746 : : // triangle control points
1747 : :
1748 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1749 [ # # ]: 0 : if (1 != nread) {
1750 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1751 [ # # ]: 0 : fclose(fp);
1752 : 0 : return;
1753 : : }
1754 [ # # ]: 0 : assert(strcmp(array_name, "triCtrlPts") == 0); // check start of array
1755 [ # # ]: 0 : for(ii=0; ii<num_tri*NUM_TRI_CPTS; ii++)
1756 : : {
1757 : : nread = fscanf(fp, "%lf %lf %lf\n", &tri_ctrl_pts[ii*3],
1758 : 0 : &tri_ctrl_pts[ii*3+1],
1759 [ # # ]: 0 : &tri_ctrl_pts[ii*3+2]);
1760 [ # # ]: 0 : if (3 != nread) {
1761 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1762 [ # # ]: 0 : fclose(fp);
1763 : 0 : return;
1764 : : }
1765 : : }
1766 : :
1767 : : // quad control points
1768 : :
1769 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1770 [ # # ]: 0 : if (1 != nread) {
1771 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1772 [ # # ]: 0 : fclose(fp);
1773 : 0 : return;
1774 : : }
1775 [ # # ]: 0 : assert(strcmp(array_name, "quadCtrlPts") == 0); // check start of array
1776 [ # # ]: 0 : for(ii=0; ii<num_quad*NUM_QUAD_CPTS; ii++)
1777 : : {
1778 : : nread = fscanf(fp, "%lf %lf %lf\n", &quad_ctrl_pts[ii*3],
1779 : 0 : &quad_ctrl_pts[ii*3+1],
1780 [ # # ]: 0 : &quad_ctrl_pts[ii*3+2]);
1781 [ # # ]: 0 : if (3 != nread) {
1782 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1783 [ # # ]: 0 : fclose(fp);
1784 : 0 : return;
1785 : : }
1786 : : }
1787 : :
1788 : : }
1789 [ # # ]: 0 : fclose(fp);
1790 : : }
1791 : :
1792 : : //===========================================================================
1793 : : // Function: importResults
1794 : : // Purpose: import only the results to from a file
1795 : : // Note: allocates the normal and tangent arrays
1796 : : // Date: 10/28/2002
1797 : : // Author: sjowen
1798 : : //===========================================================================
1799 : 0 : void importResults(const char *fileName, int *numEdge, int *numTri,
1800 : : int *numQuad, double** edgeCtrlPts, double** triCtrlPts,
1801 : : double** quadCtrlPts )
1802 : : {
1803 [ # # ]: 0 : FILE *fp = fopen(fileName, "r");
1804 [ # # ]: 0 : assert(fp != NULL); // couldn't open file for reading
1805 : :
1806 : : // read the header info
1807 : :
1808 : : char version[128];
1809 [ # # ]: 0 : char* s = fgets(version, 128, fp);
1810 [ # # ]: 0 : if (NULL == s) {
1811 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1812 [ # # ]: 0 : fclose(fp);
1813 : 0 : return;
1814 : : }
1815 [ # # ]: 0 : assert(strcmp(version, "Cholla version 1.0 Results\n") == 0); // don't recognize file
1816 : : char filetime[128];
1817 [ # # ]: 0 : s = fgets(filetime, 128, fp);
1818 [ # # ]: 0 : if (NULL == s) {
1819 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in mesh file %s\n", fileName);
[ # # ][ # # ]
1820 [ # # ]: 0 : fclose(fp);
1821 : 0 : return;
1822 : : }
1823 [ # # ]: 0 : int nread = fscanf(fp, "%d %d %d\n", numEdge, numTri, numQuad);
1824 [ # # ]: 0 : if (3 != nread) {
1825 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1826 [ # # ]: 0 : fclose(fp);
1827 : 0 : return;
1828 : : }
1829 : 0 : int num_edge = *numEdge;
1830 : 0 : int num_tri = *numTri;
1831 : 0 : int num_quad = *numQuad;
1832 : :
1833 [ # # ][ # # ]: 0 : *edgeCtrlPts = new double [3*num_edge*NUM_EDGE_CPTS];
1834 [ # # ]: 0 : if (numTri != 0)
1835 [ # # ][ # # ]: 0 : *triCtrlPts = new double [3*num_tri*NUM_TRI_CPTS];
1836 [ # # ]: 0 : if (numQuad != 0)
1837 [ # # ][ # # ]: 0 : *quadCtrlPts = new double [3*num_quad*NUM_QUAD_CPTS];
1838 : :
1839 : 0 : double *edge_ctrl_pts = *edgeCtrlPts;
1840 : 0 : double *tri_ctrl_pts = *triCtrlPts;
1841 : 0 : double *quad_ctrl_pts = *quadCtrlPts;
1842 : :
1843 : : // edge control points
1844 : :
1845 : : int ii;
1846 : : char array_name[128];
1847 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1848 [ # # ]: 0 : if (1 != nread) {
1849 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1850 [ # # ]: 0 : fclose(fp);
1851 : 0 : return;
1852 : : }
1853 [ # # ]: 0 : assert(strcmp(array_name, "edgeCtrlPts") == 0); // check start of array
1854 [ # # ]: 0 : for(ii=0; ii<num_edge*NUM_EDGE_CPTS; ii++)
1855 : : {
1856 : : nread = fscanf(fp, "%lf %lf %lf\n", &edge_ctrl_pts[ii*3],
1857 : 0 : &edge_ctrl_pts[ii*3+1],
1858 [ # # ]: 0 : &edge_ctrl_pts[ii*3+2]);
1859 [ # # ]: 0 : if (3 != nread) {
1860 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1861 [ # # ]: 0 : fclose(fp);
1862 : 0 : return;
1863 : : }
1864 : : }
1865 : :
1866 : : // triangle control points
1867 : :
1868 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1869 [ # # ]: 0 : if (1 != nread) {
1870 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1871 [ # # ]: 0 : fclose(fp);
1872 : 0 : return;
1873 : : }
1874 [ # # ]: 0 : assert(strcmp(array_name, "triCtrlPts") == 0); // check start of array
1875 [ # # ]: 0 : for(ii=0; ii<num_tri*NUM_TRI_CPTS; ii++)
1876 : : {
1877 : : nread = fscanf(fp, "%lf %lf %lf\n", &tri_ctrl_pts[ii*3],
1878 : 0 : &tri_ctrl_pts[ii*3+1],
1879 [ # # ]: 0 : &tri_ctrl_pts[ii*3+2]);
1880 [ # # ]: 0 : if (3 != nread) {
1881 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1882 [ # # ]: 0 : fclose(fp);
1883 : 0 : return;
1884 : : }
1885 : : }
1886 : :
1887 : : // quad control points
1888 : :
1889 [ # # ]: 0 : nread = fscanf(fp, "%s\n", array_name);
1890 [ # # ]: 0 : if (1 != nread) {
1891 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1892 [ # # ]: 0 : fclose(fp);
1893 : 0 : return;
1894 : : }
1895 [ # # ]: 0 : assert(strcmp(array_name, "quadCtrlPts") == 0); // check start of array
1896 [ # # ]: 0 : for(ii=0; ii<num_tri*NUM_QUAD_CPTS; ii++)
1897 : : {
1898 : : nread = fscanf(fp, "%lf %lf %lf\n", &quad_ctrl_pts[ii*3],
1899 : 0 : &quad_ctrl_pts[ii*3+1],
1900 [ # # ]: 0 : &quad_ctrl_pts[ii*3+2]);
1901 [ # # ]: 0 : if (3 != nread) {
1902 [ # # ][ # # ]: 0 : PRINT_ERROR("Format error in file %s\n", fileName);
[ # # ][ # # ]
1903 [ # # ]: 0 : fclose(fp);
1904 : 0 : return;
1905 : : }
1906 : : }
1907 : :
1908 [ # # ]: 0 : fclose(fp);
1909 : :
1910 : : }
1911 : :
1912 : : //===========================================================================
1913 : : // Function: time_stamp
1914 : : // Purpose: write the current time to a file
1915 : : // Date: 11/28/2002
1916 : : // Author: sjowen
1917 : : //===========================================================================
1918 : 0 : static void time_stamp( FILE *fp )
1919 : : {
1920 : : struct tm *newtime;
1921 : 0 : bool am = true;
1922 : : time_t long_time;
1923 : :
1924 : 0 : time( &long_time ); /* Get time as long integer. */
1925 : 0 : newtime = localtime( &long_time ); /* Convert to local time. */
1926 : :
1927 [ # # ]: 0 : if( newtime->tm_hour > 12 ) /* Set up extension. */
1928 : 0 : am = false;
1929 [ # # ]: 0 : if( newtime->tm_hour > 12 ) /* Convert from 24-hour */
1930 : 0 : newtime->tm_hour -= 12; /* to 12-hour clock. */
1931 [ # # ]: 0 : if( newtime->tm_hour == 0 ) /*Set hour to 12 if midnight. */
1932 : 0 : newtime->tm_hour = 12;
1933 : :
1934 [ # # ][ # # ]: 0 : fprintf( fp, "%.19s %s\n", asctime( newtime ), am ? "AM" : "PM" );
1935 : 0 : }
1936 : :
1937 : : //===========================================================================
1938 : : // Function: evalBezierEdge
1939 : : // Purpose: Evaluate a quartic Bezier curve at a specified parametric
1940 : : // location given its control points. Evaluates location, and/or tangent.
1941 : :
1942 : : // Date: 01/20/2004
1943 : : // Author: jfowler
1944 : : //===========================================================================
1945 : 0 : void evalBezierEdge( int numEdge, int numVert, int* edgeVert,
1946 : : double* vert, double* edgeCtrlPts,
1947 : : int numLocs, double* paramLocation,
1948 : : double* location, double* tangent )
1949 : : {
1950 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
1951 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
1952 : 0 : CubitFacet **facet_array = NULL;
1953 : 0 : CubitQuadFacet **qfacet_array = NULL;
1954 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
1955 : : int ii, index, numTri, numQuad;
1956 : : int *triEdge, *quadEdge;
1957 : 0 : triEdge = quadEdge = NULL;
1958 : 0 : numTri = numQuad = 0;
1959 : 0 : CubitStatus stat = CUBIT_SUCCESS;
1960 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
1961 : : edgeVert, vert, point_array, edge_array, facet_array,
1962 [ # # ]: 0 : qfacet_array, facet_list);
1963 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
1964 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
1965 : 0 : return;
1966 : : }
1967 : :
1968 : : CubitFacetEdge *edge_ptr;
1969 : : // Set the control points on the edges.
1970 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
1971 : : {
1972 : 0 : edge_ptr = edge_array[ii];
1973 : 0 : index = 3 * (ii * NUM_EDGE_CPTS);
1974 [ # # ]: 0 : edge_ptr->set_control_points( &(edgeCtrlPts[index]) );
1975 : : }
1976 : :
1977 : : // do the evaluations
1978 : : int jj, tindex;
1979 : 0 : index = tindex = 0;
1980 : :
1981 [ # # ]: 0 : for ( ii = 0; ii < numEdge; ii++ ) {
1982 [ # # ]: 0 : for ( jj = 0; jj < numLocs; jj++ ) {
1983 : :
1984 [ # # ]: 0 : CubitVector eval_point;
1985 : 0 : CubitVector *eval_point_ptr = &eval_point;
1986 [ # # ]: 0 : CubitVector eval_tangent;
1987 : 0 : CubitVector *eval_tangent_ptr = NULL;
1988 [ # # ]: 0 : if ( tangent != NULL )
1989 : 0 : eval_tangent_ptr = &eval_tangent;
1990 : 0 : edge_ptr = edge_array[ii];
1991 [ # # ]: 0 : edge_ptr->evaluate_single(paramLocation[jj],eval_point_ptr);
1992 [ # # ]: 0 : location[index++] = eval_point.x();
1993 [ # # ]: 0 : location[index++] = eval_point.y();
1994 [ # # ]: 0 : location[index++] = eval_point.z();
1995 [ # # ]: 0 : if ( tangent != NULL ) {
1996 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(paramLocation[jj],eval_tangent_ptr);
1997 [ # # ]: 0 : eval_tangent.normalize();
1998 [ # # ]: 0 : tangent[tindex++] = eval_tangent.x();
1999 [ # # ]: 0 : tangent[tindex++] = eval_tangent.y();
2000 [ # # ]: 0 : tangent[tindex++] = eval_tangent.z();
2001 : : }
2002 : : }
2003 : : }
2004 : : // delete the temp arrays
2005 : :
2006 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2007 [ # # ][ # # ]: 0 : delete edge_array[ii];
2008 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
2009 [ # # ][ # # ]: 0 : delete point_array[ii];
2010 [ # # ]: 0 : delete [] point_array;
2011 [ # # ][ # # ]: 0 : delete [] edge_array;
[ # # ]
2012 : :
2013 : : }
2014 : :
2015 : : //===========================================================================
2016 : : // Function: evalBezierEdgeFromTans
2017 : : // Purpose: Evaluate a quartic Bezier curve at a specified parametric
2018 : : // location given its tangents at the end-points. Evaluates location,
2019 : : // and/or tangent.
2020 : : // Date: 01/20/2004
2021 : : // Author: jfowler
2022 : : //===========================================================================
2023 : 0 : void evalBezierEdgeFromTans( int numEdge, int numVert, int* edgeVert,
2024 : : double* vert, double* edgeVertTang,
2025 : : int numLocs, double* paramLocation,
2026 : : double* location, double* tangent )
2027 : : {
2028 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
2029 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
2030 : 0 : CubitFacet **facet_array = NULL;
2031 : 0 : CubitQuadFacet **qfacet_array = NULL;
2032 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
2033 : : int ii, index, numTri, numQuad;
2034 : : int *triEdge, *quadEdge;
2035 : 0 : triEdge = quadEdge = NULL;
2036 : 0 : numTri = numQuad = 0;
2037 : 0 : CubitStatus stat = CUBIT_SUCCESS;
2038 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
2039 : : edgeVert, vert, point_array, edge_array, facet_array,
2040 [ # # ]: 0 : qfacet_array, facet_list);
2041 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
2042 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
2043 : 0 : return;
2044 : : }
2045 : :
2046 : : // set the control points on edges
2047 : :
2048 : 0 : int mydebug = 0;
2049 : 0 : FILE *fp = NULL;
2050 [ # # ]: 0 : if (mydebug)
2051 [ # # ]: 0 : fp = fopen("edges.txt", "w");
2052 : : CubitFacetEdge *edge_ptr;
2053 [ # # ][ # # ]: 0 : CubitVector T0, T1, P0, P1;
[ # # ][ # # ]
2054 [ # # ][ # # ]: 0 : CubitVector Pi[3];
2055 : : CubitPoint *pt0, *pt1;
2056 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2057 : : {
2058 : 0 : index = 6 * ii;
2059 : 0 : edge_ptr = edge_array[ii];
2060 [ # # ]: 0 : pt0 = edge_ptr->point(0);
2061 [ # # ]: 0 : pt1 = edge_ptr->point(1);
2062 [ # # ][ # # ]: 0 : P0 = pt0->coordinates();
2063 [ # # ][ # # ]: 0 : P1 = pt1->coordinates();
2064 [ # # ]: 0 : T0.x( edgeVertTang[ index ] );
2065 [ # # ]: 0 : T0.y( edgeVertTang[ index + 1] );
2066 [ # # ]: 0 : T0.z( edgeVertTang[ index + 2] );
2067 [ # # ]: 0 : T1.x( edgeVertTang[ index + 3] );
2068 [ # # ]: 0 : T1.y( edgeVertTang[ index + 4] );
2069 [ # # ]: 0 : T1.z( edgeVertTang[ index + 5] );
2070 [ # # ]: 0 : FacetEvalTool::init_edge_control_points_single(P0, P1, T0, T1, Pi);
2071 [ # # ]: 0 : if (mydebug)
2072 : : {
2073 : : int kk;
2074 [ # # ]: 0 : for(kk=0; kk<3; kk++);
2075 [ # # ][ # # ]: 0 : fprintf(fp, "(%d) %.6f %.6f %.6f\n", edge_ptr->id(), P0.x(), P0.y(), P0.z());
[ # # ][ # # ]
[ # # ]
2076 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[0].x(), Pi[0].y(), Pi[0].z());
[ # # ][ # # ]
2077 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[1].x(), Pi[1].y(), Pi[1].z());
[ # # ][ # # ]
2078 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", Pi[2].x(), Pi[2].y(), Pi[2].z());
[ # # ][ # # ]
2079 [ # # ][ # # ]: 0 : fprintf(fp, " %.6f %.6f %.6f\n", P1.x(), P1.y(), P1.z());
[ # # ][ # # ]
2080 : : }
2081 [ # # ]: 0 : edge_ptr->control_points( Pi, 4 );
2082 : : }
2083 [ # # ]: 0 : if (mydebug)
2084 [ # # ]: 0 : fclose(fp);
2085 : :
2086 : : // do the evaluations
2087 : : int jj, tindex;
2088 : 0 : index = tindex = 0;
2089 : :
2090 [ # # ]: 0 : for ( ii = 0; ii < numEdge; ii++ ) {
2091 [ # # ]: 0 : for ( jj = 0; jj < numLocs; jj++ ) {
2092 : :
2093 [ # # ]: 0 : CubitVector eval_point;
2094 : 0 : CubitVector *eval_point_ptr = &eval_point;
2095 [ # # ]: 0 : CubitVector eval_tangent;
2096 : 0 : CubitVector *eval_tangent_ptr = NULL;
2097 [ # # ]: 0 : if ( tangent != NULL )
2098 : 0 : eval_tangent_ptr = &eval_tangent;
2099 : 0 : edge_ptr = edge_array[ii];
2100 [ # # ]: 0 : edge_ptr->evaluate_single(paramLocation[jj],eval_point_ptr);
2101 [ # # ]: 0 : location[index++] = eval_point.x();
2102 [ # # ]: 0 : location[index++] = eval_point.y();
2103 [ # # ]: 0 : location[index++] = eval_point.z();
2104 [ # # ]: 0 : if ( tangent != NULL ) {
2105 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(paramLocation[jj],eval_tangent_ptr);
2106 [ # # ]: 0 : eval_tangent.normalize();
2107 [ # # ]: 0 : tangent[tindex++] = eval_tangent.x();
2108 [ # # ]: 0 : tangent[tindex++] = eval_tangent.y();
2109 [ # # ]: 0 : tangent[tindex++] = eval_tangent.z();
2110 : : }
2111 : : }
2112 : : }
2113 : : // delete the temp arrays
2114 : :
2115 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2116 [ # # ][ # # ]: 0 : delete edge_array[ii];
2117 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
2118 [ # # ][ # # ]: 0 : delete point_array[ii];
2119 [ # # ]: 0 : delete [] point_array;
2120 [ # # ][ # # ]: 0 : delete [] edge_array;
[ # # ]
2121 : : }
2122 : :
2123 : : //===========================================================================
2124 : : // Function: projBezierEdgeFromTans
2125 : : // Purpose: Project a point to an edge given tangents at the endpoints.
2126 : : // The routine only finds one point. Rarely there might be
2127 : : // two or more closest points.
2128 : : // Date: 01/20/2004
2129 : : // Author: jfowler
2130 : : //===========================================================================
2131 : 0 : void projToBezierEdgeFromTans( int numEdge, int numVert, int* edgeVert,
2132 : : double* vert, double* edgeVertTang,
2133 : : int numLocs, double* xyz,
2134 : : double* xyzOnEdge, double* tangent )
2135 : : {
2136 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
2137 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
2138 : 0 : CubitFacet **facet_array = NULL;
2139 : 0 : CubitQuadFacet **qfacet_array = NULL;
2140 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
2141 : : int ii, index, numTri, numQuad;
2142 : : int *triEdge, *quadEdge;
2143 : 0 : triEdge = quadEdge = NULL;
2144 : 0 : numTri = numQuad = 0;
2145 : 0 : CubitStatus stat = CUBIT_SUCCESS;
2146 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
2147 : : edgeVert, vert, point_array, edge_array, facet_array,
2148 [ # # ]: 0 : qfacet_array, facet_list);
2149 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
2150 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
2151 : 0 : return;
2152 : : }
2153 : :
2154 : : // set the control points on edges
2155 : :
2156 : : CubitFacetEdge *edge_ptr;
2157 [ # # ][ # # ]: 0 : CubitVector T0, T1, P0, P1;
[ # # ][ # # ]
2158 [ # # ][ # # ]: 0 : CubitVector Pi[3];
2159 : : CubitPoint *pt0, *pt1;
2160 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2161 : : {
2162 : 0 : index = 6 * ii;
2163 : 0 : edge_ptr = edge_array[ii];
2164 [ # # ]: 0 : pt0 = edge_ptr->point(0);
2165 [ # # ]: 0 : pt1 = edge_ptr->point(1);
2166 [ # # ][ # # ]: 0 : P0 = pt0->coordinates();
2167 [ # # ][ # # ]: 0 : P1 = pt1->coordinates();
2168 [ # # ]: 0 : T0.x( edgeVertTang[ index ] );
2169 [ # # ]: 0 : T0.y( edgeVertTang[ index + 1] );
2170 [ # # ]: 0 : T0.z( edgeVertTang[ index + 2] );
2171 [ # # ]: 0 : T1.x( edgeVertTang[ index + 3] );
2172 [ # # ]: 0 : T1.y( edgeVertTang[ index + 4] );
2173 [ # # ]: 0 : T1.z( edgeVertTang[ index + 5] );
2174 [ # # ]: 0 : FacetEvalTool::init_edge_control_points_single(P0, P1, T0, T1, Pi);
2175 : :
2176 [ # # ]: 0 : edge_ptr->control_points( Pi, 4 );
2177 : : }
2178 : :
2179 : : // do the projections
2180 : : int jj, kk, xyzindex, tindex;
2181 : 0 : index = tindex = 0;
2182 : 0 : double dsquared, dsquared_test, dderiv, dderiv2, t, t_min = 0.0;
2183 [ # # ]: 0 : CubitVector xyz_pt;//, second_d;
2184 : 0 : t = 0.0;
2185 : 0 : dsquared_test = CUBIT_DBL_MAX;
2186 [ # # ]: 0 : for ( ii = 0; ii < numEdge; ii++ ) {
2187 : 0 : edge_ptr = edge_array[ii];
2188 [ # # ]: 0 : for ( jj = 0; jj < numLocs; jj++ ) {
2189 : 0 : xyzindex = 3*jj;
2190 [ # # ]: 0 : xyz_pt.x(xyz[xyzindex]);
2191 [ # # ]: 0 : xyz_pt.y(xyz[xyzindex+1]);
2192 [ # # ]: 0 : xyz_pt.z(xyz[xyzindex+2]);
2193 : 0 : t = -1.0;
2194 : 0 : dsquared_test = CUBIT_DBL_MAX;
2195 : : // Get square of distance from xyz to curve at intervals of 0.2
2196 : : // in the parameter t. Save the closest distance found as a
2197 : : // starting point for Newton-Raphson iterations. It was found
2198 : : // that an interval of 0.5 was too coarse and that the N-R
2199 : : // sometimes converged on a false value.
2200 [ # # ]: 0 : CubitVector eval_point;
2201 : 0 : CubitVector *eval_point_ptr = &eval_point;
2202 [ # # ]: 0 : CubitVector eval_tangent;
2203 : 0 : CubitVector *eval_tangent_ptr = &eval_tangent;
2204 [ # # ]: 0 : CubitVector second_d;
2205 : 0 : CubitVector *second_d_ptr = &second_d;
2206 [ # # ]: 0 : for ( kk = 0; kk < 11; kk++ ) {
2207 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2208 [ # # ][ # # ]: 0 : dsquared = ( (xyz_pt.x()-eval_point.x())*(xyz_pt.x()-eval_point.x()) +
[ # # ][ # # ]
2209 [ # # ][ # # ]: 0 : (xyz_pt.y()-eval_point.y())*(xyz_pt.y()-eval_point.y()) +
[ # # ][ # # ]
2210 [ # # ][ # # ]: 0 : (xyz_pt.z()-eval_point.z())*(xyz_pt.z()-eval_point.z()) );
[ # # ][ # # ]
2211 [ # # ]: 0 : if ( fabs(dsquared) < fabs(dsquared_test) ) {
2212 : 0 : t_min = t;
2213 : 0 : dsquared_test = dsquared;
2214 : : }
2215 : 0 : t += 0.2;
2216 : : }
2217 : : double dderiva, dderivb;
2218 [ # # ]: 0 : if ( t_min == -1.00 ) {
2219 : : // Check whether the slope changed signs between -1.0 and -0.8.
2220 : : // If so, the min must have occurred in this interval -- set
2221 : : // t_min to -0.9.
2222 : 0 : t = -1.0;
2223 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2224 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2225 [ # # ][ # # ]: 0 : dderiva = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2226 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2227 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2228 : 0 : t = -0.8;
2229 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2230 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2231 [ # # ][ # # ]: 0 : dderivb = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2232 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2233 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2234 [ # # ]: 0 : if ( dderiva*dderivb < 0.0 ) t_min = -0.9;
2235 : :
2236 [ # # ]: 0 : } else if ( t_min == 1.00 ) {
2237 : : // Check the other end of the range, similarly.
2238 : 0 : t = 1.0;
2239 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2240 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2241 [ # # ][ # # ]: 0 : dderiva = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2242 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2243 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2244 : 0 : t = 0.8;
2245 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2246 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2247 [ # # ][ # # ]: 0 : dderivb = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2248 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2249 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2250 [ # # ]: 0 : if ( dderiva*dderivb < 0.0 ) t_min = 0.9;
2251 : : }
2252 : 0 : t = t_min;
2253 [ # # ][ # # ]: 0 : if ( (t > -1.0) && (t < 1.0) ) {
2254 : : int mm;
2255 : 0 : mm = 0;
2256 [ # # ]: 0 : while ( mm < 10 ) { // to avoid possible infinite loop
2257 : 0 : mm++;
2258 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2259 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2260 [ # # ][ # # ]: 0 : dderiv = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2261 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2262 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2263 : :
2264 [ # # ]: 0 : edge_ptr->evaluate_2nd_derivative(t, second_d_ptr);
2265 : :
2266 [ # # ][ # # ]: 0 : dderiv2 = (eval_point.x()-xyz_pt.x())*second_d.x() +
[ # # ]
2267 [ # # ][ # # ]: 0 : eval_tangent.x()*eval_tangent.x() +
2268 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*second_d.y() +
[ # # ]
2269 [ # # ][ # # ]: 0 : eval_tangent.y()*eval_tangent.y() +
2270 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*second_d.z() +
[ # # ]
2271 [ # # ][ # # ]: 0 : eval_tangent.z()*eval_tangent.z();
2272 : :
2273 : 0 : t = t - dderiv/dderiv2; // Newton-Raphson
2274 : :
2275 [ # # ]: 0 : if ( t < -1.0 ) {
2276 : 0 : t = -1.0;
2277 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2278 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2279 : 0 : break;
2280 : : }
2281 [ # # ]: 0 : if ( t > 1.0 ) {
2282 : 0 : t = 1.0;
2283 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2284 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2285 : 0 : break;
2286 : : }
2287 [ # # ]: 0 : if ( fabs(dderiv) < 1.e-11 ) break;
2288 : 0 : }
2289 : : } else { // At an endpoint of the paramaterization.
2290 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2291 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2292 : : }
2293 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.x();
2294 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.y();
2295 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.z();
2296 [ # # ]: 0 : if ( tangent != NULL ) {
2297 [ # # ]: 0 : eval_tangent.normalize();
2298 [ # # ]: 0 : tangent[tindex++] = eval_tangent.x();
2299 [ # # ]: 0 : tangent[tindex++] = eval_tangent.y();
2300 [ # # ]: 0 : tangent[tindex++] = eval_tangent.z();
2301 : : }
2302 : : }
2303 : : }
2304 : :
2305 : : // delete the temp arrays
2306 : :
2307 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2308 [ # # ][ # # ]: 0 : delete edge_array[ii];
2309 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
2310 [ # # ][ # # ]: 0 : delete point_array[ii];
2311 [ # # ]: 0 : delete [] point_array;
2312 [ # # ][ # # ]: 0 : delete [] edge_array;
[ # # ]
2313 : :
2314 : : }
2315 : :
2316 : :
2317 : : //===========================================================================
2318 : : // Function: projBezierEdge
2319 : : // Purpose: Project a point to an edge.
2320 : : // The routine only finds one point. Rarely there might be
2321 : : // two or more closest points.
2322 : : // Date: 01/20/2004
2323 : : // Author: jfowler
2324 : : //===========================================================================
2325 : 0 : void projToBezierEdge( int numEdge, int numVert, int* edgeVert,
2326 : : double* vert, double* edgeCtrlPts,
2327 : : int numLocs, double* xyz,
2328 : : double* xyzOnEdge, double* tangent, double* t_value )
2329 : : {
2330 [ # # ][ # # ]: 0 : CubitPoint **point_array = new CubitPoint * [numVert];
2331 [ # # ][ # # ]: 0 : CubitFacetEdge **edge_array = new CubitFacetEdge * [numEdge];
2332 : 0 : CubitFacet **facet_array = NULL;
2333 : 0 : CubitQuadFacet **qfacet_array = NULL;
2334 [ # # ]: 0 : DLIList<FacetEntity *> facet_list;
2335 : : int ii, index, numTri, numQuad;
2336 : : int *triEdge, *quadEdge;
2337 : 0 : triEdge = quadEdge = NULL;
2338 : 0 : numTri = numQuad = 0;
2339 : 0 : CubitStatus stat = CUBIT_SUCCESS;
2340 : : stat = build_facets(numTri, numQuad, numEdge, numVert, triEdge, quadEdge,
2341 : : edgeVert, vert, point_array, edge_array, facet_array,
2342 [ # # ]: 0 : qfacet_array, facet_list);
2343 [ # # ]: 0 : if (CUBIT_SUCCESS != stat) {
2344 [ # # ][ # # ]: 0 : PRINT_ERROR("Failed to build facets.\n");
[ # # ][ # # ]
2345 : 0 : return;
2346 : : }
2347 : :
2348 : : // set the control points on edges
2349 : :
2350 : : CubitFacetEdge *edge_ptr;
2351 : : //CubitVector T0, T1, P0, P1;
2352 : : //CubitVector Pi[3];
2353 : :
2354 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2355 : : {
2356 : 0 : edge_ptr = edge_array[ii];
2357 : 0 : index = 3 * (ii * NUM_EDGE_CPTS);
2358 [ # # ]: 0 : edge_ptr->set_control_points( &(edgeCtrlPts[index]) );
2359 : : }
2360 : :
2361 : : // do the projections
2362 : : int jj, kk, xyzindex, tindex;
2363 : 0 : index = tindex = 0;
2364 : 0 : double dsquared, dsquared_test, dderiv, dderiv2, t, t_min = 0.0;
2365 [ # # ]: 0 : CubitVector xyz_pt;//, second2142_d;
2366 : 0 : t = 0.0;
2367 : 0 : dsquared_test = CUBIT_DBL_MAX;
2368 [ # # ]: 0 : for ( ii = 0; ii < numEdge; ii++ ) {
2369 : 0 : edge_ptr = edge_array[ii];
2370 [ # # ]: 0 : for ( jj = 0; jj < numLocs; jj++ ) {
2371 : 0 : xyzindex = 3*jj;
2372 [ # # ]: 0 : xyz_pt.x(xyz[xyzindex]);
2373 [ # # ]: 0 : xyz_pt.y(xyz[xyzindex+1]);
2374 [ # # ]: 0 : xyz_pt.z(xyz[xyzindex+2]);
2375 : 0 : t = -1.0;
2376 : 0 : dsquared_test = CUBIT_DBL_MAX;
2377 : : // Get square of distance from xyz to curve at intervals of 0.2
2378 : : // in the parameter t. Save the closest distance found as a
2379 : : // starting point for Newton-Raphson iterations. It was found
2380 : : // that an interval of 0.5 was too coarse and that the N-R
2381 : : // sometimes converged on a false value.
2382 [ # # ]: 0 : CubitVector eval_point;
2383 : 0 : CubitVector *eval_point_ptr = &eval_point;
2384 [ # # ]: 0 : CubitVector eval_tangent;
2385 : 0 : CubitVector *eval_tangent_ptr = &eval_tangent;
2386 [ # # ]: 0 : CubitVector second_d;
2387 : 0 : CubitVector *second_d_ptr = &second_d;
2388 [ # # ]: 0 : for ( kk = 0; kk < 11; kk++ ) {
2389 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2390 [ # # ][ # # ]: 0 : dsquared = ( (xyz_pt.x()-eval_point.x())*(xyz_pt.x()-eval_point.x()) +
[ # # ][ # # ]
2391 [ # # ][ # # ]: 0 : (xyz_pt.y()-eval_point.y())*(xyz_pt.y()-eval_point.y()) +
[ # # ][ # # ]
2392 [ # # ][ # # ]: 0 : (xyz_pt.z()-eval_point.z())*(xyz_pt.z()-eval_point.z()) );
[ # # ][ # # ]
2393 [ # # ]: 0 : if ( fabs(dsquared) < fabs(dsquared_test) ) {
2394 : 0 : t_min = t;
2395 : 0 : dsquared_test = dsquared;
2396 : : }
2397 : 0 : t += 0.2;
2398 : : }
2399 : : double dderiva, dderivb;
2400 [ # # ]: 0 : if ( t_min == -1.00 ) {
2401 : : // Check whether the slope changed signs between -1.0 and -0.8.
2402 : : // If so, the min must have occurred in this interval -- set
2403 : : // t_min to -0.9.
2404 : 0 : t = -1.0;
2405 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2406 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2407 [ # # ][ # # ]: 0 : dderiva = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2408 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2409 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2410 : 0 : t = -0.8;
2411 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2412 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2413 [ # # ][ # # ]: 0 : dderivb = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2414 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2415 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2416 [ # # ]: 0 : if ( dderiva*dderivb < 0.0 ) t_min = -0.9;
2417 : :
2418 [ # # ]: 0 : } else if ( t_min == 1.00 ) {
2419 : : // Check the other end of the range, similarly.
2420 : 0 : t = 1.0;
2421 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2422 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2423 [ # # ][ # # ]: 0 : dderiva = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2424 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2425 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2426 : 0 : t = 0.8;
2427 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2428 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2429 [ # # ][ # # ]: 0 : dderivb = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2430 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2431 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2432 [ # # ]: 0 : if ( dderiva*dderivb < 0.0 ) t_min = 0.9;
2433 : : }
2434 : 0 : t = t_min;
2435 [ # # ][ # # ]: 0 : if ( (t > -1.0) && (t < 1.0) ) {
2436 : : int mm;
2437 : 0 : mm = 0;
2438 [ # # ]: 0 : while ( mm < 10 ) { // to avoid possible infinite loop
2439 : 0 : mm++;
2440 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2441 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2442 [ # # ][ # # ]: 0 : dderiv = 2.0*( (eval_point.x()-xyz_pt.x())*eval_tangent.x() +
[ # # ]
2443 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*eval_tangent.y() +
[ # # ]
2444 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*eval_tangent.z() );
[ # # ]
2445 : :
2446 [ # # ]: 0 : edge_ptr->evaluate_2nd_derivative(t, second_d_ptr);
2447 : :
2448 [ # # ][ # # ]: 0 : dderiv2 = (eval_point.x()-xyz_pt.x())*second_d.x() +
[ # # ]
2449 [ # # ][ # # ]: 0 : eval_tangent.x()*eval_tangent.x() +
2450 [ # # ][ # # ]: 0 : (eval_point.y()-xyz_pt.y())*second_d.y() +
[ # # ]
2451 [ # # ][ # # ]: 0 : eval_tangent.y()*eval_tangent.y() +
2452 [ # # ][ # # ]: 0 : (eval_point.z()-xyz_pt.z())*second_d.z() +
[ # # ]
2453 [ # # ][ # # ]: 0 : eval_tangent.z()*eval_tangent.z();
2454 : :
2455 : 0 : t = t - dderiv/dderiv2; // Newton-Raphson
2456 : :
2457 [ # # ]: 0 : if ( t < -1.0 ) {
2458 : 0 : t = -1.0;
2459 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2460 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2461 : 0 : break;
2462 : : }
2463 [ # # ]: 0 : if ( t > 1.0 ) {
2464 : 0 : t = 1.0;
2465 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2466 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2467 : 0 : break;
2468 : : }
2469 [ # # ]: 0 : if ( fabs(dderiv) < 1.e-11 ) break;
2470 : 0 : }
2471 : : } else { // At an endpoint of the paramaterization.
2472 [ # # ]: 0 : edge_ptr->evaluate_single(t,eval_point_ptr);
2473 [ # # ]: 0 : edge_ptr->evaluate_single_tangent(t,eval_tangent_ptr);
2474 : : }
2475 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.x();
2476 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.y();
2477 [ # # ]: 0 : xyzOnEdge[index++] = eval_point.z();
2478 [ # # ]: 0 : if ( tangent != NULL ) {
2479 [ # # ]: 0 : eval_tangent.normalize();
2480 [ # # ]: 0 : tangent[tindex++] = eval_tangent.x();
2481 [ # # ]: 0 : tangent[tindex++] = eval_tangent.y();
2482 [ # # ]: 0 : tangent[tindex++] = eval_tangent.z();
2483 : : }
2484 : : }
2485 : : }
2486 : :
2487 : 0 : *t_value = t;
2488 : :
2489 : : // delete the temp arrays
2490 : :
2491 [ # # ]: 0 : for (ii=0; ii<numEdge; ii++)
2492 [ # # ][ # # ]: 0 : delete edge_array[ii];
2493 [ # # ]: 0 : for (ii=0; ii<numVert; ii++)
2494 [ # # ][ # # ]: 0 : delete point_array[ii];
2495 [ # # ]: 0 : delete [] point_array;
2496 [ # # ][ # # ]: 0 : delete [] edge_array;
[ # # ]
2497 : :
2498 : : }
2499 : :
2500 : : //===========================================================================
2501 : : // Function: constructQuadNormalsFromFile
2502 : : // Purpose: Similar to constructQuadNormals, except that instead of
2503 : : // approximating the surface normals and tangents from the input
2504 : : // facets, the normals and derivatives are extracted from a file
2505 : : // that was exported from Cubit along with the exodus file. The
2506 : : // purpose for this is that Cubit will have the exact surface
2507 : : // normals and tangents, instead of the approximations from the
2508 : : // geometry.
2509 : : // Date: 08/31/2004
2510 : : // Author: mlstate
2511 : : //===========================================================================
2512 : 0 : CubitStatus constructQuadNormalsFromFile( const char *filename,
2513 : : int numQuad,
2514 : : int numEdge,
2515 : : int* quadEdge,
2516 : : int* edgeVert,
2517 : : double* quadVertNorms,
2518 : : double* edgeVertTang )
2519 : : {
2520 : 0 : CubitStatus status = CUBIT_SUCCESS;
2521 : 0 : int num_file_quads = 0,
2522 : 0 : num_file_edges = 0,
2523 : 0 : *edge_nodes_file = NULL,
2524 : 0 : *quad_nodes_file = NULL;
2525 : 0 : double *edge_node_tangs_file = NULL,
2526 : 0 : *quad_node_norms_file = NULL;
2527 : :
2528 : : status = readMBGNormalFile( filename,
2529 : : NULL,
2530 : : &num_file_quads,
2531 : : &num_file_edges,
2532 : : &edge_nodes_file,
2533 : : &edge_node_tangs_file,
2534 : : NULL,
2535 : : NULL,
2536 : : &quad_nodes_file,
2537 [ # # ]: 0 : &quad_node_norms_file );
2538 : :
2539 [ # # ][ # # ]: 0 : if ( status == CUBIT_SUCCESS && numEdge > 0 )
2540 : : {
2541 : : status = extractEdgeTangsFromFile( num_file_edges,
2542 : : edge_nodes_file,
2543 : : edge_node_tangs_file,
2544 : : numEdge,
2545 : : edgeVert,
2546 [ # # ]: 0 : edgeVertTang );
2547 : : }
2548 [ # # ][ # # ]: 0 : if ( status == CUBIT_SUCCESS && numQuad > 0 )
2549 : : {
2550 : : status = extractQuadNormalsFromFile( num_file_quads,
2551 : : quad_nodes_file,
2552 : : quad_node_norms_file,
2553 : : numQuad,
2554 : : edgeVert,
2555 : : quadEdge,
2556 [ # # ]: 0 : quadVertNorms );
2557 : : }
2558 : :
2559 [ # # ]: 0 : delete [] edge_nodes_file;
2560 [ # # ]: 0 : delete [] quad_nodes_file;
2561 [ # # ]: 0 : delete [] edge_node_tangs_file;
2562 [ # # ]: 0 : delete [] quad_node_norms_file;
2563 : 0 : return status;
2564 : : }
2565 : :
2566 : : //===========================================================================
2567 : : // Function: constructTriNormalsFromFile
2568 : : // Purpose: Similar to constructTriNormals, except that instead of
2569 : : // approximating the surface normals and tangents from the input
2570 : : // facets, the normals and derivatives are extracted from a file
2571 : : // that was exported from Cubit along with the exodus file. The
2572 : : // purpose for this is that Cubit will have the exact surface
2573 : : // normals and tangents, instead of the approximations from the
2574 : : // geometry.
2575 : : // Date: 08/31/2004
2576 : : // Author: mlstate
2577 : : //===========================================================================
2578 : 0 : CubitStatus constructTriNormalsFromFile( const char *filename,
2579 : : int numTri,
2580 : : int numEdge,
2581 : : int* triEdge,
2582 : : int* edgeVert,
2583 : : double* triVertNorms,
2584 : : double* edgeVertTang )
2585 : : {
2586 : 0 : CubitStatus status = CUBIT_SUCCESS;
2587 : 0 : int num_file_tris = 0,
2588 : 0 : num_file_edges = 0,
2589 : 0 : *edge_nodes_file = NULL,
2590 : 0 : *tri_nodes_file = NULL;
2591 : 0 : double *edge_node_tangs_file = NULL,
2592 : 0 : *tri_node_norms_file = NULL;
2593 : :
2594 : : status = readMBGNormalFile( filename,
2595 : : &num_file_tris,
2596 : : NULL,
2597 : : &num_file_edges,
2598 : : &edge_nodes_file,
2599 : : &edge_node_tangs_file,
2600 : : &tri_nodes_file,
2601 : : &tri_node_norms_file,
2602 : : NULL,
2603 [ # # ]: 0 : NULL );
2604 : :
2605 [ # # ][ # # ]: 0 : if ( status == CUBIT_SUCCESS && numEdge > 0 )
2606 : : {
2607 : : status = extractEdgeTangsFromFile( num_file_edges,
2608 : : edge_nodes_file,
2609 : : edge_node_tangs_file,
2610 : : numEdge,
2611 : : edgeVert,
2612 [ # # ]: 0 : edgeVertTang );
2613 : : }
2614 [ # # ][ # # ]: 0 : if ( status == CUBIT_SUCCESS && numTri > 0 )
2615 : : {
2616 : : status = extractTriNormalsFromFile( num_file_tris,
2617 : : tri_nodes_file,
2618 : : tri_node_norms_file,
2619 : : numTri,
2620 : : edgeVert,
2621 : : triEdge,
2622 [ # # ]: 0 : triVertNorms );
2623 : : }
2624 : :
2625 [ # # ]: 0 : delete [] edge_nodes_file;
2626 [ # # ]: 0 : delete [] tri_nodes_file;
2627 [ # # ]: 0 : delete [] edge_node_tangs_file;
2628 [ # # ]: 0 : delete [] tri_node_norms_file;
2629 : 0 : return status;
2630 : : }
2631 : :
2632 : : //===========================================================================
2633 : : // Function: extractEdgeTangsFromFile
2634 : : // Purpose: extract from the file the edge tangents.
2635 : : // Date: 08/31/2004
2636 : : // Author: mlstate
2637 : : //===========================================================================
2638 : 0 : static CubitStatus extractEdgeTangsFromFile
2639 : : (
2640 : : int num_file_edges,
2641 : : int *edge_nodes_file,
2642 : : double *edge_node_tangs_file,
2643 : : int numEdge,
2644 : : int *edgeVert,
2645 : : double *edgeVertTang
2646 : : )
2647 : : {
2648 : : #define MIN_EDGE_NODE_ID( n1, n2 ) n1 < n2 ? n1 : n2
2649 : 0 : CubitStatus status = CUBIT_SUCCESS;
2650 [ # # ][ # # ]: 0 : DLIList<int *> *edge_hash = new DLIList<int*> [num_file_edges];
[ # # # # ]
2651 [ # # ]: 0 : int * edge_data = new int [num_file_edges * 3 ];
2652 : :
2653 : : int i;
2654 : :
2655 [ # # ]: 0 : for ( i = 0; i < num_file_edges; i++ )
2656 : : {
2657 : 0 : int *data = &(edge_data[i*3]);
2658 : :
2659 : 0 : data[0] = edge_nodes_file[i*2];
2660 : 0 : data[1] = edge_nodes_file[i*2 + 1];
2661 : 0 : data[2] = i;
2662 : :
2663 [ # # ]: 0 : int key = MIN_EDGE_NODE_ID( data[0], data[1] );
2664 : 0 : key %= num_file_edges;
2665 [ # # ]: 0 : edge_hash[key].append( data );
2666 : : }
2667 : :
2668 [ # # ]: 0 : for ( i = 0; i < numEdge; i++ )
2669 : : {
2670 : 0 : CubitBoolean found = CUBIT_FALSE;
2671 : 0 : int *this_edge = &(edgeVert[i*2]);
2672 [ # # ]: 0 : int key = MIN_EDGE_NODE_ID( this_edge[0], this_edge[1] );
2673 : 0 : key %= num_file_edges;
2674 [ # # ]: 0 : for ( int j = 0; j < edge_hash[key].size(); j++ )
2675 : : {
2676 : 0 : int *data = edge_hash[key].get_and_step();
2677 [ # # ][ # # ]: 0 : if ( data[0] == this_edge[0] &&
2678 : 0 : data[1] == this_edge[1] )
2679 : : {
2680 : 0 : edgeVertTang[ i*2*3 ] = edge_node_tangs_file[ data[2] * 2 * 3 ];
2681 : 0 : edgeVertTang[ i*2*3 + 1 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 1 ];
2682 : 0 : edgeVertTang[ i*2*3 + 2 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 2 ];
2683 : 0 : edgeVertTang[ i*2*3 + 3 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 3 ];
2684 : 0 : edgeVertTang[ i*2*3 + 4 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 4 ];
2685 : 0 : edgeVertTang[ i*2*3 + 5 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 5 ];
2686 : 0 : found = CUBIT_TRUE;
2687 : 0 : break;
2688 : : }
2689 [ # # ][ # # ]: 0 : else if ( data[0] == this_edge[1] &&
2690 : 0 : data[1] == this_edge[0] )
2691 : : {
2692 : 0 : edgeVertTang[ i*2*3 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 3 ];
2693 : 0 : edgeVertTang[ i*2*3 + 1 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 4 ];
2694 : 0 : edgeVertTang[ i*2*3 + 2 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 5 ];
2695 : 0 : edgeVertTang[ i*2*3 + 3 ] = edge_node_tangs_file[ data[2] * 2 * 3 ];
2696 : 0 : edgeVertTang[ i*2*3 + 4 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 1 ];
2697 : 0 : edgeVertTang[ i*2*3 + 5 ] = edge_node_tangs_file[ data[2] * 2 * 3 + 2 ];
2698 : 0 : found = CUBIT_TRUE;
2699 : 0 : break;
2700 : : }
2701 : : }
2702 [ # # ]: 0 : if ( found == CUBIT_FALSE )
2703 : : {
2704 : 0 : status = CUBIT_FAILURE;
2705 : 0 : break;
2706 : : }
2707 : : }
2708 [ # # ]: 0 : delete [] edge_data;
2709 [ # # ][ # # ]: 0 : delete [] edge_hash;
2710 : :
2711 [ # # ]: 0 : return status;
2712 : : }
2713 : :
2714 : : //===========================================================================
2715 : : // Function: extractTriNormalsFromFile
2716 : : // Purpose: extract from the file the normals for a block of tri elements.
2717 : : // Date: 08/31/2004
2718 : : // Author: mlstate
2719 : : //===========================================================================
2720 : 0 : static CubitStatus extractTriNormalsFromFile
2721 : : (
2722 : : int num_file_tris,
2723 : : int *tri_nodes_file,
2724 : : double *tri_node_norms_file,
2725 : : int numTri,
2726 : : int *edgeVert,
2727 : : int *triEdge,
2728 : : double *triVertNorms
2729 : : )
2730 : : {
2731 : : #define MIN_TRI_NODE_ID( n1, n2, n3 ) n1 < n2 ? n1 < n3 ? n1 : n3 : n2 < n3 ? n2 : n3
2732 : :
2733 : 0 : CubitStatus status = CUBIT_SUCCESS;
2734 [ # # ][ # # ]: 0 : DLIList< int* > *tri_hash = new DLIList< int*>[num_file_tris];
[ # # # # ]
2735 [ # # ]: 0 : int * tri_data = new int [num_file_tris * 4 ];
2736 : :
2737 : : int i;
2738 : :
2739 [ # # ]: 0 : for ( i = 0; i < num_file_tris; i++ )
2740 : : {
2741 : 0 : int *data = &(tri_data[i*4]);
2742 : :
2743 : 0 : data[0] = tri_nodes_file[i*3];
2744 : 0 : data[1] = tri_nodes_file[i*3 + 1];
2745 : 0 : data[2] = tri_nodes_file[i*3 + 2];
2746 : 0 : data[3] = i;
2747 : :
2748 [ # # ][ # # ]: 0 : int key = MIN_TRI_NODE_ID( data[0], data[1], data[2] );
[ # # ]
2749 : 0 : key %= num_file_tris;
2750 [ # # ]: 0 : tri_hash[key].append( data );
2751 : : }
2752 : :
2753 [ # # ]: 0 : for ( i = 0; i < numTri; i++ )
2754 : : {
2755 : 0 : CubitBoolean found = CUBIT_FALSE;
2756 : : int this_tri[3];
2757 : :
2758 : 0 : constructTriVerts( &triEdge[i*3], edgeVert, this_tri );
2759 : :
2760 [ # # ][ # # ]: 0 : int key = MIN_TRI_NODE_ID( this_tri[0], this_tri[1], this_tri[2] );
[ # # ]
2761 : 0 : key %= num_file_tris;
2762 [ # # ][ # # ]: 0 : for ( int j = 0; j < tri_hash[key].size(); j++ )
2763 : : {
2764 [ # # ]: 0 : int *data = tri_hash[key].get_and_step();
2765 : :
2766 [ # # ][ # # ]: 0 : if ( this_tri[0] == data[0] &&
2767 [ # # ]: 0 : this_tri[1] == data[1] &&
2768 : 0 : this_tri[2] == data[2] )
2769 : : {
2770 : 0 : triVertNorms[ i*3*3 ] = tri_node_norms_file[ data[3]*3*3 ];
2771 : 0 : triVertNorms[ i*3*3 + 1 ] = tri_node_norms_file[ data[3]*3*3 + 1 ];
2772 : 0 : triVertNorms[ i*3*3 + 2 ] = tri_node_norms_file[ data[3]*3*3 + 2 ];
2773 : 0 : triVertNorms[ i*3*3 + 3 ] = tri_node_norms_file[ data[3]*3*3 + 3 ];
2774 : 0 : triVertNorms[ i*3*3 + 4 ] = tri_node_norms_file[ data[3]*3*3 + 4 ];
2775 : 0 : triVertNorms[ i*3*3 + 5 ] = tri_node_norms_file[ data[3]*3*3 + 5 ];
2776 : 0 : triVertNorms[ i*3*3 + 6 ] = tri_node_norms_file[ data[3]*3*3 + 6 ];
2777 : 0 : triVertNorms[ i*3*3 + 7 ] = tri_node_norms_file[ data[3]*3*3 + 7 ];
2778 : 0 : triVertNorms[ i*3*3 + 8 ] = tri_node_norms_file[ data[3]*3*3 + 8 ];
2779 : 0 : found = CUBIT_TRUE;
2780 : 0 : break;
2781 : : }
2782 [ # # ][ # # ]: 0 : else if ( this_tri[0] == data[1] &&
2783 [ # # ]: 0 : this_tri[1] == data[2] &&
2784 : 0 : this_tri[2] == data[0] )
2785 : : {
2786 : 0 : triVertNorms[ i*3*3 ] = tri_node_norms_file[ data[3]*3*3 + 3 ];
2787 : 0 : triVertNorms[ i*3*3 + 1 ] = tri_node_norms_file[ data[3]*3*3 + 4 ];
2788 : 0 : triVertNorms[ i*3*3 + 2 ] = tri_node_norms_file[ data[3]*3*3 + 5 ];
2789 : 0 : triVertNorms[ i*3*3 + 3 ] = tri_node_norms_file[ data[3]*3*3 + 6 ];
2790 : 0 : triVertNorms[ i*3*3 + 4 ] = tri_node_norms_file[ data[3]*3*3 + 7 ];
2791 : 0 : triVertNorms[ i*3*3 + 5 ] = tri_node_norms_file[ data[3]*3*3 + 8 ];
2792 : 0 : triVertNorms[ i*3*3 + 6 ] = tri_node_norms_file[ data[3]*3*3 ];
2793 : 0 : triVertNorms[ i*3*3 + 7 ] = tri_node_norms_file[ data[3]*3*3 + 1 ];
2794 : 0 : triVertNorms[ i*3*3 + 8 ] = tri_node_norms_file[ data[3]*3*3 + 2 ];
2795 : 0 : found = CUBIT_TRUE;
2796 : 0 : break;
2797 : : }
2798 [ # # ][ # # ]: 0 : else if ( this_tri[0] == data[2] &&
2799 [ # # ]: 0 : this_tri[1] == data[0] &&
2800 : 0 : this_tri[2] == data[1] )
2801 : : {
2802 : 0 : triVertNorms[ i*3*3 ] = tri_node_norms_file[ data[3]*3*3 + 6 ];
2803 : 0 : triVertNorms[ i*3*3 + 1 ] = tri_node_norms_file[ data[3]*3*3 + 7 ];
2804 : 0 : triVertNorms[ i*3*3 + 2 ] = tri_node_norms_file[ data[3]*3*3 + 8 ];
2805 : 0 : triVertNorms[ i*3*3 + 3 ] = tri_node_norms_file[ data[3]*3*3 ];
2806 : 0 : triVertNorms[ i*3*3 + 4 ] = tri_node_norms_file[ data[3]*3*3 + 1 ];
2807 : 0 : triVertNorms[ i*3*3 + 5 ] = tri_node_norms_file[ data[3]*3*3 + 2 ];
2808 : 0 : triVertNorms[ i*3*3 + 6 ] = tri_node_norms_file[ data[3]*3*3 + 3 ];
2809 : 0 : triVertNorms[ i*3*3 + 7 ] = tri_node_norms_file[ data[3]*3*3 + 4 ];
2810 : 0 : triVertNorms[ i*3*3 + 8 ] = tri_node_norms_file[ data[3]*3*3 + 5 ];
2811 : 0 : found = CUBIT_TRUE;
2812 : 0 : break;
2813 : : }
2814 : : }
2815 [ # # ]: 0 : if ( found == CUBIT_FALSE )
2816 : : {
2817 : 0 : status = CUBIT_FAILURE;
2818 : 0 : break;
2819 : : }
2820 : : }
2821 : :
2822 [ # # ]: 0 : delete [] tri_data;
2823 [ # # ][ # # ]: 0 : delete [] tri_hash;
2824 [ # # ]: 0 : return status;
2825 : : }
2826 : :
2827 : : //===========================================================================
2828 : : // Function: extractQuadNormalsFromFile
2829 : : // Purpose: extract from the file the normals for a block of quad elements.
2830 : : // Date: 08/31/2004
2831 : : // Author: mlstate
2832 : : //===========================================================================
2833 : 0 : static CubitStatus extractQuadNormalsFromFile
2834 : : (
2835 : : int num_file_quads,
2836 : : int *quad_nodes_file,
2837 : : double *quad_node_norms_file,
2838 : : int numQuad,
2839 : : int *edgeVert,
2840 : : int *quadEdge,
2841 : : double *quadVertNorms
2842 : : )
2843 : : {
2844 : : #define MIN_QUAD_NODE_ID( n1, n2, n3, n4 ) MIN_EDGE_NODE_ID( (MIN_EDGE_NODE_ID( n1, n2 )), (MIN_EDGE_NODE_ID( n3, n4 )) )
2845 : :
2846 : 0 : CubitStatus status = CUBIT_SUCCESS;
2847 [ # # ][ # # ]: 0 : DLIList< int* > *quad_hash = new DLIList< int*>[num_file_quads];
[ # # # # ]
2848 [ # # ]: 0 : int * quad_data = new int [num_file_quads * 5 ];
2849 : :
2850 : : int i;
2851 [ # # ]: 0 : for ( i = 0; i < num_file_quads; i++ )
2852 : : {
2853 : 0 : int *data = &(quad_data[i*5]);
2854 : :
2855 : 0 : data[0] = quad_nodes_file[i*4];
2856 : 0 : data[1] = quad_nodes_file[i*4 + 1];
2857 : 0 : data[2] = quad_nodes_file[i*4 + 2];
2858 : 0 : data[3] = quad_nodes_file[i*4 + 3];
2859 : 0 : data[4] = i;
2860 : :
2861 [ # # ][ # # ]: 0 : int key = MIN_QUAD_NODE_ID( data[0], data[1], data[2], data[3] );
[ # # ][ # # ]
[ # # ]
2862 : 0 : key %= num_file_quads;
2863 [ # # ]: 0 : quad_hash[key].append( data );
2864 : : }
2865 : :
2866 [ # # ]: 0 : for ( i = 0; i < numQuad; i++ )
2867 : : {
2868 : 0 : CubitBoolean found = CUBIT_FALSE;
2869 : : int this_quad[4];
2870 : :
2871 : 0 : constructQuadVerts( &quadEdge[i*4], edgeVert, this_quad );
2872 : :
2873 [ # # ][ # # ]: 0 : int key = MIN_QUAD_NODE_ID( this_quad[0], this_quad[1], this_quad[2], this_quad[3] );
[ # # ][ # # ]
[ # # ]
2874 : 0 : key %= num_file_quads;
2875 [ # # ][ # # ]: 0 : for ( int j = 0; j < quad_hash[key].size(); j++ )
2876 : : {
2877 [ # # ]: 0 : int *data = quad_hash[key].get_and_step();
2878 : :
2879 [ # # ][ # # ]: 0 : if ( this_quad[0] == data[0] &&
2880 [ # # ]: 0 : this_quad[1] == data[1] &&
2881 [ # # ]: 0 : this_quad[2] == data[2] &&
2882 : 0 : this_quad[3] == data[3] )
2883 : : {
2884 : 0 : quadVertNorms[ i*4*3 ] = quad_node_norms_file[ data[4]*4*3 ];
2885 : 0 : quadVertNorms[ i*4*3 + 1 ] = quad_node_norms_file[ data[4]*4*3 + 1 ];
2886 : 0 : quadVertNorms[ i*4*3 + 2 ] = quad_node_norms_file[ data[4]*4*3 + 2 ];
2887 : 0 : quadVertNorms[ i*4*3 + 3 ] = quad_node_norms_file[ data[4]*4*3 + 3 ];
2888 : 0 : quadVertNorms[ i*4*3 + 4 ] = quad_node_norms_file[ data[4]*4*3 + 4 ];
2889 : 0 : quadVertNorms[ i*4*3 + 5 ] = quad_node_norms_file[ data[4]*4*3 + 5 ];
2890 : 0 : quadVertNorms[ i*4*3 + 6 ] = quad_node_norms_file[ data[4]*4*3 + 6 ];
2891 : 0 : quadVertNorms[ i*4*3 + 7 ] = quad_node_norms_file[ data[4]*4*3 + 7 ];
2892 : 0 : quadVertNorms[ i*4*3 + 8 ] = quad_node_norms_file[ data[4]*4*3 + 8 ];
2893 : 0 : quadVertNorms[ i*4*3 + 9 ] = quad_node_norms_file[ data[4]*4*3 + 9 ];
2894 : 0 : quadVertNorms[ i*4*3 + 10] = quad_node_norms_file[ data[4]*4*3 + 10];
2895 : 0 : quadVertNorms[ i*4*3 + 11] = quad_node_norms_file[ data[4]*4*3 + 11];
2896 : 0 : found = CUBIT_TRUE;
2897 : 0 : break;
2898 : : }
2899 [ # # ][ # # ]: 0 : else if ( this_quad[0] == data[1] &&
2900 [ # # ]: 0 : this_quad[1] == data[2] &&
2901 [ # # ]: 0 : this_quad[2] == data[3] &&
2902 : 0 : this_quad[3] == data[0] )
2903 : : {
2904 : 0 : quadVertNorms[ i*4*3 ] = quad_node_norms_file[ data[4]*4*3 + 3 ];
2905 : 0 : quadVertNorms[ i*4*3 + 1 ] = quad_node_norms_file[ data[4]*4*3 + 4 ];
2906 : 0 : quadVertNorms[ i*4*3 + 2 ] = quad_node_norms_file[ data[4]*4*3 + 5 ];
2907 : 0 : quadVertNorms[ i*4*3 + 3 ] = quad_node_norms_file[ data[4]*4*3 + 6 ];
2908 : 0 : quadVertNorms[ i*4*3 + 4 ] = quad_node_norms_file[ data[4]*4*3 + 7 ];
2909 : 0 : quadVertNorms[ i*4*3 + 5 ] = quad_node_norms_file[ data[4]*4*3 + 8 ];
2910 : 0 : quadVertNorms[ i*4*3 + 6 ] = quad_node_norms_file[ data[4]*4*3 + 9 ];
2911 : 0 : quadVertNorms[ i*4*3 + 7 ] = quad_node_norms_file[ data[4]*4*3 + 10 ];
2912 : 0 : quadVertNorms[ i*4*3 + 8 ] = quad_node_norms_file[ data[4]*4*3 + 11 ];
2913 : 0 : quadVertNorms[ i*4*3 + 9 ] = quad_node_norms_file[ data[4]*4*3 + 0 ];
2914 : 0 : quadVertNorms[ i*4*3 + 10] = quad_node_norms_file[ data[4]*4*3 + 1 ];
2915 : 0 : quadVertNorms[ i*4*3 + 11] = quad_node_norms_file[ data[4]*4*3 + 2 ];
2916 : 0 : found = CUBIT_TRUE;
2917 : 0 : break;
2918 : : }
2919 [ # # ][ # # ]: 0 : else if ( this_quad[0] == data[2] &&
2920 [ # # ]: 0 : this_quad[1] == data[3] &&
2921 [ # # ]: 0 : this_quad[2] == data[0] &&
2922 : 0 : this_quad[3] == data[1] )
2923 : : {
2924 : 0 : quadVertNorms[ i*4*3 ] = quad_node_norms_file[ data[4]*4*3 + 6 ];
2925 : 0 : quadVertNorms[ i*4*3 + 1 ] = quad_node_norms_file[ data[4]*4*3 + 7 ];
2926 : 0 : quadVertNorms[ i*4*3 + 2 ] = quad_node_norms_file[ data[4]*4*3 + 8 ];
2927 : 0 : quadVertNorms[ i*4*3 + 3 ] = quad_node_norms_file[ data[4]*4*3 + 9 ];
2928 : 0 : quadVertNorms[ i*4*3 + 4 ] = quad_node_norms_file[ data[4]*4*3 + 10 ];
2929 : 0 : quadVertNorms[ i*4*3 + 5 ] = quad_node_norms_file[ data[4]*4*3 + 11 ];
2930 : 0 : quadVertNorms[ i*4*3 + 6 ] = quad_node_norms_file[ data[4]*4*3 + 0 ];
2931 : 0 : quadVertNorms[ i*4*3 + 7 ] = quad_node_norms_file[ data[4]*4*3 + 1 ];
2932 : 0 : quadVertNorms[ i*4*3 + 8 ] = quad_node_norms_file[ data[4]*4*3 + 2 ];
2933 : 0 : quadVertNorms[ i*4*3 + 9 ] = quad_node_norms_file[ data[4]*4*3 + 3 ];
2934 : 0 : quadVertNorms[ i*4*3 + 10] = quad_node_norms_file[ data[4]*4*3 + 4 ];
2935 : 0 : quadVertNorms[ i*4*3 + 11] = quad_node_norms_file[ data[4]*4*3 + 5 ];
2936 : 0 : found = CUBIT_TRUE;
2937 : 0 : break;
2938 : : }
2939 [ # # ][ # # ]: 0 : else if ( this_quad[0] == data[3] &&
2940 [ # # ]: 0 : this_quad[1] == data[0] &&
2941 [ # # ]: 0 : this_quad[2] == data[1] &&
2942 : 0 : this_quad[3] == data[2] )
2943 : : {
2944 : 0 : quadVertNorms[ i*4*3 ] = quad_node_norms_file[ data[4]*4*3 + 9 ];
2945 : 0 : quadVertNorms[ i*4*3 + 1 ] = quad_node_norms_file[ data[4]*4*3 + 10 ];
2946 : 0 : quadVertNorms[ i*4*3 + 2 ] = quad_node_norms_file[ data[4]*4*3 + 11 ];
2947 : 0 : quadVertNorms[ i*4*3 + 3 ] = quad_node_norms_file[ data[4]*4*3 + 0 ];
2948 : 0 : quadVertNorms[ i*4*3 + 4 ] = quad_node_norms_file[ data[4]*4*3 + 1 ];
2949 : 0 : quadVertNorms[ i*4*3 + 5 ] = quad_node_norms_file[ data[4]*4*3 + 2 ];
2950 : 0 : quadVertNorms[ i*4*3 + 6 ] = quad_node_norms_file[ data[4]*4*3 + 3 ];
2951 : 0 : quadVertNorms[ i*4*3 + 7 ] = quad_node_norms_file[ data[4]*4*3 + 4 ];
2952 : 0 : quadVertNorms[ i*4*3 + 8 ] = quad_node_norms_file[ data[4]*4*3 + 5 ];
2953 : 0 : quadVertNorms[ i*4*3 + 9 ] = quad_node_norms_file[ data[4]*4*3 + 6 ];
2954 : 0 : quadVertNorms[ i*4*3 + 10] = quad_node_norms_file[ data[4]*4*3 + 7 ];
2955 : 0 : quadVertNorms[ i*4*3 + 11] = quad_node_norms_file[ data[4]*4*3 + 8 ];
2956 : 0 : found = CUBIT_TRUE;
2957 : 0 : break;
2958 : : }
2959 : : }
2960 [ # # ]: 0 : if ( found == CUBIT_FALSE )
2961 : : {
2962 : 0 : status = CUBIT_FAILURE;
2963 : 0 : break;
2964 : : }
2965 : : }
2966 : :
2967 [ # # ]: 0 : delete [] quad_data;
2968 [ # # ][ # # ]: 0 : delete [] quad_hash;
2969 [ # # ]: 0 : return status;
2970 : : }
2971 : :
2972 : : //===========================================================================
2973 : : // Function: constructTriVerts
2974 : : // Purpose: Given an array of edge indices for a tri element and an array
2975 : : // of edge to vertex info, construct an array containing the vertices
2976 : : // on the input tri.
2977 : : // Date: 08/31/2004
2978 : : // Author: mlstate
2979 : : //===========================================================================
2980 : 0 : static void constructTriVerts
2981 : : (
2982 : : int triEdge[3], // <I> Array of edges on this tri.
2983 : : int *edgeVert, // <I> Array of edge to vertex info.
2984 : : int triVert[3] // <O> The vertices on the input tri.
2985 : : )
2986 : : {
2987 : 0 : int *verts_edge1 = &(edgeVert[triEdge[0]*2]);
2988 : 0 : int *verts_edge2 = &(edgeVert[triEdge[1]*2]);
2989 : 0 : int *verts_edge3 = &(edgeVert[triEdge[2]*2]);
2990 : :
2991 [ # # ][ # # ]: 0 : if ( verts_edge1[0] != verts_edge2[0] &&
2992 : 0 : verts_edge1[0] != verts_edge2[1] )
2993 : : {
2994 : 0 : triVert[0] = verts_edge1[0];
2995 : 0 : triVert[1] = verts_edge1[1];
2996 : : }
2997 : : else
2998 : : {
2999 : 0 : triVert[0] = verts_edge1[1];
3000 : 0 : triVert[1] = verts_edge1[0];
3001 : : }
3002 [ # # ][ # # ]: 0 : if ( verts_edge3[0] != verts_edge2[0] &&
3003 : 0 : verts_edge3[0] != verts_edge2[1] )
3004 : : {
3005 : 0 : triVert[2] = verts_edge3[1];
3006 : : }
3007 : : else
3008 : : {
3009 : 0 : triVert[2] = verts_edge3[0];
3010 : : }
3011 : 0 : }
3012 : :
3013 : : //===========================================================================
3014 : : // Function: constructQuadVerts
3015 : : // Purpose: Given an array of edge indices for a quad element and an array
3016 : : // of edge to vertex info, construct an array containing the vertices
3017 : : // on the input quad.
3018 : : // Date: 08/31/2004
3019 : : // Author: mlstate
3020 : : //===========================================================================
3021 : 0 : static void constructQuadVerts
3022 : : (
3023 : : int QuadEdge[4], // <I> Array of edges on this quad.
3024 : : int *edgeVert, // <I> Array of edge to vertex info.
3025 : : int QuadVert[4] // <O> The vertices on the input quad.
3026 : : )
3027 : : {
3028 : 0 : int *verts_edge1 = &(edgeVert[QuadEdge[0]*2]);
3029 : 0 : int *verts_edge2 = &(edgeVert[QuadEdge[1]*2]);
3030 : 0 : int *verts_edge3 = &(edgeVert[QuadEdge[2]*2]);
3031 : 0 : int *verts_edge4 = &(edgeVert[QuadEdge[3]*2]);
3032 : :
3033 [ # # ][ # # ]: 0 : if ( verts_edge1[0] != verts_edge2[0] &&
3034 : 0 : verts_edge1[0] != verts_edge2[1] )
3035 : : {
3036 : 0 : QuadVert[0] = verts_edge1[0];
3037 : 0 : QuadVert[1] = verts_edge1[1];
3038 : : }
3039 : : else
3040 : : {
3041 : 0 : QuadVert[0] = verts_edge1[1];
3042 : 0 : QuadVert[1] = verts_edge1[0];
3043 : : }
3044 : :
3045 [ # # ][ # # ]: 0 : if ( verts_edge3[0] != verts_edge4[0] &&
3046 : 0 : verts_edge3[0] != verts_edge4[1] )
3047 : : {
3048 : 0 : QuadVert[2] = verts_edge3[0];
3049 : 0 : QuadVert[3] = verts_edge3[1];
3050 : : }
3051 : : else
3052 : : {
3053 : 0 : QuadVert[2] = verts_edge3[1];
3054 : 0 : QuadVert[3] = verts_edge3[0];
3055 : : }
3056 : 0 : }
3057 : :
3058 : : //===========================================================================
3059 : : // Function: readMBGNormalFile
3060 : : // Purpose: Read the raw data out of the Normal/Tangents file.
3061 : : // Date: 08/31/2004
3062 : : // Author: mlstate
3063 : : //===========================================================================
3064 : 0 : CubitStatus readMBGNormalFile
3065 : : (
3066 : : const char *filename,
3067 : : int *num_tris, // <O>
3068 : : int *num_quads, // <O>
3069 : : int *num_edges, // <O>
3070 : : int **edge_nodes, // <OF> len = 2 * num_edges
3071 : : double **edge_node_tangents, // <OF> len = 3 * 2 * num_edges
3072 : : int **tri_nodes, // <OF> len = 3 * num_tris
3073 : : double **tri_node_normals, // <OF> len = 3 * 3 * num_tris
3074 : : int **quad_nodes, // <OF> len = 4 * num_quads
3075 : : double **quad_node_normals // <OF> len = 3 * 4 * num_quads
3076 : : )
3077 : : {
3078 : : #define MAX_FILE_LINE 512
3079 : :
3080 : 0 : int num_t = 0;
3081 : 0 : int num_q = 0;
3082 : 0 : int curr_quad = 0;
3083 : 0 : int curr_tri = 0;
3084 : 0 : int num_e = 0;
3085 : : int n[4];
3086 : :
3087 [ # # ]: 0 : if ( num_tris ) *num_tris = 0;
3088 [ # # ]: 0 : if ( num_quads ) *num_quads = 0;
3089 [ # # ]: 0 : if ( num_edges ) *num_edges = 0;
3090 [ # # ]: 0 : if ( edge_nodes ) *edge_nodes = NULL;
3091 [ # # ]: 0 : if ( edge_node_tangents ) *edge_node_tangents = NULL;
3092 [ # # ]: 0 : if ( tri_nodes ) *tri_nodes = NULL;
3093 [ # # ]: 0 : if ( quad_nodes ) *quad_nodes = NULL;
3094 [ # # ]: 0 : if ( tri_node_normals ) *tri_node_normals = NULL;
3095 [ # # ]: 0 : if ( quad_node_normals ) *quad_node_normals = NULL;
3096 : :
3097 [ # # ]: 0 : if ( strlen( filename ) == 0 )
3098 : : {
3099 [ # # ][ # # ]: 0 : PRINT_ERROR(" No filename specified\n");
[ # # ][ # # ]
3100 : 0 : return CUBIT_FAILURE;
3101 : : }
3102 : :
3103 [ # # ]: 0 : ifstream geom_file( filename );
3104 : :
3105 : : char fileline[MAX_FILE_LINE];
3106 : :
3107 [ # # ][ # # ]: 0 : while ( geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3108 : : {
3109 : : int block_id,
3110 : : num_corners,
3111 : : num_elems;
3112 : :
3113 [ # # ]: 0 : if ( !strcmp( fileline, END_OF_BLOCKS ) ) break;
3114 [ # # ]: 0 : if ( 3 != sscanf( fileline, "BLK %d %d %d",
3115 : 0 : &block_id, &num_corners, &num_elems ) )
3116 : : {
3117 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3118 [ # # ]: 0 : geom_file.close();
3119 : 0 : return CUBIT_FAILURE;
3120 : : }
3121 : :
3122 [ # # ]: 0 : if ( num_corners == 3 )
3123 : : {
3124 : : checkMemoryAllocations( num_corners, num_elems, &num_t, tri_nodes,
3125 [ # # ]: 0 : tri_node_normals );
3126 : : }
3127 [ # # ]: 0 : else if ( num_corners == 4 )
3128 : : {
3129 : : checkMemoryAllocations( num_corners, num_elems, &num_q, quad_nodes,
3130 [ # # ]: 0 : quad_node_normals );
3131 : : }
3132 : : else
3133 : : {
3134 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file can only support quads and tris.\n" );
[ # # ][ # # ]
3135 [ # # ]: 0 : geom_file.close();
3136 : 0 : return CUBIT_FAILURE;
3137 : : }
3138 : :
3139 [ # # ]: 0 : for ( int ielem = 0; ielem < num_elems; ielem++ )
3140 : : {
3141 : 0 : int *this_nodes = NULL;
3142 : 0 : double *this_norms = NULL;
3143 : :
3144 [ # # ][ # # ]: 0 : if ( !geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3145 : : {
3146 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3147 [ # # ]: 0 : geom_file.close();
3148 : 0 : return CUBIT_FAILURE;
3149 : : }
3150 : :
3151 [ # # ]: 0 : if ( num_corners == 4 )
3152 : : {
3153 [ # # ]: 0 : if ( quad_nodes ) this_nodes = &((*quad_nodes )[4 * curr_quad]);
3154 [ # # ]: 0 : if ( quad_node_normals ) this_norms = &((*quad_node_normals)[4 * 3 * curr_quad]);
3155 : 0 : curr_quad++;
3156 : : }
3157 [ # # ]: 0 : else if ( num_corners == 3 )
3158 : : {
3159 [ # # ]: 0 : if ( tri_nodes ) this_nodes = &((*tri_nodes )[3 * curr_tri]);
3160 [ # # ]: 0 : if ( tri_node_normals ) this_norms = &((*tri_node_normals)[3 * 3 * curr_tri]);
3161 : 0 : curr_tri++;
3162 : : }
3163 [ # # ][ # # ]: 0 : if ( ( num_corners == 4 && 4 != sscanf( fileline, "%d %d %d %d",
3164 [ # # ][ # # ]: 0 : &n[0], &n[1], &n[2], &n[3] ) ) ||
3165 [ # # ]: 0 : ( num_corners == 3 && 3 != sscanf( fileline, "%d %d %d",
3166 : 0 : &n[0], &n[1], &n[2] ) ) )
3167 : : {
3168 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3169 [ # # ]: 0 : geom_file.close();
3170 : 0 : return CUBIT_FAILURE;
3171 : : }
3172 : :
3173 [ # # ]: 0 : if ( this_nodes )
3174 : 0 : memcpy( this_nodes, n, num_corners * sizeof( int ) );
3175 : :
3176 : : // Read in the surface normals at the vertices of this element.
3177 : :
3178 [ # # ]: 0 : for ( int icorner = 0; icorner < num_corners; icorner++ )
3179 : : {
3180 [ # # ][ # # ]: 0 : if ( !geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3181 : : {
3182 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3183 [ # # ]: 0 : geom_file.close();
3184 : 0 : return CUBIT_FAILURE;
3185 : : }
3186 [ # # ]: 0 : if ( this_norms )
3187 : : {
3188 : : double norm[3];
3189 : :
3190 [ # # ]: 0 : if ( 3 != sscanf( fileline, "%lf %lf %lf", &norm[0], &norm[1], &norm[2] ) )
3191 : : {
3192 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3193 [ # # ]: 0 : geom_file.close();
3194 : 0 : return CUBIT_FAILURE;
3195 : : }
3196 : 0 : memcpy( &(this_norms[icorner * 3]), norm, 3 * sizeof( double ) );
3197 : : }
3198 : : }
3199 : : }
3200 : : }
3201 : :
3202 [ # # ][ # # ]: 0 : if ( !geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3203 : : {
3204 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3205 [ # # ]: 0 : geom_file.close();
3206 : 0 : return CUBIT_FAILURE;
3207 : : }
3208 : :
3209 [ # # ]: 0 : if ( 1 != sscanf( fileline, "%d", &num_e ) )
3210 : : {
3211 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3212 [ # # ]: 0 : geom_file.close();
3213 : 0 : return CUBIT_FAILURE;
3214 : : }
3215 : :
3216 [ # # ][ # # ]: 0 : if ( edge_nodes || edge_node_tangents )
3217 : : {
3218 [ # # ][ # # ]: 0 : if ( edge_nodes ) *edge_nodes = new int [ 2 * num_e ];
[ # # ]
3219 [ # # ][ # # ]: 0 : if ( edge_node_tangents ) *edge_node_tangents = new double [3 * 2 * num_e];
[ # # ]
3220 : :
3221 [ # # ]: 0 : for ( int iedge = 0; iedge < num_e; iedge++ )
3222 : : {
3223 [ # # ][ # # ]: 0 : if ( !geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3224 : : {
3225 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3226 [ # # ]: 0 : geom_file.close();
3227 : 0 : return CUBIT_FAILURE;
3228 : : }
3229 [ # # ]: 0 : if ( 2 != sscanf( fileline, "%d %d", &n[0], &n[1] ) )
3230 : : {
3231 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3232 [ # # ]: 0 : geom_file.close();
3233 : 0 : return CUBIT_FAILURE;
3234 : : }
3235 : :
3236 [ # # ]: 0 : for ( int iend = 0; iend < 2; iend++ )
3237 : : {
3238 : : double d[3];
3239 : :
3240 [ # # ]: 0 : if ( edge_nodes )
3241 : : {
3242 : 0 : (*edge_nodes)[iedge * 2 + iend ] = n[iend];
3243 : : }
3244 : :
3245 [ # # ][ # # ]: 0 : if ( !geom_file.getline( fileline, MAX_FILE_LINE ) )
[ # # ]
3246 : : {
3247 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3248 [ # # ]: 0 : geom_file.close();
3249 : 0 : return CUBIT_FAILURE;
3250 : : }
3251 [ # # ]: 0 : if ( 3 != sscanf( fileline, "%lf %lf %lf", &d[0], &d[1], &d[2] ) )
3252 : : {
3253 [ # # ][ # # ]: 0 : PRINT_ERROR( "MBG Normal file has incorrect file format\n" );
[ # # ][ # # ]
3254 [ # # ]: 0 : geom_file.close();
3255 : 0 : return CUBIT_FAILURE;
3256 : : }
3257 : :
3258 : 0 : int index = ((iedge*2) + iend) * 3;
3259 [ # # ]: 0 : if ( edge_node_tangents )
3260 : : {
3261 : 0 : (*edge_node_tangents)[ index ] = d[0];
3262 : 0 : (*edge_node_tangents)[ index + 1 ] = d[1];
3263 : 0 : (*edge_node_tangents)[ index + 2 ] = d[2];
3264 : : }
3265 : : }
3266 : : }
3267 : : }
3268 : :
3269 [ # # ]: 0 : if ( num_quads ) *num_quads = num_q;
3270 [ # # ]: 0 : if ( num_tris ) *num_tris = num_t;
3271 [ # # ]: 0 : if ( num_edges ) *num_edges = num_e;
3272 : :
3273 [ # # ]: 0 : geom_file.close();
3274 [ # # ]: 0 : return CUBIT_SUCCESS;
3275 : : }
3276 : :
3277 : : //===========================================================================
3278 : : // Function: checkMemoryAllocations
3279 : : // Purpose: While reading data out of the Normal/Tangents file, reallocate
3280 : : // Arrays to ensure that we have enough memory.
3281 : : // Date: 08/31/2004
3282 : : // Author: mlstate
3283 : : //===========================================================================
3284 : 0 : static void checkMemoryAllocations
3285 : : (
3286 : : int num_nodes_per_elem,
3287 : : int additional_num_elems,
3288 : : int *num_elems,
3289 : : int **nodes,
3290 : : double **node_normals
3291 : : )
3292 : : {
3293 : 0 : int old_num = *num_elems;
3294 : :
3295 : 0 : *num_elems += additional_num_elems;
3296 [ # # ][ # # ]: 0 : if ( node_normals && *node_normals )
3297 : : {
3298 [ # # ]: 0 : double *new_mem = new double [ *num_elems * 3 * num_nodes_per_elem ];
3299 : 0 : memcpy( new_mem, *node_normals, old_num * 3 * num_nodes_per_elem * sizeof( double ) );
3300 : 0 : delete *node_normals;
3301 : 0 : *node_normals = new_mem;
3302 : : }
3303 [ # # ]: 0 : else if ( node_normals )
3304 : : {
3305 [ # # ]: 0 : *node_normals = new double [ *num_elems * num_nodes_per_elem * 3 ];
3306 : : }
3307 : :
3308 [ # # ][ # # ]: 0 : if ( nodes && *nodes )
3309 : : {
3310 [ # # ]: 0 : int *new_mem = new int [ *num_elems * num_nodes_per_elem ];
3311 : 0 : memcpy( new_mem, *nodes, old_num * num_nodes_per_elem * sizeof( int ) );
3312 : 0 : delete *nodes;
3313 : 0 : *nodes = new_mem;
3314 : : }
3315 [ # # ]: 0 : else if ( nodes )
3316 : : {
3317 [ # # ]: 0 : *nodes = new int [ *num_elems * num_nodes_per_elem ];
3318 : : }
3319 : 0 : }
3320 : :
3321 : : // This code classify the type of edge with respect to convexity,
3322 : : // C0 or C1 continuity, and topology. The returned int uses the following
3323 : : // number code:
3324 : : // -1 = non-convex, manifold edge
3325 : : // 0 = non-feature, manifold edge (ie, C1 edge)
3326 : : // 1 = convex, manifold edge
3327 : : // 2 = non-manifold edge
3328 : : // Numbers higher than two are error codes
3329 : : // 3 = edge only attached to a single facet (boundary of 2-d sheet).
3330 : 0 : static int classify_local_convexity_at_edge(CubitFacetEdge *edge_ptr)
3331 : : {
3332 : : //easy check... if not a feature return 0
3333 [ # # ][ # # ]: 0 : if(!edge_ptr->is_feature())
3334 : 0 : return 0;
3335 : :
3336 : : //else we know that we are at feature... determine the type
3337 : 0 : CubitFacet *adj_facet_1 = NULL, *adj_facet_2 = NULL;
3338 : 0 : CubitPoint *point_1 = NULL, *point_2 = NULL, *point_3 = NULL;
3339 [ # # ][ # # ]: 0 : CubitVector v2, v3;
3340 [ # # ][ # # ]: 0 : CubitVector v0 = edge_ptr->point(0)->coordinates();
3341 [ # # ][ # # ]: 0 : CubitVector v1 = edge_ptr->point(1)->coordinates();
3342 [ # # ]: 0 : DLIList<CubitFacet*> adj_facets;
3343 [ # # ]: 0 : edge_ptr->facets(adj_facets);
3344 : :
3345 : : //if non-manifold edge, mark as two
3346 [ # # ][ # # ]: 0 : if(adj_facets.size() > 2){
3347 : 0 : return 2;
3348 : : }
3349 : : //if boundary on 2-d (shouldn't happen) mark as 3
3350 : : //anything above 2 is essentially an error code...
3351 [ # # ][ # # ]: 0 : else if(adj_facets.size() < 2){
3352 : 0 : return 3;
3353 : : }
3354 : : //otherwise, determine whether it is 1 or -1
3355 : : //this else should be a separate function: classify_edge_convexity()
3356 [ # # ]: 0 : adj_facet_1 = adj_facets.get_and_step();
3357 [ # # ]: 0 : int edge_index = adj_facet_1->edge_index(edge_ptr);
3358 : : //order such that *_1 is the facet that uses
3359 : : // this edge as the forward edge...
3360 [ # # ][ # # ]: 0 : if(adj_facet_1->edge_use(edge_index) == 1)
3361 [ # # ]: 0 : adj_facet_2 = adj_facets.get();
3362 : : else{
3363 [ # # ]: 0 : adj_facet_1 = adj_facets.get();
3364 [ # # ]: 0 : adj_facets.reset();
3365 [ # # ]: 0 : adj_facet_2 = adj_facets.get();
3366 : : }
3367 : : // order the points so that the Jacobian of the matrix
3368 : : // will tell use whether we are convex or not at this edge.
3369 [ # # ]: 0 : adj_facet_1->points(point_1, point_2, point_3);
3370 [ # # ][ # # ]: 0 : if(point_1 == edge_ptr->point(0)){
3371 [ # # ][ # # ]: 0 : v2 = point_3->coordinates();
3372 : : }
3373 [ # # ][ # # ]: 0 : else if(point_1 == edge_ptr->point(1)){
3374 [ # # ][ # # ]: 0 : v2 = point_2->coordinates();
3375 : : }
3376 : : else{
3377 [ # # ][ # # ]: 0 : v2 = point_1->coordinates();
3378 : : }
3379 [ # # ]: 0 : adj_facet_2->points(point_1, point_2, point_3);
3380 [ # # ][ # # ]: 0 : if(point_1 == edge_ptr->point(0)){
3381 [ # # ][ # # ]: 0 : v3 = point_2->coordinates();
3382 : : }
3383 [ # # ][ # # ]: 0 : else if(point_1 == edge_ptr->point(1)){
3384 [ # # ][ # # ]: 0 : v3 = point_3->coordinates();
3385 : : }
3386 : : else{
3387 [ # # ][ # # ]: 0 : v3 = point_1->coordinates();
3388 : : }
3389 : :
3390 [ # # ]: 0 : v1-=v0;
3391 [ # # ]: 0 : v2-=v0;
3392 [ # # ]: 0 : v3-=v0;
3393 : : //if Jacobian is negative, we are convex here.
3394 [ # # ][ # # ]: 0 : if( ( v1 % (v2* v3) ) < 0.0){
[ # # ]
3395 : 0 : return 1;
3396 : : }
3397 : : else{
3398 : 0 : return -1;
3399 [ # # ]: 0 : }
3400 : :
3401 [ + - ][ + - ]: 6540 : }
3402 : : /*
3403 : : int test_edges(int num_edge, int* edge_vert, double* my_vert, CubitFacetEdge** edge_array){
3404 : : int ii;
3405 : : int return_val = 0.0;
3406 : : CubitFacetEdgeData *cfed_ptr = NULL;
3407 : : int my_v1, my_v2;
3408 : : CubitVector vert1, vert2;
3409 : : CubitPoint *poi1, *poi2;
3410 : : CubitVector coord1, coord2;
3411 : : for(ii=0; ii<num_edge; ++ii){
3412 : : cfed_ptr = CAST_TO(edge_array[ii], CubitFacetEdgeData);
3413 : : poi1 = cfed_ptr->start_node();
3414 : : poi2 = cfed_ptr->end_node();
3415 : : coord1.set(poi1->x(),poi1->y(),poi1->z());
3416 : : coord2.set(poi2->x(),poi2->y(),poi2->z());
3417 : : my_v1 = edge_vert[ii*2];
3418 : : my_v2 = edge_vert[ii*2+1];
3419 : : vert1.set(my_vert[my_v1*3],my_vert[my_v1*3+1],my_vert[my_v1*3+2]);
3420 : : vert2.set(my_vert[my_v2*3],my_vert[my_v2*3+1],my_vert[my_v2*3+2]);
3421 : : double my_dist = vert1.distance_between(coord1);
3422 : : if(my_dist > .01){
3423 : : PRINT_INFO("Edge (%i)(1), distance (%f)\n",ii,my_dist);
3424 : : return_val++;
3425 : : }
3426 : : my_dist = vert2.distance_between(coord2);
3427 : : if(my_dist > .01){
3428 : : PRINT_INFO("Edge (%i)(2), distance (%f)\n",ii,my_dist);
3429 : : return_val++;
3430 : : }
3431 : : }
3432 : : return return_val;
3433 : : }
3434 : :
3435 : :
3436 : : */
|