[Data assimilation] Interpolation for the base background state and ensemble simulations

.Exec_dev/DataAssimilation/ERF_Prob.H

@@ -4,7 +4,6 @@
 #include <string>
 
 #include "AMReX_REAL.H"
-
 #include "ERF_ProbCommon.H"
 
 struct ProbParm : ProbParmDefaults {
@@ -29,6 +28,8 @@ class Problem : public ProblemBase
 public:
     Problem(const amrex_real* problo, const amrex_real* probhi);
 
+#include "Prob/ERF_InitDensityHSE.H"
+#include "Prob/ERF_InitRayleighDamping.H"
 #include "Prob/ERF_InitCustomPert.H"
 
 protected:

.Exec_dev/DataAssimilation/ERF_Prob.cpp

@@ -129,20 +129,24 @@ Problem::init_custom_pert_vels (
     Array4<Real      > const& x_vel_pert,
     Array4<Real      > const& y_vel_pert,
     Array4<Real      > const& z_vel_pert,
-    Array4<Real      > const& r_hse,
+    amrex::Array4<amrex::Real const> const& z_nd,
     GeometryData const& geomdata,
     Array4<Real const> const& /*mf_u*/,
     Array4<Real const> const& /*mf_v*/,
     const SolverChoice& /*sc*/,
     const int /*lev*/)
 {
+    ignore_unused(z_nd);
     // --------------------------------------------------------
     // Per-ensemble perturbation controls
     // --------------------------------------------------------
     Real ens_pert_amplitude = 0.0;
     ParmParse pp_ens("ensemble_pert");
     pp_ens.query("amplitude", ens_pert_amplitude);
 
+    Real xc = parms.xc; Real yc = parms.yc;
+    Real R  = parms.R ; Real beta = parms.beta;
+    Real sigma = parms.sigma;
 
   // Set the x-velocity
   ParallelFor(xbx, [=, parms_d=parms] AMREX_GPU_DEVICE(int i, int j, int k) noexcept

.Exec_dev/DataAssimilation/GNUmakefile

@@ -28,6 +28,6 @@ USE_ASSERTION = TRUE
 Bpack := ./Make.package
 Blocs := .
 
-ERF_HOME := ../../..
-ERF_PROBLEM_DIR = $(ERF_HOME)/.Exec_dev/Assimilation
+ERF_HOME := ../..
+ERF_PROBLEM_DIR = $(ERF_HOME)/.Exec_dev/DataAssimilation
 include $(ERF_HOME)/Exec/Make.ERF

.Exec_dev/DataAssimilation/ReadPltFileAndInterpolate/Makefile

@@ -0,0 +1,58 @@
+
+
+#---------------------------
+# Makefile for AMReX project
+#---------------------------
+CXX = mpicxx
+CXXFLAGS = -O2 -std=c++17 -fopenmp -Wall
+
+# Replace with your AMReX paths
+AMREX_INC = /pscratch/sd/n/nataraj2/AIModCon/GitHub/amrex_repos/amrex/include
+AMREX_LIB = /pscratch/sd/n/nataraj2/AIModCon/GitHub/amrex_repos/amrex/lib
+
+HDRS = ReadPlotFile.H
+LIBS = -lamrex -lpthread -lgfortran -ldl -lm -lstdc++
+
+# Default: do nothing if 'make' is called without target
+all:
+	@echo "Use 'make serial' or 'make parallel'"
+
+# -------------------------
+# Serial build
+# -------------------------
+SRCS_SERIAL = main_serial.cpp ReadPlotFile.cpp
+OBJS_SERIAL = $(SRCS_SERIAL:.cpp=.o)
+EXE_SERIAL  = out_serial
+
+serial: $(EXE_SERIAL)
+
+$(EXE_SERIAL): $(OBJS_SERIAL)
+	$(CXX) $(CXXFLAGS) -o $@ $^ -L$(AMREX_LIB) $(LIBS)
+
+# -------------------------
+# Parallel build
+# -------------------------
+SRCS_PARALLEL = main_parallel.cpp ReadPlotFile.cpp
+OBJS_PARALLEL = $(SRCS_PARALLEL:.cpp=.o)
+EXE_PARALLEL  = out_parallel
+
+parallel: $(EXE_PARALLEL)
+
+$(EXE_PARALLEL): $(OBJS_PARALLEL)
+	$(CXX) $(CXXFLAGS) -o $@ $^ -L$(AMREX_LIB) $(LIBS)
+
+# -------------------------
+# Compile rules for object files
+# -------------------------
+%.o: %.cpp $(HDRS)
+	$(CXX) $(CXXFLAGS) -I$(AMREX_INC) -c $< -o $@
+
+# -------------------------
+# Clean
+# -------------------------
+clean:
+	rm -f *.o out_serial out_parallel
+
+.PHONY: all serial parallel clean
+
+

.Exec_dev/DataAssimilation/ReadPltFileAndInterpolate/ReadPlotFile.H

@@ -0,0 +1,60 @@
+#ifndef READ_PLOTFILE_H_
+#define READ_PLOTFILE_H_
+
+#include <string>
+#include <AMReX.H>
+#include <AMReX_Vector.H>
+#include <AMReX_MultiFab.H>
+#include <AMReX_PlotFileUtil.H>
+
+amrex::Vector<std::string> ReadVarNames(const std::string& filename);
+
+// ------------------------------------------------------------
+// ReadPlotFile
+//
+// var_filename  : text file containing variable names (one per line)
+// plot_filename : AMReX plotfile directory
+// mf            : MultiFab that will be defined and filled
+// ------------------------------------------------------------
+void ReadPlotFile(const std::string& var_filename,
+                  amrex::PlotFileData& pf,
+                  amrex::MultiFab& mf);
+
+void
+GetCoarseMultiFabOnFineDMap(const amrex::Geometry& geom_coarse,
+                            const amrex::Geometry& geom_fine,
+                            const amrex::MultiFab& multifab_coarse,
+                            const amrex::MultiFab& multifab_fine,
+                            amrex::MultiFab& coarse_multifab_on_fine_dmap);
+
+void
+CreateNodalMultiFabFromCellCenteredMultiFab (amrex::MultiFab& mf_nc,      // output nodal MF
+                    amrex::MultiFab& mf_cc,      // input cell-centered MF (coarse)
+                    const amrex::Geometry& geom);
+
+void
+CreateCellCenteredMultiFabFromNodalMultiFab(amrex::MultiFab& tmp_cc,   // cell-centered output
+                                            amrex::MultiFab& mf_nc);   // node-centered input
+
+// Populate a fine, cell-centered MultiFab from coarse nodal data using
+// trilinear interpolation at fine cell centers.
+// mf_cc_fine provides the fine BoxArray + DistributionMapping template.
+void
+PopulateFineCellCenteredFromCoarseNodal(const amrex::Geometry& geom_coarse,
+                                        const amrex::Geometry& geom_fine,
+                                        const amrex::MultiFab& coarse_multifab_on_fine_dmap,
+                                        const amrex::MultiFab& mf_cc_fine,
+                                        amrex::MultiFab& mf_cc_from_coarse);
+
+void
+check_large_values(const amrex::MultiFab& mf);
+
+void
+ApplyNeumannBCs(const amrex::Geometry& geom,
+                     amrex::MultiFab& mf_cc);
+
+void WriteCustomDataFile(const amrex::Geometry& geom,
+                         const amrex::MultiFab& mf_cc,
+                         const std::string& filename_custom);
+
+#endif