From 62ab2748216181b2d782274bed9f75b29437785d Mon Sep 17 00:00:00 2001 From: Andrew Nolan Date: Mon, 10 Aug 2026 12:46:50 -0400 Subject: [PATCH 1/7] Add initial work on MOAB mesh initialization --- .../src/drivers/coupled/MoabInterface.cpp | 46 +++++++++++++++++++ .../omega/src/drivers/coupled/MoabInterface.h | 30 ++++++++++++ 2 files changed, 76 insertions(+) create mode 100644 components/omega/src/drivers/coupled/MoabInterface.cpp create mode 100644 components/omega/src/drivers/coupled/MoabInterface.h diff --git a/components/omega/src/drivers/coupled/MoabInterface.cpp b/components/omega/src/drivers/coupled/MoabInterface.cpp new file mode 100644 index 000000000000..a4072e9e05d3 --- /dev/null +++ b/components/omega/src/drivers/coupled/MoabInterface.cpp @@ -0,0 +1,46 @@ +//===-- drivers/coupled/MoabInterface.cpp - MOAB coupling bridge +//-----------===// +// +//===---------------------------------------------------------------------===// + +#include "MoabInterface.h" + +#include "Decomp.h" +#include "HorzMesh.h" +#include "Logging.h" +#include "MachEnv.h" + +#include "moab/iMOAB.h" + +#include + +namespace OMEGA { + +namespace { + +int Pid = -1; + +// Registers this ocean instance as a MOAB application and returns its pid +int registerApplication(MPI_Comm Comm, int OcnID) { + + ErrCode Err; + int LocalPid; + int CompID = OcnID; + + std::string AppName = "OMEGA_MB_" + std::format("{:02}", OcnID); + + Err = iMOAB_RegisterApplication(AppName.c_str(), &Comm, &CompID, &LocalPid); + + if (Err != moab::MB_SUCCESS) + ABORT_ERROR("iMOAB_RegisterApplication: FAIL"); + + return LocalPid; +} +} // namespace + +int moabInit(MPI_Comm Comm, int OcnID) { + Pid = registerApplication(Comm, OcnID); + + return Pid; +} +} // namespace OMEGA diff --git a/components/omega/src/drivers/coupled/MoabInterface.h b/components/omega/src/drivers/coupled/MoabInterface.h new file mode 100644 index 000000000000..d2866dff17bd --- /dev/null +++ b/components/omega/src/drivers/coupled/MoabInterface.h @@ -0,0 +1,30 @@ +#ifndef OMEGA_MOABINTERFACE_H +#define OMEGA_MOABINTERFACE_H +//===-- drivers/coupled/MoabInterface.h - MOAB coupling bridge -*- C++ -*-===// +// +// \file +// \brief +// +// +//===---------------------------------------------------------------------===// + +#include "DataTypes.h" + +#include +#include + +namespace OMEGA { + +int moabInit(MPI_Comm Comm, int OcnID); + +void moabDefineTagStoragee(int Pid, const std::string &Cpl2OcnFieldNames, + const std::string &Ocn2CplFieldNames) + + // Read x2o tag storage (NFields, NCellsOwned) into buffer + void moabImportTagStorgae(Real *Buffer, int Len); + +// Write Buffer (NFields, NCellsOwned) into the o2x tag storage +void moabExportTagStorage(const Real *Buffer, int Len) + +} // namespace OMEGA +#endif // !OMEGA_MOABINTERFACE_H From 68fa1ff277e69dc9b3833d20afe8f299698f16fe Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Sat, 15 Aug 2026 23:42:40 -0700 Subject: [PATCH 2/7] Wire MOAB into Omega's CMake build Detect MOAB via OMEGA_MOAB_ROOT when COMP_INTERFACE=moab, define HAVE_MOAB, and conditionally add/link MoabInterface.cpp. No effect when MOAB is absent. --- components/omega/cime_config/buildlib_cmake | 5 +++++ components/omega/src/CMakeLists.txt | 16 ++++++++++++++++ 2 files changed, 21 insertions(+) diff --git a/components/omega/cime_config/buildlib_cmake b/components/omega/cime_config/buildlib_cmake index 4d2c76655009..8591e76aa578 100755 --- a/components/omega/cime_config/buildlib_cmake +++ b/components/omega/cime_config/buildlib_cmake @@ -31,6 +31,11 @@ def buildlib(bldroot, installpath, case): f' -DOMEGA_PARMETIS_ROOT="{parmetis_root}"' ) + comp_interface = case.get_value("COMP_INTERFACE") + moab_root = os.environ.get("MOAB_ROOT") + if comp_interface == "moab" and moab_root: + cmake_args += f' -DOMEGA_MOAB_ROOT="{moab_root}"' + print(f"omega cmake options: '{cmake_args}'") return cmake_args diff --git a/components/omega/src/CMakeLists.txt b/components/omega/src/CMakeLists.txt index dc7e764da3ab..745e81e334d1 100644 --- a/components/omega/src/CMakeLists.txt +++ b/components/omega/src/CMakeLists.txt @@ -117,6 +117,16 @@ if("${OMEGA_BUILD_MODE}" STREQUAL "E3SM") drivers/coupled/omega_cpl_indices.F90 ) + # MOAB is optional: only present when COMP_INTERFACE=moab (buildlib_cmake + # passes OMEGA_MOAB_ROOT in that case). The non-MOAB build must be + # unaffected when it isn't found. + if(OMEGA_MOAB_ROOT) + find_package(MOAB CONFIG PATHS ${OMEGA_MOAB_ROOT} NO_DEFAULT_PATH) + endif() + if(MOAB_FOUND) + list(APPEND OCN_SRC drivers/coupled/MoabInterface.cpp) + endif() + # create ocn lib add_library(ocn ${OCN_SRC}) @@ -138,6 +148,12 @@ if("${OMEGA_BUILD_MODE}" STREQUAL "E3SM") csm_share ) + if(MOAB_FOUND) + target_compile_definitions(ocn PRIVATE HAVE_MOAB) + target_include_directories(ocn PRIVATE ${MOAB_INCLUDE_DIRS}) + target_link_libraries(ocn PRIVATE ${MOAB_LIBRARIES}) + endif() + endif() # build Omega executable From 40b561ca448e793dc03ae063560add64f8ab1556 Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Sat, 15 Aug 2026 23:43:05 -0700 Subject: [PATCH 3/7] Implement MOAB mesh and tag-storage builder registerApplication only registered an app. Now builds an owned-cells mesh (vertex compaction, polygon elements, GLOBAL_ID tags), domain tags, and x2o/o2x tag get/set. Also fixes header typos that blocked compilation. --- .../src/drivers/coupled/MoabInterface.cpp | 268 +++++++++++++++++- .../omega/src/drivers/coupled/MoabInterface.h | 25 +- 2 files changed, 282 insertions(+), 11 deletions(-) diff --git a/components/omega/src/drivers/coupled/MoabInterface.cpp b/components/omega/src/drivers/coupled/MoabInterface.cpp index a4072e9e05d3..1455de2fe52f 100644 --- a/components/omega/src/drivers/coupled/MoabInterface.cpp +++ b/components/omega/src/drivers/coupled/MoabInterface.cpp @@ -6,12 +6,18 @@ #include "MoabInterface.h" #include "Decomp.h" +#include "GlobalConstants.h" #include "HorzMesh.h" #include "Logging.h" #include "MachEnv.h" +#include "OmegaKokkos.h" #include "moab/iMOAB.h" +#include +#include +#include +#include #include namespace OMEGA { @@ -20,6 +26,62 @@ namespace { int Pid = -1; +// Full CIME x2o/o2x field-name lists (colon-separated, null-terminated), +// captured by moabDefineTagStorage and reused by moabImportTagStorage/ +// moabExportTagStorage. +std::string X2oTagNames; +std::string O2xTagNames; + +void checkMoabErr(ErrCode Err, const std::string &Msg) { + if (Err != moab::MB_SUCCESS) + LOG_CRITICAL("MOAB error in {}: code {}", Msg, static_cast(Err)); +} + +// iMOAB's tag storage API is double-only. These are templates (not plain +// `if constexpr` in a non-template function) specifically so the untaken +// branch is genuinely discarded-and-not-instantiated for the concrete T, +// rather than needing to type-check for every T in every build: when T is +// double (the common case) the buffer is read/written directly with no +// intermediate copy; only a build where Real != double needs the scratch +// conversion. +template +void getDoubleTagStorage(const std::string &TagNames, T *Buffer, int Len) { + int EntType = 1; // elements + if constexpr (std::is_same_v) { + ErrCode Err = iMOAB_GetDoubleTagStorage(&Pid, TagNames.c_str(), &Len, + &EntType, Buffer); + checkMoabErr(Err, "iMOAB_GetDoubleTagStorage"); + } else { + std::vector Scratch(Len); + ErrCode Err = iMOAB_GetDoubleTagStorage(&Pid, TagNames.c_str(), &Len, + &EntType, Scratch.data()); + checkMoabErr(Err, "iMOAB_GetDoubleTagStorage"); + for (int I = 0; I < Len; ++I) + Buffer[I] = static_cast(Scratch[I]); + } +} + +template +void setDoubleTagStorage(const std::string &TagNames, const T *Buffer, + int Len) { + int EntType = 1; // elements + if constexpr (std::is_same_v) { + // iMOAB_SetDoubleTagStorage takes a non-const double* even though it + // only reads from it; the const_cast is safe, not a real mutation. + ErrCode Err = iMOAB_SetDoubleTagStorage(&Pid, TagNames.c_str(), &Len, + &EntType, + const_cast(Buffer)); + checkMoabErr(Err, "iMOAB_SetDoubleTagStorage"); + } else { + std::vector Scratch(Len); + for (int I = 0; I < Len; ++I) + Scratch[I] = static_cast(Buffer[I]); + ErrCode Err = iMOAB_SetDoubleTagStorage(&Pid, TagNames.c_str(), &Len, + &EntType, Scratch.data()); + checkMoabErr(Err, "iMOAB_SetDoubleTagStorage"); + } +} + // Registers this ocean instance as a MOAB application and returns its pid int registerApplication(MPI_Comm Comm, int OcnID) { @@ -27,20 +89,218 @@ int registerApplication(MPI_Comm Comm, int OcnID) { int LocalPid; int CompID = OcnID; - std::string AppName = "OMEGA_MB_" + std::format("{:02}", OcnID); + char OcnIDStr[8]; + std::snprintf(OcnIDStr, sizeof(OcnIDStr), "%02d", OcnID); + std::string AppName = std::string("OMEGA_MB_") + OcnIDStr; Err = iMOAB_RegisterApplication(AppName.c_str(), &Comm, &CompID, &LocalPid); - - if (Err != moab::MB_SUCCESS) - ABORT_ERROR("iMOAB_RegisterApplication: FAIL"); + checkMoabErr(Err, "iMOAB_RegisterApplication"); return LocalPid; } + +// Defines a dense double tag with 1 component per entity on elements and +// sets it from Data (size NCellsOwned). +void defineAndSetElemDoubleTag(int LocalPid, const std::string &TagName, + int NCellsOwned, const double *Data) { + ErrCode Err; + int TagType = DENSE_DOUBLE; + int NumCo = 1; + int TagIndex; + + Err = iMOAB_DefineTagStorage(&LocalPid, TagName.c_str(), &TagType, &NumCo, + &TagIndex); + checkMoabErr(Err, "iMOAB_DefineTagStorage(" + TagName + ")"); + + int EntType = 1; // elements + int Len = NCellsOwned; + Err = iMOAB_SetDoubleTagStorage(&LocalPid, TagName.c_str(), &Len, &EntType, + const_cast(Data)); + checkMoabErr(Err, "iMOAB_SetDoubleTagStorage(" + TagName + ")"); +} + +// Builds the MOAB mesh for this ocean instance from the default +// Decomp/HorzMesh, scoped to owned cells only (matches MPAS-Ocean's +// init_moab_mpas, which never registers halo elements or ghost layers for +// its component-side MOAB application). +void createMesh(int LocalPid, int OcnID) { + + Decomp *DefDecomp = Decomp::getDefault(); + HorzMesh *DefMesh = HorzMesh::getDefault(); + const I4 NCellsOwned = DefDecomp->NCellsOwned; + const I4 MaxEdges = DefDecomp->MaxEdges; + + deepCopy(DefDecomp->CellIDH, DefDecomp->CellID); + deepCopy(DefDecomp->VertexIDH, DefDecomp->VertexID); + deepCopy(DefDecomp->VerticesOnCellH, DefDecomp->VerticesOnCell); + deepCopy(DefDecomp->NEdgesOnCellH, DefDecomp->NEdgesOnCell); + deepCopy(DefMesh->XVertexH, DefMesh->XVertex); + deepCopy(DefMesh->YVertexH, DefMesh->YVertex); + deepCopy(DefMesh->ZVertexH, DefMesh->ZVertex); + + // Single-pass compaction: remap Decomp's local vertex indices (which + // range over owned+halo vertices) to a compact, first-encountered + // numbering covering only the vertices actually referenced by owned + // cells, and build the (NCellsOwned, MaxEdges) connectivity array in + // that compact numbering. Cells with fewer than MaxEdges sides are + // padded by repeating the cell's last real vertex, matching + // mpas_moabmesh.F's convention. + std::unordered_map CompactIndex; + std::vector CompactToLocal; + std::vector Connectivity(static_cast(NCellsOwned) * MaxEdges); + + I4 Offset = 0; + for (I4 Cell = 0; Cell < NCellsOwned; ++Cell) { + const I4 NEdges = DefDecomp->NEdgesOnCellH(Cell); + I4 LastCompact = -1; + for (I4 Edge = 0; Edge < NEdges; ++Edge) { + const I4 LocalVertex = DefDecomp->VerticesOnCellH(Cell, Edge); + auto It = CompactIndex.find(LocalVertex); + I4 Compact; + if (It == CompactIndex.end()) { + Compact = static_cast(CompactToLocal.size()); + CompactIndex.emplace(LocalVertex, Compact); + CompactToLocal.push_back(LocalVertex); + } else { + Compact = It->second; + } + Connectivity[Offset + Edge] = Compact; + LastCompact = Compact; + } + for (I4 Edge = NEdges; Edge < MaxEdges; ++Edge) + Connectivity[Offset + Edge] = LastCompact; + Offset += MaxEdges; + } + + int NCompactVerts = static_cast(CompactToLocal.size()); + + // Vertex coordinates, unit-sphere Cartesian, in compact order. + std::vector Coords(3 * static_cast(NCompactVerts)); + for (int V = 0; V < NCompactVerts; ++V) { + const I4 LocalVertex = CompactToLocal[V]; + Coords[3 * V + 0] = DefMesh->XVertexH[LocalVertex] / DefMesh->SphereRadius; + Coords[3 * V + 1] = DefMesh->YVertexH[LocalVertex] / DefMesh->SphereRadius; + Coords[3 * V + 2] = DefMesh->ZVertexH[LocalVertex] / DefMesh->SphereRadius; + } + + ErrCode Err; + int CoordsLen = 3 * NCompactVerts; + int Dim = 3; + Err = iMOAB_CreateVertices(&LocalPid, &CoordsLen, &Dim, Coords.data()); + checkMoabErr(Err, "iMOAB_CreateVertices"); + + int NElem = static_cast(NCellsOwned); + int MBType = 4; // MBPOLYGON + int NNodesPerElem = static_cast(MaxEdges); + int BlockID = 100 * OcnID + MachEnv::getDefault()->getMyTask(); + Err = iMOAB_CreateElements(&LocalPid, &NElem, &MBType, &NNodesPerElem, + Connectivity.data(), &BlockID); + checkMoabErr(Err, "iMOAB_CreateElements"); + + // GLOBAL_ID tags on vertices and elements, needed both for + // ResolveSharedEntities and so the coupler-side offline weight file's + // row/column numbering lines up with this mesh. + std::string TagName = "GLOBAL_ID"; + int TagType = DENSE_INTEGER; + int NumCo = 1; + int TagIndex; + Err = iMOAB_DefineTagStorage(&LocalPid, TagName.c_str(), &TagType, &NumCo, + &TagIndex); + checkMoabErr(Err, "iMOAB_DefineTagStorage(GLOBAL_ID)"); + + std::vector VertexGlobalIDs(NCompactVerts); + for (int V = 0; V < NCompactVerts; ++V) + VertexGlobalIDs[V] = DefDecomp->VertexIDH[CompactToLocal[V]]; + + int VertEntType = 0; // vertices + Err = iMOAB_SetIntTagStorage(&LocalPid, TagName.c_str(), &NCompactVerts, + &VertEntType, VertexGlobalIDs.data()); + checkMoabErr(Err, "iMOAB_SetIntTagStorage(GLOBAL_ID, vertices)"); + + std::vector CellGlobalIDs(NElem); + for (I4 Cell = 0; Cell < NCellsOwned; ++Cell) + CellGlobalIDs[Cell] = DefDecomp->CellIDH[Cell]; + + int ElemEntType = 1; // elements + Err = iMOAB_SetIntTagStorage(&LocalPid, TagName.c_str(), &NElem, + &ElemEntType, CellGlobalIDs.data()); + checkMoabErr(Err, "iMOAB_SetIntTagStorage(GLOBAL_ID, elements)"); + + Err = iMOAB_ResolveSharedEntities(&LocalPid, &NCompactVerts, + VertexGlobalIDs.data()); + checkMoabErr(Err, "iMOAB_ResolveSharedEntities"); + + Err = iMOAB_UpdateMeshInfo(&LocalPid); + checkMoabErr(Err, "iMOAB_UpdateMeshInfo"); +} + +// Defines and sets the domain tags (lon/lat/area/mask/frac) that the +// coupler-side mapper reads for the offline weight file's area-weighted +// interpolation. Values are computed the same way as the MCT path's +// ocn_set_domain_mct (via omega_get_lonlat_cell/omega_get_area_cell in +// omega_cxx2f_interface.cpp) so both drivers see identical numbers. +void setDomainTags(int LocalPid) { + + HorzMesh *DefMesh = HorzMesh::getDefault(); + const I4 NCellsOwned = DefMesh->NCellsOwned; + + deepCopy(DefMesh->LonCellH, DefMesh->LonCell); + deepCopy(DefMesh->LatCellH, DefMesh->LatCell); + deepCopy(DefMesh->AreaCellH, DefMesh->AreaCell); + + const Real SphereRadius2 = DefMesh->SphereRadius * DefMesh->SphereRadius; + + std::vector Lon(NCellsOwned), Lat(NCellsOwned), Area(NCellsOwned); + std::vector MaskFrac(NCellsOwned, 1.0); + for (I4 Cell = 0; Cell < NCellsOwned; ++Cell) { + Lon[Cell] = static_cast(DefMesh->LonCellH[Cell] * Rad2Deg); + Lat[Cell] = static_cast(DefMesh->LatCellH[Cell] * Rad2Deg); + Area[Cell] = static_cast(DefMesh->AreaCellH[Cell] / SphereRadius2); + } + + defineAndSetElemDoubleTag(LocalPid, "lon", NCellsOwned, Lon.data()); + defineAndSetElemDoubleTag(LocalPid, "lat", NCellsOwned, Lat.data()); + defineAndSetElemDoubleTag(LocalPid, "area", NCellsOwned, Area.data()); + defineAndSetElemDoubleTag(LocalPid, "mask", NCellsOwned, MaskFrac.data()); + defineAndSetElemDoubleTag(LocalPid, "frac", NCellsOwned, MaskFrac.data()); +} + } // namespace int moabInit(MPI_Comm Comm, int OcnID) { Pid = registerApplication(Comm, OcnID); + createMesh(Pid, OcnID); + setDomainTags(Pid); return Pid; } + +void moabDefineTagStorage(int LocalPid, const std::string &Cpl2OcnFieldNames, + const std::string &Ocn2CplFieldNames) { + + X2oTagNames = Cpl2OcnFieldNames; + O2xTagNames = Ocn2CplFieldNames; + + ErrCode Err; + int TagType = DENSE_DOUBLE; + int NumCo = 1; + int TagIndex; + + Err = iMOAB_DefineTagStorage(&LocalPid, X2oTagNames.c_str(), &TagType, + &NumCo, &TagIndex); + checkMoabErr(Err, "iMOAB_DefineTagStorage(x2o fields)"); + + Err = iMOAB_DefineTagStorage(&LocalPid, O2xTagNames.c_str(), &TagType, + &NumCo, &TagIndex); + checkMoabErr(Err, "iMOAB_DefineTagStorage(o2x fields)"); +} + +void moabImportTagStorage(Real *Buffer, int Len) { + getDoubleTagStorage(X2oTagNames, Buffer, Len); +} + +void moabExportTagStorage(const Real *Buffer, int Len) { + setDoubleTagStorage(O2xTagNames, Buffer, Len); +} + } // namespace OMEGA diff --git a/components/omega/src/drivers/coupled/MoabInterface.h b/components/omega/src/drivers/coupled/MoabInterface.h index d2866dff17bd..720861f623ae 100644 --- a/components/omega/src/drivers/coupled/MoabInterface.h +++ b/components/omega/src/drivers/coupled/MoabInterface.h @@ -3,8 +3,10 @@ //===-- drivers/coupled/MoabInterface.h - MOAB coupling bridge -*- C++ -*-===// // // \file -// \brief +// \brief Bridge between Omega's ocean core and the MOAB coupling driver. // +// MoabInterface.h/.cpp are the only Omega source files that include +// moab/iMOAB.h. // //===---------------------------------------------------------------------===// @@ -15,16 +17,25 @@ namespace OMEGA { +// Registers this ocean instance as a MOAB application, builds its mesh from +// the default Decomp/HorzMesh (owned cells only), and defines the domain +// tags (GLOBAL_ID, lon, lat, area, mask, frac). Returns the MOAB +// application id. int moabInit(MPI_Comm Comm, int OcnID); -void moabDefineTagStoragee(int Pid, const std::string &Cpl2OcnFieldNames, - const std::string &Ocn2CplFieldNames) +// Defines dense, double, single-component tags for the coupler's full +// x2o/o2x field lists (colon-separated, null-terminated) on the ocean's +// MOAB application. +void moabDefineTagStorage(int Pid, const std::string &Cpl2OcnFieldNames, + const std::string &Ocn2CplFieldNames); - // Read x2o tag storage (NFields, NCellsOwned) into buffer - void moabImportTagStorgae(Real *Buffer, int Len); +// Reads the x2o tag storage (NFields, NCellsOwned, unrolled by tag) into +// Buffer. +void moabImportTagStorage(Real *Buffer, int Len); -// Write Buffer (NFields, NCellsOwned) into the o2x tag storage -void moabExportTagStorage(const Real *Buffer, int Len) +// Writes Buffer (NFields, NCellsOwned, unrolled by tag) into the o2x tag +// storage. +void moabExportTagStorage(const Real *Buffer, int Len); } // namespace OMEGA #endif // !OMEGA_MOABINTERFACE_H From b853cae8f006be65e4b7e7bc0847006605347a35 Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Sat, 15 Aug 2026 23:43:26 -0700 Subject: [PATCH 4/7] Wire Omega's ocean cap to MOAB through the coupled bridge Pass full CIME field-name lists to the C++ side, pick MOAB vs MCT layout at compile time, and call the mesh/tag-storage functions around ocean init and run. Drops the unused runtime layout selection it replaces. --- .../src/drivers/coupled/ocn_comp_mct.F90 | 22 +++--- .../src/drivers/coupled/omega_cpl_indices.F90 | 12 ++++ .../drivers/coupled/omega_cxx2f_interface.cpp | 68 ++++++++++++++++--- .../drivers/coupled/omega_f2cxx_interface.F90 | 28 ++------ 4 files changed, 86 insertions(+), 44 deletions(-) diff --git a/components/omega/src/drivers/coupled/ocn_comp_mct.F90 b/components/omega/src/drivers/coupled/ocn_comp_mct.F90 index 7f5e5e7fc21a..0d85e810ac8f 100644 --- a/components/omega/src/drivers/coupled/ocn_comp_mct.F90 +++ b/components/omega/src/drivers/coupled/ocn_comp_mct.F90 @@ -45,9 +45,7 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) use omega_f2cxx_mod, only: & omega_ocn_init1, & - omega_ocn_init2, & - omega_get_layout_mct, & - omega_get_layout_moab + omega_ocn_init2 use omega_cpl_indices, only: & num_coupler_imports, & @@ -56,6 +54,8 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) export_field_names, & import_field_indices, & export_field_indices, & + cpl_x2o_field_names, & + cpl_o2x_field_names, & omega_set_cpl_indices use mct_mod, only: mct_gsMap_lsize @@ -96,7 +96,6 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) integer(IN) :: & coupling_time_step, case_start_tod, case_start_ymd, cur_tod, cur_ymd integer(kind=c_int) :: start_type_c - integer(kind=c_int) :: layout character(kind=c_char, len=CL), target :: calendar_c character(kind=c_char, len=CL), target :: ocn_log_fname_c @@ -186,12 +185,6 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) ! populate the import/export field name and index arrays call omega_set_cpl_indices() -#ifdef HAVE_MOAB - layout = omega_get_layout_moab() -#else - layout = omega_get_layout_mct() -#endif - call omega_ocn_init1( & mpicom_ocn, & OCN_ID, & @@ -209,7 +202,9 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) c_loc(import_field_names), & c_loc(export_field_names), & c_loc(import_field_indices), & - c_loc(export_field_indices) & + c_loc(export_field_indices), & + c_loc(cpl_x2o_field_names), & + c_loc(cpl_o2x_field_names) & ) !------------------------------------------------------------------------- @@ -236,7 +231,10 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) ! TODO: Get case config info and add as MetaData to Omega - ! TODO: ifdef HAVE_MOAB + ! Under HAVE_MOAB, omega_ocn_init2 ignores these MCT attribute-vector + ! pointers and attaches its own MOAB-backed buffers instead (see + ! omega_cxx2f_interface.cpp); they're still passed here unconditionally + ! since x2o/o2x are always allocated above regardless of driver. call omega_ocn_init2(c_loc(x2o%rAttr), c_loc(o2x%rAttr)) end subroutine ocn_init_mct diff --git a/components/omega/src/drivers/coupled/omega_cpl_indices.F90 b/components/omega/src/drivers/coupled/omega_cpl_indices.F90 index 0b22be5918c6..c591b762daf9 100644 --- a/components/omega/src/drivers/coupled/omega_cpl_indices.F90 +++ b/components/omega/src/drivers/coupled/omega_cpl_indices.F90 @@ -19,12 +19,20 @@ module omega_cpl_indices import_field_indices(num_omega_imports), & export_field_indices(num_omega_exports) + ! Full CIME x2o/o2x field-name lists (colon-separated), needed by the + ! MOAB bridge to define coupler tags covering every field CIME expects + ! on this app, not just the fields Omega currently imports/exports. + character(len=:, kind=c_char), public, allocatable, target :: & + cpl_x2o_field_names, & + cpl_o2x_field_names + public :: omega_set_cpl_indices contains subroutine omega_set_cpl_indices() + use, intrinsic :: iso_c_binding, only: c_null_char use mct_mod, only: mct_aVect, mct_aVect_init, mct_aVect_clean use seq_flds_mod, only: seq_flds_o2x_fields, seq_flds_x2o_fields @@ -36,6 +44,10 @@ subroutine omega_set_cpl_indices() call mct_aVect_init(x2o, rList=seq_flds_x2o_fields, lsize=1) call mct_aVect_init(o2x, rList=seq_flds_o2x_fields, lsize=1) + ! full CIME field lists, for the MOAB bridge's tag definitions + cpl_x2o_field_names = trim(seq_flds_x2o_fields)//c_null_char + cpl_o2x_field_names = trim(seq_flds_o2x_fields)//c_null_char + ! total number of import/export fields in coupler data arrays num_coupler_imports = size(x2o%rAttr, 1) num_coupler_exports = size(o2x%rAttr, 1) diff --git a/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp b/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp index 2347521336f6..5d80e7986600 100644 --- a/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp +++ b/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp @@ -3,6 +3,7 @@ // //===----------------------------------------------------------------------===// #include "DataTypes.h" +#include "Decomp.h" #include "Logging.h" #include "MachEnv.h" #include "OceanDriver.h" @@ -12,6 +13,11 @@ #include "TimeMgr.h" #include "TimeStepper.h" #include +#include + +#ifdef HAVE_MOAB +#include "MoabInterface.h" +#endif // helper C++ functions namespace { @@ -39,6 +45,17 @@ std::map buildFieldIndexMap(const char *FieldNames, } return FieldIdx; } + +#ifdef HAVE_MOAB +// Coupling buffers for the MOAB path: attached once (by address) to +// SfcCoupling in omega_ocn_init2, refilled from/drained to MOAB tag storage +// around each omega_ocn_run call. +int MoabNCouplerImports = 0; +int MoabNCouplerExports = 0; +std::vector MoabCplToOcn; +std::vector MoabOcnToCpl; +#endif + } // namespace extern "C" { @@ -60,7 +77,9 @@ void omega_ocn_init1( const char *ImportFieldNames, // [in] array of import field names const char *ExportFieldNames, // [in] array of export field names const int *ImportFieldIndices, // [in] array of import field indices - const int *ExportFieldIndices // [in] array of export field indices + const int *ExportFieldIndices, // [in] array of export field indices + const char *Cpl2OcnFieldNames, // [in] full CIME x2o field list (MOAB) + const char *Ocn2CplFieldNames // [in] full CIME o2x field list (MOAB) ) { // Create the C MPI_Comm from the Fortran one @@ -101,7 +120,12 @@ void omega_ocn_init1( OMEGA::TimeInitParams TimeParams{StartTime, std::nullopt}; OMEGA::CouplingInitParams CouplingParams{ NCouplerImports, NCouplerExports, ImportIdxMap, - ExportIdxMap, CouplingInterval, OMEGA::CouplingLayout::MCT}; + ExportIdxMap, CouplingInterval, +#ifdef HAVE_MOAB + OMEGA::CouplingLayout::MOAB}; +#else + OMEGA::CouplingLayout::MCT}; +#endif OMEGA::ocnInit1(Comm, OcnID, YamlConfigFile, OcnLogFile, StartTypeEnum, TimeParams, CouplingParams); @@ -110,11 +134,35 @@ void omega_ocn_init1( LOG_INFO("ocnInit: Finished initializing ocean model"); int ErrAll; + +#ifdef HAVE_MOAB + // Decomp/HorzMesh exist by now (built inside OMEGA::ocnInit1 above), so + // the MOAB mesh can be constructed from them. + MoabNCouplerImports = NCouplerImports; + MoabNCouplerExports = NCouplerExports; + int MoabPid = OMEGA::moabInit(Comm, OcnID); + OMEGA::moabDefineTagStorage(MoabPid, Cpl2OcnFieldNames, Ocn2CplFieldNames); +#endif } void omega_ocn_init2(const double *cpl_to_ocn_data, double *ocn_to_cpl_data) { Pacer::start("Init2", 0); +#ifdef HAVE_MOAB + // The MCT attribute-vector pointers above are meaningless under MOAB; + // attach Omega's own buffers instead, filled from MOAB tag storage. + // These are refilled/drained in place (same memory address) by + // omega_ocn_run on every subsequent coupling interval. + const int NCellsOwned = static_cast(OMEGA::Decomp::getDefault()->NCellsOwned); + MoabCplToOcn.assign(static_cast(MoabNCouplerImports) * NCellsOwned, + 0); + MoabOcnToCpl.assign(static_cast(MoabNCouplerExports) * NCellsOwned, + 0); + OMEGA::moabImportTagStorage(MoabCplToOcn.data(), + static_cast(MoabCplToOcn.size())); + OMEGA::ocnInit2(MoabCplToOcn.data(), MoabOcnToCpl.data()); +#else OMEGA::ocnInit2(cpl_to_ocn_data, ocn_to_cpl_data); +#endif Pacer::stop("Init2", 0); } @@ -127,7 +175,15 @@ int omega_ocn_run(bool WriteRestart) { OMEGA::TimeInstant CurrTime = ModelClock->getCurrentTime(); Pacer::start("Run", 0); +#ifdef HAVE_MOAB + OMEGA::moabImportTagStorage(MoabCplToOcn.data(), + static_cast(MoabCplToOcn.size())); +#endif ErrRun = OMEGA::ocnRun(CurrTime, WriteRestart); +#ifdef HAVE_MOAB + OMEGA::moabExportTagStorage(MoabOcnToCpl.data(), + static_cast(MoabOcnToCpl.size())); +#endif Pacer::stop("Run", 0); return ErrRun; @@ -158,14 +214,6 @@ int omega_ocn_finalize() { return ErrFinalize; } -int omega_get_layout_mct() { - return static_cast(OMEGA::CouplingLayout::MCT); -} - -int omega_get_layout_moab() { - return static_cast(OMEGA::CouplingLayout::MOAB); -} - int omega_get_ncells_local() { OMEGA::Decomp *OcnDecomp = OMEGA::Decomp::getDefault(); diff --git a/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 b/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 index 7f00263fffa7..a998a4138a3d 100644 --- a/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 +++ b/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 @@ -27,7 +27,9 @@ subroutine omega_ocn_init1( & import_field_names, & export_field_names, & import_field_indices, & - export_field_indices) bind(c) + export_field_indices, & + cpl_x2o_field_names, & + cpl_o2x_field_names) bind(c) use, intrinsic :: iso_c_binding, only: c_int, c_char, c_ptr @@ -52,7 +54,9 @@ subroutine omega_ocn_init1( & import_field_names, & export_field_names, & import_field_indices, & - export_field_indices + export_field_indices, & + cpl_x2o_field_names, & + cpl_o2x_field_names end subroutine omega_ocn_init1 @@ -84,26 +88,6 @@ subroutine omega_ocn_finalize() bind(c) end subroutine omega_ocn_finalize - function omega_get_layout_mct() result(layout_mct) bind(c) - - use, intrinsic :: iso_c_binding, only: c_int - - implicit none - - integer(kind=c_int) :: layout_mct - - end function omega_get_layout_mct - - function omega_get_layout_moab() result(layout_moab) bind(c) - - use, intrinsic :: iso_c_binding, only: c_int - - implicit none - - integer(kind=c_int) :: layout_moab - - end function omega_get_layout_moab - function omega_get_ncells_local() result(ncells_local) bind(c) use, intrinsic :: iso_c_binding, only: c_int From 07cf0c8fc10a2c9178c106e0d8ba984c9e206530 Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Sun, 16 Aug 2026 19:09:24 -0700 Subject: [PATCH 5/7] Set mpoid so the coupler can find the ocean's MOAB mesh moabInit's app id was staying local to the C++ side. The coupler-side migration/mapping code reads seq_comm_mct's mpoid to locate the ocean mesh, so leaving it at -1 meant it never found it, causing the crash and memory corruption seen at runtime. Fortran now copies the id into mpoid. --- .../omega/src/drivers/coupled/ocn_comp_mct.F90 | 10 ++++++++++ .../src/drivers/coupled/omega_cxx2f_interface.cpp | 12 +++++++++++- .../src/drivers/coupled/omega_f2cxx_interface.F90 | 10 ++++++++++ 3 files changed, 31 insertions(+), 1 deletion(-) diff --git a/components/omega/src/drivers/coupled/ocn_comp_mct.F90 b/components/omega/src/drivers/coupled/ocn_comp_mct.F90 index 0d85e810ac8f..dcf009836078 100644 --- a/components/omega/src/drivers/coupled/ocn_comp_mct.F90 +++ b/components/omega/src/drivers/coupled/ocn_comp_mct.F90 @@ -46,6 +46,10 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) use omega_f2cxx_mod, only: & omega_ocn_init1, & omega_ocn_init2 +#ifdef HAVE_MOAB + use omega_f2cxx_mod, only: omega_get_moab_pid + use seq_comm_mct, only: mpoid +#endif use omega_cpl_indices, only: & num_coupler_imports, & @@ -207,6 +211,12 @@ subroutine ocn_init_mct(EClock, cdata, x2o, o2x, NLFilename) c_loc(cpl_o2x_field_names) & ) +#ifdef HAVE_MOAB + ! tell the coupler-side migration/mapping code (cplcomp_exchange_mod, + ! prep_ocn_mod) which MOAB app id is this ocean instance's own mesh + mpoid = omega_get_moab_pid() +#endif + !------------------------------------------------------------------------- ! initialize MCT gsmap, domain, and attribute vectors !------------------------------------------------------------------------- diff --git a/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp b/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp index 5d80e7986600..ab25ac1e57b1 100644 --- a/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp +++ b/components/omega/src/drivers/coupled/omega_cxx2f_interface.cpp @@ -54,6 +54,12 @@ int MoabNCouplerImports = 0; int MoabNCouplerExports = 0; std::vector MoabCplToOcn; std::vector MoabOcnToCpl; + +// The ocean's own MOAB application id, returned by moabInit. The Fortran +// cap must copy this into seq_comm_mct's mpoid module variable, since that +// is what the coupler-side migration/mapping code reads to find the +// ocean's mesh. +int MoabPid = -1; #endif } // namespace @@ -140,11 +146,15 @@ void omega_ocn_init1( // the MOAB mesh can be constructed from them. MoabNCouplerImports = NCouplerImports; MoabNCouplerExports = NCouplerExports; - int MoabPid = OMEGA::moabInit(Comm, OcnID); + MoabPid = OMEGA::moabInit(Comm, OcnID); OMEGA::moabDefineTagStorage(MoabPid, Cpl2OcnFieldNames, Ocn2CplFieldNames); #endif } +#ifdef HAVE_MOAB +int omega_get_moab_pid() { return MoabPid; } +#endif + void omega_ocn_init2(const double *cpl_to_ocn_data, double *ocn_to_cpl_data) { Pacer::start("Init2", 0); #ifdef HAVE_MOAB diff --git a/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 b/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 index a998a4138a3d..42c7554f0219 100644 --- a/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 +++ b/components/omega/src/drivers/coupled/omega_f2cxx_interface.F90 @@ -88,6 +88,16 @@ subroutine omega_ocn_finalize() bind(c) end subroutine omega_ocn_finalize + function omega_get_moab_pid() result(moab_pid) bind(c) + + use, intrinsic :: iso_c_binding, only: c_int + + implicit none + + integer(kind=c_int) :: moab_pid + + end function omega_get_moab_pid + function omega_get_ncells_local() result(ncells_local) bind(c) use, intrinsic :: iso_c_binding, only: c_int From f6d95ad227c9acd11af7e12093f1dcbede63005f Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Sun, 16 Aug 2026 19:40:19 -0700 Subject: [PATCH 6/7] Define an aream domain tag alongside area/mask component_init_areacor_moab reads area:aream:mask as one combined tag; a missing tag fails the whole read. aream gets a -9999 placeholder, same convention MCT uses, since the coupler fills in the real value later. --- .../omega/src/drivers/coupled/MoabInterface.cpp | 11 +++++++++-- 1 file changed, 9 insertions(+), 2 deletions(-) diff --git a/components/omega/src/drivers/coupled/MoabInterface.cpp b/components/omega/src/drivers/coupled/MoabInterface.cpp index 1455de2fe52f..b3d2f798e8ec 100644 --- a/components/omega/src/drivers/coupled/MoabInterface.cpp +++ b/components/omega/src/drivers/coupled/MoabInterface.cpp @@ -234,11 +234,16 @@ void createMesh(int LocalPid, int OcnID) { checkMoabErr(Err, "iMOAB_UpdateMeshInfo"); } -// Defines and sets the domain tags (lon/lat/area/mask/frac) that the +// Defines and sets the domain tags (lon/lat/area/aream/mask/frac) that the // coupler-side mapper reads for the offline weight file's area-weighted // interpolation. Values are computed the same way as the MCT path's // ocn_set_domain_mct (via omega_get_lonlat_cell/omega_get_area_cell in -// omega_cxx2f_interface.cpp) so both drivers see identical numbers. +// omega_cxx2f_interface.cpp) so both drivers see identical numbers. aream +// is given a placeholder, same as ocn_set_domain_mct does for MCT: the +// coupler computes the real value from the mapping file and pushes it back +// down onto this tag (component_init_areacor_moab's 'x2c' exchange), but +// the tag must already exist here for that exchange to have somewhere to +// write it, and for its combined "area:aream:mask" read to succeed. void setDomainTags(int LocalPid) { HorzMesh *DefMesh = HorzMesh::getDefault(); @@ -252,6 +257,7 @@ void setDomainTags(int LocalPid) { std::vector Lon(NCellsOwned), Lat(NCellsOwned), Area(NCellsOwned); std::vector MaskFrac(NCellsOwned, 1.0); + std::vector AreaM(NCellsOwned, -9999.0); for (I4 Cell = 0; Cell < NCellsOwned; ++Cell) { Lon[Cell] = static_cast(DefMesh->LonCellH[Cell] * Rad2Deg); Lat[Cell] = static_cast(DefMesh->LatCellH[Cell] * Rad2Deg); @@ -261,6 +267,7 @@ void setDomainTags(int LocalPid) { defineAndSetElemDoubleTag(LocalPid, "lon", NCellsOwned, Lon.data()); defineAndSetElemDoubleTag(LocalPid, "lat", NCellsOwned, Lat.data()); defineAndSetElemDoubleTag(LocalPid, "area", NCellsOwned, Area.data()); + defineAndSetElemDoubleTag(LocalPid, "aream", NCellsOwned, AreaM.data()); defineAndSetElemDoubleTag(LocalPid, "mask", NCellsOwned, MaskFrac.data()); defineAndSetElemDoubleTag(LocalPid, "frac", NCellsOwned, MaskFrac.data()); } From abdb6e1089e65f52acbdb57702dfef916f2397cc Mon Sep 17 00:00:00 2001 From: Vijay Mahadevan Date: Mon, 17 Aug 2026 08:49:31 -0700 Subject: [PATCH 7/7] Fix off-by-one in polygon connectivity indices iMOAB_CreateElements expects connectivity entries to be 1-based (it computes firstVertex + connectivity[j] - 1). Compact vertex indices were written 0-based, so every vertex reference pointed one slot too early. --- components/omega/src/drivers/coupled/MoabInterface.cpp | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/components/omega/src/drivers/coupled/MoabInterface.cpp b/components/omega/src/drivers/coupled/MoabInterface.cpp index b3d2f798e8ec..b582fae5b78f 100644 --- a/components/omega/src/drivers/coupled/MoabInterface.cpp +++ b/components/omega/src/drivers/coupled/MoabInterface.cpp @@ -164,11 +164,15 @@ void createMesh(int LocalPid, int OcnID) { } else { Compact = It->second; } - Connectivity[Offset + Edge] = Compact; + // iMOAB_CreateElements indexes into the just-created vertex range + // as connectivity[j] + firstVertex - 1, i.e. it expects 1-based + // indices; Compact itself stays 0-based everywhere else (it also + // indexes Coords/CompactToLocal directly). + Connectivity[Offset + Edge] = Compact + 1; LastCompact = Compact; } for (I4 Edge = NEdges; Edge < MaxEdges; ++Edge) - Connectivity[Offset + Edge] = LastCompact; + Connectivity[Offset + Edge] = LastCompact + 1; Offset += MaxEdges; }