Merge pull request seahorce-scidac#71 from seahorce-scidac/remove_rho…

…_comp

Remove Rho_comp and Theta_comp and associated ...

Exec/ParticlesOverSeaMount/prob.cpp

@@ -139,15 +139,14 @@ init_custom_prob(
         const Real z = z_r(i,j,k);
 
         state(i, j, k, Temp_comp) = 1.;
-        state(i, j, k, Rho_comp) = 1.;
 
         state(i,j,k,Temp_comp)=parms.T0+8.0*std::exp(z/50.0_rt);
 #ifdef ROMSX_USE_SALINITY
         state(i,j,k,Salt_comp)=parms.S0;
 #endif
 
         // Set scalar = 0 everywhere
-        state(i, j, k, RhoScalar_comp) = parms.rho0;
+        state(i, j, k, Scalar_comp) = 0.0;
     });
 
   // Construct a box that is on x-faces

Exec/Seamount/prob.cpp

@@ -160,15 +160,14 @@ init_custom_prob(
         const Real z = z_r(i,j,k);
 
         state(i, j, k, Temp_comp) = 1.;
-        state(i, j, k, Rho_comp) = 1.;
 
         state(i,j,k,Temp_comp)=parms.T0+7.5_rt*std::exp(z/1000.0_rt);
 #ifdef ROMSX_USE_SALINITY
         state(i,j,k,Salt_comp)=parms.S0;
 #endif
 
         // Set scalar = 0 everywhere
-        state(i, j, k, RhoScalar_comp) = parms.rho0;
+        state(i, j, k, Scalar_comp) = 0.0;
     });
 
   // Construct a box that is on x-faces

Exec/Upwelling/prob.cpp

@@ -147,15 +147,14 @@ init_custom_prob(
         const Real z = z_r(i,j,k);
 
         state(i, j, k, Temp_comp) = 1.;
-        state(i, j, k, Rho_comp) = 1.;
 
         state(i,j,k,Temp_comp)=parms.T0+8.0*std::exp(z/50.0_rt);
 #ifdef ROMSX_USE_SALINITY
         state(i,j,k,Salt_comp)=parms.S0;
 #endif
 
         // Set scalar = 0 everywhere
-        state(i, j, k, RhoScalar_comp) = parms.rho0;
+        state(i, j, k, Scalar_comp) = 0.0;
     });
 
   // Construct a box that is on x-faces

Source/DataStruct.H

@@ -18,18 +18,6 @@ enum class Coord {
     x, y, z
 };
 
-enum class Stagger {
-    None, x, y, z
-};
-
-enum class AdvectedQuantity {
-    unity, u, v, w, theta, scalar
-};
-
-enum class AdvectingQuantity {
-    rho_u, rho_v, rho_w
-};
-
 enum class AdvectionDir {
     x, y, z
 };
@@ -38,18 +26,6 @@ enum class AdvectionScheme {
     centered4, upstream3
 };
 
-enum class DiffusionDir {
-    x, y, z
-};
-
-enum class MomentumEqn {
-    x, y, z
-};
-
-enum class MolecDiffType {
-    None, Constant, ConstantAlpha
-};
-
 enum class IC_BC_Type {
     Ideal, Real
 };

Source/Derive.H

