Create alpha and beta that depend on the state. This gets the 3D

sphere const vel, rho case working. For the nbhd in question, the MEB state redist equations create a quadratic equation for beta. The code checks for this type of nbhd and then solved the quadratic. It does not treat another possible challenging nbhd that could exist but doesn't seem to arise in our 3D sphere set up (at least for a while). This more complex nbhd would require solving a cubic equation. Code needs a lot of clean up.
AMReX-Fluids · Jun 20, 2024 · afcc5c8 · afcc5c8
1 parent edbc077
commit afcc5c8
Show file tree

Hide file tree

Showing 3 changed files with 378 additions and 190 deletions.
diff --git a/src/incflo_advance.cpp b/src/incflo_advance.cpp
@@ -108,15 +108,15 @@ void incflo::Advance(Real orig_mass, Real& prev_mass, Real& prev_vol)
 
     amrex::Print() << "Pred:Change in domain volume  " << (vol - prev_vol) << std::endl;
     // For plane geoms:
-    // amrex::Print() << "Pred:Change in domain mass   " << (sum - (prev_mass - 1.*dx[1]*16*m_dt*2.)) << std::endl;
-    // amrex::Print() << "Pred:Change speed   " << (sum - prev_mass) / (1.* dx[1]*16 *m_dt) << std::endl;
-    // amrex::Print() << "Pred:Change error   " << (sum - prev_mass - vol + prev_vol)/ (vol-prev_vol) << std::endl;
+    amrex::Print() << "Pred:Change in domain mass   " << (sum - (prev_mass - 1.*dx[1]*16*m_dt*2.)) << std::endl;
+    amrex::Print() << "Pred:Change speed   " << (sum - prev_mass) / (1.* dx[1]*16 *m_dt) << std::endl;
+    amrex::Print() << "Pred:Change error   " << (sum - prev_mass - vol + prev_vol)/ (vol-prev_vol) << std::endl;
 
     // amrex::Print() << "Pred:Change over time in last time step divided by time and area " << (sum - prev_mass) / m_dt / 1.6 <<
     //                   " using " << sum << " " << prev_mass << " " << m_dt << std::endl;
 
     // For cyl :
-    amrex::Print() << "Pred:Change in domain mass 2D " << (sum - (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.)) << std::endl;
+//    amrex::Print() << "Pred:Change in domain mass 2D " << (sum - (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.)) << std::endl;
     // amrex::Print() << "Pred:Change speed   " << (sum - prev_mass - (vol-prev_vol)*(2.45+1.625)/2.) / ((2.45-1.)* dx[1]*64 *m_dt) << std::endl;
 
     // For sphere :
@@ -157,30 +157,31 @@ void incflo::Advance(Real orig_mass, Real& prev_mass, Real& prev_vol)
     amrex::Print() << "Corr:Change in domain volume " << (vol - prev_vol) <<std::endl;
     static Real verr = 0.;
     verr += (vol-prev_vol);
-    amrex::Print() << "Corr:Sum volume error " << verr << std::endl;
+    amrex::Print() << "Corr:Sum (volume error) " << verr << std::endl;
     //amrex::Print() << "Corr:Change in domain mass   " << (sum - prev_mass) << std::endl;
-    // amrex::Print() << "Corr:Change in domain mass   " << (sum - (prev_mass - 1.*dx[1]*16*m_dt*2.)) << std::endl;
-    // amrex::Print() << "Corr:Change speed   " << (sum - prev_mass) / (1.* dx[1]*16 *m_dt) << std::endl;
-
-
-    amrex::Print() << "Corr:Change in domain mass 2D " << (sum - (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.)) << std::endl;
-    // This only works if we start at time 0
-    // Real m_calc = (orig_mass + (2.45-1.)*dx[1]*ncelly*m_t_new[0]*1.);
-    // amrex::Print() << "Corr:Mass conservation error 2D " << (sum - m_calc) << std::endl;
-    Real m_calc = (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.);
-    static Real err = 0.;
-    err += (sum - m_calc);
-    amrex::Print() << "Corr:Sum conservation error " << err << std::endl;
-    //if ( std::abs( (vol - prev_vol) -  (sum - prev_mass) ) > 1e-14 ) Abort("Vol diff and dens diff are off");
-    // amrex::Print() << "Corr:Sum of mass at time = " << m_cur_time+m_dt << " " << sum << " " << std::endl;
-//  amrex::Print() << "Change over time divided by time and area " << (sum - orig_mass) / (m_cur_time+m_dt) / 1.6 <<
-//                    " using " << sum << " " << orig_mass << " " << m_cur_time+m_dt << std::endl;
-    // amrex::Print() << "Corr:Change over time in last time step divided by time and area " << (sum - prev_mass) / m_dt / 1.6 <<
-    //                   " using " << sum << " " << prev_mass << " " << m_dt << std::endl;
+    amrex::Print() << "Corr:Change in domain mass   " << (sum - (prev_mass - 1.*dx[1]*16*m_dt*2.)) << std::endl;
+    amrex::Print() << "Corr:Change speed   " << (sum - prev_mass) / (1.* dx[1]*16 *m_dt) << std::endl;
+
+
+//     amrex::Print() << "Corr:Change in domain mass 2D " << (sum - (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.)) << std::endl;
+//     // This only works if we start at time 0
+//     // Real m_calc = (orig_mass + (2.45-1.)*dx[1]*ncelly*m_t_new[0]*1.);
+//     // amrex::Print() << "Corr:Mass conservation error 2D " << (sum - m_calc) << std::endl;
+//     Real m_calc = (prev_mass + (2.45-1.)*dx[1]*ncelly*m_dt*1.);
+//     static Real err = 0.;
+//     err += (sum - m_calc);
+//     amrex::Print() << "Corr:Sum conservation error " << err << std::endl;
+//     //if ( std::abs( (vol - prev_vol) -  (sum - prev_mass) ) > 1e-14 ) Abort("Vol diff and dens diff are off");
+//     // amrex::Print() << "Corr:Sum of mass at time = " << m_cur_time+m_dt << " " << sum << " " << std::endl;
+// //  amrex::Print() << "Change over time divided by time and area " << (sum - orig_mass) / (m_cur_time+m_dt) / 1.6 <<
+// //                    " using " << sum << " " << orig_mass << " " << m_cur_time+m_dt << std::endl;
+//     // amrex::Print() << "Corr:Change over time in last time step divided by time and area " << (sum - prev_mass) / m_dt / 1.6 <<
+//     //                   " using " << sum << " " << prev_mass << " " << m_dt << std::endl;
+
+#endif
 
     prev_mass = sum;
     prev_vol = vol;
-#endif
 
     // Stop timing current time step
     Real end_step = static_cast<Real>(ParallelDescriptor::second()) - strt_step;