DUNE PDELab (git)

dunefunctionslfsindexcache.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_PDELAB_GRIDFUNCTIONSPACE_DUNEFUNCTIONSLFSINDEXCACHE_HH
4#define DUNE_PDELAB_GRIDFUNCTIONSPACE_DUNEFUNCTIONSLFSINDEXCACHE_HH
5
6#include <dune/pdelab/gridfunctionspace/lfsindexcache.hh>
7#include <dune/pdelab/gridfunctionspace/dunefunctionsgridfunctionspace.hh>
8
9namespace Dune {
10 namespace PDELab {
11
12 template<typename LFS, typename C>
13 class LFSIndexCacheBase<LFS,C,Experimental::DuneFunctionsCacheTag, false>
14 {
15
16 enum DOFFlags
17 {
18 DOF_NONCONSTRAINED = 0,
19 DOF_CONSTRAINED = 1<<0,
20 DOF_DIRICHLET = 1<<1
21 };
22
23 public:
24
25 using LocalFunctionSpace = LFS;
26 using GFS = typename LFS::Traits::GridFunctionSpace;
27 using Ordering = typename GFS::Ordering;
28 using DOFIndex = typename Ordering::Traits::DOFIndex;
29 using DI = DOFIndex;
30 using ContainerIndex = typename Ordering::Traits::ContainerIndex;
31 using CI = ContainerIndex;
32 using size_type = std::size_t;
33
34 typedef std::vector<CI> CIVector;
35 typedef std::unordered_map<DI,CI> CIMap;
36
37 typedef std::unordered_map<const CI*,std::pair<size_type,bool> > InverseMap;
38
39 struct ConstraintsEntry
40 : public std::pair<const CI*,typename C::mapped_type::mapped_type>
41 {
42 typedef CI ContainerIndex;
43 typedef typename C::mapped_type::mapped_type Weight;
44
45 const ContainerIndex& containerIndex() const
46 {
47 return *(this->first);
48 }
49
50 const Weight& weight() const
51 {
52 return this->second;
53 }
54 };
55
56 //typedef std::pair<CI,typename C::mapped_type::mapped_type> ConstraintsEntry;
57
58 typedef std::vector<ConstraintsEntry> ConstraintsVector;
59 typedef typename ConstraintsVector::const_iterator ConstraintsIterator;
60
61 LFSIndexCacheBase(const LFS& lfs, const C& constraints, bool enable_constraints_caching)
62 : _lfs(lfs)
63 , _enable_constraints_caching(enable_constraints_caching)
64 , _container_indices(lfs.maxSize())
65 , _dof_flags(lfs.maxSize(),0)
66 , _constraints_iterators(lfs.maxSize())
67 , _inverse_cache_built(false)
68 , _gfs_constraints(constraints)
69 {
70 }
71
72 void update()
73 {}
74
75 DI dofIndex(size_type i) const
76 {
77 return _lfs.dofIndex(i);
78 }
79
80 CI containerIndex(size_type i) const
81 {
82 return _lfs.containerIndex(i);
83 }
84
85 CI containerIndex(const DI& i) const
86 {
87 return i;
88 }
89
90 bool isConstrained(size_type i) const
91 {
92 return _dof_flags[i] & DOF_CONSTRAINED;
93 }
94
95 bool isDirichletConstraint(size_type i) const
96 {
97 return _dof_flags[i] & DOF_DIRICHLET;
98 }
99
100 ConstraintsIterator constraintsBegin(size_type i) const
101 {
102 assert(isConstrained(i));
103 return _constraints_iterators[i].first;
104 }
105
106 ConstraintsIterator constraintsEnd(size_type i) const
107 {
108 assert(isConstrained(i));
109 return _constraints_iterators[i].second;
110 }
111
112 const LocalFunctionSpace& localFunctionSpace() const
113 {
114 return _lfs;
115 }
116
117 size_type size() const
118 {
119 return _lfs.size();
120 }
121
122 bool constraintsCachingEnabled() const
123 {
124 return _enable_constraints_caching;
125 }
126
127 private:
128
129 const LFS& _lfs;
130 const bool _enable_constraints_caching;
131 CIVector _container_indices;
132 std::vector<unsigned char> _dof_flags;
133 std::vector<std::pair<ConstraintsIterator,ConstraintsIterator> > _constraints_iterators;
134 mutable CIMap _container_index_map;
135 ConstraintsVector _constraints;
136 mutable bool _inverse_cache_built;
137 mutable InverseMap _inverse_map;
138
139 const C& _gfs_constraints;
140
141 };
142
143
144 template<typename LFS>
145 class LFSIndexCacheBase<LFS,EmptyTransformation,Experimental::DuneFunctionsCacheTag, false>
146 {
147
148 public:
149
150 using LocalFunctionSpace = LFS;
151 using GFS = typename LFS::Traits::GridFunctionSpace;
152 using Ordering = typename GFS::Ordering;
153 using DOFIndex = typename Ordering::Traits::DOFIndex;
154 using DI = DOFIndex;
155 using ContainerIndex = typename Ordering::Traits::ContainerIndex;
156 using CI = ContainerIndex;
157 using size_type = std::size_t;
158
159 struct ConstraintsEntry
160 : public std::pair<const CI*,double>
161 {
162 typedef CI ContainerIndex;
163 typedef double Weight;
164
165 const ContainerIndex& containerIndex() const
166 {
167 return *(this->first);
168 }
169
170 const Weight& weight() const
171 {
172 return this->second;
173 }
174 };
175
176 typedef std::vector<ConstraintsEntry> ConstraintsVector;
177 typedef typename ConstraintsVector::const_iterator ConstraintsIterator;
178
179 typedef std::vector<CI> CIVector;
180 typedef std::unordered_map<DI,CI> CIMap;
181
182 explicit LFSIndexCacheBase(const LFS& lfs)
183 : _lfs(lfs)
184 {}
185
186 template<typename C>
187 LFSIndexCacheBase(const LFS& lfs, const C& c, bool enable_constraints_caching)
188 : _lfs(lfs)
189 {}
190
191 DI dofIndex(size_type i) const
192 {
193 return _lfs.dofIndex(i);
194 }
195
196 CI containerIndex(size_type i) const
197 {
198 return _lfs.containerIndex(i);
199 }
200
201 CI containerIndex(const DI& i) const
202 {
203 return _lfs.containerIndex(i);
204 }
205
206 bool isConstrained(size_type i) const
207 {
208 return false;
209 }
210
211 bool isDirichletConstraint(size_type i) const
212 {
213 return false;
214 }
215
216 ConstraintsIterator constraintsBegin(size_type i) const
217 {
218 return _constraints.begin();
219 }
220
221 ConstraintsIterator constraintsEnd(size_type i) const
222 {
223 return _constraints.end();
224 }
225
226 const LocalFunctionSpace& localFunctionSpace() const
227 {
228 return _lfs;
229 }
230
231 size_type size() const
232 {
233 return _lfs.size();
234 }
235
236 bool constraintsCachingEnabled() const
237 {
238 return false;
239 }
240
241 void update()
242 {}
243
244 private:
245
246 const LFS& _lfs;
247 CIVector _container_indices;
248 mutable CIMap _container_index_map;
249 const ConstraintsVector _constraints;
250
251 };
252
253 } // namespace PDELab
254} // namespace Dune
255
256#endif // DUNE_PDELAB_GRIDFUNCTIONSPACE_DUNEFUNCTIONSLFSINDEXCACHE_HH
void constraints(const GFS &gfs, CG &cg, const bool verbose=false)
construct constraints
Definition: constraints.hh:749
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)