omega for self-gravity can be set in the input file

PrincetonUniversity · Mar 2, 2024 · 640a7b4 · 640a7b4
1 parent dbfad59
commit 640a7b4
Show file tree

Hide file tree

Showing 4 changed files with 4 additions and 13 deletions.
diff --git a/src/crdiffusion/crdiffusion.cpp b/src/crdiffusion/crdiffusion.cpp
@@ -45,7 +45,6 @@ CRDiffusion::CRDiffusion(MeshBlock *pmb, ParameterInput *pin) :
   Dpara_ = pin->GetReal("crdiffusion", "Dpara");
   Dperp_ = pin->GetReal("crdiffusion", "Dperp");
   Lambda_ = pin->GetReal("crdiffusion", "Lambda");
-  zeta_factor_ = pin->GetReal("crdiffusion", "zeta_factor");
 
   output_defect = pin->GetOrAddBoolean("crdiffusion", "output_defect", false);
   if (output_defect)

diff --git a/src/crdiffusion/mg_crdiffusion.cpp b/src/crdiffusion/mg_crdiffusion.cpp
@@ -265,13 +265,6 @@ void MGCRDiffusionDriver::Solve(int stage, Real dt) {
 
   crtlist_->DoTaskListOneStage(pmy_mesh_, stage);
 
-  for (auto itr = vmg_.begin(); itr < vmg_.end(); itr++) {
-    Multigrid *pmg = *itr;
-    CRDiffusion *pcrdiff = pmg->pmy_block_->pcrdiff;
-    Hydro *phydro = pmg->pmy_block_->phydro;
-    pcrdiff->CalculateIonizationRate(phydro->w);
-  }
-
   return;
 }
 

diff --git a/src/gravity/mg_gravity.cpp b/src/gravity/mg_gravity.cpp
@@ -43,6 +43,7 @@ MGGravityDriver::MGGravityDriver(Mesh *pm, ParameterInput *pin)
                       pm->MGGravitySourceMaskFunction_, pm->MGGravitySourceMaskFunction_,
                       1, 0, 0) {
   four_pi_G_ = pmy_mesh_->four_pi_G_;
+  omega_ = pin->GetOrAddReal("gravity", "omega", 1.15);
   eps_ = pin->GetOrAddReal("gravity", "threshold", -1.0);
   niter_ = pin->GetOrAddInteger("gravity", "niteration", -1);
   ffas_ = pin->GetOrAddBoolean("gravity", "fas", ffas_);
@@ -226,7 +227,7 @@ void MGGravity::Smooth(AthenaArray<Real> &u, const AthenaArray<Real> &src,
   if (rlev <= 0) dx = rdx_*static_cast<Real>(1<<(-rlev));
   else           dx = rdx_/static_cast<Real>(1<<rlev);
   Real dx2 = SQR(dx);
-  Real isix = omega_/6.0;
+  Real isix = static_cast<MGGravityDriver*>(pmy_driver_)->omega_/6.0;
   color ^= pmy_driver_->coffset_;
 
 #pragma omp parallel for num_threads(pmy_driver_->nthreads_) if (th && (ku-kl) >= minth_)

diff --git a/src/gravity/mg_gravity.hpp b/src/gravity/mg_gravity.hpp
@@ -41,9 +41,6 @@ class MGGravity : public Multigrid {
   void CalculateFASRHS(AthenaArray<Real> &def, const AthenaArray<Real> &src,
                 const AthenaArray<Real> &coeff, const AthenaArray<Real> &matrix,
                 int rlev, int il, int iu, int jl, int ju, int kl, int ku, bool th) final;
-
- private:
-  static constexpr Real omega_ = 1.15;
 };
 
 
@@ -57,8 +54,9 @@ class MGGravityDriver : public MultigridDriver {
   void Solve(int stage, Real dt = 0.0) final;
   void ProlongateOctetBoundariesFluxCons(AthenaArray<Real> &dst,
                  AthenaArray<Real> &cbuf, const AthenaArray<bool> &ncoarse) final;
+  friend class MGGravity;
  private:
-  Real four_pi_G_;
+  Real four_pi_G_, omega_;
   GravityBoundaryTaskList *gtlist_;
 };