Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .github/workflows/clang-tidy.yml
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,7 @@ jobs:
with:
repo-name: 'meshFields-openmpi'
repo-path: 'SCOREC/meshFields'
repo-ref: 'b1482bbba288df210784b2345eae08e34faabdc4'
repo-ref: 'b9984520c20f73955a3b4e25e4907780b91d430f'
cache: true
options: '-DCMAKE_CXX_COMPILER=`which mpicxx`
-DCMAKE_C_COMPILER=`which mpicc`
Expand Down
2 changes: 1 addition & 1 deletion .github/workflows/cmake-test.yml
Original file line number Diff line number Diff line change
Expand Up @@ -132,7 +132,7 @@ jobs:
with:
repo-name: 'meshFields'
repo-path: 'SCOREC/meshFields'
repo-ref: 'b1482bbba288df210784b2345eae08e34faabdc4'
repo-ref: 'b9984520c20f73955a3b4e25e4907780b91d430f'
cache: true
cache-suffix: ${{ matrix.python_api == 'ON' && '-shared' || '' }}
options: '-DCMAKE_CXX_COMPILER=`which mpicxx`
Expand Down
28 changes: 15 additions & 13 deletions src/pcms/field/evaluator/mesh_fields_backend.h
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,7 @@

#include <Kokkos_Core.hpp>
#include <MeshField.hpp>
#include <MeshField_Config.hpp>
#include <memory>

#include "pcms/field/layout/mesh_fields.h"
Expand Down Expand Up @@ -53,7 +54,7 @@ class MeshFieldBackendImpl : public MeshFieldBackend<T>
{
auto self = const_cast<MeshFieldBackendImpl<T, Dim, Order>*>(this);
return self->mesh_field_.triangleLocalPointEval(localCoords, offsets,
shape_field_);
shape_field_.field);
}

void SetData(Rank1View<const T, DeviceMemorySpace> data, size_t num_nodes,
Expand All @@ -65,7 +66,7 @@ class MeshFieldBackendImpl : public MeshFieldBackend<T>
mesh_.nents(dim), KOKKOS_CLASS_LAMBDA(size_t ent) {
for (size_t n = 0; n < num_nodes; ++n) {
for (size_t c = 0; c < num_components; ++c) {
shape_field_(ent, n, c, topo) =
shape_field_.field(ent, n, c, topo) =
data[ent * stride + n * num_components + c];
}
}
Expand All @@ -82,7 +83,7 @@ class MeshFieldBackendImpl : public MeshFieldBackend<T>
for (size_t n = 0; n < num_nodes; ++n) {
for (size_t c = 0; c < num_components; ++c) {
data[ent * stride + n * num_components + c] =
shape_field_(ent, n, c, topo);
shape_field_.field(ent, n, c, topo);
}
}
});
Expand All @@ -91,9 +92,8 @@ class MeshFieldBackendImpl : public MeshFieldBackend<T>
private:
Omega_h::Mesh& mesh_;
MeshField::OmegahMeshField<DefaultExecutionSpace, Dim> mesh_field_;
using ShapeField =
decltype(mesh_field_.template CreateLagrangeField<T, Order, 1>());
ShapeField shape_field_;
using FWC = decltype(mesh_field_.template CreateLagrangeField<T, Order, 1>());
FWC shape_field_; // FWC = FieldWithController; keeps ctrlr alive
};

// ---------------------------------------------------------------------------
Expand Down Expand Up @@ -240,7 +240,7 @@ struct FillCoordinatesAndIndicesFunctor
const auto owner_idx = owning_elem_ids_(i);
LO count = Kokkos::atomic_sub_fetch(&elem_counts_(owner_idx), 1);
LO index = offsets_(owner_idx) + count - 1;
for (int j = 0; j < (dim_ + 1); ++j) {
for (int j = 0; j < dim_; ++j) {
coordinates_(index, j) = coord[j];
}
indices_(index) = i;
Expand Down Expand Up @@ -404,8 +404,9 @@ struct FillCoordinatesDeviceFunctor
global_coords_(orig_idx, 1)};
const auto local =
Omega_h::barycentric_from_global<2, 2>(point, vertex_coords);
for (int j = 0; j < (dim_ + 1); ++j)
for (int j = 0; j < dim_; ++j) {

Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@cwsmith @Joshua-Kloepfer This code is required due to the assert statements in MeshFields. I'm wondering if those asserts should have tolerances to allow for floating point computations, or if there is a better way to handle this.

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A PR with a looser tolerance is here: SCOREC/meshFields#97

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I can remove the check here in pcms, but this seems to create a few incompatible versions between meshfields 1.0.0 and latest. Would that be a problem? If so, we'd probably need a patch on the meshfields side.

coordinates_(index, j) = local[j];
}
indices_(index) = orig_idx;
}
};
Expand Down Expand Up @@ -495,7 +496,7 @@ struct MeshFieldsAdapter2LocalizationHint

