added test case for low Mach averaging (including new IC for Tombouli…

…des)

.gitattributes

@@ -26,3 +26,4 @@ test/ref_solns/plate150K_1step/restart_output.sol.h5 filter=lfs diff=lfs merge=l
 test/ref_solns/interpInlet/restart_output.sol.h5 filter=lfs diff=lfs merge=lfs -text
 test/meshes/spongeBox.msh filter=lfs diff=lfs merge=lfs -text
 test/ref_solns/sgsLoMach/restart_output.sol.h5 filter=lfs diff=lfs merge=lfs -text
+test/ref_solns/aveLoMach/restart_output.sol.h5 filter=lfs diff=lfs merge=lfs -text

src/tomboulides.cpp

@@ -46,7 +46,9 @@ using namespace mfem;
 
 /// forward declarations
 void vel_exact_tgv2d(const Vector &x, double t, Vector &u);
+void vel_tgv2d_uniform(const Vector &x, double t, Vector &u);
 void vel_exact_pipe(const Vector &x, double t, Vector &u);
+void vel_channel(const Vector &x, double t, Vector &u);
 
 static double radius(const Vector &pos) { return pos[0]; }
 FunctionCoefficient radius_coeff(radius);
@@ -386,10 +388,20 @@ void Tomboulides::initializeSelf() {
   // set IC if we have one at this point
   if (!ic_string_.empty()) {
     if (ic_string_ == "tgv2d") {
-      std::cout << "Setting tgv2d IC..." << std::endl;
+      if (rank0_) std::cout << "Setting tgv2d IC..." << std::endl;
       VectorFunctionCoefficient u_excoeff(2, vel_exact_tgv2d);
       u_excoeff.SetTime(0.0);
       u_curr_gf_->ProjectCoefficient(u_excoeff);
+    } else if (ic_string_ == "tgv2d_uniform") {
+      if (rank0_) std::cout << "Setting tgv2d+uniform IC..." << std::endl;
+      VectorFunctionCoefficient u_excoeff(2, vel_tgv2d_uniform);
+      u_excoeff.SetTime(0.0);
+      u_curr_gf_->ProjectCoefficient(u_excoeff);
+    } else if (ic_string_ == "channel") {
+      if (rank0_) std::cout << "Setting channel IC..." << std::endl;
+      VectorFunctionCoefficient u_excoeff(3, vel_channel);
+      u_excoeff.SetTime(0.0);
+      u_curr_gf_->ProjectCoefficient(u_excoeff);
     }
   }
 
@@ -1672,3 +1684,57 @@ void vel_exact_pipe(const Vector &x, double t, Vector &u) {
   u(0) = 0.0;
   u(1) = 2.0 * (1 - x[0] * x[0]);
 }
+
+/// Used to set the velocity IC with TG field and uniform
+void vel_tgv2d_uniform(const Vector &x, double t, Vector &u) {
+  const double u0 = 1.0;
+  const double F = 0.1;
+  const double PI = 3.14159265359;
+  double twoPi = 2.0 * PI;
+
+  u(0) = u0;
+  u(1) = 0.0;
+
+  u(0) += +F * std::sin(twoPi * x[0]) * std::cos(twoPi * x[1]);
+  u(1) += -F * std::cos(twoPi * x[0]) * std::sin(twoPi * x[1]);
+}
+
+/// Used to set the channel IC
+void vel_channel(const Vector &x, double t, Vector &u) {
+  double PI = 3.14159265359;
+  double Lx = 25.0;
+  double Ly = 2.0;
+  double Lz = 9.4;
+  double Umean = 1.0;
+  double uInt = 0.1;
+  int nModes = 4;
+  double uM;
+  double ax, by, cz;
+  double AA, BB, CC;
+  double wall;
+
+  // expects channel height (-1,1)
+  wall = (1.0 - std::pow(x(1), 8.0));
+  u(0) = Umean * wall;
+  u(1) = 0.0;
+  u(2) = 0.0;
+
+  for (int n = 1; n <= nModes; n++) {
+    ax = 4.0 * PI / Lx * (double)n;
+    by = 2.0 * PI / Ly * (double)n;
+    cz = 2.0 * PI / Lz * (double)n;
+
+    AA = 1.0;
+    BB = 1.0;
+    CC = -(AA * ax + BB * by) / cz;
+
+    uM = uInt / (double)n;
+
+    u(0) += uM * AA * cos(ax * (x(0) + (double)(n - 1) * Umean)) * sin(by * x(1)) *
+            sin(cz * (x(2) + 0.5 * (double)(n - 1) * Umean)) * wall;
+    u(1) += uM * BB * sin(ax * (x(0) + (double)(n - 1) * Umean)) * cos(by * x(1)) *
+            sin(cz * (x(2) + 0.5 * (double)(n - 1) * Umean)) * wall;
+    u(2) += uM * CC * sin(ax * (x(0) + (double)(n - 1) * Umean)) * sin(by * x(1)) *
+            cos(cz * (x(2) + 0.5 * (double)(n - 1) * Umean)) * wall;
+  }
+}

test/Makefile.am

@@ -31,6 +31,7 @@ EXTRA_DIST  = tap-driver.sh test_tps_splitcomm.py soln_differ inputs meshes lte-
 	      ref_solns/lequere-varmu/*.h5 \
 	      ref_solns/taylor-couette/*.h5 \
 	      ref_solns/pipe/*.h5 \
+	      ref_solns/aveLoMach/*.h5 \
 		      vpath.sh die.sh count_gpus.sh sniff_mpirun.sh \
 		      cyl3d.gpu.test cyl3d.mflow.gpu.test wedge.gpu.test \
 	          averaging.gpu.test cyl3d.test cyl3d.gpu.python.test cyl3d.mflow.test cyl3d.dtconst.test \
@@ -45,7 +46,7 @@ EXTRA_DIST  = tap-driver.sh test_tps_splitcomm.py soln_differ inputs meshes lte-
 		      tabulated.test lte_mixture.test distance_fcn.test \
               sgsSmag.test sgsSigma.test heatEq.test sponge.test plate.test pipe.mix.test lte2noneq-restart.test \
               coupled-3d.interface.test plasma.axisym.test plasma.axisym.lte1d.test \
-	      lomach-flow.test lomach-lequere.test interpInlet.test sgsLoMach.test autoPeriodic.test
+	      lomach-flow.test lomach-lequere.test interpInlet.test sgsLoMach.test autoPeriodic.test aveLoMach.test
 
 TESTS = vpath.sh
 XFAIL_TESTS =
@@ -110,7 +111,8 @@ TESTS += cyl3d.test \
 	 lomach-flow.test \
 	 lomach-lequere.test \
          interpInlet.test \
-         autoPeriodic.test
+         autoPeriodic.test \
+         aveLoMach.test
 
 if PYTHON_ENABLED
 TESTS += cyl3d.python.test \

test/aveLoMach.test

@@ -0,0 +1,30 @@
+#!./bats
+# -*- mode: sh -*-
+
+TEST="aveLoMach"
+RUNFILE="inputs/input.aveLoMach.ini"
+EXE="../src/tps"
+RESTART="ref_solns/aveLoMach/restart_output.sol.h5"
+
+setup() {
+    SOLN_FILE=restart_output.sol.h5
+    REF_FILE=ref_solns/aveLoMach/restart_output.sol.h5
+}
+
+@test "[$TEST] check for input file $RUNFILE" {
+    test -s $RUNFILE
+}
+
+@test "[$TEST] run tps with input -> $RUNFILE" {
+    rm -rf output/*
+    rm -f $SOLN_FILE
+    run $EXE --runFile $RUNFILE
+    [[ ${status} -eq 0 ]]
+    test -s $SOLN_FILE
+}
+
+@test "[$TEST] verify tps output with input -> $RUNFILE" {
+    test -s $SOLN_FILE
+    test -s $REF_FILE
+    h5diff -r --delta=1e-12 $SOLN_FILE $REF_FILE meanVel
+}

test/inputs/input.aveLoMach.ini

@@ -0,0 +1,83 @@
+[solver]
+type = loMach
+
+[loMach]
+flow-solver = tomboulides
+thermo-solver = calorically-perfect
+mesh = meshes/flatBox.msh
+order = 1
+maxIters = 100
+outputFreq = 100
+fluid = dry_air
+refLength = 1.0
+equation_system = navier-stokes
+enablePressureForcing = False
+enableGravity = False
+openSystem = True
+ambientPressure = 101326.
+sgsModel = none
+#sgsModel = smagorinsky
+sgsModelConstant = 0.09
+
+[loMach/calperfect]
+#viscosity-model = constant
+#density-model = constant
+#constant-visc/mu0 = 1.552e-5
+#constant-density = 1.1839
+linear-solver-rtol = 1.0e-12
+linear-solver-max-iter = 2000
+viscosity-model = sutherland
+sutherland/mu0 = 1.68e-5
+sutherland/T0 = 273.0
+sutherland/S0 = 110.4
+numerical-integ = false
+#Prandtl = 0.72
+Prandtl = 0.01
+
+[loMach/tomboulides]
+linear-solver-rtol = 1.0e-12
+linear-solver-max-iter = 2000
+numerical-integ = false
+ic = tgv2d_uniform
+
+[io]
+outdirBase = output
+#enableRestart = True                                                                                                                                                         
+#restartMode = variableP
+
+[time]
+integrator = curlcurl
+cfl = 0.5
+#dt_fixed = 1.0e-4
+dt_initial = 1.0e-4
+dtFactor = 0.01
+
+[initialConditions]
+rho = 1.1839
+rhoU = 0.
+rhoV = 0.
+rhoW = 0.
+temperature = 298.15
+pressure = 101325.0
+
+[averaging]
+enableContinuation = False
+#enableContinuation = True
+#restartRMS = True
+saveMeanHist = True
+startIter = 1
+sampleFreq = 10
+
+[boundaryConditions]
+numWalls = 0
+numInlets = 0
+numOutlets = 0
+
+[periodicity]
+enablePeriodic = True
+periodicX = True
+periodicY = True
+periodicZ = True
+#xTrans = 1.0
+#yTrans = 1.0
+#zTrans = 0.1