Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
moab::LloydSmoother Class Reference

#include <LloydSmoother.hpp>

+ Collaboration diagram for moab::LloydSmoother:

Public Member Functions

 LloydSmoother (Interface *impl)
 
 LloydSmoother (Interface *impl, ParallelComm *pc, Range &elems, Tag cds_tag=0, Tag fixed_tag=0, double abs_tol=-1.0, double rel_tol=1.0e-6)
 
 ~LloydSmoother ()
 
ErrorCode perform_smooth ()
 
Interfacemb_impl ()
 
ParallelCommpcomm ()
 
void pcomm (ParallelComm *pc)
 
Rangeelems ()
 
const Rangeelems () const
 
Tag fixed_tag ()
 
void fixed_tag (Tag fixed)
 
Tag coords_tag ()
 
void coords_tag (Tag coords)
 
double abs_tol ()
 
void abs_tol (double tol)
 
double rel_tol ()
 
void rel_tol (double tol)
 
int num_its ()
 
void num_its (int num)
 
int report_its ()
 
void report_its (int num)
 

Private Member Functions

ErrorCode initialize ()
 

Private Attributes

InterfacembImpl
 
ParallelCommmyPcomm
 
Range myElems
 
Tag coordsTag
 
Tag fixedTag
 
double absTol
 
double relTol
 
int reportIts
 
int numIts
 
bool iCreatedTag
 

Detailed Description

Definition at line 26 of file LloydSmoother.hpp.

Constructor & Destructor Documentation

◆ LloydSmoother() [1/2]

moab::LloydSmoother::LloydSmoother ( Interface impl)

◆ LloydSmoother() [2/2]

moab::LloydSmoother::LloydSmoother ( Interface impl,
ParallelComm pc,
Range elems,
Tag  cds_tag = 0,
Tag  fixed_tag = 0,
double  abs_tol = -1.0,
double  rel_tol = 1.0e-6 
)

Definition at line 17 of file LloydSmoother.cpp.

18  : mbImpl( impl ), myPcomm( pc ), myElems( elms ), coordsTag( ctag ), fixedTag( ftag ), absTol( at ), relTol( rt ),
19  reportIts( 0 ), numIts( 0 ), iCreatedTag( false )
20 {
21 }

◆ ~LloydSmoother()

moab::LloydSmoother::~LloydSmoother ( )

Definition at line 23 of file LloydSmoother.cpp.

24 {
25  if( iCreatedTag && fixedTag )
26  {
27  ErrorCode rval = mbImpl->tag_delete( fixedTag );MB_CHK_SET_ERR_RET( rval, "Failed to delete the fixed tag" );
28  }
29 }

References ErrorCode, fixedTag, iCreatedTag, MB_CHK_SET_ERR_RET, mbImpl, and moab::Interface::tag_delete().

Member Function Documentation

◆ abs_tol() [1/2]

double moab::LloydSmoother::abs_tol ( )
inline

Definition at line 127 of file LloydSmoother.hpp.

128  {
129  return absTol;
130  }

References absTol.

◆ abs_tol() [2/2]

void moab::LloydSmoother::abs_tol ( double  tol)
inline

Definition at line 134 of file LloydSmoother.hpp.

135  {
136  absTol = tol;
137  }

References absTol.

◆ coords_tag() [1/2]

Tag moab::LloydSmoother::coords_tag ( )
inline

Definition at line 113 of file LloydSmoother.hpp.

114  {
115  return coordsTag;
116  }

References coordsTag.

◆ coords_tag() [2/2]

void moab::LloydSmoother::coords_tag ( Tag  coords)
inline

Definition at line 120 of file LloydSmoother.hpp.

121  {
122  coordsTag = coords;
123  }

References coordsTag.

◆ elems() [1/2]

Range& moab::LloydSmoother::elems ( )
inline

Definition at line 85 of file LloydSmoother.hpp.

86  {
87  return myElems;
88  }

References myElems.

◆ elems() [2/2]

const Range& moab::LloydSmoother::elems ( ) const
inline

Definition at line 92 of file LloydSmoother.hpp.

93  {
94  return myElems;
95  }

References myElems.

◆ fixed_tag() [1/2]

Tag moab::LloydSmoother::fixed_tag ( )
inline

