1 /** @file gsBasis.hpp
2
3 @brief Provides implementation of Basis default operations.
4
5 This file is part of the G+Smo library.
6
7 This Source Code Form is subject to the terms of the Mozilla Public
8 License, v. 2.0. If a copy of the MPL was not distributed with this
9 file, You can obtain one at http://mozilla.org/MPL/2.0/.
10
11 Author(s): A. Mantzaflaris
12 */
13
14 #pragma once
15
16 #include <gsCore/gsBasisFun.h>
17 #include <gsCore/gsDomainIterator.h>
18 #include <gsCore/gsBoundary.h>
19 #include <gsCore/gsGeometry.h>
20
21 namespace gismo
22 {
23
24 template<class T>
gsBasis()25 gsBasis<T>::gsBasis()
26 { }
27
28 template<class T>
gsBasis(const gsBasis & other)29 gsBasis<T>::gsBasis(const gsBasis& other) : Base(other)
30 { }
31
32 template<class T>
~gsBasis()33 gsBasis<T>::~gsBasis()
34 { }
35
36 template<class T>
function(index_t i) const37 gsBasisFun<T> gsBasis<T>::function(index_t i) const
38 {
39 return gsBasisFun<T>(*this,i);
40 }
41
42
43 // Evaluates a linear combination of basis functions (default implementation)
44 template<class T>
evalFunc_into(const gsMatrix<T> & u,const gsMatrix<T> & coefs,gsMatrix<T> & result) const45 void gsBasis<T>::evalFunc_into(const gsMatrix<T> &u,
46 const gsMatrix<T> & coefs,
47 gsMatrix<T>& result) const
48 {
49 gsMatrix<T> B ;
50 gsMatrix<index_t> actives;
51
52 // compute function values
53 this->eval_into(u,B);
54 // compute active functions
55 this->active_into(u,actives);
56
57 // compute result as linear combination of
58 // "coefs(actives)" and B
59 linearCombination_into( coefs, actives, B, result );
60 }
61
62
63 // Evaluates the Jacobian of the function given by coefs (default implementation)
64 // For each point, result contains a geomDim x parDim matrix block containing the Jacobian matrix
65 template<class T>
jacobianFunc_into(const gsMatrix<T> & u,const gsMatrix<T> & coefs,gsMatrix<T> & result) const66 void gsBasis<T>::jacobianFunc_into(const gsMatrix<T> &u, const gsMatrix<T> & coefs, gsMatrix<T>& result) const
67 {
68 const index_t n = coefs.cols();
69 const index_t numPts = u.cols(); // at how many points to evaluate the gradients
70 const index_t pardim = this->dim();
71
72 result.setZero( n, pardim * numPts );
73 gsMatrix<T> B;
74 gsMatrix<index_t> ind;
75
76 this->deriv_into(u,B);
77 this->active_into(u,ind);
78 const index_t numAct=ind.rows();
79
80 for (index_t p = 0; p < numPts; ++p) // p = point
81 for (index_t c=0; c<n; ++c ) // c = component
82 for ( index_t a=0; a< numAct ; ++a ) // a = active function
83 {
84 result.block(c,p*pardim, 1, pardim).noalias() +=
85 coefs(ind(a,p), c) * B.block(a*pardim, p, pardim, 1).transpose();
86 }
87 }
88
89 // Evaluates the partial derivatives of the function given by coefs (default implementation)
90 // For each point, result contains one column with stacked gradients of the components
91 template<class T>
derivFunc_into(const gsMatrix<T> & u,const gsMatrix<T> & coefs,gsMatrix<T> & result) const92 void gsBasis<T>::derivFunc_into(const gsMatrix<T> &u, const gsMatrix<T> & coefs, gsMatrix<T>& result) const
93 {
94 gsMatrix<T> B;
95 gsMatrix<index_t> actives;
96
97 // compute first derivatives
98 this->deriv_into(u,B);
99 // compute active functions
100 this->active_into(u,actives);
101
102 // compute result as linear combination of
103 // "coefs(actives)" and B
104 linearCombination_into( coefs, actives, B, result );
105 }
106
107 // Evaluates the second derivatives of the function given by coefs (default implementation)
108 template<class T>
deriv2Func_into(const gsMatrix<T> & u,const gsMatrix<T> & coefs,gsMatrix<T> & result) const109 void gsBasis<T>::deriv2Func_into(const gsMatrix<T> & u,
110 const gsMatrix<T> & coefs,
111 gsMatrix<T>& result) const
112 {
113 gsMatrix<T> B;
114 gsMatrix<index_t> actives;
115
116 // compute second derivatives
117 this->deriv2_into(u,B);
118 // compute active functions
119 this->active_into(u,actives);
120
121 // compute result as linear combination of
122 // "coefs(actives)" and B
123 linearCombination_into( coefs, actives, B, result );
124 }
125
126 template<class T>
evalAllDersFunc_into(const gsMatrix<T> & u,const gsMatrix<T> & coefs,const unsigned n,std::vector<gsMatrix<T>> & result) const127 void gsBasis<T>::evalAllDersFunc_into(const gsMatrix<T> &u,
128 const gsMatrix<T> & coefs,
129 const unsigned n,
130 std::vector< gsMatrix<T> >& result) const
131 {
132 // resize result so that it will hold
133 // function values and up to the n-th derivatives
134 result.resize(n+1);
135
136 // B will contain the derivatives up to order n
137 std::vector< gsMatrix<T> >B;
138 // actives will contain the indices of the basis functions
139 // which are active at the evaluation points
140 gsMatrix<index_t> actives;
141
142 this->evalAllDers_into(u,n,B);
143 this->active_into(u,actives);
144
145 // for derivatives 0 to n, evaluate the function by linear combination
146 // of coefficients with the respective function values/derivatives
147 for( unsigned i = 0; i <= n; i++)
148 linearCombination_into( coefs, actives, B[i], result[i] );
149 }
150
151
152 template<class T>
linearCombination_into(const gsMatrix<T> & coefs,const gsMatrix<index_t> & actives,const gsMatrix<T> & values,gsMatrix<T> & result)153 void gsBasis<T>::linearCombination_into(const gsMatrix<T> & coefs,
154 const gsMatrix<index_t> & actives,
155 const gsMatrix<T> & values,
156 gsMatrix<T> & result)
157 {
158 const index_t numPts = values.cols() ;
159 const index_t tarDim = coefs.cols() ;
160 const index_t stride = values.rows() / actives.rows();
161
162 GISMO_ASSERT( actives.rows() * stride == values.rows(),
163 "Number of values and actives does not fit together");
164
165 result.resize( tarDim * stride, numPts );
166 result.setZero();
167
168 for ( index_t pt = 0; pt < numPts; ++pt ) // For pt, i.e., for every column of u
169 for ( index_t i = 0; i < actives.rows(); ++i ) // for all nonzero basis functions
170 for ( index_t c = 0; c < tarDim; ++c ) // for all components of the geometry
171 {
172 result.block( stride * c, pt, stride, 1).noalias() +=
173 coefs( actives(i,pt), c) * values.block( stride * i, pt, stride, 1);
174 }
175 }
176
177
178 template<class T>
laplacian(const gsMatrix<T> & u) const179 inline gsMatrix<T> gsBasis<T>::laplacian(const gsMatrix<T> & u ) const
180 {
181 gsMatrix<T> tmp;
182 this->deriv2_into(u,tmp);
183 return tmp.colwise().sum();
184 }
185
186 template<class T> inline
collocationMatrix(const gsMatrix<T> & u,gsSparseMatrix<T> & result) const187 void gsBasis<T>::collocationMatrix(const gsMatrix<T> & u, gsSparseMatrix<T> & result) const
188 {
189 result.resize( u.cols(), this->size() );
190 gsMatrix<T> ev;
191 gsMatrix<index_t> act;
192
193 eval_into (u.col(0), ev);
194 active_into(u.col(0), act);
195 result.reservePerColumn( act.rows() );
196 for (index_t i=0; i!=act.rows(); ++i)
197 result.insert(0, act.at(i) ) = ev.at(i);
198
199 for (index_t k=1; k!=u.cols(); ++k)
200 {
201 eval_into (u.col(k), ev );
202 active_into(u.col(k), act);
203 for (index_t i=0; i!=act.rows(); ++i)
204 result.insert(k, act.at(i) ) = ev.at(i);
205 }
206
207 result.makeCompressed();
208 }
209
210 template<class T> inline
interpolateData(gsMatrix<T> const & vals,gsMatrix<T> const & pts) const211 memory::unique_ptr<gsGeometry<T> > gsBasis<T>::interpolateData( gsMatrix<T> const& vals,
212 gsMatrix<T> const& pts) const
213 {
214 GISMO_ASSERT (dim() == pts.rows() , "Wrong dimension of the points("<<
215 pts.rows()<<", expected "<<dim() <<").");
216 GISMO_ASSERT (this->size() == pts.cols() , "Expecting as many points as the basis functions." );
217 GISMO_ASSERT (this->size() == vals.cols(), "Expecting as many values as the number of points." );
218
219 gsSparseMatrix<T> Cmat;
220 collocationMatrix(pts, Cmat);
221 gsMatrix<T> x ( this->size(), vals.rows());
222
223 // typename gsSparseSolver<T>::BiCGSTABIdentity solver( Cmat );
224 // typename gsSparseSolver<T>::BiCGSTABDiagonal solver( Cmat );
225 // typename gsSparseSolver<T>::QR solver( Cmat );
226 typename gsSparseSolver<T>::BiCGSTABILUT solver( Cmat );
227
228 // Solves for many right hand side columns
229 x = solver.solve( vals.transpose() );
230
231 // gsDebug <<"gs Interpolate error : " << solver.error() << std::"\n";
232 // gsDebug <<"gs Interpolate iters : " << solver.iterations() << std::"\n";
233 // gsDebug <<"intpl sol : " << x.transpose() << std::"\n";
234
235 return makeGeometry( give(x) );
236 }
237
238 template<class T> inline
interpolateAtAnchors(gsMatrix<T> const & vals) const239 memory::unique_ptr<gsGeometry<T> > gsBasis<T>::interpolateAtAnchors(gsMatrix<T> const & vals) const
240 {
241 GISMO_ASSERT (this->size() == vals.cols(),
242 "Expecting as many values as the number of basis functions." );
243 gsMatrix<T> pts;
244 anchors_into(pts);
245 return interpolateData(vals, pts);
246 }
247
248 /*
249 template<class T> inline
250 gsGeometry<T> * gsBasis<T>::interpolateAtAnchors(gsFunction<T> const & func) const
251 {
252 gsMatrix<T> pts, vals;
253 anchors_into(pts);
254 func.eval_into(pts, vals);
255 return interpolateData(vals, pts);
256 }
257
258 template<class T>
259 gsGeometry<T> * gsBasis<T>::projectL2(gsFunction<T> const & func) const
260 {
261
262 return NULL;
263 }
264 */
265
266 template<class T> inline
anchors_into(gsMatrix<T> &) const267 void gsBasis<T>::anchors_into(gsMatrix<T>&) const
268 { GISMO_NO_IMPLEMENTATION }
269
270 template<class T> inline
anchor_into(index_t,gsMatrix<T> &) const271 void gsBasis<T>::anchor_into(index_t, gsMatrix<T>&) const
272 { GISMO_NO_IMPLEMENTATION }
273
274
275 template<class T> inline
connectivityAtAnchors(gsMesh<T> & mesh) const276 void gsBasis<T>::connectivityAtAnchors(gsMesh<T> & mesh) const
277 {
278 gsMatrix<T> nodes = anchors();
279 nodes.transposeInPlace();// coefficient vectors have ctrl points at rows
280 connectivity(nodes, mesh);
281 }
282
283 template<class T>
connectivity(const gsMatrix<T> &,gsMesh<T> &) const284 void gsBasis<T>::connectivity(const gsMatrix<T> &, gsMesh<T> &) const
285 { GISMO_NO_IMPLEMENTATION }
286
287 template<class T>
active_into(const gsMatrix<T> &,gsMatrix<index_t> &) const288 void gsBasis<T>::active_into(const gsMatrix<T> &, gsMatrix<index_t>&) const
289 { GISMO_NO_IMPLEMENTATION }
290
291 template <class T>
isActive(const index_t,const gsVector<T> &) const292 bool gsBasis<T>::isActive(const index_t, const gsVector<T> &) const
293 { GISMO_NO_IMPLEMENTATION }
294
295 template<class T>
numActive_into(const gsMatrix<T> &,gsVector<index_t> &) const296 void gsBasis<T>::numActive_into(const gsMatrix<T> &, gsVector<index_t>&) const
297 { GISMO_NO_IMPLEMENTATION }
298
299 template<class T>
activeCoefs_into(const gsVector<T> &,const gsMatrix<T> &,gsMatrix<T> &) const300 void gsBasis<T>::activeCoefs_into(const gsVector<T> &, const gsMatrix<T> &,
301 gsMatrix<T>&) const
302 { GISMO_NO_IMPLEMENTATION }
303
304 template<class T>
305 gsMatrix<index_t>
allBoundary() const306 gsBasis<T>::allBoundary() const
307 { GISMO_NO_IMPLEMENTATION }
308
309 template<class T>
310 gsMatrix<index_t>
boundaryOffset(boxSide const &,index_t) const311 gsBasis<T>::boundaryOffset(boxSide const &,index_t) const
312 { GISMO_NO_IMPLEMENTATION }
313
314 template<class T>
315 index_t
functionAtCorner(boxCorner const &) const316 gsBasis<T>::functionAtCorner(boxCorner const &) const
317 { GISMO_NO_IMPLEMENTATION }
318
319 /// @cond
320 template<class T>
boundaryBasis_impl(boxSide const &) const321 gsBasis<T>* gsBasis<T>::boundaryBasis_impl(boxSide const &) const
322 { GISMO_NO_IMPLEMENTATION }
323 /// @endcond
324
325 template<class T>
componentBasis(boxComponent b) const326 typename gsBasis<T>::uPtr gsBasis<T>::componentBasis(boxComponent b) const
327 {
328 GISMO_ASSERT( b.totalDim() == this->dim(), "The dimensions do not agree." );
329
330 const short_t dim = this->dim();
331
332 uPtr result;
333 short_t d=0;
334 for (short_t i=0; i<dim; ++i)
335 {
336 boxComponent::location loc = b.locationForDirection(i);
337 if (loc)
338 {
339 if (result)
340 result = result->boundaryBasis( boxSide(loc+2*d) );
341 else
342 result = this->boundaryBasis( boxSide(loc+2*d) );
343 }
344 else
345 ++d;
346 }
347
348 if (!result)
349 result = clone();
350
351 return result;
352 }
353
354 template<class T>
componentBasis_withIndices(boxComponent b,gsMatrix<index_t> & indices,bool noBoundary) const355 typename gsBasis<T>::uPtr gsBasis<T>::componentBasis_withIndices(boxComponent b, gsMatrix<index_t>& indices, bool noBoundary) const
356 {
357 GISMO_ASSERT( b.totalDim() == this->dim(), "The dimensions do not agree." );
358 const short_t dim = this->dim();
359
360 uPtr result;
361 short_t d=0;
362 for (short_t i=0; i<dim; ++i)
363 {
364 boxComponent::location loc = b.locationForDirection(i);
365 if (loc)
366 {
367 if (result)
368 {
369 gsMatrix<index_t> tmp = result->boundary( boxSide(loc+2*d) );
370 for (index_t j=0; j<tmp.size(); ++j)
371 tmp(j,0) = indices(tmp(j,0),0);
372 tmp.swap(indices);
373 result = result->boundaryBasis( boxSide(loc+2*d) );
374 }
375 else
376 {
377 indices = this->boundary( boxSide(loc+2*d) );
378 result = this->boundaryBasis( boxSide(loc+2*d) );
379 }
380 }
381 else
382 ++d;
383 }
384
385 if (!result)
386 {
387 result = clone();
388 const index_t sz = this->size();
389 indices.resize(sz,1);
390 for (index_t i=0;i<sz;++i)
391 indices(i,0) = i;
392 }
393
394 if (noBoundary && d > 0)
395 {
396
397 gsMatrix<index_t> bdy_indices = result->allBoundary();
398
399 const index_t indices_sz = indices.rows();
400 const index_t bdy_indices_sz = bdy_indices.rows();
401
402 // Copy all entries from indices to indices_cleaned except
403 // those with indices in bdy_indices
404
405 gsMatrix<index_t> indices_cleaned(indices_sz - bdy_indices_sz, 1);
406 index_t j = 0, t = 0;
407 for (index_t i = 0; i < indices_sz; ++i)
408 {
409 if (util::greater(i, bdy_indices(j, 0)) && j < bdy_indices_sz)
410 ++j;
411 if (util::less(i, bdy_indices(j, 0)) || j == bdy_indices_sz)
412 {
413 indices_cleaned(t, 0) = indices(i, 0);
414 ++t;
415 }
416 }
417 GISMO_ASSERT(t == indices_cleaned.rows(), "Internal error.");
418 indices.swap(indices_cleaned);
419 }
420
421 return result;
422
423 }
424
425 template<class T>
support() const426 gsMatrix<T> gsBasis<T>::support() const
427 { GISMO_NO_IMPLEMENTATION }
428
429 template<class T>
support(const index_t &) const430 gsMatrix<T> gsBasis<T>::support(const index_t &) const
431 { GISMO_NO_IMPLEMENTATION }
432
433 template<class T>
supportInterval(index_t dir) const434 gsMatrix<T> gsBasis<T>::supportInterval(index_t dir) const
435 { return support().row(dir); }
436
437 template<class T>
eval_into(const gsMatrix<T> &,gsMatrix<T> &) const438 void gsBasis<T>::eval_into(const gsMatrix<T> &, gsMatrix<T>&) const
439 { GISMO_NO_IMPLEMENTATION }
440
441 template<class T>
evalSingle_into(index_t,const gsMatrix<T> &,gsMatrix<T> &) const442 void gsBasis<T>::evalSingle_into(index_t, const gsMatrix<T> &, gsMatrix<T>&) const
443 { GISMO_NO_IMPLEMENTATION }
444
445 template<class T>
deriv_into(const gsMatrix<T> &,gsMatrix<T> &) const446 void gsBasis<T>::deriv_into(const gsMatrix<T> &, gsMatrix<T>&) const
447 { GISMO_NO_IMPLEMENTATION }
448
449 template<class T>
derivSingle_into(index_t,const gsMatrix<T> &,gsMatrix<T> &) const450 void gsBasis<T>::derivSingle_into(index_t,
451 const gsMatrix<T> &,
452 gsMatrix<T>&) const
453 { GISMO_NO_IMPLEMENTATION }
454
455 template<class T>
deriv2_into(const gsMatrix<T> &,gsMatrix<T> &) const456 void gsBasis<T>::deriv2_into(const gsMatrix<T> &, gsMatrix<T>&) const
457 { GISMO_NO_IMPLEMENTATION }
458
459 template<class T>
deriv2Single_into(index_t,const gsMatrix<T> &,gsMatrix<T> &) const460 void gsBasis<T>::deriv2Single_into(index_t,
461 const gsMatrix<T> &,
462 gsMatrix<T>&) const
463 { GISMO_NO_IMPLEMENTATION }
464
465 template<class T>
evalAllDers_into(const gsMatrix<T> & u,int n,std::vector<gsMatrix<T>> & result) const466 void gsBasis<T>::evalAllDers_into(const gsMatrix<T> & u, int n,
467 std::vector<gsMatrix<T> >& result) const
468 {
469 result.resize(n+1);
470
471 switch(n)
472 {
473 case 0:
474 eval_into(u, result[0]);
475 break;
476 case 1:
477 eval_into (u, result[0]);
478 deriv_into(u, result[1]);
479 break;
480 case 2:
481 eval_into (u, result[0]);
482 deriv_into (u, result[1]);
483 deriv2_into(u, result[2]);
484 break;
485 default:
486 GISMO_ERROR("evalAllDers implemented for order up to 2<"<<n<< " for "<<*this);
487 break;
488 }
489 }
490
491 template<class T>
evalAllDersSingle_into(index_t,const gsMatrix<T> &,int,gsMatrix<T> &) const492 void gsBasis<T>::evalAllDersSingle_into(index_t, const gsMatrix<T> &,
493 int, gsMatrix<T>&) const
494 { GISMO_NO_IMPLEMENTATION }
495
496 template<class T>
evalDerSingle_into(index_t,const gsMatrix<T> &,int,gsMatrix<T> &) const497 void gsBasis<T>::evalDerSingle_into(index_t, const
498 gsMatrix<T> &, int,
499 gsMatrix<T>&) const
500 { GISMO_NO_IMPLEMENTATION }
501
502
503 template<class T>
create() const504 typename gsBasis<T>::uPtr gsBasis<T>::create() const
505 { GISMO_NO_IMPLEMENTATION }
506
507
508 template<class T>
tensorize(const gsBasis &) const509 typename gsBasis<T>::uPtr gsBasis<T>::tensorize(const gsBasis &) const
510 { GISMO_NO_IMPLEMENTATION }
511
512 template<class T>
513 typename gsBasis<T>::domainIter
makeDomainIterator() const514 gsBasis<T>::makeDomainIterator() const
515 { GISMO_NO_IMPLEMENTATION }
516
517 template<class T>
518 typename gsBasis<T>::domainIter
makeDomainIterator(const boxSide &) const519 gsBasis<T>::makeDomainIterator(const boxSide &) const
520 { GISMO_NO_IMPLEMENTATION }
521
522 template<class T>
numElements() const523 size_t gsBasis<T>::numElements() const
524 { GISMO_NO_IMPLEMENTATION }
525
526 template<class T>
numElements(boxSide const &) const527 size_t gsBasis<T>::numElements(boxSide const &) const
528 { GISMO_NO_IMPLEMENTATION }
529
530 template<class T>
elementIndex(const gsVector<T> &) const531 size_t gsBasis<T>::elementIndex(const gsVector<T> &) const
532 { GISMO_NO_IMPLEMENTATION }
533
534 template<class T>
elementInSupportOf(index_t j) const535 gsMatrix<T> gsBasis<T>::elementInSupportOf(index_t j) const
536 { GISMO_NO_IMPLEMENTATION }
537
538 template<class T>
component(short_t) const539 const gsBasis<T>& gsBasis<T>::component(short_t) const
540 { GISMO_NO_IMPLEMENTATION }
541
542 template<class T>
component(short_t i)543 gsBasis<T>& gsBasis<T>::component(short_t i)
544 { return const_cast<gsBasis<T>&>(const_cast<const gsBasis<T>*>(this)->component(i));}
545
546 template<class T>
asElements(gsMatrix<T> const &,int) const547 std::vector<index_t> gsBasis<T>::asElements(gsMatrix<T> const &, int) const
548 { GISMO_NO_IMPLEMENTATION }
549
550 template<class T>
refine(gsMatrix<T> const &,int)551 void gsBasis<T>::refine(gsMatrix<T> const &, int)
552 { GISMO_NO_IMPLEMENTATION }
553
554 template<class T>
refineElements(std::vector<index_t> const &)555 void gsBasis<T>::refineElements(std::vector<index_t> const &)
556 { GISMO_NO_IMPLEMENTATION }
557
558 template<class T>
refineElements_withCoefs(gsMatrix<T> &,std::vector<index_t> const &)559 void gsBasis<T>::refineElements_withCoefs(gsMatrix<T> &,std::vector<index_t> const &)
560 { GISMO_NO_IMPLEMENTATION }
561
562 template<class T>
uniformRefine(int,int)563 void gsBasis<T>::uniformRefine(int, int)
564 { GISMO_NO_IMPLEMENTATION }
565
566 template<class T>
uniformRefine_withCoefs(gsMatrix<T> &,int,int)567 void gsBasis<T>::uniformRefine_withCoefs(gsMatrix<T>& , int , int )
568 { GISMO_NO_IMPLEMENTATION }
569
570 template<class T>
uniformRefine_withTransfer(gsSparseMatrix<T,RowMajor> &,int,int)571 void gsBasis<T>::uniformRefine_withTransfer(gsSparseMatrix<T,RowMajor> &,
572 int, int)
573 { GISMO_NO_IMPLEMENTATION }
574
575 template<class T>
uniformCoarsen(int)576 void gsBasis<T>::uniformCoarsen(int)
577 { GISMO_NO_IMPLEMENTATION }
578
579 template<class T>
uniformCoarsen_withTransfer(gsSparseMatrix<T,RowMajor> &,int)580 void gsBasis<T>::uniformCoarsen_withTransfer(gsSparseMatrix<T,RowMajor> &,
581 int)
582 { GISMO_NO_IMPLEMENTATION }
583
584 template<class T>
degreeElevate(short_t const &,short_t const)585 void gsBasis<T>::degreeElevate(short_t const &, short_t const)
586 { GISMO_NO_IMPLEMENTATION }
587
588 template<class T>
degreeReduce(short_t const &,short_t const)589 void gsBasis<T>::degreeReduce(short_t const &, short_t const)
590 { GISMO_NO_IMPLEMENTATION }
591
592 template<class T>
degreeIncrease(short_t const &,short_t const)593 void gsBasis<T>::degreeIncrease(short_t const &, short_t const)
594 { GISMO_NO_IMPLEMENTATION }
595
596 template<class T>
degreeDecrease(short_t const &,short_t const)597 void gsBasis<T>::degreeDecrease(short_t const &, short_t const)
598 { GISMO_NO_IMPLEMENTATION }
599
600 template<class T>
setDegree(short_t const & i)601 void gsBasis<T>::setDegree(short_t const& i)
602 {
603 const short_t dm = this->dim();
604 for (short_t k = 0; k!=dm; ++k)
605 {
606 const short_t p = this->degree(k);
607 if ( i > p )
608 {
609 this->degreeElevate(i-p, k);
610 }
611 else if ( i < p )
612 {
613 this->degreeReduce(p-i, k);
614 }
615 }
616 }
617
618 template<class T>
setDegreePreservingMultiplicity(short_t const & i)619 void gsBasis<T>::setDegreePreservingMultiplicity(short_t const& i)
620 {
621 for ( short_t d = 0; d < dim(); ++ d )
622 {
623 if ( i > degree(d) )
624 degreeIncrease(i-degree(d),d);
625 else if ( i < degree(d) )
626 degreeDecrease(-i+degree(d),d);
627 }
628 }
629
630 template<class T>
elevateContinuity(int const &)631 void gsBasis<T>::elevateContinuity(int const &)
632 { GISMO_NO_IMPLEMENTATION }
633
634 template<class T>
reduceContinuity(int const &)635 void gsBasis<T>::reduceContinuity(int const &)
636 { GISMO_NO_IMPLEMENTATION }
637
638 template<class T>
domain() const639 gsDomain<T> * gsBasis<T>::domain() const
640 { GISMO_NO_IMPLEMENTATION }
641
642 template<class T>
maxDegree() const643 short_t gsBasis<T>::maxDegree() const
644 { GISMO_NO_IMPLEMENTATION }
645
646 template<class T>
minDegree() const647 short_t gsBasis<T>::minDegree() const
648 { GISMO_NO_IMPLEMENTATION }
649
650 template<class T>
totalDegree() const651 short_t gsBasis<T>::totalDegree() const
652 { GISMO_NO_IMPLEMENTATION }
653
654 template<class T>
degree(short_t) const655 short_t gsBasis<T>::degree(short_t) const
656 { GISMO_NO_IMPLEMENTATION }
657
658 template<class T>
reverse()659 void gsBasis<T>::reverse()
660 { GISMO_NO_IMPLEMENTATION }
661
662 template<class T>
matchWith(const boundaryInterface &,const gsBasis<T> &,gsMatrix<index_t> &,gsMatrix<index_t> &) const663 void gsBasis<T>::matchWith(const boundaryInterface &, const gsBasis<T> &,
664 gsMatrix<index_t> &, gsMatrix<index_t> &) const
665 { GISMO_NO_IMPLEMENTATION }
666
667 template<class T>
getMinCellLength() const668 T gsBasis<T>::getMinCellLength() const
669 {
670 const domainIter it = this->makeDomainIterator();
671 T h = 0;
672 for (; it->good(); it->next() )
673 {
674 const T sz = it->getMinCellLength();
675 if ( sz < h || h == 0 ) h = sz;
676 }
677 return h;
678 }
679
680 template<class T>
getMaxCellLength() const681 T gsBasis<T>::getMaxCellLength() const
682 {
683 const domainIter it = this->makeDomainIterator();
684 T h = 0;
685 for (; it->good(); it->next() )
686 {
687 const T sz = it->getMaxCellLength();
688 if ( sz > h ) h = sz;
689 }
690 return h;
691 }
692
693 // gsBasis<T>::linearComb(active, evals, m_tmpCoefs, result);
694 // gsBasis<T>::jacobianFromGradients(active, grads, m_tmpCoefs, result);
695
696 /*
697 template<class T>
698 void gsBasis<T>::linearComb(const gsMatrix<index_t> & actives,
699 const gsMatrix<T> & basisVals,
700 const gsMatrix<T> & coefs,
701 gsMatrix<T>& result )
702 {
703 // basisVals.rows()==1 (or else basisVals.rows() == coefs.cols() and .cwiseProd)
704 result.resize(coefs.cols(), basisVals.cols()) ;
705
706 for ( index_t j=0; j!=basisVals.cols() ; j++ ) // for all basis function values
707 {
708 //todo grab result.col(j)
709 result.col(j) = basisVals(0,j) * coefs.row( actives(0,j) ) ;//transpose ?
710 for ( index_t i=1; i< actives.rows() ; i++ )
711 result.col(j) += basisVals(i,j) * coefs.row( actives(i,j) ) ;
712 }
713 }
714 */
715
716 }; // namespace gismo
717