Remove uses of global 'mesh' from Laplacian solvers

include/invert_laplace.hxx

@@ -105,22 +105,22 @@ const int INVERT_OUT_RHS = 32768; ///< Use input value in RHS at outer boundary
 /// Base class for Laplacian inversion
 class Laplacian {
 public:
-  Laplacian(Options *options = nullptr, const CELL_LOC loc = CELL_CENTRE);
+  Laplacian(Options *options = nullptr, const CELL_LOC loc = CELL_CENTRE, Mesh* mesh_in = mesh);
   virtual ~Laplacian() {}
 
   /// Set coefficients for inversion. Re-builds matrices if necessary
   virtual void setCoefA(const Field2D &val) = 0;
   virtual void setCoefA(const Field3D &val) { setCoefA(DC(val)); }
   virtual void setCoefA(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefA(f);
   }
 
   virtual void setCoefC(const Field2D &val) = 0;
   virtual void setCoefC(const Field3D &val) { setCoefC(DC(val)); }
   virtual void setCoefC(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefC(f);
   }
@@ -130,7 +130,7 @@ public:
   }
   virtual void setCoefC1(const Field3D &val) { setCoefC1(DC(val)); }
   virtual void setCoefC1(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefC1(f);
   }
@@ -140,31 +140,31 @@ public:
   }
   virtual void setCoefC2(const Field3D &val) { setCoefC2(DC(val)); }
   virtual void setCoefC2(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefC2(f);
   }
 
   virtual void setCoefD(const Field2D &val) = 0;
   virtual void setCoefD(const Field3D &val) { setCoefD(DC(val)); }
   virtual void setCoefD(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefD(f);
   }
 
   virtual void setCoefEx(const Field2D &val) = 0;
   virtual void setCoefEx(const Field3D &val) { setCoefEx(DC(val)); }
   virtual void setCoefEx(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefEx(f);
   }
 
   virtual void setCoefEz(const Field2D &val) = 0;
   virtual void setCoefEz(const Field3D &val) { setCoefEz(DC(val)); }
   virtual void setCoefEz(BoutReal r) {
-    Field2D f(r);
+    Field2D f(r, localmesh);
     f.setLocation(location);
     setCoefEz(f);
   }
@@ -192,7 +192,7 @@ public:
    * 
    * @param[in] opt  The options section to use. By default "laplace" will be used
    */
-  static Laplacian *create(Options *opt = nullptr, const CELL_LOC loc = CELL_CENTRE);
+  static Laplacian *create(Options *opt = nullptr, const CELL_LOC loc = CELL_CENTRE, Mesh *mesh_in = mesh);
   static Laplacian* defaultInstance(); ///< Return pointer to global singleton
 
   static void cleanup(); ///< Frees all memory
@@ -227,7 +227,10 @@ protected:
                     const Field2D *a, const Field2D *ccoef, 
                     const Field2D *d,
                     bool includeguards=true);
-  CELL_LOC location;
+  CELL_LOC location;   ///< staggered grid location of this solver
+  Mesh* localmesh;     ///< Mesh object for this solver
+  Coordinates* coords; ///< Coordinates object, so we only have to call
+                       ///  localmesh->getCoordinates(location) once
 private:
   /// Singleton instance
   static Laplacian *instance;

include/invert_parderiv.hxx

@@ -64,15 +64,16 @@ public:
    * with pure virtual members, so can't be created directly.
    * To create an InvertPar object call the create() static function.
    */ 
-  InvertPar(Options *UNUSED(opt)) {}
+  InvertPar(Options *UNUSED(opt), Mesh *mesh_in = mesh)
+    : localmesh(mesh_in) {}
   virtual ~InvertPar() {}
 
   /*!
    * Create an instance of InvertPar
    * 
    * Note: For consistency this should be renamed "create" and take an Options* argument
    */
-  static InvertPar* Create();
+  static InvertPar* Create(Mesh *mesh_in = mesh);
 
   /*!
    * Solve the system of equations
@@ -100,35 +101,38 @@ public:
    */
   virtual void setCoefA(const Field2D &f) = 0;
   virtual void setCoefA(const Field3D &f) {setCoefA(DC(f));}
-  virtual void setCoefA(BoutReal f) {setCoefA(Field2D(f));}
+  virtual void setCoefA(BoutReal f) {setCoefA(Field2D(f, localmesh));}
 
   /*!
    * Set the Grad2_par2 coefficient B
    */ 
   virtual void setCoefB(const Field2D &f) = 0;
   virtual void setCoefB(const Field3D &f) {setCoefB(DC(f));}
-  virtual void setCoefB(BoutReal f) {setCoefB(Field2D(f));}
+  virtual void setCoefB(BoutReal f) {setCoefB(Field2D(f, localmesh));}
 
   /*!
    * Set the D2DYDZ coefficient C
    */
   virtual void setCoefC(const Field2D &f) = 0;
   virtual void setCoefC(const Field3D &f) {setCoefB(DC(f));}
-  virtual void setCoefC(BoutReal f) {setCoefB(Field2D(f));}
+  virtual void setCoefC(BoutReal f) {setCoefB(Field2D(f, localmesh));}
 
   /*!
    * Set the D2DZ2 coefficient D
    */ 
   virtual void setCoefD(const Field2D &f) = 0;
   virtual void setCoefD(const Field3D &f) {setCoefB(DC(f));}
-  virtual void setCoefD(BoutReal f) {setCoefB(Field2D(f));}
+  virtual void setCoefD(BoutReal f) {setCoefB(Field2D(f, localmesh));}
 
   /*!
    * Set the DDY coefficient E
    */
   virtual void setCoefE(const Field2D &f) = 0;
   virtual void setCoefE(const Field3D &f) {setCoefB(DC(f));}
-  virtual void setCoefE(BoutReal f) {setCoefB(Field2D(f));}
+  virtual void setCoefE(BoutReal f) {setCoefB(Field2D(f, localmesh));}
+
+protected:
+  Mesh* localmesh; ///< Mesh object for this solver
 
 private:
 };