Definition at line 99 of file LloydSmoother.hpp.

100  {
101  return fixedTag;
102  }

References fixedTag.

◆ fixed_tag() [2/2]

void moab::LloydSmoother::fixed_tag ( Tag  fixed)
inline

Definition at line 106 of file LloydSmoother.hpp.

107  {
108  fixedTag = fixed;
109  }

References fixedTag.

◆ initialize()

ErrorCode moab::LloydSmoother::initialize ( )
private

Definition at line 207 of file LloydSmoother.cpp.

208 {
209  ErrorCode rval = MB_SUCCESS;
210 
211  // first, check for tag; if we don't have it, make one and mark skin as fixed
212  if( !fixedTag )
213  {
214  unsigned char fixed = 0x0;
215  rval = mbImpl->tag_get_handle( "", 1, MB_TYPE_OPAQUE, fixedTag, MB_TAG_DENSE | MB_TAG_CREAT, &fixed );MB_CHK_SET_ERR( rval, "Trouble making fixed tag" );
216  iCreatedTag = true;
217 
218  // get the skin; get facets, because we might need to filter on shared entities
219  Skinner skinner( mbImpl );
220  Range skin, skin_verts;
221  rval = skinner.find_skin( 0, myElems, false, skin );MB_CHK_SET_ERR( rval, "Unable to find skin" );
222 
223 #ifdef MOAB_HAVE_MPI
224  // need to do a little extra if we're working in parallel
225  if( myPcomm )
226  {
227  // filter out ghost and interface facets
228  rval = myPcomm->filter_pstatus( skin, PSTATUS_GHOST | PSTATUS_INTERFACE, PSTATUS_NOT );MB_CHK_SET_ERR( rval, "Failed to filter on shared status" );
229  }
230 #endif
231  // get the vertices from those entities
232  rval = mbImpl->get_adjacencies( skin, 0, false, skin_verts, Interface::UNION );MB_CHK_SET_ERR( rval, "Trouble getting vertices" );
233 
234  // mark them fixed
235  std::vector< unsigned char > marks( skin_verts.size(), 0x1 );
236  rval = mbImpl->tag_set_data( fixedTag, skin_verts, &marks[0] );MB_CHK_SET_ERR( rval, "Unable to set tag on skin" );
237  }
238 
239  return MB_SUCCESS;
240 }

References ErrorCode, moab::ParallelComm::filter_pstatus(), moab::Skinner::find_skin(), fixedTag, moab::Interface::get_adjacencies(), iCreatedTag, MB_CHK_SET_ERR, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TYPE_OPAQUE, mbImpl, myElems, myPcomm, PSTATUS_GHOST, PSTATUS_INTERFACE, PSTATUS_NOT, moab::Range::size(), moab::Interface::tag_get_handle(), moab::Interface::tag_set_data(), and moab::Interface::UNION.

Referenced by perform_smooth().

◆ mb_impl()

Interface* moab::LloydSmoother::mb_impl ( )
inline

Definition at line 64 of file LloydSmoother.hpp.

65  {
66  return mbImpl;
67  }

References mbImpl.

◆ num_its() [1/2]

int moab::LloydSmoother::num_its ( )
inline

Definition at line 155 of file LloydSmoother.hpp.

156  {
157  return numIts;
158  }

References numIts.

Referenced by DeformMeshRemap::execute().

◆ num_its() [2/2]

void moab::LloydSmoother::num_its ( int  num)
inline

Definition at line 159 of file LloydSmoother.hpp.

160  {
161  numIts = num;
162  }

References numIts.

◆ pcomm() [1/2]

ParallelComm* moab::LloydSmoother::pcomm ( )
inline

Definition at line 71 of file LloydSmoother.hpp.

72  {
73  return myPcomm;
74  }

References myPcomm.

◆ pcomm() [2/2]

void moab::LloydSmoother::pcomm ( ParallelComm pc)
inline

Definition at line 78 of file LloydSmoother.hpp.

79  {
80  myPcomm = pc;
81  }

References myPcomm.

◆ perform_smooth()

ErrorCode moab::LloydSmoother::perform_smooth ( )

Definition at line 31 of file LloydSmoother.cpp.

