root/src/secOrd_op_brusselator.cpp

DEFINITIONS

1. Brusselator
2. Brusselator
3. solve
4. mass_matrix
5. massMatrix_constOne
6. add_massMatrixConstOne_to
7. assign_firstArg_quadMassMatrix
8. assign_secondArg_quadMassMatrix
9. Data
10. Data
11. Data

   1 /*! \file secOrd_op_brusselator.cpp
   2     \brief Implementation of secOrd_op_brusselator.h
   3 
   4      Revision history
   5      --------------------------------------------------
   6 
   7           Revised by Christian Power May 2016
   8           Originally written by Christian Power
   9                (power22c@gmail.com) Februar 2016
  10 
  11      \author Christian Power 
  12      \date 18. May 2016
  13      \copyright Copyright (c) 2016 Christian Power. All rights reserved.
  14  */
  15 
  16 #include <stdexcept>
  17 #include <config.h>
  18 #include <dune/fem/solver/cginverseoperator.hh>
  19 #include <dune/fem/solver/oemsolver.hh>
  20 #include <dune/fem/operator/linear/istloperator.hh>
  21 #include <dune/fem/solver/istlsolver.hh>
  22 #include "secOrd_op_brusselator.h"
  23 #include "secOrd_op_brusselator_impl.h"
  24 #include "io_parameter.h"
  25 #include "grid.h"
  26 
  27 using namespace std;
  28 using FE_function = Esfem::Grid::Scal_FEfun::Dune_FEfun;
  29 using CG_method = Dune::Fem::CGInverseOperator<FE_function>;
  30 
  31 using Linear_operator = Dune::Fem::ISTLLinearOperator<FE_function, FE_function>;
  32 using GMRes
  33 = Dune::Fem::ISTLGMResOp<FE_function, Linear_operator>;
  34 // = Dune::Fem::ISTLCGOp<FE_function, LinearOperatorType>;
  35 
  36 struct Esfem::SecOrd_op::Brusselator::Data{
  37   bool owns;
  38   FE_function* fef1_ptr {nullptr};
  39   FE_function* fef2_ptr {nullptr};
  40   const FE_function& fef1_ref;
  41   const FE_function& fef2_ref;
  42   FE_function tmp_var;
  43   Brusselator_op bruss_op;
  44   Linear_operator bruss_matrix;
  45   CG_method bruss_cg;
  46   GMRes bruss_gmres;
  47   explicit Data(const Io::Parameter&, const Grid::Grid_and_time&,
  48                 const Growth);
  49   explicit Data(const Io::Parameter&, const Grid::Grid_and_time&,
  50                 const Growth, const FE_function&, const FE_function&);
  51   ~Data();
  52   Data(const Data&) = delete;
  53   Data& operator=(const Data&) = delete;
  54 };
  55 
  56 Esfem::SecOrd_op::Brusselator::
  57 Brusselator(const Io::Parameter& p, const Grid::Grid_and_time& gt, const Growth type)
  58 try : d_ptr {make_unique<Data>(p, gt, type)}
  59 {}
  60 catch(const std::exception&){
  61   std::throw_with_nested(std::logic_error
  62                          {"Error in constructor of Brusselator."});
  63  }
  64  catch(...){
  65    throw std::logic_error {"Unkown error in constructor of Brusselator."};
  66  }
  67 
  68 Esfem::SecOrd_op::Brusselator::
  69 Brusselator(const Io::Parameter& p, const Grid::Grid_and_time& gt, const Growth type,
  70             const Grid::Scal_FEfun& fef1, const Grid::Scal_FEfun& fef2)
  71 try :d_ptr {make_unique<Data>(p, gt, type, fef1, fef2)}
  72 {}
  73  catch(...){
  74    throw_with_nested(std::logic_error {"Error in constructor of Brusselator."});
  75  }
  76 Esfem::SecOrd_op::Brusselator::~Brusselator() = default;
  77 
  78 void Esfem::SecOrd_op::Brusselator::
  79 solve(const Grid::Scal_FEfun& rhs, Grid::Scal_FEfun& lhs) const{
  80   const FE_function& rhs_ref = rhs;
  81   FE_function& lhs_ref = lhs;
  82   d_ptr -> bruss_cg(rhs_ref, lhs_ref);
  83   // d_ptr -> bruss_op.jacobian(rhs_ref, d_ptr -> bruss_matrix);
  84   // d_ptr -> bruss_gmres(rhs_ref, lhs_ref);
  85 }
  86 void Esfem::SecOrd_op::Brusselator::
  87 mass_matrix(const Grid::Scal_FEfun& rhs, Grid::Scal_FEfun& lhs) const{
  88   const FE_function& rhs_ref = rhs;
  89   FE_function& lhs_ref = lhs;
  90   d_ptr -> bruss_op.mass_matrix(rhs_ref, lhs_ref);
  91 }
  92 void Esfem::SecOrd_op::Brusselator::massMatrix_constOne(Grid::Scal_FEfun& fef) const{
  93   FE_function& fef_ref = fef;
  94   d_ptr -> bruss_op.massMatrix_constOne(fef_ref);
  95 }
  96 void Esfem::SecOrd_op::Brusselator::
  97 add_massMatrixConstOne_to(Grid::Scal_FEfun& fef) const{
  98   FE_function& fef_ref = fef;
  99   d_ptr -> bruss_op.massMatrix_constOne(d_ptr -> tmp_var);
 100   fef_ref.axpy(1., d_ptr -> tmp_var);
 101 }
 102 void Esfem::SecOrd_op::Brusselator::
 103 assign_firstArg_quadMassMatrix(const Grid::Scal_FEfun& fef){
 104   if(d_ptr -> owns){
 105     const FE_function& fef_ref = fef;
 106     d_ptr -> fef1_ptr -> assign(fef_ref);
 107   }
 108   else
 109     throw std::logic_error
 110     {"Error in Brusselator::assign_firstArg_quadMassMatrix().  "
 111         "*this has used a wrong constructor."};
 112 }
 113 void Esfem::SecOrd_op::Brusselator::
 114 assign_secondArg_quadMassMatrix(const Grid::Scal_FEfun& fef){
 115   if(d_ptr -> owns){
 116     const FE_function& fef_ref = fef;
 117     d_ptr -> fef2_ptr -> assign(fef_ref);
 118   }
 119   else
 120     throw std::logic_error
 121     {"Error in Brusselator::assign_secondArg_quadMassMatrix().  "
 122         "*this has used a wrong constructor."};
 123 }
 124 
 125 // testing
 126 void Esfem::SecOrd_op::Brusselator::
 127 operator()(const Grid::Scal_FEfun& rhs, Grid::Scal_FEfun& lhs) const{
 128   const FE_function& rhs_ref = rhs;
 129   FE_function& lhs_ref = lhs;
 130   d_ptr -> bruss_op(rhs,lhs);
 131   // d_ptr -> bruss_op.jacobian(rhs_ref, d_ptr -> bruss_matrix);
 132   // d_ptr -> bruss_matrix(rhs_ref, lhs_ref);
 133 }
 134 
 135 // ----------------------------------------------------------------------
 136 // Internal implementation
 137 
 138 Esfem::SecOrd_op::Brusselator::Data::
 139 Data(const Io::Parameter& p, const Grid::Grid_and_time& gt, const Growth type)
 140   :owns {true},
 141    fef1_ptr {new FE_function {"fef1", gt.fe_space()}},
 142    fef2_ptr {new FE_function {"fef2", gt.fe_space()}},
 143    fef1_ref {*fef1_ptr},
 144    fef2_ref {*fef2_ptr},
 145    tmp_var {"tmp_var", gt.fe_space()},
 146    bruss_op {p, gt, type, fef1_ref, fef2_ref},
 147    bruss_matrix {"assempled elliptic operator", gt.fe_space(), gt.fe_space()},
 148    bruss_cg {bruss_op, p.eps(), p.eps()},
 149    bruss_gmres {bruss_matrix, p.eps(), p.eps()}
 150   {}
 151 Esfem::SecOrd_op::Brusselator::Data::
 152 Data(const Io::Parameter& p, const Grid::Grid_and_time& gt,
 153      const Growth type, const FE_function& fef1,
 154      const FE_function& fef2)
 155   :owns {false}, fef1_ref {fef1}, fef2_ref {fef2},
 156    tmp_var {"tmp_var", gt.fe_space()},
 157    bruss_op {p, gt, type, fef1_ref, fef2_ref},
 158    bruss_matrix {"assempled elliptic operator", gt.fe_space(), gt.fe_space()},
 159    bruss_cg {bruss_op, p.eps(), p.eps()},
 160    bruss_gmres {bruss_matrix, p.eps(), p.eps()}
 161   {}
 162 Esfem::SecOrd_op::Brusselator::Data::
 163 ~Data(){
 164   if(owns){
 165     delete fef1_ptr;
 166     delete fef2_ptr;
 167     fef1_ptr = fef2_ptr = nullptr;
 168   }
 169 }
 170 
 171 /*! Log:
 172  */