diff --git a/amr-wind/utilities/sampling/CMakeLists.txt b/amr-wind/utilities/sampling/CMakeLists.txt index c815cc8644..e1f7e540ff 100644 --- a/amr-wind/utilities/sampling/CMakeLists.txt +++ b/amr-wind/utilities/sampling/CMakeLists.txt @@ -12,5 +12,6 @@ target_sources(${amr_wind_lib_name} KineticEnergy.cpp Enstrophy.cpp FreeSurface.cpp + VolumeSampler.cpp WaveEnergy.cpp ) diff --git a/amr-wind/utilities/sampling/VolumeSampler.H b/amr-wind/utilities/sampling/VolumeSampler.H new file mode 100644 index 0000000000..3d5fdab88b --- /dev/null +++ b/amr-wind/utilities/sampling/VolumeSampler.H @@ -0,0 +1,63 @@ +#ifndef VOLUMESAMPLER_H +#define VOLUMESAMPLER_H + +#include "amr-wind/utilities/sampling/SamplerBase.H" + +namespace amr_wind::sampling { + +/** Sample data on uniform grids + * \ingroup sampling + * + * Defines probe locations on a volume defined by an + * ` volume that originates from the point `lo` and extends to `hi`. + * The resolution of the volume is defined by the vector npts_dir. + * + */ +class VolumeSampler : public SamplerBase::Register +{ +public: + static std::string identifier() { return "VolumeSampler"; } + + explicit VolumeSampler(const CFDSim& /*unused*/); + + ~VolumeSampler() override; + + /** Read user inputs and initialize the sampling object + * + * \param key Prefix used to parse inputs from file + */ + void initialize(const std::string& key) override; + + //! Populate and return a vector of probe locations to be sampled + void sampling_locations(SampleLocType& /*locs*/) const override; + + void + define_netcdf_metadata(const ncutils::NCGroup& /*unused*/) const override; + void + populate_netcdf_metadata(const ncutils::NCGroup& /*unused*/) const override; + + //! Name of this sampling object + std::string label() const override { return m_label; } + std::string& label() override { return m_label; } + + //! Unique identifier for this set of probe locations + int id() const override { return m_id; } + int& id() override { return m_id; } + + //! Number of probe locations along the line + int num_points() const override { return m_npts; } + +private: + amrex::Vector m_hi; + amrex::Vector m_lo; + amrex::Vector m_npts_dir; + + std::string m_label; + + int m_id{-1}; + int m_npts{0}; +}; + +} // namespace amr_wind::sampling + +#endif /* VOLUMESAMPLER_H */ diff --git a/amr-wind/utilities/sampling/VolumeSampler.cpp b/amr-wind/utilities/sampling/VolumeSampler.cpp new file mode 100644 index 0000000000..5b182ffa45 --- /dev/null +++ b/amr-wind/utilities/sampling/VolumeSampler.cpp @@ -0,0 +1,77 @@ +#include + +#include "amr-wind/utilities/sampling/VolumeSampler.H" + +#include "AMReX_ParmParse.H" + +namespace amr_wind::sampling { + +VolumeSampler::VolumeSampler(const CFDSim& /*unused*/) {} + +VolumeSampler::~VolumeSampler() = default; + +void VolumeSampler::initialize(const std::string& key) +{ + amrex::ParmParse pp(key); + + pp.getarr("hi", m_hi); + pp.getarr("lo", m_lo); + pp.getarr("num_points", m_npts_dir); + AMREX_ALWAYS_ASSERT(static_cast(m_hi.size()) == AMREX_SPACEDIM); + AMREX_ALWAYS_ASSERT(static_cast(m_lo.size()) == AMREX_SPACEDIM); + AMREX_ALWAYS_ASSERT(static_cast(m_npts_dir.size()) == AMREX_SPACEDIM); + + // Update total number of points + const size_t tmp = + static_cast(m_npts_dir[0]) * m_npts_dir[1] * m_npts_dir[2]; + if (tmp > static_cast(std::numeric_limits::max())) { + amrex::Abort( + "VolumeSampler: Volume definition " + key + + " exceeds 32-bit integer limits"); + } + m_npts = static_cast(tmp); +} + +void VolumeSampler::sampling_locations(SampleLocType& locs) const +{ + locs.resize(m_npts); + + const amrex::Array dx = { + (m_hi[0] - m_lo[0]) / m_npts_dir[0], + (m_hi[1] - m_lo[1]) / m_npts_dir[1], + (m_hi[2] - m_lo[2]) / m_npts_dir[2]}; + + int idx = 0; + for (int k = 0; k < m_npts_dir[2]; ++k) { + for (int j = 0; j < m_npts_dir[1]; ++j) { + for (int i = 0; i < m_npts_dir[0]; ++i) { + locs[idx][0] = m_lo[0] + dx[0] * i; + locs[idx][1] = m_lo[1] + dx[1] * j; + locs[idx][2] = m_lo[2] + dx[2] * k; + ++idx; + } + } + } +} + +#ifdef AMR_WIND_USE_NETCDF +void VolumeSampler::define_netcdf_metadata(const ncutils::NCGroup& grp) const +{ + const std::vector ijk{m_npts_dir[0], m_npts_dir[1], m_npts_dir[2]}; + grp.put_attr("sampling_type", identifier()); + grp.put_attr("ijk_dims", ijk); + grp.put_attr("lo", m_lo); + grp.put_attr("hi", m_hi); +} + +void VolumeSampler::populate_netcdf_metadata(const ncutils::NCGroup&) const {} +#else +void VolumeSampler::define_netcdf_metadata( + const ncutils::NCGroup& /*unused*/) const +{} +void VolumeSampler::populate_netcdf_metadata( + const ncutils::NCGroup& /*unused*/) const +{} +#endif + +} // namespace amr_wind::sampling diff --git a/docs/sphinx/user/inputs_Sampling.rst b/docs/sphinx/user/inputs_Sampling.rst index 6d9e2f110c..34a1465360 100644 --- a/docs/sphinx/user/inputs_Sampling.rst +++ b/docs/sphinx/user/inputs_Sampling.rst @@ -139,4 +139,16 @@ Example:: The first line of the file contains the total number of probes for this set. This is followed by the coordinates (three real numbers), one line per probe. +Sampling on a volume +````````````````````` +The ``VolumeSampler`` samples a 3D volume that starts at ``lo`` and +extends to ``hi``. The resolution in all directions is specified by +``num_points``. + +Example:: + + sampling.volume1.type = VolumeSampler + sampling.volume1.hi = 3.0 1.0 0.5 + sampling.volume1.lo = 0.0 0.0 -0.5 + sampling.volume1.num_points = 30 10 10 diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 4b46029be1..73484d199c 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -176,6 +176,7 @@ add_test_re(abl_godunov_timetable) add_test_re(abl_ksgsm84_godunov) add_test_re(abl_mol_cn) add_test_re(abl_mol_explicit) +add_test_re(abl_sampling) add_test_re(abl_stable) add_test_re(abl_unstable) add_test_re(abl_unstable_constant_wall_model) diff --git a/test/test_files/abl_sampling/abl_sampling.inp b/test/test_files/abl_sampling/abl_sampling.inp new file mode 100644 index 0000000000..b79db92da1 --- /dev/null +++ b/test/test_files/abl_sampling/abl_sampling.inp @@ -0,0 +1,85 @@ +#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨# +# 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.5 # Use this constant dt if > 0 +time.cfl = 0.95 # CFL factor + +#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨# +# INPUT AND OUTPUT # +#.......................................# +time.plot_interval = 10 # Steps between plot files +time.checkpoint_interval = 5 # 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 +turbulence.model = Smagorinsky +Smagorinsky_coeffs.Cs = 0.135 + + +incflo.physics = ABL +ICNS.source_terms = BoussinesqBuoyancy CoriolisForcing ABLForcing +BoussinesqBuoyancy.reference_temperature = 300.0 +ABL.reference_temperature = 300.0 +CoriolisForcing.latitude = 41.3 +ABLForcing.abl_forcing_height = 90 + +incflo.velocity = 6.128355544951824 5.142300877492314 0.0 + +ABL.temperature_heights = 650.0 750.0 1000.0 +ABL.temperature_values = 300.0 308.0 308.75 + +ABL.kappa = .41 +ABL.surface_roughness_z0 = 0.15 + +#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨# +# 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 = 1 1 0 # Periodicity x y z (0/1) + +# Boundary conditions +zlo.type = "wall_model" + +zhi.type = "slip_wall" +zhi.temperature_type = "fixed_gradient" +zhi.temperature = 0.003 # tracer is used to specify potential temperature gradient + +#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨# +# VERBOSITY # +#.......................................# +incflo.verbose = 0 # incflo_level + +#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨# +# SAMPLING # +#.......................................# +incflo.post_processing = sampling +sampling.output_frequency = 1 +sampling.fields = velocity +sampling.labels = volume1 +sampling.volume1.type = VolumeSampler +sampling.volume1.lo = 200.0 200.0 200.0 +sampling.volume1.hi = 500.0 500.0 500.0 +sampling.volume1.num_points = 30 10 10 + diff --git a/unit_tests/utilities/test_sampling.cpp b/unit_tests/utilities/test_sampling.cpp index e1d9cfa0d2..26d57262c4 100644 --- a/unit_tests/utilities/test_sampling.cpp +++ b/unit_tests/utilities/test_sampling.cpp @@ -4,6 +4,7 @@ #include "amr-wind/utilities/sampling/Sampling.H" #include "amr-wind/utilities/sampling/SamplingContainer.H" #include "amr-wind/utilities/sampling/PlaneSampler.H" +#include "amr-wind/utilities/sampling/VolumeSampler.H" namespace amr_wind_tests { @@ -292,4 +293,20 @@ TEST_F(SamplingTest, plane_sampler) #endif } +TEST_F(SamplingTest, volume_sampler) +{ + initialize_mesh(); + amrex::ParmParse pp("volume"); + pp.addarr("hi", amrex::Vector{1.0, 1.0, 1.0}); + pp.addarr("lo", amrex::Vector{0.0, 0.0, 0.0}); + pp.addarr("num_points", amrex::Vector{3, 5, 5}); + + amr_wind::sampling::VolumeSampler volume(sim()); + volume.initialize("volume"); + amr_wind::sampling::VolumeSampler::SampleLocType locs; + volume.sampling_locations(locs); + + ASSERT_EQ(locs.size(), 3 * 5 * 5); +} + } // namespace amr_wind_tests