Single Boundary Plane

amr-wind/wind_energy/ABLBoundaryPlane.H

@@ -60,7 +60,8 @@ public:
         const int lev,
         const Field* /*fld*/,
         const amrex::Real time,
-        const amrex::Vector<amrex::Real>& /*times*/);
+        const amrex::Vector<amrex::Real>& /*times*/,
+        const bool frozen_data);
 
     void interpolate(const amrex::Real /*time*/);
     bool is_populated(amrex::Orientation /*ori*/) const;
@@ -139,6 +140,7 @@ public:
         const std::string& chkname, const Field& field) const;
 
     void read_file(const bool /* nph_target_time*/);
+    void read_frozen_file();
 
     void populate_data(
         const int /*lev*/,
@@ -256,6 +258,9 @@ private:
 
     //! output format for bndry output
     std::string m_out_fmt{"native"};
+
+    //! Read frozen data
+    bool m_read_frozen_data{false};
 };
 
 } // namespace amr_wind

amr-wind/wind_energy/ABLBoundaryPlane.cpp

@@ -162,37 +162,43 @@ void InletData::read_data_native(
     const int lev,
     const Field* fld,
     const amrex::Real time,
-    const amrex::Vector<amrex::Real>& times)
+    const amrex::Vector<amrex::Real>& times,
+    const bool frozen_data)
 {
     const size_t nc = fld->num_comp();
     const int nstart =
         static_cast<int>(m_components[static_cast<int>(fld->id())]);
 
     const int idx = utils::closest_index(times, time, constants::LOOSE_TOL);
     const int idxp1 = idx + 1;
-
-    m_tn = times[idx];
-    m_tnp1 = times[idxp1];
-
     auto ori = oit();
-
-    if (!(m_tn <= time + constants::LOOSE_TOL) ||
-        !(time <= m_tnp1 + constants::LOOSE_TOL)) {
-        amrex::Abort(
-            "ABLBoundaryPlane.cpp InletData::read_data_native() check "
-            "failed\n"
-            "Left time quantities should be <= right time quantities. Indices "
-            "supplied for debugging.\n"
-            "m_tn = " +
-            std::to_string(m_tn) + ", time + LOOSE_TOL = " +
-            std::to_string(time + constants::LOOSE_TOL) +
-            "\n"
-            "time = " +
-            std::to_string(time) + ", m_tnp1 + LOOSE_TOL = " +
-            std::to_string(m_tnp1 + constants::LOOSE_TOL) +
-            "\n"
-            "idx = " +
-            std::to_string(idx) + ", idxp1 = " + std::to_string(idxp1));
+    if (frozen_data) {
+        // For frozen data, both tn and tnp1 are set to the same time
+        m_tn = times[0];
+        m_tnp1 = times[idx];
+    } else {
+        m_tn = times[idx];
+        m_tnp1 = times[idxp1];
+
+        if (!(m_tn <= time + constants::LOOSE_TOL) ||
+            !(time <= m_tnp1 + constants::LOOSE_TOL)) {
+            amrex::Abort(
+                "ABLBoundaryPlane.cpp InletData::read_data_native() check "
+                "failed\n"
+                "Left time quantities should be <= right time quantities. "
+                "Indices "
+                "supplied for debugging.\n"
+                "m_tn = " +
+                std::to_string(m_tn) + ", time + LOOSE_TOL = " +
+                std::to_string(time + constants::LOOSE_TOL) +
+                "\n"
+                "time = " +
+                std::to_string(time) + ", m_tnp1 + LOOSE_TOL = " +
+                std::to_string(m_tnp1 + constants::LOOSE_TOL) +
+                "\n"
+                "idx = " +
+                std::to_string(idx) + ", idxp1 = " + std::to_string(idxp1));
+        }
     }
     AMREX_ALWAYS_ASSERT(fld->num_comp() == bndry_n[ori].nComp());
     AMREX_ASSERT(bndry_n[ori].boxArray() == bndry_np1[ori].boxArray());
@@ -315,6 +321,7 @@ ABLBoundaryPlane::ABLBoundaryPlane(CFDSim& sim)
     pp.queryarr("bndry_var_names", m_var_names);
     pp.get("bndry_file", m_filename);
     pp.query("bndry_output_format", m_out_fmt);
+    pp.query("bndry_read_frozen_data", m_read_frozen_data);
 
 #ifndef AMR_WIND_USE_NETCDF
     if (m_out_fmt == "netcdf") {
@@ -343,27 +350,38 @@ void ABLBoundaryPlane::post_init_actions()
     if (!m_is_initialized) {
         return;
     }
+    amrex::Print() << "Initializing ABL Boundary Plane IO\n"
+                   << m_read_frozen_data << std::endl;
     initialize_data();
     write_header();
     write_file();
     read_header();
-    read_file(false);
+    if (m_read_frozen_data) {
+        amrex::Print() << "Reading frozen ABL Boundary Plane Inflow Data\n";
+        read_frozen_file();
+    } else {
+        read_file(false);
+    }
 }
 
 void ABLBoundaryPlane::pre_advance_work()
 {
     if (!m_is_initialized) {
         return;
     }
-    read_file(true);
+    if (!m_read_frozen_data) {
+        read_file(true);
+    }
 }
 
 void ABLBoundaryPlane::pre_predictor_work()
 {
     if (!m_is_initialized) {
         return;
     }
-    read_file(false);
+    if (!m_read_frozen_data) {
+        read_file(false);
+    }
 }
 
 void ABLBoundaryPlane::post_advance_work()
@@ -1293,7 +1311,8 @@ void ABLBoundaryPlane::read_file(const bool nph_target_time)
                     bndry2[ori].read(facename2);
 
                     m_in_data.read_data_native(
-                        oit, bndry1, bndry2, lev, fld, time, m_in_times);
+                        oit, bndry1, bndry2, lev, fld, time, m_in_times,
+                        m_read_frozen_data);
                 }
             }
         }
