Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
constraint_interface.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2010-2015: G-CSC, Goethe University Frankfurt
3 * Author: Andreas Vogel
4 *
5 * This file is part of UG4.
6 *
7 * UG4 is free software: you can redistribute it and/or modify it under the
8 * terms of the GNU Lesser General Public License version 3 (as published by the
9 * Free Software Foundation) with the following additional attribution
10 * requirements (according to LGPL/GPL v3 §7):
11 *
12 * (1) The following notice must be displayed in the Appropriate Legal Notices
13 * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 *
15 * (2) The following notice must be displayed at a prominent place in the
16 * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 *
18 * (3) The following bibliography is recommended for citation and must be
19 * preserved in all covered files:
20 * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 * parallel geometric multigrid solver on hierarchically distributed grids.
22 * Computing and visualization in science 16, 4 (2013), 151-164"
23 * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 * flexible software system for simulating pde based models on high performance
25 * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 *
27 * This program is distributed in the hope that it will be useful,
28 * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 * GNU Lesser General Public License for more details.
31 */
32
33#ifndef __H__UG__LIB_DISC__SPATIAL_DISC__CONSTRAINTS__CONSTRAINTS__INTERFACE__
34#define __H__UG__LIB_DISC__SPATIAL_DISC__CONSTRAINTS__CONSTRAINTS__INTERFACE__
35
36// extern headers
37#include <vector>
38
39// intern headers
44
45namespace ug{
46
48
61template <typename TAlgebra>
63{
64 public:
66 typedef TAlgebra algebra_type;
67
69 typedef typename algebra_type::matrix_type matrix_type;
70
72 typedef typename algebra_type::vector_type vector_type;
73
74 public:
76 virtual void adjust_jacobian(matrix_type& J, const vector_type& u,
77 ConstSmartPtr<DoFDistribution> dd, int type, number time = 0.0,
79 const number s_a0 = 1.0) = 0;
80
82 virtual void adjust_defect(vector_type& d, const vector_type& u,
83 ConstSmartPtr<DoFDistribution> dd, int type, number time = 0.0,
85 const std::vector<number>* vScaleMass = NULL,
86 const std::vector<number>* vScaleStiff = NULL) = 0;
87
89 virtual void adjust_linear(matrix_type& mat, vector_type& rhs,
90 ConstSmartPtr<DoFDistribution> dd, int type, number time = 0.0) = 0;
91
93 virtual void adjust_rhs(vector_type& rhs, const vector_type& u,
94 ConstSmartPtr<DoFDistribution> dd, int type, number time = 0.0) = 0;
95
98 number time = 0.0) = 0;
99
104 number time = 0.0) {adjust_solution(c, dd, type, time);}
105
108 ConstSmartPtr<DoFDistribution> dd, int type, number time = 0.0)
109 {adjust_correction(d, dd, type, time);}
110
113 number time = 0.0,
115 const std::vector<number>* vScaleMass = NULL,
116 const std::vector<number>* vScaleStiff = NULL) {};
117
122 int type,
123 number time = 0.0) {};
124
129 int type,
130 number time = 0.0) {};
131
133 virtual void adjust_restriction(vector_type& uCoarse, GridLevel coarseLvl,
134 const vector_type& uFine, GridLevel fineLvl,
135 int type) {};
136
138 virtual void adjust_prolongation(vector_type& uFine, GridLevel fineLvl,
139 const vector_type& uCoarse, GridLevel coarseLvl,
140 int type) {};
141
143 virtual void modify_solution(vector_type& uMod, const vector_type& u,
145
151
153 virtual int type() const = 0;
154
156 virtual ~IConstraint() {};
157
158};
159
160template <typename TDomain, typename TAlgebra>
161class IDomainConstraint : public IConstraint<TAlgebra>
162{
163 public:
165 typedef TDomain domain_type;
166
168 typedef TAlgebra algebra_type;
169
171 typedef typename algebra_type::matrix_type matrix_type;
172
174 typedef typename algebra_type::vector_type vector_type;
175
176 public:
180
183 {
184 m_spApproxSpace = approxSpace;
185 }
186
192
198
200 virtual int type() const = 0;
201
204 {
205 m_spAssTuner = spAssemblingTuner;
206 }
207
208 // //////////////////////////
209 // Error estimator
210 // //////////////////////////
211 public:
213
222
224 bool err_est_enabled() const {return m_bDoErrEst;}
225
228
229 private:
232
233 protected:
236
237 // ///////////////////////////
238
239
240 protected:
243 {
244 return m_spApproxSpace->dof_distribution(gl);
245 }
246
247 protected:
250
253
254};
255
256} // end namespace ug
257
258#endif /* __H__UG__LIB_DISC__SPATIAL_DISC__CONSTRAINTS__CONSTRAINTS__INTERFACE__ */
Definition smart_pointer.h:296
Definition smart_pointer.h:108
base class for approximation spaces without type of algebra or dof distribution
Definition approximation_space.h:279
The AssemblingTuner class combines tools to adapt the assembling routine.
Definition ass_tuner.h:90
Definition grid_level.h:42
interface for adjustment of constraints
Definition constraint_interface.h:63
virtual void adjust_restriction(matrix_type &R, ConstSmartPtr< DoFDistribution > ddCoarse, ConstSmartPtr< DoFDistribution > ddFine, int type, number time=0.0)
sets constraints in restriction
Definition constraint_interface.h:126
virtual void adjust_correction(vector_type &c, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0)
Definition constraint_interface.h:103
virtual void modify_solution(vector_type &uMod, const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type)
modifies solution vector before calling the assembling routine
Definition constraint_interface.h:143
virtual void adjust_linear(matrix_type &mat, vector_type &rhs, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0)=0
adapts matrix and rhs (linear case) to enforce constraints
virtual void adjust_defect(vector_type &d, const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol=SPNULL, const std::vector< number > *vScaleMass=NULL, const std::vector< number > *vScaleStiff=NULL)=0
adapts defect to enforce constraints
algebra_type::vector_type vector_type
Type of algebra vector.
Definition constraint_interface.h:72
virtual void modify_solution(SmartPtr< VectorTimeSeries< vector_type > > vSolMod, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, ConstSmartPtr< DoFDistribution > dd, int type)
modify_solution for instationary case
Definition constraint_interface.h:147
virtual void adjust_solution(vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0)=0
sets the constraints in a solution vector
virtual ~IConstraint()
virtual destructor
Definition constraint_interface.h:156
virtual int type() const =0
returns the type of constraints
virtual void adjust_rhs(vector_type &rhs, const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0)=0
adapts a rhs to enforce constraints
virtual void adjust_restriction(vector_type &uCoarse, GridLevel coarseLvl, const vector_type &uFine, GridLevel fineLvl, int type)
sets the constraints in a solution vector
Definition constraint_interface.h:133
virtual void adjust_prolongation(matrix_type &P, ConstSmartPtr< DoFDistribution > ddFine, ConstSmartPtr< DoFDistribution > ddCoarse, int type, number time=0.0)
sets constraints in prolongation
Definition constraint_interface.h:119
algebra_type::matrix_type matrix_type
Type of algebra matrix.
Definition constraint_interface.h:69
virtual void adjust_error(const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol=SPNULL, const std::vector< number > *vScaleMass=NULL, const std::vector< number > *vScaleStiff=NULL)
adjusts the assembled error estimator values in the attachments according to the constraint
Definition constraint_interface.h:112
virtual void adjust_jacobian(matrix_type &J, const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol=SPNULL, const number s_a0=1.0)=0
adapts jacobian to enforce constraints
TAlgebra algebra_type
Algebra type.
Definition constraint_interface.h:66
virtual void adjust_prolongation(vector_type &uFine, GridLevel fineLvl, const vector_type &uCoarse, GridLevel coarseLvl, int type)
sets the constraints in a solution vector
Definition constraint_interface.h:138
virtual void adjust_linear_residual(vector_type &d, const vector_type &u, ConstSmartPtr< DoFDistribution > dd, int type, number time=0.0)
adjust linear residual
Definition constraint_interface.h:107
Definition constraint_interface.h:162
SmartPtr< ApproximationSpace< TDomain > > m_spApproxSpace
Approximation Space.
Definition constraint_interface.h:249
algebra_type::vector_type vector_type
Type of algebra vector.
Definition constraint_interface.h:174
SmartPtr< IErrEstData< TDomain > > m_spErrEstData
error estimation object associated to the element discretization
Definition constraint_interface.h:235
algebra_type::matrix_type matrix_type
Type of algebra matrix.
Definition constraint_interface.h:171
bool err_est_enabled() const
find out whether or not a posteriori error estimation is to be performed for this disc
Definition constraint_interface.h:224
virtual int type() const =0
returns the type of constraints
virtual SmartPtr< IErrEstData< TDomain > > err_est_data()
returns the pointer to the error estimator data object (or NULL)
Definition constraint_interface.h:227
IDomainConstraint()
constructor
Definition constraint_interface.h:178
void set_error_estimator(SmartPtr< IErrEstData< TDomain > > ee)
sets the pointer to an error estimator data object (or NULL)
Definition constraint_interface.h:221
SmartPtr< ApproximationSpace< TDomain > > approximation_space()
returns approximation space
Definition constraint_interface.h:188
ConstSmartPtr< ApproximationSpace< TDomain > > approximation_space() const
returns approximation space
Definition constraint_interface.h:194
ConstSmartPtr< DoFDistribution > dd(const GridLevel &gl) const
returns the level dof distribution
Definition constraint_interface.h:242
ConstSmartPtr< AssemblingTuner< TAlgebra > > m_spAssTuner
Assemble adapter.
Definition constraint_interface.h:252
void set_ass_tuner(ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssemblingTuner=NULL)
sets the assemble adapter for the constraints
Definition constraint_interface.h:203
virtual void set_approximation_space(SmartPtr< ApproximationSpace< TDomain > > approxSpace)
sets the approximation space
Definition constraint_interface.h:182
TDomain domain_type
Domain Type.
Definition constraint_interface.h:165
bool m_bDoErrEst
flag indicating whether or not a posteriori error estimation is to be performed for this disc
Definition constraint_interface.h:231
TAlgebra algebra_type
Algebra type.
Definition constraint_interface.h:168
Base class for error estimator data.
Definition err_est_data.h:71
time series of solutions and corresponding time point
Definition solution_time_series.h:59
const NullSmartPtr SPNULL
The equivalent to NULL for smart pointers.
Definition smart_pointer.h:90
double number
Definition types.h:124
the ug namespace