@@ -6,12 +6,6 @@
 
 namespace derived {
 
-void romsx_derrhodivide(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  const amrex::FArrayBox& datfab,
-  const int scalar_index);
-
 void romsx_dernull(
   const amrex::Box& bx,
   amrex::FArrayBox& dromsxab,
@@ -22,40 +16,6 @@ void romsx_dernull(
   amrex::Real time,
   const int* bcrec,
   const int level);
-
-void romsx_derscalar(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int dcomp,
-  int ncomp,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& geomdata,
-  amrex::Real time,
-  const int* bcrec,
-  const int level);
-void romsx_deromega(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int dcomp,
-  int ncomp,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& geomdata,
-  amrex::Real time,
-  const int* bcrec,
-  const int level);
-
-#ifdef ROMSX_USE_SALINITY
-void romsx_dersalt(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int dcomp,
-  int ncomp,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& geomdata,
-  amrex::Real time,
-  const int* bcrec,
-  const int level);
-#endif
 }
 #endif
 

Source/Derive.cpp

@@ -4,23 +4,6 @@
 
 namespace derived {
 
-void romsx_derrhodivide(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  const amrex::FArrayBox& datfab,
-  const int scalar_index)
-{
-  // This routine divides any cell-centered conserved quantity by density
-  auto const dat = datfab.array();
-  auto primitive  = dromsxab.array();
-
-  amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
-    const amrex::Real rho       = dat(i, j, k, Rho_comp);
-    const amrex::Real conserved = dat(i, j, k, scalar_index);
-    primitive(i,j,k) = conserved / rho;
-  });
-}
-
 void
 romsx_dernull(
   const amrex::Box& /*bx*/,
@@ -35,51 +18,5 @@ romsx_dernull(
 {
   // This routine does nothing -- we use it as a placeholder.
 }
-
-void
-romsx_derscalar(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int /*dcomp*/,
-  int /*ncomp*/,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& /*geomdata*/,
-  amrex::Real /*time*/,
-  const int* /*bcrec*/,
-  const int /*level*/)
-{
-  romsx_derrhodivide(bx, dromsxab, datfab, RhoScalar_comp);
-}
-
-void
-romsx_deromega(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int /*dcomp*/,
-  int /*ncomp*/,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& /*geomdata*/,
-  amrex::Real /*time*/,
-  const int* /*bcrec*/,
-  const int /*level*/)
-{
-  romsx_derrhodivide(bx, dromsxab, datfab, Omega_comp);
-}
-#ifdef ROMSX_USE_SALINITY
-void
-romsx_dersalt(
-  const amrex::Box& bx,
-  amrex::FArrayBox& dromsxab,
-  int /*dcomp*/,
-  int /*ncomp*/,
-  const amrex::FArrayBox& datfab,
-  const amrex::Geometry& /*geomdata*/,
-  amrex::Real /*time*/,
-  const int* /*bcrec*/,
-  const int /*level*/)
-{
-  romsx_derrhodivide(bx, dromsxab, datfab, Salt_comp);
-}
-#endif
 }
 

Source/IO/Plotfile.cpp

@@ -171,9 +171,6 @@ ROMSX::WritePlotFile (int which, Vector<std::string> plot_var_names)
             }
         };
 