@@ -1302,6 +1321,82 @@ void ABLBoundaryPlane::read_file(const bool nph_target_time)
     m_in_data.interpolate(time);
 }
 
+void ABLBoundaryPlane::read_frozen_file()
+{
+    BL_PROFILE("amr-wind::ABLBoundaryPlane::read_frozen_file");
+    if (m_io_mode != io_mode::input) {
+        return;
+    }
+    const amrex::Real time = m_time.current_time();
+    if (m_out_fmt == "native") {
+
+        const int t_step1 = m_in_timesteps[0];
+        const int t_step2 = m_in_timesteps[0];
+
+        const std::string chkname1 =
+            m_filename + amrex::Concatenate("/bndry_output", t_step1);
+        const std::string chkname2 =
+            m_filename + amrex::Concatenate("/bndry_output", t_step2);
+
+        const std::string level_prefix = "Level_";
+
+        const int nlevels = boundary_native_file_levels();
+        const auto bndry_bas =
+            read_bndry_native_boxarrays(chkname1, *(m_fields[0]));
+        for (int lev = 0; lev < nlevels; ++lev) {
+            for (auto* fld : m_fields) {
+                auto& field = *fld;
+
+                const auto& ba = bndry_bas[lev];
+                amrex::DistributionMapping dm{ba};
+
+                amrex::BndryRegister bndry1(
+                    ba, dm, m_in_rad, m_out_rad, m_extent_rad,
+                    field.num_comp());
+                amrex::BndryRegister bndry2(
+                    ba, dm, m_in_rad, m_out_rad, m_extent_rad,
+                    field.num_comp());
+
+                bndry1.setVal(1.0e13);
+                bndry2.setVal(1.0e13);
+
+                std::string filename1 = amrex::MultiFabFileFullPrefix(
+                    lev, chkname1, level_prefix, field.name());
+                std::string filename2 = amrex::MultiFabFileFullPrefix(
+                    lev, chkname2, level_prefix, field.name());
+
+                for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
+                    auto ori = oit();
+
+                    if ((!m_in_data.is_populated(ori)) ||
+                        ((field.bc_type()[ori] != BC::mass_inflow) &&
+                         (field.bc_type()[ori] != BC::mass_inflow_outflow))) {
+                        continue;
+                    }
+
+                    std::string facename1 =
+                        amrex::Concatenate(filename1 + '_', ori, 1);
+                    std::string facename2 =
+                        amrex::Concatenate(filename2 + '_', ori, 1);
+
+                    bndry1[ori].read(facename1);
+                    bndry2[ori].read(facename2);
+
+                    m_in_data.read_data_native(
+                        oit, bndry1, bndry2, lev, fld, time, m_in_times,
+                        m_read_frozen_data);
+                }
+            }
+        }
+    } else {
+        amrex::Abort(
+            "ABLBoundaryPlane::read_frozen_file() is only implemented for "
+            "native format.");
+    }
+
+    m_in_data.interpolate(time);
+}
+
 // NOLINTNEXTLINE(readability-convert-member-functions-to-static)
 void ABLBoundaryPlane::populate_data(
     const int lev,
@@ -1318,35 +1413,36 @@ void ABLBoundaryPlane::populate_data(
         return;
     }
 
-    if (!(m_in_data.tn() <= time + constants::LOOSE_TOL) &&
-        !(time <= m_in_data.tnp1() + constants::LOOSE_TOL)) {
-        amrex::Abort(
-            "ABLBoundaryPlane.cpp ABLBoundaryPlane::populate_data() check 1"
-            "failed\n"
-            "Left time quantities should be <= right time quantities\n"
-            "m_in_data.tn() = " +
-            std::to_string(m_in_data.tn()) + ", time + LOOSE_TOL = " +
-            std::to_string(time + constants::LOOSE_TOL) +
-            "\n"
-            "time = " +
-            std::to_string(time) + ", m_in_data.tnp1() + LOOSE_TOL = " +
-            std::to_string(m_in_data.tnp1() + constants::LOOSE_TOL));
-    }
-    if (!(std::abs(time - m_in_data.tinterp()) < constants::LOOSE_TOL)) {
-        amrex::Abort(
-            "ABLBoundaryPlane.cpp ABLBoundaryPlane::populate_data() check 2"
-            "failed\n"
-            "Left time quantities should be < right time quantities. "
-            "Additional quantities supplied on second line for debugging.\n"
-            "std::abs(time - m_in_data.tinterp()) = " +
-            std::to_string(std::abs(time - m_in_data.tinterp())) +
-            ", LOOSE_TOL = " + std::to_string(constants::LOOSE_TOL) +
-            "\n"
-            "time = " +
-            std::to_string(time) +
-            ", m_in_data.tinterp() = " + std::to_string(m_in_data.tinterp()));
+    if (!m_read_frozen_data) {
+        if (!(m_in_data.tn() <= time + constants::LOOSE_TOL) &&
+            !(time <= m_in_data.tnp1() + constants::LOOSE_TOL)) {
+            amrex::Abort(
+                "ABLBoundaryPlane.cpp ABLBoundaryPlane::populate_data() check 1"
+                "failed\n"
+                "Left time quantities should be <= right time quantities\n"
+                "m_in_data.tn() = " +
+                std::to_string(m_in_data.tn()) + ", time + LOOSE_TOL = " +
+                std::to_string(time + constants::LOOSE_TOL) +
+                "\n"
+                "time = " +
+                std::to_string(time) + ", m_in_data.tnp1() + LOOSE_TOL = " +
+                std::to_string(m_in_data.tnp1() + constants::LOOSE_TOL));
+        }
+        if (!(std::abs(time - m_in_data.tinterp()) < constants::LOOSE_TOL)) {
+            amrex::Abort(
+                "ABLBoundaryPlane.cpp ABLBoundaryPlane::populate_data() check 2"
+                "failed\n"
+                "Left time quantities should be < right time quantities. "
+                "Additional quantities supplied on second line for debugging.\n"
+                "std::abs(time - m_in_data.tinterp()) = " +
+                std::to_string(std::abs(time - m_in_data.tinterp())) +
+                ", LOOSE_TOL = " + std::to_string(constants::LOOSE_TOL) +
+                "\n"
+                "time = " +
+                std::to_string(time) + ", m_in_data.tinterp() = " +
+                std::to_string(m_in_data.tinterp()));
+        }
     }
-
     for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
         auto ori = oit();
         if ((!m_in_data.is_populated(ori)) ||

test/test_files/abl_rans_frozen_input/abl_rans_frozen_input.inp

@@ -0,0 +1,92 @@
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#            SIMULATION STOP            #
+#.......................................#
+time.stop_time               =   22000.0     # Max (simulated) time to evolve
+time.max_step                =   10          # Max number of time steps
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#         TIME STEP COMPUTATION         #
+#.......................................#
+time.fixed_dt         =   0.4        # Use this constant dt if > 0
+time.cfl              =   0.95         # CFL factor
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#            INPUT AND OUTPUT           #
+#.......................................#
+io.restart_file = "../abl_bndry_output_native/chk00005"
+time.plot_interval            =  1      # Steps between plot files
+time.checkpoint_interval      =  -1       # Steps between checkpoint files
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#               PHYSICS                 #
+#.......................................#
+incflo.gravity          =   0.  0. -9.81  # Gravitational force (3D)
+incflo.density          = 1.0          # Reference density
+incflo.use_godunov = 1
+incflo.diffusion_type = 2
+transport.viscosity = 1.0e-5
+transport.laminar_prandtl = 0.7
+transport.turbulent_prandtl = 0.3333
+transport.reference_temperature = 290.0
+turbulence.model                            = KLAxell
+TKE.source_terms                            = KransAxell
+incflo.physics = ABL
+ICNS.source_terms = CoriolisForcing GeostrophicForcing
+CoriolisForcing.east_vector = 1.0 0.0 0.0
+CoriolisForcing.north_vector = 0.0 1.0 0.0
+CoriolisForcing.latitude = 90.0
+CoriolisForcing.rotational_time_period = 125663.706143592
+GeostrophicForcing.geostrophic_wind = 10.0 0.0 0.0
+incflo.velocity = 10.0 0.0 0.0
+ABL.temperature_heights = 0.0 2000.0
+ABL.temperature_values = 290.0 290.0
+ABL.rans_1dprofile_file                     = "rans_1d.info"
+ABL.surface_temp_flux = 0.0 
+ABL.perturb_temperature = false
+ABL.cutoff_height = 50.0
+ABL.perturb_velocity = true
+ABL.perturb_ref_height = 50.0
+ABL.Uperiods = 4.0
+ABL.Vperiods = 4.0
+ABL.deltaU = 1.0
+ABL.deltaV = 1.0
+ABL.kappa = .41
+ABL.surface_roughness_z0 = 0.01
+ABL.bndry_file = "../abl_bndry_output_native/bndry_files"
+ABL.bndry_io_mode = 1
+ABL.bndry_var_names = velocity temperature tke
+ABL.bndry_output_format = native
+ABL.bndry_read_frozen_data = true
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#        ADAPTIVE MESH REFINEMENT       #
+#.......................................#
+amr.n_cell              = 48 48 48    # Grid cells at coarsest AMRlevel
+amr.max_level           = 0           # Max AMR level in hierarchy 
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#              GEOMETRY                 #
+#.......................................#
+geometry.prob_lo        =   0.       0.     0.  # Lo corner coordinates
+geometry.prob_hi        =   1000.  1000.  1000.  # Hi corner coordinates
+geometry.is_periodic    =   0   0   0   # Periodicity x y z (0/1)
+# Boundary conditions
+xlo.type = "mass_inflow_outflow"
+xlo.density = 1.0
+xlo.temperature = 0.0
+xlo.tke = 1.0 
+xhi.type = "mass_inflow_outflow"
+xhi.density = 1.0
+xhi.temperature = 0.0
+xhi.tke = 1.0 
+ylo.type = "mass_inflow_outflow"
+ylo.density = 1.0
+ylo.temperature = 0.0
+ylo.tke = 1.0 
+yhi.type = "mass_inflow_outflow"
+yhi.density = 1.0
+yhi.temperature = 0.0
+yhi.tke = 1.0 
+zlo.type =   "wall_model"
+zhi.type =   "slip_wall"
+zhi.temperature_type = "fixed_gradient"
+zhi.temperature = 0.0
+#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
+#              VERBOSITY                #
+#.......................................#
+incflo.verbose          =   0          # incflo_level