32 {
33  ErrorCode rval;
34 
35  if( myElems.empty() )
36  {
37  MB_SET_ERR( MB_FAILURE, "No elements specified to Lloyd smoother" );
38  }
40  {
41  MB_SET_ERR( MB_FAILURE, "Elements of unequal dimension specified to Lloyd smoother" );
42  }
43 
45 
46  // first figure out tolerance to use
47  if( 0 > absTol )
48  {
49  // no tolerance set - get one relative to bounding box around elements
50  BoundBox bb;
51  rval = bb.update( *mbImpl, myElems );MB_CHK_SET_ERR( rval, "Failed to compute bounding box around elements" );
52  absTol = relTol * bb.diagonal_length();
53  }
54 
55  // initialize if we need to
56  rval = initialize();MB_CHK_SET_ERR( rval, "Failed to initialize" );
57 
58  // get all vertices
59  Range verts;
60  rval = mbImpl->get_adjacencies( myElems, 0, false, verts, Interface::UNION );MB_CHK_SET_ERR( rval, "Failed to get all vertices" );
61 
62  // perform Lloyd relaxation:
63  // 1. setup: set vertex centroids from vertex coords; filter to owned verts; get fixed tags
64 
65  // get all verts coords into tag; don't need to worry about filtering out fixed verts,
66  // we'll just be setting to their fixed coords
67  std::vector< double > vcentroids( 3 * verts.size() );
68  if( !coordsTag )
69  {
70  rval = mbImpl->get_coords( verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to get vert coords" );
71  }
72  else
73  {
74  rval = mbImpl->tag_get_data( coordsTag, verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to get vert coords tag values" );
75  }
76 
77  Tag centroid;
78  rval = mbImpl->tag_get_handle( "", 3, MB_TYPE_DOUBLE, centroid, MB_TAG_CREAT | MB_TAG_DENSE );MB_CHK_SET_ERR( rval, "Couldn't get tag handle" );
79  rval = mbImpl->tag_set_data( centroid, verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed setting centroid tag" );
80 
81  Range owned_verts, shared_owned_verts;
82 #ifdef MOAB_HAVE_MPI
83  // filter verts down to owned ones and get fixed tag for them
84  if( myPcomm && myPcomm->size() > 1 )
85  {
86  rval = myPcomm->filter_pstatus( verts, PSTATUS_NOT_OWNED, PSTATUS_NOT, -1, &owned_verts );MB_CHK_SET_ERR( rval, "Failed to filter on pstatus" );
87  // get shared owned verts, for exchanging tags
88  rval = myPcomm->filter_pstatus( owned_verts, PSTATUS_SHARED, PSTATUS_AND, -1, &shared_owned_verts );MB_CHK_SET_ERR( rval, "Failed to filter for shared owned" );
89  // workaround: if no shared owned verts, put a non-shared one in the list, to prevent
90  // exchanging tags for all shared entities
91  if( shared_owned_verts.empty() ) shared_owned_verts.insert( *verts.begin() );
92  }
93  else
94  owned_verts = verts;
95 #else
96  owned_verts = verts;
97 #endif
98 
99  std::vector< unsigned char > fix_tag( owned_verts.size() );
100  rval = mbImpl->tag_get_data( fixedTag, owned_verts, &fix_tag[0] );MB_CHK_SET_ERR( rval, "Failed to get fixed tag" );
101 
102  // now fill vcentroids array with positions of just owned vertices, since those are the ones
103  // we're actually computing
104  vcentroids.resize( 3 * owned_verts.size() );
105  rval = mbImpl->tag_get_data( centroid, owned_verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to get centroid tag" );
106 
107  // some declarations for later iterations
108  std::vector< double > fcentroids( 3 * myElems.size() ); // fcentroids for element centroids
109  std::vector< double > ctag(
110  3 * CN::MAX_NODES_PER_ELEMENT ); // temporary coordinate storage for verts bounding an element
111  const EntityHandle* conn; // const ptr & size to elem connectivity
112  int nconn;
113  Range::iterator eit, vit; // for iterating over elems, verts
114  int e, v; // for indexing into centroid vectors
115  std::vector< EntityHandle > adj_elems; // used in vertex iteration
116 
117  // 2. while !converged
118  double resid = DBL_MAX;
119  numIts = 0;
120  while( resid > absTol )
121  {
122  numIts++;
123  resid = 0.0;
124 
125  // 2a. foreach elem: centroid = sum(vertex centroids)/num_verts_in_cell
126  for( eit = myElems.begin(), e = 0; eit != myElems.end(); ++eit, e++ )
127  {
128  // get verts for this elem
129  rval = mbImpl->get_connectivity( *eit, conn, nconn );MB_CHK_SET_ERR( rval, "Failed to get connectivity" );
130  // get centroid tags for those verts
131  rval = mbImpl->tag_get_data( centroid, conn, nconn, &ctag[0] );MB_CHK_SET_ERR( rval, "Failed to get centroid" );
132  fcentroids[3 * e + 0] = fcentroids[3 * e + 1] = fcentroids[3 * e + 2] = 0.0;
133  for( v = 0; v < nconn; v++ )
134  {
135  fcentroids[3 * e + 0] += ctag[3 * v + 0];
136  fcentroids[3 * e + 1] += ctag[3 * v + 1];
137  fcentroids[3 * e + 2] += ctag[3 * v + 2];
138  }
139  for( v = 0; v < 3; v++ )
140  fcentroids[3 * e + v] /= nconn;
141  }
142  rval = mbImpl->tag_set_data( centroid, myElems, &fcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to set elem centroid" );
143 
144  // 2b. foreach owned vertex:
145  for( vit = owned_verts.begin(), v = 0; vit != owned_verts.end(); ++vit, v++ )
146  {
147  // if !fixed
148  if( fix_tag[v] ) continue;
149  // vertex centroid = sum(cell centroids)/ncells
150  adj_elems.clear();
151  rval = mbImpl->get_adjacencies( &( *vit ), 1, dim, false, adj_elems );MB_CHK_SET_ERR( rval, "Failed getting adjs" );
152  rval = mbImpl->tag_get_data( centroid, &adj_elems[0], adj_elems.size(), &fcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to get elem centroid" );
153  double vnew[] = { 0.0, 0.0, 0.0 };
154  for( e = 0; e < (int)adj_elems.size(); e++ )
155  {
156  vnew[0] += fcentroids[3 * e + 0];
157  vnew[1] += fcentroids[3 * e + 1];
158  vnew[2] += fcentroids[3 * e + 2];
159  }
160  for( e = 0; e < 3; e++ )
161  vnew[e] /= adj_elems.size();
162  double delta = ( CartVect( vnew ) - CartVect( &vcentroids[3 * v] ) ).length();
163  resid = ( v ? std::max( resid, delta ) : delta );
164  for( e = 0; e < 3; e++ )
165  vcentroids[3 * v + e] = vnew[e];
166  }
167 
168  // set the centroid tag; having them only in vcentroids array isn't enough, as vertex
169  // centroids are accessed randomly in loop over faces
170  rval = mbImpl->tag_set_data( centroid, owned_verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to set vertex centroid" );
171 
172 #ifdef MOAB_HAVE_MPI
173  // 2c. exchange tags on owned verts
174  if( myPcomm && myPcomm->size() > 1 )
175  {
176  rval = myPcomm->exchange_tags( centroid, shared_owned_verts );MB_CHK_SET_ERR( rval, "Failed to exchange tags" );
177  }
178 #endif
179 
180  if( reportIts && !( numIts % reportIts ) )
181  {
182  double global_max = resid;
183 #ifdef MOAB_HAVE_MPI
184  // global reduce for maximum delta, then report it
185  if( myPcomm && myPcomm->size() > 1 )
186  MPI_Reduce( &resid, &global_max, 1, MPI_DOUBLE, MPI_MAX, 0, myPcomm->comm() );
187  if( !myPcomm || !myPcomm->rank() )
188 #endif
189  std::cout << "Max residual = " << global_max << std::endl;
190  }
191 
192  } // end while
193 
194  // write the tag back onto vertex coordinates
195  if( !coordsTag )
196  {
197  rval = mbImpl->set_coords( owned_verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to set vertex coords" );
198  }
199  else
200  {
201  rval = mbImpl->tag_set_data( coordsTag, owned_verts, &vcentroids[0] );MB_CHK_SET_ERR( rval, "Failed to set vert coords tag values" );
202  }
203 
204  return MB_SUCCESS;
205 }

References absTol, moab::Range::begin(), moab::ParallelComm::comm(), coordsTag, moab::BoundBox::diagonal_length(), dim, moab::Interface::dimension_from_handle(), moab::Range::empty(), moab::Range::end(), ErrorCode, moab::ParallelComm::exchange_tags(), moab::ParallelComm::filter_pstatus(), fixedTag, moab::Interface::get_adjacencies(), moab::Interface::get_connectivity(), moab::Interface::get_coords(), initialize(), moab::Range::insert(), length(), moab::CN::MAX_NODES_PER_ELEMENT, MB_CHK_SET_ERR, MB_SET_ERR, MB_SUCCESS, MB_TAG_CREAT, MB_TAG_DENSE, MB_TYPE_DOUBLE, mbImpl, myElems, myPcomm, numIts, PSTATUS_AND, PSTATUS_NOT, PSTATUS_NOT_OWNED, PSTATUS_SHARED, moab::ParallelComm::rank(), moab::Range::rbegin(), relTol, reportIts, moab::Interface::set_coords(), moab::Range::size(), moab::ParallelComm::size(), moab::Interface::tag_get_data(), moab::Interface::tag_get_handle(), moab::Interface::tag_set_data(), moab::Interface::UNION, and moab::BoundBox::update().

Referenced by DeformMeshRemap::execute(), and run_local_smoother().

◆ rel_tol() [1/2]

double moab::LloydSmoother::rel_tol ( )
inline

Definition at line 141 of file LloydSmoother.hpp.

142  {
143  return relTol;
144  }

References relTol.

◆ rel_tol() [2/2]

void moab::LloydSmoother::rel_tol ( double  tol)
inline

Definition at line 148 of file LloydSmoother.hpp.

149  {
150  relTol = tol;
151  }

References relTol.

◆ report_its() [1/2]

int moab::LloydSmoother::report_its ( )
inline

Definition at line 166 of file LloydSmoother.hpp.

167  {
168  return reportIts;
169  }

References reportIts.

◆ report_its() [2/2]

void moab::LloydSmoother::report_its ( int  num)
inline

Definition at line 170 of file LloydSmoother.hpp.

171  {
172  reportIts = num;
173  }

References reportIts.

Member Data Documentation

◆ absTol

double moab::LloydSmoother::absTol
private

Definition at line 195 of file LloydSmoother.hpp.

Referenced by abs_tol(), and perform_smooth().

◆ coordsTag

Tag moab::LloydSmoother::coordsTag
private

Definition at line 189 of file LloydSmoother.hpp.

Referenced by coords_tag(), and perform_smooth().

◆ fixedTag

Tag moab::LloydSmoother::fixedTag
private

Definition at line 192 of file LloydSmoother.hpp.

Referenced by fixed_tag(), initialize(), perform_smooth(), and ~LloydSmoother().

◆ iCreatedTag

bool moab::LloydSmoother::iCreatedTag
private

Definition at line 204 of file LloydSmoother.hpp.

Referenced by initialize(), and ~LloydSmoother().

◆ mbImpl

Interface* moab::LloydSmoother::mbImpl
private

Definition at line 180 of file LloydSmoother.hpp.

Referenced by initialize(), mb_impl(), perform_smooth(), and ~LloydSmoother().

◆ myElems

Range moab::LloydSmoother::myElems
private

Definition at line 186 of file LloydSmoother.hpp.

Referenced by elems(), initialize(), and perform_smooth().

◆ myPcomm

ParallelComm* moab::LloydSmoother::myPcomm
private

Definition at line 183 of file LloydSmoother.hpp.

Referenced by initialize(), pcomm(), and perform_smooth().

◆ numIts

int moab::LloydSmoother::numIts
private

Definition at line 201 of file LloydSmoother.hpp.

Referenced by num_its(), and perform_smooth().

◆ relTol

double moab::LloydSmoother::relTol
private

Definition at line 195 of file LloydSmoother.hpp.

Referenced by perform_smooth(), and rel_tol().

◆ reportIts

int moab::LloydSmoother::reportIts
private

Definition at line 198 of file LloydSmoother.hpp.

Referenced by perform_smooth(), and report_its().


The documentation for this class was generated from the following files: