Update BOUT++

hermes-3.cxx

@@ -121,8 +121,8 @@ class DecayLengthBoundary : public BoundaryOp {
 
     // Get cell radial length
     Coordinates* coord = mesh->getCoordinates();
-    auto dx = coord->dx;
-    auto g11 = coord->g11;
+    auto dx = coord->dx();
+    auto g11 = coord->g11();
     Coordinates::FieldMetric dr =
         dx / sqrt(g11); // cell radial length. dr = dx/(Bpol * R) and g11 = (Bpol*R)**2
 
@@ -254,28 +254,7 @@ int Hermes::init(bool restarting) {
             .doc("Normalise input metric tensor? (assumes input is in SI units)")
             .withDefault<bool>(true)) {
       Coordinates* coord = mesh->getCoordinates();
-      // To use non-orthogonal metric
-      // Normalise
-      coord->dx /= rho_s0 * rho_s0 * Bnorm;
-      coord->Bxy /= Bnorm;
-      // Metric is in grid file - just need to normalise
-      coord->g11 /= SQ(Bnorm * rho_s0);
-      coord->g22 *= SQ(rho_s0);
-      coord->g33 *= SQ(rho_s0);
-      coord->g12 /= Bnorm;
-      coord->g13 /= Bnorm;
-      coord->g23 *= SQ(rho_s0);
-
-      coord->J *= Bnorm / rho_s0;
-
-      coord->g_11 *= SQ(Bnorm * rho_s0);
-      coord->g_22 /= SQ(rho_s0);
-      coord->g_33 /= SQ(rho_s0);
-      coord->g_12 *= Bnorm;
-      coord->g_13 *= Bnorm;
-      coord->g_23 /= SQ(rho_s0);
-
-      coord->geometry(); // Calculate other metrics
+      coord->normaliseMetric(rho_s0, Bnorm);
     }
   }
 

include/braginskii_conduction.hxx

@@ -79,7 +79,7 @@ private:
   ///
   void transform_impl(GuardedOptions& state) override;
 
-  std::string conduction_method{"orginal"};
+  bout::ConductionMethod conduction_method{bout::ConductionMethod::Original};
 };
 
 namespace {

include/div_ops.hxx

@@ -43,7 +43,7 @@
 /*!
  * Diffusion in index space
  *
- * Similar to using Div_par_diffusion(SQ(mesh->dy)*mesh->g_22, f)
+ * Similar to using Div_par_diffusion(SQ(mesh->dy())*mesh->g_22(), f)
  *
  * @param[in] The field to be differentiated
  * @param[in] bndry_flux  Are fluxes through the boundary calculated?

include/integrate.hxx

@@ -75,7 +75,7 @@ auto cellAverage(Function func, const RegionType& region) {
     result.allocate();
 
     // Get the coordinate Jacobian
-    auto J = result.getCoordinates()->J;
+    auto J = result.getCoordinates()->J();
     BOUT_FOR(i, region) {
       // Offset indices
       auto yp = i.yp();

include/neutral_mixed.hxx

@@ -82,7 +82,7 @@ private:
   bool output_ddt; ///< Save time derivatives?
   bool diagnose;   ///< Save additional diagnostics?
 
-  std::string conduction_method{"orginal"};
+  bout::ConductionMethod conduction_method{bout::ConductionMethod::Original};
 
   // Flow diagnostics
   Field3D pf_adv_perp_xlow, pf_adv_perp_ylow, pf_adv_par_ylow;

include/polarisation_drift.hxx

@@ -69,7 +69,7 @@ struct PolarisationDrift : public NamedComponent<PolarisationDrift> {
 private:
   std::unique_ptr<Laplacian> phiSolver; // Laplacian solver in X-Z
 
-  Coordinates::FieldMetric Bsq; // Cached SQ(coord->Bxy)
+  Coordinates::FieldMetric Bsq; // Cached SQ(coord->Bxy())
 
   // Diagnostic outputs
   bool diagnose;       ///< Save diagnostic outputs?

include/rate_helper.hxx

@@ -156,7 +156,7 @@ struct RateHelper {
     // Populate cell_props differently according to the rate function type
     std::visit(
         [this, &cell_props, &do_averaging, &result](auto&& rate_calc_func) {
-          auto J = result.rate.getCoordinates()->J;
+          auto J = result.rate.getCoordinates()->J();
           BOUT_FOR(i, region) {
             // Get densities for this cell
             collect_densities(i, cell_props, do_averaging);

include/relax_potential.hxx

@@ -116,7 +116,7 @@ private:
   BoutReal phi_boundary_last_update; ///< Time when last updated
   bool phi_core_averagey;            ///< Average phi core boundary in Y?
 
-  Coordinates::FieldMetric Bsq;           ///< SQ(coord->Bxy)
+  Coordinates::FieldMetric Bsq;           ///< SQ(coord->Bxy())
   VectorMetric Curlb_B;                   ///< Curvature vector Curl(b/B)
   BoutReal hyper_z;                       ///< Hyper-viscosity in Z
   Coordinates::FieldMetric viscosity;     ///< Perpendicular Kinematic viscosity

include/scale_timederivs.hxx

@@ -40,7 +40,7 @@ private:
   void transform_impl(GuardedOptions& state) override {
 
     auto* coord = bout::globals::mesh->getCoordinates();
-    auto dl2 = coord->g_22 * SQ(coord->dy);
+    auto dl2 = coord->g_22() * SQ(coord->dy());
 
     // Scale by parallel heat conduction CFL timescale
     auto Te = get<Field3D>(state["species"]["e"]["temperature"]);

include/vorticity.hxx

@@ -127,7 +127,7 @@ private:
   bool split_n0;                        // Split phi into n=0 and n!=0 components
   std::unique_ptr<LaplaceXY> laplacexy; // Laplacian solver in X-Y (n=0)
 
-  Coordinates::FieldMetric Bsq;       // SQ(coord->Bxy)
+  Coordinates::FieldMetric Bsq;       // SQ(coord->Bxy())
   VectorMetric Curlb_B;               // Curvature vector Curl(b/B)
   BoutReal hyper_z;                   ///< Hyper-viscosity in Z
   Coordinates::FieldMetric viscosity; ///< Kinematic viscosity

src/braginskii_conduction.cxx

@@ -114,8 +114,7 @@ BraginskiiConduction::BraginskiiConduction(const std::string& name, Options& all
   }
   substitutePermissions("sp", species);
 
-  conduction_method =
-      alloptions["conduction_method"].withDefault<std::string>(conduction_method);
+  conduction_method = alloptions["conduction_method"].withDefault(conduction_method);
 }
 
 void BraginskiiConduction::transform_impl(GuardedOptions& state) {

src/braginskii_electron_viscosity.cxx

@@ -50,7 +50,7 @@ void BraginskiiElectronViscosity::transform_impl(GuardedOptions& state) {
   const Field3D V = get<Field3D>(species["velocity"]);
 
   Coordinates* coord = P.getCoordinates();
-  const Field3D Bxy = coord->Bxy;
+  const Field3D Bxy = coord->Bxy();
   const Field3D sqrtB = sqrt(Bxy);
 
   // Parallel electron viscosity

src/braginskii_ion_viscosity.cxx

@@ -120,7 +120,7 @@ BraginskiiIonViscosity::BraginskiiIonViscosity(const std::string& name,
     Curlb_B.y *= SQ(Lnorm);
     Curlb_B.z *= SQ(Lnorm);
 
-    Curlb_B *= 2. / coord->Bxy;
+    Curlb_B *= 2. / coord->Bxy();
   }
   if (bounce_frequency) {
     const Options& units = alloptions["units"];
@@ -146,7 +146,7 @@ void BraginskiiIonViscosity::transform_impl(GuardedOptions& state) {
   GuardedOptions allspecies = state["species"];
 
   auto coord = mesh->getCoordinates();
-  const auto Bxy = coord->Bxy;
+  const auto Bxy = coord->Bxy();
   const auto sqrtB = sqrt(Bxy);
   const auto Grad_par_logB = Grad_par(log(Bxy));
 

src/classical_diffusion.cxx

@@ -12,7 +12,7 @@ ClassicalDiffusion::ClassicalDiffusion(std::string name, Options& alloptions, So
            readWrite("species:{all_species}:{output}")}) {
   Options& options = alloptions[name];
 
-  Bsq = SQ(bout::globals::mesh->getCoordinates()->Bxy);
+  Bsq = SQ(bout::globals::mesh->getCoordinates()->Bxy());
 
   diagnose =
       options["diagnose"].doc("Output additional diagnostics?").withDefault<bool>(false);

src/diamagnetic_drift.cxx

@@ -47,7 +47,7 @@ DiamagneticDrift::DiamagneticDrift(std::string name, Options& alloptions,
   Curlb_B.y *= SQ(Lnorm);
   Curlb_B.z *= SQ(Lnorm);
 
-  Curlb_B *= 2. / mesh->getCoordinates()->Bxy;
+  Curlb_B *= 2. / mesh->getCoordinates()->Bxy();
 
   // Set drift to zero through sheath boundaries.
   // Flux through those cell faces should be set by sheath.