dune-fem 2.8.0
Loading...
Searching...
No Matches
femscheme.hh
Go to the documentation of this file.
1/**************************************************************************
2
3 The dune-fem module is a module of DUNE (see www.dune-project.org).
4 It is based on the dune-grid interface library
5 extending the grid interface by a number of discretization algorithms
6 for solving non-linear systems of partial differential equations.
7
8 Copyright (C) 2003 - 2015 Robert Kloefkorn
9 Copyright (C) 2003 - 2010 Mario Ohlberger
10 Copyright (C) 2004 - 2015 Andreas Dedner
11 Copyright (C) 2005 Adrian Burri
12 Copyright (C) 2005 - 2015 Mirko Kraenkel
13 Copyright (C) 2006 - 2015 Christoph Gersbacher
14 Copyright (C) 2006 - 2015 Martin Nolte
15 Copyright (C) 2011 - 2015 Tobias Malkmus
16 Copyright (C) 2012 - 2015 Stefan Girke
17 Copyright (C) 2013 - 2015 Claus-Justus Heine
18 Copyright (C) 2013 - 2014 Janick Gerstenberger
19 Copyright (C) 2013 Sven Kaulman
20 Copyright (C) 2013 Tom Ranner
21 Copyright (C) 2015 Marco Agnese
22 Copyright (C) 2015 Martin Alkaemper
23
24
25 The dune-fem module is free software; you can redistribute it and/or
26 modify it under the terms of the GNU General Public License as
27 published by the Free Software Foundation; either version 2 of
28 the License, or (at your option) any later version.
29
30 The dune-fem module is distributed in the hope that it will be useful,
31 but WITHOUT ANY WARRANTY; without even the implied warranty of
32 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
33 GNU General Public License for more details.
34
35 You should have received a copy of the GNU General Public License along
36 with this program; if not, write to the Free Software Foundation, Inc.,
37 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
38
39**************************************************************************/
40#ifndef DUNE_FEM_SCHEMES_FEMSCHEME_HH
41#define DUNE_FEM_SCHEMES_FEMSCHEME_HH
42
43#include <type_traits>
44#include <utility>
45#include <iostream>
46#include <memory>
47#include <dune/common/typeutilities.hh>
48
49// include discrete function space
51
53
57
58// FemScheme
59//----------
60
61template < class Op, class DF, typename = void >
63{
64 static const bool value = false;
66};
67template < class Op, class DF>
68struct AddDirichletBC<Op,DF,std::enable_if_t<std::is_void< decltype( std::declval<const Op>().
69 setConstraints( std::declval<DF&>() ) )>::value > >
70{
71 static const bool value = true;
72 using DirichletBlockVector = typename Op::DirichletBlockVector;
73};
74
75
76template< class Operator, class LinearInverseOperator >
78{
79public:
81 typedef typename Operator::ModelType ModelType;
82 typedef typename Operator::DomainFunctionType DomainFunctionType;
83 typedef typename Operator::RangeFunctionType RangeFunctionType;
84 typedef typename Operator::RangeFunctionType DiscreteFunctionType;
86 typedef typename DiscreteFunctionType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
87 typedef LinearInverseOperator LinearInverseOperatorType;
88
90 typedef typename ModelType::GridPartType GridPartType;
91 static_assert( std::is_same< typename DiscreteFunctionSpaceType::GridPartType, GridPartType >::value,
92 "GridPart of Space has to be identical to GridPart of Model class" );
93
95 typedef typename GridPartType::GridType GridType;
96
98 typedef typename DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType;
99
100 typedef typename Operator::JacobianOperatorType JacobianOperatorType;
101 typedef typename Operator::JacobianOperatorType LinearOperatorType;
104
105 typedef typename FunctionSpaceType::RangeType RangeType;
106 static const int dimRange = FunctionSpaceType::dimRange;
109 /*********************************************************/
110
112 : space_( space ),
113 // the elliptic operator (implicit)
114 implicitOperator_( space, space, model, parameter ),
115 // create linear operator (domainSpace,rangeSpace)
116 invOp_( parameter )
117 {}
118
121
122 template <typename O = DifferentiableOperatorType>
123 auto setQuadratureOrders(unsigned int interior, unsigned int surface)
124 -> Dune::void_t< decltype( std::declval< O >().setQuadratureOrders(0,0) ) >
125 {
126 fullOperator().setQuadratureOrders(interior,surface);
127 }
128
129 template <typename O = Operator>
130 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
132 {
133 implicitOperator_.setConstraints( u );
134 }
135 template <typename O = Operator>
136 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
138 {
139 implicitOperator_.setConstraints( u, v );
140 }
141 template <typename O = Operator>
142 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
144 {
145 implicitOperator_.subConstraints( u, v );
146 }
147 template <typename O = Operator>
148 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
150 {
151 implicitOperator_.setConstraints( value, u );
152 }
153 // template <typename O = Operator>
154 // std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,
155 // const decltype(implicitOperator_.dirichletBlocks())&>
156 const auto &dirichletBlocks() const
157 {
159 return implicitOperator_.dirichletBlocks();
160 }
161
163 {
164 implicitOperator_( arg, dest );
165 }
166 template <class GridFunction>
167 auto operator() ( const GridFunction &arg, DiscreteFunctionType &dest ) const
168 -> Dune::void_t<decltype(std::declval<const Operator&>()(arg,dest))>
169 {
170 implicitOperator_( arg, dest );
171 }
172
174 {
175 SolverInfo(bool pconverged,int plinearIterations,int pnonlinearIterations)
176 : converged(pconverged), linearIterations(plinearIterations), nonlinearIterations(pnonlinearIterations)
177 {}
181 };
182 void setErrorMeasure(ErrorMeasureType &errorMeasure) const
183 {
184 invOp_.setErrorMeasure(errorMeasure);
185 }
187 {
189 DiscreteFunctionType rhs0 = rhs;
190 setZeroConstraints( rhs0 );
191 setModelConstraints( solution );
192 invOp_( rhs0, solution );
193 invOp_.unbind();
195 }
197 {
198 DiscreteFunctionType bnd( solution );
199 bnd.clear();
200 setModelConstraints( solution );
202 invOp_( bnd, solution );
203 invOp_.unbind();
205 }
206
207 template< class GridFunction, std::enable_if_t<
208 std::is_same< decltype(
209 std::declval< const DifferentiableOperatorType >().jacobian(
210 std::declval< const GridFunction& >(), std::declval< JacobianOperatorType& >()
211 )
212 ), void >::value, int> i = 0
213 >
214 void jacobian( const GridFunction &ubar, JacobianOperatorType &linOp ) const
215 {
216 implicitOperator_.jacobian(ubar, linOp);
217 }
218
219 const GridPartType &gridPart () const { return space().gridPart(); }
220 const DiscreteFunctionSpaceType &space( ) const { return space_; }
221
222 const ModelType &model() const
223 {
224 return implicitOperator_.model();
225 }
227 {
228 return implicitOperator_.model();
229 }
230protected:
231 template <typename O = Operator>
232 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
234 template<class...Args>
235 void setZeroConstraints(Args&&...) const { }
236 template <typename O = Operator>
237 std::enable_if_t<AddDirichletBC<O,DomainFunctionType>::value,void>
238 setModelConstraints( DiscreteFunctionType &u ) const { fullOperator().setConstraints( u ); }
239 template<class...Args>
240 void setModelConstraints(Args&&... ) const { }
241 const DiscreteFunctionSpaceType &space_; // discrete function space
244};
245
246#endif // #ifndef DUNE_FEM_SCHEMES_FEMSCHEME_HH
STL namespace.
static ParameterContainer & container()
Definition: io/parameter.hh:193
Definition: femscheme.hh:63
void DirichletBlockVector
Definition: femscheme.hh:65
static const bool value
Definition: femscheme.hh:64
Definition: femscheme.hh:78
Operator::ModelType ModelType
type of the mathematical model
Definition: femscheme.hh:81
Operator::JacobianOperatorType LinearOperatorType
Definition: femscheme.hh:101
const DifferentiableOperatorType & fullOperator() const
Definition: femscheme.hh:119
SolverInfo solve(const DiscreteFunctionType &rhs, DiscreteFunctionType &solution) const
Definition: femscheme.hh:186
DifferentiableOperatorType & fullOperator()
Definition: femscheme.hh:120
DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType
type of function space (scalar functions,
Definition: femscheme.hh:98
Operator::JacobianOperatorType JacobianOperatorType
Definition: femscheme.hh:100
Dune::Fem::NewtonInverseOperator< LinearOperatorType, LinearInverseOperatorType > InverseOperatorType
Definition: femscheme.hh:102
GridPartType::GridType GridType
type of underyling hierarchical grid needed for data output
Definition: femscheme.hh:95
FemScheme(const DiscreteFunctionSpaceType &space, ModelType &model, const Dune::Fem::ParameterReader &parameter=Dune::Fem::Parameter::container())
Definition: femscheme.hh:111
void setModelConstraints(Args &&...) const
Definition: femscheme.hh:240
const GridPartType & gridPart() const
Definition: femscheme.hh:219
Operator::RangeFunctionType RangeFunctionType
Definition: femscheme.hh:83
const auto & dirichletBlocks() const
Definition: femscheme.hh:156
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > setConstraints(const DiscreteFunctionType &u, DiscreteFunctionType &v) const
Definition: femscheme.hh:137
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > setModelConstraints(DiscreteFunctionType &u) const
Definition: femscheme.hh:238
static const int dimRange
Definition: femscheme.hh:106
ModelType & model()
Definition: femscheme.hh:226
Operator::RangeFunctionType DiscreteFunctionType
Definition: femscheme.hh:84
const DiscreteFunctionSpaceType & space() const
Definition: femscheme.hh:220
SolverInfo solve(DiscreteFunctionType &solution) const
Definition: femscheme.hh:196
typename AddDirichletBC< Operator, DomainFunctionType >::DirichletBlockVector DirichletBlockVector
Definition: femscheme.hh:108
Operator::DomainFunctionType DomainFunctionType
Definition: femscheme.hh:82
const ModelType & model() const
Definition: femscheme.hh:222
void setErrorMeasure(ErrorMeasureType &errorMeasure) const
Definition: femscheme.hh:182
auto setQuadratureOrders(unsigned int interior, unsigned int surface) -> Dune::void_t< decltype(std::declval< O >().setQuadratureOrders(0, 0)) >
Definition: femscheme.hh:123
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > setConstraints(DomainFunctionType &u) const
Definition: femscheme.hh:131
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > subConstraints(const DiscreteFunctionType &u, DiscreteFunctionType &v) const
Definition: femscheme.hh:143
void jacobian(const GridFunction &ubar, JacobianOperatorType &linOp) const
Definition: femscheme.hh:214
static constexpr bool addDirichletBC
Definition: femscheme.hh:107
DifferentiableOperatorType implicitOperator_
Definition: femscheme.hh:242
LinearInverseOperator LinearInverseOperatorType
Definition: femscheme.hh:87
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > setConstraints(const RangeType &value, DiscreteFunctionType &u) const
Definition: femscheme.hh:149
void setZeroConstraints(Args &&...) const
Definition: femscheme.hh:235
InverseOperatorType::ErrorMeasureType ErrorMeasureType
Definition: femscheme.hh:103
ModelType::GridPartType GridPartType
grid view (e.g. leaf grid view) provided in the template argument list
Definition: femscheme.hh:90
std::enable_if_t< AddDirichletBC< O, DomainFunctionType >::value, void > setZeroConstraints(DiscreteFunctionType &u) const
Definition: femscheme.hh:233
DiscreteFunctionType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
Definition: femscheme.hh:86
const DiscreteFunctionSpaceType & space_
Definition: femscheme.hh:241
Operator DifferentiableOperatorType
Definition: femscheme.hh:85
InverseOperatorType invOp_
Definition: femscheme.hh:243
void operator()(const DiscreteFunctionType &arg, DiscreteFunctionType &dest) const
Definition: femscheme.hh:162
FunctionSpaceType::RangeType RangeType
Definition: femscheme.hh:105
Definition: femscheme.hh:174
int nonlinearIterations
Definition: femscheme.hh:180
int linearIterations
Definition: femscheme.hh:179
SolverInfo(bool pconverged, int plinearIterations, int pnonlinearIterations)
Definition: femscheme.hh:175
bool converged
Definition: femscheme.hh:178
std::function< bool(const RangeFunctionType &w, const RangeFunctionType &dw, double residualNorm) > ErrorMeasureType
Definition: newtoninverseoperator.hh:230
void unbind()
Definition: newtoninverseoperator.hh:305
int linearIterations() const
Definition: newtoninverseoperator.hh:311
bool converged() const
Definition: newtoninverseoperator.hh:333
int iterations() const
Definition: newtoninverseoperator.hh:309
void bind(const OperatorType &op)
Definition: newtoninverseoperator.hh:303
void setErrorMeasure(ErrorMeasureType finished)
Definition: newtoninverseoperator.hh:301