MOAB: Mesh Oriented datABase  (version 5.2.1)
TempestRemapper.hpp
Go to the documentation of this file.
00001 /*
00002  * =====================================================================================
00003  *
00004  *       Filename:  TempestRemapper.hpp
00005  *
00006  *    Description:  Interface to the TempestRemap library to enable intersection and
00007  *                  high-order conservative remapping of climate solution from
00008  *                  arbitrary resolution of source and target grids on the sphere.
00009  *
00010  *         Author:  Vijay S. Mahadevan (vijaysm), mahadevan@anl.gov
00011  *
00012  * =====================================================================================
00013  */
00014 
00015 #ifndef MB_TEMPESTREMAPPER_HPP
00016 #define MB_TEMPESTREMAPPER_HPP
00017 
00018 #include "moab/Remapping/Remapper.hpp"
00019 #include "moab/IntxMesh/Intx2MeshOnSphere.hpp"
00020 #include "moab/IntxMesh/IntxUtils.hpp"
00021 
00022 // Tempest includes
00023 #ifdef MOAB_HAVE_TEMPESTREMAP
00024 #include "netcdfcpp.h"
00025 #include "TempestRemapAPI.h"
00026 #else
00027 #error "This tool depends on TempestRemap library. Reconfigure using --with-tempestremap"
00028 #endif
00029 
00030 namespace moab
00031 {
00032 
00033 // Forward declare our friend, the mapper
00034 class TempestOnlineMap;
00035 
00036 class TempestRemapper : public Remapper
00037 {
00038   public:
00039 #ifdef MOAB_HAVE_MPI
00040     TempestRemapper( moab::Interface* mbInt, moab::ParallelComm* pcomm = NULL )
00041         : Remapper( mbInt, pcomm ),
00042 #else
00043     TempestRemapper( moab::Interface* mbInt )
00044         : Remapper( mbInt ),
00045 #endif
00046           meshValidate( false ), constructEdgeMap( false ), m_source_type( DEFAULT ), m_target_type( DEFAULT )
00047     {
00048     }
00049 
00050     virtual ~TempestRemapper();
00051 
00052     // Mesh type with a correspondence to Tempest/Climate formats
00053     enum TempestMeshType
00054     {
00055         DEFAULT        = -1,
00056         CS             = 0,
00057         RLL            = 1,
00058         ICO            = 2,
00059         ICOD           = 3,
00060         OVERLAP_FILES  = 4,
00061         OVERLAP_MEMORY = 5,
00062         OVERLAP_MOAB   = 6
00063     };
00064 
00065     friend class TempestOnlineMap;
00066 
00067   public:
00068     /// <summary>
00069     ///     Initialize the TempestRemapper object internal datastructures including the mesh sets
00070     ///     and TempestRemap mesh references.
00071     /// </summary>
00072     virtual ErrorCode initialize( bool initialize_fsets = true );
00073 
00074     /// <summary>
00075     ///     Deallocate and clear any memory initialized in the TempestRemapper object
00076     /// </summary>
00077     virtual ErrorCode clear();
00078 
00079     /// <summary>
00080     ///     Generate a mesh in memory of given type (CS/RLL/ICO/MPAS(structured)) and store it
00081     ///     under the context specified by the user.
00082     /// </summary>
00083     moab::ErrorCode GenerateMesh( Remapper::IntersectionContext ctx, TempestMeshType type );
00084 
00085     /// <summary>
00086     ///     Load a mesh from disk of given type and store it under the context specified by the
00087     ///     user.
00088     /// </summary>
00089     moab::ErrorCode LoadMesh( Remapper::IntersectionContext ctx, std::string inputFilename, TempestMeshType type );
00090 
00091     /// <summary>
00092     ///     Construct a source covering mesh such that it completely encompasses the target grid in
00093     ///     parallel. This operation is critical to ensure that the parallel advancing-front
00094     ///     intersection algorithm can the intersection mesh only locally without any process
00095     ///     communication.
00096     /// </summary>
00097     moab::ErrorCode ConstructCoveringSet( double tolerance = 1e-8, double radius_src = 1.0, double radius_tgt = 1.0,
00098                                           double boxeps = 0.1, bool regional_mesh = false );
00099 
00100     /// <summary>
00101     ///     Compute the intersection mesh between the source and target grids that have been
00102     ///     instantiated in the Remapper. This function invokes the parallel advancing-front
00103     ///     intersection algorithm internally for spherical meshes and can handle arbitrary
00104     ///     unstructured grids (CS, RLL, ICO, MPAS) with and without holes.
00105     /// </summary>
00106     moab::ErrorCode ComputeOverlapMesh( bool kdtree_search = true, bool use_tempest = false );
00107 
00108     /* Converters between MOAB and Tempest representations */
00109 
00110     /// <summary>
00111     ///     Convert the TempestRemap mesh object to a corresponding MOAB mesh representation
00112     ///     according to the intersection context.
00113     /// </summary>
00114     moab::ErrorCode ConvertTempestMesh( Remapper::IntersectionContext ctx );
00115 
00116     /// <summary>
00117     ///     Convert the MOAB mesh representation to a corresponding TempestRemap mesh object
00118     ///     according to the intersection context.
00119     /// </summary>
00120     moab::ErrorCode ConvertMeshToTempest( Remapper::IntersectionContext ctx );
00121 
00122     /// <summary>
00123     ///     Get the TempestRemap mesh object according to the intersection context.
00124     /// </summary>
00125     Mesh* GetMesh( Remapper::IntersectionContext ctx );
00126 
00127     /// <summary>
00128     ///     Set the TempestRemap mesh object according to the intersection context.
00129     /// </summary>
00130     void SetMesh( Remapper::IntersectionContext ctx, Mesh* mesh, bool overwrite = true );
00131 
00132     /// <summary>
00133     ///     Get the covering mesh (TempestRemap) object.
00134     /// </summary>
00135     Mesh* GetCoveringMesh();
00136 
00137     /// <summary>
00138     ///     Get the MOAB mesh set corresponding to the intersection context.
00139     /// </summary>
00140     moab::EntityHandle& GetMeshSet( Remapper::IntersectionContext ctx );
00141 
00142     /// <summary>
00143     ///     Const overload. Get the MOAB mesh set corresponding to the intersection context.
00144     /// </summary>
00145     moab::EntityHandle GetMeshSet( Remapper::IntersectionContext ctx ) const;
00146 
00147     /// <summary>
00148     ///     Get the mesh element entities corresponding to the intersection context.
00149     /// </summary>
00150     moab::Range& GetMeshEntities( Remapper::IntersectionContext ctx );
00151 
00152     /// <summary>
00153     ///     Const overload. Get the mesh element entities corresponding to the intersection context.
00154     /// </summary>
00155     const moab::Range& GetMeshEntities( Remapper::IntersectionContext ctx ) const;
00156 
00157     /// <summary>
00158     ///     Get the mesh vertices corresponding to the intersection context. Useful for point-cloud
00159     ///     meshes
00160     /// </summary>
00161     moab::Range& GetMeshVertices( Remapper::IntersectionContext ctx );
00162 
00163     /// <summary>
00164     ///     Const overload. Get the mesh vertices corresponding to the intersection context. Useful
00165     ///     for point-cloud meshes.
00166     /// </summary>
00167     const moab::Range& GetMeshVertices( Remapper::IntersectionContext ctx ) const;
00168 
00169     /// <summary>
00170     ///     Get access to the underlying source covering set if available. Else return the source
00171     ///     set.
00172     /// </summary>
00173     moab::EntityHandle& GetCoveringSet();
00174 
00175     /// <summary>
00176     ///     Set the mesh type corresponding to the intersection context
00177     /// </summary>
00178     void SetMeshType( Remapper::IntersectionContext ctx, TempestMeshType type );
00179 
00180     /// <summary>
00181     ///     Get the mesh type corresponding to the intersection context
00182     /// </summary>
00183     TempestMeshType GetMeshType( Remapper::IntersectionContext ctx ) const;
00184 
00185     /// <summary>
00186     ///     Get the global ID corresponding to the local entity ID according to the context (source,
00187     ///     target, intersection)
00188     /// </summary>
00189     int GetGlobalID( Remapper::IntersectionContext ctx, int localID );
00190 
00191     /// <summary>
00192     ///     Get the local ID corresponding to the global entity ID according to the context (source,
00193     ///     target, intersection)
00194     /// </summary>
00195     int GetLocalID( Remapper::IntersectionContext ctx, int globalID );
00196 
00197     /// <summary>
00198     ///     Gather the overlap mesh and asssociated source/target data and write it out to disk
00199     ///     using the TempestRemap output interface. This information can then be used with the
00200     ///     "GenerateOfflineMap" tool in TempestRemap as needed.
00201     /// </summary>
00202     moab::ErrorCode WriteTempestIntersectionMesh( std::string strOutputFileName, const bool fAllParallel,
00203                                                   const bool fInputConcave, const bool fOutputConcave );
00204 
00205     /// <summary>
00206     ///     Generate the necessary metadata and specifically the GLL node numbering for DoFs for a
00207     ///     CS mesh. This negates the need for running external code like HOMME to output the
00208     ///     numbering needed for computing maps. The functionality is used through the `mbconvert`
00209     ///     tool to compute processor-invariant Global DoF IDs at GLL nodes.
00210     /// </summary>
00211     moab::ErrorCode GenerateCSMeshMetadata( const int ntot_elements, moab::Range& entities,
00212                                             moab::Range* secondary_entities, const std::string dofTagName, int nP );
00213 
00214     /// <summary>
00215     ///     Generate the necessary metadata for DoF node numbering in a given mesh.
00216     ///     Currently, only the functionality to generate numbering on CS grids is supported.
00217     /// </summary>
00218     moab::ErrorCode GenerateMeshMetadata( Mesh& mesh, const int ntot_elements, moab::Range& entities,
00219                                           moab::Range* secondary_entities, const std::string dofTagName, int nP );
00220 
00221     /// <summary>
00222     ///     Compute the local and global IDs for elements in source/target/coverage meshes.
00223     /// </summary>
00224     moab::ErrorCode ComputeGlobalLocalMaps();
00225 
00226     /// <summary>
00227     ///     Get all the ghosted overlap entities that were accumulated to enable conservation in
00228     /// parallel
00229     /// </summary>
00230     moab::ErrorCode GetOverlapAugmentedEntities( moab::Range& sharedGhostEntities );
00231 
00232 #ifndef MOAB_HAVE_MPI
00233     /// <summary>
00234     ///    Internal method to assign vertex and element global IDs if one does not exist already
00235     /// </summary>
00236     moab::ErrorCode assign_vertex_element_IDs( Tag idtag, EntityHandle this_set,
00237                                             const int dimension = 2, const int start_id = 1 );
00238 #endif
00239 
00240   public:               // public members
00241     bool meshValidate;  // Validate the mesh after loading from file
00242 
00243     bool constructEdgeMap;  //  Construct the edge map within the TempestRemap datastructures
00244 
00245     static const bool verbose = true;
00246 
00247   private:
00248     moab::ErrorCode convert_overlap_mesh_sorted_by_source();
00249 
00250     // private methods
00251     moab::ErrorCode load_tempest_mesh_private( std::string inputFilename, Mesh** tempest_mesh );
00252 
00253     moab::ErrorCode convert_mesh_to_tempest_private( Mesh* mesh, moab::EntityHandle meshset, moab::Range& entities,
00254                                                      moab::Range* pverts );
00255 
00256     moab::ErrorCode convert_tempest_mesh_private( TempestMeshType type, Mesh* mesh, moab::EntityHandle& meshset,
00257                                                   moab::Range& entities, moab::Range* vertices );
00258 
00259     moab::ErrorCode augment_overlap_set();
00260 
00261     /* Source meshset, mesh and entity references */
00262     Mesh* m_source;
00263     TempestMeshType m_source_type;
00264     moab::Range m_source_entities;
00265     moab::Range m_source_vertices;
00266     moab::EntityHandle m_source_set;
00267     int max_source_edges;
00268     bool point_cloud_source;
00269 
00270     /* Target meshset, mesh and entity references */
00271     Mesh* m_target;
00272     TempestMeshType m_target_type;
00273     moab::Range m_target_entities;
00274     moab::Range m_target_vertices;
00275     moab::EntityHandle m_target_set;
00276     int max_target_edges;
00277     bool point_cloud_target;
00278 
00279     /* Overlap meshset, mesh and entity references */
00280     Mesh* m_overlap;
00281     TempestMeshType m_overlap_type;
00282     moab::Range m_overlap_entities;
00283     moab::EntityHandle m_overlap_set;
00284     std::vector< std::pair< int, int > > m_sorted_overlap_order;
00285 
00286     /* Intersection context on a sphere */
00287     moab::Intx2MeshOnSphere* mbintx;
00288 
00289     /* Parallel - migrated mesh that is in the local view */
00290     Mesh* m_covering_source;
00291     moab::EntityHandle m_covering_source_set;
00292     moab::Range m_covering_source_entities;
00293     moab::Range m_covering_source_vertices;
00294 
00295     /* local to glboal and global to local ID maps */
00296     std::map< int, int > gid_to_lid_src, gid_to_lid_covsrc, gid_to_lid_tgt;
00297     std::map< int, int > lid_to_gid_src, lid_to_gid_covsrc, lid_to_gid_tgt;
00298 
00299     IntxAreaUtils::AreaMethod m_area_method;
00300 
00301     bool rrmgrids;
00302     bool is_parallel, is_root;
00303     int rank, size;
00304 };
00305 
00306 // Inline functions
00307 inline Mesh* TempestRemapper::GetMesh( Remapper::IntersectionContext ctx )
00308 {
00309     switch( ctx )
00310     {
00311         case Remapper::SourceMesh:
00312             return m_source;
00313         case Remapper::TargetMesh:
00314             return m_target;
00315         case Remapper::OverlapMesh:
00316             return m_overlap;
00317         case Remapper::CoveringMesh:
00318             return m_covering_source;
00319         case Remapper::DEFAULT:
00320         default:
00321             return NULL;
00322     }
00323 }
00324 
00325 inline void TempestRemapper::SetMesh( Remapper::IntersectionContext ctx, Mesh* mesh, bool overwrite )
00326 {
00327     switch( ctx )
00328     {
00329         case Remapper::SourceMesh:
00330             if( !overwrite && m_source ) return;
00331             if( overwrite && m_source ) delete m_source;
00332             m_source = mesh;
00333             break;
00334         case Remapper::TargetMesh:
00335             if( !overwrite && m_target ) return;
00336             if( overwrite && m_target ) delete m_target;
00337             m_target = mesh;
00338             break;
00339         case Remapper::OverlapMesh:
00340             if( !overwrite && m_overlap ) return;
00341             if( overwrite && m_overlap ) delete m_overlap;
00342             m_overlap = mesh;
00343             break;
00344         case Remapper::CoveringMesh:
00345             if( !overwrite && m_covering_source ) return;
00346             if( overwrite && m_covering_source ) delete m_covering_source;
00347             m_covering_source = mesh;
00348             break;
00349         case Remapper::DEFAULT:
00350         default:
00351             break;
00352     }
00353 }
00354 
00355 inline moab::EntityHandle& TempestRemapper::GetMeshSet( Remapper::IntersectionContext ctx )
00356 {
00357     switch( ctx )
00358     {
00359         case Remapper::SourceMesh:
00360             return m_source_set;
00361         case Remapper::TargetMesh:
00362             return m_target_set;
00363         case Remapper::OverlapMesh:
00364             return m_overlap_set;
00365         case Remapper::CoveringMesh:
00366             return m_covering_source_set;
00367         case Remapper::DEFAULT:
00368         default:
00369             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_set );
00370     }
00371 }
00372 
00373 inline moab::EntityHandle TempestRemapper::GetMeshSet( Remapper::IntersectionContext ctx ) const
00374 {
00375     switch( ctx )
00376     {
00377         case Remapper::SourceMesh:
00378             return m_source_set;
00379         case Remapper::TargetMesh:
00380             return m_target_set;
00381         case Remapper::OverlapMesh:
00382             return m_overlap_set;
00383         case Remapper::CoveringMesh:
00384             return m_covering_source_set;
00385         case Remapper::DEFAULT:
00386         default:
00387             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_set );
00388     }
00389 }
00390 
00391 inline moab::Range& TempestRemapper::GetMeshEntities( Remapper::IntersectionContext ctx )
00392 {
00393     switch( ctx )
00394     {
00395         case Remapper::SourceMesh:
00396             return m_source_entities;
00397         case Remapper::TargetMesh:
00398             return m_target_entities;
00399         case Remapper::OverlapMesh:
00400             return m_overlap_entities;
00401         case Remapper::CoveringMesh:
00402             return m_covering_source_entities;
00403         case Remapper::DEFAULT:
00404         default:
00405             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_entities );
00406     }
00407 }
00408 
00409 inline const moab::Range& TempestRemapper::GetMeshEntities( Remapper::IntersectionContext ctx ) const
00410 {
00411     switch( ctx )
00412     {
00413         case Remapper::SourceMesh:
00414             return m_source_entities;
00415         case Remapper::TargetMesh:
00416             return m_target_entities;
00417         case Remapper::OverlapMesh:
00418             return m_overlap_entities;
00419         case Remapper::CoveringMesh:
00420             return m_covering_source_entities;
00421         case Remapper::DEFAULT:
00422         default:
00423             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_overlap_entities );
00424     }
00425 }
00426 
00427 inline moab::Range& TempestRemapper::GetMeshVertices( Remapper::IntersectionContext ctx )
00428 {
00429     switch( ctx )
00430     {
00431         case Remapper::SourceMesh:
00432             return m_source_vertices;
00433         case Remapper::TargetMesh:
00434             return m_target_vertices;
00435         case Remapper::CoveringMesh:
00436             return m_covering_source_vertices;
00437         case Remapper::DEFAULT:
00438         default:
00439             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_source_vertices );
00440     }
00441 }
00442 
00443 inline const moab::Range& TempestRemapper::GetMeshVertices( Remapper::IntersectionContext ctx ) const
00444 {
00445     switch( ctx )
00446     {
00447         case Remapper::SourceMesh:
00448             return m_source_vertices;
00449         case Remapper::TargetMesh:
00450             return m_target_vertices;
00451         case Remapper::CoveringMesh:
00452             return m_covering_source_vertices;
00453         case Remapper::DEFAULT:
00454         default:
00455             MB_SET_ERR_RET_VAL( "Invalid context passed to GetMeshSet", m_source_vertices );
00456     }
00457 }
00458 
00459 inline void TempestRemapper::SetMeshType( Remapper::IntersectionContext ctx, TempestRemapper::TempestMeshType type )
00460 {
00461     switch( ctx )
00462     {
00463         case Remapper::SourceMesh:
00464             m_source_type = type;
00465             break;
00466         case Remapper::TargetMesh:
00467             m_target_type = type;
00468             break;
00469         case Remapper::OverlapMesh:
00470             m_overlap_type = type;
00471             break;
00472         case Remapper::DEFAULT:
00473         default:
00474             break;
00475     }
00476 }
00477 
00478 inline TempestRemapper::TempestMeshType TempestRemapper::GetMeshType( Remapper::IntersectionContext ctx ) const
00479 {
00480     switch( ctx )
00481     {
00482         case Remapper::SourceMesh:
00483             return m_source_type;
00484         case Remapper::TargetMesh:
00485             return m_target_type;
00486         case Remapper::OverlapMesh:
00487             return m_overlap_type;
00488         case Remapper::DEFAULT:
00489         default:
00490             return TempestRemapper::DEFAULT;
00491     }
00492 }
00493 
00494 inline Mesh* TempestRemapper::GetCoveringMesh()
00495 {
00496     return m_covering_source;
00497 }
00498 
00499 inline moab::EntityHandle& TempestRemapper::GetCoveringSet()
00500 {
00501     return m_covering_source_set;
00502 }
00503 
00504 inline int TempestRemapper::GetGlobalID( Remapper::IntersectionContext ctx, int localID )
00505 {
00506     switch( ctx )
00507     {
00508         case Remapper::SourceMesh:
00509             return lid_to_gid_src[localID];
00510         case Remapper::TargetMesh:
00511             return lid_to_gid_tgt[localID];
00512         case Remapper::CoveringMesh:
00513             return lid_to_gid_covsrc[localID];
00514         case Remapper::OverlapMesh:
00515         case Remapper::DEFAULT:
00516         default:
00517             return -1;
00518     }
00519 }
00520 
00521 inline int TempestRemapper::GetLocalID( Remapper::IntersectionContext ctx, int globalID )
00522 {
00523     switch( ctx )
00524     {
00525         case Remapper::SourceMesh:
00526             return gid_to_lid_src[globalID];
00527         case Remapper::TargetMesh:
00528             return gid_to_lid_tgt[globalID];
00529         case Remapper::CoveringMesh:
00530             return gid_to_lid_covsrc[globalID];
00531         case Remapper::DEFAULT:
00532         case Remapper::OverlapMesh:
00533         default:
00534             return -1;
00535     }
00536 }
00537 
00538 }  // namespace moab
00539 
00540 #endif  // MB_TEMPESTREMAPPER_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines