1 /* *****************************************************************
2     MESQUITE -- The Mesh Quality Improvement Toolkit
3 
4     Copyright 2004 Sandia Corporation and Argonne National
5     Laboratory.  Under the terms of Contract DE-AC04-94AL85000
6     with Sandia Corporation, the U.S. Government retains certain
7     rights in this software.
8 
9     This library is free software; you can redistribute it and/or
10     modify it under the terms of the GNU Lesser General Public
11     License as published by the Free Software Foundation; either
12     version 2.1 of the License, or (at your option) any later version.
13 
14     This library is distributed in the hope that it will be useful,
15     but WITHOUT ANY WARRANTY; without even the implied warranty of
16     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17     Lesser General Public License for more details.
18 
19     You should have received a copy of the GNU Lesser General Public License
20     (lgpl.txt) along with this library; if not, write to the Free Software
21     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
22 
23     diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov,
24     pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov
25 
26   ***************************************************************** */
27 #ifndef MESQUITE_PATCHDATA_HPP
28 #define MESQUITE_PATCHDATA_HPP
29 /*!
30   \file   PatchData.hpp
31   \brief    This file contains the PatchData class and its associated mementos.
32 
33 
34   The PatchData class provides the mesh information and functionality to Mesquite.
35   The PatchDataVerticesMemento class allows the state of a PatchData object to be saved
36   in order to later restore that object to its previous state.
37 
38   \author Thomas Leurent
39   \author Michael Brewer
40   \date   2002-01-17
41 */
42 
43 #include "Mesquite.hpp"
44 #include "MsqVertex.hpp"
45 #include "MsqMeshEntity.hpp"
46 #include "MsqVertex.hpp"
47 #include "MeshInterface.hpp"
48 #include "MsqError.hpp"
49 #include "Settings.hpp"
50 #include "NodeSet.hpp"
51 #include "MappingFunction.hpp"
52 
53 #include <cstddef>
54 #include <cstdlib>
55 #include <map>
56 #include <vector>
57 #include <iosfwd>
58 
59 namespace MBMesquite
60 {
61   class ExtraData;
62   class PatchDataVerticesMemento;
63   class Mesh;
64   class Settings;
65 
66   /*!
67     Contains all the mesh information necessary for
68     one iteration of the optimization algorithms over a
69     local mesh patch. */
70   class PatchData
71   {
72   public:
73       // Constructor/Destructor
74     MESQUITE_EXPORT PatchData();
75     MESQUITE_EXPORT ~PatchData();
76 
attach_settings(const Settings * p_settings)77     MESQUITE_EXPORT void attach_settings( const Settings* p_settings )
78       { mSettings = p_settings; }
settings() const79     MESQUITE_EXPORT const Settings* settings() const
80       { return mSettings; }
81 
82       /**\brief For use by testing code -- create patch explicitly
83        *
84        * Create a patch containing elements of the same type and
85        * without any higher-order nodes.
86        *
87        *\param num_vertex   Number of vertices in patch
88        *\param vtx_coords   Array of vertex coords.  Length must be 3*num_vertex
89        *\param type         Element type
90        *\param connectivity Element connectivity, specified as a list
91        *                    of vertex numbers, beginning with zero.
92        *\param vertex_fixed_flags Optional array to specify which vertices
93        *                    are to be marked as fixed.  If not specified,
94        *                    no vertices are fixed.
95        */
96 	MESQUITE_EXPORT
97     void fill( size_t num_vertex, const double* vtx_coords,
98                size_t num_elem, EntityTopology type,
99                const size_t* connectivity,
100                const bool* vertex_fixed_flags,
101                MsqError& err );
102 
103       /**\brief For use by testing code -- create patch explicitly
104        *
105        * Create a patch containing elements without any higher-order nodes.
106        *
107        *\param num_vertex   Number of vertices in patch
108        *\param vtx_coords   Array of vertex coords.  Length must be 3*num_vertex
109        *\param elem_types   The type of each element
110        *\param connectivity Element connectivity, specified as a list
111        *                    of vertex numbers, beginning with zero.
112        *\param vertex_fixed_flags Optional array to specify which vertices
113        *                    are to be marked as fixed.  If NULL,
114        *                    no vertices are fixed.
115        */
116 	MESQUITE_EXPORT
117     void fill( size_t num_vertex, const double* vtx_coords,
118                size_t num_elem, const EntityTopology* elem_types,
119                const size_t* connectivity,
120                const bool* vertex_fixed_flags,
121                MsqError& err );
122 
123       /**\brief For use by testing code -- create patch explicitly
124        *
125        * Most general form of fill function.  Works for polygons,
126        * elements with higher-order nodes, etc.
127        *
128        *\param num_vertex   Number of vertices in patch
129        *\param vtx_coords   Array of vertex coords.  Length must be 3*num_vertex
130        *\param elem_types   The type of each element
131        *\param vertex_per_elem The length of the connectivity list for each element.
132        *\param connectivity Element connectivity, specified as a list
133        *                    of vertex numbers, beginning with zero.
134        *\param vertex_fixed_flags Optional array to specify which vertices
135        *                    are to be marked as fixed.  If NULL,
136        *                    no vertices are fixed.
137        */
138 	MESQUITE_EXPORT
139     void fill( size_t num_vertex, const double* vtx_coords,
140                size_t num_elem, const EntityTopology* elem_types,
141                const size_t* vertex_per_elem,
142                const size_t* elem_connectivity,
143                const bool* vertex_fixed_flags,
144                MsqError& err );
145 
146       /**\brief Create global patch
147        *
148        * Create a global patch - mesh should be initialized first.
149        */
150 	MESQUITE_EXPORT
151     void fill_global_patch( MsqError& err );
152 
153 
154  	MESQUITE_EXPORT
155    void set_mesh_entities(
156                    std::vector<Mesh::ElementHandle>& patch_elems,
157                    std::vector<Mesh::VertexHandle>& free_vertices,
158                    MsqError& err );
159 
160   private:
161       //! Doesn't allow PatchData to be copied implicitly.
162       //! Mementos such as PatchDataVerticesMemento should be used when necessary.
163     PatchData(const PatchData &pd);
164       //! Doesn't allow a PatchData object to be assigned to another.
165       //! Mementos such as PatchDataVerticesMemento should be used when necessary.
166     PatchData& operator=(const PatchData &pd);
167 
168   public:
169 
170     enum ComputedInfo {
171       MIN_UNSIGNED_AREA = 0, //!< minimum volume or area out of all elements in the patch
172       MAX_UNSIGNED_AREA, //!< maximum volume or area out of all elements in the patch
173       MIN_EDGE_LENGTH, //!< minimum edge length in the patch
174       MAX_EDGE_LENGTH, //!< maximum edge length in the patch
175       MINMAX_SIGNED_DET2D, //!< minimum and maximum corner area out of all elements in the patch
176       MINMAX_SIGNED_DET3D, //!< minimum and maximum corner volume out of all elements in the patch
177       AVERAGE_DET3D, //!< average corner determinant out of all elements in the patch
178       MAX_COMPUTED_INFO_ENUM
179     };
180 
181     //! This function clears the patch information such as maximum volume, etc ...
182 	MESQUITE_EXPORT
clear_computed_info()183     void clear_computed_info() { haveComputedInfos = 0; }
184 
185 	MESQUITE_EXPORT
have_computed_info(ComputedInfo info) const186     bool have_computed_info( ComputedInfo info ) const
187       { return 0 != (haveComputedInfos&(1<<info)); }
188 
189     //! Returns the maximum volume or area out of all the elements in the patch
190     //! This information is stored in the patch and should not decrease performance
191     //! when used properly. See also PatchData::clear_computed_info() .
192 	MESQUITE_EXPORT
193     void get_minmax_element_unsigned_area(double& min, double& max, MsqError &err);
194 
195 	MESQUITE_EXPORT
196     void get_minmax_edge_length(double& min, double& max) const;
197 
198     //! Returns average corner determinant over all corners in the patch
199     //! This information is stored in the patch and should not decrease performance
200     //! when used properly. See also PatchData::clear_computed_info() .
201 //    double get_average_Lambda_3d(MsqError &err);
202 
203       //! Removes data
204  	MESQUITE_EXPORT
205    void clear();
206       //! Reorders the mesh data
207  	MESQUITE_EXPORT
208    void reorder();
209 
210       //! number of vertices in the patch.
num_nodes() const211     MESQUITE_EXPORT size_t num_nodes() const
212       { return vertexArray.size();}
num_free_vertices() const213     MESQUITE_EXPORT size_t num_free_vertices() const
214       { return numFreeVertices; }
num_slave_vertices() const215     MESQUITE_EXPORT size_t num_slave_vertices() const
216       { return numSlaveVertices; }
num_fixed_vertices() const217     MESQUITE_EXPORT size_t num_fixed_vertices() const
218       { return num_nodes() - num_free_vertices() - num_slave_vertices(); }
219       //! number of elements in the Patch.
num_elements() const220     MESQUITE_EXPORT size_t num_elements() const
221       { return elementArray.size(); }
222 
is_vertex_free(size_t index) const223     MESQUITE_EXPORT bool is_vertex_free(size_t index) const
224       { return index < numFreeVertices; }
is_vertex_not_free(size_t index) const225     MESQUITE_EXPORT bool is_vertex_not_free( size_t index ) const
226       { return index >= numFreeVertices; }
is_vertex_slave(size_t index) const227     MESQUITE_EXPORT bool is_vertex_slave(size_t index) const
228       { return index >= numFreeVertices && (index - numFreeVertices) < numSlaveVertices; }
is_vertex_fixed(size_t index) const229     MESQUITE_EXPORT bool is_vertex_fixed(size_t index) const
230       { return index >= numFreeVertices + numSlaveVertices; }
231 
232       //! number of element corners (number of vertex uses) in patch
233     MESQUITE_EXPORT size_t num_corners() const;
234 
235       //! Returns a pointer to the start of the vertex array.
236     MESQUITE_EXPORT const MsqVertex* get_vertex_array( MsqError& err ) const;
237     //MsqVertex* get_vertex_array(MsqError &err);
get_vertex_array() const238     MESQUITE_EXPORT const MsqVertex* get_vertex_array() const { return arrptr(vertexArray); }
239     //MsqVertex* get_vertex_array()             { return arrptr(vertexArray); }
240 
241       //! Returns a pointer to the start of the element array.
242     MESQUITE_EXPORT const MsqMeshEntity* get_element_array( MsqError& err ) const;
243     MESQUITE_EXPORT MsqMeshEntity* get_element_array(MsqError &err);
244 
get_connectivity_array()245     MESQUITE_EXPORT size_t* get_connectivity_array( )
246       { return arrptr(elemConnectivityArray); }
247 
get_element_handles_array()248     MESQUITE_EXPORT Mesh::ElementHandle* get_element_handles_array( )
249       { return arrptr(elementHandlesArray); }
250 
get_vertex_handles_array()251     MESQUITE_EXPORT Mesh::VertexHandle* get_vertex_handles_array()
252       { return arrptr(vertexHandlesArray); }
253 
254       //! Returns the start of the vertex->element array.
255       //! For each vertex in the patch, this array holds
256       //! the number of elements the vertex is attached to,
257       //! followed by the indices of those elements.
258     //const size_t* get_vertex_to_elem_array(MsqError &err);
259       //! Returns the start of the vertex->element offset
260       //! array (v2e_o).  For vertex i, v2e_o[i] is the
261       //! index into the vertex->element array (v2e) where
262       //! vertex i's data begins.  So, v2e[v2e_o[i]] gives
263       //! you the number of elements vertex i is attached
264       //! to, and v2e[v2e_o[i]+1] gives you the index of
265       //! the first element attached to vertex i.
266     //const size_t* get_vertex_to_elem_offset(MsqError &err);
267 
268     //MsqVertex& vertex_by_index(size_t index);
269     MESQUITE_EXPORT const MsqVertex& vertex_by_index(size_t index) const;
270     MESQUITE_EXPORT MsqMeshEntity& element_by_index(size_t index);
271     MESQUITE_EXPORT const MsqMeshEntity& element_by_index(size_t index) const;
272     MESQUITE_EXPORT size_t get_vertex_index(MsqVertex* vertex);
273     MESQUITE_EXPORT size_t get_element_index(MsqMeshEntity* element);
274 
275       //! Get the coordinates of vertices attached to the specified element
276 	MESQUITE_EXPORT
277     void get_element_vertex_coordinates(size_t elem_index,
278                                         std::vector<Vector3D> &coords,
279                                         MsqError &err);
280       /*! Get the indices of vertices of specified element. !inefficient!*/
281 	MESQUITE_EXPORT
282     void get_element_vertex_indices(size_t elem_index,
283                                     std::vector<size_t> &vertex_indices,
284                                     MsqError &err);
285       /*! Get the indices of the elements attached to the specified vertex. */
286 	MESQUITE_EXPORT
287     void get_vertex_element_indices(size_t vertex_index,
288                                     std::vector<size_t> &elem_indices,
289                                     MsqError &err);
290 
291       /** Get the indices of elements adjacent to the specified vertex,
292        *  and having the specified dimension */
293 	MESQUITE_EXPORT
294     void get_vertex_element_indices( size_t vertex_index,
295                                      unsigned element_dimension,
296                                      std::vector<size_t>& elem_indices,
297                                      MsqError& err );
298 
299       /*! Get indices of elements attached to specified vertex */
300 	MESQUITE_EXPORT
301     const size_t* get_vertex_element_adjacencies( size_t vertex_index,
302                                                   size_t& array_len_out,
303                                                   MsqError& err );
304 
305       /*! Get the indices of vertices that are attached to vertex (given by
306         vertex_index) by an element edge.
307       */
308 	MESQUITE_EXPORT
309     void get_adjacent_vertex_indices(size_t vertex_index,
310                                      std::vector<size_t> &vert_indices,
311                                      MsqError &err) const;
312 
313 
314       /*! \brief Get the indices of entities attached to entity
315 	(given by ent_ind).
316         adj_ents is filled with the indices into the entity array of elements
317         adjacent to the given element via an n-dimensional entity.
318 
319       */
320 	MESQUITE_EXPORT
321     void get_adjacent_entities_via_n_dim(int n, size_t ent_ind,
322                                          std::vector<size_t> &adj_ents,
323                                          MsqError &err);
324 
325       /*! Create the arrays that store which elements are attached
326         to each node.  If you know how many total vertex uses there are,
327         pass it in.  Otherwise the PatchData will calculate that number.
328       */
329 	MESQUITE_EXPORT
330     void generate_vertex_to_element_data();
331 
332       /*!
333       */
334 	MESQUITE_EXPORT
335     void set_vertex_coordinates(const Vector3D &coords,
336                                 size_t index,
337                                 MsqError &err);
338       /*! Add delta to the index-th free vertex in the patch
339       */
340 	MESQUITE_EXPORT
341     void move_vertex( const Vector3D &delta, size_t index, MsqError &err);
342 
343       /*! Adjust the position of the specified vertex so that it
344           lies on its constraining domain.  The actual domain constraint
345           is managed by the MeshSet's MeshDomain object.
346       */
347 	MESQUITE_EXPORT
348     void snap_vertex_to_domain(size_t vertex_index, MsqError &err);
349 
350     /*! Returns whether a domain is associated with the MeshSet from which
351         the Patch originates.
352         If false, you cannot ask for a surface normal. */
353 	MESQUITE_EXPORT
domain_set() const354     bool domain_set() const
355     { return 0 != myDomain; }
356 
357       /*\brief Get domain normal at vertex location
358        *
359        * Get the normal of the domain associated with the passed
360        * element handle at the location of the specified vertex.
361        *\param vert_index  The index of the vertex in this PatchData
362        *\param elem_handle The handle of the element passed to the domain
363        *\param normal_out  The resulting domain normal
364        *\param err         Error flag.  Possible error conditions include:
365        *                   invalid input data, no domain associated with
366        *                   element, no domain at all, etc.
367        */
368 	MESQUITE_EXPORT
369     void get_domain_normal_at_vertex( size_t vert_index,
370                                       Mesh::ElementHandle element,
371                                       Vector3D& normal_out,
372                                       MsqError& err );
373 
374       /*! Get the normal to the domain at the centroid (projected to the
375           domain) of a given element.
376           Normal is returned in Vector3D &surf_norm.  If the normal cannot
377           be determined, or if the underlying domain is not a surface,
378           the normal will be set to (0,0,0).
379           Check PatchData::domain_set() is not false first.
380       */
381 	MESQUITE_EXPORT
382     void get_domain_normal_at_element(size_t elem_index,
383                                       Vector3D &surf_norm,
384                                       MsqError &err);
385 
386       /** Get surface normals at element corners.
387        *  normals_out must be of sufficient size to hold
388        *  the normals of all the corners.
389        **/
390 	MESQUITE_EXPORT
391     void get_domain_normals_at_corners( size_t element_index,
392                                         Vector3D normals_out[],
393                                         MsqError& err ) ;
394 
395  	MESQUITE_EXPORT
396     void get_domain_normal_at_corner( size_t elemen_index,
397                                       unsigned corner,
398                                       Vector3D& normal,
399                                       MsqError& err );
400 
401 	MESQUITE_EXPORT
402     void get_domain_normal_at_mid_edge( size_t element_index,
403                                         unsigned edge_number,
404                                         Vector3D& normal,
405                                         MsqError& err );
406 
407       //! Alternative signature. Same functionality.
408 	MESQUITE_EXPORT
get_domain_normal_at_element(const MsqMeshEntity * elem_ptr,Vector3D & surf_norm,MsqError & err)409     void get_domain_normal_at_element(const MsqMeshEntity* elem_ptr,
410                                       Vector3D &surf_norm, MsqError &err)
411     { get_domain_normal_at_element(size_t(elem_ptr-&(elementArray[0])), surf_norm, err); }
412 
413 	MESQUITE_EXPORT
get_domain_normal_at_sample(size_t element_index,Sample location,Vector3D & surf_norm,MsqError & err)414     void get_domain_normal_at_sample( size_t element_index,
415                                       Sample location,
416                                       Vector3D &surf_norm, MsqError &err)
417     {
418       switch(location.dimension) {
419         case 0:
420           get_domain_normal_at_corner( element_index, location.number, surf_norm, err );
421           break;
422         case 1:
423           get_domain_normal_at_mid_edge( element_index, location.number, surf_norm, err );
424           break;
425         case 2:
426           assert(location.number == 0);
427           get_domain_normal_at_element( element_index, surf_norm, err );
428           break;
429         default:
430           MSQ_SETERR(err)("Invalid dimension for surface element subentity.\n", MsqError::INVALID_ARG );
431       }
432     }
433 
434       //! Moves free vertices and then snaps the free vertices to the domain.
435       /*\param dk an array of directions, ordered like the vertices in
436         the PatchData.
437         \param nb_vtx number of vertices.
438         \param step_size a scalar that multiplies the vectors given in dk.
439       */
440 	MESQUITE_EXPORT
441     void move_free_vertices_constrained(Vector3D dk[], size_t nb_vtx,
442                                         double step_size, MsqError &err);
443 
444     /*! Moves free vertices from a memento position along a certain direction
445       and then snaps the free vertices to the domain.
446       \param dk an array of directions, ordered like the vertices in
447       the PatchData.
448       \param nb_vtx number of vertices.
449       \param step_size a scalar that multiplies the vectors given in dk.
450     */
451 	MESQUITE_EXPORT
452     void set_free_vertices_constrained(PatchDataVerticesMemento* memento,
453                                        Vector3D dk[], size_t nb_vtx,
454                                        double step_size, MsqError &err);
455 
456     //! Project gradient vector terms onto geometric domain
457 	MESQUITE_EXPORT
458     void project_gradient( std::vector<Vector3D>& gradient, MsqError& err );
459 
460       //!Calculates the distance each vertex has moved from its original
461       //!position as defined by the PatchDataVerticesMememnto.
462 	MESQUITE_EXPORT
463     double get_max_vertex_movement_squared(PatchDataVerticesMemento* memento,
464                                            MsqError &err);
465 
466       //! Updates the underlying mesh (the MBMesquite::Mesh implementation) with
467       //! new node coordinates and flag values.
468       //!\param tag If non-null, store vertex coords in tag rather than
469       //!           updating the coords in the mesh database.  Used for
470       //!           Jacobi optimizations.
471 	MESQUITE_EXPORT
472     void update_mesh(MsqError &err, const TagHandle* tag = 0);
473 
474       //! Calculate new location for all slave higher-order nodes using
475       //! mapping function.  Called by update_mesh().
476 	MESQUITE_EXPORT
477     void update_slave_node_coordinates( MsqError& err );
478 	MESQUITE_EXPORT
479     void update_slave_node_coordinates( const size_t* elem_indices,
480                                         size_t num_elem,
481                                         MsqError& err );
482 
483       //!Remove the soft_fixed flag from all vertices in the patch.
484 	MESQUITE_EXPORT
485     void set_all_vertices_soft_free(MsqError &err);
486       //!Add a soft_fixed flag to all vertices in the patch.
487 	MESQUITE_EXPORT
488     void set_all_vertices_soft_fixed(MsqError &err);
489       //!Add a soft_fixed flag to all free vertices in the patch.
490 	MESQUITE_EXPORT
491     void set_free_vertices_soft_fixed(MsqError &err);
492 
493       //!Mark vertex as culled (soft fixed)
494 	MESQUITE_EXPORT
set_vertex_culled(size_t vtx_index)495     void set_vertex_culled( size_t vtx_index )
496       { vertexArray[vtx_index].flags() |= MsqVertex::MSQ_CULLED; }
497       //!Mark vertex as culled (soft fixed)
498 	MESQUITE_EXPORT
clear_vertex_culled(size_t vtx_index)499     void clear_vertex_culled( size_t vtx_index )
500       { vertexArray[vtx_index].flags() &= ~MsqVertex::MSQ_CULLED; }
501       //! check if vertex is culled
502 	MESQUITE_EXPORT
check_vertex_culled(size_t vtx_index) const503     int check_vertex_culled( size_t vtx_index ) const
504       { return vertexArray[vtx_index].get_flags() | MsqVertex::MSQ_CULLED; }
505 
506       //! Fills a PatchData with the elements attached to a center vertex.
507       //! Note that all entities in the sub-patch are copies of the entities
508       //! in 'this' patch.  As such, moving a vertex in the sub-patch
509       //! won't move the corresponding vertex in the source patch.  Also,
510       //! calling 'update_mesh()' on the sub-patch WILL modify the TSTT
511       //! mesh, but the source patch won't see the changes.
512 	MESQUITE_EXPORT
513     void get_subpatch(size_t center_vertex_index,
514                       unsigned num_adj_elem_layers,
515                       PatchData &pd_to_fill,
516                       MsqError &err);
517 
518 	MESQUITE_EXPORT
519     void get_free_vertex_coordinates( std::vector<Vector3D>& coords_out ) const;
520 
521       //! Creates a memento that holds the current
522       //! state of the PatchData coordinates.
523 	MESQUITE_EXPORT
524     PatchDataVerticesMemento* create_vertices_memento( MsqError &err );
525 
526       //! reinstantiates a memento to holds the current
527       //! state of the PatchData coordinates. Improves memory management.
528 	MESQUITE_EXPORT
529     void recreate_vertices_memento( PatchDataVerticesMemento* memento,
530                                     MsqError &err );
531 
532     //! Restore the PatchData coordinates to the state
533     //! contained in the memento.
534 	MESQUITE_EXPORT
535     void set_to_vertices_memento(PatchDataVerticesMemento* memento,
536                                  MsqError &err);
537 
538     //! Sets the originating meshSet. This is normally done in MeshSet::get_next_patch().
539     //! This function is only for tests purposes.
540 	MESQUITE_EXPORT
541     void set_mesh(Mesh* ms);
542 
543     //! Returns the originating meshSet.
544 	MESQUITE_EXPORT
get_mesh() const545     Mesh* get_mesh() const
546       { return myMesh; }
547 
548 	MESQUITE_EXPORT
549     void set_domain( MeshDomain* dm );
550 
551 	MESQUITE_EXPORT
get_domain() const552     MeshDomain* get_domain() const
553       { return myDomain; }
554 
555 
556 	MESQUITE_EXPORT
get_settings() const557     const Settings* get_settings() const
558       { return mSettings; }
559 	MESQUITE_EXPORT
get_mapping_function(EntityTopology type) const560     const MappingFunction* get_mapping_function( EntityTopology type ) const
561       { return mSettings->get_mapping_function(type); }
562 	MESQUITE_EXPORT
get_mapping_function_2D(EntityTopology type) const563     const MappingFunction2D* get_mapping_function_2D( EntityTopology type ) const
564       { return mSettings->get_mapping_function_2D(type); }
565 	MESQUITE_EXPORT
get_mapping_function_3D(EntityTopology type) const566     const MappingFunction3D* get_mapping_function_3D( EntityTopology type ) const
567       { return mSettings->get_mapping_function_3D(type); }
568 
569     //! Get R^3 coordinates for logical sample location.
570 	MESQUITE_EXPORT
571     void get_sample_location( size_t element_index,
572                               Sample sample,
573                               Vector3D& result,
574                               MsqError& err ) const;
575 
576     //! This function returns a NodeSet indicating which
577     //! nodes in the specified element are not slaved.
578 	MESQUITE_EXPORT
579     NodeSet non_slave_node_set( size_t elem_idx ) const;
580 
581 	MESQUITE_EXPORT
get_samples(size_t element,NodeSet non_slave_nodes) const582     NodeSet get_samples( size_t element, NodeSet non_slave_nodes ) const
583     {
584         // If we have a mapping function, use it
585       const EntityTopology type = element_by_index(element).get_element_type();
586       const MappingFunction* f;
587       if (mSettings && (f = mSettings->get_mapping_function(type)))
588         return f->sample_points( non_slave_nodes );
589         // Otherwise default to sampling at all non-slave nodes
590       non_slave_nodes.set_all_corner_nodes(type);
591       return non_slave_nodes;
592     }
593 
594 	MESQUITE_EXPORT
get_samples(size_t element) const595     NodeSet get_samples( size_t element ) const
596       { return get_samples( element, non_slave_node_set( element ) ); }
597 
598 	MESQUITE_EXPORT
599     void get_samples( size_t element, std::vector<Sample>& samples_out, MsqError& err ) const;
600 
601     //! Display the coordinates and connectivity information
602     friend std::ostream& operator<<( std::ostream&, const PatchData& );
603 
604    private:
605 
606     /* allow access to the following to functions */
607     friend class MBMesquite::ExtraData;
608     /**\brief Attach an ExtraData object to this PatchData */
609     bool attach_extra_data( ExtraData* data );
610     /**\brief Remove an ExtraData object from this PatchData */
611     bool remove_extra_data( ExtraData* data );
612 
613     /**\brief notify all attached ExtraData instances that patch contents have changed */
614     void notify_new_patch( );
615     /**\brief notify all attached ExtraData instances that a subpatch is being initalized
616      *\param sub_patch  The new, already populated subpatch
617      *\param vertex_index_map For the i-th vertex in the subpatch, the
618      *                  i-th entry in this list is the corresponding index in
619      *                  this patch.
620      *\param element_index_map For the i-th element in the subpatch, the
621      *                  i-th entry in this list is the corresponding index in
622      *                  this patch.
623      */
624     void notify_sub_patch( PatchData& sub_patch,
625                            const size_t* vertex_index_map,
626                            const size_t* element_index_map,
627                            MsqError& err );
628     /**\brief notify all attached ExtraData instances that this patch is being destroyed */
629     void notify_patch_destroyed();
630 
631       /** Call before initialize_data to change vertex_flags for
632        *  higher-order nodes to MSQ_DEPENDENT.
633        */
634     void enslave_higher_order_nodes( const size_t* element_offset_array,
635                                      unsigned char* vertex_flags,
636                                      MsqError& err ) const;
637 
638       /** Call after filling vertex handle and connectivity arrays to
639        * finish initializing the PatchData.  Reorders vertex handles array
640        * such that all higher-order nodes are at end of array, updates
641        * element connectivity array appropriately, initalizes numCornerVertices,
642        * and per-element vertex and node counts.
643        *
644        * NOTE:  If the patch contains higher-order elements, this function
645        *        will re-order the nodes in the vertex array. Do *NOT* assume
646        *        vertex indices are the same after calling this function!
647        *
648        * NOTE:  This function expects the following data to be initalized:
649        *         vertexHandlesArray
650        *         elemConnectivityArray
651        *         the topology type for all elements in elementArray
652        *        The function assumes the following data has not been
653        *        initialized and therefore does not need to be updated:
654        *         vertexArray
655        *
656        * \param elem_offset_array Offset into connectivity array for each element
657        */
658     void initialize_data( size_t* elem_offset_array, unsigned char* vertex_flags, MsqError& err );
659 
660       /** Code common to misc. methods for populating patch data.
661        *  Remove duplicates from an array of handles.
662        *\param handles    The array of handles to uniquify.
663        *\param count      As input, the lenght of the #handles
664        *                  array.  As output, the number of unique
665        *                  handles remaining in the array.
666        *\param index_map  If non-null, this must be an array of the
667        *                  same length as the handles array.  If this
668        *                  array is passed, the entry cooresponding
669        *                  to each handle in the input #handles array will
670        *                  be set to the index of that handle in the output
671        *                  array.
672        */
673     static void make_handles_unique( Mesh::EntityHandle* handles,
674                                      size_t& count,
675                                      size_t* index_map = 0 );
676 
677       /*\brief Note that the passed info has been calculated and stored */
note_have_info(ComputedInfo info)678     void note_have_info( ComputedInfo info )
679       { haveComputedInfos |= (1<<info); }
680 
681       /*\brief Update cached domain normal data */
682     void update_cached_normals( MsqError& );
683 
684     Mesh* myMesh;              //!< The Mesh used to fill this PatchData [may be NULL]
685     MeshDomain* myDomain;      //!< The geometric domain of the mesh [may be NULL]
686 
687       //! Cached data for vertices in PatchData::vertexHandlesArray,
688       //! or vertex data for a temporary patch.
689     std::vector<MsqVertex> vertexArray;
690       //! The list of handles for the vertices in this patch
691       //! May be empty if PatchData::myMesh is NULL
692     std::vector<Mesh::VertexHandle> vertexHandlesArray;
693       //! The number of vertices in PatchData::vertexArray that are
694       //! free vertices.  The vertex array is sorted such that
695       //! free vertices are first in the array.  This value
696       //! is therefore the offset at which slave vertices begin
697       //! in the array.
698     size_t numFreeVertices;
699       //! The number of slave vertices in vertexArray.  The
700       //! vertices are ordered such that all fixed vertices occur
701       //! after the slave vertices.  Therefore the offset at which
702       //! the fixed vertices begin is numFreeVertices+numSlaveVertices
703     size_t numSlaveVertices;
704       //! Cached data for elements in PatchData::elementHandlesArray
705       //! or element data for a temporary patch.
706     std::vector<MsqMeshEntity> elementArray;
707       //! The hist of handles for elements in this patch.
708       //! May be empty if PatchData::myMesh is NULL
709     std::vector<Mesh::ElementHandle> elementHandlesArray;
710       //! Element connectivity data.  The concatenation of the
711       //! connectivity list of each element in PatchData::elementArray.
712       //! Each element in PatchData::elementArray has a pointer into
713       //! this array at the correct offset for that element's connectivity data.
714     std::vector<size_t> elemConnectivityArray;
715       //! The concatenation of the adjacency lists of all the vertices
716       //! in PatchData::vertexArray.  Each value in the array is an index into
717       //! PatchData::elementArray indicating that the corresponding element uses
718       //! the vertex.  May be empty if vertex adjacency data has not been
719       //! requested.
720     std::vector<size_t> vertAdjacencyArray;
721       //! This array is indexed by vertex indices and specifies the
722       //! offset in \vertAdjacencyArray at which the adjacency list
723       //! for the corresponding vertex begins.  May be empty if vertex
724       //! adjacency data has not been requested.
725     std::vector<size_t> vertAdjacencyOffsets;
726       //! Index into normalData at which the normal for the corresponding
727       //! vertex index is located. Only vertices constrained to a single
728       //! domain with a topological dimension of 2 have a unique domain
729       //! normal.
730     std::vector<unsigned> vertexNormalIndices;
731       //! Storage space for cached domain normal data.  Pointers in
732       //! PatchData::vertexNormalPointers point into this list.
733     std::vector<Vector3D> normalData;
734       //! Storage space for cached domain DOF for vertices.  IF
735       //! a domain exists and PatchData::normalData is not empty, but
736       //! this array is, it may be assumed that all vertices have
737       //! have a DOF == 2.
738     std::vector<unsigned short> vertexDomainDOF;
739 
740       // Arrays in which to store temporary data
741       // (avoids reallocation of temp space)
742     std::vector<size_t> offsetArray;
743     std::vector<unsigned char> byteArray;
744     mutable std::vector<bool> bitMap;
745 
746       // Patch Computed Information (maxs, mins, etc ... )
747     double computedInfos[MAX_COMPUTED_INFO_ENUM];
748       // Bit map indicating which values in PatchData::computedInfos
749       // are valud (which values have been calculated.)
750     unsigned haveComputedInfos;
751 
752     ExtraData* dataList;
753 
754     const Settings* mSettings;
755     static const Settings defaultSettings;
756 
757   };
758 
759   void print_patch_data( const PatchData& pd );
760 
761 
762   /*! \brief Contains a copy of the coordinates of a PatchData.
763 
764     Use PatchDataVerticesMemento when you want to change the coordinates
765     of a PatchData object but also have the option to restore them.
766     This class can only be instantiated through PatchData::create_vertices_memento().
767   */
768   class PatchDataVerticesMemento
769   {
770   public:
clear()771     void clear()
772     {
773       originator = 0;
774       vertices.clear();
775       normalData.clear();
776     }
777   private:
778     // Constructor accessible only to originator (i.e. PatchData)
779     friend class PatchData;
PatchDataVerticesMemento()780     PatchDataVerticesMemento()
781       : originator(0)
782     {}
783 
784     PatchData* originator; //!< PatchData whose state is kept
785     std::vector<MsqVertex> vertices;
786     std::vector<Vector3D> normalData;
787   };
788 
clear()789   inline void PatchData::clear()
790   {
791     vertexArray.clear();
792     vertexHandlesArray.clear();
793     elementArray.clear();
794     elementHandlesArray.clear();
795     elemConnectivityArray.clear();
796     vertAdjacencyArray.clear();
797     vertAdjacencyOffsets.clear();
798     vertexNormalIndices.clear();
799     normalData.clear();
800     //vertexDomainDOF.clear();
801     numFreeVertices = 0;
802     numSlaveVertices = 0;
803     haveComputedInfos = 0;
804     myMesh = 0;
805     myDomain = 0;
806   }
807 
808   /*! \brief Returns an array of all vertices in the PatchData.
809   */
get_vertex_array(MsqError & err) const810   inline const MsqVertex* PatchData::get_vertex_array(MsqError &err) const
811   {
812     if (vertexArray.empty())
813       MSQ_SETERR(err)( "No vertex array defined", MsqError::INVALID_STATE );
814     return arrptr(vertexArray);
815   }
816 
817   /*! \brief Returns the PatchData element array.
818   */
get_element_array(MsqError & err) const819   inline const MsqMeshEntity* PatchData::get_element_array(MsqError &err) const
820   {
821     if (elementArray.empty())
822       MSQ_SETERR(err)( "No element array defined", MsqError::INVALID_STATE );
823     return arrptr(elementArray);
824   }
get_element_array(MsqError & err)825   inline MsqMeshEntity* PatchData::get_element_array(MsqError &err)
826   {
827     if (elementArray.empty())
828       MSQ_SETERR(err)( "No element array defined", MsqError::INVALID_STATE );
829     return arrptr(elementArray);
830   }
831 
832   /*! \brief set the coordinates of the index-th vertex in the raw array
833   */
set_vertex_coordinates(const Vector3D & coords,size_t index,MsqError & err)834   inline void PatchData::set_vertex_coordinates(const Vector3D &coords,
835                                                 size_t index,
836                                                 MsqError &err)
837   {
838     if (index >= vertexArray.size()) {
839       MSQ_SETERR(err)( "Index bigger than numVertices.", MsqError::INVALID_ARG );
840       return;
841     }
842 
843     vertexArray[index] = coords;
844 
845     if (numSlaveVertices) {
846       size_t num_elem;
847       const size_t *indices;
848       indices = get_vertex_element_adjacencies( index, num_elem, err ); MSQ_ERRRTN(err);
849       update_slave_node_coordinates( indices, num_elem, err ); MSQ_ERRRTN(err);
850     }
851   }
852   /*! \brief increment the coordinates of the index-th vertex in the raw array
853   */
move_vertex(const Vector3D & delta,size_t index,MsqError & err)854   inline void PatchData::move_vertex( const Vector3D &delta,
855                                       size_t index,
856                                       MsqError &err)
857   {
858     if (index >= vertexArray.size()) {
859       MSQ_SETERR(err)( "Index bigger than numVertices.", MsqError::INVALID_ARG );
860       return;
861     }
862 
863     vertexArray[index] += delta;
864 
865     if (numSlaveVertices) {
866       size_t num_elem;
867       const size_t *indices;
868       indices = get_vertex_element_adjacencies( index, num_elem, err ); MSQ_ERRRTN(err);
869       update_slave_node_coordinates( indices, num_elem, err ); MSQ_ERRRTN(err);
870     }
871   }
872 
873   //inline MsqVertex& PatchData::vertex_by_index(size_t index)
874   //{
875   //  return vertexArray[index];
876   //}
877 
vertex_by_index(size_t index) const878   inline const MsqVertex& PatchData::vertex_by_index( size_t index ) const
879   {
880     assert(index < vertexArray.size());
881     return vertexArray[index];
882   }
883 
element_by_index(size_t index)884   inline MsqMeshEntity& PatchData::element_by_index(size_t index)
885   {
886      assert(index < elementArray.size());
887      return elementArray[index];
888   }
889 
element_by_index(size_t index) const890   inline const MsqMeshEntity& PatchData::element_by_index( size_t index ) const
891   {
892      assert(index < elementArray.size());
893      return elementArray[index];
894   }
895 
896   /*! gets the index of a vertex in the PatchData vertex array,
897     given a pointer to the vertex. */
get_vertex_index(MsqVertex * vertex)898   inline size_t PatchData::get_vertex_index(MsqVertex* vertex)
899   {
900     return vertex - arrptr(vertexArray);
901   }
902 
get_element_index(MsqMeshEntity * element)903   inline size_t PatchData::get_element_index(MsqMeshEntity* element)
904   {
905     return element - arrptr(elementArray);
906   }
907 
get_free_vertex_coordinates(std::vector<Vector3D> & coords_out) const908   inline void PatchData::get_free_vertex_coordinates( std::vector<Vector3D>& coords_out ) const
909   {
910     coords_out.resize( num_free_vertices() );
911     std::copy( vertexArray.begin(), vertexArray.begin()+num_free_vertices(),
912                    coords_out.begin() );
913   }
914 
915 
916 
917   /*!
918     This function instantiate PatchDataVerticesMemento object and returns a pointer to it.
919     The PatchDataVerticesMemento contains the current state of the PatchData coordinates.
920     It can be used to restore the same PatchData object to those coordinates.
921 
922     It is the responsibility of the caller to discard the PatchDataVerticesMemento
923     when not needed any more.
924   */
create_vertices_memento(MsqError & err)925   inline PatchDataVerticesMemento* PatchData::create_vertices_memento(MsqError& err)
926   {
927     PatchDataVerticesMemento* memento = new PatchDataVerticesMemento;
928     recreate_vertices_memento( memento, err );
929     if (MSQ_CHKERR(err)) {
930       delete memento;
931       return 0;
932     }
933     return memento;
934   }
935 
936   /*!
937     This function reuses an existing PatchDataVerticesMemento object.
938     The PatchDataVerticesMemento contains the current state of the PatchData coordinates.
939     It can be used to restore the same PatchData object to those coordinates.
940 
941     It is the responsibility of the caller to delete the PatchDataVerticesMemento
942     when it is no longer needed.
943   */
recreate_vertices_memento(PatchDataVerticesMemento * memento,MsqError &)944   inline void PatchData::recreate_vertices_memento(PatchDataVerticesMemento* memento,
945                                                    MsqError& /*err*/)
946   {
947     memento->originator = this;
948 
949     size_t num_vtx = num_free_vertices() + num_slave_vertices();
950 
951     memento->vertices.resize( num_vtx );
952     std::copy( vertexArray.begin(), vertexArray.begin()+num_vtx, memento->vertices.begin() );
953 
954     int num_normal;
955     if (normalData.empty())
956       num_normal = 0;
957     else if (vertexNormalIndices.empty())
958       num_normal = num_vtx;
959     else {
960       num_normal = num_vtx;
961       while (num_normal != 0 && vertexNormalIndices[--num_normal] >= normalData.size());
962       if (num_normal == 0) {
963         if (vertexNormalIndices[0] < normalData.size())
964           num_normal = vertexNormalIndices[0] + 1;
965       }
966       else
967         num_normal = vertexNormalIndices[num_normal] + 1;
968     }
969 
970     memento->normalData.resize( num_normal );
971     std::copy( normalData.begin(), normalData.begin()+num_normal, memento->normalData.begin() );
972   }
973 
974   /*!
975     This function restores a PatchData object coordinates to a previous state hold in
976     a PatchDataVerticesMemento object (see create_vertices_memento() ).
977 
978     The function checks whether the memento originates from this particular PatchData object.
979     The function does not destroy the memento object: this is the caller responsibility.
980   */
set_to_vertices_memento(PatchDataVerticesMemento * memento,MsqError & err)981   inline void PatchData::set_to_vertices_memento(PatchDataVerticesMemento* memento,
982                                                  MsqError &err)
983   {
984     if (memento->originator != this)
985     {
986       MSQ_SETERR(err)("Memento may only be used to restore the PatchData "
987                       "object from which it was created.",
988                       MsqError::INVALID_ARG);
989       return;
990     }
991 
992     if (memento->vertices.size() != num_free_vertices()+num_slave_vertices())
993     {
994       MSQ_SETERR(err)("Unable to restore patch coordinates.  Number of "
995                       "vertices in PatchData has changed.",
996                       MsqError::INVALID_STATE);
997       return;
998     }
999 
1000       // copies the memento array into the PatchData array.
1001     std::copy( memento->vertices.begin(), memento->vertices.end(), vertexArray.begin() );
1002     std::copy( memento->normalData.begin(), memento->normalData.end(), normalData.begin() );
1003   }
1004 
1005 } // namespace
1006 
1007 
1008 #endif
1009