Line data    Source code

       1              : //
       2              : // Class Poisson
       3              : //   Base class for solvers for the Poisson problem
       4              : //
       5              : 
       6              : #ifndef IPPL_POISSON_H
       7              : #define IPPL_POISSON_H
       8              : 
       9              : #include "Utility/ParameterList.h"
      10              : 
      11              : #include "Field/Field.h"
      12              : 
      13              : namespace ippl {
      14              : 
      15              :     template <typename FieldLHS, typename FieldRHS>
      16              :     class Poisson {
      17              :         constexpr static unsigned Dim = FieldLHS::dim;
      18              :         typedef typename FieldLHS::Mesh_t Mesh;
      19              :         typedef typename FieldLHS::Centering_t Centering;
      20              :         typedef typename Mesh::matrix_type Matrix_t;
      21              :         typedef Field<Matrix_t, Dim, Mesh, Centering> MField_t;
      22              : 
      23              :     public:
      24              :         using lhs_type = FieldLHS;
      25              :         using rhs_type = FieldRHS;
      26              :         using Tlhs     = typename FieldLHS::value_type;
      27              :         using Trhs     = typename FieldRHS::value_type;
      28              : 
      29              :         using grad_type = Field<Vector<Tlhs, Dim>, Dim, Mesh, Centering>;
      30              : 
      31              :         /*!
      32              :          * Represents the types of fields that should
      33              :          * be output by the solver
      34              :          */
      35              :         enum OutputType {
      36              :             SOL          = 0b01,
      37              :             GRAD         = 0b10,
      38              :             SOL_AND_GRAD = 0b11
      39              :         };
      40              : 
      41              :         /*!
      42              :          * Default constructor for electrostatic solvers;
      43              :          * desired output type defaults to solution only
      44              :          */
      45            0 :         Poisson()
      46            0 :             : grad_mp(nullptr) {
      47              :             static_assert(std::is_floating_point<Trhs>::value, "Not a floating point type");
      48            0 :             setDefaultParameters();
      49            0 :         }
      50              : 
      51            0 :         Poisson(lhs_type& lhs, rhs_type& rhs)
      52            0 :             : grad_mp(nullptr) {
      53            0 :             setLhs(lhs);
      54            0 :             setRhs(rhs);
      55              : 
      56              :             static_assert(std::is_floating_point<Trhs>::value, "Not a floating point type");
      57            0 :             setDefaultParameters();
      58            0 :         }
      59              : 
      60              :         /*!
      61              :          * Update one of the solver's parameters
      62              :          * @param key The parameter key
      63              :          * @param value The new value
      64              :          * @throw IpplException Fails if there is no existing parameter with the given key
      65              :          */
      66              :         template <typename T>
      67              :         void updateParameter(const std::string& key, const T& value) {
      68              :             params_m.update<T>(key, value);
      69              :         }
      70              : 
      71              :         /*!
      72              :          * Updates all solver parameters based on values in another parameter set
      73              :          * @param params Parameter list with updated values
      74              :          * @throw IpplException Fails if the provided parameter list includes keys not already
      75              :          * present
      76              :          */
      77              :         void updateParameters(const ParameterList& params) { params_m.update(params); }
      78              : 
      79              :         /*!
      80              :          * Merges another parameter set into the solver's parameters, overwriting
      81              :          * existing parameters in case of conflict
      82              :          * @param params Parameter list with desired values
      83              :          */
      84            0 :         void mergeParameters(const ParameterList& params) { params_m.merge(params); }
      85              : 
      86              :         /*!
      87              :          * Set the problem LHS
      88              :          * @param lhs Reference to problem LHS field
      89              :          */
      90            0 :         void setLhs(lhs_type& lhs) { lhs_mp = &lhs; }
      91              : 
      92              :         /*!
      93              :          * Set the problem RHS
      94              :          * @param rhs Reference to problem RHS field
      95              :          */
      96            0 :         virtual void setRhs(rhs_type& rhs) { rhs_mp = &rhs; }
      97              : 
      98              :         /*!
      99              :          * Get the Hessian matrix of the solution
     100              :          * @return Matrix field containing the Hessian of the lhs
     101              :          */
     102            0 :         virtual MField_t* getHessian() { return nullptr; }
     103              : 
     104              :         /*!
     105              :          * Set the field in which the gradient of the computed potential
     106              :          * should be stored
     107              :          * @param grad Reference to field in which to store the gradient
     108              :          */
     109            0 :         void setGradient(grad_type& grad) { grad_mp = &grad; }
     110              : 
     111              :         /*!
     112              :          * Solve the Poisson problem described by
     113              :          * -laplace(lhs) = rhs
     114              :          */
     115              :         virtual void solve() = 0;
     116              : 
     117            0 :         virtual ~Poisson() {}
     118              : 
     119              :     protected:
     120              :         ParameterList params_m;
     121              : 
     122              :         rhs_type* rhs_mp = nullptr;
     123              :         lhs_type* lhs_mp = nullptr;
     124              : 
     125              :         grad_type* grad_mp;
     126              : 
     127              :         /*!
     128              :          * Utility function for initializing a solver's default
     129              :          * parameters (to be overridden for each base class)
     130              :          */
     131            0 :         virtual void setDefaultParameters() { this->params_m.add("output_type", SOL); }
     132              :     };
     133              : }  // namespace ippl
     134              : 
     135              : #endif