MOAB: Mesh Oriented datABase
(version 5.2.1)
|
00001 /** 00002 * \class moab::Coupler 00003 * \author Tim Tautges 00004 * 00005 * \brief This class couples data between meshes. 00006 * 00007 * The coupler interpolates solution data at a set of points. Data 00008 * being interpolated resides on a source mesh, in a tag. 00009 * Applications calling this coupler send in entities, usually points 00010 * or vertices, and receive back the tag value interpolated at those 00011 * points. Entities in the source mesh containing those points 00012 * do not have to reside on the same processor. 00013 * 00014 * To use, an application should: 00015 * - instantiate this coupler by calling the constructor collectively 00016 * on all processors in the communicator 00017 * - call locate_points, which locates the points to be interpolated and 00018 * (optionally) caches the results in this object 00019 * - call interpolate, which does the interpolation 00020 * 00021 * Multiple interpolations can be done after locating the points. 00022 * 00023 */ 00024 #ifndef COUPLER_HPP 00025 #define COUPLER_HPP 00026 00027 #include "moab/Range.hpp" 00028 #include "moab/Interface.hpp" 00029 #include "moab/CartVect.hpp" 00030 #include "moab/TupleList.hpp" 00031 00032 #include <sstream> 00033 00034 namespace moab 00035 { 00036 00037 class ParallelComm; 00038 00039 class AdaptiveKDTree; 00040 00041 class TupleList; 00042 00043 class Coupler 00044 { 00045 public: 00046 enum Method 00047 { 00048 CONSTANT, 00049 LINEAR_FE, 00050 QUADRATIC_FE, 00051 SPECTRAL, 00052 SPHERICAL 00053 }; 00054 00055 enum IntegType 00056 { 00057 VOLUME 00058 }; 00059 00060 /* Constructor 00061 * Constructor, which also optionally initializes the coupler 00062 * \param pc ParallelComm object to be used with this coupler, representing the union 00063 * of processors containing source and target meshes 00064 * \param local_elems Local elements in the source mesh 00065 * \param coupler_id Id of this coupler, should be the same over all procs 00066 * \param init_tree If true, initializes kdtree inside the constructor 00067 */ 00068 Coupler( Interface* impl, ParallelComm* pc, Range& local_elems, int coupler_id, bool init_tree = true, 00069 int max_ent_dim = 3 ); 00070 00071 /* Destructor 00072 */ 00073 virtual ~Coupler(); 00074 00075 /** 00076 * \brief Initialize the kdtree, locally and across communicator 00077 */ 00078 ErrorCode initialize_tree(); 00079 00080 /* \brief Locate points on the source mesh 00081 * This function finds the element/processor/natural coordinates for the 00082 * source mesh element containing each point, optionally storing the results 00083 * on the target mesh processor. Relative tolerance is compared to bounding 00084 * box diagonal length. Tolerance is compared to [-1,1] parametric extent 00085 * in the reference element. 00086 * \param xyz Point locations (interleaved) being located 00087 * \param num_points Number of points in xyz 00088 * \param rel_eps Relative tolerance for the non-linear iteration inside a given element 00089 * \param abs_eps Absolute tolerance for the non-linear iteration inside a given element 00090 * \param tl Tuple list containing the results, with each tuple 00091 * consisting of (p, i), p = proc, i = index on that proc 00092 * \param store_local If true, stores the tuple list in targetPts 00093 */ 00094 ErrorCode locate_points( double* xyz, unsigned int num_points, double rel_eps = 0.0, double abs_eps = 0.0, 00095 TupleList* tl = NULL, bool store_local = true ); 00096 00097 /* \brief Locate entities on the source mesh 00098 * This function finds the element/processor/natural coordinates for the 00099 * source mesh element containing each entity, optionally storing the results 00100 * on the target mesh processor. Location of each target mesh entity passed in 00101 * is its centroid (for non-vertices, the avg of its vertex positions). 00102 * Relative tolerance is compared to bounding 00103 * box diagonal length. Tolerance is compared to [-1,1] parametric extent 00104 * in the reference element. 00105 * \param ents Entities being located 00106 * \param rel_eps Relative tolerance for the non-linear iteration inside a given element 00107 * \param abs_eps Absolute tolerance for the non-linear iteration inside a given element 00108 * \param tl Tuple list containing the results, with each tuple 00109 * consisting of (p, i), p = proc, i = index on that proc 00110 * \param store_local If true, stores the tuple list in targetPts 00111 * 00112 */ 00113 ErrorCode locate_points( Range& ents, double rel_eps = 0.0, double abs_eps = 0.0, TupleList* tl = NULL, 00114 bool store_local = true ); 00115 00116 /* \brief Interpolate data from the source mesh onto points 00117 * All entities/points or, if tuple_list is input, only those points 00118 * are interpolated from the source mesh. Application should 00119 * allocate enough memory in interp_vals to hold interpolation results. 00120 * 00121 * If normalization is requested, technique used depends on the coupling 00122 * method. 00123 * 00124 * \param method Interpolation/normalization method 00125 * \param tag Tag on source mesh holding data to be interpolated 00126 * \param interp_vals Memory holding interpolated data 00127 * \param tl Tuple list of points to be interpolated, in format used by targetPts 00128 * (see documentation for targetPts below); if NULL, all locations 00129 * in targetPts are interpolated 00130 * \param normalize If true, normalization is done according to method 00131 */ 00132 ErrorCode interpolate( Coupler::Method method, Tag tag, double* interp_vals, TupleList* tl = NULL, 00133 bool normalize = true ); 00134 00135 /* \brief Interpolate data from the source mesh onto points 00136 * All entities/points or, if tuple_list is input, only those points 00137 * are interpolated from the source mesh. Application should 00138 * allocate enough memory in interp_vals to hold interpolation results. 00139 * 00140 * If normalization is requested, technique used depends on the coupling 00141 * method. 00142 * 00143 * \param methods Interpolation/normalization method 00144 * \param tag_name Name of tag on source mesh holding data to be interpolated 00145 * \param interp_vals Memory holding interpolated data 00146 * \param tl Tuple list of points to be interpolated, in format used by targetPts 00147 * (see documentation for targetPts below); if NULL, all locations 00148 * in targetPts are interpolated 00149 * \param normalize If true, normalization is done according to method 00150 */ 00151 ErrorCode interpolate( Coupler::Method method, const std::string& tag_name, double* interp_vals, 00152 TupleList* tl = NULL, bool normalize = true ); 00153 00154 /* \brief Interpolate data from multiple tags 00155 * All entities/points or, if tuple_list is input, only those points 00156 * are interpolated from the source mesh. Application should 00157 * allocate enough memory in interp_vals to hold interpolation results. 00158 * 00159 * In this variant, multiple tags, possibly with multiple interpolation 00160 * methods, are specified. Sum of values in points_per_method should be 00161 * the number of points in tl or, if NULL, targetPts. 00162 * 00163 * If normalization is requested, technique used depends on the coupling 00164 * method. 00165 * 00166 * \param methods Vector of Interpolation/normalization methods 00167 * \param tag_names Names of tag being interpolated for each method 00168 * \param points_per_method Number of points for each method 00169 * \param num_methods Length of vectors in previous 3 arguments 00170 * \param interp_vals Memory holding interpolated data 00171 * \param tl Tuple list of points to be interpolated; if NULL, all locations 00172 * stored in this object are interpolated 00173 * \param normalize If true, normalization is done according to method 00174 */ 00175 ErrorCode interpolate( Coupler::Method* methods, const std::string* tag_names, int* points_per_method, 00176 int num_methods, double* interp_vals, TupleList* tl = NULL, bool normalize = true ); 00177 00178 /* \brief Interpolate data from multiple tags 00179 * All entities/points or, if tuple_list is input, only those points 00180 * are interpolated from the source mesh. Application should 00181 * allocate enough memory in interp_vals to hold interpolation results. 00182 * 00183 * In this variant, multiple tags, possibly with multiple interpolation 00184 * methods, are specified. Sum of values in points_per_method should be 00185 * the number of points in tl or, if NULL, targetPts. 00186 * 00187 * If normalization is requested, technique used depends on the coupling 00188 * method. 00189 * 00190 * \param methods Vector of Interpolation/normalization methods 00191 * \param tag_names Names of tag being interpolated for each method 00192 * \param points_per_method Number of points for each method 00193 * \param num_methods Length of vectors in previous 3 arguments 00194 * \param interp_vals Memory holding interpolated data 00195 * \param tl Tuple list of points to be interpolated; if NULL, all locations 00196 * stored in this object are interpolated 00197 * \param normalize If true, normalization is done according to method 00198 */ 00199 ErrorCode interpolate( Coupler::Method* methods, Tag* tag_names, int* points_per_method, int num_methods, 00200 double* interp_vals, TupleList* tl = NULL, bool normalize = true ); 00201 00202 /* \brief Normalize a field over an entire mesh 00203 * A field existing on the vertices of elements of a mesh is integrated 00204 * over all elements in the mesh. The integrated value is normalized 00205 * and the normalization factor is saved to a new tag 00206 * on the mesh entity set. 00207 * 00208 * \param root_set Entity Set representing the entire mesh 00209 * \param norm_tag Tag containing field data to integrate 00210 * \param integ_type Type of integration to perform 00211 * \param num_integ_pts The number of Gaussian integration points to use in each dimension 00212 */ 00213 ErrorCode normalize_mesh( EntityHandle root_set, const char* norm_tag, Coupler::IntegType integ_type, 00214 int num_integ_pts ); 00215 00216 /* \brief Normalize a field over subsets of entities 00217 * A field existing on the vertices of elements of a mesh is integrated 00218 * over subsets of elements identified by the tags and values. The integrated 00219 * values are normalized and the normalization factor is saved to a new tag 00220 * on the entity sets which contain the elements of a subset. 00221 * 00222 * \param root_set Entity Set from the mesh from which to select subsets 00223 * \param norm_tag Tag containing field data to integrate 00224 * \param tag_names Array of tag names used for selecting element subsets 00225 * \param num_tags Number of tag names 00226 * \param tag_values Array of tag values passed as strings; the array will be 00227 * the same length as that for tag names however some entries may be 00228 * NULL indicating that tag should be matched for existence and not value 00229 * \param integ_type Type of integration to perform 00230 * \param num_integ_pts The number of Gaussian integration points to use in each dimension 00231 */ 00232 ErrorCode normalize_subset( EntityHandle root_set, const char* norm_tag, const char** tag_names, int num_tags, 00233 const char** tag_values, Coupler::IntegType integ_type, int num_integ_pts ); 00234 00235 /* \brief Normalize a field over subsets of entities 00236 * A field existing on the vertices of elements of a mesh is integrated 00237 * over subsets of elements identified by the tags and values. The integrated 00238 * values are normalized and the normalization factor is saved to a new tag 00239 * on the entity sets which contain the elements of a subset. 00240 * 00241 * \param root_set Entity Set from the mesh from which to select subsets 00242 * \param norm_tag Tag containing field data to integrate 00243 * \param tag_handles Array of tag handles used for selecting element subsets 00244 * \param num_tags Number of tag handles 00245 * \param tag_values Array of tag values passed as strings; the array will be 00246 * the same length as that for tag handles however some entries may be 00247 * NULL indicating that tag should be matched for existence and not value 00248 * \param integ_type Type of integration to perform 00249 * \param num_integ_pts The number of Gaussian integration points to use in each dimension 00250 */ 00251 ErrorCode normalize_subset( EntityHandle root_set, const char* norm_tag, Tag* tag_handles, int num_tags, 00252 const char** tag_values, Coupler::IntegType integ_type, int num_integ_pts ); 00253 00254 /* \brief Retrieve groups of entities matching tags and values(if present) 00255 * Retrieve a vector of vectors of entity handles matching the 00256 * tags and values. The entity set passed is used as the search domain. 00257 * 00258 * \param norm_tag Tag containing field data to integrate 00259 * \param entity_sets Pointer to vector of vectors of entity set handles 00260 * \param entity_groups Pointer to vector of vectors of entity handles from each entity set 00261 * \param integ_type Type of integration to perform 00262 * \param num_integ_pts The number of Gaussian integration points to use in each dimension 00263 */ 00264 ErrorCode do_normalization( const char* norm_tag, std::vector< std::vector< EntityHandle > >& entity_sets, 00265 std::vector< std::vector< EntityHandle > >& entity_groups, 00266 Coupler::IntegType integ_type, int num_integ_pts ); 00267 00268 /* \brief Retrieve groups of entities matching tags and values(if present) 00269 * Retrieve a vector of vectors of entity handles matching the 00270 * tags and values. The entity set passed is used as the search domain. 00271 * 00272 * \param root_set Set from which to search for matching entities 00273 * \param tag_names Array of tag names used to select entities 00274 * \param tag_values Array of tag values used to select entities 00275 * \param num_tags Number of tag names 00276 * \param entity_sets Pointer to vector of vectors of entity set handles found in the search 00277 * \param entity_groups Pointer to vector of vectors of entity handles from each entity set 00278 */ 00279 ErrorCode get_matching_entities( EntityHandle root_set, const char** tag_names, const char** tag_values, 00280 int num_tags, std::vector< std::vector< EntityHandle > >* entity_sets, 00281 std::vector< std::vector< EntityHandle > >* entity_groups ); 00282 00283 /* \brief Retrieve groups of entities matching tags and values(if present) 00284 * Retrieve a vector of vectors of entity handles matching the 00285 * tags and values. The entity set passed is used as the search domain. 00286 * 00287 * \param root_set Set from which to search for matching entities 00288 * \param tag_handles Array of tag handles used to select entities 00289 * \param tag_values Array of tag values used to select entities 00290 * \param num_tags Number of tag handles 00291 * \param entity_sets Pointer to vector of vectors of entity set handles found in the search 00292 * \param entity_groups Pointer to vector of vectors of entity handles from each entity set 00293 */ 00294 ErrorCode get_matching_entities( EntityHandle root_set, Tag* tag_handles, const char** tag_values, int num_tags, 00295 std::vector< std::vector< EntityHandle > >* entity_sets, 00296 std::vector< std::vector< EntityHandle > >* entity_groups ); 00297 00298 /* \brief Return an array of tuples of tag values for each Entity Set 00299 * A list of n-tuples will be constructed with 1 n-tuple for each Entity Set. 00300 * The n-tuple will have an component for each tag given. It is assumed all 00301 * of the tags are integer tags. 00302 * 00303 * \param ent_sets Array of Entity Set handles to use for retrieving tag data 00304 * \param num_sets Number of Entity Sets 00305 * \param tag_names Array of tag names 00306 * \param num_tags Number of tag names 00307 * \param tuples The returned tuple_list structure 00308 */ 00309 ErrorCode create_tuples( Range& ent_sets, const char** tag_names, unsigned int num_tags, TupleList** tuples ); 00310 00311 /* \brief Return an array of tuples of tag values for each Entity Set 00312 * A list of n-tuples will be constructed with 1 n-tuple for each Entity Set. 00313 * The n-tuple will have an component for each tag given. It is assumed all 00314 * of the tags are integer tags. 00315 * 00316 * \param ent_sets Array of Entity Set handles to use for retrieving tag data 00317 * \param num_sets Number of Entity Sets 00318 * \param tag_handles Array of tag handles 00319 * \param num_tags Number of tag handles 00320 * \param tuples The returned tuple_list structure 00321 */ 00322 ErrorCode create_tuples( Range& ent_sets, Tag* tag_handles, unsigned int num_tags, TupleList** tuples ); 00323 00324 /* \brief Consolidate an array of n-tuples lists into one n-tuple list with no duplicates 00325 * An array of list of n-tuples are consolidated into a single list of n-tuples 00326 * with all duplicates removed. Only integer columns in the tuple_list are assumed to 00327 * be used. 00328 * 00329 * \param all_tuples Array of tuple_lists to consolidate to one 00330 * \param num_tuples Number of tuple_lists 00331 * \param unique_tuples The consolidated tuple_list with no duplicates 00332 */ 00333 ErrorCode consolidate_tuples( TupleList** all_tuples, unsigned int num_tuples, TupleList** unique_tuples ); 00334 00335 /* \brief Calculate integrated field values for groups of entities 00336 * An integrated field value, as defined by the field function, 00337 * is calculated for each group of entities passed in. 00338 * 00339 * \param groups The vector contains vectors of entity handles, each representing a group 00340 * \param integ_vals The integrated field values for each group 00341 * \param norm_tag The tag name of the vertex-based field to be integrated 00342 * \param num_integ_pts The number of Gaussian integration points to use in each dimension 00343 * \param integ_type Type of integration to perform 00344 */ 00345 ErrorCode get_group_integ_vals( std::vector< std::vector< EntityHandle > >& groups, 00346 std::vector< double >& integ_vals, const char* norm_tag, int num_integ_pts, 00347 Coupler::IntegType integ_type ); 00348 00349 /* \brief Apply a normalization factor to group of entities 00350 * Multiply a normalization factor with the value of norm_tag for each vertex 00351 * of each entity in a group. Save the value back to norm_tag on each vertex. 00352 * 00353 * \param entity_sets The vector contains vectors of entity set handles, each containing the 00354 * members of a group \param norm_factors The normalization factors for each group \param 00355 * norm_tag The tag to be normalized on each group \param integ_type Type of integration to 00356 * perform 00357 */ 00358 ErrorCode apply_group_norm_factor( std::vector< std::vector< EntityHandle > >& entity_sets, 00359 std::vector< double >& norm_factors, const char* norm_tag, 00360 Coupler::IntegType integ_type ); 00361 00362 /* 00363 * this method will look at source (and target sets?) sets, and look for the SEM_DIMS tag 00364 * if it exists, it will trigger a spectral element caching, with the order specified 00365 */ 00366 ErrorCode initialize_spectral_elements( EntityHandle rootSource, EntityHandle rootTarget, bool& specSou, 00367 bool& specTar ); 00368 00369 /* 00370 * this method will put in an array, interleaved, the points of interest for coupling 00371 * with a target mesh (so where do we need to compute the field of interest) 00372 */ 00373 ErrorCode get_gl_points_on_elements( Range& targ_elems, std::vector< double >& vpos, int& numPointsOfInterest ); 00374 00375 /* Get functions */ 00376 00377 inline Interface* mb_impl() const 00378 { 00379 return mbImpl; 00380 } 00381 inline AdaptiveKDTree* my_tree() const 00382 { 00383 return myTree; 00384 } 00385 inline EntityHandle local_root() const 00386 { 00387 return localRoot; 00388 } 00389 inline const std::vector< double >& all_boxes() const 00390 { 00391 return allBoxes; 00392 } 00393 inline ParallelComm* my_pc() const 00394 { 00395 return myPc; 00396 } 00397 inline const Range& target_ents() const 00398 { 00399 return targetEnts; 00400 } 00401 inline int my_id() const 00402 { 00403 return myId; 00404 } 00405 inline const Range& my_range() const 00406 { 00407 return myRange; 00408 } 00409 inline TupleList* mapped_pts() const 00410 { 00411 return mappedPts; 00412 } 00413 inline int num_its() const 00414 { 00415 return numIts; 00416 } 00417 // used for spherical tests 00418 inline void set_spherical( bool arg1 = true ) 00419 { 00420 spherical = arg1; 00421 } 00422 00423 private: 00424 // Given a coordinate position, find all entities containing 00425 // the point and the natural coords in those ents 00426 ErrorCode nat_param( double xyz[3], std::vector< EntityHandle >& entities, std::vector< CartVect >& nat_coords, 00427 double epsilon = 0.0 ); 00428 00429 ErrorCode interp_field( EntityHandle elem, CartVect nat_coord, Tag tag, double& field ); 00430 00431 ErrorCode constant_interp( EntityHandle elem, Tag tag, double& field ); 00432 00433 ErrorCode test_local_box( double* xyz, int from_proc, int remote_index, int index, bool& point_located, 00434 double rel_eps = 0.0, double abs_eps = 0.0, TupleList* tl = NULL ); 00435 00436 /* \brief MOAB instance 00437 */ 00438 Interface* mbImpl; 00439 00440 /* \brief Kdtree for local mesh 00441 */ 00442 AdaptiveKDTree* myTree; 00443 00444 /* \brief Local root of the kdtree 00445 */ 00446 EntityHandle localRoot; 00447 00448 /* \brief Min/max bounding boxes for all proc tree roots 00449 */ 00450 std::vector< double > allBoxes; 00451 00452 /* \brief ParallelComm object for this coupler 00453 */ 00454 ParallelComm* myPc; 00455 00456 /* \brief Id of this coupler 00457 */ 00458 int myId; 00459 00460 /* \brief Range of source elements 00461 */ 00462 Range myRange; 00463 00464 /* \brief Range of target entities 00465 */ 00466 Range targetEnts; 00467 00468 /* \brief List of locally mapped tuples 00469 * Tuples contain the following: 00470 * n = # mapped points 00471 * vul[i] = local handle of mapped entity 00472 * vr[3*i..3*i+2] = natural coordinates in mapped entity 00473 */ 00474 TupleList* mappedPts; 00475 00476 /* \brief Tuple list of target points and interpolated data 00477 * Tuples contain the following: 00478 * n = # target points 00479 * vi[3*i] = remote proc mapping target point 00480 * vi[3*i+1] = local index of target point 00481 * vi[3*i+2] = remote index of target point 00482 */ 00483 TupleList* targetPts; 00484 00485 /* \brief Number of iterations of tree building before failing 00486 * 00487 */ 00488 int numIts; 00489 00490 // Entity dimension 00491 int max_dim; 00492 00493 // A cached spectral element for source and target, separate 00494 // Assume that their number of GL points (order + 1) does not change 00495 // If it does change, we need to reinitialize it 00496 void* _spectralSource; 00497 void* _spectralTarget; 00498 moab::Tag _xm1Tag, _ym1Tag, _zm1Tag; 00499 int _ntot; 00500 00501 // spherical coupling 00502 bool spherical; 00503 }; 00504 00505 inline ErrorCode Coupler::interpolate( Coupler::Method method, Tag tag, double* interp_vals, TupleList* tl, 00506 bool normalize ) 00507 { 00508 int num_pts = ( tl ? tl->get_n() : targetPts->get_n() ); 00509 return interpolate( &method, &tag, &num_pts, 1, interp_vals, tl, normalize ); 00510 } 00511 00512 } // namespace moab 00513 00514 #endif