Branch data Line data Source code
1 : : #include "moab/SpectralMeshTool.hpp"
2 : : #include "moab/Interface.hpp"
3 : : #include "Internals.hpp"
4 : : #include "moab/ReadUtilIface.hpp"
5 : : #include "moab/CN.hpp"
6 : :
7 : : #include <cmath>
8 : : #include <assert.h>
9 : :
10 : : namespace moab
11 : : {
12 : : const short int SpectralMeshTool::permute_array[] = { 0, 1, 13, 25, 3, 2, 14, 26, 7, 6, 18, 30, 11, 10, 22, 34 };
13 : :
14 : : // lin_permute_array does the same to get the linear vertices of the coarse quad
15 : : const short int SpectralMeshTool::lin_permute_array[] = { 0, 25, 34, 11 };
16 : :
17 : 0 : Tag SpectralMeshTool::spectral_vertices_tag( const bool create_if_missing )
18 : : {
19 : 0 : ErrorCode rval = MB_SUCCESS;
20 [ # # ][ # # ]: 0 : if( !svTag && create_if_missing )
21 : : {
22 [ # # ]: 0 : if( !spectralOrder )
23 : : {
24 : : // should already be a spectral order tag...
25 [ # # ][ # # ]: 0 : MB_SET_ERR_RET_VAL( "Spectral order must be set before creating spectral vertices tag", 0 );
[ # # ][ # # ]
[ # # ]
26 : : }
27 : :
28 : : // create it
29 [ # # ]: 0 : std::vector< EntityHandle > dum_val( spectralOrderp1 * spectralOrderp1, 0 );
30 : : rval = mbImpl->tag_get_handle( "SPECTRAL_VERTICES", spectralOrderp1 * spectralOrderp1, MB_TYPE_HANDLE, svTag,
31 [ # # ][ # # ]: 0 : MB_TAG_DENSE | MB_TAG_CREAT, &dum_val[0] );
32 : : }
33 : :
34 [ # # ]: 0 : return ( rval == MB_SUCCESS ? svTag : 0 );
35 : : }
36 : :
37 : : /** \brief Return tag used to store spectral order
38 : : * \param create_if_missing If true, will create this tag if it doesn't exist already
39 : : */
40 : 0 : Tag SpectralMeshTool::spectral_order_tag( const bool create_if_missing )
41 : : {
42 : 0 : ErrorCode rval = MB_SUCCESS;
43 [ # # ][ # # ]: 0 : if( !soTag && create_if_missing )
44 : : {
45 : :
46 : : // create it
47 : 0 : int dum = 0;
48 [ # # ]: 0 : rval = mbImpl->tag_get_handle( "SPECTRAL_ORDER", 1, MB_TYPE_INTEGER, soTag, MB_TAG_DENSE | MB_TAG_CREAT, &dum );
49 : : }
50 : :
51 [ # # ]: 0 : return ( rval == MB_SUCCESS ? soTag : 0 );
52 : : }
53 : :
54 : : /** \brief Convert representation from coarse to fine
55 : : * Each element in set, or in interface if set is not input, is converted to fine elements, using
56 : : * vertices in SPECTRAL_VERTICES tagged array
57 : : * \param spectral_set Set containing spectral elements
58 : : */
59 : 0 : ErrorCode SpectralMeshTool::convert_to_fine( EntityHandle /* spectral_set */ )
60 : : {
61 : 0 : return MB_NOT_IMPLEMENTED;
62 : : }
63 : :
64 : : /** \brief Convert representation from fine to coarse
65 : : * Each element in set, or in interface if set is not input, is converted to coarse elements, with
66 : : * fine vertices put into SPECTRAL_VERTICES tagged array. NOTE: This function assumes that each
67 : : * order^d (fine) elements comprise each coarse element, and are in order of fine elements in each
68 : : * coarse element. If order is input as 0, looks for a SPECTRAL_ORDER tag on the mesh.
69 : : * \param order Order of the spectral mesh
70 : : * \param spectral_set Set containing spectral elements
71 : : */
72 : 0 : ErrorCode SpectralMeshTool::convert_to_coarse( int order, EntityHandle spectral_set )
73 : : {
74 [ # # ]: 0 : if( order ) spectralOrder = order;
75 [ # # ]: 0 : if( !spectralOrder )
76 [ # # ][ # # ]: 0 : { MB_SET_ERR( MB_FAILURE, "Spectral order must be set or input before converting to spectral mesh" ); }
[ # # ][ # # ]
[ # # ]
77 : :
78 [ # # ][ # # ]: 0 : Range tmp_ents, ents;
79 [ # # ]: 0 : ErrorCode rval = mbImpl->get_entities_by_handle( spectral_set, tmp_ents );
80 [ # # ][ # # ]: 0 : if( MB_SUCCESS != rval || ents.empty() ) return rval;
[ # # ][ # # ]
81 : :
82 : : // get the max-dimensional elements from it
83 [ # # ][ # # ]: 0 : ents = tmp_ents.subset_by_dimension( 3 );
84 [ # # ][ # # ]: 0 : if( ents.empty() ) ents = tmp_ents.subset_by_dimension( 2 );
[ # # ][ # # ]
85 [ # # ][ # # ]: 0 : if( ents.empty() ) ents = tmp_ents.subset_by_dimension( 1 );
[ # # ][ # # ]
86 [ # # ][ # # ]: 0 : if( ents.empty() ) { MB_SET_ERR( MB_FAILURE, "Can't find any entities for conversion" ); }
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
87 : :
88 : : // get a ptr to connectivity
89 [ # # ][ # # ]: 0 : if( ents.psize() != 1 ) { MB_SET_ERR( MB_FAILURE, "Entities must be in one chunk for conversion" ); }
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ]
90 : : EntityHandle* conn;
91 : : int count, verts_per_e;
92 [ # # ][ # # ]: 0 : rval = mbImpl->connect_iterate( ents.begin(), ents.end(), conn, verts_per_e, count );
[ # # ]
93 [ # # ][ # # ]: 0 : if( MB_SUCCESS != rval || count != (int)ents.size() ) return rval;
[ # # ][ # # ]
94 : :
95 [ # # ]: 0 : Range tmp_range;
96 [ # # ][ # # ]: 0 : return create_spectral_elems( conn, ents.size(), CN::Dimension( TYPE_FROM_HANDLE( *ents.begin() ) ), tmp_range );
[ # # ][ # # ]
[ # # ][ # # ]
97 : : }
98 : :
99 : : template < class T >
100 : 0 : ErrorCode SpectralMeshTool::create_spectral_elems( const T* conn, int num_fine_elems, int dim, Range& output_range,
101 : : int start_idx, Range* local_gids )
102 : : {
103 [ # # ][ # # ]: 0 : assert( spectralOrder && num_fine_elems );
[ # # ][ # # ]
104 : :
105 : : // now create num_coarse_elems
106 : : // compute the number of local elems, accounting for coarse or fine representation
107 : : // spectral_unit is the # fine elems per coarse elem, or spectralOrder^2
108 : 0 : int spectral_unit = spectralOrder * spectralOrder;
109 : 0 : int num_coarse_elems = num_fine_elems / spectral_unit;
110 : :
111 : : EntityHandle* new_conn;
112 : : EntityHandle start_elem;
113 : : ReadUtilIface* rmi;
114 [ # # ][ # # ]: 0 : ErrorCode rval = mbImpl->query_interface( rmi );
115 [ # # ][ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
116 : :
117 [ # # ][ # # ]: 0 : int verts_per_felem = spectralOrderp1 * spectralOrderp1, verts_per_celem = std::pow( (double)2.0, dim );
118 : :
119 [ # # ][ # # ]: 0 : rval = rmi->get_element_connect( num_coarse_elems, verts_per_celem, ( 2 == dim ? MBQUAD : MBHEX ), 0, start_elem,
120 [ # # ][ # # ]: 0 : new_conn );MB_CHK_SET_ERR( rval, "Failed to create elems" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
121 : :
122 [ # # ][ # # ]: 0 : output_range.insert( start_elem, start_elem + num_coarse_elems - 1 );
123 : :
124 : : // read connectivity into that space
125 : :
126 : : // permute_array takes a 4*order^2-long vector of integers, representing the connectivity of
127 : : // order^2 elems (fine elems in a coarse elem), and picks out the ids of the vertices necessary
128 : : // to get a lexicographically-ordered array of vertices in a spectral element of that order
129 : : // assert(verts_per_felem == (sizeof(permute_array)/sizeof(unsigned int)));
130 : :
131 : : // we're assuming here that elems was empty on input
132 : : int count;
133 : 0 : EntityHandle* sv_ptr = NULL;
134 [ # # ][ # # ]: 0 : rval = mbImpl->tag_iterate( spectral_vertices_tag( true ), output_range.begin(), output_range.end(), count,
135 [ # # ][ # # ]: 0 : (void*&)sv_ptr );MB_CHK_SET_ERR( rval, "Failed to get SPECTRAL_VERTICES ptr" );
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
[ # # ][ # # ]
136 [ # # ][ # # ]: 0 : assert( count == num_coarse_elems );
137 : 0 : int f = start_idx, fs = 0, fl = 0;
138 [ # # ][ # # ]: 0 : for( int c = 0; c < num_coarse_elems; c++ )
139 : : {
140 [ # # ][ # # ]: 0 : for( int i = 0; i < verts_per_celem; i++ )
141 : 0 : new_conn[fl + i] = conn[f + lin_permute_array[i]];
142 : 0 : fl += verts_per_celem;
143 [ # # ][ # # ]: 0 : for( int i = 0; i < verts_per_felem; i++ )
144 : 0 : sv_ptr[fs + i] = conn[f + permute_array[i]];
145 : 0 : f += verts_per_celem * spectral_unit;
146 : 0 : fs += verts_per_felem;
147 : : }
148 [ # # ][ # # ]: 0 : if( local_gids ) std::copy( sv_ptr, sv_ptr + verts_per_felem * num_coarse_elems, range_inserter( *local_gids ) );
[ # # ][ # # ]
[ # # ][ # # ]
149 : :
150 : 0 : return MB_SUCCESS;
151 : : }
152 : :
153 : : // force instantiation of a few specific types
154 : : template ErrorCode SpectralMeshTool::create_spectral_elems< int >( const int* conn, int num_fine_elems, int dim,
155 : : Range& output_range, int start_idx,
156 : : Range* local_gids );
157 : : template ErrorCode SpectralMeshTool::create_spectral_elems< EntityHandle >( const EntityHandle* conn,
158 : : int num_fine_elems, int dim,
159 : : Range& output_range, int start_idx,
160 : : Range* local_gids );
161 [ + - ][ + - ]: 228 : } // namespace moab
|