offsets_ = Kokkos::View<LO*, HostMemorySpace>("offsets", mesh.nelems() + 1);
coordinates_ = Kokkos::View<Real**, HostMemorySpace>(
"coordinates", num_valid_, mesh.dim() + 1);
"coordinates", num_valid_, mesh.dim());
indices_ = Kokkos::View<LO*, HostMemorySpace>("indices", num_valid_);

if (num_missing_ > 0) {
Expand Down Expand Up @@ -537,8 +538,9 @@ struct MeshFieldsAdapter2LocalizationHint
global_coords(orig_idx, 1)};
const auto local =
Omega_h::barycentric_from_global<2, 2>(point, vertex_coords);
for (int j = 0; j < (mesh.dim() + 1); ++j)
for (int j = 0; j < mesh.dim(); ++j) {
coordinates_(index, j) = local[j];
}
indices_(index) = static_cast<LO>(orig_idx);
}

Expand Down Expand Up @@ -632,7 +634,7 @@ struct MeshFieldsAdapter2LocalizationHint

// Step 6: Fill coordinates and indices on device
coordinates_d_ =
Kokkos::View<Real**>("coordinates_d", num_valid_, mesh.dim() + 1);
Kokkos::View<Real**>("coordinates_d", num_valid_, mesh.dim());
indices_d_ = Kokkos::View<LO*>("indices_d", num_valid_);
const auto tris2verts = mesh.ask_elem_verts();
const auto mesh_coords = mesh.coords();
Expand All @@ -654,7 +656,7 @@ struct MeshFieldsAdapter2LocalizationHint
// Create host mirrors for compatibility (lazy copy - only if needed)
offsets_ = Kokkos::create_mirror_view(offsets_d_);
coordinates_ = Kokkos::View<Real**, HostMemorySpace>(
"coordinates_", num_valid_, mesh.dim() + 1);
"coordinates_", num_valid_, mesh.dim());
DeepCopyMismatchLayouts(coordinates_, coordinates_d_);
indices_ = Kokkos::create_mirror_view(indices_d_);
if (num_missing_ > 0) {
Expand All @@ -681,4 +683,4 @@ struct MeshFieldsAdapter2LocalizationHint

} // namespace pcms

#endif // PCMS_ADAPTER_MESHFIELDS_MESH_FIELDS_BACKEND_H
#endif // PCMS_ADAPTER_MESHFIELDS_MESH_FIELDS_BACKEND_H
6 changes: 3 additions & 3 deletions src/pcms/transfer/mass_matrix_integrator.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@

