Add a CustomVelocity/Scalar UDF and a Rankine UDF (Exawind#940)

amr-wind/boundary_conditions/BCInterface.cpp

@@ -127,8 +127,7 @@ std::pair<const std::string, const std::string> BCIface::get_dirichlet_udfs()
 
                 if (has_inflow_udf && (inflow_udf != val)) {
                     amrex::Abort(
-                        "BCVelocity: Inflow UDF must be same for all inflow "
-                        "faces");
+                        "BC: Inflow UDF must be same for all inflow faces");
                 } else {
                     inflow_udf = val;
                     has_inflow_udf = true;
@@ -143,8 +142,7 @@ std::pair<const std::string, const std::string> BCIface::get_dirichlet_udfs()
 
                 if (has_wall_udf && (wall_udf != val)) {
                     amrex::Abort(
-                        "BCVelocity: Wall UDF must be same for all wall "
-                        "faces");
+                        "BC: Wall UDF must be same for all wall faces");
                 } else {
                     wall_udf = val;
                     has_wall_udf = true;
@@ -311,13 +309,16 @@ void BCScalar::read_values()
     const auto& bctype = m_field.bc_type();
     auto& bcval = m_field.bc_values();
     const int ndim = m_field.num_comp();
+    const auto udfs = get_dirichlet_udfs();
+    const auto const_dirichlet_inflow = udfs.first == "ConstDirichlet";
+
     for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
         auto ori = oit();
         const auto& bcid = bcnames[ori];
         const auto bct = bctype[ori];
 
         amrex::ParmParse pp(bcid);
-        if (bct == BC::mass_inflow) {
+        if ((bct == BC::mass_inflow) && (const_dirichlet_inflow)) {
             pp.getarr(fname.c_str(), bcval[ori], 0, ndim);
         } else {
             pp.queryarr(fname.c_str(), bcval[ori], 0, ndim);

amr-wind/boundary_conditions/CMakeLists.txt

@@ -3,6 +3,7 @@ target_sources(${amr_wind_lib_name}
       #C++
       BCInterface.cpp
       FixedGradientBC.cpp
+      scalar_bcs.cpp
    )
 
 add_subdirectory(wall_models)

amr-wind/boundary_conditions/scalar_bcs.H

@@ -0,0 +1,36 @@
+#ifndef SCALAR_BCS_H
+#define SCALAR_BCS_H
+
+#include "amr-wind/core/FieldBCOps.H"
+#include "amr-wind/core/FieldFillPatchOps.H"
+#include "amr-wind/core/FieldRepo.H"
+#include "amr-wind/core/SimTime.H"
+#include "amr-wind/physics/udfs/CustomScalar.H"
+
+namespace amr_wind::scalar_bc {
+
+template <typename WallOp>
+void register_inflow_scalar_dirichlet(
+    Field& field,
+    const std::string& inflow_udf,
+    const amrex::AmrCore& mesh,
+    const SimTime& time)
+{
+    if (inflow_udf == "CustomScalar") {
+        using InflowOp = BCOpCreator<udf::CustomScalar, WallOp>;
+        field.register_fill_patch_op<FieldFillPatchOps<InflowOp>>(
+            mesh, time, InflowOp(field));
+    } else {
+        amrex::Abort("Scalar BC: Invalid dirichlet BC type = " + inflow_udf);
+    }
+}
+
+void register_scalar_dirichlet(
+    Field& field,
+    const amrex::AmrCore& mesh,
+    const SimTime& time,
+    std::pair<const std::string, const std::string> udfs);
+
+} // namespace amr_wind::scalar_bc
+
+#endif /* SCALAR_BCS_H */

amr-wind/boundary_conditions/scalar_bcs.cpp

@@ -0,0 +1,25 @@
+#include "amr-wind/boundary_conditions/scalar_bcs.H"
+
+namespace amr_wind::scalar_bc {
+void register_scalar_dirichlet(
+    Field& field,
+    const amrex::AmrCore& mesh,
+    const SimTime& time,
+    std::pair<const std::string, const std::string> udfs)
+{
+    const std::string inflow_udf = udfs.first;
+    const std::string wall_udf = udfs.second;
+
+    if ((inflow_udf == "ConstDirichlet") && (wall_udf == "ConstDirichlet")) {
+        return;
+    }
+
+    if (wall_udf != "ConstDirichlet") {
+        amrex::Abort(
+            "Scalar BC: Only constant dirichlet supported for Wall BC");
+    }
+
+    register_inflow_scalar_dirichlet<ConstDirichlet>(
+        field, inflow_udf, mesh, time);
+}
+} // namespace amr_wind::scalar_bc

amr-wind/boundary_conditions/velocity_bcs.H

@@ -8,6 +8,8 @@
 #include "amr-wind/physics/udfs/LinearProfile.H"
 #include "amr-wind/physics/udfs/PowerLawProfile.H"
 #include "amr-wind/physics/udfs/BurggrafLid.H"
+#include "amr-wind/physics/udfs/Rankine.H"
+#include "amr-wind/physics/udfs/CustomVelocity.H"
 
 namespace amr_wind::vel_bc {
 
@@ -30,6 +32,14 @@ void register_inflow_vel_dirichlet(
         using InflowOp = BCOpCreator<udf::BurggrafLid, WallOp>;
         field.register_fill_patch_op<FieldFillPatchOps<InflowOp>>(
             mesh, time, InflowOp(field));
+    } else if (inflow_udf == "Rankine") {
+        using InflowOp = BCOpCreator<udf::Rankine, WallOp>;
+        field.register_fill_patch_op<FieldFillPatchOps<InflowOp>>(
+            mesh, time, InflowOp(field));
+    } else if (inflow_udf == "CustomVelocity") {
+        using InflowOp = BCOpCreator<udf::CustomVelocity, WallOp>;
+        field.register_fill_patch_op<FieldFillPatchOps<InflowOp>>(
+            mesh, time, InflowOp(field));
     } else {
         amrex::Abort("Velocity BC: Invalid dirichlet BC type = " + inflow_udf);
     }

amr-wind/equation_systems/BCOps.H

@@ -4,6 +4,7 @@
 #include "amr-wind/equation_systems/PDETraits.H"
 #include "amr-wind/equation_systems/PDEOps.H"
 #include "amr-wind/boundary_conditions/BCInterface.H"
+#include "amr-wind/boundary_conditions/scalar_bcs.H"
 
 namespace amr_wind::pde {
 
@@ -25,6 +26,9 @@ struct BCOp<
     {
         BCScalar bc(m_fields.field);
         bc(PDE::default_bc_value);
+        const auto udfs = bc.get_dirichlet_udfs();
+        scalar_bc::register_scalar_dirichlet(
+            m_fields.field, m_fields.repo.mesh(), m_time, udfs);
 
         // Used for fillpatch operation on the source term
         BCSrcTerm bc_src(m_fields.src_term);

amr-wind/physics/udfs/CMakeLists.txt

@@ -5,4 +5,7 @@ target_sources(${amr_wind_lib_name}
   LinearProfile.cpp
   PowerLawProfile.cpp
   BurggrafLid.cpp
+  Rankine.cpp
+  CustomVelocity.cpp
+  CustomScalar.cpp
   )

amr-wind/physics/udfs/CustomScalar.H

@@ -0,0 +1,61 @@
+#ifndef CUSTOM_SCALAR_H
+#define CUSTOM_SCALAR_H
+
+#include "AMReX_Geometry.H"
+#include "AMReX_Gpu.H"
+
+namespace amr_wind {
+
+class Field;
+
+namespace udf {
+
+struct CustomScalar
+{
+    struct DeviceOp
+    {
+        // clang-format off
+        // Declare parameters here if needed. For example:
+        // amrex::Real foo{1.0};
+        // clang-format on
+
+        AMREX_GPU_DEVICE
+        inline void operator()(
+            const amrex::IntVect& /*iv*/,
+            amrex::Array4<amrex::Real> const& /*field*/,
+            amrex::GeometryData const& /*geom*/,
+            const amrex::Real /*time*/,
+            amrex::Orientation /*ori*/,
+            const int /*comp*/,
+            const int /*dcomp*/,
+            const int /*orig_comp*/) const
+        {
+            // Compute quantities to set the field values. For example:
+            // clang-format off
+            // const auto* problo = geom.ProbLo();
+            // const auto* dx = geom.CellSize();
+            // const auto x = problo[0] + (iv[0] + 0.5) * dx[0];
+            // const auto y = problo[1] + (iv[1] + 0.5) * dx[1];
+            // const auto z = problo[2] + (iv[2] + 0.5) * dx[2];
+            // const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> vel = {1.0, 0.0, 0.0};
+
+            // Once the above is done, fill the field as:
+            // field(iv[0], iv[1], iv[2], dcomp + comp) = vel[orig_comp + comp];
+            // clang-format on
+        }
+    };
+    using DeviceType = DeviceOp;
+
+    static std::string identifier() { return "CustomScalar"; }
+
+    explicit CustomScalar(const Field& fld);
+
+    DeviceType device_instance() const { return m_op; }
+
+    DeviceOp m_op;
+};
+
+} // namespace udf
+} // namespace amr_wind
+
+#endif /* CUSTOM_SCALAR_H */

amr-wind/physics/udfs/CustomScalar.cpp

@@ -0,0 +1,34 @@
+#include "amr-wind/physics/udfs/CustomScalar.H"
+#include "amr-wind/core/Field.H"
+#include "amr-wind/core/FieldRepo.H"
+#include "amr-wind/core/vs/vector.H"
+#include "amr-wind/equation_systems/icns/icns.H"
+
+#include "AMReX_ParmParse.H"
+
+namespace amr_wind::udf {
+
+CustomScalar::CustomScalar(const Field& fld)
+{
+    // This is a where the user can set some user defined variables
+    // This capability can be activated with the following in the input file:
+    // xlo.type = "mass_inflow"
+    // xlo.temperature.inflow_type = CustomScalar
+    // CustomScalar.foo = 1.0
+
+    // clang-format off
+    //{
+    //    amrex::ParmParse pp("CustomScalar");
+    //    pp.query("foo", m_op.foo);
+    //}
+    // clang-format on
+    const int ncomp = fld.num_comp();
+    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(
+        (ncomp == 1), "CustomScalar requires field with 1 component");
+
+    amrex::Abort(
+        "Please define the body of this function and the corresponding struct "
+        "in the header file before using it. Then remove this message");
+}
+
+} // namespace amr_wind::udf

amr-wind/physics/udfs/CustomVelocity.H

@@ -0,0 +1,62 @@
+#ifndef CUSTOM_VELOCITY_H
+#define CUSTOM_VELOCITY_H
+
+#include "AMReX_Geometry.H"
+#include "AMReX_Gpu.H"
+
+namespace amr_wind {
+
+class Field;
+
+namespace udf {
+
+struct CustomVelocity
+{
+    struct DeviceOp
+    {
+        // clang-format off
+        // Declare parameters here if needed. For example:
+        // amrex::Real foo{1.0};
+        // amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> bar = {0.0};
+        // clang-format on
+
+        AMREX_GPU_DEVICE
+        inline void operator()(
+            const amrex::IntVect& /*iv*/,
+            amrex::Array4<amrex::Real> const& /*field*/,
+            amrex::GeometryData const& /*geom*/,
+            const amrex::Real /*time*/,
+            amrex::Orientation /*ori*/,
+            const int /*comp*/,
+            const int /*dcomp*/,
+            const int /*orig_comp*/) const
+        {
+            // Compute quantities to set the field values. For example:
+            // clang-format off
+            // const auto* problo = geom.ProbLo();
+            // const auto* dx = geom.CellSize();
+            // const auto x = problo[0] + (iv[0] + 0.5) * dx[0];
+            // const auto y = problo[1] + (iv[1] + 0.5) * dx[1];
+            // const auto z = problo[2] + (iv[2] + 0.5) * dx[2];
+            // const amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> vel = {1.0, 0.0, 0.0};
+
+            // Once the above is done, fill the field as:
+            // field(iv[0], iv[1], iv[2], dcomp + comp) = vel[orig_comp + comp];
+            // clang-format on
+        }
+    };
+    using DeviceType = DeviceOp;
+
+    static std::string identifier() { return "CustomVelocity"; }
+
+    explicit CustomVelocity(const Field& fld);
+
+    DeviceType device_instance() const { return m_op; }
+
+    DeviceOp m_op;
+};
+
+} // namespace udf
+} // namespace amr_wind
+
+#endif /* CUSTOM_VELOCITY_H */

amr-wind/physics/udfs/CustomVelocity.cpp

@@ -0,0 +1,37 @@
+#include "amr-wind/physics/udfs/CustomVelocity.H"
+#include "amr-wind/core/Field.H"
+#include "amr-wind/core/FieldRepo.H"
+#include "amr-wind/core/vs/vector.H"
+#include "amr-wind/equation_systems/icns/icns.H"
+
+#include "AMReX_ParmParse.H"
+
+namespace amr_wind::udf {
+
+CustomVelocity::CustomVelocity(const Field& /*fld*/)
+{
+    // This is a where the user can set some user defined variables
+    // This capability can be activated with the following in the input file:
+    // xlo.type = "mass_inflow"
+    // xlo.velocity.inflow_type = CustomVelocity
+    // CustomVelocity.foo = 1.0
+
+    // clang-format off
+    //{
+    //    const int ncomp = fld.num_comp();
+    //    amrex::ParmParse pp("CustomVelocity");
+    //    pp.query("foo", m_op.foo);
+    //    amrex::Vector<amrex::Real> vel(0.0, ncomp);
+    //    pp.getarr("velocity", vel);
+    //    AMREX_ALWAYS_ASSERT(vel.size() == ncomp);
+    //    for (int i = 0; i < ncomp; ++i) {
+    //        m_op.bar[i] = vel[i];
+    //    }
+    //}
+    // clang-format on
+    amrex::Abort(
+        "Please define the body of this function and the corresponding struct "
+        "in the header file before using it. Then remove this message");
+}
+
+} // namespace amr_wind::udf