Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
TempestRemapper.hpp
Go to the documentation of this file.
1 /*
2  * =====================================================================================
3  *
4  * Filename: TempestRemapper.hpp
5  *
6  * Description: Interface to the TempestRemap library to enable intersection and
7  * high-order conservative remapping of climate solution from
8  * arbitrary resolution of source and target grids on the sphere.
9  *
10  * Author: Vijay S. Mahadevan (vijaysm), [email protected]
11  *
12  * =====================================================================================
13  */
14 
15 #ifndef MB_TEMPESTREMAPPER_HPP
16 #define MB_TEMPESTREMAPPER_HPP
17 
21 
22 // Tempest includes
23 #ifdef MOAB_HAVE_TEMPESTREMAP
24 #include "netcdfcpp.h"
25 #include "TempestRemapAPI.h"
26 #else
27 #error "This tool depends on TempestRemap library. Reconfigure using --with-tempestremap"
28 #endif
29 
30 namespace moab
31 {
32 
33 // Forward declare our friend, the mapper
34 class TempestOnlineMap;
35 
36 class TempestRemapper : public Remapper
37 {
38  public:
39 #ifdef MOAB_HAVE_MPI
40  TempestRemapper( moab::Interface* mbInt, moab::ParallelComm* pcomm = NULL, bool offlineMode = false )
41  : Remapper( mbInt, pcomm ),
42 #else
43  TempestRemapper( moab::Interface* mbInt, bool offlineMode = false )
44  : Remapper( mbInt ),
45 #endif
46  offlineWorkflow( offlineMode ), meshValidate( false ), constructEdgeMap( false ), m_source_type( DEFAULT ),
48  {
49  }
50 
51  virtual ~TempestRemapper();
52 
53  // Mesh type with a correspondence to Tempest/Climate formats
55  {
56  DEFAULT = -1,
57  CS = 0,
58  RLL = 1,
59  ICO = 2,
60  ICOD = 3,
63  OVERLAP_MOAB = 6
64  };
65 
66  friend class TempestOnlineMap;
67 
68  public:
69  /// <summary>
70  /// Initialize the TempestRemapper object internal datastructures including the mesh sets
71  /// and TempestRemap mesh references.
72  /// </summary>
73  virtual ErrorCode initialize( bool initialize_fsets = true );
74 
75  /// <summary>
76  /// Deallocate and clear any memory initialized in the TempestRemapper object
77  /// </summary>
78  virtual ErrorCode clear();
79 
80  /// <summary>
81  /// Generate a mesh in memory of given type (CS/RLL/ICO/MPAS(structured)) and store it
82  /// under the context specified by the user.
83  /// </summary>
85 
86  /// <summary>
87  /// Load a mesh from disk of given type and store it under the context specified by the
88  /// user.
89  /// </summary>
90  moab::ErrorCode LoadMesh( Remapper::IntersectionContext ctx, std::string inputFilename, TempestMeshType type );
91 
92  /// <summary>
93  /// Construct a source covering mesh such that it completely encompasses the target grid in
94  /// parallel. This operation is critical to ensure that the parallel advancing-front
95  /// intersection algorithm can the intersection mesh only locally without any process
96  /// communication.
97  /// </summary>
99  double radius_src = 1.0,
100  double radius_tgt = 1.0,
101  double boxeps = 0.1,
102  bool regional_mesh = false,
103  bool gnomonic = true,
104  int order = 1 );
105 
106  /// <summary>
107  /// Compute the intersection mesh between the source and target grids that have been
108  /// instantiated in the Remapper. This function invokes the parallel advancing-front
109  /// intersection algorithm internally for spherical meshes and can handle arbitrary
110  /// unstructured grids (CS, RLL, ICO, MPAS) with and without holes.
111  /// </summary>
112  moab::ErrorCode ComputeOverlapMesh( bool kdtree_search = true, bool use_tempest = false, int nLayers = 0 );
113 
114  /* Converters between MOAB and Tempest representations */
115 
116  /// <summary>
117  /// Convert the TempestRemap mesh object to a corresponding MOAB mesh representation
118  /// according to the intersection context.
119  /// </summary>
121 
122  /// <summary>
123  /// Convert the MOAB mesh representation to a corresponding TempestRemap mesh object
124  /// according to the intersection context.
125  /// </summary>
127 
128  /// <summary>
129  /// Get the TempestRemap mesh object according to the intersection context.
130  /// </summary>
132 
133  /// <summary>
134  /// Set the TempestRemap mesh object according to the intersection context.
135  /// </summary>
136  void SetMesh( Remapper::IntersectionContext ctx, Mesh* mesh, bool overwrite = true );
137 
138  void SetMeshSet( Remapper::IntersectionContext ctx /* Remapper::CoveringMesh*/,
139  moab::EntityHandle mset,
140  moab::Range* entities = nullptr );
141  /// <summary>
142  /// Get the covering mesh (TempestRemap) object.
143  /// </summary>
144  Mesh* GetCoveringMesh();
145 
146  /// <summary>
147  /// Get the MOAB mesh set corresponding to the intersection context.
148  /// </summary>
150 
151  /// <summary>
152  /// Const overload. Get the MOAB mesh set corresponding to the intersection context.
153  /// </summary>
155 
156  /// <summary>
157  /// Get the mesh element entities corresponding to the intersection context.
158  /// </summary>
160 
161  /// <summary>
162  /// Const overload. Get the mesh element entities corresponding to the intersection context.
163  /// </summary>
165 
166  /// <summary>
167  /// Get the mesh vertices corresponding to the intersection context. Useful for point-cloud
168  /// meshes
169  /// </summary>
171 
172  /// <summary>
173  /// Const overload. Get the mesh vertices corresponding to the intersection context. Useful
174  /// for point-cloud meshes.
175  /// </summary>
177 
178  /// <summary>
179  /// Get access to the underlying source covering set if available. Else return the source
180  /// set.
181  /// </summary>
183 
184  /// <summary>
185  /// Set the mesh type corresponding to the intersection context
186  /// </summary>
187  void SetMeshType( Remapper::IntersectionContext ctx, const std::vector< int >& metadata );
188 
189  /// <summary>
190  /// reconstruct mesh, used now only for IO; need a better solution maybe
191  /// </summary>
193 
194  /// <summary>
195  /// Get the mesh type corresponding to the intersection context
196  /// </summary>
198 
199  /// <summary>
200  /// Get the global ID corresponding to the local entity ID according to the context (source,
201  /// target, intersection)
202  /// </summary>
203  int GetGlobalID( Remapper::IntersectionContext ctx, int localID );
204 
205  /// <summary>
206  /// Get the local ID corresponding to the global entity ID according to the context (source,
207  /// target, intersection)
208  /// </summary>
209  int GetLocalID( Remapper::IntersectionContext ctx, int globalID );
210 
211  /// <summary>
212  /// Gather the overlap mesh and asssociated source/target data and write it out to disk
213  /// using the TempestRemap output interface. This information can then be used with the
214  /// "GenerateOfflineMap" tool in TempestRemap as needed.
215  /// </summary>
216  moab::ErrorCode WriteTempestIntersectionMesh( std::string strOutputFileName,
217  const bool fAllParallel,
218  const bool fInputConcave,
219  const bool fOutputConcave );
220 
221  /// <summary>
222  /// Generate the necessary metadata and specifically the GLL node numbering for DoFs for a
223  /// CS mesh. This negates the need for running external code like HOMME to output the
224  /// numbering needed for computing maps. The functionality is used through the `mbconvert`
225  /// tool to compute processor-invariant Global DoF IDs at GLL nodes.
226  /// </summary>
227  moab::ErrorCode GenerateCSMeshMetadata( const int ntot_elements,
229  moab::Range* secondary_entities,
230  const std::string& dofTagName,
231  int nP );
232 
233  /// <summary>
234  /// Generate the necessary metadata for DoF node numbering in a given mesh.
235  /// Currently, only the functionality to generate numbering on CS grids is supported.
236  /// </summary>
238  const int ntot_elements,
240  moab::Range* secondary_entities,
241  const std::string dofTagName,
242  int nP );
243 
244  /// <summary>
245  /// Compute the local and global IDs for elements in source/target/coverage meshes.
246  /// </summary>
248 
249  /// <summary>
250  /// Get all the ghosted overlap entities that were accumulated to enable conservation in
251  /// parallel
252  /// </summary>
254 
255 #ifndef MOAB_HAVE_MPI
256  /// <summary>
257  /// Internal method to assign vertex and element global IDs if one does not exist already
258  /// </summary>
260  EntityHandle this_set,
261  const int dimension = 2,
262  const int start_id = 1 );
263 #endif
264  // <summary>
265  /// Get the masks that could have been defined
266  /// </summary>
267  ErrorCode GetIMasks( Remapper::IntersectionContext ctx, std::vector< int >& masks );
268 
269  public: // public members
270  const bool offlineWorkflow; // check whether we are in an offline workflow context (mbtempest)
271  bool meshValidate; // Validate the mesh after loading from file
272 
273  bool constructEdgeMap; // Construct the edge map within the TempestRemap datastructures
274 
275  static const bool verbose = true;
276 
277  private:
279 
280  // private methods
281  moab::ErrorCode load_tempest_mesh_private( std::string inputFilename, Mesh** tempest_mesh );
282 
284  moab::EntityHandle meshset,
286  moab::Range* pverts );
287 
289  Mesh* mesh,
290  moab::EntityHandle& meshset,
292  moab::Range* vertices );
293 
295 
296  /* Source meshset, mesh and entity references */
297  Mesh* m_source;
304  std::vector< int > m_source_metadata;
305 
306  /* Target meshset, mesh and entity references */
307  Mesh* m_target;
314  std::vector< int > m_target_metadata;
315 
316  /* Overlap meshset, mesh and entity references */
317  Mesh* m_overlap;
321  std::vector< std::pair< int, int > > m_sorted_overlap_order;
322 
323  /* Intersection context on a sphere */
325 
326  /* Parallel - migrated mesh that is in the local view */
331 
332  /* local to glboal and global to local ID maps */
335 
337 
338  // this is a copy of the original source/target set,
339  // with some ghosts created for parallel cases
342 
343  bool rrmgrids;
345  int rank, size;
346 };
347 
348 // Inline functions
350 {
351  switch( ctx )
352  {
354  return m_source;
356  return m_target;
358  return m_overlap;
360  return m_covering_source;
361  case Remapper::DEFAULT: // Remapper::SourceMeshWithGhosts does not need one
362  default:
363  return NULL;
364  }
365 }
366 
367 inline void TempestRemapper::SetMesh( Remapper::IntersectionContext ctx, Mesh* mesh, bool overwrite )
368 {
369  switch( ctx )
370  {
372  if( !overwrite && m_source ) return;
373  if( overwrite && m_source ) delete m_source;
374  m_source = mesh;
375  break;
377  if( !overwrite && m_target ) return;
378  if( overwrite && m_target ) delete m_target;
379  m_target = mesh;
380  break;
382  if( !overwrite && m_overlap ) return;
383  if( overwrite && m_overlap ) delete m_overlap;
384  m_overlap = mesh;
385  break;
387  if( !overwrite && m_covering_source ) return;
388  if( overwrite && m_covering_source ) delete m_covering_source;
389  m_covering_source = mesh;
390  break;
391  case Remapper::DEFAULT: // Remapper::SourceMeshWithGhosts does not need one
392  default:
393  break;
394  }
395 }
396 
398 {
399  switch( ctx )
400  {
402  delete m_source;
403  m_source = new Mesh;
404  m_source_set = meshSet;
406  m_source->CalculateFaceAreas( false ); // fInputConcave is false ?
407  break;
409  // not needed yet
410  break;
412  // not needed yet
413  break;
415  // not needed yet
416  break;
417  case Remapper::DEFAULT:
418  default:
419  break;
420  }
421 }
422 
423 
425 {
426  switch( ctx )
427  {
429  return m_source_set;
431  return m_target_set;
433  return m_overlap_set;
435  return m_covering_source_set;
440  case Remapper::DEFAULT:
441  default:
442  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_set );
443  }
444 }
445 
447 {
448  switch( ctx )
449  {
451  return m_source_set;
453  return m_target_set;
455  return m_overlap_set;
457  return m_covering_source_set;
462  case Remapper::DEFAULT:
463  default:
464  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_set );
465  }
466 }
467 
469 {
470  switch( ctx )
471  {
473  return m_source_entities;
475  return m_target_entities;
477  return m_overlap_entities;
480  case Remapper::DEFAULT:
481  default:
482  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_entities );
483  }
484 }
485 
487 {
488  switch( ctx )
489  {
491  return m_source_entities;
493  return m_target_entities;
495  return m_overlap_entities;
498  case Remapper::DEFAULT:
499  default:
500  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_entities );
501  }
502 }
503 
505 {
506  switch( ctx )
507  {
509  return m_source_vertices;
511  return m_target_vertices;
514  case Remapper::DEFAULT:
515  default:
516  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_source_vertices );
517  }
518 }
519 
521 {
522  switch( ctx )
523  {
525  return m_source_vertices;
527  return m_target_vertices;
530  case Remapper::DEFAULT:
531  default:
532  MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_source_vertices );
533  }
534 }
535 
536 inline void TempestRemapper::SetMeshType( Remapper::IntersectionContext ctx, const std::vector< int >& metadata )
537 {
538  switch( ctx )
539  {
541  m_source_type = static_cast< moab::TempestRemapper::TempestMeshType >( metadata[0] );
542  if( metadata[0] == 1 ) // RLL mesh
543  {
544  m_source_metadata.resize( 2 );
545  m_source_metadata[0] = metadata[1];
546  m_source_metadata[1] = metadata[2];
547  }
548  else
549  {
550  m_source_metadata.resize( 1 );
551  m_source_metadata[0] = metadata[1];
552  }
553  break;
555  m_target_type = static_cast< moab::TempestRemapper::TempestMeshType >( metadata[0] );
556  if( metadata[0] == 1 ) // RLL mesh
557  {
558  m_target_metadata.resize( 2 );
559  m_target_metadata[0] = metadata[1];
560  m_target_metadata[1] = metadata[2];
561  }
562  else
563  {
564  m_target_metadata.resize( 1 );
565  m_target_metadata[0] = metadata[1];
566  }
567  break;
570  default:
571  break;
572  }
573 }
574 
576 {
577  switch( ctx )
578  {
580  return m_source_type;
582  return m_target_type;
584  return m_overlap_type;
585  case Remapper::DEFAULT: // not need yet case Remapper::SourceMeshWithGhosts
586  default:
588  }
589 }
590 
592 {
593  return m_covering_source;
594 }
595 
597 {
598  return m_covering_source_set;
599 }
600 
602 {
603  switch( ctx )
604  {
606  return lid_to_gid_src[localID];
608  return lid_to_gid_tgt[localID];
610  return lid_to_gid_covsrc[localID];
613  case Remapper::DEFAULT:
614  default:
615  return -1;
616  }
617 }
618 
620 {
621  switch( ctx )
622  {
624  return gid_to_lid_src[globalID];
626  return gid_to_lid_tgt[globalID];
628  return gid_to_lid_covsrc[globalID];
629  case Remapper::DEFAULT:
632  default:
633  return -1;
634  }
635 }
636 
637 } // namespace moab
638 
639 #endif // MB_TEMPESTREMAPPER_HPP