#include <KokkosController.hpp>
#include <MeshField.hpp>
#include <MeshField_Config.hpp>
#include <MeshField_Element.hpp>
#include <MeshField_Fail.hpp>
#include <MeshField_For.hpp>
Expand Down Expand Up @@ -38,7 +39,7 @@ class MassMatrixIntegrator : public MeshField::Integrator
// std::cerr << " Number points per Elem : " << numPtsPerElem << "\n";
assert(numPtsPerElem >= 1);
const size_t ptDim = p.extent(1);
assert(ptDim == fe.MeshEntDim + 1);
assert(ptDim == fe.MeshEntDim);
// Copy values needed in the kernel to avoid capturing host references
// (mesh and fe are host objects and cannot be dereferenced on the device)
const auto numElems = mesh.nelems();
Expand All @@ -51,8 +52,7 @@ class MassMatrixIntegrator : public MeshField::Integrator
const auto last = first + numPtsPerElem;
for (auto pt = first; pt < last; pt++) {
// FIXME better way to fill? pass kokkos::subview to getValues?
Kokkos::Array<MeshField::Real, FieldElement::MeshEntDim + 1>
localCoord;
Kokkos::Array<MeshField::Real, FieldElement::MeshEntDim> localCoord;
for (auto i = 0; i < localCoord.size(); i++) {
localCoord[i] = p(pt, i);
}
Expand Down
13 changes: 8 additions & 5 deletions src/pcms/transfer/omega_h_form_integrator_utils.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@
#include "pcms/utility/assert.h"
#include <MeshField_Integrate.hpp>
#include <MeshField_Shape.hpp>
#include <MeshField_Config.hpp>
#include <Omega_h_shape.hpp>
#include <Kokkos_MathematicalFunctions.hpp>

Expand Down Expand Up @@ -72,13 +73,15 @@ inline void CheckOmegaHScalarLagrangeLayout(
}

[[nodiscard]] OMEGA_H_INLINE Omega_h::Vector<2> GlobalFromBarycentric(
const MeshField::Vector3& barycentric_coord,
const MeshField::Vector2& barycentric_coord,
const Omega_h::Few<Omega_h::Vector<2>, 3>& verts_coord)
{
Omega_h::Vector<2> real_coords = {0.0, 0.0};
const Omega_h::Real xi3 = 1.0 - barycentric_coord[0] - barycentric_coord[1];
const Omega_h::Real xi[3] = {barycentric_coord[0], barycentric_coord[1], xi3};
for (int i = 0; i < 3; ++i) {
real_coords[0] += barycentric_coord[i] * verts_coord[i][0];
real_coords[1] += barycentric_coord[i] * verts_coord[i][1];
real_coords[0] += xi[i] * verts_coord[i][0];
real_coords[1] += xi[i] * verts_coord[i][1];
}
return real_coords;
}
Expand Down Expand Up @@ -250,15 +253,15 @@ OMEGA_H_INLINE void ForEachIntersectionSubtriangle(
// source_order + target_order), which are only known at construction.
struct IntegrationData
{
Kokkos::View<MeshField::Vector3*> bary_coords; // barycentric coordinates
Kokkos::View<MeshField::Vector2*> bary_coords; // barycentric coordinates
Kokkos::View<Omega_h::Real*> weights; // quadrature weights

explicit IntegrationData(int order)
{
auto ip_vec = MeshField::getIntegrationPoints<MeshField::Triangle>(order);
const std::size_t num_ip = ip_vec.size();

bary_coords = Kokkos::View<MeshField::Vector3*>("bary_coords", num_ip);
bary_coords = Kokkos::View<MeshField::Vector2*>("bary_coords", num_ip);
weights = Kokkos::View<Omega_h::Real*>("weights", num_ip);

auto bary_coords_host = Kokkos::create_mirror_view(bary_coords);
Expand Down
2 changes: 1 addition & 1 deletion src/pcms/transfer/omega_h_mass_integrator.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -123,7 +123,7 @@ OmegaHMassIntegrator::OmegaHMassIntegrator(
omf(mesh);
auto coordField = omf.getCoordField();
const auto [shp, map] = MeshField::Omegah::getTriangleElement<1>(mesh);
MeshField::FieldElement coordFe(mesh.nelems(), coordField, shp, map);
MeshField::FieldElement coordFe(mesh.nelems(), coordField.field, shp, map);
auto elm_mass_dev = buildElementMassMatrix(mesh, coordFe);

// Build COO sparsity pattern on device: each element contributes a 3x3 block.
Expand Down
5 changes: 5 additions & 0 deletions test/test_omega_h_mass_integrator.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@
#include <pcms/transfer/omega_h_mass_integrator.hpp>
#include <KokkosController.hpp>
#include <MeshField.hpp>
#include <MeshField_Config.hpp>
#include <MeshField_Element.hpp>
#include <petscksp.h>
#include <map>
Expand Down Expand Up @@ -60,7 +61,11 @@ std::map<std::pair<pcms::GO, pcms::GO>, pcms::Real> BuildReferenceMassMap(
omf(mesh);
auto coordField = omf.getCoordField();
const auto [shp, map] = MeshField::Omegah::getTriangleElement<1>(mesh);
#if MeshFields_VERSION < 10000
MeshField::FieldElement coordFe(mesh.nelems(), coordField, shp, map);
#else
MeshField::FieldElement coordFe(mesh.nelems(), coordField.field, shp, map);
#endif
auto elm_mass_dev = buildElementMassMatrix(mesh, coordFe);
auto elm_mass_host =
Kokkos::create_mirror_view_and_copy(Kokkos::HostSpace{}, elm_mass_dev);
Expand Down
Loading