G+Smo  25.01.0
Geometry + Simulation Modules
 
Loading...
Searching...
No Matches
gsContainerBasis.h
Go to the documentation of this file.
1
16#pragma once
17
18#include<gsCore/gsBasis.h>
19
20namespace gismo
21{
22
23 template<short_t d,class T>
24 class gsContainerBasis : public gsBasis<T>
25 {
26
27 public:
29 typedef memory::shared_ptr< gsContainerBasis > Ptr;
30
32 typedef memory::unique_ptr< gsContainerBasis > uPtr;
33
34 gsContainerBasis(index_t numSubspaces)
35 {
36 basisContainer.clear();
37 for (index_t i = 0; i != numSubspaces; i++)
38 basisContainer.push_back(gsTensorBSplineBasis<d,T>());
39
40 }
41
42 ~gsContainerBasis() {};
43
44 // boundary(boxSide const & side)
45
46// implementations of gsBasis
47 public:
48
49 static uPtr make( const gsContainerBasis& other)
50 { return uPtr( new gsContainerBasis( other ) ); }
51
52
53 public:
54
55 GISMO_CLONE_FUNCTION(gsContainerBasis)
56
57 short_t domainDim() const
58 {
59 return d;
60 }
61
62 void connectivity(const gsMatrix<T> & nodes,
63 gsMesh<T> & mesh) const
64 {
65 GISMO_UNUSED(nodes); GISMO_UNUSED(mesh);
66 GISMO_NO_IMPLEMENTATION;
67 }
68
69 memory::unique_ptr<gsGeometry<T> > makeGeometry( gsMatrix<T> coefs ) const
70 {
71 GISMO_UNUSED(coefs);
72 GISMO_NO_IMPLEMENTATION;
73 }
74
75 std::ostream &print(std::ostream &os) const
76 {
77 GISMO_UNUSED(os);
78 GISMO_NO_IMPLEMENTATION;
79 }
80
81 void uniformRefine(int numKnots = 1, int mul=1, int dir=-1)
82 {
83 for (size_t i=0; i< basisContainer.size(); ++i)
84 basisContainer[i].uniformRefine(numKnots,mul,dir);
85 }
86
87 // Returm max degree of all the spaces, otherwise i =
88 short_t degree(short_t dir) const
89 {
90 short_t deg = 0;
91 for (size_t i=0; i< basisContainer.size(); ++i)
92 if (basisContainer[i].degree(dir) > deg)
93 deg = basisContainer[i].degree(dir);
94
95 return deg;
96 }
97
98 index_t size() const {
99 index_t sz = 0;
100 for (size_t i=0; i< basisContainer.size(); ++i)
101 sz += basisContainer[i].size();
102 return sz;
103 }
104
105 void reverse()
106 {
107 for (size_t i=0; i< basisContainer.size(); ++i)
108 {
109 if (gsTensorBSplineBasis<d, T> * mb =
110 dynamic_cast<gsTensorBSplineBasis<d, T> *>(&basisContainer[i]) )
111 {
112 gsTensorBSplineBasis<d, T> basis_temp = dynamic_cast<gsTensorBSplineBasis<d,T>&>(basisContainer[i]);
113 gsTensorBSplineBasis<d, T> newTensorBasis(basis_temp.knots(1),basis_temp.knots(0));
114 basisContainer[i] = newTensorBasis;
115 }
116 /*
117 else if (gsTensorNurbsBasis<d, T> * mb =
118 dynamic_cast<gsTensorNurbsBasis<d, T> *>(&basisContainer[i]) )
119 {
120 gsTensorNurbsBasis<d, T> basis_temp = dynamic_cast<gsTensorNurbsBasis<d,T>&>(basisContainer[i]);
121 basis_temp.swapDirections(1,0);
122 basisContainer[i] = basis_temp;
123 }
124 */
125 else
126 gsInfo << "Works for now just with gsTensorBSplineBasis<d, T> \n";
127
128 }
129 }
130
131 index_t nPieces() const {return basisContainer.size();}
132
133 typename gsBasis<T>::domainIter makeDomainIterator(const boxSide & side) const override
134 {
135 // Using the inner basis for iterating
136 return basisContainer[0].makeDomainIterator(side);
137 }
138
139 typename gsBasis<T>::domainIter makeDomainIterator() const override
140 {
141 // Using the inner basis for iterating
142 return basisContainer[0].makeDomainIterator();
143 }
144
145 gsBasis<T> & component(short_t i) override
146 {
147 return basisContainer[0].component(i);
148 }
149
150 const gsBasis<T> & component(short_t i) const override
151 {
152 return basisContainer[0].component(i);
153 }
154
155/* void matchWith(const boundaryInterface & bi, const gsBasis<T> & other,
156 gsMatrix<index_t> & bndThis, gsMatrix<index_t> & bndOther) const;
157*/
158
159 gsMatrix<T> support() const override
160 {
161 return basisContainer[0].support();
162 }
163
164 gsMatrix<T> support(const index_t & i) const override
165 {
166 return basisContainer[0].support(i);
167 }
168
169 gsBasis<T>* boundaryBasis_impl(boxSide const & s) const
170 {
171 /*
172 if (basisContainer.size() != 1)
173 {
174 typename gsBSplineTraits<d-1, T>::Basis* bBasis = new typename gsBSplineTraits<d-1, T>::Basis(*basisContainer[s.index()].boundaryBasis(s));
175 return bBasis;
176 }
177 else {
178 typename gsBSplineTraits<d-1, T>::Basis* bBasis = new typename gsBSplineTraits<d-1, T>::Basis(*basisContainer[0].boundaryBasis(s));
179 return bBasis;
180 }
181 */
182
183 gsDebugVar(s.index());
184 // Maybe not working for approx C1 Basis functions
185 typename gsBSplineTraits<d-1, T>::Basis* bBasis = new typename gsBSplineTraits<d-1, T>::Basis(*basisContainer[0].boundaryBasis(s));
186 gsDebugVar(*bBasis);
187 return bBasis;
188 }
189
190 void active_into(const gsMatrix<T> & u, gsMatrix<index_t> & result) const
191 {
192 GISMO_ASSERT(u.rows() == d, "Dimension of the points in active_into is wrong");
193 //GISMO_ASSERT(u.cols() == 1, "Active_into is wrong");
194
195 index_t nr = 0;
196 std::vector<gsMatrix<index_t>> result_temp;
197 result_temp.resize(basisContainer.size());
198 for (size_t i=0; i< basisContainer.size(); ++i)
199 {
200 basisContainer[i].active_into(u, result_temp[i]);
201 nr += result_temp[i].rows();
202 }
203
204 result.resize(nr,u.cols());
205
206 index_t shift = 0;
207 index_t shift_rows = 0;
208 for (size_t i=0; i< basisContainer.size(); ++i)
209 {
210 result_temp[i].array() += shift;
211 result.block(shift_rows, 0, result_temp[i].rows(), u.cols()) = result_temp[i];
212 shift += basisContainer[i].size();
213 shift_rows += result_temp[i].rows();
214
215 }
216 }
217
218 void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
219 {
220 result.resize(0, u.cols());
221 for (size_t i=0; i< basisContainer.size(); ++i)
222 {
223 gsMatrix<T> result_temp;
224 basisContainer[i].eval_into(u, result_temp);
225 result.conservativeResize(result.rows()+result_temp.rows(), result.cols());
226 result.bottomRows(result_temp.rows()) = result_temp;
227
228 }
229 }
230
231 void deriv_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
232 {
233 result.resize(0, u.cols());
234 for (size_t i=0; i< basisContainer.size(); ++i)
235 {
236 gsMatrix<T> result_temp;
237 basisContainer[i].deriv_into(u, result_temp);
238 result.conservativeResize(result.rows()+result_temp.rows(), result.cols());
239 result.bottomRows(result_temp.rows()) = result_temp;
240
241 }
242 }
243
244 void deriv2_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
245 {
246 result.resize(0, u.cols());
247 for (size_t i=0; i< basisContainer.size(); ++i)
248 {
249 gsMatrix<T> result_temp;
250 basisContainer[i].deriv2_into(u, result_temp);
251 result.conservativeResize(result.rows()+result_temp.rows(), result.cols());
252 result.bottomRows(result_temp.rows()) = result_temp;
253
254 }
255 }
256
257 void matchWith(const boundaryInterface &bi, const gsBasis<T> &other, gsMatrix<index_t> &bndThis,
258 gsMatrix<index_t> &bndOther, int) const override
259 {
260 GISMO_UNUSED(bi); GISMO_UNUSED(other); GISMO_UNUSED(bndThis); GISMO_UNUSED(bndOther);
261 // First side
262 // edge
263 gsDebug << "matchWith() function is not implemented yet \n";
264
265 /*
266 gsMatrix<index_t> indizes1, indizes2;
267 short_t side_id = bi.first().side().index();
268 indizes1 = basisG1Container[side_id].boundaryOffset(boxSide(bi.first().side()), 0);
269 indizes2 = basisG1Container[side_id].boundaryOffset(boxSide(bi.first().side()), 1);
270
271 bndThis = indizes1;
272 bndThis.resize(indizes1.rows() + indizes2.rows(), indizes1.cols());
273 bndThis.block(0, 0, indizes1.rows(), indizes1.cols()) = indizes1;
274 bndThis.block(indizes1.rows(), 0, indizes2.rows(), indizes1.cols()) = indizes2;
275
276 index_t shift = 0;
277 for (index_t i = 0; i < side_id; ++i)
278 shift += basisG1Container[side_id].size();
279
280 bndOther.array() += shift;
281
282 bndOther = bndThis;
283 */
284 }
285
286
287 gsMatrix<index_t> boundaryOffset(const boxSide & bside, index_t offset) const
288 {
289 gsMatrix<index_t> result;
290
291 // TODO TEST WITH:
292 /*
293 index_t shift = 0;
294 result.resize(0, 1);
295 for (size_t i=0; i< basisContainer.size(); ++i)
296 {
297 gsMatrix<index_t> result_temp;
298 result_temp = basisContainer[i].boundaryOffset(bside, offset);
299 result_temp.array() += shift;
300 result.conservativeResize(result.rows()+result_temp.rows(), 1);
301 result.bottomRows(result_temp.rows()) = result_temp;
302 shift += basisContainer[i].size();
303 }
304 */
305 // Edges
306 short_t side_id = bside.index();
307 if (basisContainer.size() != 1)
308 {
309 index_t shift = 0;
310 for (index_t i=0; i< side_id; ++i)
311 shift += basisContainer[i].size();
312 result = basisContainer[side_id].boundaryOffset(bside, offset);
313 result.array() += shift;
314 }
315 else if (basisContainer.size() == 1)
316 result = basisContainer[0].boundaryOffset(bside, offset);
317
318
319 // Vertices:
320 if (basisContainer.size() != 1)
321 {
322 std::vector<boxCorner> containedCorners;
323 bside.getContainedCorners(d, containedCorners);
324
325 GISMO_ASSERT(containedCorners.size() != 0, "No contained corner");
326
327
328 for (size_t nc = 0; nc < containedCorners.size(); nc++)
329 {
330 index_t corner_id = containedCorners[nc].m_index + 4; // + 4 included because of 4 sides!
331
332 index_t shift = 0;
333 for (index_t i = 0; i < corner_id; ++i)
334 shift += basisContainer[i].size();
335
336 gsMatrix<index_t> result_temp;
337 result_temp = basisContainer[corner_id].boundaryOffset(bside, offset);
338 result_temp.array() += shift;
339
340 index_t r_rows = result.rows() + result_temp.rows();
341 result.conservativeResize(r_rows, 1);
342 result.bottomRows(result_temp.rows()) = result_temp;
343 }
344 }
345 return result;
346 }
347
348 index_t functionAtCorner(const boxCorner & corner) const
349 {
350 if (basisContainer.size() != 1)
351 {
352 index_t corner_id = corner.m_index + 4; // + 4 included because of 4 sides!
353 index_t shift = 0;
354 for (index_t i = 0; i < corner_id; ++i)
355 shift += basisContainer[i].size();
356
357 return basisContainer[corner_id].functionAtCorner(corner) + shift;
358 }
359 else
360 return basisContainer[0].functionAtCorner(corner);
361
362 }
363
364// implementations of gsContainerBasis
365 public:
366
367 // basisContainer:
368 void setBasis(index_t row, gsTensorBSplineBasis<d,T> basis) { basisContainer[row] = basis; }
369
370 const gsBasis<T> & piece(const index_t k) const { return basisContainer[k]; }
371 // basisContainer END
372
373 //std::vector<gsTensorBSplineBasis<d,T>> & getBasisContainer() { return basisContainer; }
374
375 // Data members
376 protected:
377
378 // Collection of the subspaces
379 std::vector<gsTensorBSplineBasis<d, T>> basisContainer;
380 };
381
382}
gismo::gsBasis::uniformRefine
virtual void uniformRefine(int numKnots=1, int mul=1, short_t dir=-1)
Refine the basis uniformly by inserting numKnots new knots with multiplicity mul on each knot span.
Definition gsBasis.hpp:603
gismo::gsBasis::boundaryBasis
gsBasis< T >::uPtr boundaryBasis(boxSide const &s)
Returns the boundary basis for side s.
gismo::gsBasis::matchWith
virtual void matchWith(const boundaryInterface &bi, const gsBasis< T > &other, gsMatrix< index_t > &bndThis, gsMatrix< index_t > &bndOther, index_t offset=0) const
Computes the indices of DoFs that match on the interface bi. The interface is assumed to be a common ...
Definition gsBasis.hpp:707
gismo::gsFunctionSet::basis
const gsBasis< T > & basis(const index_t k) const
Helper which casts and returns the k-th piece of this function set as a gsBasis.
Definition gsFunctionSet.hpp:33
gsBasis.h
Provides declaration of Basis abstract interface.
short_t
#define short_t
Definition gsConfig.h:35
index_t
#define index_t
Definition gsConfig.h:32
gsDebug
#define gsDebug
Definition gsDebug.h:61
GISMO_NO_IMPLEMENTATION
#define GISMO_NO_IMPLEMENTATION
Definition gsDebug.h:129
GISMO_UNUSED
#define GISMO_UNUSED(x)
Definition gsDebug.h:112
gsInfo
#define gsInfo
Definition gsDebug.h:43
GISMO_ASSERT
#define GISMO_ASSERT(cond, message)
Definition gsDebug.h:89
gismo
The G+Smo namespace, containing all definitions for the library.