Various fixes

src/pcms/coupler/field_communicator.hpp

@@ -35,7 +35,11 @@ class FieldCommunicator
       serializer_(std::move(serializer))
   {
     PCMS_ALWAYS_ASSERT(&layout_comm.GetLayout() == &field.GetLayout());
-    comm_buffer_.resize(layout_comm.GetMsgSize());
+    // The exchange plan is per DOF holder; the field payload carries
+    // num_components values per holder, so the buffer is scaled accordingly.
+    comm_buffer_.resize(
+      layout_comm.GetMsgSize() *
+      static_cast<size_t>(field.GetLayout().GetNumComponents()));
     comm_ =
       layout_comm_.GetChannel().CreateComm<T>(name, layout_comm.GetMPIComm());
     layout_comm_.SetOutMessageLayout<T>(comm_);

src/pcms/coupler/field_layout_communicator.h

@@ -41,7 +41,19 @@ class FieldLayoutCommunicator
   template <typename T>
   void SetOutMessageLayout(redev::BidirectionalComm<T>& comm)
   {
-    comm.SetOutMessageLayout(plan_.dest_ranks, plan_.offsets);
+    // plan_.offsets is a per-holder CSR layout. Each holder carries
+    // num_components field values, so the field message's offsets are scaled
+    // by num_components while the destination ranks are unchanged.
+    const int num_comp = layout_.GetNumComponents();
+    if (num_comp == 1) {
+      comm.SetOutMessageLayout(plan_.dest_ranks, plan_.offsets);
+      return;
+    }
+    redev::LOs scaled_offsets(plan_.offsets.size());
+    for (size_t i = 0; i < plan_.offsets.size(); ++i) {
+      scaled_offsets[i] = plan_.offsets[i] * num_comp;
+    }
+    comm.SetOutMessageLayout(plan_.dest_ranks, scaled_offsets);
   }
 
   void UpdateLayout();

src/pcms/coupler/field_serializer.h

@@ -20,36 +20,42 @@ class FieldSerializer
   {
     auto data = field.GetDOFHolderDataHost();
     auto owned = layout.GetOwnedHost();
-    // owned/permutation/buffer index the flat node-major wire order
-    // (dof * num_components + component); read the field values through the
-    // 2D [dof][comp] view.
+    // The exchange plan is per DOF holder: owned[i] and permutation[i] are
+    // indexed by holder. All num_components components of a holder share its
+    // location, so they occupy one contiguous block permutation[i]*num_comp in
+    // the wire buffer.
     if (buffer.size() > 0) {
       const LO num_dof = static_cast<LO>(data.extent(0));
       const LO num_comp = static_cast<LO>(data.extent(1));
       for (LO i = 0; i < num_dof; ++i) {
-        for (LO c = 0; c < num_comp; ++c) {
-          const LO flat = i * num_comp + c;
-          if (owned[flat])
-            buffer[permutation[flat]] = data(i, c);
+        if (owned[i]) {
+          for (LO c = 0; c < num_comp; ++c) {
+            buffer[permutation[i] * num_comp + c] = data(i, c);
+          }
         }
       }
     }
-    return static_cast<int>(data.size());
+    return static_cast<int>(buffer.size());
   }
 
   virtual void Deserialize(
     FieldData<T>& field, const FieldLayout& layout,
     Rank1View<const T, HostMemorySpace> buffer,
     Rank1View<const LO, HostMemorySpace> permutation) const
   {
-    Kokkos::View<T*, HostMemorySpace> sorted("sorted", permutation.size());
+    const LO num_dof = layout.GetNumOwnedDofHolder();
+    const LO num_comp = layout.GetNumComponents();
+    Kokkos::View<T*, HostMemorySpace> sorted("sorted", layout.OwnedSize());
     auto owned = layout.GetOwnedHost();
-    for (LO i = 0; i < static_cast<LO>(sorted.size()); ++i) {
-      if (owned[i])
-        sorted[i] = buffer[permutation[i]];
+    for (LO i = 0; i < num_dof; ++i) {
+      if (owned[i]) {
+        for (LO c = 0; c < num_comp; ++c) {
+          sorted[i * num_comp + c] = buffer[permutation[i] * num_comp + c];
+        }
+      }
     }
-    field.SetDOFHolderDataHost(Rank2View<const T, HostMemorySpace>(
-      sorted.data(), layout.GetNumOwnedDofHolder(), layout.GetNumComponents()));
+    field.SetDOFHolderDataHost(
+      Rank2View<const T, HostMemorySpace>(sorted.data(), num_dof, num_comp));
   }
 
   virtual ~FieldSerializer() noexcept = default;

src/pcms/coupler/serializer/xgc.h

@@ -35,14 +35,15 @@ class XGCFieldSerializer : public FieldSerializer<T>
 
     auto data = xgc_field->GetDOFHolderDataHost();
     auto owned = layout.GetOwnedHost();
+    // Per-holder plan: owned[i]/permutation[i] index holders; a holder's
+    // num_components values form one contiguous block in the wire buffer.
     if (buffer.size() > 0) {
       const LO num_dof = static_cast<LO>(data.extent(0));
       const LO num_comp = static_cast<LO>(data.extent(1));
       for (LO i = 0; i < num_dof; ++i) {
-        for (LO c = 0; c < num_comp; ++c) {
-          const LO flat = i * num_comp + c;
-          if (owned[flat]) {
-            buffer[permutation[flat]] = data(i, c);
+        if (owned[i]) {
+          for (LO c = 0; c < num_comp; ++c) {
+            buffer[permutation[i] * num_comp + c] = data(i, c);
           }
         }
       }
@@ -72,9 +73,11 @@ class XGCFieldSerializer : public FieldSerializer<T>
       }
     }
     if (rank_participates_) {
-      for (LO i = 0; i < static_cast<LO>(current.size()); ++i) {
+      for (LO i = 0; i < num_dof; ++i) {
         if (owned[i]) {
-          full_data[i] = buffer[permutation[i]];
+          for (LO c = 0; c < num_comp; ++c) {
+            full_data[i * num_comp + c] = buffer[permutation[i] * num_comp + c];
+          }
         }
       }
     }

src/pcms/discretization/discretization/omega_h.hpp

@@ -28,6 +28,9 @@ class OmegaHDiscretization : public Discretization
   Rank1View<const ClassificationId, DeviceMemorySpace>
   GetEntityClassificationIds(int entity_dim) const override;
 
+  // The Omega_h mesh this discretization is defined on.
+  Omega_h::Mesh& GetMesh() const noexcept { return mesh_; }
+
 private:
   static constexpr int max_entity_dim_ = Region;
 

src/pcms/field/evaluator/uniform_grid.h

@@ -186,7 +186,7 @@ class UniformGridPointEvaluator : public PointEvaluator<Real, LayoutPolicy>
 
 private:
   std::shared_ptr<const UniformGridFieldLayout<Dim>> layout_;
-  const UniformGrid<Dim>& grid_;
+  UniformGrid<Dim> grid_;
   UniformGridFieldLocalizationHint<Dim> hint_;
   Real fill_value_;
 };
@@ -304,7 +304,7 @@ class UniformGridEvaluatorFactory : public FieldEvaluatorFactory<Real>
 
 private:
   std::shared_ptr<const UniformGridFieldLayout<Dim>> layout_;
-  const UniformGrid<Dim>& grid_;
+  UniformGrid<Dim> grid_;
 };
 
 using UniformGridEvaluatorFactory2D = UniformGridEvaluatorFactory<2>;

src/pcms/field/evaluator/uniform_grid_spline.h

@@ -171,7 +171,7 @@ class UniformGridSplinePointEvaluator2D
 
 private:
   std::shared_ptr<const UniformGridFieldLayout<2>> layout_;
-  const UniformGrid<2>& grid_;
+  UniformGrid<2> grid_;
   UniformGridFieldLocalizationHint<2> hint_;
   Real fill_value_;
   Kokkos::View<Real*, DeviceMemorySpace> x_coords_;
@@ -287,7 +287,7 @@ class UniformGridSplineEvaluatorFactory2D : public FieldEvaluatorFactory<Real>
 
 private:
   std::shared_ptr<const UniformGridFieldLayout<2>> layout_;
-  const UniformGrid<2>& grid_;
+  UniformGrid<2> grid_;
 };
 
 } // namespace pcms

src/pcms/field/layout/omega_h_lagrange.cpp

@@ -3,6 +3,7 @@
 #include "pcms/utility/mesh_geometry.h"
 #include "pcms/utility/omega_h_array_utils.h"
 #include "pcms/utility/profile.h"
+#include <Omega_h_class.hpp>
 #include <stdexcept>
 
 namespace pcms
@@ -48,6 +49,35 @@ Omega_h::Write<Omega_h::GO> BuildGids(Omega_h::Mesh& mesh, int entity_dim,
   return gids;
 }
 
+// Some meshes (e.g. certain XGC meshes) are stored without any geometric
+// classification (no class_id/class_dim tags). The layout and its
+// discretization read these tags on every mesh dimension, so derive a default
+// classification when it is missing. classify_elements +
+// finalize_classification populate a geometrically meaningful class_dim on all
+// dimensions but do not create class_id (that normally comes from a geometric
+// model we don't have here), so fill any missing class_id with 0.
+// Classification is only consulted for coupling DOF matching
+// (field_exchange_planner); standalone field transfer never reads it, so a
+// default id is sufficient to keep such meshes usable.
+void EnsureClassification(Omega_h::Mesh& mesh)
+{
+  if (mesh.has_tag(mesh.dim(), "class_id") &&
+      mesh.has_tag(mesh.dim(), "class_dim")) {
+    return;
+  }
+  if (!mesh.has_tag(mesh.dim(), "class_dim")) {
+    Omega_h::classify_elements(&mesh);
+    Omega_h::finalize_classification(&mesh);
+  }
+  for (int d = 0; d <= mesh.dim(); ++d) {
+    if (!mesh.has_tag(d, "class_id")) {
+      mesh.add_tag<Omega_h::ClassId>(
+        d, "class_id", 1,
+        Omega_h::Read<Omega_h::ClassId>(mesh.nents(d), 0, "class_id"));
+    }
+  }
+}
+
 Kokkos::View<bool*, DeviceMemorySpace> BuildOwned(Omega_h::Mesh& mesh,
                                                   int entity_dim)
 {
@@ -87,6 +117,7 @@ OmegaHLagrangeLayout::OmegaHLagrangeLayout(Omega_h::Mesh& mesh, int order,
   PCMS_FUNCTION_TIMER;
   int entity_dim = EntityDimForOrder(order_, mesh_.dim());
 
+  EnsureClassification(mesh_);
   gids_ = BuildGids(mesh_, entity_dim, global_id_name_);
   gids_host_ = Omega_h::HostWrite<Omega_h::GO>(gids_);
   coords_ = get_entity_centroids(mesh_, entity_dim);
@@ -130,6 +161,7 @@ OmegaHLagrangeLayout::OmegaHLagrangeLayout(
   PCMS_FUNCTION_TIMER;
   int entity_dim = EntityDimForOrder(order_, mesh_.dim());
 
+  EnsureClassification(mesh_);
   gids_ = BuildGids(mesh_, entity_dim, global_id_name_);
   gids_host_ = Omega_h::HostWrite<Omega_h::GO>(gids_);
   coords_ = get_entity_centroids(mesh_, entity_dim);

src/pcms/localization/adj_search.cpp

@@ -160,7 +160,10 @@ OMEGA_H_INLINE void add_support_if_unvisited_and_within_cutoff(
   if (!visited.notVisited(candidate_id))
     return;
 
-  visited.push_back(candidate_id);
+  // If the visited buffer is full, stop here so the BFS terminates. Without
+  // this, unrecorded vertices are treated as unvisited and re-queued forever.
+  if (!visited.push_back(candidate_id))
+    return;
 
   Omega_h::Real candidate_coords[max_dim];
   for (Omega_h::LO k = 0; k < dim; ++k) {

src/pcms/localization/queue_visited.hpp

@@ -56,7 +56,7 @@ class Track
   ~Track() {}
 
   OMEGA_H_INLINE
-  void push_back(const int& item);
+  bool push_back(const int& item);
 
   OMEGA_H_INLINE
   int size();
@@ -107,15 +107,18 @@ bool Queue::isFull() const
 }
 
 OMEGA_H_INLINE
-void Track::push_back(const int& item)
+bool Track::push_back(const int& item)
 {
   if (count == MAX_SIZE_TRACK) {
-    printf("track is full %d\n", count);
-    return;
+    // Visited buffer is full. Report failure so callers can stop expanding the
+    // search; otherwise new vertices can never be marked visited and the BFS
+    // re-queues them forever (infinite loop).
+    return false;
   }
   last = (last + 1) % MAX_SIZE_TRACK;
   tracking_array[last] = item;
   count++;
+  return true;
 }
 
 OMEGA_H_INLINE

test/test_omega_h_lagrange_field.cpp

@@ -210,6 +210,30 @@ TEST_CASE("OmegaHLagrangeField order-1: serialize / deserialize round-trip")
   pcms::test::CheckSerializeDeserialize(*factory.GetLayout(), field.GetData());
 }
 
+TEST_CASE(
+  "OmegaHLagrangeField order-1: multi-component serialize / deserialize "
+  "round-trip")
+{
+  auto lib = Omega_h::Library{};
+  auto mesh = MakeBox2D(lib.world());
+  const int nc = 3;
+  // The MeshFields backend is scalar-only; multi-component needs OmegaH.
+  auto factory = pcms::LagrangeFunctionSpace::FromMesh(
+    mesh, 1, nc, pcms::CoordinateSystem::Cartesian, "global",
+    pcms::LagrangeFunctionSpace::Backend::OmegaH);
+  auto field = factory.CreateField<Real>(pcms::FieldMetadata{});
+
+  const int n = factory.GetLayout()->GetNumOwnedDofHolder();
+  std::vector<Real> data(static_cast<size_t>(n) * nc);
+  for (int i = 0; i < n; ++i)
+    for (int c = 0; c < nc; ++c)
+      data[static_cast<size_t>(i) * nc + c] = i + 0.25 * c;
+  field.GetData().SetDOFHolderDataHost(
+    pcms::Rank2View<const Real, pcms::HostMemorySpace>(data.data(), n, nc));
+
+  pcms::test::CheckSerializeDeserialize(*factory.GetLayout(), field.GetData());
+}
+
 // ---- Order-0 field tests ----------------------------------------------------
 
 TEST_CASE("OmegaHLagrangeField order-0: set/get DOF data round-trip")