-        // Note: All derived variables must be computed in order of "derived_names" defined in ROMSX.H
-        calculate_derived("scalar",      derived::romsx_derscalar);
-
 #ifdef ROMSX_USE_PARTICLES
         if (containerHasElement(plot_var_names, "tracer_particle_count"))
         {

Source/IO/ReadFromWRFBdy.cpp

@@ -557,18 +557,12 @@ convert_wrfbdy_data(int which, const Box& domain, Vector<Vector<FArrayBox>>& bdy
 
     Array4<Real const> u_arr   = NC_xvel_fab.const_array();
     Array4<Real const> v_arr   = NC_yvel_fab.const_array();
-    Array4<Real const> r_arr   = NC_rho_fab.const_array();
-    Array4<Real const> rth_arr = NC_temp_fab.const_array();
 
     int ntimes = bdy_data.size();
     for (int nt = 0; nt < ntimes; nt++)
     {
         Array4<Real> bdy_u_arr  = bdy_data[nt][WRFBdyVars::U].array();  // This is face-centered
         Array4<Real> bdy_v_arr  = bdy_data[nt][WRFBdyVars::V].array();
-        Array4<Real> bdy_r_arr  = bdy_data[nt][WRFBdyVars::R].array();
-        Array4<Real> bdy_t_arr  = bdy_data[nt][WRFBdyVars::T].array();
-        Array4<Real> bdy_qv_arr = bdy_data[nt][WRFBdyVars::QV].array();
-        Array4<Real> mu_arr     = bdy_data[nt][WRFBdyVars::MU].array(); // This is cell-centered
 
         int ilo  = domain.smallEnd()[0];
         int ihi  = domain.bigEnd()[0];
@@ -591,22 +585,6 @@ convert_wrfbdy_data(int which, const Box& domain, Vector<Vector<FArrayBox>>& bdy
             bdy_u_arr(i,j,k) = new_bdy;
         });
 
-#ifndef AMREX_USE_GPU
-        if (nt == 0) {
-            FArrayBox diff(bx_u,1);
-            diff.template copy<RunOn::Device>(bdy_data[0][WRFBdyVars::U]);
-            diff.template minus<RunOn::Device>(NC_xvel_fab);
-            if (which == 0)
-                amrex::Print() << "Max norm of diff between initial U and bdy U on lo x face: " << diff.norm(0) << std::endl;
-            if (which == 1)
-                amrex::Print() << "Max norm of diff between initial U and bdy U on hi x face: " << diff.norm(0) << std::endl;
-            if (which == 2)
-                amrex::Print() << "Max norm of diff between initial U and bdy U on lo y face: " << diff.norm(0) << std::endl;
-            if (which == 3)
-                amrex::Print() << "Max norm of diff between initial U and bdy U on hi y face: " << diff.norm(0) << std::endl;
-        }
-#endif
-
         auto& bx_v  = bdy_data[0][WRFBdyVars::V].box();
         amrex::ParallelFor(bx_v, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
             Real xmu;
@@ -623,87 +601,6 @@ convert_wrfbdy_data(int which, const Box& domain, Vector<Vector<FArrayBox>>& bdy
             bdy_v_arr(i,j,k) = new_bdy;
         });
 
-#ifndef AMREX_USE_GPU
-        if (nt == 0) {
-            FArrayBox diff(bx_v,1);
-            diff.template copy<RunOn::Device>(bdy_data[0][WRFBdyVars::V]);
-            diff.template minus<RunOn::Device>(NC_yvel_fab);
-            if (which == 0)
-                amrex::Print() << "Max norm of diff between initial V and bdy V on lo x face: " << diff.norm(0) << std::endl;
-            if (which == 1)
-                amrex::Print() << "Max norm of diff between initial V and bdy V on hi x face: " << diff.norm(0) << std::endl;
-            if (which == 2)
-                amrex::Print() << "Max norm of diff between initial V and bdy V on lo y face: " << diff.norm(0) << std::endl;
-            if (which == 3)
-                amrex::Print() << "Max norm of diff between initial V and bdy V on hi y face: " << diff.norm(0) << std::endl;
-        }
-#endif
-
-        auto& bx_t = bdy_data[0][WRFBdyVars::T].box(); // Note this is currently "THM" aka the perturbational moist pot. temp.
-
-        // Define density
-        amrex::ParallelFor(bx_t, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
-
-            Real xmu = c1h_arr(0,0,k) * (mu_arr(i,j,0) + mub_arr(i,j,0)) + c2h_arr(0,0,k);
-
-            Real dpht =  (ph_arr(i,j,k+1) + phb_arr(i,j,k+1)) - (ph_arr(i,j,k) + phb_arr(i,j,k));
-
-            bdy_r_arr(i,j,k) = -xmu / ( dpht * rdnw_arr(0,0,k) );
-
-            //if (nt == 0 and std::abs(r_arr(i,j,k) - bdy_r_arr(i,j,k)) > 0.) {
-            //    amrex::Print() << "INIT VS BDY DEN " << IntVect(i,j,k) << " " << r_arr(i,j,k) << " " << bdy_r_arr(i,j,k) <<
-            //                    " " << std::abs(r_arr(i,j,k) - bdy_r_arr(i,j,k)) << std::endl;
-            //}
-        });
-
-        // Define theta
-        amrex::Real theta_ref = 300.;
-        amrex::ParallelFor(bx_t, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
-
-            Real xmu  = (mu_arr(i,j,0) + mub_arr(i,j,0));
-            Real xmu_mult = c1h_arr(0,0,k) * xmu + c2h_arr(0,0,k);
-
-            Real new_bdy_Th = bdy_t_arr(i,j,k) / xmu_mult + theta_ref;
-
-            Real qv_fac = (1. + bdy_qv_arr(i,j,k) / 0.622 / xmu_mult);
-
-            new_bdy_Th /= qv_fac;
-
-            bdy_t_arr(i,j,k) = new_bdy_Th * bdy_r_arr(i,j,k);
-
-            //if (nt == 0 and std::abs(rth_arr(i,j,k) - bdy_t_arr(i,j,k)) > 0.) {
-            //    amrex::Print() << "INIT VS BDY TH " << IntVect(i,j,k) << " " << rth_arr(i,j,k) << " " << bdy_t_arr(i,j,k) <<
-            //                    " " << std::abs(th_arr(i,j,k) - bdy_t_arr(i,j,k)) << std::endl;
-            //}
-        });
-
-#ifndef AMREX_USE_GPU
-        if (nt == 0) {
-            FArrayBox diff(bx_t,1);
-            diff.template copy<RunOn::Device>(bdy_data[0][WRFBdyVars::R]);
-            //diff.template mult<RunOn::Device>(NC_rho_fab);
-            diff.template minus<RunOn::Device>(NC_rho_fab);
-            if (which == 0)
-                amrex::Print() << "Max norm of diff between initial r and bdy r on lo x face: " << diff.norm(0) << std::endl;
-            if (which == 1)
-                amrex::Print() << "Max norm of diff between initial r and bdy r on hi x face: " << diff.norm(0) << std::endl;
-            if (which == 2)
-                amrex::Print() << "Max norm of diff between initial r and bdy r on lo y face: " << diff.norm(0) << std::endl;
-            if (which == 3)
-                amrex::Print() << "Max norm of diff between initial r and bdy r on hi y face: " << diff.norm(0) << std::endl;
-
-            diff.template copy<RunOn::Device>(bdy_data[0][WRFBdyVars::T]);
-            diff.template minus<RunOn::Device>(NC_temp_fab);
-            if (which == 0)
-                amrex::Print() << "Max norm of diff between initial rTh and bdy rTh on lo x face: " << diff.norm(0) << std::endl;
-            if (which == 1)
-                amrex::Print() << "Max norm of diff between initial rTh and bdy rTh on hi x face: " << diff.norm(0) << std::endl;
-            if (which == 2)
-                amrex::Print() << "Max norm of diff between initial rTh and bdy rTh on lo y face: " << diff.norm(0) << std::endl;
-            if (which == 3)
-                amrex::Print() << "Max norm of diff between initial rTh and bdy rTh on hi y face: " << diff.norm(0) << std::endl;
-        }
-#endif
     } // ntimes
 }
 #endif // ROMSX_USE_NETCDF

Source/IO/ReadFromWRFInput.cpp

@@ -235,17 +235,5 @@ ROMSX::read_from_wrfinput(int lev, int idx,
         {
             w_arr(i,j,k) /= msfw_arr(i,j,0);
         });
-
-        //
-        // WRF decomposes (1/rho) rather than rho so rho = 1/(ALB + AL)
-        //
-        NC_rho_fab[idx].template plus<RunOn::Device>(NC_rhop_fab[idx], 0, 0, 1);
-        NC_rho_fab[idx].template invert<RunOn::Device>(1.0);
-
-        const Real theta_ref = 300.0;
-        NC_temp_fab[idx].template plus<RunOn::Device>(theta_ref);
-
-        // Now multiply by rho to get (rho theta) instead of theta
-        NC_temp_fab[idx].template mult<RunOn::Device>(NC_rho_fab[idx],0,0,1);
 }
 #endif // ROMSX_USE_NETCDF