Branch data Line data Source code
1 : : #include "CubitFacet.hpp"
2 : : #include "CubitPoint.hpp"
3 : : #include "CubitVector.hpp"
4 : : #include "GeometryDefines.h"
5 : : #include "CubitPlane.hpp"
6 : : #include "CubitFacetEdge.hpp"
7 : : #include "CubitFacetData.hpp"
8 : : #include "FacetEvalTool.hpp"
9 : : #include "CastTo.hpp"
10 : : #include "FacetEntity.hpp"
11 : : #include "CubitPointData.hpp"
12 : : #include "CubitFacetEdgeData.hpp"
13 : : #include "TDFacetBoundaryPoint.hpp"
14 : : #include "GfxDebug.hpp"
15 : : #include "CubitMessage.hpp"
16 : :
17 : : #ifndef CUBIT_MAX
18 : : #define CUBIT_MAX(a,b) (((a) > (b)) ? (a) : (b))
19 : : #endif
20 : : #ifndef CUBIT_MIN
21 : : #define CUBIT_MIN(a,b) (((a) < (b)) ? (a) : (b))
22 : : #endif
23 : : #define min3(a,b,c) CUBIT_MIN((CUBIT_MIN((a),(b))),(c))
24 : : #define max3(a,b,c) CUBIT_MAX((CUBIT_MAX((a),(b))),(c))
25 : :
26 : : //===========================================================================
27 : : //Function Name: CubitFacet
28 : : //
29 : : //Member Type: PUBLIC
30 : : //Description: constructor
31 : : //===========================================================================
32 : 1584 : CubitFacet::CubitFacet( )
33 : : : cachedPlane(NULL), facetWeight(0.0), patchCtrlPts(NULL),
34 [ + - ]: 1584 : markedFlag(0), isFlat(999), isBackwards(0), toolID(0)
35 : :
36 : : {
37 : :
38 : 1584 : }
39 : :
40 : : //===========================================================================
41 : : //Function Name: ~CubitFacetData
42 : : //
43 : : //Member Type: PUBLIC
44 : : //Description: destructor
45 : : //===========================================================================
46 [ # # ]: 0 : CubitFacet::~CubitFacet()
47 : : {
48 [ # # ]: 0 : if (cachedPlane )
49 : 0 : delete cachedPlane;
50 : 0 : cachedPlane = NULL;
51 [ # # ]: 0 : if (patchCtrlPts)
52 [ # # ]: 0 : delete [] patchCtrlPts;
53 [ # # ]: 0 : }
54 : :
55 : :
56 : : //===========================================================================
57 : : //Function Name: normal
58 : : //
59 : : //Member Type: PUBLIC
60 : : //Description: return the facet normal
61 : : //===========================================================================
62 : 14102 : CubitVector CubitFacet::normal()
63 : : {
64 [ + - ][ + - ]: 14102 : CubitPlane fac_plane = plane();
65 [ + - ][ + - ]: 14102 : CubitVector the_normal = fac_plane.normal();
66 [ - + ]: 14102 : if (isBackwards)
67 [ # # ][ # # ]: 0 : the_normal = -the_normal;
68 : 14102 : return the_normal;
69 : : }
70 : :
71 : : //-------------------------------------------------------------------------
72 : : // Purpose : set/get the bezier patch internal control points.
73 : : //
74 : : // Special Notes : Allocate if necessary
75 : : //
76 : : // Creator : Steve Owen
77 : : //
78 : : // Creation Date : 06/28/00
79 : : //-------------------------------------------------------------------------
80 : 264 : void CubitFacet::set_control_points( CubitVector points[6] )
81 : : {
82 [ + + ]: 264 : if (!patchCtrlPts) {
83 [ + - ][ + + ]: 924 : patchCtrlPts = new CubitVector [6];
84 : : }
85 : 264 : memcpy( patchCtrlPts, points, 6 *sizeof(CubitVector) );
86 : 264 : }
87 : :
88 : 0 : void CubitFacet::set_control_points( const double *pt_array )
89 : : {
90 [ # # ]: 0 : if (!patchCtrlPts) {
91 [ # # ][ # # ]: 0 : patchCtrlPts = new CubitVector [6];
92 : : }
93 : : int ii;
94 [ # # ]: 0 : for (ii=0; ii<6; ii++)
95 : : {
96 : 0 : patchCtrlPts[ii].x(pt_array[3*ii]);
97 : 0 : patchCtrlPts[ii].y(pt_array[3*ii+1]);
98 : 0 : patchCtrlPts[ii].z(pt_array[3*ii+2]);
99 : : }
100 : 0 : }
101 : :
102 : 396 : void CubitFacet::get_control_points( CubitVector points[6] )
103 : : {
104 [ - + ]: 396 : assert(patchCtrlPts != 0);
105 : 396 : memcpy( points, patchCtrlPts, 6 *sizeof(CubitVector) );
106 : 396 : }
107 : :
108 : 15686 : const CubitPlane &CubitFacet::plane()
109 : : {
110 [ + + ]: 15686 : if( ! cachedPlane )
111 : : {
112 : : CubitPoint *p0, *p1, *p2;
113 [ + - ]: 1584 : points(p0, p1, p2);
114 [ + - ][ + - ]: 1584 : CubitVector v1 = p1->coordinates() - p0->coordinates();
[ + - ]
115 [ + - ][ + - ]: 1584 : CubitVector v2 = p2->coordinates() - p0->coordinates();
[ + - ]
116 [ + - ][ + - ]: 1584 : cachedPlane = new CubitPlane( v1 * v2, p0->coordinates() );
[ + - ][ + - ]
117 : : }
118 : 15686 : return *cachedPlane;
119 : : }
120 : :
121 : 264 : void CubitFacet::update_plane()
122 : : {
123 [ - + ]: 264 : if ( ! cachedPlane )
124 : 264 : return;
125 : :
126 : : CubitPoint *p0, *p1, *p2;
127 [ + - ]: 264 : points(p0, p1, p2);
128 [ + - ][ + - ]: 264 : CubitVector v1 = p1->coordinates() - p0->coordinates();
[ + - ]
129 [ + - ][ + - ]: 264 : CubitVector v2 = p2->coordinates() - p0->coordinates();
[ + - ]
130 [ + - ]: 264 : CubitVector normal = v1 * v2;
131 [ + - ][ - + ]: 264 : if (is_backwards()) normal = -normal;
[ # # ][ # # ]
132 [ + - ][ + - ]: 264 : cachedPlane->set(normal, p0->coordinates() );
133 : : }
134 : :
135 : : //-------------------------------------------------------------------------
136 : : // Purpose : determine if the facet is flat (all normal are the same)
137 : : //
138 : : // Special Notes :
139 : : //
140 : : // Creator : Steve Owen
141 : : //
142 : : // Creation Date : 5/4/01
143 : : //-------------------------------------------------------------------------
144 : 924 : int CubitFacet::is_flat()
145 : : {
146 [ + + ]: 924 : if (isFlat != 999)
147 : : {
148 : 132 : return isFlat;
149 : : }
150 : :
151 : : // check the edges first. If on a boundary then assume not flat for now
152 : : // This is to account for any in-plane curvature in the boundary edges
153 : :
154 : : int ii;
155 : 792 : CubitBoolean on_boundary = CUBIT_FALSE;
156 [ + + ][ + + ]: 1672 : for (ii=0; ii<3 && !on_boundary; ii++)
157 : : {
158 : 880 : CubitPoint *point_ptr = point(ii);
159 : 880 : CubitFacetEdge *edge_ptr = edge(ii);
160 : 880 : TDFacetBoundaryPoint *td = TDFacetBoundaryPoint::get_facet_boundary_point(point_ptr);
161 [ + + ][ + - ]: 880 : if (td != NULL || (edge_ptr && edge_ptr->num_adj_facets() <= 1))
[ - + ][ + + ]
162 : : {
163 : 748 : on_boundary = CUBIT_TRUE;
164 : : }
165 : : }
166 [ + + ]: 792 : if (on_boundary)
167 : : {
168 : 748 : isFlat = 0;
169 : : }
170 : : else
171 : : {
172 : :
173 : : CubitPoint *p0, *p1, *p2;
174 [ + - ]: 44 : points(p0, p1, p2);
175 : :
176 [ + - ]: 44 : CubitVector n0 = p0->normal( this );
177 [ + - ]: 44 : CubitVector n1 = p1->normal( this );
178 [ + - ]: 44 : CubitVector n2 = p2->normal( this );
179 [ + - ]: 44 : double dot0 = n0 % n1;
180 [ + - ]: 44 : double dot1 = n1 % n2;
181 [ + - ]: 44 : double dot2 = n2 % n0;
182 : :
183 : 44 : double tol = 1.0 - GEOMETRY_RESABS;
184 [ - + ][ # # ]: 44 : if (fabs(dot0) > tol && fabs(dot1) > tol && fabs(dot2) > tol)
[ # # ]
185 : 0 : isFlat = 1;
186 : : else
187 : 44 : isFlat = 0;
188 : : }
189 : :
190 : 924 : return isFlat;
191 : : }
192 : :
193 : : //===========================================================================
194 : : //Function Name: evaluate_position
195 : : //
196 : : //Member Type: PUBLIC
197 : : //Description: evaluate the facet at a position (projects to facet)
198 : : // eval_normal is NULL if normal not needed
199 : : //===========================================================================
200 : 0 : CubitStatus CubitFacet::evaluate_position( const CubitVector &start_position,
201 : : CubitVector *eval_point,
202 : : CubitVector *eval_normal)
203 : : {
204 : 0 : CubitStatus stat = CUBIT_SUCCESS;
205 : :
206 [ # # ]: 0 : if (is_flat())
207 : : {
208 [ # # ]: 0 : if (eval_point != NULL)
209 : 0 : closest_point(start_position, *eval_point);
210 [ # # ]: 0 : if (eval_normal != NULL)
211 [ # # ]: 0 : *eval_normal = normal();
212 : : }
213 : : else // project to the smooth facet
214 : : {
215 : : // project the position to the planar facet
216 [ # # ]: 0 : CubitVector close_point;
217 [ # # ]: 0 : stat = closest_point(start_position, close_point);
218 : :
219 : : // get the area coordinates of this point as a starting guess
220 [ # # ]: 0 : CubitVector area_coordinates;
221 [ # # ]: 0 : FacetEvalTool::facet_area_coordinate(this, close_point, area_coordinates);
222 : :
223 : : // now evaluate the smooth facet (this may alter the area coords)
224 [ # # ]: 0 : CubitVector proj_point;
225 : : CubitBoolean outside;
226 [ # # ]: 0 : double tol = sqrt(area()) * 1.0e-3;
227 : : stat = FacetEvalTool::project_to_facet( this, close_point, area_coordinates,
228 [ # # ]: 0 : proj_point, outside, tol );
229 : :
230 [ # # ]: 0 : if (eval_point != NULL)
231 : : {
232 [ # # ]: 0 : *eval_point = proj_point;
233 : : }
234 : : // compute the smooth normal if required
235 [ # # ]: 0 : if (eval_normal != NULL)
236 : : {
237 [ # # ]: 0 : FacetEvalTool::eval_facet_normal(this, area_coordinates, *eval_normal);
238 : : }
239 : : }
240 : :
241 : 0 : return stat;
242 : : }
243 : :
244 : : //===========================================================================
245 : : //Function Name: evaluate
246 : : //
247 : : //Member Type: PUBLIC
248 : : //Description: evaluate the facet at area coordinates
249 : : // eval_normal is NULL if normal not needed
250 : : //===========================================================================
251 : 0 : CubitStatus CubitFacet::evaluate( CubitVector &areacoord,
252 : : CubitVector *eval_point,
253 : : CubitVector *eval_normal )
254 : : {
255 [ # # ]: 0 : if (isBackwards)
256 : : {
257 : 0 : double temp = areacoord.y();
258 : 0 : areacoord.y( areacoord.z() );
259 : 0 : areacoord.z( temp );
260 : : }
261 : 0 : return FacetEvalTool::eval_facet( this, areacoord, eval_point, eval_normal );
262 : : }
263 : :
264 : : //===========================================================================
265 : : //Function Name: closest_point
266 : : //
267 : : //Member Type: PUBLIC
268 : : //Description: return the closest point on plane defined by the facet
269 : : //===========================================================================
270 : 0 : CubitStatus CubitFacet::closest_point(const CubitVector &point,
271 : : CubitVector &closest_point )
272 : : {
273 [ # # ][ # # ]: 0 : CubitPlane fac_plane = plane();
274 [ # # ][ # # ]: 0 : closest_point = fac_plane.project( point );
275 : 0 : return CUBIT_SUCCESS;
276 : : }
277 : :
278 : : //===========================================================================
279 : : //Function Name: closest_point_trimmed
280 : : //
281 : : //Member Type: PUBLIC
282 : : //Description: return the closest point on the facet to the point (trimmed
283 : : // to its boundary)
284 : : //===========================================================================
285 : 0 : CubitStatus CubitFacet::closest_point_trimmed( const CubitVector &mypoint,
286 : : CubitVector &closest_point,
287 : : CubitPoint *&next_edge_p1,
288 : : CubitPoint *&next_edge_p2)
289 : : {
290 [ # # ][ # # ]: 0 : CubitVector p1 = point(0)->coordinates();
291 [ # # ][ # # ]: 0 : CubitVector p2 = point(1)->coordinates();
292 [ # # ][ # # ]: 0 : CubitVector p3 = point(2)->coordinates();
293 : : //First get the edge vectors.
294 [ # # ]: 0 : p1.x();
295 [ # # ]: 0 : p2.x();
296 : :
297 [ # # ]: 0 : CubitVector y1 = p2 - p1;
298 [ # # ]: 0 : CubitVector y2 = p3 - p2;
299 [ # # ]: 0 : CubitVector y3 = p1 - p3;
300 : : //Now get the vectors from the point to the vertices of the facet.
301 [ # # ]: 0 : CubitVector w1 = mypoint - p1;
302 [ # # ]: 0 : CubitVector w2 = mypoint - p2;
303 [ # # ]: 0 : CubitVector w3 = mypoint - p3;
304 : : //Now cross the edge vectors with the vectors to the point. If the point
305 : : //in questionis inside the facet then each of these vectors will be in in
306 : : //the same direction as the normal of this facet.
307 [ # # ]: 0 : CubitVector x1 = y1*w1;
308 [ # # ]: 0 : CubitVector x2 = y2*w2;
309 [ # # ]: 0 : CubitVector x3 = y3*w3;
310 : :
311 : : //Now take the dot products to help us determine if it is in the triangle.
312 [ # # ]: 0 : CubitVector n = normal();
313 [ # # ]: 0 : double d1 = x1%n;
314 [ # # ]: 0 : double d2 = x2%n;
315 [ # # ]: 0 : double d3 = x3%n;
316 : :
317 : : //If this is true, then we just take the closest point to this facet.
318 [ # # ][ # # ]: 0 : if ( d1 >= -GEOMETRY_RESABS && d2 >= -GEOMETRY_RESABS && d3 >= -GEOMETRY_RESABS )
[ # # ]
319 : : {
320 [ # # ]: 0 : CubitStatus rv = this->closest_point(mypoint, closest_point);
321 [ # # ]: 0 : if ( rv != CUBIT_SUCCESS )
322 : : {
323 [ # # ][ # # ]: 0 : PRINT_ERROR("Closest Point Trimmed Error. Point in facet but can't"
[ # # ]
324 [ # # ]: 0 : " calc. point.\n");
325 : 0 : return CUBIT_FAILURE;
326 : : }
327 : 0 : next_edge_p1 = NULL;
328 : 0 : next_edge_p2 = NULL;
329 : 0 : return CUBIT_SUCCESS;
330 : : }
331 : 0 : CubitBoolean close_p1 = CUBIT_FALSE;
332 : 0 : CubitBoolean close_p2 = CUBIT_FALSE;
333 : 0 : CubitBoolean close_p3 = CUBIT_FALSE;
334 : : double k1, k2, k3;
335 : : //Now with the "d" values, determine which point or edge we
336 : : //should project to. We have to go through each of these and
337 : : //determine which is the closest.
338 [ # # ]: 0 : if ( d1 < 0.0 )
339 : : {
340 [ # # ]: 0 : double w1_dot_y1 = w1%y1;
341 [ # # ]: 0 : double y1_squared = y1.length_squared();
342 [ # # ][ # # ]: 0 : if ( y1_squared <= CUBIT_DBL_MIN && y1_squared >= -CUBIT_DBL_MIN )
343 : : {
344 [ # # ][ # # ]: 0 : PRINT_ERROR("Length of facet edge too small.\n");
[ # # ][ # # ]
345 : 0 : return CUBIT_FAILURE;
346 : : }
347 : 0 : k1 = w1_dot_y1/y1_squared;
348 [ # # ]: 0 : if ( k1 < 0.0 )
349 : 0 : close_p1 = CUBIT_TRUE;
350 [ # # ]: 0 : else if ( k1 > 1.0 )
351 : 0 : close_p2 = CUBIT_TRUE;
352 [ # # ][ # # ]: 0 : else if ( k1 >= 0.0 && k1 <= 1.0 + GEOMETRY_RESABS )
353 : : {
354 : : //So we know that y1, is the closest edge.
355 [ # # ][ # # ]: 0 : closest_point = p1 + k1*y1;
[ # # ]
356 [ # # ]: 0 : next_edge_p1 = point(0);
357 [ # # ]: 0 : next_edge_p2 = point(1);
358 : 0 : return CUBIT_SUCCESS;
359 : : }
360 : : }
361 [ # # ]: 0 : if ( d2 < 0.0 )
362 : : {
363 [ # # ]: 0 : double w2_dot_y2 = w2%y2;
364 [ # # ]: 0 : double y2_squared = y2.length_squared();
365 [ # # ][ # # ]: 0 : if ( y2_squared <= CUBIT_DBL_MIN && y2_squared >= -CUBIT_DBL_MIN )
366 : : {
367 [ # # ][ # # ]: 0 : PRINT_ERROR("Length of facet edge too small.\n");
[ # # ][ # # ]
368 : 0 : return CUBIT_FAILURE;
369 : : }
370 : 0 : k2 = w2_dot_y2/y2_squared;
371 [ # # ]: 0 : if ( k2 < 0.0 )
372 : : {
373 : 0 : close_p2 = CUBIT_TRUE;
374 : : }
375 [ # # ]: 0 : else if ( k2 > 1.0 )
376 : : {
377 : 0 : close_p3 = CUBIT_TRUE;
378 : : }
379 [ # # ][ # # ]: 0 : else if ( k2 >= 0.0 && k2 <= 1.0 + GEOMETRY_RESABS )
380 : : {
381 : : //So we know that y2, is the closest edge.
382 [ # # ][ # # ]: 0 : closest_point = p2 + k2*y2;
[ # # ]
383 [ # # ]: 0 : next_edge_p1 = point(1);
384 [ # # ]: 0 : next_edge_p2 = point(2);
385 : 0 : return CUBIT_SUCCESS;
386 : : }
387 : : }
388 [ # # ]: 0 : if ( d3 < 0.0 )
389 : : {
390 [ # # ]: 0 : double w3_dot_y3 = w3%y3;
391 [ # # ]: 0 : double y3_squared = y3.length_squared();
392 [ # # ][ # # ]: 0 : if ( y3_squared <= CUBIT_DBL_MIN && y3_squared >= -CUBIT_DBL_MIN )
393 : : {
394 [ # # ][ # # ]: 0 : PRINT_ERROR("Length of facet edge too small.\n");
[ # # ][ # # ]
395 : 0 : return CUBIT_FAILURE;
396 : : }
397 : 0 : k3 = w3_dot_y3/y3_squared;
398 [ # # ]: 0 : if ( k3 < 0.0 )
399 : : {
400 : 0 : close_p3 = CUBIT_TRUE;
401 : : }
402 [ # # ]: 0 : else if ( k3 > 1.0 )
403 : : {
404 : 0 : close_p1 = CUBIT_TRUE;
405 : : }
406 [ # # ][ # # ]: 0 : else if ( k3 >= 0.0 && k3 <= 1.0 + GEOMETRY_RESABS )
407 : : {
408 : : //So we know that y3, is the closest edge.
409 [ # # ][ # # ]: 0 : closest_point = p3 + k3*y3;
[ # # ]
410 [ # # ]: 0 : next_edge_p1 = point(2);
411 [ # # ]: 0 : next_edge_p2 = point(0);
412 : 0 : return CUBIT_SUCCESS;
413 : : }
414 : : }
415 : : //Now we have the distances, and which edges the point is closest to.
416 [ # # ][ # # ]: 0 : if ( close_p1 && !close_p2 && !close_p3 )
[ # # ]
417 : : {
418 [ # # ]: 0 : closest_point = p1 ;
419 [ # # ]: 0 : next_edge_p1 = point(0);
420 : 0 : next_edge_p2 = NULL;
421 : 0 : return CUBIT_SUCCESS;
422 : : }
423 [ # # ][ # # ]: 0 : else if ( close_p2 && !close_p1 && !close_p3 )
[ # # ]
424 : : {
425 [ # # ]: 0 : closest_point = p2;
426 [ # # ]: 0 : next_edge_p1 = point(1);
427 : 0 : next_edge_p2 = NULL;
428 : 0 : return CUBIT_SUCCESS;
429 : : }
430 [ # # ][ # # ]: 0 : else if ( close_p3 && !close_p1 && !close_p2 )
[ # # ]
431 : : {
432 [ # # ]: 0 : closest_point = p3;
433 [ # # ]: 0 : next_edge_p1 = point(2);
434 : 0 : next_edge_p2 = NULL;
435 : 0 : return CUBIT_SUCCESS;
436 : : }
437 [ # # ][ # # ]: 0 : else if( close_p1 && close_p2 && !close_p3 )
[ # # ]
438 : : {
439 [ # # ][ # # ]: 0 : if( w1.length_squared() < w2.length_squared() )
[ # # ]
440 [ # # ]: 0 : closest_point = p1;
441 : : else
442 [ # # ]: 0 : closest_point = p2;
443 : 0 : return CUBIT_SUCCESS;
444 : : }
445 [ # # ][ # # ]: 0 : else if( close_p2 && close_p3 && !close_p1 )
[ # # ]
446 : : {
447 [ # # ][ # # ]: 0 : if( w2.length_squared() < w3.length_squared() )
[ # # ]
448 [ # # ]: 0 : closest_point = p2;
449 : : else
450 [ # # ]: 0 : closest_point = p3;
451 : 0 : return CUBIT_SUCCESS;
452 : : }
453 [ # # ][ # # ]: 0 : else if( close_p1 && close_p3 && !close_p2 )
[ # # ]
454 : : {
455 [ # # ][ # # ]: 0 : if( w1.length_squared() < w3.length_squared() )
[ # # ]
456 [ # # ]: 0 : closest_point = p1;
457 : : else
458 [ # # ]: 0 : closest_point = p3;
459 : 0 : return CUBIT_SUCCESS;
460 : : }
461 : :
462 [ # # ][ # # ]: 0 : PRINT_ERROR("Problems finding the closest point to a facet.\n");
[ # # ][ # # ]
463 : 0 : return CUBIT_FAILURE;
464 : : }
465 : :
466 : : //===========================================================================
467 : : //Function Name: contains
468 : : //
469 : : //Member Type: PUBLIC
470 : : //Description: determine if point is contained in facet
471 : : //===========================================================================
472 : 2486 : CubitBoolean CubitFacet::contains( CubitPoint *p1 )
473 : : {
474 [ + + ]: 7964 : for ( int i = 0; i < 3; i++ )
475 [ + + ]: 6424 : if ( point(i) == p1 )
476 : 946 : return CUBIT_TRUE;
477 : 1540 : return CUBIT_FALSE;
478 : : }
479 : :
480 : : //===========================================================================
481 : : //Function Name: shared_facet
482 : : //
483 : : //Member Type: PUBLIC
484 : : //Description: Find the other facet that shares these two points.
485 : : // Note: assumes max two facets per edge
486 : : //===========================================================================
487 : 0 : CubitFacet* CubitFacet::shared_facet( CubitPoint *p1, CubitPoint *p2 )
488 : : {
489 : : //Find the other facet that shares these two points.
490 : : int ii;
491 [ # # ]: 0 : DLIList<CubitFacet*> facet_list;
492 [ # # ]: 0 : p1->facets(facet_list);
493 : :
494 [ # # ][ # # ]: 0 : for ( ii = facet_list.size(); ii > 0; ii-- )
495 : : {
496 [ # # ]: 0 : CubitFacet *t_facet = facet_list.get_and_step();
497 [ # # ]: 0 : if ( t_facet == this )
498 : 0 : continue;
499 [ # # ][ # # ]: 0 : if ( t_facet->contains(p2) )
500 : : {
501 [ # # ][ # # ]: 0 : assert( t_facet->contains(p1) );
502 : 0 : return t_facet;
503 : : }
504 : : }
505 [ # # ]: 0 : return (CubitFacet*)NULL;
506 : : }
507 : :
508 : : //===========================================================================
509 : : //Function Name: shared_facets
510 : : //
511 : : //Member Type: PUBLIC
512 : : //Description: Make a list of all facets adjacent this facet at the edge
513 : : // defined by p1, p2
514 : : //===========================================================================
515 : 396 : void CubitFacet::shared_facets( CubitPoint *p1, CubitPoint *p2,
516 : : DLIList <CubitFacet *> &adj_facet_list)
517 : : {
518 : : //Find the other facets that share these two points.
519 : : int ii;
520 [ + - ]: 396 : DLIList<CubitFacet*> facet_list;
521 [ + - ]: 396 : p1->facets(facet_list);
522 : :
523 [ + - ][ + + ]: 2244 : for ( ii = facet_list.size(); ii > 0; ii-- )
524 : : {
525 [ + - ]: 1848 : CubitFacet *t_facet = facet_list.get_and_step();
526 [ + + ]: 1848 : if ( t_facet == this )
527 : 396 : continue;
528 [ + - ][ + + ]: 1452 : if ( t_facet->contains(p2) )
529 : : {
530 [ + - ][ - + ]: 396 : assert( t_facet->contains(p1) );
531 [ + - ]: 1452 : adj_facet_list.append( t_facet );
532 : : }
533 [ + - ]: 396 : }
534 : 396 : }
535 : :
536 : :
537 : : //===========================================================================
538 : : //Function Name: shared_facet_on_surf
539 : : //
540 : : //Member Type: PUBLIC
541 : : //Description: Find the other facet that shares these two points and that
542 : : // is has the same tool id.
543 : : // Note: assumes max two facets per edge
544 : : //===========================================================================
545 : 132 : CubitFacet* CubitFacet::shared_facet_on_surf( CubitPoint *p1, CubitPoint *p2,
546 : : int tool_id )
547 : : {
548 : : //Find the other facet that shares these two points and that
549 : : // is marked with flag.
550 : : int ii;
551 [ + - ]: 132 : DLIList<CubitFacet*> facet_list;
552 [ + - ]: 132 : p1->facets(facet_list);
553 : :
554 [ + - ][ + + ]: 858 : for ( ii = facet_list.size(); ii > 0; ii-- )
555 : : {
556 [ + - ]: 748 : CubitFacet *t_facet = facet_list.get_and_step();
557 [ + + ]: 748 : if ( t_facet == this )
558 : 132 : continue;
559 [ + - ][ + + ]: 616 : if ( t_facet->contains(p2) )
560 : : {
561 [ + - ][ + + ]: 132 : if (tool_id == t_facet->tool_id())
562 : : {
563 [ + - ][ - + ]: 22 : assert( t_facet->contains(p1) );
564 : 22 : return t_facet;
565 : : }
566 : : }
567 : : }
568 [ + - ]: 132 : return (CubitFacet*)NULL;
569 : : }
570 : :
571 : : //===========================================================================
572 : : //Function Name: center
573 : : //
574 : : //Member Type: PUBLIC
575 : : //Description: return the centroid of this facet
576 : : //===========================================================================
577 : 0 : CubitVector CubitFacet::center()
578 : : {
579 : 0 : CubitVector vec( point(0)->coordinates() );
580 [ # # ]: 0 : vec += point(1)->coordinates();
581 [ # # ]: 0 : vec += point(2)->coordinates();
582 : 0 : vec /= 3.0;
583 : 0 : return vec;
584 : : }
585 : :
586 : : //===========================================================================
587 : : //Function Name: draw
588 : : //
589 : : //Member Type: PUBLIC
590 : : //Description: draw the facet
591 : : //===========================================================================
592 : 0 : void CubitFacet::debug_draw(int color, int flush_it, int draw_uv )
593 : : {
594 [ # # ]: 0 : if ( color == -1 )
595 : 0 : color = CUBIT_RED_INDEX;
596 : 0 : GfxDebug::draw_facet(this, color, draw_uv);
597 [ # # ]: 0 : if ( flush_it )
598 : 0 : GfxDebug::flush();
599 : 0 : }
600 : :
601 : : const int CubitFacet::point_edge_conn[30][2] = { {4, 8}, {8, 11}, {11, 13}, {13, 14},
602 : : {3, 7}, {7, 10}, {10, 12},
603 : : {2, 6}, {6, 9},
604 : : {1, 5},
605 : : {14, 12}, {12, 9}, {9, 5}, {5, 0},
606 : : {13, 10}, {10, 6}, {6, 1},
607 : : {11, 7}, {7, 2},
608 : : {8, 3},
609 : : {0, 1}, {1, 2}, {2, 3}, {3, 4},
610 : : {5, 6}, {6, 7}, {7, 8},
611 : : {9, 10}, {10, 11},
612 : : {12, 13} };
613 : : const int CubitFacet::point_facet_conn[16][3] = {
614 : : {0, 1, 5}, {1, 6, 5}, {1, 2, 6}, {2, 7, 6}, {2, 3, 7}, {3, 8, 7}, {3, 4, 8},
615 : : {5, 6, 9}, {6, 10, 9}, {6, 7, 10}, {7, 11, 10}, {7, 8, 11},
616 : : {9, 10, 12}, {10, 13, 12}, {10, 11, 13},
617 : : {12, 13, 14} };
618 : : const double CubitFacet::my_points[15][3] = {
619 : : {1.00, 0.00, 0.00},{0.75, 0.25, 0.00},{0.50, 0.50, 0.00},
620 : : {0.25, 0.75, 0.00},{0.00, 1.00, 0.00},{0.75, 0.00, 0.25},
621 : : {0.50, 0.25, 0.25},{0.25, 0.50, 0.25},{0.00, 0.75, 0.25},
622 : : {0.50, 0.00, 0.50},{0.25, 0.25, 0.50},{0.00, 0.50, 0.50},
623 : : {0.25, 0.00, 0.75},{0.00, 0.25, 0.75},{0.00, 0.00, 1.00} };
624 : :
625 : :
626 : :
627 : : //===========================================================================
628 : : //Function Name: min_diagonal
629 : : //
630 : : //Member Type: PUBLIC
631 : : //Description: from the three points find the minimum diagonal.
632 : : //===========================================================================
633 : 0 : double CubitFacet::min_diagonal()
634 : : {
635 : : //from the three points find the minimum diagonal.
636 [ # # ][ # # ]: 0 : CubitVector p1 = point(0)->coordinates();
637 [ # # ][ # # ]: 0 : CubitVector p2 = point(1)->coordinates();
638 [ # # ][ # # ]: 0 : CubitVector p3 = point(2)->coordinates();
639 [ # # ]: 0 : CubitVector temp1 = p2-p1;
640 [ # # ][ # # ]: 0 : CubitVector mid_side_1 = p1 + temp1/2.0;
641 [ # # ]: 0 : CubitVector temp2 = p3-p2;
642 [ # # ][ # # ]: 0 : CubitVector mid_side_2 = p2 + temp2/2.0;
643 [ # # ]: 0 : CubitVector temp3 = p1-p3;
644 [ # # ][ # # ]: 0 : CubitVector mid_side_3 = p3 + temp3/2.0;
645 [ # # ]: 0 : CubitVector diagv_1 = p3 - mid_side_1;
646 [ # # ]: 0 : CubitVector diagv_2 = p2 - mid_side_3;
647 [ # # ]: 0 : CubitVector diagv_3 = p1 - mid_side_2;
648 : :
649 [ # # ]: 0 : double diag_1 = diagv_1.length();
650 [ # # ]: 0 : double diag_2 = diagv_2.length();
651 [ # # ]: 0 : double diag_3 = diagv_3.length();
652 [ # # ][ # # ]: 0 : if ( diag_1 >= diag_2 && diag_1 >= diag_3 )
653 : 0 : return diag_1;
654 [ # # ][ # # ]: 0 : else if ( diag_2 > diag_1 && diag_2 > diag_3 )
655 : 0 : return diag_2;
656 : 0 : return diag_3;
657 : : }
658 : :
659 : : //-------------------------------------------------------------------------
660 : : // Purpose : Get the two points of the edge opposite the
661 : : // passed point.
662 : : //
663 : : // Special Notes :
664 : : //
665 : : // Creator : Jason Kraftcheck
666 : : //
667 : : // Creation Date : 03/07/00
668 : : //-------------------------------------------------------------------------
669 : 0 : void CubitFacet::opposite_edge( CubitPoint* thepoint,
670 : : CubitPoint*& p1, CubitPoint *& p2 )
671 : : {
672 : 0 : p1 = p2 = 0;
673 : : int i;
674 : :
675 [ # # ]: 0 : for( i = 0; i < 3; i++ )
676 [ # # ]: 0 : if( point(i) == thepoint )
677 : 0 : break;
678 : :
679 [ # # ]: 0 : if( i < 3 ) //point is on this triangle
680 : : {
681 : 0 : p1 = point((i+1)%3);
682 : 0 : p2 = point((i+2)%3);
683 : : }
684 : 0 : }
685 : :
686 : : //-------------------------------------------------------------------------
687 : : // Purpose : Local modification functions.
688 : : //
689 : : // Special Notes :
690 : : //
691 : : // Creator :
692 : : //
693 : : // Creation Date : 03/25/00
694 : : //-------------------------------------------------------------------------
695 : 0 : CubitPoint* CubitFacet::split_edge( CubitPoint* /* edge_pt1 */,
696 : : CubitPoint* /* edge_pt2 */,
697 : : const CubitVector& /* position */ )
698 : : {
699 : : // this function should be implemented in child class since it creates
700 : : // new facet and point entities
701 : : assert(1);
702 : 0 : CubitPoint* new_point = NULL;
703 : 0 : return new_point;
704 : : }
705 : :
706 : : //-------------------------------------------------------------------------
707 : : // Purpose : insert a point into the facet
708 : : //
709 : : // Special Notes : create two new facets and return them
710 : : //
711 : : // Creator :
712 : : //
713 : : // Creation Date :
714 : : //-------------------------------------------------------------------------
715 : 0 : CubitPoint* CubitFacet::insert_point( const CubitVector& /* position */,
716 : : CubitFacet*& /* new_tri1 */,
717 : : CubitFacet*& /* new_tri2 */)
718 : : {
719 : : // this function should be implemented in child class since it creates
720 : : // new facet and point entities
721 : : assert(1);
722 : 0 : CubitPoint* new_point = NULL;
723 : 0 : return new_point;
724 : : }
725 : :
726 : :
727 : : //-------------------------------------------------------------------------
728 : : // Purpose : return the index of the other index not used by
729 : : // points p1 and p2
730 : : //
731 : : // Special Notes :
732 : : //
733 : : // Creator :
734 : : //
735 : : // Creation Date :
736 : : //-------------------------------------------------------------------------
737 : 0 : int CubitFacet::other_index( CubitPoint* pt1, CubitPoint* pt2 )
738 : : {
739 : : int i;
740 [ # # ][ # # ]: 0 : for( i = 0; (point(i) != pt1) && (i < 3); i++ );
[ # # ]
741 [ # # ]: 0 : if( i == 3 ) return -1;
742 : :
743 : 0 : int i1 = (i+1)%3;
744 : 0 : int i2 = (i+2)%3;
745 : :
746 [ # # ]: 0 : if( point(i1) == pt2 ) return i2;
747 [ # # ]: 0 : else if( point(i2) == pt2 ) return i1;
748 : 0 : else return -1;
749 : : }
750 : :
751 : : //-------------------------------------------------------------------------
752 : : // Purpose : compute the angle at one of the points on the facet
753 : : //
754 : : // Special Notes : angle is in radians
755 : : //
756 : : // Creator : Steve Owen
757 : : //
758 : : // Creation Date : 06/28/00
759 : : //-------------------------------------------------------------------------
760 : 3564 : double CubitFacet::angle( CubitPoint *pt )
761 : : {
762 : 3564 : int ii, index = -1;
763 : 3564 : CubitBoolean found=CUBIT_FALSE;
764 [ + + ][ + + ]: 9900 : for (ii=0; ii<3 && !found; ii++) {
765 [ + - ][ + + ]: 6336 : if (pt==this->point(ii)) {
766 : 3564 : index = ii;
767 : 3564 : found = CUBIT_TRUE;
768 : : }
769 : : }
770 [ - + ]: 3564 : if (!found) {
771 : 0 : return 0.0e0;
772 : : }
773 [ + - ][ + - ]: 3564 : CubitVector e0 = this->point((index+1)%3)->coordinates() - pt->coordinates();
[ + - ][ + - ]
774 [ + - ][ + - ]: 3564 : CubitVector e1 = this->point((index+2)%3)->coordinates() - pt->coordinates();
[ + - ][ + - ]
775 [ + - ]: 3564 : e0.normalize();
776 [ + - ]: 3564 : e1.normalize();
777 : : double myangle;
778 [ + - ]: 3564 : double cosangle = e0%e1;
779 [ - + ]: 3564 : if (cosangle >= 1.0) {
780 : 0 : myangle = 0.0e0;
781 : : }
782 [ - + ]: 3564 : else if( cosangle <= -1.0 ) {
783 : 0 : myangle = CUBIT_PI;
784 : : }
785 : : else {
786 : 3564 : myangle = acos(cosangle);
787 : : }
788 : 3564 : return myangle;
789 : : }
790 : :
791 : : //-------------------------------------------------------------------------
792 : : // Purpose : return the edge on a facet and its corresponding index.
793 : : //
794 : : // Special Notes : The index of the edge corresponds to the point index on
795 : : // the facet opposite the edge.
796 : : //
797 : : // Creator : Steve Owen
798 : : //
799 : : // Creation Date : 06/28/00
800 : : //-------------------------------------------------------------------------
801 : 0 : CubitFacetEdge *CubitFacet::edge_from_pts( CubitPoint *p1,
802 : : CubitPoint *p2,
803 : : int &edge_index )
804 : : {
805 [ # # ]: 0 : if ((point(0) == p1 && point(1) == p2) ||
[ # # # # ]
[ # # ]
806 [ # # ]: 0 : (point(0) == p2 && point(1) == p1)) {
807 : 0 : edge_index = 2;
808 : 0 : return edge(2);
809 : : }
810 [ # # ]: 0 : if ((point(1) == p1 && point(2) == p2) ||
[ # # # # ]
[ # # ]
811 [ # # ]: 0 : (point(1) == p2 && point(2) == p1)) {
812 : 0 : edge_index = 0;
813 : 0 : return edge(0);
814 : : }
815 [ # # ]: 0 : if ((point(2) == p1 && point(0) == p2) ||
[ # # # # ]
[ # # ]
816 [ # # ]: 0 : (point(2) == p2 && point(0) == p1)) {
817 : 0 : edge_index = 1;
818 : 0 : return edge(1);
819 : : }
820 : 0 : edge_index = -1;
821 : 0 : return NULL;
822 : : }
823 : :
824 : : //-------------------------------------------------------------------------
825 : : // Purpose : return the index of an edge on a facet given two
826 : : // points on the facet
827 : : //
828 : : // Special Notes : The index of the edge corresponds to the point index on
829 : : // the facet opposite the edge.
830 : : //
831 : : // Creator : Steve Owen
832 : : //
833 : : // Creation Date : 06/28/00
834 : : //-------------------------------------------------------------------------
835 : 3168 : int CubitFacet::edge_index( CubitPoint *p1, CubitPoint *p2, int &sense )
836 : : {
837 [ + + ][ + + ]: 3168 : if (point(0) == p1 && point(1) == p2)
[ + + ]
838 : : {
839 : 572 : sense = 1;
840 : 572 : return 2;
841 : : }
842 [ + + ][ + + ]: 2596 : if (point(0) == p2 && point(1) == p1)
[ + + ]
843 : : {
844 : 484 : sense = -1;
845 : 484 : return 2;
846 : : }
847 [ + + ][ + - ]: 2112 : if (point(1) == p1 && point(2) == p2)
[ + + ]
848 : : {
849 : 726 : sense = 1;
850 : 726 : return 0;
851 : : }
852 [ + + ][ + - ]: 1386 : if (point(1) == p2 && point(2) == p1)
[ + + ]
853 : : {
854 : 330 : sense = -1;
855 : 330 : return 0;
856 : : }
857 [ + + ][ + - ]: 1056 : if (point(2) == p1 && point(0) == p2)
[ + + ]
858 : : {
859 : 396 : sense = 1;
860 : 396 : return 1;
861 : : }
862 [ + - ][ + - ]: 660 : if (point(2) == p2 && point(0) == p1)
[ + - ]
863 : : {
864 : 660 : sense = -1;
865 : 660 : return 1;
866 : : }
867 : 0 : sense = 0;
868 : 0 : return -1;
869 : : }
870 : :
871 : : // overloaded function - based on edge pointer
872 : 44 : int CubitFacet::edge_index( CubitFacetEdge *theedge )
873 : : {
874 [ + + ]: 44 : if (edge(0) == theedge) return 0;
875 [ + + ]: 22 : if (edge(1) == theedge) return 1;
876 [ + - ]: 11 : if (edge(2) == theedge) return 2;
877 : 0 : return -1;
878 : : }
879 : :
880 : :
881 : : //-------------------------------------------------------------------------
882 : : // Purpose : return the index of a point on a facet
883 : : //
884 : : // Creator : Steve Owen
885 : : //
886 : : // Creation Date : 06/28/00
887 : : //-------------------------------------------------------------------------
888 : 121 : int CubitFacet::point_index( CubitPoint *pt )
889 : : {
890 [ + + ]: 121 : if (point(0) == pt) return 0;
891 [ + + ]: 77 : if (point(1) == pt) return 1;
892 [ + - ]: 66 : if (point(2) == pt) return 2;
893 : 0 : return -1;
894 : : }
895 : :
896 : : //-------------------------------------------------------------------------
897 : : // Purpose : get the points that define an edge of a facet.
898 : : //
899 : : // Special Notes : The index of the edge corresponds to the point index on
900 : : // the facet opposite the edge. The points are returned
901 : : // based on the edgeUse order
902 : : //
903 : : // Creator : Steve Owen
904 : : //
905 : : // Creation Date : 06/28/00
906 : : //-------------------------------------------------------------------------
907 : 1694 : void CubitFacet::get_edge_pts( int index, CubitPoint *&p1, CubitPoint *&p2 )
908 : : {
909 [ + - ][ - + ]: 1694 : assert(index >=0 && index <=2);
910 : :
911 : 1694 : p1 = point((index+1)%3);
912 : 1694 : p2 = point((index+2)%3);
913 : 1694 : }
914 : :
915 : : //-------------------------------------------------------------------------
916 : : // Purpose : update the bounding box of the facet based on the control
917 : : // polygon of the facet.
918 : : //
919 : : // Creator : Steve Owen
920 : : //
921 : : // Creation Date : 06/28/00
922 : : //-------------------------------------------------------------------------
923 : 264 : void CubitFacet::update_bezier_bounding_box( )
924 : : {
925 : : int i,j;
926 [ + - ][ + + ]: 1584 : CubitVector ctrl_pts[5], min, max;
[ + - ][ + - ]
927 : : CubitFacetEdge *myedge;
928 [ + - ][ + - ]: 264 : CubitBox bbox = bounding_box();
929 [ + - ][ + - ]: 264 : min = bbox.minimum();
930 [ + - ][ + - ]: 264 : max = bbox.maximum();
931 [ + + ]: 1056 : for (i=0; i<3; i++) {
932 [ + - ]: 792 : myedge = edge(i);
933 [ - + ]: 792 : assert(myedge != 0);
934 [ + - ]: 792 : myedge->control_points( this, ctrl_pts );
935 [ + + ]: 3168 : for (j=1; j<4; j++) {
936 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].x() < min.x()) min.x( ctrl_pts[j].x() );
[ - + ][ # # ]
[ # # ]
937 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].y() < min.y()) min.y( ctrl_pts[j].y() );
[ - + ][ # # ]
[ # # ]
938 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].z() < min.z()) min.z( ctrl_pts[j].z() );
[ - + ][ # # ]
[ # # ]
939 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].x() > max.x()) max.x( ctrl_pts[j].x() );
[ - + ][ # # ]
[ # # ]
940 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].y() > max.y()) max.y( ctrl_pts[j].y() );
[ - + ][ # # ]
[ # # ]
941 [ + - ][ + - ]: 2376 : if (ctrl_pts[j].z() > max.z()) max.z( ctrl_pts[j].z() );
[ - + ][ # # ]
[ # # ]
942 : : }
943 : : }
944 [ + - ][ + + ]: 1848 : CubitVector patch_ctrl_pts[6];
945 [ + - ]: 264 : get_control_points( patch_ctrl_pts );
946 : : assert(patch_ctrl_pts != 0);
947 [ + + ]: 1848 : for (j=0; j<6; j++) {
948 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].x() < min.x()) min.x( patch_ctrl_pts[j].x() );
[ - + ][ # # ]
[ # # ]
949 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].y() < min.y()) min.y( patch_ctrl_pts[j].y() );
[ - + ][ # # ]
[ # # ]
950 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].z() < min.z()) min.z( patch_ctrl_pts[j].z() );
[ - + ][ # # ]
[ # # ]
951 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].x() > max.x()) max.x( patch_ctrl_pts[j].x() );
[ - + ][ # # ]
[ # # ]
952 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].y() > max.y()) max.y( patch_ctrl_pts[j].y() );
[ - + ][ # # ]
[ # # ]
953 [ + - ][ + - ]: 1584 : if (patch_ctrl_pts[j].z() > max.z()) max.z( patch_ctrl_pts[j].z() );
[ - + ][ # # ]
[ # # ]
954 : : }
955 [ + - ][ + - ]: 264 : bBox.reset(min,max);
956 : 264 : }
957 : :
958 : : //-------------------------------------------------------------------------
959 : : // Purpose : reset the bounding box of the facet based on the control
960 : : // polygon of the facet.
961 : : //
962 : : // Creator : Steve Owen
963 : : //
964 : : // Creation Date : 3/14/02
965 : : //-------------------------------------------------------------------------
966 : 264 : void CubitFacet::reset_bounding_box( )
967 : : {
968 : :
969 : 264 : CubitPoint *p1 = point (0);
970 : 264 : CubitPoint *p2 = point (1);
971 : 264 : CubitPoint *p3 = point (2);
972 : :
973 : : // define the bounding box
974 : :
975 [ - + ][ # # ]: 264 : if (!patchCtrlPts || is_flat())
[ + - ]
976 : : {
977 [ + - ][ + - ]: 264 : CubitVector bbox_min, bbox_max;
978 [ + - ][ + - ]: 264 : bbox_min.x(min3(p1->x(),p2->x(),p3->x()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
979 [ + - ][ + - ]: 264 : bbox_min.y(min3(p1->y(),p2->y(),p3->y()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
980 [ + - ][ + - ]: 264 : bbox_min.z(min3(p1->z(),p2->z(),p3->z()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
981 [ + - ][ + - ]: 264 : bbox_max.x(max3(p1->x(),p2->x(),p3->x()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
982 [ + - ][ + - ]: 264 : bbox_max.y(max3(p1->y(),p2->y(),p3->y()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
983 [ + - ][ + - ]: 264 : bbox_max.z(max3(p1->z(),p2->z(),p3->z()));
[ + + ][ + - ]
[ + - ][ + - ]
[ + + ][ + - ]
[ + - ][ + + ]
[ + - ][ + - ]
[ + - ][ + - ]
984 [ + - ]: 264 : bBox.reset(bbox_min,bbox_max);
985 : : }
986 : : else
987 : : {
988 : 0 : update_bezier_bounding_box( );
989 : : }
990 : 264 : }
991 : :
992 : : //-------------------------------------------------------------------------
993 : : // Purpose : get the next CCW (rhr) facetedge on the facet
994 : : //
995 : : // Special Notes :
996 : : //
997 : : // Creator : Steve Owen
998 : : //
999 : : // Creation Date : 06/28/00
1000 : : //-------------------------------------------------------------------------
1001 : 2244 : CubitFacetEdge *CubitFacet::next_edge( CubitFacetEdge *theedge )
1002 : : {
1003 : : int i;
1004 [ + - ]: 5280 : for (i=0; i<3; i++) {
1005 [ + + ]: 5280 : if (theedge == edge(i)) {
1006 : 2244 : return edge((i+1)%3);
1007 : : }
1008 : : }
1009 : 0 : return NULL;
1010 : : }
1011 : :
1012 : : //-------------------------------------------------------------------------
1013 : : // Purpose : get the prev CCW (rhr) facetedge on the facet
1014 : : //
1015 : : // Special Notes :
1016 : : //
1017 : : // Creator : Steve Owen
1018 : : //
1019 : : // Creation Date : 06/28/00
1020 : : //-------------------------------------------------------------------------
1021 : 0 : CubitFacetEdge *CubitFacet::prev_edge( CubitFacetEdge *theedge )
1022 : : {
1023 : : int i;
1024 [ # # ]: 0 : for (i=0; i<3; i++) {
1025 [ # # ]: 0 : if (theedge == edge(i)) {
1026 : 0 : return edge((i+2)%3);
1027 : : }
1028 : : }
1029 : 0 : return NULL;
1030 : : }
1031 : :
1032 : : //-------------------------------------------------------------------------
1033 : : // Purpose : reorient the facet
1034 : : //
1035 : : // Special Notes :
1036 : : //
1037 : : // Creator : Steve Owen
1038 : : //
1039 : : // Creation Date : 06/28/00
1040 : : //-------------------------------------------------------------------------
1041 : 0 : void CubitFacet::flip()
1042 : : {
1043 : : // must be implemented in child class
1044 : 0 : assert(0);
1045 : :
1046 : : }
1047 : :
1048 : : //-------------------------------------------------------------------------
1049 : : // Purpose : return the next edge on the triangle at the point
1050 : : //
1051 : : // Special Notes :
1052 : : //
1053 : : // Creator : Steve Owen
1054 : : //
1055 : : // Creation Date : 4/31/01
1056 : : //-------------------------------------------------------------------------
1057 : 0 : CubitFacetEdge *CubitFacet::next_edge_at_point( CubitFacetEdge *edge_ptr,
1058 : : CubitPoint *point_ptr )
1059 : : {
1060 : 0 : int eidx = edge_index( edge_ptr );
1061 : 0 : int pidx = point_index( point_ptr );
1062 [ # # # # ]: 0 : switch(eidx)
1063 : : {
1064 : : case 0:
1065 [ # # # ]: 0 : switch(pidx)
1066 : : {
1067 : 0 : case 1: eidx = 2; break;
1068 : 0 : case 2: eidx = 1; break;
1069 : 0 : default: eidx = -1; break;
1070 : : }
1071 : 0 : break;
1072 : : case 1:
1073 [ # # # ]: 0 : switch(pidx)
1074 : : {
1075 : 0 : case 0: eidx = 2; break;
1076 : 0 : case 2: eidx = 0; break;
1077 : 0 : default: eidx = -1; break;
1078 : : }
1079 : 0 : break;
1080 : : case 2:
1081 [ # # # ]: 0 : switch(pidx)
1082 : : {
1083 : 0 : case 0: eidx = 1; break;
1084 : 0 : case 1: eidx = 0; break;
1085 : 0 : default: eidx = -1; break;
1086 : : }
1087 : 0 : break;
1088 : : default:
1089 : 0 : eidx = -1;
1090 : 0 : break;
1091 : : }
1092 [ # # ]: 0 : if (eidx == -1)
1093 : 0 : return (CubitFacetEdge *)NULL;
1094 : :
1095 : 0 : return edge( eidx );
1096 : : }
1097 : :
1098 : : //======================================================================
1099 : : // Function: get_edge control_points (PUBLIC)
1100 : : // Description: return the control points on the facet edges based
1101 : : // on its edge use directions
1102 : : // Author: sjowen
1103 : : // Date: 5/01
1104 : : //======================================================================
1105 : 264 : CubitStatus CubitFacet::get_edge_control_points( CubitVector P[3][5] )
1106 : : {
1107 : : CubitStatus stat;
1108 : : CubitFacetEdge *edge_ptr;
1109 [ + - ][ + + ]: 1584 : CubitVector ctrl_pts[5];
1110 : : int ii, jj;
1111 [ + + ]: 1056 : for (ii=0; ii<3; ii++) {
1112 [ + - ]: 792 : edge_ptr = edge(ii);
1113 [ + - ]: 792 : stat = edge_ptr->control_points(this, ctrl_pts);
1114 [ - + ]: 792 : if (stat!= CUBIT_SUCCESS)
1115 : 0 : return stat;
1116 [ + + ]: 4752 : for (jj=0; jj<5; jj++)
1117 : : {
1118 [ + - ]: 3960 : P[ii][jj] = ctrl_pts[jj];
1119 : : }
1120 : : }
1121 : 264 : return stat;
1122 : : }
1123 : :
1124 : : //======================================================================
1125 : : // Function: area (PUBLIC)
1126 : : // Description: compute the area of the facet
1127 : : // Author: sjowen
1128 : : // Date: 3/02
1129 : : //======================================================================
1130 : 8316 : double CubitFacet::area( )
1131 : : {
1132 [ + - ][ + - ]: 8316 : CubitVector e0(point(0)->coordinates(), point(1)->coordinates());
[ + - ][ + - ]
[ + - ]
1133 [ + - ][ + - ]: 8316 : CubitVector e1(point(0)->coordinates(), point(2)->coordinates());
[ + - ][ + - ]
[ + - ]
1134 [ + - ][ + - ]: 8316 : double area = (e0 * e1).length() * 0.5;
1135 : 8316 : return area;
1136 : : }
1137 : :
1138 : 0 : double CubitFacet::aspect_ratio()
1139 : : {
1140 : : static const double normal_coeff = sqrt( 3. ) / 6.;
1141 : :
1142 : : // three vectors for each side
1143 [ # # ][ # # ]: 0 : CubitVector a = point(1)->coordinates() - point(0)->coordinates();
[ # # ][ # # ]
[ # # ]
1144 [ # # ][ # # ]: 0 : CubitVector b = point(2)->coordinates() - point(1)->coordinates();
[ # # ][ # # ]
[ # # ]
1145 [ # # ][ # # ]: 0 : CubitVector c = point(0)->coordinates() - point(2)->coordinates();
[ # # ][ # # ]
[ # # ]
1146 : :
1147 [ # # ]: 0 : double a1 = a.length();
1148 [ # # ]: 0 : double b1 = b.length();
1149 [ # # ]: 0 : double c1 = c.length();
1150 : :
1151 [ # # ]: 0 : double hm = a1 > b1 ? a1 : b1;
1152 [ # # ]: 0 : hm = hm > c1 ? hm : c1;
1153 : :
1154 [ # # ]: 0 : CubitVector ab = a * b;
1155 [ # # ]: 0 : double denominator = ab.length();
1156 : :
1157 [ # # ]: 0 : if( denominator < CUBIT_DBL_MIN )
1158 : 0 : return (double)CUBIT_DBL_MAX;
1159 : : else
1160 : : {
1161 : : double aspect_ratio;
1162 : 0 : aspect_ratio = normal_coeff * hm * (a1 + b1 + c1) / denominator;
1163 : :
1164 [ # # ]: 0 : if( aspect_ratio > 0 )
1165 [ # # ]: 0 : return (double) CUBIT_MIN( aspect_ratio, CUBIT_DBL_MAX );
1166 [ # # ]: 0 : return (double) CUBIT_MAX( aspect_ratio, -CUBIT_DBL_MAX );
1167 : : }
1168 : : }
1169 : :
1170 : :
1171 : : //======================================================================
1172 : : // Function: init_patch (PUBLIC)
1173 : : // Description: computes the interior control points for the facet and
1174 : : // stores them with the facet. Assumes edge control points
1175 : : // have already been computed
1176 : : // Author: sjowen
1177 : : // Date: 10/03
1178 : : //======================================================================
1179 : 0 : CubitStatus CubitFacet::init_patch( )
1180 : : {
1181 : 0 : return FacetEvalTool::init_bezier_facet( this );
1182 : : }
1183 : :
1184 : : //===========================================================================
1185 : : //Function Name: add_edge
1186 : : //
1187 : : //Member Type: PUBLIC
1188 : : //Description: add an existing edge to a facet. The edge must contain
1189 : : // points that are already part of this facet
1190 : : //===========================================================================
1191 : 0 : void CubitFacet::add_edge(CubitFacetEdge *edge)
1192 : : {
1193 : 0 : CubitPoint *p0 = edge->point(0);
1194 : 0 : CubitPoint *p1 = edge->point(1);
1195 : :
1196 : : // find the points on this facet
1197 : :
1198 : 0 : int p0_index = point_index( p0 );
1199 : 0 : int p1_index = point_index( p1 );
1200 [ # # ][ # # ]: 0 : assert(p0_index >=0 && p1_index >= 0);
1201 [ # # ]: 0 : assert(p0_index != p1_index);
1202 : :
1203 : : // add the edge based on the relative index of the points
1204 : : // set the edge use based on their order
1205 : :
1206 : : int edge_index;
1207 [ # # ]: 0 : if ((p0_index+1)%3 == p1_index)
1208 : : {
1209 : 0 : edge_index = (p1_index + 1)%3;
1210 [ # # ]: 0 : assert(this->edge(edge_index) == NULL);
1211 : 0 : this->edge( edge, edge_index );
1212 : 0 : this->edge_use(1, edge_index);
1213 : : }
1214 [ # # ]: 0 : else if((p0_index+2)%3 == p1_index)
1215 : : {
1216 : 0 : edge_index = (p0_index + 1)%3;
1217 [ # # ]: 0 : assert(this->edge(edge_index) == NULL);
1218 : 0 : this->edge( edge, edge_index );
1219 : 0 : this->edge_use(-1, edge_index);
1220 : : }
1221 : : else
1222 : : {
1223 : 0 : assert(0); // shouldn't happen
1224 : : }
1225 : :
1226 : : // add the facet to the edge
1227 : :
1228 : 0 : edge->add_facet( this );
1229 : 0 : }
1230 : :
1231 : :
1232 : 0 : CubitPoint *CubitFacet::opposite_point( CubitFacetEdge *edge )
1233 : : {
1234 : : int i;
1235 [ # # ]: 0 : for( i = 0; i < 3; i++ )
1236 [ # # ][ # # ]: 0 : if( point(i) != edge->point(0) && point(i) != edge->point(1) )
[ # # ]
1237 : 0 : return point(i);
1238 : :
1239 : 0 : return NULL;
1240 : : }
1241 : :
1242 : 0 : CubitVector CubitFacet::update_normal( void )
1243 : : {
1244 : 0 : this->update_plane();
1245 : 0 : return normal();
1246 : : }
1247 : :
1248 : 0 : void CubitFacet::unlink_from_children( void )
1249 : : {
1250 : :
1251 : 0 : this->point(0)->remove_facet(this);
1252 : 0 : this->point(1)->remove_facet(this);
1253 : 0 : this->point(2)->remove_facet(this);
1254 : 0 : this->edge(0)->remove_facet(this);
1255 : 0 : this->edge(1)->remove_facet(this);
1256 : 0 : this->edge(2)->remove_facet(this);
1257 : :
1258 : 0 : }
1259 : :
1260 : 0 : CubitFacetEdge* CubitFacet::shared_edge( CubitFacet *cubit_facet )
1261 : : {
1262 [ # # ]: 0 : for( int i=0; i<3; i++ )
1263 [ # # ]: 0 : for( int j=0; j<3; j++ )
1264 [ # # ]: 0 : if( edge(i) == cubit_facet->edge(j) )
1265 : 0 : return edge(i);
1266 : :
1267 : 0 : return NULL;
1268 [ + - ][ + - ]: 6540 : }
|
