KASKADE 7 development version
errorEstimator.hh
Go to the documentation of this file.
1/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2/* */
3/* This file is part of the library KASKADE 7 */
4/* https://www.zib.de/research/projects/kaskade7-finite-element-toolbox */
5/* */
6/* Copyright (C) 2002-2012 Zuse Institute Berlin */
7/* */
8/* KASKADE 7 is distributed under the terms of the ZIB Academic License. */
9/* see $KASKADE/academic.txt */
10/* */
11/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
12
13#ifndef ALGORITHM_3_ERROR_ESTIMATOR_EE
14#define ALGORITHM_3_ERROR_ESTIMATOR_EE
15
16#include <algorithm>
17#include <iostream>
18#include <utility>
19#include <vector>
20
21#include <boost/timer/timer.hpp>
22#include <boost/fusion/include/at_c.hpp>
23
24#include "algorithm/newton_bridge.hh"
25#include "algorithm/adaptationStrategy.hh"
26//#include "algorithm/errorDistribution.hh"
27//#include "algorithm/mynorm.hh" // file does not exist
28//#include "algorithm/errorDistribution3.hh"
29#include "linalg/jacobiPreconditioner.hh"
30// #include "linalg/minresSolver.hh" // file does not exist
31#include "fem/forEach.hh"
32#include "fem/hierarchicErrorEstimator.hh"
33#include "fem/variables.hh"
34#include "linalg/direct.hh"
35//#include "linalg/blockDiagonalSchurPreconditioner.hh"
36#include "linalg/tcg.hh"
37#include "utilities/enums.hh"
38
39// forward declarations
40namespace Dune
41{
42 struct InverseOperatorResult;
43 template <class,int> class FieldVector;
44}
45
46namespace Kaskade
47{
48 namespace ErrorEstimator_Detail
49 {
50
51 template <class Functional>
52 std::integral_constant<int,Functional::stateId> getStateIdImpl(decltype(Functional::stateId)*);
53
54 template <class Functional>
55 std::integral_constant<int,Functional::yIdx> getStateIdImpl(decltype(Functional::yIdx)*);
56
57 template <class Functional>
58 std::integral_constant<int,Functional::controlId> getControlIdImpl(decltype(Functional::controlId)*);
59
60 template <class Functional>
61 std::integral_constant<int,Functional::uIdx> getControlIdImpl(decltype(Functional::uIdx)*);
62
63 template <class Functional>
64 std::integral_constant<int,Functional::adjointId> getAdjointIdImpl(decltype(Functional::adjointId)*);
65
66 template <class Functional>
67 std::integral_constant<int,Functional::pIdx> getAdjointIdImpl(decltype(Functional::pIdx)*);
68
69 template <class Functional>
70 std::integral_constant<int,Functional::lIdx> getAdjointIdImpl(decltype(Functional::lIdx)*);
71
72 template <class Functional>
73 constexpr int getStateId()
74 {
75 typedef decltype(getStateIdImpl<Functional>(nullptr)) type;
76 return type::value;
77 }
78
79 template <class Functional>
80 constexpr int getControlId()
81 {
82 typedef decltype(getControlIdImpl<Functional>(nullptr)) type;
83 return type::value;
84 }
85
86 template <class Functional>
87 constexpr int getAdjointId()
88 {
89 typedef decltype(getAdjointIdImpl<Functional>(nullptr)) type;
90 return type::value;
91 }
92
93 template < template <class,class,class> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription>
94 struct Traits
95 {
96 typedef VariableSet<VariableSetDescription> VarSet;
97 typedef VariableSet<ExtensionVariableSetDescription> ExtensionVarSet;
98
99 typedef Functional<VariableSetDescription,VariableSetDescription,VariableSetDescription> Functional_HH;
100 typedef Functional<VariableSetDescription,ExtensionVariableSetDescription,VariableSetDescription> Functional_HE;
101 typedef Functional<ExtensionVariableSetDescription,VariableSetDescription,ExtensionVariableSetDescription> Functional_EH;
102 typedef Functional<ExtensionVariableSetDescription,ExtensionVariableSetDescription,ExtensionVariableSetDescription> Functional_EE;
103
104 typedef typename Functional_HH::Scalar Scalar;
105
106 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HH> > Assembler_HH;
107 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HE> > Assembler_HE;
108 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EH> > Assembler_EH;
109 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EE> > Assembler_EE;
110
111 static constexpr int stateId = getStateId<Functional_HH>();
112 static constexpr int controlId = getControlId<Functional_HH>();
113 static constexpr int adjointId = getAdjointId<Functional_HH>();
114
115 static constexpr int dim = VariableSetDescription::Grid::dimension;
116
117 /*
118 * operators:
119 * the state operator is called A
120 * the control enters the state equation via operator B
121 * the derivatives of the lagrange functional are called Lyy, Lyu, Luy, Luu
122 *
123 * the endings indicate which ansatz and test spaces are used (L=ansatz space H=extension space), i.e. Luy_HE maps from the
124 * ansatz space to the extension space
125 */
126 typedef AssembledGalerkinOperator<Assembler_HH> Operator_HH;
127 typedef AssembledGalerkinOperator<Assembler_EH> Operator_EH;
128 typedef AssembledGalerkinOperator<Assembler_HE> Operator_HE;
129 typedef AssembledGalerkinOperator<Assembler_EE> Operator_EE;
130
131 typedef AssembledGalerkinOperator<Assembler_HH,stateId,stateId+1,stateId,stateId+1> Lyy_HH;
132 typedef AssembledGalerkinOperator<Assembler_EH,stateId,stateId+1,stateId,stateId+1> Lyy_EH;
133 typedef AssembledGalerkinOperator<Assembler_HE,stateId,stateId+1,stateId,stateId+1> Lyy_HE;
134 typedef AssembledGalerkinOperator<Assembler_EE,stateId,stateId+1,stateId,stateId+1> Lyy_EE;
135
136 typedef AssembledGalerkinOperator<Assembler_HH,stateId,stateId+1,controlId,controlId+1> Lyu_HH;
137 typedef AssembledGalerkinOperator<Assembler_EH,stateId,stateId+1,controlId,controlId+1> Lyu_EH;
138 typedef AssembledGalerkinOperator<Assembler_HE,stateId,stateId+1,controlId,controlId+1> Lyu_HE;
139 typedef AssembledGalerkinOperator<Assembler_EE,stateId,stateId+1,controlId,controlId+1> Lyu_EE;
140
141 typedef AssembledGalerkinOperator<Assembler_HH,controlId,controlId+1,stateId,stateId+1> Luy_HH;
142 typedef AssembledGalerkinOperator<Assembler_EH,controlId,controlId+1,stateId,stateId+1> Luy_EH;
143 typedef AssembledGalerkinOperator<Assembler_HE,controlId,controlId+1,stateId,stateId+1> Luy_HE;
144 typedef AssembledGalerkinOperator<Assembler_EE,controlId,controlId+1,stateId,stateId+1> Luy_EE;
145
146 typedef AssembledGalerkinOperator<Assembler_HH,controlId,controlId+1,controlId,controlId+1> Luu_HH;
147 typedef AssembledGalerkinOperator<Assembler_EH,controlId,controlId+1,controlId,controlId+1> Luu_EH;
148 typedef AssembledGalerkinOperator<Assembler_HE,controlId,controlId+1,controlId,controlId+1> Luu_HE;
149 typedef AssembledGalerkinOperator<Assembler_EE,controlId,controlId+1,controlId,controlId+1> Luu_EE;
150
151 typedef AssembledGalerkinOperator<Assembler_HH,adjointId,adjointId+1,controlId,controlId+1> B_HH;
152 typedef AssembledGalerkinOperator<Assembler_EH,adjointId,adjointId+1,controlId,controlId+1> B_EH;
153 typedef AssembledGalerkinOperator<Assembler_HE,adjointId,adjointId+1,controlId,controlId+1> B_HE;
154 typedef AssembledGalerkinOperator<Assembler_EE,adjointId,adjointId+1,controlId,controlId+1> B_EE;
155
156 typedef AssembledGalerkinOperator<Assembler_HH,adjointId,adjointId+1,stateId,stateId+1> A_HH;
157 typedef AssembledGalerkinOperator<Assembler_EH,adjointId,adjointId+1,stateId,stateId+1> A_EH;
158 typedef AssembledGalerkinOperator<Assembler_HE,adjointId,adjointId+1,stateId,stateId+1> A_HE;
159 typedef AssembledGalerkinOperator<Assembler_EE,adjointId,adjointId+1,stateId,stateId+1> A_EE;
160
161 typedef AssembledGalerkinOperator<Assembler_HH,controlId,controlId+1,adjointId,adjointId+1> BT_HH;
162 typedef AssembledGalerkinOperator<Assembler_EH,controlId,controlId+1,adjointId,adjointId+1> BT_EH;
163 typedef AssembledGalerkinOperator<Assembler_HE,controlId,controlId+1,adjointId,adjointId+1> BT_HE;
164 typedef AssembledGalerkinOperator<Assembler_EE,controlId,controlId+1,adjointId,adjointId+1> BT_EE;
165
166 typedef AssembledGalerkinOperator<Assembler_HH,stateId,stateId+1,adjointId,adjointId+1> AT_HH;
167 typedef AssembledGalerkinOperator<Assembler_EH,stateId,stateId+1,adjointId,adjointId+1> AT_EH;
168 typedef AssembledGalerkinOperator<Assembler_HE,stateId,stateId+1,adjointId,adjointId+1> AT_HE;
169 typedef AssembledGalerkinOperator<Assembler_EE,stateId,stateId+1,adjointId,adjointId+1> AT_EE;
170
174 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<>::type Vector;
175 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<stateId,stateId+1>::type VectorY;
176 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<controlId,controlId+1>::type VectorU;
177 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<adjointId,adjointId+1>::type VectorP;
178 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<>::type ExtensionVector;
179 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<stateId,stateId+1>::type ExtensionVectorY;
180 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<controlId,controlId+1>::type ExtensionVectorU;
181 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<adjointId,adjointId+1>::type ExtensionVectorP;
182
183 /*
184 * Coefficient vector initializers
185 */
186 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<> Vector_Initializer;
187 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<stateId,stateId+1> VectorY_Initializer;
188 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<controlId,controlId+1> VectorU_Initializer;
189 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<adjointId,adjointId+1> VectorP_Initializer;
190 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<> ExtensionVector_Initializer;
191 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<stateId,stateId+1> ExtensionVectorY_Initializer;
192 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<controlId,controlId+1> ExtensionVectorU_Initializer;
193 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<adjointId,adjointId+1> ExtensionVectorP_Initializer;
194 };
195 }
196
197 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
198 class NormFunctional, class LinearSolverLA, class LinearSolverHA, class LinearSolverLM=LinearSolverLA, class LinearSolverHM=LinearSolverHA,
199 bool lumpM=false, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
200 class HierarchicalBasisErrorEstimator2 : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
201 {
202 static constexpr int noOfVariables = ExtensionVariableSetDescription::noOfVariables;
203 template <class AnsatzVars, class TestVars, class OriginVars> using ErrorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
204 typedef ErrorEstimator_Detail::Traits<ErrorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
205 // typedef typename LF::AnsatzVars LVars;
206 static constexpr int yId = 1, uId = 0, pId = 2;
207 public:
208 typedef Functional<VariableSetDescription,VariableSetDescription,VariableSetDescription,lumpM> Functional_HH;
209 typedef Functional<ExtensionVariableSetDescription,VariableSetDescription,ExtensionVariableSetDescription,lumpM> Functional_EH;
210 typedef Functional<VariableSetDescription,ExtensionVariableSetDescription,VariableSetDescription,lumpM> Functional_HE;
211 typedef Functional<ExtensionVariableSetDescription,ExtensionVariableSetDescription,ExtensionVariableSetDescription,lumpM> Functional_EE;
212 typedef HierarchicErrorEstimator<LinearizationAt<Functional_HH>,ExtensionVariableSetDescription,ExtensionVariableSetDescription> ErrorEstimator;
213 // typedef HierarchicErrorEstimator<LinearizationAt<Functional_HH>,ExtensionVariableSetDescription,ExtensionVariableSetDescription,HierarchicErrorEstimatorDetail::TakeAllD2<LinearizationAt<Functional_HH> > > ErrorEstimator;
214 typedef VariationalFunctionalAssembler<ErrorEstimator> Assembler;
215 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HH> > Ass_HH;
216 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HE> > Ass_HE;
217 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EH> > Ass_EH;
218 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EE> > Ass_EE;
219
220 typedef typename Traits::AT_HH AT_HH;
221 typedef typename Traits::AT_EE AT_EE;
222 typedef typename Traits::AT_EH AT_EH;
223
224 typedef AssembledGalerkinOperator<Ass_EE,pId,pId+1,yId,yId+1> Ahh;
225 typedef AssembledGalerkinOperator<Ass_HH,pId,pId+1,yId,yId+1> All;
226 typedef AssembledGalerkinOperator<Ass_HE,pId,pId+1,yId,yId+1> Alh;
227 typedef AssembledGalerkinOperator<Ass_EH,pId,pId+1,yId,yId+1> Ahl;
228 typedef AssembledGalerkinOperator<Ass_EE,yId,yId+1,yId,yId+1> Myyhh;
229
230 typedef AssembledGalerkinOperator<Ass_HH,pId,pId+1,uId,uId+1> Bll;
231 typedef AssembledGalerkinOperator<Ass_HE,pId,pId+1,uId,uId+1> Blh;
232 typedef AssembledGalerkinOperator<Ass_EH,pId,pId+1,uId,uId+1> Bhl;
233 typedef AssembledGalerkinOperator<Ass_EE,pId,pId+1,uId,uId+1> Bhh;
234
235 typedef AssembledGalerkinOperator<Ass_HH,uId,uId+1,pId,pId+1> BTll;
236 typedef AssembledGalerkinOperator<Ass_HE,uId,uId+1,pId,pId+1> BTlh;
237 typedef AssembledGalerkinOperator<Ass_EH,uId,uId+1,pId,pId+1> BThl;
238 typedef AssembledGalerkinOperator<Ass_EE,uId,uId+1,pId,pId+1> BThh;
239
240 typedef AssembledGalerkinOperator<Ass_HH,uId,uId+1,uId,uId+1> Muull;
241 typedef AssembledGalerkinOperator<Ass_HE,uId,uId+1,uId,uId+1> Muulh;
242 typedef AssembledGalerkinOperator<Ass_EH,uId,uId+1,uId,uId+1> Muuhl;
243 typedef AssembledGalerkinOperator<Ass_EE,uId,uId+1,uId,uId+1> Muuhh;
244
245 typedef typename Functional_HH::Scalar Scalar;
246 static constexpr int dim = VariableSetDescription::Grid::dimension;
247 typedef VariableSet<VariableSetDescription> VarSet;
248 typedef VariableSet<ExtensionVariableSetDescription> ExtVarSet;
249 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<> CVL;
250 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<yId,yId+1> CVLY;
251 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<uId,uId+1> CVLU;
252 typedef typename VariableSetDescription::template CoefficientVectorRepresentation<pId,pId+1> CVLP;
253 typedef typename CVL::type CoefficientVectorL;
254 typedef typename CVLY::type CoefficientVectorLY;
255 typedef typename CVLU::type CoefficientVectorLU;
256 typedef typename CVLP::type CoefficientVectorLP;
257 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<0,1>::type CoefficientVectorH01;
258 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<0,2>::type CoefficientVectorH02;
259 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<2,3> CVH23;
260 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<2,3>::type CoefficientVectorH23;
261 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<> CVH;
262 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<yId,yId+1> CVHY;
263 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<uId,uId+1> CVHU;
264 typedef typename ExtensionVariableSetDescription::template CoefficientVectorRepresentation<pId,pId+1> CVHP;
265 typedef typename CVH::type CoefficientVectorH;
266 typedef typename CVHY::type CoefficientVectorHY;
267 typedef typename CVHU::type CoefficientVectorHU;
268 typedef typename CVHP::type CoefficientVectorHP;
269 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
270 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
271
272 typedef MatrixAsTriplet<Scalar> Matrix;
273
274
275 HierarchicalBasisErrorEstimator2(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
276 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
277 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
278 extensionSpace(extensionSpace_), assembler(extensionVariableSetDescription.spaces),
279 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
280 {
281 initAssemblers();
282 }
283
284 virtual ~HierarchicalBasisErrorEstimator2(){}
285
286 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
287 {
288 double ilutTol = 0.01;
289 int ilutLfil = 100;
290 double minresTol = 1e-3;
291
292 boost::timer::cpu_timer overallTimer;
293 using boost::fusion::at_c;
294 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
295 Bridge::Vector<VarSet> const& xl = dynamic_cast<const Bridge::Vector<VarSet>&>(x_);
296 Bridge::Vector<ExtVarSet> xe(extensionVariableSetDescription);
297 assembleAt(xl,xe);
298 Bridge::Vector<VarSet> const& dxl = dynamic_cast<const Bridge::Vector<VarSet>&>(dx_);
299
300 /***************************************************************************************/
301 // estimate error in dual variable
302 typename ExtensionVariableSetDescription::VariableSet tmpRep(extensionVariableSetDescription), errorEstimate_E(extensionVariableSetDescription);
303 typename VariableSetDescription::VariableSet tmpRepL(variableSetDescription), errorEstimate_H(variableSetDescription);
304
305 All all(*ass_HH,false);
306 Ahh ahh(*ass_EE,false);
307 Alh alh(*ass_HE,false);
308 Ahl ahl(*ass_EH,false);
309
310 Matrix mAll = all.template get<Matrix>();
311 mAll.transpose();
312 Matrix mAhh = ahh.template get<Matrix>();
313 mAhh.transpose();
314 Matrix mAlh = alh.template get<Matrix>();
315 mAlh.transpose();
316 Matrix mAhl = ahl.template get<Matrix>();
317 mAhl.transpose();
318
319 boost::timer::cpu_timer atimer;
320 assembler.assemble(ErrorEstimator(LinearizationAt<Functional_HH>(*F_HH,xl.get()),dxl.get()));
321 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << "s" << std::endl;
322 typedef AssembledGalerkinOperator<Assembler,yId,yId+1,pId,pId+1> StateOperator;
323 StateOperator stateOp(assembler,false);
324 CoefficientVectorHY rhs0(assembler.template rhs<yId,yId+1>());
325 CoefficientVectorHP sol0(CVHP::init(extensionVariableSetDescription));
326
327
328 tmpRep *= 0;
329 at_c<yId>(tmpRep.data) = at_c<0>(rhs0.data);
330 CoefficientVectorLP xlp(CVLP::init(variableSetDescription));
331 at_c<0>(xlp.data) = at_c<pId>(dxl.get().data).coefficients();
332 std::vector<Scalar> tmpRepVec, xlpVec;
333 IstlInterfaceDetail::toVector(rhs0,tmpRepVec);
334 if(verbose){
335 std::cout << "initial rhs entries:" << std::endl;
336 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
337 for(size_t i=0; i<tmpRepVec.size(); ++i){
338 if(std::isnan(tmpRepVec[i])) std::cout << "nan at " << i << std::endl;
339
340 Scalar val = tmpRepVec[i];
341 if(val < minRhs) minRhs = val;
342 if(val > maxRhs) maxRhs = val;
343 rhsSum += val*val;
344 }
345 std::cout << "max rhs: " << maxRhs << std::endl;
346 std::cout << "min rhs: " << minRhs << std::endl;
347 std::cout << "sum: " << rhsSum << std::endl;
348 }
349 IstlInterfaceDetail::toVector(xlp,xlpVec);
350 mAhl.usmv(-1.0,xlpVec,tmpRepVec);
351
352 if(verbose)
353 {
354 std::cout << "corrected rhs entries:" << std::endl;
355 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
356 for(size_t i=0; i<tmpRepVec.size(); ++i){
357 if(std::isnan(tmpRepVec[i])) std::cout << "nan at " << i << std::endl;
358
359 Scalar val = tmpRepVec[i];
360 if(val < minRhs) minRhs = val;
361 if(val > maxRhs) maxRhs = val;
362 rhsSum += val*val;
363 }
364 std::cout << "max rhs: " << maxRhs << std::endl;
365 std::cout << "min rhs: " << minRhs << std::endl;
366 std::cout << "sum: " << rhsSum << std::endl;
367 }
368 IstlInterfaceDetail::fromVector(tmpRepVec,rhs0);
369 Dune::InverseOperatorResult res;
370 boost::timer::cpu_timer timer;
371 if(verbose) std::cout << "first solve: ";
372 //LinearSolverHA(mAhh,DirectType::MUMPS,MatrixProperties::GENERAL).apply(sol0,rhs0);
373 if(fast)
374 {
375 MatrixRepresentedOperator<Matrix,CoefficientVectorH23,CoefficientVectorH23> mAhhOp(mAhh);
376 ILUTPreconditioner<MatrixRepresentedOperator<Matrix,CoefficientVectorH23,CoefficientVectorH23> > prec(mAhhOp,ilutLfil,ilutTol,1);
377
378 MINRESSolverK<CoefficientVectorH23> minres(mAhhOp, prec, minresTol, 100, 1);
379 minres.apply(sol0,rhs0,res);
380 }
381 else InverseLinearOperator<DirectSolver<CoefficientVectorHP,CoefficientVectorHP> >(DirectSolver<CoefficientVectorHP,CoefficientVectorHP>(mAhh,DirectType::UMFPACK3264,MatrixProperties::GENERAL)).apply(rhs0,sol0);
382 if(verbose) std::cout << boost::timer::format(timer.elapsed()) << std::endl;
383 // InverseLinearOperator<DirectSolver<CoefficientVectorHP,CoefficientVectorHP> >(DirectSolver<CoefficientVectorHP,CoefficientVectorHP>(mAhh,DirectType::MUMPS,MatrixProperties::GENERAL)).apply(tmpRepVec);
384 // IstlInterfaceDetail::fromVector(tmpRepVec,sol0);
385
386 // //MyH1SemiNorm myNorm;
387 // //std::cout << "||rhs||=" << myNorm(tmpRep) << std::endl;
388 // LinearSolverHA(stateOp,DirectType::MUMPS,MatrixProperties::GENERAL).apply(sol0,rhs0);
389 // //directInverseOperator(stateOp,DirectType::MUMPS,MatrixProperties::GENERAL).apply(rhs0,sol0);
390 // tmpRep *= 0;
391 at_c<pId>(tmpRep.data) = at_c<0>(sol0.data);
392 at_c<pId>(errorEstimate_E.data) = at_c<0>(sol0.data);
393 //std::cout << "||sol||=" << myNorm(tmpRep) << std::endl;
394 Scalar maxVal = std::numeric_limits<Scalar>::min(),
395 minVal = std::numeric_limits<Scalar>::max(),
396 sum = 0;
397 if(verbose){
398 for(int i=0; i<at_c<pId>(errorEstimate_E.data).size(); ++i)
399 {
400 for(int j=0; j<at_c<pId>(errorEstimate_E.data).coefficients()[i].size(); ++j)
401 {
402 Scalar val = at_c<pId>(errorEstimate_E.data).coefficients()[i][j];
403 sum += val*val;
404 if(val > maxVal) maxVal = val;
405 if(val < minVal) minVal = val;
406 }
407 }
408 std::cout << "p: minimal value: " << minVal << std::endl;
409 std::cout << "p: maximal value: " << maxVal << std::endl;
410 std::cout << "p: l2-sum: " << sum << std::endl;
411 }
412 /***************************************************************************************/
413 // estimate error in control variable
414 // rhs vectors
415 CoefficientVectorLU rhs1L(ass_EH->template rhs<uId,uId+1>());
416 CoefficientVectorHU rhs1H(ass_EE->template rhs<uId,uId+1>());
417 Bhl bhl(*ass_EH,false);
418 Blh blh(*ass_HE,false);
419 Bll bll(*ass_HH,false);
420 Bhh bhh(*ass_EE,false);
421
422 MatrixAsTriplet<Scalar> tmpmat = blh.template get<MatrixAsTriplet<Scalar> >();
423 std::vector<Scalar> tmpx, tmpy;
424 IstlInterfaceDetail::toVector(sol0,tmpy);
425 IstlInterfaceDetail::toVector(rhs1L,tmpx);
426 tmpmat.usmtv(-1.0,tmpy,tmpx);
427 IstlInterfaceDetail::fromVector(tmpx,rhs1L);
428 //std::cout << "rhs1L=" << at_c<0>(rhs1L.data)[0] << std::endl;
429
430 tmpx.clear(); tmpy.clear();
431 IstlInterfaceDetail::toVector(sol0,tmpy);
432 IstlInterfaceDetail::toVector(rhs1H,tmpx);
433 tmpmat = bhh.template get<MatrixAsTriplet<Scalar> >();
434 tmpmat.usmtv(-1.0,tmpy,tmpx);
435 IstlInterfaceDetail::fromVector(tmpx,rhs1H);
436 // bophh.applyscaleadd(-1.0,sol0,rhs1H);
437 //std::cout << "rhs1H=" << at_c<0>(rhs1H.data)[0] << std::endl;
438
439 Muuhh muuhh(*ass_EE,true);
440 Muull muull(*ass_HH,true);
441 Muulh muulh(*ass_HE,true);
442 CoefficientVectorHU solhu(CVHU::init(extensionVariableSetDescription));
443 CoefficientVectorLU sollu(CVLU::init(variableSetDescription));
444 //LinearSolverLM(muull,DirectType::MUMPS,MatrixProperties::GENERAL).apply(sollu,rhs1L);
445 timer.start();
446 if(verbose) std::cout << "second solve: ";
447 JacobiPreconditioner<Muull>(muull).apply(sollu,rhs1L);
448 // directInverseOperator(muull,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(rhs1L,sollu);
449 if(verbose) std::cout << boost::timer::format(timer.elapsed()) << std::endl;
450 muulh.applyscaleadd(-1.0,sollu,rhs1H);
451 JacobiPreconditioner<Muuhh>(muuhh).apply(solhu,rhs1H);
452 // LinearSolverHM(muuhh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solhu,rhs1H);
453 // directInverseOperator(muuhh,DirectType::MUMPS,MatrixProperties::GENERAL).apply(rhs1H,solhu);
454
455 //std::cout << "sollu=" << at_c<0>(sollu.data)[0] << std::endl;
456 //std::cout << "solhu=" << at_c<0>(solhu.data)[0] << std::endl;
457 at_c<uId>(errorEstimate_E.data) = at_c<0>(solhu.data);
458 at_c<uId>(errorEstimate_H.data) = at_c<0>(sollu.data);
459 /***************************************************************************************/
460 // estimate error in state variable
461 CoefficientVectorLY solly(CVLP::init(variableSetDescription));
462 CoefficientVectorHY solhy(CVHP::init(extensionVariableSetDescription));
463 CoefficientVectorLP rhslp(ass_HH->template rhs<pId,pId+1>());
464 CoefficientVectorHP rhshp(ass_EE->template rhs<pId,pId+1>());
465
466 std::vector<Scalar> tmprhslp, tmprhshp;
467 IstlInterfaceDetail::toVector(rhslp,tmprhslp);
468 IstlInterfaceDetail::toVector(rhshp,tmprhshp);
469 if(verbose){
470 std::cout << "initial rhsl entries:" << std::endl;
471 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
472 for(size_t i=0; i<tmprhslp.size(); ++i){
473 if(std::isnan(tmprhslp[i])) std::cout << "nan at " << i << std::endl;
474
475 Scalar val = tmprhslp[i];
476 if(val < minRhs) minRhs = val;
477 if(val > maxRhs) maxRhs = val;
478 rhsSum += val*val;
479 }
480 std::cout << "max rhs: " << maxRhs << std::endl;
481 std::cout << "min rhs: " << minRhs << std::endl;
482 std::cout << "sum: " << rhsSum << std::endl;
483 }
484 if(verbose){
485 std::cout << "initial rhsh entries:" << std::endl;
486 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
487 for(size_t i=0; i<tmprhshp.size(); ++i){
488 if(std::isnan(tmprhshp[i])) std::cout << "nan at " << i << std::endl;
489
490 Scalar val = tmprhshp[i];
491 if(val < minRhs) minRhs = val;
492 if(val > maxRhs) maxRhs = val;
493 rhsSum += val*val;
494 }
495 std::cout << "max rhs: " << maxRhs << std::endl;
496 std::cout << "min rhs: " << minRhs << std::endl;
497 std::cout << "sum: " << rhsSum << std::endl;
498 }
499
500 bll.applyscaleadd(-1.0,sollu,rhslp);
501 if(verbose){
502 std::cout << "1: rhsl entries:" << std::endl;
503 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
504 for(size_t i=0; i<tmprhslp.size(); ++i){
505 if(std::isnan(tmprhslp[i])) std::cout << "nan at " << i << std::endl;
506
507 Scalar val = tmprhslp[i];
508 if(val < minRhs) minRhs = val;
509 if(val > maxRhs) maxRhs = val;
510 rhsSum += val*val;
511 }
512 std::cout << "max rhs: " << maxRhs << std::endl;
513 std::cout << "min rhs: " << minRhs << std::endl;
514 std::cout << "sum: " << rhsSum << std::endl;
515 }
516 bhl.applyscaleadd(-1.0,solhu,rhslp);
517 if(verbose){
518 std::cout << "2: rhsl entries:" << std::endl;
519 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
520 for(size_t i=0; i<tmprhslp.size(); ++i){
521 if(std::isnan(tmprhslp[i])) std::cout << "nan at " << i << std::endl;
522
523 Scalar val = tmprhslp[i];
524 if(val < minRhs) minRhs = val;
525 if(val > maxRhs) maxRhs = val;
526 rhsSum += val*val;
527 }
528 std::cout << "max rhs: " << maxRhs << std::endl;
529 std::cout << "min rhs: " << minRhs << std::endl;
530 std::cout << "sum: " << rhsSum << std::endl;
531 }
532
533 blh.applyscaleadd(-1.0,sollu,rhshp);
534 if(verbose){
535 std::cout << "1: rhsh entries:" << std::endl;
536 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
537 for(size_t i=0; i<tmprhshp.size(); ++i){
538 if(std::isnan(tmprhshp[i])) std::cout << "nan at " << i << std::endl;
539
540 Scalar val = tmprhshp[i];
541 if(val < minRhs) minRhs = val;
542 if(val > maxRhs) maxRhs = val;
543 rhsSum += val*val;
544 }
545 std::cout << "max rhs: " << maxRhs << std::endl;
546 std::cout << "min rhs: " << minRhs << std::endl;
547 std::cout << "sum: " << rhsSum << std::endl;
548 }
549 bhh.applyscaleadd(-1.0,solhu,rhshp);
550 if(verbose) {
551 std::cout << "2: rhsh entries:" << std::endl;
552 Scalar maxRhs=std::numeric_limits<Scalar>::min(), minRhs=std::numeric_limits<Scalar>::max(), rhsSum=0;
553 for(size_t i=0; i<tmprhshp.size(); ++i){
554 if(std::isnan(tmprhshp[i])) std::cout << "nan at " << i << std::endl;
555
556 Scalar val = tmprhshp[i];
557 if(val < minRhs) minRhs = val;
558 if(val > maxRhs) maxRhs = val;
559 rhsSum += val*val;
560 }
561 std::cout << "max rhs: " << maxRhs << std::endl;
562 std::cout << "min rhs: " << minRhs << std::endl;
563 std::cout << "sum: " << rhsSum << std::endl;
564 }
565
566 // LinearSolverLA(ahh,DirectType::MUMPS,MatrixProperties::GENERAL).apply(solhy,rhshp);
567 // ahl.applyscaleadd(-1.0,solhy,rhslp);
568 // LinearSolverLA(all,DirectType::MUMPS,MatrixProperties::GENERAL).apply(solly,rhslp);
569 timer.start();
570 if(verbose) std::cout << "third solve: ";
571 LinearSolverLA(all,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solly,rhslp);
572 if(verbose) std::cout << boost::timer::format(timer.elapsed()) << std::endl;
573 if(verbose)std::cout << "first solve for last variable finished." << std::endl;
574 alh.applyscaleadd(-1.0,solly,rhshp);
575 timer.start();
576 if(verbose) std::cout << "fourth solve: ";
577 if(fast)
578 {
579 ILUTPreconditioner<Ahh> prec2(ahh,ilutLfil,ilutTol,1);
580
581 MINRESSolverK<CoefficientVectorH23> minres2(ahh, prec2, minresTol, 100, 1);
582 minres2.apply(solhy,rhshp,res);
583 }
584 else LinearSolverHA(ahh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solhy,rhshp);
585 if(verbose) std::cout << boost::timer::format(timer.elapsed()) << std::endl;
586 at_c<yId>(errorEstimate_E.data) = at_c<0>(solhy.data);
587 at_c<yId>(errorEstimate_H.data) = at_c<0>(solly.data);
588
589 if(verbose){
590 Scalar tmpMax = maxVal, tmpMin = minVal, tmpSum = sum;
591 maxVal = std::numeric_limits<Scalar>::min();
592 minVal = std::numeric_limits<Scalar>::max();
593 sum = 0;
594
595 for(int i=0; i<at_c<uId>(errorEstimate_E.data).size(); ++i)
596 {
597 for(int j=0; j<at_c<uId>(errorEstimate_E.data).coefficients()[i].size(); ++j)
598 {
599 Scalar val = at_c<uId>(errorEstimate_E.data).coefficients()[i][j];
600 sum += val*val;
601 if(val > maxVal) maxVal = val;
602 if(val < minVal) minVal = val;
603 }
604 }
605 std::cout << "u: minimal value: " << minVal << std::endl;
606 std::cout << "u: maximal value: " << maxVal << std::endl;
607 std::cout << "u: l2-sum: " << sum << std::endl;
608 tmpMax = std::max(tmpMax,maxVal), tmpMin = std::min(tmpMin,minVal), tmpSum += sum;
609 maxVal = std::numeric_limits<Scalar>::min();
610 minVal = std::numeric_limits<Scalar>::max();
611 sum = 0;
612
613 for(int i=0; i<at_c<yId>(errorEstimate_E.data).size(); ++i)
614 {
615 for(int j=0; j<at_c<yId>(errorEstimate_E.data).coefficients()[i].size(); ++j)
616 {
617 Scalar val = at_c<yId>(errorEstimate_E.data).coefficients()[i][j];
618 sum += val*val;
619 if(val > maxVal) maxVal = val;
620 if(val < minVal) minVal = val;
621 }
622 }
623 std::cout << "y: minimal value: " << minVal << std::endl;
624 std::cout << "y: maximal value: " << maxVal << std::endl;
625 std::cout << "y: l2-sum: " << sum << std::endl;
626
627 tmpSum += sum;
628 std::cout << "minimal value: " << std::min(tmpMin,minVal) << std::endl;
629 std::cout << "maximal value: " << std::max(tmpMax,maxVal) << std::endl;
630 std::cout << "l2-sum: " << tmpSum << std::endl;
631 }
632 /***************************************************************************************/
633 // Transfer error indicators to cells.
634 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
635 errorDistribution.clear();
636 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
637 Scalar errLevel(0.0), minRefine(0.2);
638
639 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
640 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
641 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
642 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
643 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
644 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
645
646 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
647 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
648 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
649 mde = distError;
650 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
651
652 if(verbose)
653 {
654 std::string name = "errorDistribution_";
655 name += std::to_string(step);
656 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
657 }
658 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
659 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
660 // Aufteilen von v_i auf die einzelnen Zellen
661 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
662 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
663 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
664
665 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
666
667 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
668 }
669
670 void refineGrid()
671 {
672 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
673// Scalar tol = squaredFraction*totalErrorSquared;
674// Scalar sum = 0;
675// std::sort(errorDistribution.begin(), errorDistribution.end(), [](std::pair<double,int> const& a, std::pair<double,int> const& b) -> bool
676// {
677// return a.first > b.first;
678// });
679// std::vector<size_t> cellIds;
680// for(size_t i=0; i<errorDistribution.size(); ++i)
681// {
682// sum += errorDistribution[i].first;
683// cellIds.push_back(errorDistribution[i].second);
684// if(sum > tol) break;
685// }
686//
687// forEachCell(variableSetDescription.gridView, [&](Cell const& cell)
688// {
689// if(std::find(cellIds.begin(),cellIds.end(),variableSetDescription.gridView.indexSet().index(cell))!=cellIds.end()) extensionSpace.gridManager().mark(1,cell);
690// });
691// extensionSpace.gridManager().adaptAtOnce();
692 }
693
694 double estimatedAbsoluteError() const final
695 {
696 return sqrt(fabs(totalErrorSquared));
697 }
698
699 size_t gridSize() const final
700 {
701 return extensionSpace.gridManager().grid().size(0);
702 }
703
704 template <typename... Args>
705 void initFunctionals(const Args&... args)
706 {
707 F_HH.reset(new Functional_HH(args...));
708 F_HE.reset(new Functional_HE(args...));
709 F_EH.reset(new Functional_EH(args...));
710 F_EE.reset(new Functional_EE(args...));
711 }
712
713 private:
714 void initAssemblers()
715 {
716 ass_HH.reset(new Ass_HH(extensionVariableSetDescription.spaces));
717 ass_HE.reset(new Ass_HE(extensionVariableSetDescription.spaces));
718 ass_EH.reset(new Ass_EH(extensionVariableSetDescription.spaces));
719 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
720 }
721
722 void assembleAt(const Bridge::Vector<VarSet>& xl, const Bridge::Vector<ExtVarSet>& xe)
723 {
724 ass_HH->assemble(linearization(*F_HH,xl.get()));
725 ass_HE->assemble(linearization(*F_HE,xl.get()));
726 ass_EH->assemble(linearization(*F_EH,xe.get()));
727 ass_EE->assemble(linearization(*F_EE,xe.get()));
728 }
729
730 std::unique_ptr<Functional_HH> F_HH;
731 std::unique_ptr<Functional_HE> F_HE;
732 std::unique_ptr<Functional_EH> F_EH;
733 std::unique_ptr<Functional_EE> F_EE;
734 std::unique_ptr<Ass_HH> ass_HH;
735 std::unique_ptr<Ass_HE> ass_HE;
736 std::unique_ptr<Ass_EH> ass_EH;
737 std::unique_ptr<Ass_EE> ass_EE;
738 // LF& lf;
739 // LVars& lvars;
740 NormFunctional& normFunctional;
741 VariableSetDescription& variableSetDescription;
742 ExtensionVariableSetDescription& extensionVariableSetDescription;
743 ExtensionSpace& extensionSpace;
744 Assembler assembler;
745 Scalar squaredFraction;
746 Scalar totalErrorSquared;
747 std::vector<std::pair<double,int> > errorDistribution;
748 bool verbose, fast;
749 // AdjustRHS<CoefficientVector> adjustRHS;
750 };
751
752
753 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
754 class NormFunctional, class LinearSolverHA, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
755 class StateEquationHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
756 {
757 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
758 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
759 public:
760 // Functionals
761 typedef typename Traits::Functional_HH Functional_HH;
762 typedef typename Traits::Functional_HE Functional_HE;
763 typedef typename Traits::Functional_EE Functional_EE;
764
765 // assemblers
766 typedef typename Traits::Assembler_HE Ass_HE;
767 typedef typename Traits::Assembler_EE Ass_EE;
768
769 // operators
770 typedef typename Traits::B_HE Blh;
771 typedef typename Traits::A_HE Alh;
772 typedef typename Traits::A_EE Ahh;
773
774 typedef typename Traits::Scalar Scalar;
775 static constexpr int dim = VariableSetDescription::Grid::dimension;
776
777 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
778 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
779
780
781
782 StateEquationHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
783 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
784 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
785 extensionSpace(extensionSpace_),
786 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
787 {
788 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
789 ass_HE.reset(new Ass_HE(extensionVariableSetDescription.spaces));
790 }
791
792 virtual ~StateEquationHBErrorEstimator(){}
793
794 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
795 {
796 boost::timer::cpu_timer overallTimer;
797 using boost::fusion::at_c;
798 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
799 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
800 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
801 boost::timer::cpu_timer atimer;
802 assembleAt(xl,xe);
803 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
804 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
805
806 /***************************************************************************************/
807 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
808 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
809 Blh blh(*ass_HE,false);
810 Alh alh(*ass_HE,false);
811 Ahh ahh(*ass_EE,false);
812 typename Traits::ExtensionVectorP rhs0(ass_EE->template rhs<Traits::adjointId,Traits::adjointId+1>());
813 rhs0 *= -1.0;
814 typename Traits::ExtensionVectorY sol0(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
815 typename Traits::VectorU sollu(Traits::VectorU_Initializer::init(variableSetDescription));
816 typename Traits::VectorY solly(Traits::VectorY_Initializer::init(variableSetDescription));
817 at_c<0>(solly.data) = at_c<Traits::stateId>(dxl.get().data).coefficients();
818 at_c<0>(sollu.data) = at_c<Traits::controlId>(dxl.get().data).coefficients();
819 blh.applyscaleadd(-1.0,sollu,rhs0);
820 alh.applyscaleadd(-1.0,solly,rhs0);
821 boost::timer::cpu_timer ctimer;
822 LinearSolverHA(ahh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol0,rhs0);
823 if(verbose) std::cout << "ERROR ESTIMATOR: computation time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
824 at_c<Traits::stateId>(errorEstimate_E.data) = at_c<0>(sol0.data);
825
826 /***************************************************************************************/
827 // Transfer error indicators to cells.
828 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
829 errorDistribution.clear();
830 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
831 Scalar errLevel(0.0), minRefine(0.2);
832
833 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
834 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
835 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
836 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
837 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
838 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
839
840 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
841 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
842 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
843 mde = distError;
844 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
845
846 if(verbose)
847 {
848 std::string name = "errorDistribution_";
849 name += std::to_string(step);
850 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
851 }
852 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
853 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
854 // Aufteilen von v_i auf die einzelnen Zellen
855 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
856 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
857 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
858
859 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
860
861 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
862 }
863
864 void refineGrid()
865 {
866 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
867 }
868
869 double estimatedAbsoluteError() const final
870 {
871 return sqrt(fabs(totalErrorSquared));
872 }
873
874 size_t gridSize() const final
875 {
876 return extensionSpace.gridManager().grid().size(0);
877 }
878
879 template <typename... Args>
880 void initFunctionals(const Args&... args)
881 {
882 F_HE.reset(new Functional_HE(args...));
883 F_EE.reset(new Functional_EE(args...));
884 }
885
886 private:
887 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
888 {
889 ass_HE->assemble(linearization(*F_HE,xl.get()));
890 ass_EE->assemble(linearization(*F_EE,xe.get()));
891 }
892
893 std::unique_ptr<Functional_HE> F_HE;
894 std::unique_ptr<Functional_EE> F_EE;
895 std::unique_ptr<Ass_HE> ass_HE;
896 std::unique_ptr<Ass_EE> ass_EE;
897 NormFunctional& normFunctional;
898 VariableSetDescription& variableSetDescription;
899 ExtensionVariableSetDescription& extensionVariableSetDescription;
900 ExtensionSpace& extensionSpace;
901 Scalar squaredFraction;
902 Scalar totalErrorSquared;
903 std::vector<std::pair<double,int> > errorDistribution;
904 bool verbose, fast;
905 };
906
907
908 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
909 class NormFunctional, class LinearSolverHA, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
910 class VariationalEquationHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
911 {
912 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
913 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
914 public:
915 // Functionals
916 typedef typename Traits::Functional_HE Functional_HE;
917 typedef typename Traits::Functional_EE Functional_EE;
918
919 // assemblers
920 typedef typename Traits::Assembler_HE Ass_HE;
921 typedef typename Traits::Assembler_EE Ass_EE;
922
923 // operators
924 typedef typename Traits::Luu_HE MLH;
925 typedef typename Traits::Luu_EE MHH;
926 typedef typename Traits::BT_HE BTLH;
927
928 typedef typename Traits::Scalar Scalar;
929 static constexpr int dim = VariableSetDescription::Grid::dimension;
930
931 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
932 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
933
934 VariationalEquationHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
935 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
936 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
937 extensionSpace(extensionSpace_),
938 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
939 {
940 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
941 ass_HE.reset(new Ass_HE(extensionVariableSetDescription.spaces));
942 }
943
944 virtual ~VariationalEquationHBErrorEstimator(){}
945
946 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
947 {
948 boost::timer::cpu_timer overallTimer;
949 using boost::fusion::at_c;
950 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
951 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
952 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
953 boost::timer::cpu_timer atimer;
954 assembleAt(xl,xe);
955 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
956 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
957
958 /***************************************************************************************/
959 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
960 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
961
962 MLH mlh(*ass_HE,false);
963 MHH mhh(*ass_EE,false);
964 BTLH btlh(*ass_HE,false);
965
966 typename Traits::ExtensionVectorU sol(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
967 typename Traits::VectorU sollu(at_c<Traits::controlId>(dxl.get().data.coefficients()));
968 typename Traits::VectorP sollp(at_c<Traits::adjointId>(dxl.get().data.coefficients()));
969 typename Traits::ExtensionVectorU rhs(ass_EE->template rhs<Traits::controlId,Traits::controlId+1>());
970 rhs *= -1.0;
971 mlh.applyscaleadd(-1.0,sollu,rhs);
972 btlh.applyscaleadd(-1.0,sollp,rhs);
973
974 LinearSolverHA(mhh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol,rhs);
975 at_c<Traits::control>(errorEstimate_E.data) = at_c<0>(sol.data);
976
977 /***************************************************************************************/
978 // Transfer error indicators to cells.
979 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
980 errorDistribution.clear();
981 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
982 Scalar errLevel(0.0), minRefine(0.2);
983
984 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
985 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
986 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
987 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
988 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
989 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
990
991 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
992 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
993 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
994 mde = distError;
995 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
996
997 if(verbose)
998 {
999 std::string name = "errorDistribution_";
1000 name += std::to_string(step);
1001 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1002 }
1003 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1004 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1005 // Aufteilen von v_i auf die einzelnen Zellen
1006 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
1007 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1008 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1009
1010 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1011
1012 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1013 }
1014
1015 void refineGrid()
1016 {
1017 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
1018 }
1019
1020 double estimatedAbsoluteError() const final
1021 {
1022 return sqrt(fabs(totalErrorSquared));
1023 }
1024
1025 size_t gridSize() const final
1026 {
1027 return extensionSpace.gridManager().grid().size(0);
1028 }
1029
1030 template <typename... Args>
1031 void initFunctionals(const Args&... args)
1032 {
1033 F_HE.reset(new Functional_HE(args...));
1034 F_EE.reset(new Functional_EE(args...));
1035 }
1036
1037 private:
1038 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
1039 {
1040 ass_HE->assemble(linearization(*F_HE,xl.get()));
1041 ass_EE->assemble(linearization(*F_EE,xe.get()));
1042 }
1043
1044 std::unique_ptr<Functional_HE> F_HE;
1045 std::unique_ptr<Functional_EE> F_EE;
1046 std::unique_ptr<Ass_HE> ass_HE;
1047 std::unique_ptr<Ass_EE> ass_EE;
1048 NormFunctional& normFunctional;
1049 VariableSetDescription& variableSetDescription;
1050 ExtensionVariableSetDescription& extensionVariableSetDescription;
1051 ExtensionSpace& extensionSpace;
1052 Scalar squaredFraction;
1053 Scalar totalErrorSquared;
1054 std::vector<std::pair<double,int> > errorDistribution;
1055 bool verbose, fast;
1056 };
1057
1058
1059 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
1060 class NormFunctional, class LinearSolverHA, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
1061 class AdjointEquationHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
1062 {
1063 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
1064 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
1065 public:
1066 // assemblers
1067 typedef typename Traits::Assembler_HE Ass_HE;
1068 typedef typename Traits::Assembler_EE Ass_EE;
1069
1070 // operators
1071 typedef typename Traits::Lyy_HE Lyy_HE;
1072 typedef typename Traits::AT_HE AT_HE;
1073 typedef typename Traits::AT_EE AT_EE;
1074
1075 typedef typename Traits::Scalar Scalar;
1076 static constexpr int dim = VariableSetDescription::Grid::dimension;
1077
1078 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
1079 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
1080
1081
1082
1083 AdjointEquationHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
1084 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
1085 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
1086 extensionSpace(extensionSpace_),
1087 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
1088 {
1089 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
1090 ass_HE.reset(new Ass_HE(extensionVariableSetDescription.spaces));
1091 }
1092
1093 virtual ~AdjointEquationHBErrorEstimator(){}
1094
1095 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
1096 {
1097 boost::timer::cpu_timer overallTimer;
1098 using boost::fusion::at_c;
1099 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
1100 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
1101 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
1102 boost::timer::cpu_timer atimer;
1103 assembleAt(xl,xe);
1104 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1105 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
1106
1107 /***************************************************************************************/
1108 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
1109 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
1110 Lyy_HE Lyy_HE(*ass_HE,false);
1111 AT_HE at_HE(*ass_HE,false);
1112 AT_EE at_EE(*ass_EE,false);
1113 typename Traits::ExtensionVectorY rhs0(ass_EE->template rhs<Traits::stateId,Traits::stateId+1>());
1114 rhs0 *= -1.0;
1115 typename Traits::ExtensionVectorP sol0(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
1116 typename Traits::VectorP sollp(Traits::VectorP_Initializer::init(variableSetDescription));
1117 typename Traits::VectorY solly(Traits::VectorY_Initializer::init(variableSetDescription));
1118 at_c<0>(solly.data) = at_c<Traits::stateId>(dxl.get().data).coefficients();
1119 at_c<0>(sollp.data) = at_c<Traits::adjointId>(dxl.get().data).coefficients();
1120 Lyy_HE.applyscaleadd(-1.0,solly,rhs0);
1121 at_HE.applyscaleadd(-1.0,sollp,rhs0);
1122 boost::timer::cpu_timer ctimer;
1123 LinearSolverHA(at_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol0,rhs0);
1124 if(verbose) std::cout << "ERROR ESTIMATOR: computation time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1125 at_c<Traits::stateId>(errorEstimate_E.data) = at_c<0>(sol0.data);
1126
1127 /***************************************************************************************/
1128 // Transfer error indicators to cells.
1129 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
1130 errorDistribution.clear();
1131 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
1132 Scalar errLevel(0.0), minRefine(0.2);
1133
1134 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1135 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1136 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
1137 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
1138 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
1139 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
1140
1141 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
1142 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
1143 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
1144 mde = distError;
1145 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
1146
1147 if(verbose)
1148 {
1149 std::string name = "errorDistribution_";
1150 name += std::to_string(step);
1151 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1152 }
1153 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1154 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1155 // Aufteilen von v_i auf die einzelnen Zellen
1156 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
1157 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1158 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1159
1160 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1161
1162 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1163 }
1164
1165 void refineGrid()
1166 {
1167 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
1168 }
1169
1170 double estimatedAbsoluteError() const final
1171 {
1172 return sqrt(fabs(totalErrorSquared));
1173 }
1174
1175 size_t gridSize() const final
1176 {
1177 return extensionSpace.gridManager().grid().size(0);
1178 }
1179
1180 template <typename... Args>
1181 void initFunctionals(const Args&... args)
1182 {
1183 F_HE.reset(new typename Traits::Functional_HE(args...));
1184 F_EE.reset(new typename Traits::Functional_EE(args...));
1185 }
1186
1187 private:
1188 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
1189 {
1190 ass_HE->assemble(linearization(*F_HE,xl.get()));
1191 ass_EE->assemble(linearization(*F_EE,xe.get()));
1192 }
1193
1194 std::unique_ptr<typename Traits::Functional_HE> F_HE;
1195 std::unique_ptr<typename Traits::Functional_EE> F_EE;
1196 std::unique_ptr<Ass_HE> ass_HE;
1197 std::unique_ptr<Ass_EE> ass_EE;
1198 NormFunctional& normFunctional;
1199 VariableSetDescription& variableSetDescription;
1200 ExtensionVariableSetDescription& extensionVariableSetDescription;
1201 ExtensionSpace& extensionSpace;
1202 Scalar squaredFraction;
1203 Scalar totalErrorSquared;
1204 std::vector<std::pair<double,int> > errorDistribution;
1205 bool verbose, fast;
1206 };
1207
1208
1209 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
1210 class NormFunctional, class LinearSolverLA, class LinearSolverHA, class LinearSolverLU, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
1211 class AdjointEquationLinearPropagationHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
1212 {
1213 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
1214 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
1215 public:
1216 // assemblers
1217
1218 // operators
1219 typedef typename Traits::Lyy_HE Lyy_HE;
1220 typedef typename Traits::AT_HE AT_HE;
1221 typedef typename Traits::AT_EE AT_EE;
1222
1223 typedef typename Traits::Scalar Scalar;
1224 static constexpr int dim = VariableSetDescription::Grid::dimension;
1225
1226 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
1227 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
1228
1229
1230
1231 AdjointEquationLinearPropagationHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
1232 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
1233 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
1234 extensionSpace(extensionSpace_),
1235 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
1236 {
1237 ass_EE.reset(new typename Traits::Assembler_EE(extensionVariableSetDescription.spaces));
1238 ass_HE.reset(new typename Traits::Assembler_HE(extensionVariableSetDescription.spaces));
1239 ass_EH.reset(new typename Traits::Assembler_EH(extensionVariableSetDescription.spaces));
1240 ass_HH.reset(new typename Traits::Assembler_HH(extensionVariableSetDescription.spaces));
1241 }
1242
1243 virtual ~AdjointEquationLinearPropagationHBErrorEstimator(){}
1244
1245 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
1246 {
1247 boost::timer::cpu_timer overallTimer;
1248 using boost::fusion::at_c;
1249 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
1250 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
1251 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
1252 boost::timer::cpu_timer atimer;
1253 assembleAt(xl,xe);
1254 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1255 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
1256
1257 /***************************************************************************************/
1261 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
1262 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
1263 Lyy_HE Lyy_HE(*ass_HE,false);
1264 AT_HE at_HE(*ass_HE,false);
1265 AT_EE at_EE(*ass_EE,false);
1266 typename Traits::ExtensionVectorY rhs0(ass_EE->template rhs<Traits::stateId,Traits::stateId+1>());
1267 rhs0 *= -1.0;
1268 typename Traits::ExtensionVectorP sol0(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
1269 typename Traits::VectorP sollp(Traits::VectorP_Initializer::init(variableSetDescription));
1270 typename Traits::VectorY solly(Traits::VectorY_Initializer::init(variableSetDescription));
1271 at_c<0>(solly.data) = at_c<Traits::stateId>(dxl.get().data).coefficients();
1272 at_c<0>(sollp.data) = at_c<Traits::adjointId>(dxl.get().data).coefficients();
1273 Lyy_HE.applyscaleadd(-1.0,solly,rhs0);
1274 at_HE.applyscaleadd(-1.0,sollp,rhs0);
1275 boost::timer::cpu_timer ctimer;
1276 LinearSolverHA(at_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol0,rhs0);
1277 if(verbose) std::cout << "ERROR ESTIMATOR: computation time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1278
1279 std::cout << "ERR4 sol0 " << (sol0*sol0) << std::endl;
1283 // variational equality
1284 typename Traits::VectorU rhsU(Traits::VectorU_Initializer::init(variableSetDescription));
1285 typename Traits::BT_EH bt_EH(*ass_EH,false);
1286 std::cout << "ERR4 rhsU0 " << (rhsU*rhsU) << std::endl;
1287 bt_EH.applyscaleadd(-1.0,sol0,rhsU);
1288 std::cout << "ERR4 rhsU1 " << (rhsU*rhsU) << std::endl;
1289
1290 typename Traits::Luu_HH luu(*ass_HH,false);
1291 typename Traits::VectorU solU(Traits::VectorU_Initializer::init(variableSetDescription));
1292 std::cout << "ERR4 solU0 " << (solU*solU) << std::endl;
1293 LinearSolverLU(luu,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solU,rhsU);
1294
1295 std::cout << "ERR4 solU " << (solU*solU) << std::endl;
1296 // state equation
1297 typename Traits::VectorP rhsP(Traits::VectorP_Initializer::init(variableSetDescription));
1298 typename Traits::B_HH b_HH(*ass_HH,false);
1299 b_HH.applyscaleadd(-1.0,solU,rhsP);
1300 typename Traits::A_HH a_HH(*ass_HH,false);
1301 typename Traits::VectorY solY(Traits::VectorY_Initializer::init(variableSetDescription));
1302 LinearSolverLA(a_HH,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solY,rhsP);
1303 std::cout << "ERR4 solY " << (solY*solY) << std::endl;
1304 at_c<Traits::stateId>(errorEstimate_H.data) = at_c<0>(solY.data);
1305 at_c<Traits::controlId>(errorEstimate_H.data) = at_c<0>(solU.data);
1306 /***************************************************************************************/
1307 // Transfer error indicators to cells.
1308 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
1309 errorDistribution.clear();
1310 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
1311 Scalar errLevel(0.0), minRefine(0.2);
1312
1313 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1314 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1315 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
1316 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
1317 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
1318 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
1319
1320 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
1321 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
1322 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
1323 mde = distError;
1324 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
1325
1326 if(verbose)
1327 {
1328 std::string name = "errorDistribution_";
1329 name += std::to_string(step);
1330 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1331 }
1332 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1333 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1334 // Aufteilen von v_i auf die einzelnen Zellen
1335 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
1336 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1337 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1338
1339 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1340
1341 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1342 }
1343
1344 void refineGrid()
1345 {
1346 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
1347 }
1348
1349 double estimatedAbsoluteError() const final
1350 {
1351 return sqrt(fabs(totalErrorSquared));
1352 }
1353
1354 size_t gridSize() const final
1355 {
1356 return extensionSpace.gridManager().grid().size(0);
1357 }
1358
1359 template <typename... Args>
1360 void initFunctionals(const Args&... args)
1361 {
1362 F_HH.reset(new typename Traits::Functional_HH(args...));
1363 F_HE.reset(new typename Traits::Functional_HE(args...));
1364 F_EH.reset(new typename Traits::Functional_EH(args...));
1365 F_EE.reset(new typename Traits::Functional_EE(args...));
1366 }
1367
1368 private:
1369 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
1370 {
1371 ass_HH->assemble(linearization(*F_HH,xl.get()));
1372 ass_HE->assemble(linearization(*F_HE,xl.get()));
1373 ass_EH->assemble(linearization(*F_EH,xe.get()));
1374 ass_EE->assemble(linearization(*F_EE,xe.get()));
1375 }
1376
1377 std::unique_ptr<typename Traits::Functional_HH> F_HH;
1378 std::unique_ptr<typename Traits::Functional_HE> F_HE;
1379 std::unique_ptr<typename Traits::Functional_EH> F_EH;
1380 std::unique_ptr<typename Traits::Functional_EE> F_EE;
1381 std::unique_ptr<typename Traits::Assembler_HH> ass_HH;
1382 std::unique_ptr<typename Traits::Assembler_HE> ass_HE;
1383 std::unique_ptr<typename Traits::Assembler_EH> ass_EH;
1384 std::unique_ptr<typename Traits::Assembler_EE> ass_EE;
1385 NormFunctional& normFunctional;
1386 VariableSetDescription& variableSetDescription;
1387 ExtensionVariableSetDescription& extensionVariableSetDescription;
1388 ExtensionSpace& extensionSpace;
1389 Scalar squaredFraction;
1390 Scalar totalErrorSquared;
1391 std::vector<std::pair<double,int> > errorDistribution;
1392 bool verbose, fast;
1393 };
1394
1395 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
1396 class NormFunctional, class LinearSolverLA, class LinearSolverHA, class LinearSolverHU, class LinearSolverLU, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration, bool lump=false>
1397 class AnotherHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
1398 {
1399 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,lump>;
1400 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
1401 public:
1402 // assemblers
1403
1404 // operators
1405 typedef typename Traits::Lyy_HE Lyy_HE;
1406 typedef typename Traits::AT_HE AT_HE;
1407 typedef typename Traits::AT_EE AT_EE;
1408
1409 typedef typename Traits::Scalar Scalar;
1410 static constexpr int dim = VariableSetDescription::Grid::dimension;
1411
1412 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
1413 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
1414
1415
1416
1417 AnotherHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
1418 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
1419 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
1420 extensionSpace(extensionSpace_),
1421 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
1422 {
1423 ass_EE.reset(new typename Traits::Assembler_EE(extensionVariableSetDescription.spaces));
1424 ass_HE.reset(new typename Traits::Assembler_HE(variableSetDescription.spaces));
1425 ass_EH.reset(new typename Traits::Assembler_EH(extensionVariableSetDescription.spaces));
1426 ass_HH.reset(new typename Traits::Assembler_HH(variableSetDescription.spaces));
1427 }
1428
1429 virtual ~AnotherHBErrorEstimator(){}
1430
1431 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
1432 {
1433 boost::timer::cpu_timer overallTimer;
1434 using boost::fusion::at_c;
1435 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
1436 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
1437 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
1438 boost::timer::cpu_timer atimer;
1439 assembleAt(xl,xe);
1440 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1441 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
1442
1443 std::cout << "ids: " << Traits::stateId << " " << Traits::controlId << " " << Traits::adjointId << std::endl;
1444 /***************************************************************************************/
1448 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
1449 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
1450 Lyy_HE Lyy_HE(*ass_HE,false);
1451 AT_HE at_HE(*ass_HE,false);
1452 AT_EE at_EE(*ass_EE,false);
1453 typename Traits::ExtensionVectorY rhs0(ass_EE->template rhs<Traits::stateId,Traits::stateId+1>());
1454 rhs0 *= -1.0;
1455 typename Traits::ExtensionVectorP sol0(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
1456 typename Traits::VectorP sollp(Traits::VectorP_Initializer::init(variableSetDescription));
1457 typename Traits::VectorY solly(Traits::VectorY_Initializer::init(variableSetDescription));
1458 at_c<0>(solly.data) = at_c<Traits::stateId>(dxl.get().data).coefficients();
1459 at_c<0>(sollp.data) = at_c<Traits::adjointId>(dxl.get().data).coefficients();
1460 Lyy_HE.applyscaleadd(-1.0,solly,rhs0);
1461 at_HE.applyscaleadd(-1.0,sollp,rhs0);
1462 boost::timer::cpu_timer ctimer;
1463 //JacobiPreconditioner<AT_EE>(at_EE).apply(sol0,rhs0);
1464 LinearSolverHA(at_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol0,rhs0);
1465 if(verbose) std::cout << "ERROR ESTIMATOR: computation time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1466
1467 std::cout << "|solPH|^2=" << (sol0*sol0) << std::endl;
1468
1472 // variational equality
1473 typename Traits::VectorU rhsUL(ass_HH->template rhs<Traits::controlId,Traits::controlId+1>());
1474 rhsUL *= -1.0;
1475 typename Traits::ExtensionVectorU rhsUH(ass_EE->template rhs<Traits::controlId,Traits::controlId+1>());
1476 rhsUH *= -1.0;
1477 std::cout << "|rhsUH|^2=" << (rhsUH*rhsUH) << std::endl;
1478 std::cout << "|rhsUL|^2=" << (rhsUL*rhsUL) << std::endl;
1479 typename Traits::BT_EH bt_EH(*ass_EH,false);
1480 typename Traits::BT_EE bt_EE(*ass_EE,false);
1481 bt_EH.applyscaleadd(-1.0,sol0,rhsUL);
1482 bt_EE.applyscaleadd(-1.0,sol0,rhsUH);
1483 std::cout << "|rhsUH|^2=" << (rhsUH*rhsUH) << std::endl;
1484 std::cout << "|rhsUL|^2=" << (rhsUL*rhsUL) << std::endl;
1485 typename Traits::Luu_EE luu_EE(*ass_EE,true);
1486 typename Traits::Luu_EH luu_EH(*ass_EH,false);
1487 typename Traits::Luu_HH luu_HH(*ass_HH,false);
1488 typename Traits::VectorU solUL(Traits::VectorU_Initializer::init(variableSetDescription));
1489 typename Traits::ExtensionVectorU solUH(Traits::ExtensionVectorU_Initializer::init(extensionVariableSetDescription));
1490 std::cout << "computing u_h" << std::endl;
1491 JacobiPreconditioner<typename Traits::Luu_EE>(luu_EE).apply(solUH,rhsUH);//LinearSolverHU(luu_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solUH,rhsUH);
1492 std::cout << "|solUH|^2=" << (solUH*solUH) << std::endl;
1493 //luu_EH.applyscaleadd(-1.0,solUH,rhsUL);
1494 std::cout << "|rhsUL|^2=" << (rhsUL*rhsUL) << std::endl;
1495 std::cout << "computing u_l" << std::endl;
1496 JacobiPreconditioner<typename Traits::Luu_HH>(luu_HH).apply(solUL,rhsUL);//LinearSolverLU(luu_HH,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solUL,rhsUL);
1497 std::cout << "|solUH|^2=" << (solUH*solUH) << std::endl;
1498 std::cout << "|solUL|^2=" << (solUL*solUL) << std::endl;
1499
1500 // state equation
1501 typename Traits::VectorP rhsPL(ass_HH->template rhs<Traits::adjointId,Traits::adjointId+1>());
1502 rhsPL *= -1.0;
1503 typename Traits::ExtensionVectorP rhsPH(ass_EE->template rhs<Traits::adjointId,Traits::adjointId+1>());
1504 rhsPH *= -1.0;
1505 typename Traits::B_HH b_HH(*ass_HH,false);
1506 typename Traits::B_HE b_HE(*ass_HE,false);
1507 typename Traits::B_EH b_EH(*ass_EH,false);
1508 typename Traits::B_EE b_EE(*ass_EE,false);
1509 b_HH.applyscaleadd(-1.0,solUL,rhsPL);
1510 b_HE.applyscaleadd(-1.0,solUL,rhsPH);
1511 b_EH.applyscaleadd(-1.0,solUH,rhsPL);
1512 b_HH.applyscaleadd(-1.0,solUH,rhsPH);
1513
1514 typename Traits::A_HH a_HH(*ass_HH,false);
1515 typename Traits::A_EH a_EH(*ass_EH,false);
1516 typename Traits::A_EE a_EE(*ass_EE,false);
1517
1518 typename Traits::VectorY solYL(Traits::VectorY_Initializer::init(variableSetDescription));
1519 typename Traits::ExtensionVectorY solYH(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
1520 JacobiPreconditioner<typename Traits::A_EE>(a_EE).apply(solYH,rhsPH);//LinearSolverHA(a_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solYH,rhsPH);
1521
1522 a_EH.applyscaleadd(-1.0,solYH,rhsPL);
1523 JacobiPreconditioner<typename Traits::A_HH>(a_HH).apply(solYL,rhsPL);//LinearSolverLA(a_HH,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(solYL,rhsPL);
1524 std::cout << "|solYH|^2=" << (solYH*solYH) << std::endl;
1525 std::cout << "|solYL|^2=" << (solYL*solYL) << std::endl;
1526 at_c<Traits::stateId>(errorEstimate_H.data) = at_c<0>(solYL.data);
1527 at_c<Traits::controlId>(errorEstimate_H.data) = at_c<0>(solUL.data);
1528 at_c<Traits::stateId>(errorEstimate_E.data) = at_c<0>(solYH.data);
1529 at_c<Traits::controlId>(errorEstimate_E.data) = at_c<0>(solUH.data);
1530 /***************************************************************************************/
1531 // Transfer error indicators to cells.
1532 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
1533 errorDistribution.clear();
1534 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
1535 Scalar errLevel(0.0), minRefine(0.2);
1536
1537 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1538 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1539 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
1540 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
1541 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
1542 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
1543
1544 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
1545 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
1546 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
1547 mde = distError;
1548 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
1549
1550 if(verbose)
1551 {
1552 std::string name = "errorDistribution_";
1553 name += std::to_string(step);
1554 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1555 }
1556 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1557 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1558 // Aufteilen von v_i auf die einzelnen Zellen
1559 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
1560 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1561 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1562
1563 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1564
1565 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1566 }
1567
1568 void refineGrid()
1569 {
1570 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,squaredFraction*totalErrorSquared);
1571 }
1572
1573 double estimatedAbsoluteError() const final
1574 {
1575 return sqrt(fabs(totalErrorSquared));
1576 }
1577
1578 size_t gridSize() const final
1579 {
1580 return extensionSpace.gridManager().grid().size(0);
1581 }
1582
1583 template <typename... Args>
1584 void initFunctionals(const Args&... args)
1585 {
1586 F_HH.reset(new typename Traits::Functional_HH(args...));
1587 F_HE.reset(new typename Traits::Functional_HE(args...));
1588 F_EH.reset(new typename Traits::Functional_EH(args...));
1589 F_EE.reset(new typename Traits::Functional_EE(args...));
1590 }
1591
1592 private:
1593 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
1594 {
1595 ass_HH->assemble(linearization(*F_HH,xl.get()));
1596 ass_HE->assemble(linearization(*F_HE,xl.get()));
1597 ass_EH->assemble(linearization(*F_EH,xe.get()));
1598 ass_EE->assemble(linearization(*F_EE,xe.get()));
1599 }
1600
1601 std::unique_ptr<typename Traits::Functional_HH> F_HH;
1602 std::unique_ptr<typename Traits::Functional_HE> F_HE;
1603 std::unique_ptr<typename Traits::Functional_EH> F_EH;
1604 std::unique_ptr<typename Traits::Functional_EE> F_EE;
1605 std::unique_ptr<typename Traits::Assembler_HH> ass_HH;
1606 std::unique_ptr<typename Traits::Assembler_HE> ass_HE;
1607 std::unique_ptr<typename Traits::Assembler_EH> ass_EH;
1608 std::unique_ptr<typename Traits::Assembler_EE> ass_EE;
1609 NormFunctional& normFunctional;
1610 VariableSetDescription& variableSetDescription;
1611 ExtensionVariableSetDescription& extensionVariableSetDescription;
1612 ExtensionSpace& extensionSpace;
1613 Scalar squaredFraction;
1614 Scalar totalErrorSquared;
1615 std::vector<std::pair<double,int> > errorDistribution;
1616 bool verbose, fast;
1617 };
1618
1619
1620 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
1621 class NormFunctional, class LinearSolverHA, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
1622 class StupidHBErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
1623 {
1624 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
1625 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
1626 public:
1627 // Functionals
1628 typedef typename Traits::Functional_HH Functional_HH;
1629 typedef typename Traits::Functional_HE Functional_HE;
1630 typedef typename Traits::Functional_EE Functional_EE;
1631
1632 // assemblers
1633 typedef typename Traits::Assembler_HE Ass_HE;
1634 typedef typename Traits::Assembler_EE Ass_EE;
1635
1636 // operators
1637 typedef typename Traits::Operator_HE A_HE;
1638 typedef typename Traits::Operator_EE A_EE;
1639
1640 typedef typename Traits::Scalar Scalar;
1641 static constexpr int dim = VariableSetDescription::Grid::dimension;
1642
1643 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
1644 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
1645
1646
1647
1648 StupidHBErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
1649 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
1650 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
1651 extensionSpace(extensionSpace_),
1652 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
1653 {
1654 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
1655 ass_HE.reset(new Ass_HE(variableSetDescription.spaces));
1656 }
1657
1658 virtual ~StupidHBErrorEstimator(){}
1659
1660 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
1661 {
1662 boost::timer::cpu_timer overallTimer;
1663 using boost::fusion::at_c;
1664 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
1665 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
1666 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
1667 boost::timer::cpu_timer atimer;
1668 assembleAt(xl,xe);
1669 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1670 Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
1671
1672 /***************************************************************************************/
1673 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
1674 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
1675 A_HE alh(*ass_HE,false);
1676 A_EE ahh(*ass_EE,false);
1677 typename Traits::ExtensionVector rhs_H(ass_EE->rhs());
1678 rhs_H *= -1.0;
1679 typename Traits::ExtensionVector sol_H(Traits::ExtensionVector_Initializer::init(extensionVariableSetDescription));
1680 typename Traits::Vector sol_L(Traits::Vector_Initializer::init(variableSetDescription));
1681 sol_L = dxl.get();
1682 alh.applyscaleadd(-1.0,sol_L,rhs_H);
1683 boost::timer::cpu_timer ctimer;
1684 LinearSolverHA(ahh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol_H,rhs_H);
1685 if(verbose) std::cout << "ERROR ESTIMATOR: computation time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1686 errorEstimate_E = sol_H;
1687
1688 /***************************************************************************************/
1689 // Transfer error indicators to cells.
1690 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
1691 errorDistribution.clear();
1692 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
1693 Scalar errLevel(0.0), minRefine(0.2);
1694
1695 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1696 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1697 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
1698 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
1699 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
1700 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
1701
1702 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
1703 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
1704 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
1705 mde = distError;
1706 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
1707
1708 if(verbose)
1709 {
1710 std::string name = "errorDistribution_";
1711 name += std::to_string(step);
1712 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1713 }
1714 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1715 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1716 // Aufteilen von v_i auf die einzelnen Zellen
1717 CellIterator cend = extensionVariableSetDescription.gridView.template end<0>();
1718 for (CellIterator ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1719 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1720
1721 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1722
1723 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1724 }
1725
1726 void refineGrid()
1727 {
1728 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
1729 }
1730
1731 double estimatedAbsoluteError() const final
1732 {
1733 return sqrt(fabs(totalErrorSquared));
1734 }
1735
1736 size_t gridSize() const final
1737 {
1738 return extensionSpace.gridManager().grid().size(0);
1739 }
1740
1741 template <typename... Args>
1742 void initFunctionals(const Args&... args)
1743 {
1744 F_HE.reset(new Functional_HE(args...));
1745 F_EE.reset(new Functional_EE(args...));
1746 }
1747
1748 private:
1749 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
1750 {
1751 ass_HE->assemble(linearization(*F_HE,xl.get()));
1752 ass_EE->assemble(linearization(*F_EE,xe.get()));
1753 }
1754
1755 std::unique_ptr<Functional_HE> F_HE;
1756 std::unique_ptr<Functional_EE> F_EE;
1757 std::unique_ptr<Ass_HE> ass_HE;
1758 std::unique_ptr<Ass_EE> ass_EE;
1759 NormFunctional& normFunctional;
1760 VariableSetDescription& variableSetDescription;
1761 ExtensionVariableSetDescription& extensionVariableSetDescription;
1762 ExtensionSpace& extensionSpace;
1763 Scalar squaredFraction;
1764 Scalar totalErrorSquared;
1765 std::vector<std::pair<double,int> > errorDistribution;
1766 bool verbose, fast;
1767 };
1768
1769
1770
1771 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
1772 class NormFunctional, class LinearSystemSolver_H, class LinearSystemSolver_L, class LinearSolverA_H, class LinearSolverA_L, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
1773 class MartinsErrorEstimator : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
1774 {
1775 template <class AnsatzVars, class TestVars, class OriginVars> using EstimatorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
1776 template <class AnsatzVars, class TestVars, class OriginVars> using LumpedFunctional = Functional<AnsatzVars,TestVars,OriginVars,true>;
1777 typedef ErrorEstimator_Detail::Traits<EstimatorFunctional,VariableSetDescription,ExtensionVariableSetDescription> FullTraits;
1778 typedef ErrorEstimator_Detail::Traits<LumpedFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
1779 public:
1780 // Functionals
1781 typedef typename FullTraits::Functional_HH Functional_HH;
1782 typedef typename Traits::Functional_HE Functional_HE;
1783 typedef typename Traits::Functional_EH Functional_EH;
1784 typedef typename Traits::Functional_EE Functional_EE;
1785
1786 // assemblers
1787 typedef typename Traits::Assembler_EE Ass_EE;
1788 typedef typename Traits::Assembler_HE Ass_HE;
1789 typedef typename Traits::Assembler_EH Ass_EH;
1790 typedef typename FullTraits::Assembler_HH Ass_HH;
1791
1792 // operators
1793 typedef typename FullTraits::Operator_HH H_HH;
1794 typedef typename Traits::Operator_HE H_HE;
1795 typedef typename Traits::Operator_EH H_EH;
1796 typedef typename Traits::Operator_EE H_EE;
1797
1798 typedef typename FullTraits::B_HH B_HH;
1799 typedef typename Traits::B_HE B_HE;
1800 typedef typename Traits::B_EH B_EH;
1801 typedef typename Traits::B_EE B_EE;
1802
1803 typedef typename FullTraits::Lyy_HH Hyy_HH;
1804 typedef typename Traits::Lyy_EH Hyy_EH;
1805 typedef typename Traits::Lyy_HE Hyy_HE;
1806 typedef typename Traits::Lyy_EE Hyy_EE;
1807
1808 typedef typename FullTraits::A_HH A_HH;
1809 typedef typename Traits::A_HE A_HE;
1810 typedef typename Traits::A_EE A_EE;
1811
1812 typedef typename FullTraits::AT_HH AT_HH;
1813 typedef typename Traits::AT_HE AT_HE;
1814 typedef typename Traits::AT_EE AT_EE;
1815
1816 typedef typename Traits::Scalar Scalar;
1817 static constexpr int dim = VariableSetDescription::Grid::dimension;
1818
1819 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Iterator CellIterator;
1820 typedef typename ExtensionVariableSetDescription::GridView::template Codim<0>::Entity Cell;
1821
1822
1823
1824 MartinsErrorEstimator(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
1825 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
1826 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
1827 extensionSpace(extensionSpace_),
1828 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
1829 {
1830 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
1831 ass_HE.reset(new Ass_HE(variableSetDescription.spaces));
1832 ass_EH.reset(new Ass_EH(extensionVariableSetDescription.spaces));
1833 ass_HH.reset(new Ass_HH(variableSetDescription.spaces));
1834 }
1835
1836 virtual ~MartinsErrorEstimator(){}
1837
1838 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
1839 {
1840 boost::timer::cpu_timer overallTimer;
1841 using boost::fusion::at_c;
1842 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
1843 Bridge::Vector<typename Traits::VarSet> const& xl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(x_);
1844 Bridge::Vector<typename Traits::ExtensionVarSet> xe(extensionVariableSetDescription);
1845 boost::timer::cpu_timer atimer;
1846 assembleAt(xl,xe);
1847 if(verbose) std::cout << "ERROR ESTIMATOR: assembly time: " << boost::timer::format(atimer.elapsed(), 2, "%t") << std::endl;
1848// Bridge::Vector<typename Traits::VarSet> const& dxl = dynamic_cast<const Bridge::Vector<typename Traits::VarSet>&>(dx_);
1849
1850 /***************************************************************************************/
1851 H_EE h_EE(*ass_EE,false);
1852 H_EH h_EH(*ass_EH,false);
1853 H_HH h_HH(*ass_HH,false);
1854 H_HE h_HE(*ass_HE,false);
1855
1856 typename Traits::ExtensionVector rhs_h(ass_EE->rhs());
1857 rhs_h *= -1.0;
1858 typename Traits::ExtensionVector sol_h(Traits::ExtensionVector_Initializer::init(extensionVariableSetDescription));
1859
1860 LinearSystemSolver_H(h_EE,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol_h,rhs_h);
1861
1862 typename Traits::Vector rhs_l(Traits::Vector_Initializer::init(variableSetDescription)),
1863 sol_l(Traits::Vector_Initializer::init(variableSetDescription));
1864 h_EH.applyscaleadd(-1.0,sol_h,rhs_l);
1865
1866 LinearSystemSolver_L(h_HH,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(sol_l,rhs_l);
1867
1868 // compute w_y
1869 B_HH bll(*ass_HH,false);
1870 B_HE blh(*ass_HE,false);
1871 B_EH bhl(*ass_EH,false);
1872 B_EE bhh(*ass_EE,false);
1873 A_HH all(*ass_HH,false);
1874 A_HE alh(*ass_HE,false);
1875 A_EE ahh(*ass_EE,false);
1876
1877 typename Traits::VectorP rhsP_l(Traits::VectorP_Initializer::init(variableSetDescription));
1878 typename Traits::ExtensionVectorP rhsP_h(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
1879
1880 typename Traits::VectorU wu_l(Traits::VectorU_Initializer::init(variableSetDescription));
1881 typename Traits::ExtensionVectorU wu_h(Traits::ExtensionVectorU_Initializer::init(extensionVariableSetDescription));
1882
1883 at_c<0>(wu_l.data) = at_c<Traits::controlId>(sol_l.data);
1884 at_c<0>(wu_h.data) = at_c<Traits::controlId>(sol_h.data);
1885 bhh.applyscaleadd(-1.0,wu_h,rhsP_h);
1886 blh.applyscaleadd(-1.0,wu_l,rhsP_h);
1887 bhl.applyscaleadd(-1.0,wu_h,rhsP_l);
1888 bll.applyscaleadd(-1.0,wu_l,rhsP_l);
1889
1890 typename Traits::VectorY wy_l(Traits::VectorY_Initializer::init(variableSetDescription));
1891 typename Traits::ExtensionVectorY wy_h(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
1892
1893 LinearSolverA_L(all,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(wy_l,rhsP_l);
1894 alh.applyscaleadd(-1.0,wy_l,rhsP_h);
1895 LinearSolverA_H(ahh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(wy_h,rhsP_h);
1896
1897 // compute w_y
1898 typename Traits::VectorY rhsY_l(Traits::VectorY_Initializer::init(variableSetDescription));
1899 typename Traits::ExtensionVectorY rhsY_h(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
1900
1901 typename Traits::VectorP wp_l(Traits::VectorP_Initializer::init(variableSetDescription));
1902 typename Traits::ExtensionVectorP wp_h(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
1903
1904 Hyy_HH hyyll(*ass_HH,false);
1905 Hyy_EH hyyhl(*ass_EH,false);
1906 Hyy_HE hyylh(*ass_HE,false);
1907 Hyy_EE hyyhh(*ass_EE,false);
1908
1909 AT_HH atll(*ass_HH,false);
1910 AT_HE atlh(*ass_HE,false);
1911 AT_EE athh(*ass_EE,false);
1912
1913 hyyhh.applyscaleadd(-1.0,wy_h,rhsY_h);
1914 hyylh.applyscaleadd(-1.0,wy_l,rhsY_h);
1915 hyyhl.applyscaleadd(-1.0,wy_h,rhsY_l);
1916 hyyll.applyscaleadd(-1.0,wy_l,rhsY_l);
1917
1918 LinearSolverA_L(atll,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(wp_l,rhsY_l);
1919 atlh.applyscaleadd(-1.0,wp_l,rhsY_h);
1920 LinearSolverA_H(athh,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(wp_h,rhsY_h);
1921
1922 // compute error indicator
1923 typename Traits::ExtensionVector weights(Traits::ExtensionVector_Initializer::init(extensionVariableSetDescription));
1924 at_c<Traits::controlId>(weights.data) = at_c<0>(wu_h.data);
1925 at_c<Traits::stateId>(weights.data) = at_c<0>(wy_h.data);
1926 at_c<Traits::adjointId>(weights.data) = at_c<0>(wp_h.data);
1927
1928 typename Traits::ExtensionVector errorIndicator(Traits::ExtensionVector_Initializer::init(extensionVariableSetDescription));
1929 typename Traits::Vector tmp(Traits::Vector_Initializer::init(variableSetDescription)),
1930 tmp2(Traits::Vector_Initializer::init(variableSetDescription));
1931
1932 h_EE.apply(weights,errorIndicator);
1933 h_EH.apply(weights,tmp);
1934 LinearSystemSolver_L(h_HH,DirectType::UMFPACK3264,MatrixProperties::GENERAL).apply(tmp2,tmp);
1935 h_HE.applyscaleadd(-1.0,tmp2,errorIndicator);
1936
1937 typename ExtensionVariableSetDescription::VariableSet errorEstimate_E(extensionVariableSetDescription);
1938 typename VariableSetDescription::VariableSet errorEstimate_H(variableSetDescription); // estimate error in state variable
1939 errorEstimate_E = errorIndicator;
1940
1941 /***************************************************************************************/
1942 // Transfer error indicators to cells.
1943 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
1944 errorDistribution.clear();
1945 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
1946 Scalar errLevel(0.0), minRefine(0.2);
1947
1948 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1949 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
1950 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
1951 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
1952 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
1953 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
1954
1955 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
1956 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
1957 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
1958 mde = distError;
1959 std::cout << "nError: " << distError.dim() << std::endl;
1960 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
1961
1962 if(verbose)
1963 {
1964 std::string name = "errorDistribution_";
1965 name += std::to_string(step);
1966 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
1967 }
1968 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
1969 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
1970 // Aufteilen von v_i auf die einzelnen Zellen
1971 CellIterator cend = variableSetDescription.gridView.template end<0>();
1972 for (CellIterator ci=variableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
1973 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
1974
1975 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
1976
1977 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
1978 }
1979
1980 void refineGrid()
1981 {
1982 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
1983 }
1984
1985 double estimatedAbsoluteError() const final
1986 {
1987 return sqrt(fabs(totalErrorSquared));
1988 }
1989
1990 size_t gridSize() const final
1991 {
1992 return extensionSpace.gridManager().grid().size(0);
1993 }
1994
1995 template <typename... Args>
1996 void initFunctionals(const Args&... args)
1997 {
1998 F_HE.reset(new Functional_HE(args...));
1999 F_EH.reset(new Functional_EH(args...));
2000 F_EE.reset(new Functional_EE(args...));
2001 F_HH.reset(new Functional_HH(args...));
2002 }
2003
2004 private:
2005 void assembleAt(const Bridge::Vector<typename Traits::VarSet>& xl, const Bridge::Vector<typename Traits::ExtensionVarSet>& xe)
2006 {
2007 ass_HE->assemble(linearization(*F_HE,xl.get()));
2008 ass_EH->assemble(linearization(*F_EH,xe.get()));
2009 ass_EE->assemble(linearization(*F_EE,xe.get()));
2010 ass_HH->assemble(linearization(*F_HH,xl.get()));
2011 }
2012
2013 std::unique_ptr<Functional_HE> F_HE;
2014 std::unique_ptr<Functional_EH> F_EH;
2015 std::unique_ptr<Functional_EE> F_EE;
2016 std::unique_ptr<Functional_HH> F_HH;
2017 std::unique_ptr<Ass_HE> ass_HE;
2018 std::unique_ptr<Ass_EH> ass_EH;
2019 std::unique_ptr<Ass_EE> ass_EE;
2020 std::unique_ptr<Ass_HH> ass_HH;
2021 NormFunctional& normFunctional;
2022 VariableSetDescription& variableSetDescription;
2023 ExtensionVariableSetDescription& extensionVariableSetDescription;
2024 ExtensionSpace& extensionSpace;
2025 Scalar squaredFraction;
2026 Scalar totalErrorSquared;
2027 std::vector<std::pair<double,int> > errorDistribution;
2028 bool verbose, fast;
2029 };
2030
2031
2032
2033 template <template <class,class,class,bool> class Functional, class VariableSetDescription, class ExtensionVariableSetDescription, class ExtensionSpace,
2034 class NormFunctional, class LinearSolverLA, class LinearSolverHA, class LinearSolverLM=LinearSolverLA, class LinearSolverHM=LinearSolverHA,
2035 bool lumpM=false, template <class> class RefinementStrategy = Adaptivity::ErrorEquilibration>
2036 class HierarchicalBasisErrorEstimator3 : public AbstractHierarchicalErrorEstimator, public RefinementStrategy<typename VariableSetDescription::Grid>
2037 {
2038 static constexpr int noOfVariables = ExtensionVariableSetDescription::noOfVariables;
2039 template <class AnsatzVars, class TestVars, class OriginVars> using ErrorFunctional = Functional<AnsatzVars,TestVars,OriginVars,false>;
2040
2041 typedef ErrorEstimator_Detail::Traits<ErrorFunctional,VariableSetDescription,ExtensionVariableSetDescription> Traits;
2042 static constexpr int controlId = Traits::controlId;
2043 static constexpr int stateId = Traits::stateId;
2044 static constexpr int adjointId = Traits::adjointId;
2045
2046 public:
2047 typedef typename Traits::Scalar Scalar;
2048 typedef Functional<VariableSetDescription,VariableSetDescription,VariableSetDescription,lumpM> Functional_HH;
2049 typedef Functional<ExtensionVariableSetDescription,VariableSetDescription,ExtensionVariableSetDescription,lumpM> Functional_EH;
2050 typedef Functional<VariableSetDescription,ExtensionVariableSetDescription,VariableSetDescription,lumpM> Functional_HE;
2051 typedef Functional<ExtensionVariableSetDescription,ExtensionVariableSetDescription,ExtensionVariableSetDescription,lumpM> Functional_EE;
2052 typedef HierarchicErrorEstimator<LinearizationAt<Functional_HH>,ExtensionVariableSetDescription,ExtensionVariableSetDescription> ErrorEstimator;
2053 // typedef HierarchicErrorEstimator<LinearizationAt<Functional_HH>,ExtensionVariableSetDescription,ExtensionVariableSetDescription,HierarchicErrorEstimatorDetail::TakeAllD2<LinearizationAt<Functional_HH> > > ErrorEstimator;
2054 typedef VariationalFunctionalAssembler<ErrorEstimator> Assembler;
2055 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HH> > Ass_HH;
2056 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_HE> > Ass_HE;
2057 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EH> > Ass_EH;
2058 typedef VariationalFunctionalAssembler<LinearizationAt<Functional_EE> > Ass_EE;
2059
2060 typedef typename Traits::AT_EE AT_EE;
2061 typedef typename Traits::AT_HE AT_HE;
2062
2063 typedef typename Traits::A_HH A_HH;
2064 typedef typename Traits::A_EE A_EE;
2065 typedef typename Traits::A_EH A_EH;
2066
2067 typedef typename Traits::Luu_HH M_HH;
2068 typedef typename Traits::Luu_EE M_EE;
2069
2070 typedef typename Traits::BT_EH BT_EH;
2071 typedef typename Traits::BT_EE BT_EE;
2072
2073 typedef typename Traits::B_HH B_HH;
2074 typedef typename Traits::B_HE B_HE;
2075 typedef typename Traits::B_EH B_EH;
2076 typedef typename Traits::B_EE B_EE;
2077
2078 static constexpr int dim = Traits::dim;
2079
2080 HierarchicalBasisErrorEstimator3(NormFunctional& normFunctional_, VariableSetDescription& variableSetDescription_, ExtensionVariableSetDescription& extensionVariableSetDescription_,
2081 ExtensionSpace& extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_ = false)
2082 : RefinementStrategy<typename VariableSetDescription::Grid>(extensionSpace_.gridManager(),fraction,verbose),/*lf(lf_), lvars(lvars_),*/ normFunctional(normFunctional_), variableSetDescription(variableSetDescription_), extensionVariableSetDescription(extensionVariableSetDescription_),
2083 extensionSpace(extensionSpace_), assembler(extensionVariableSetDescription.spaces),
2084 squaredFraction(fraction*fraction), verbose(verbose_), fast(fast_)
2085 {
2086 initAssemblers();
2087 }
2088
2089 virtual ~HierarchicalBasisErrorEstimator3(){}
2090
2091 void operator()(AbstractVector const& x_, AbstractVector const& dx_, int step, AbstractVector const&)
2092 {
2093 boost::timer::cpu_timer overallTimer;
2094 using boost::fusion::at_c;
2095 if(verbose) std::cout << "ERROR ESTIMATOR: Start." << std::endl;
2096 Bridge::Vector<typename VariableSetDescription::VariableSet> const& xl = dynamic_cast<const Bridge::Vector<typename VariableSetDescription::VariableSet>&>(x_);
2097 Bridge::Vector<typename ExtensionVariableSetDescription::VariableSet> xe(extensionVariableSetDescription);
2098 assembleAt(xl,xe);
2099 Bridge::Vector<typename VariableSetDescription::VariableSet> const& dxl = dynamic_cast<const Bridge::Vector<typename VariableSetDescription::VariableSet>&>(dx_);
2100
2101 /***************************************************************************************/
2102 // estimate error in dual variable
2103 typename ExtensionVariableSetDescription::VariableSet tmpRep(extensionVariableSetDescription), errorEstimate_E(extensionVariableSetDescription);
2104 typename VariableSetDescription::VariableSet tmpRepL(variableSetDescription), errorEstimate_H(variableSetDescription);
2105
2106 AT_EE atee(*ass_EE,false);
2107 AT_HE athe(*ass_HE,false);
2108
2109 typename Traits::ExtensionVectorP dpe(Traits::ExtensionVectorP_Initializer::init(extensionVariableSetDescription));
2110// typename Traits::VectorP dph(at_c<adjointId>(dxl.get().data));
2111 typename Traits::VectorP dph(Traits::VectorP_Initializer::init(variableSetDescription));
2112 at_c<0>(dph.data) = at_c<adjointId>(dxl.get().data).coefficients();
2113 typename Traits::ExtensionVectorY rye(ass_EE->template rhs<stateId,stateId+1>());
2114 rye *= -1.;
2115
2116 athe.applyscaleadd(-1.0,dph,rye);
2117 JacobiPreconditioner<AT_EE>(atee).apply(dpe,rye);
2118
2119 at_c<adjointId>(errorEstimate_E.data) = at_c<0>(dpe.data);
2120
2121 /***************************************************************************************/
2122 // estimate error in control variable
2123
2124 M_HH mhh(*ass_HH,false);
2125 M_EE mee(*ass_EE,false);
2126
2127 BT_EH bteh(*ass_EH,false);
2128 BT_EE btee(*ass_EE,false);
2129
2130 typename Traits::VectorU ruh(ass_HH->template rhs<controlId,controlId+1>()),
2131 duh(Traits::VectorU_Initializer::init(variableSetDescription));
2132 typename Traits::ExtensionVectorU rue(ass_EE->template rhs<controlId,controlId+1>()),
2133 due(Traits::ExtensionVectorU_Initializer::init(extensionVariableSetDescription));
2134 ruh *= -1.;
2135 rue *= -1.;
2136
2137 bteh.applyscaleadd(-1.0,dpe,ruh);
2138 btee.applyscaleadd(-1.0,dpe,rue);
2139
2140 JacobiPreconditioner<M_HH>(mhh).apply(duh,ruh);
2141 JacobiPreconditioner<M_EE>(mee).apply(due,rue);
2142
2143 at_c<controlId>(errorEstimate_E.data) = at_c<0>(due.data);
2144 at_c<controlId>(errorEstimate_H.data) = at_c<0>(duh.data);
2145 /***************************************************************************************/
2146 // estimate error in state variable
2147
2148 B_HH bhh(*ass_HH,false);
2149 B_EH beh(*ass_EH,false);
2150 B_HE bhe(*ass_HE,false);
2151 B_EE bee(*ass_EE,false);
2152
2153 A_HH ahh(*ass_HH,false);
2154 A_EE aee(*ass_EE,false);
2155 A_EH aeh(*ass_EH,false);
2156
2157 typename Traits::VectorP rph(ass_HH->template rhs<adjointId,adjointId+1>());
2158 typename Traits::ExtensionVectorP rpe(ass_EE->template rhs<adjointId,adjointId+1>());
2159 rph *= -1.;
2160 rpe *= -1.;
2161
2162 bhh.applyscaleadd(-1.0,duh,rph);
2163 bhe.applyscaleadd(-1.0,duh,rpe);
2164 beh.applyscaleadd(-1.0,due,rph);
2165 bee.applyscaleadd(-1.0,due,rpe);
2166
2167 typename Traits::VectorY dyh(Traits::VectorY_Initializer::init(variableSetDescription));
2168 typename Traits::ExtensionVectorY dye(Traits::ExtensionVectorY_Initializer::init(extensionVariableSetDescription));
2169
2170 JacobiPreconditioner<A_EE>(aee).apply(dye,rpe);
2171 aeh.applyscaleadd(-1.0,dye,rph);
2172 JacobiPreconditioner<A_HH>(ahh).apply(dyh,rph);
2173
2174 at_c<stateId>(errorEstimate_E.data) = at_c<0>(dye.data);
2175 at_c<stateId>(errorEstimate_H.data) = at_c<0>(dyh.data);
2176 /***************************************************************************************/
2177 // Transfer error indicators to cells.
2178 auto const& is = extensionSpace.gridManager().grid().leafIndexSet();
2179 errorDistribution.clear();
2180 errorDistribution.resize(is.size(0),std::make_pair(0.0,0));
2181 Scalar errLevel(0.0), minRefine(0.2);
2182
2183 // typedef ErrorDistribution<NormFunctional,ExtensionVariableSetDescription> EnergyError;
2184 typedef ErrorDistribution3<NormFunctional,ExtensionVariableSetDescription> EnergyError;
2185 // EnergyError energyError(normFunctional,xl.get(),errorEstimate_E);
2186 EnergyError energyError(normFunctional,xl.get(),errorEstimate_H,errorEstimate_E);
2187 typedef VariationalFunctionalAssembler<LinearizationAt<EnergyError> > EnergyErrorAssembler;
2188 EnergyErrorAssembler eeAssembler(extensionSpace.gridManager(), energyError.getSpaces());
2189
2190 eeAssembler.assemble(linearization(energyError,xl.get()), EnergyErrorAssembler::RHS);
2191 typename EnergyError::ErrorVector distError( eeAssembler.rhs() );
2192 typename EnergyError::AnsatzVars::VariableSet mde( energyError.getVariableSetDescription() );
2193 mde = distError;
2194 //std::cout << "||mde||=" << myNorm(mde) << std::endl;
2195
2196 if(verbose)
2197 {
2198 std::string name = "errorDistribution_";
2199 name += std::to_string(step);
2200 writeVTKFile(mde.descriptions.gridView,mde.descriptions,mde,name);
2201 }
2202 // Fehler bzgl L_xx : ( sum_i estSol_i * (L_xx * estSol)_i )^{1/2} (nur die x - Komponente, nicht p)
2203 // Lokalisierung: einfachste Idee: v_i = estSol_i * (L_xx * estSol)_i
2204 // Aufteilen von v_i auf die einzelnen Zellen
2205 auto cend = extensionVariableSetDescription.gridView.template end<0>();
2206 for (auto ci=extensionVariableSetDescription.gridView.template begin<0>(); ci!=cend; ++ci)
2207 errorDistribution[is.index(*ci)] = std::make_pair( fabs(boost::fusion::at_c<0>(mde.data).value(*ci,Dune::FieldVector<Scalar,dim>(0.3))) , is.index(*ci));
2208
2209 totalErrorSquared = std::accumulate(errorDistribution.begin(), errorDistribution.end(), 0.0, HierarchicErrorEstimator_Detail::add);
2210
2211 if(verbose) std::cout << "overall error estimation time: " << boost::timer::format(overallTimer.elapsed()) << std::endl;
2212 }
2213
2214 void refineGrid()
2215 {
2216 if(totalErrorSquared > 0) this->refineGrid_impl(errorDistribution,sqrt(totalErrorSquared));
2217// Scalar tol = squaredFraction*totalErrorSquared;
2218// Scalar sum = 0;
2219// std::sort(errorDistribution.begin(), errorDistribution.end(), [](std::pair<double,int> const& a, std::pair<double,int> const& b) -> bool
2220// {
2221// return a.first > b.first;
2222// });
2223// std::vector<size_t> cellIds;
2224// for(size_t i=0; i<errorDistribution.size(); ++i)
2225// {
2226// sum += errorDistribution[i].first;
2227// cellIds.push_back(errorDistribution[i].second);
2228// if(sum > tol) break;
2229// }
2230//
2231// forEachCell(variableSetDescription.gridView, [&](Cell const& cell)
2232// {
2233// if(std::find(cellIds.begin(),cellIds.end(),variableSetDescription.gridView.indexSet().index(cell))!=cellIds.end()) extensionSpace.gridManager().mark(1,cell);
2234// });
2235// extensionSpace.gridManager().adaptAtOnce();
2236 }
2237
2238 double estimatedAbsoluteError() const final
2239 {
2240 return sqrt(fabs(totalErrorSquared));
2241 }
2242
2243 size_t gridSize() const final
2244 {
2245 return extensionSpace.gridManager().grid().size(0);
2246 }
2247
2248 template <typename... Args>
2249 void initFunctionals(const Args&... args)
2250 {
2251 F_HH.reset(new Functional_HH(args...));
2252 F_HE.reset(new Functional_HE(args...));
2253 F_EH.reset(new Functional_EH(args...));
2254 F_EE.reset(new Functional_EE(args...));
2255 }
2256
2257 private:
2258 void initAssemblers()
2259 {
2260 ass_HH.reset(new Ass_HH(extensionVariableSetDescription.spaces));
2261 ass_HE.reset(new Ass_HE(extensionVariableSetDescription.spaces));
2262 ass_EH.reset(new Ass_EH(extensionVariableSetDescription.spaces));
2263 ass_EE.reset(new Ass_EE(extensionVariableSetDescription.spaces));
2264 }
2265
2266 void assembleAt(const Bridge::Vector<typename VariableSetDescription::VariableSet>& xh, const Bridge::Vector<typename ExtensionVariableSetDescription::VariableSet>& xe)
2267 {
2268 ass_HH->assemble(linearization(*F_HH,xh.get()));
2269 ass_HE->assemble(linearization(*F_HE,xh.get()));
2270 ass_EH->assemble(linearization(*F_EH,xe.get()));
2271 ass_EE->assemble(linearization(*F_EE,xe.get()));
2272 }
2273
2274 std::unique_ptr<Functional_HH> F_HH;
2275 std::unique_ptr<Functional_HE> F_HE;
2276 std::unique_ptr<Functional_EH> F_EH;
2277 std::unique_ptr<Functional_EE> F_EE;
2278 std::unique_ptr<Ass_HH> ass_HH;
2279 std::unique_ptr<Ass_HE> ass_HE;
2280 std::unique_ptr<Ass_EH> ass_EH;
2281 std::unique_ptr<Ass_EE> ass_EE;
2282 // LF& lf;
2283 // LVars& lvars;
2284 NormFunctional& normFunctional;
2285 VariableSetDescription& variableSetDescription;
2286 ExtensionVariableSetDescription& extensionVariableSetDescription;
2287 ExtensionSpace& extensionSpace;
2288 Assembler assembler;
2289 Scalar squaredFraction;
2290 Scalar totalErrorSquared;
2291 std::vector<std::pair<double,int> > errorDistribution;
2292 bool verbose, fast;
2293 // AdjustRHS<CoefficientVector> adjustRHS;
2294 };
2295}
2296
2297#endif
Kaskade::Adaptivity::ErrorEquilibration< VariableSetDescription::Grid >::refineGrid_impl
void refineGrid_impl(Err const &err, ErrorRepresentation &errorDistribution, Scalar tol)
Definition: adaptationStrategy.hh:263
Kaskade::AdjointEquationHBErrorEstimator::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:1181
Kaskade::AdjointEquationHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:1095
Kaskade::AdjointEquationHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:1078
Kaskade::AdjointEquationHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:1079
Kaskade::AdjointEquationHBErrorEstimator::AdjointEquationHBErrorEstimator
AdjointEquationHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:1083
Kaskade::AdjointEquationLinearPropagationHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:1226
Kaskade::AdjointEquationLinearPropagationHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:1245
Kaskade::AdjointEquationLinearPropagationHBErrorEstimator::AdjointEquationLinearPropagationHBErrorEstimator
AdjointEquationLinearPropagationHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:1231
Kaskade::AdjointEquationLinearPropagationHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:1227
Kaskade::AnotherHBErrorEstimator::estimatedAbsoluteError
double estimatedAbsoluteError() const final
Definition: errorEstimator.hh:1573
Kaskade::AnotherHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:1412
Kaskade::AnotherHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:1413
Kaskade::AnotherHBErrorEstimator::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:1584
Kaskade::AnotherHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:1431
Kaskade::AnotherHBErrorEstimator::AnotherHBErrorEstimator
AnotherHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:1417
Kaskade::AssembledGalerkinOperator
An adaptor presenting a Dune::LinearOperator <domain_type,range_type> interface to a contiguous sub-b...
Definition: istlinterface.hh:758
Kaskade::AssembledGalerkinOperator::apply
virtual void apply(Domain const &x, Range &b) const
compute
Definition: istlinterface.hh:899
Kaskade::AssembledGalerkinOperator::applyscaleadd
virtual void applyscaleadd(Scalar alpha, Domain const &x, Range &b) const
Compute .
Definition: istlinterface.hh:938
Kaskade::Bridge::Vector
Mathematical Vector that supports copy-on-write, implements AbstractFunctionSpaceElement.
Definition: newton_bridge.hh:196
Kaskade::Bridge::Vector::get
Implementation const & get() const
Access to the data.
Definition: newton_bridge.hh:216
Kaskade::HierarchicErrorEstimator
Defines assembly of hierarchically extended problems for defining DLY style error estimators.
Definition: fem/hierarchicErrorEstimator.hh:153
Kaskade::HierarchicalBasisErrorEstimator2::Muuhh
AssembledGalerkinOperator< Ass_EE, uId, uId+1, uId, uId+1 > Muuhh
Definition: errorEstimator.hh:243
Kaskade::HierarchicalBasisErrorEstimator2::Alh
AssembledGalerkinOperator< Ass_HE, pId, pId+1, yId, yId+1 > Alh
Definition: errorEstimator.hh:226
Kaskade::HierarchicalBasisErrorEstimator2::BTlh
AssembledGalerkinOperator< Ass_HE, uId, uId+1, pId, pId+1 > BTlh
Definition: errorEstimator.hh:236
Kaskade::HierarchicalBasisErrorEstimator2::CVH23
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< 2, 3 > CVH23
Definition: errorEstimator.hh:259
Kaskade::HierarchicalBasisErrorEstimator2::Functional_EE
Functional< ExtensionVariableSetDescription, ExtensionVariableSetDescription, ExtensionVariableSetDescription, lumpM > Functional_EE
Definition: errorEstimator.hh:211
Kaskade::HierarchicalBasisErrorEstimator2::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:270
Kaskade::HierarchicalBasisErrorEstimator2::Muulh
AssembledGalerkinOperator< Ass_HE, uId, uId+1, uId, uId+1 > Muulh
Definition: errorEstimator.hh:241
Kaskade::HierarchicalBasisErrorEstimator2::Functional_HH
Functional< VariableSetDescription, VariableSetDescription, VariableSetDescription, lumpM > Functional_HH
Definition: errorEstimator.hh:208
Kaskade::HierarchicalBasisErrorEstimator2::Functional_EH
Functional< ExtensionVariableSetDescription, VariableSetDescription, ExtensionVariableSetDescription, lumpM > Functional_EH
Definition: errorEstimator.hh:209
Kaskade::HierarchicalBasisErrorEstimator2::Myyhh
AssembledGalerkinOperator< Ass_EE, yId, yId+1, yId, yId+1 > Myyhh
Definition: errorEstimator.hh:228
Kaskade::HierarchicalBasisErrorEstimator2::Muuhl
AssembledGalerkinOperator< Ass_EH, uId, uId+1, uId, uId+1 > Muuhl
Definition: errorEstimator.hh:242
Kaskade::HierarchicalBasisErrorEstimator2::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:286
Kaskade::HierarchicalBasisErrorEstimator2::CVL
VariableSetDescription::template CoefficientVectorRepresentation CVL
Definition: errorEstimator.hh:249
Kaskade::HierarchicalBasisErrorEstimator2::CVLU
VariableSetDescription::template CoefficientVectorRepresentation< uId, uId+1 > CVLU
Definition: errorEstimator.hh:251
Kaskade::HierarchicalBasisErrorEstimator2::HierarchicalBasisErrorEstimator2
HierarchicalBasisErrorEstimator2(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:275
Kaskade::HierarchicalBasisErrorEstimator2::CVHU
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< uId, uId+1 > CVHU
Definition: errorEstimator.hh:263
Kaskade::HierarchicalBasisErrorEstimator2::BThh
AssembledGalerkinOperator< Ass_EE, uId, uId+1, pId, pId+1 > BThh
Definition: errorEstimator.hh:238
Kaskade::HierarchicalBasisErrorEstimator2::Ass_EE
VariationalFunctionalAssembler< LinearizationAt< Functional_EE > > Ass_EE
Definition: errorEstimator.hh:218
Kaskade::HierarchicalBasisErrorEstimator2::CVLP
VariableSetDescription::template CoefficientVectorRepresentation< pId, pId+1 > CVLP
Definition: errorEstimator.hh:252
Kaskade::HierarchicalBasisErrorEstimator2::Ass_HH
VariationalFunctionalAssembler< LinearizationAt< Functional_HH > > Ass_HH
Definition: errorEstimator.hh:215
Kaskade::HierarchicalBasisErrorEstimator2::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:705
Kaskade::HierarchicalBasisErrorEstimator2::VarSet
VariableSet< VariableSetDescription > VarSet
Definition: errorEstimator.hh:247
Kaskade::HierarchicalBasisErrorEstimator2::Blh
AssembledGalerkinOperator< Ass_HE, pId, pId+1, uId, uId+1 > Blh
Definition: errorEstimator.hh:231
Kaskade::HierarchicalBasisErrorEstimator2::CVH
ExtensionVariableSetDescription::template CoefficientVectorRepresentation CVH
Definition: errorEstimator.hh:261
Kaskade::HierarchicalBasisErrorEstimator2::Ass_HE
VariationalFunctionalAssembler< LinearizationAt< Functional_HE > > Ass_HE
Definition: errorEstimator.hh:216
Kaskade::HierarchicalBasisErrorEstimator2::CVHP
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< pId, pId+1 > CVHP
Definition: errorEstimator.hh:264
Kaskade::HierarchicalBasisErrorEstimator2::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:269
Kaskade::HierarchicalBasisErrorEstimator2::Functional_HE
Functional< VariableSetDescription, ExtensionVariableSetDescription, VariableSetDescription, lumpM > Functional_HE
Definition: errorEstimator.hh:210
Kaskade::HierarchicalBasisErrorEstimator2::Assembler
VariationalFunctionalAssembler< ErrorEstimator > Assembler
Definition: errorEstimator.hh:214
Kaskade::HierarchicalBasisErrorEstimator2::CoefficientVectorH01
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< 0, 1 >::type CoefficientVectorH01
Definition: errorEstimator.hh:257
Kaskade::HierarchicalBasisErrorEstimator2::Ass_EH
VariationalFunctionalAssembler< LinearizationAt< Functional_EH > > Ass_EH
Definition: errorEstimator.hh:217
Kaskade::HierarchicalBasisErrorEstimator2::Ahl
AssembledGalerkinOperator< Ass_EH, pId, pId+1, yId, yId+1 > Ahl
Definition: errorEstimator.hh:227
Kaskade::HierarchicalBasisErrorEstimator2::CoefficientVectorH23
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< 2, 3 >::type CoefficientVectorH23
Definition: errorEstimator.hh:260
Kaskade::HierarchicalBasisErrorEstimator2::Bhh
AssembledGalerkinOperator< Ass_EE, pId, pId+1, uId, uId+1 > Bhh
Definition: errorEstimator.hh:233
Kaskade::HierarchicalBasisErrorEstimator2::All
AssembledGalerkinOperator< Ass_HH, pId, pId+1, yId, yId+1 > All
Definition: errorEstimator.hh:225
Kaskade::HierarchicalBasisErrorEstimator2::CVLY
VariableSetDescription::template CoefficientVectorRepresentation< yId, yId+1 > CVLY
Definition: errorEstimator.hh:250
Kaskade::HierarchicalBasisErrorEstimator2::Ahh
AssembledGalerkinOperator< Ass_EE, pId, pId+1, yId, yId+1 > Ahh
Definition: errorEstimator.hh:224
Kaskade::HierarchicalBasisErrorEstimator2::Bll
AssembledGalerkinOperator< Ass_HH, pId, pId+1, uId, uId+1 > Bll
Definition: errorEstimator.hh:230
Kaskade::HierarchicalBasisErrorEstimator2::Bhl
AssembledGalerkinOperator< Ass_EH, pId, pId+1, uId, uId+1 > Bhl
Definition: errorEstimator.hh:232
Kaskade::HierarchicalBasisErrorEstimator2::ErrorEstimator
HierarchicErrorEstimator< LinearizationAt< Functional_HH >, ExtensionVariableSetDescription, ExtensionVariableSetDescription > ErrorEstimator
Definition: errorEstimator.hh:212
Kaskade::HierarchicalBasisErrorEstimator2::ExtVarSet
VariableSet< ExtensionVariableSetDescription > ExtVarSet
Definition: errorEstimator.hh:248
Kaskade::HierarchicalBasisErrorEstimator2::BThl
AssembledGalerkinOperator< Ass_EH, uId, uId+1, pId, pId+1 > BThl
Definition: errorEstimator.hh:237
Kaskade::HierarchicalBasisErrorEstimator2::CVHY
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< yId, yId+1 > CVHY
Definition: errorEstimator.hh:262
Kaskade::HierarchicalBasisErrorEstimator2::Muull
AssembledGalerkinOperator< Ass_HH, uId, uId+1, uId, uId+1 > Muull
Definition: errorEstimator.hh:240
Kaskade::HierarchicalBasisErrorEstimator2::CoefficientVectorH02
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< 0, 2 >::type CoefficientVectorH02
Definition: errorEstimator.hh:258
Kaskade::HierarchicalBasisErrorEstimator2::BTll
AssembledGalerkinOperator< Ass_HH, uId, uId+1, pId, pId+1 > BTll
Definition: errorEstimator.hh:235
Kaskade::HierarchicalBasisErrorEstimator3::Ass_HH
VariationalFunctionalAssembler< LinearizationAt< Functional_HH > > Ass_HH
Definition: errorEstimator.hh:2055
Kaskade::HierarchicalBasisErrorEstimator3::Functional_EE
Functional< ExtensionVariableSetDescription, ExtensionVariableSetDescription, ExtensionVariableSetDescription, lumpM > Functional_EE
Definition: errorEstimator.hh:2051
Kaskade::HierarchicalBasisErrorEstimator3::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:2249
Kaskade::HierarchicalBasisErrorEstimator3::Assembler
VariationalFunctionalAssembler< ErrorEstimator > Assembler
Definition: errorEstimator.hh:2054
Kaskade::HierarchicalBasisErrorEstimator3::Functional_EH
Functional< ExtensionVariableSetDescription, VariableSetDescription, ExtensionVariableSetDescription, lumpM > Functional_EH
Definition: errorEstimator.hh:2049
Kaskade::HierarchicalBasisErrorEstimator3::ErrorEstimator
HierarchicErrorEstimator< LinearizationAt< Functional_HH >, ExtensionVariableSetDescription, ExtensionVariableSetDescription > ErrorEstimator
Definition: errorEstimator.hh:2052
Kaskade::HierarchicalBasisErrorEstimator3::Ass_EH
VariationalFunctionalAssembler< LinearizationAt< Functional_EH > > Ass_EH
Definition: errorEstimator.hh:2057
Kaskade::HierarchicalBasisErrorEstimator3::HierarchicalBasisErrorEstimator3
HierarchicalBasisErrorEstimator3(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:2080
Kaskade::HierarchicalBasisErrorEstimator3::Functional_HE
Functional< VariableSetDescription, ExtensionVariableSetDescription, VariableSetDescription, lumpM > Functional_HE
Definition: errorEstimator.hh:2050
Kaskade::HierarchicalBasisErrorEstimator3::Ass_HE
VariationalFunctionalAssembler< LinearizationAt< Functional_HE > > Ass_HE
Definition: errorEstimator.hh:2056
Kaskade::HierarchicalBasisErrorEstimator3::Functional_HH
Functional< VariableSetDescription, VariableSetDescription, VariableSetDescription, lumpM > Functional_HH
Definition: errorEstimator.hh:2048
Kaskade::HierarchicalBasisErrorEstimator3::Ass_EE
VariationalFunctionalAssembler< LinearizationAt< Functional_EE > > Ass_EE
Definition: errorEstimator.hh:2058
Kaskade::HierarchicalBasisErrorEstimator3::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:2091
Kaskade::InverseLinearOperator
Dune::LinearOperator interface for inverse operators.
Definition: direct.hh:553
Kaskade::JacobiPreconditionerDetail::JacobiPreconditionerBase::apply
virtual void apply(domain_type &x, range_type const &b)
Definition: jacobiPreconditioner.hh:187
Kaskade::LinearizationAt
Proxy class for the linearization of a nonlinear functional.
Definition: functional_aux.hh:615
Kaskade::MartinsErrorEstimator::Functional_HE
Traits::Functional_HE Functional_HE
Definition: errorEstimator.hh:1782
Kaskade::MartinsErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:1819
Kaskade::MartinsErrorEstimator::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:1996
Kaskade::MartinsErrorEstimator::Ass_HH
FullTraits::Assembler_HH Ass_HH
Definition: errorEstimator.hh:1790
Kaskade::MartinsErrorEstimator::Functional_HH
FullTraits::Functional_HH Functional_HH
Definition: errorEstimator.hh:1781
Kaskade::MartinsErrorEstimator::MartinsErrorEstimator
MartinsErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:1824
Kaskade::MartinsErrorEstimator::Functional_EH
Traits::Functional_EH Functional_EH
Definition: errorEstimator.hh:1783
Kaskade::MartinsErrorEstimator::Functional_EE
Traits::Functional_EE Functional_EE
Definition: errorEstimator.hh:1784
Kaskade::MartinsErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:1838
Kaskade::MartinsErrorEstimator::H_HH
FullTraits::Operator_HH H_HH
Definition: errorEstimator.hh:1793
Kaskade::MartinsErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:1820
Kaskade::MartinsErrorEstimator::estimatedAbsoluteError
double estimatedAbsoluteError() const final
Definition: errorEstimator.hh:1985
Kaskade::MatrixAsTriplet::transpose
MatrixAsTriplet & transpose()
Transposition.
Definition: triplet.hh:453
Kaskade::MatrixAsTriplet::usmtv
void usmtv(Scalar const alpha, X const &x, Y &y) const
Matrix-vector multiplication.
Definition: triplet.hh:296
Kaskade::MatrixAsTriplet::usmv
void usmv(Scalar const alpha, X const &x, Y &y) const
scaled matrix-vector multiplication
Definition: triplet.hh:304
Kaskade::MatrixRepresentedOperator
Operator with matrix-representation and coordinate mappings.
Definition: istlinterface.hh:310
Kaskade::StateEquationHBErrorEstimator::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:880
Kaskade::StateEquationHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:778
Kaskade::StateEquationHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:794
Kaskade::StateEquationHBErrorEstimator::StateEquationHBErrorEstimator
StateEquationHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:782
Kaskade::StateEquationHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:777
Kaskade::StupidHBErrorEstimator::estimatedAbsoluteError
double estimatedAbsoluteError() const final
Definition: errorEstimator.hh:1731
Kaskade::StupidHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:1644
Kaskade::StupidHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:1660
Kaskade::StupidHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:1643
Kaskade::StupidHBErrorEstimator::Functional_HH
Traits::Functional_HH Functional_HH
Definition: errorEstimator.hh:1628
Kaskade::StupidHBErrorEstimator::initFunctionals
void initFunctionals(const Args &... args)
Definition: errorEstimator.hh:1742
Kaskade::StupidHBErrorEstimator::Functional_HE
Traits::Functional_HE Functional_HE
Definition: errorEstimator.hh:1629
Kaskade::StupidHBErrorEstimator::StupidHBErrorEstimator
StupidHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:1648
Kaskade::StupidHBErrorEstimator::Functional_EE
Traits::Functional_EE Functional_EE
Definition: errorEstimator.hh:1630
Kaskade::VariableSetDescription
A class that stores information about a set of variables.
Definition: variables.hh:537
Kaskade::VariableSetDescription::gridView
GridView const & gridView
The grid view on which the variables live.
Definition: variables.hh:812
Kaskade::VariableSetDescription::spaces
Spaces spaces
A heterogeneous sequence of pointers to (const) spaces involved.
Definition: variables.hh:807
Kaskade::VariableSet
A class for storing a heterogeneous collection of FunctionSpaceElement s.
Definition: variables.hh:341
Kaskade::VariationalEquationHBErrorEstimator::Cell
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Entity Cell
Definition: errorEstimator.hh:932
Kaskade::VariationalEquationHBErrorEstimator::operator()
void operator()(AbstractVector const &x_, AbstractVector const &dx_, int step, AbstractVector const &)
Definition: errorEstimator.hh:946
Kaskade::VariationalEquationHBErrorEstimator::VariationalEquationHBErrorEstimator
VariationalEquationHBErrorEstimator(NormFunctional &normFunctional_, VariableSetDescription &variableSetDescription_, ExtensionVariableSetDescription &extensionVariableSetDescription_, ExtensionSpace &extensionSpace_, Scalar fraction=0.7, bool verbose_=false, bool fast_=false)
Definition: errorEstimator.hh:934
Kaskade::VariationalEquationHBErrorEstimator::CellIterator
ExtensionVariableSetDescription::GridView::template Codim< 0 >::Iterator CellIterator
Definition: errorEstimator.hh:931
Kaskade::VariationalFunctionalAssembler
A class for assembling Galerkin representation matrices and right hand sides for variational function...
Definition: assemblerCore.hh:1432
Kaskade::VariationalFunctionalAssembler::assemble
void assemble(F const &f, unsigned int flags=Assembler::EVERYTHING, int nThreads=0, bool verbose=false)
Assembly without block filter or cell filter.
Definition: assemblerCore.hh:1598
hierarchicErrorEstimator.hh
Error estimation via hierachic FEM.
Kaskade::Assembler::RHS
@ RHS
Definition: assemblerCore.hh:1402
DirectType::UMFPACK3264
@ UMFPACK3264
NO LONGER SUPPORTED.
MatrixProperties::GENERAL
@ GENERAL
Dune::max
Dune::FieldVector< T, n > max(Dune::FieldVector< T, n > x, Dune::FieldVector< T, n > const &y)
Componentwise maximum.
Definition: fixdune.hh:110
Dune::min
Dune::FieldVector< T, n > min(Dune::FieldVector< T, n > x, Dune::FieldVector< T, n > const &y)
Componentwise minimum.
Definition: fixdune.hh:122
Kaskade::ErrorEstimator_Detail::getStateId
constexpr int getStateId()
Definition: errorEstimator.hh:73
Kaskade::ErrorEstimator_Detail::getControlId
constexpr int getControlId()
Definition: errorEstimator.hh:80
Kaskade::ErrorEstimator_Detail::getControlIdImpl
std::integral_constant< int, Functional::controlId > getControlIdImpl(decltype(Functional::controlId) *)
Kaskade::ErrorEstimator_Detail::getAdjointId
constexpr int getAdjointId()
Definition: errorEstimator.hh:87
Kaskade::ErrorEstimator_Detail::getStateIdImpl
std::integral_constant< int, Functional::stateId > getStateIdImpl(decltype(Functional::stateId) *)
Kaskade::ErrorEstimator_Detail::getAdjointIdImpl
std::integral_constant< int, Functional::adjointId > getAdjointIdImpl(decltype(Functional::adjointId) *)
Kaskade::HierarchicErrorEstimator_Detail::add
double add(const double &x1, const std::pair< double, int > &x2)
Definition: algorithm/hierarchicErrorEstimator.hh:47
newton_bridge.hh
Bridge classes that connect low level FE algorithms to higher level algorithms.
Kaskade::ErrorEstimator_Detail::Traits
Definition: errorEstimator.hh:95
Kaskade::ErrorEstimator_Detail::Traits::Luy_EE
AssembledGalerkinOperator< Assembler_EE, controlId, controlId+1, stateId, stateId+1 > Luy_EE
Definition: errorEstimator.hh:144
Kaskade::ErrorEstimator_Detail::Traits::Luy_HE
AssembledGalerkinOperator< Assembler_HE, controlId, controlId+1, stateId, stateId+1 > Luy_HE
Definition: errorEstimator.hh:143
Kaskade::ErrorEstimator_Detail::Traits::A_EE
AssembledGalerkinOperator< Assembler_EE, adjointId, adjointId+1, stateId, stateId+1 > A_EE
Definition: errorEstimator.hh:159
Kaskade::ErrorEstimator_Detail::Traits::VectorY_Initializer
VariableSetDescription::template CoefficientVectorRepresentation< stateId, stateId+1 > VectorY_Initializer
Definition: errorEstimator.hh:187
Kaskade::ErrorEstimator_Detail::Traits::Assembler_HH
VariationalFunctionalAssembler< LinearizationAt< Functional_HH > > Assembler_HH
Definition: errorEstimator.hh:106
Kaskade::ErrorEstimator_Detail::Traits::Lyy_HE
AssembledGalerkinOperator< Assembler_HE, stateId, stateId+1, stateId, stateId+1 > Lyy_HE
Definition: errorEstimator.hh:133
Kaskade::ErrorEstimator_Detail::Traits::Luy_EH
AssembledGalerkinOperator< Assembler_EH, controlId, controlId+1, stateId, stateId+1 > Luy_EH
Definition: errorEstimator.hh:142
Kaskade::ErrorEstimator_Detail::Traits::Lyu_HE
AssembledGalerkinOperator< Assembler_HE, stateId, stateId+1, controlId, controlId+1 > Lyu_HE
Definition: errorEstimator.hh:138
Kaskade::ErrorEstimator_Detail::Traits::Functional_HE
Functional< VariableSetDescription, ExtensionVariableSetDescription, VariableSetDescription > Functional_HE
Definition: errorEstimator.hh:100
Kaskade::ErrorEstimator_Detail::Traits::AT_EE
AssembledGalerkinOperator< Assembler_EE, stateId, stateId+1, adjointId, adjointId+1 > AT_EE
Definition: errorEstimator.hh:169
Kaskade::ErrorEstimator_Detail::Traits::BT_EE
AssembledGalerkinOperator< Assembler_EE, controlId, controlId+1, adjointId, adjointId+1 > BT_EE
Definition: errorEstimator.hh:164
Kaskade::ErrorEstimator_Detail::Traits::adjointId
static constexpr int adjointId
Definition: errorEstimator.hh:113
Kaskade::ErrorEstimator_Detail::Traits::Lyy_HH
AssembledGalerkinOperator< Assembler_HH, stateId, stateId+1, stateId, stateId+1 > Lyy_HH
Definition: errorEstimator.hh:131
Kaskade::ErrorEstimator_Detail::Traits::Operator_EE
AssembledGalerkinOperator< Assembler_EE > Operator_EE
Definition: errorEstimator.hh:129
Kaskade::ErrorEstimator_Detail::Traits::Functional_HH
Functional< VariableSetDescription, VariableSetDescription, VariableSetDescription > Functional_HH
Definition: errorEstimator.hh:99
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorY_Initializer
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< stateId, stateId+1 > ExtensionVectorY_Initializer
Definition: errorEstimator.hh:191
Kaskade::ErrorEstimator_Detail::Traits::Operator_HH
AssembledGalerkinOperator< Assembler_HH > Operator_HH
Definition: errorEstimator.hh:126
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorP
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< adjointId, adjointId+1 >::type ExtensionVectorP
Definition: errorEstimator.hh:181
Kaskade::ErrorEstimator_Detail::Traits::VectorP
VariableSetDescription::template CoefficientVectorRepresentation< adjointId, adjointId+1 >::type VectorP
Definition: errorEstimator.hh:177
Kaskade::ErrorEstimator_Detail::Traits::Lyu_EH
AssembledGalerkinOperator< Assembler_EH, stateId, stateId+1, controlId, controlId+1 > Lyu_EH
Definition: errorEstimator.hh:137
Kaskade::ErrorEstimator_Detail::Traits::AT_EH
AssembledGalerkinOperator< Assembler_EH, stateId, stateId+1, adjointId, adjointId+1 > AT_EH
Definition: errorEstimator.hh:167
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorU_Initializer
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< controlId, controlId+1 > ExtensionVectorU_Initializer
Definition: errorEstimator.hh:192
Kaskade::ErrorEstimator_Detail::Traits::VectorU
VariableSetDescription::template CoefficientVectorRepresentation< controlId, controlId+1 >::type VectorU
Definition: errorEstimator.hh:176
Kaskade::ErrorEstimator_Detail::Traits::Luu_HE
AssembledGalerkinOperator< Assembler_HE, controlId, controlId+1, controlId, controlId+1 > Luu_HE
Definition: errorEstimator.hh:148
Kaskade::ErrorEstimator_Detail::Traits::Scalar
Functional_HH::Scalar Scalar
Definition: errorEstimator.hh:104
Kaskade::ErrorEstimator_Detail::Traits::Luu_EE
AssembledGalerkinOperator< Assembler_EE, controlId, controlId+1, controlId, controlId+1 > Luu_EE
Definition: errorEstimator.hh:149
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorU
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< controlId, controlId+1 >::type ExtensionVectorU
Definition: errorEstimator.hh:180
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVector_Initializer
ExtensionVariableSetDescription::template CoefficientVectorRepresentation ExtensionVector_Initializer
Definition: errorEstimator.hh:190
Kaskade::ErrorEstimator_Detail::Traits::BT_EH
AssembledGalerkinOperator< Assembler_EH, controlId, controlId+1, adjointId, adjointId+1 > BT_EH
Definition: errorEstimator.hh:162
Kaskade::ErrorEstimator_Detail::Traits::AT_HH
AssembledGalerkinOperator< Assembler_HH, stateId, stateId+1, adjointId, adjointId+1 > AT_HH
Definition: errorEstimator.hh:166
Kaskade::ErrorEstimator_Detail::Traits::Lyy_EE
AssembledGalerkinOperator< Assembler_EE, stateId, stateId+1, stateId, stateId+1 > Lyy_EE
Definition: errorEstimator.hh:134
Kaskade::ErrorEstimator_Detail::Traits::A_HE
AssembledGalerkinOperator< Assembler_HE, adjointId, adjointId+1, stateId, stateId+1 > A_HE
Definition: errorEstimator.hh:158
Kaskade::ErrorEstimator_Detail::Traits::Lyy_EH
AssembledGalerkinOperator< Assembler_EH, stateId, stateId+1, stateId, stateId+1 > Lyy_EH
Definition: errorEstimator.hh:132
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorY
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< stateId, stateId+1 >::type ExtensionVectorY
Definition: errorEstimator.hh:179
Kaskade::ErrorEstimator_Detail::Traits::VarSet
VariableSet< VariableSetDescription > VarSet
Definition: errorEstimator.hh:96
Kaskade::ErrorEstimator_Detail::Traits::VectorY
VariableSetDescription::template CoefficientVectorRepresentation< stateId, stateId+1 >::type VectorY
Definition: errorEstimator.hh:175
Kaskade::ErrorEstimator_Detail::Traits::B_HE
AssembledGalerkinOperator< Assembler_HE, adjointId, adjointId+1, controlId, controlId+1 > B_HE
Definition: errorEstimator.hh:153
Kaskade::ErrorEstimator_Detail::Traits::B_EH
AssembledGalerkinOperator< Assembler_EH, adjointId, adjointId+1, controlId, controlId+1 > B_EH
Definition: errorEstimator.hh:152
Kaskade::ErrorEstimator_Detail::Traits::VectorP_Initializer
VariableSetDescription::template CoefficientVectorRepresentation< adjointId, adjointId+1 > VectorP_Initializer
Definition: errorEstimator.hh:189
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVectorP_Initializer
ExtensionVariableSetDescription::template CoefficientVectorRepresentation< adjointId, adjointId+1 > ExtensionVectorP_Initializer
Definition: errorEstimator.hh:193
Kaskade::ErrorEstimator_Detail::Traits::A_EH
AssembledGalerkinOperator< Assembler_EH, adjointId, adjointId+1, stateId, stateId+1 > A_EH
Definition: errorEstimator.hh:157
Kaskade::ErrorEstimator_Detail::Traits::A_HH
AssembledGalerkinOperator< Assembler_HH, adjointId, adjointId+1, stateId, stateId+1 > A_HH
Definition: errorEstimator.hh:156
Kaskade::ErrorEstimator_Detail::Traits::B_EE
AssembledGalerkinOperator< Assembler_EE, adjointId, adjointId+1, controlId, controlId+1 > B_EE
Definition: errorEstimator.hh:154
Kaskade::ErrorEstimator_Detail::Traits::Assembler_HE
VariationalFunctionalAssembler< LinearizationAt< Functional_HE > > Assembler_HE
Definition: errorEstimator.hh:107
Kaskade::ErrorEstimator_Detail::Traits::BT_HH
AssembledGalerkinOperator< Assembler_HH, controlId, controlId+1, adjointId, adjointId+1 > BT_HH
Definition: errorEstimator.hh:161
Kaskade::ErrorEstimator_Detail::Traits::Luy_HH
AssembledGalerkinOperator< Assembler_HH, controlId, controlId+1, stateId, stateId+1 > Luy_HH
Definition: errorEstimator.hh:141
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVector
ExtensionVariableSetDescription::template CoefficientVectorRepresentation ::type ExtensionVector
Definition: errorEstimator.hh:178
Kaskade::ErrorEstimator_Detail::Traits::Assembler_EH
VariationalFunctionalAssembler< LinearizationAt< Functional_EH > > Assembler_EH
Definition: errorEstimator.hh:108
Kaskade::ErrorEstimator_Detail::Traits::Vector
VariableSetDescription::template CoefficientVectorRepresentation ::type Vector
Definition: errorEstimator.hh:174
Kaskade::ErrorEstimator_Detail::Traits::B_HH
AssembledGalerkinOperator< Assembler_HH, adjointId, adjointId+1, controlId, controlId+1 > B_HH
Definition: errorEstimator.hh:151
Kaskade::ErrorEstimator_Detail::Traits::Functional_EH
Functional< ExtensionVariableSetDescription, VariableSetDescription, ExtensionVariableSetDescription > Functional_EH
Definition: errorEstimator.hh:101
Kaskade::ErrorEstimator_Detail::Traits::dim
static constexpr int dim
Definition: errorEstimator.hh:115
Kaskade::ErrorEstimator_Detail::Traits::ExtensionVarSet
VariableSet< ExtensionVariableSetDescription > ExtensionVarSet
Definition: errorEstimator.hh:97
Kaskade::ErrorEstimator_Detail::Traits::VectorU_Initializer
VariableSetDescription::template CoefficientVectorRepresentation< controlId, controlId+1 > VectorU_Initializer
Definition: errorEstimator.hh:188
Kaskade::ErrorEstimator_Detail::Traits::Operator_EH
AssembledGalerkinOperator< Assembler_EH > Operator_EH
Definition: errorEstimator.hh:127
Kaskade::ErrorEstimator_Detail::Traits::Operator_HE
AssembledGalerkinOperator< Assembler_HE > Operator_HE
Definition: errorEstimator.hh:128
Kaskade::ErrorEstimator_Detail::Traits::Lyu_EE
AssembledGalerkinOperator< Assembler_EE, stateId, stateId+1, controlId, controlId+1 > Lyu_EE
Definition: errorEstimator.hh:139
Kaskade::ErrorEstimator_Detail::Traits::Vector_Initializer
VariableSetDescription::template CoefficientVectorRepresentation Vector_Initializer
Definition: errorEstimator.hh:186
Kaskade::ErrorEstimator_Detail::Traits::Luu_HH
AssembledGalerkinOperator< Assembler_HH, controlId, controlId+1, controlId, controlId+1 > Luu_HH
Definition: errorEstimator.hh:146
Kaskade::ErrorEstimator_Detail::Traits::AT_HE
AssembledGalerkinOperator< Assembler_HE, stateId, stateId+1, adjointId, adjointId+1 > AT_HE
Definition: errorEstimator.hh:168
Kaskade::ErrorEstimator_Detail::Traits::BT_HE
AssembledGalerkinOperator< Assembler_HE, controlId, controlId+1, adjointId, adjointId+1 > BT_HE
Definition: errorEstimator.hh:163
Kaskade::ErrorEstimator_Detail::Traits::Functional_EE
Functional< ExtensionVariableSetDescription, ExtensionVariableSetDescription, ExtensionVariableSetDescription > Functional_EE
Definition: errorEstimator.hh:102
Kaskade::ErrorEstimator_Detail::Traits::Lyu_HH
AssembledGalerkinOperator< Assembler_HH, stateId, stateId+1, controlId, controlId+1 > Lyu_HH
Definition: errorEstimator.hh:136
Kaskade::ErrorEstimator_Detail::Traits::controlId
static constexpr int controlId
Definition: errorEstimator.hh:112
Kaskade::ErrorEstimator_Detail::Traits::stateId
static constexpr int stateId
Definition: errorEstimator.hh:111
Kaskade::ErrorEstimator_Detail::Traits::Assembler_EE
VariationalFunctionalAssembler< LinearizationAt< Functional_EE > > Assembler_EE
Definition: errorEstimator.hh:109
Kaskade::ErrorEstimator_Detail::Traits::Luu_EH
AssembledGalerkinOperator< Assembler_EH, controlId, controlId+1, controlId, controlId+1 > Luu_EH
Definition: errorEstimator.hh:147
variables.hh
Variables and their descriptions.
Generated on Fri Sep 6 2024 13:12:58 for KASKADE 7 development version by doxygen 1.9.4