1 /** @file gsGridHierarchy.hpp
2
3 @brief Coarsening algorithms for knot vectors and bases.
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): C. Hofreither, S. Takacs
12 */
13
14 #pragma once
15
16 #include <gsNurbs/gsKnotVector.h>
17 #include <gsNurbs/gsBSplineBasis.h>
18 #include <gsNurbs/gsTensorBSplineBasis.h>
19 #include <gsIO/gsOptionList.h>
20 #include <gsAssembler/gsAssemblerOptions.h>
21 #include <gsCore/gsMultiBasis.h>
22
23 namespace gismo
24 {
25
26 template <typename T>
buildByRefinement(gsMultiBasis<T> mBasis,const gsBoundaryConditions<T> & boundaryConditions,const gsOptionList & options,index_t levels,index_t numberOfKnotsToBeInserted,index_t multiplicityOfKnotsToBeInserted,index_t unk)27 gsGridHierarchy<T> gsGridHierarchy<T>::buildByRefinement(
28 gsMultiBasis<T> mBasis,
29 const gsBoundaryConditions<T>& boundaryConditions,
30 const gsOptionList& options,
31 index_t levels,
32 index_t numberOfKnotsToBeInserted,
33 index_t multiplicityOfKnotsToBeInserted,
34 index_t unk
35 )
36 {
37 gsGridHierarchy<T> result;
38 result.m_mBases.resize(levels);
39 result.m_transferMatrices.resize(levels-1);
40 result.m_mBases[0] = give(mBasis);
41 for ( index_t i=1; i<levels; ++i )
42 {
43 result.m_mBases[i] = result.m_mBases[i-1];
44 result.m_mBases[i].uniformRefine_withTransfer(
45 result.m_transferMatrices[i-1],
46 boundaryConditions,
47 options,
48 numberOfKnotsToBeInserted,
49 multiplicityOfKnotsToBeInserted,
50 unk
51 );
52 }
53 return result;
54 }
55
56 template <typename T>
buildByCoarsening(gsMultiBasis<T> mBasis,const gsBoundaryConditions<T> & boundaryConditions,const gsOptionList & options,index_t levels,index_t degreesOfFreedom,index_t unk)57 gsGridHierarchy<T> gsGridHierarchy<T>::buildByCoarsening(
58 gsMultiBasis<T> mBasis,
59 const gsBoundaryConditions<T>& boundaryConditions,
60 const gsOptionList& options,
61 index_t levels,
62 index_t degreesOfFreedom,
63 index_t unk
64 )
65 {
66 gsGridHierarchy<T> result;
67
68 result.m_mBases.push_back(give(mBasis));
69
70 index_t lastSize = result.m_mBases[0].totalSize();
71
72 for (index_t i = 0; i < levels-1 && lastSize > degreesOfFreedom; ++i)
73 {
74 gsSparseMatrix<T, RowMajor> transferMatrix;
75 gsMultiBasis<T> coarseMBasis = result.m_mBases[i];
76 coarseMBasis.uniformCoarsen_withTransfer(
77 transferMatrix,
78 boundaryConditions,
79 options,
80 1,
81 unk
82 );
83
84 index_t newSize = coarseMBasis.totalSize();
85 // If the number of dofs could not be decreased, then cancel. However, if only the number
86 // of levels was specified, then this should be ignored (the caller might need to have a
87 // fixed number of levels).
88 if (lastSize <= newSize && degreesOfFreedom > 0)
89 break;
90 lastSize = newSize;
91
92 result.m_mBases.push_back(give(coarseMBasis));
93 result.m_transferMatrices.push_back(give(transferMatrix));
94 }
95
96 std::reverse( result.m_mBases.begin(), result.m_mBases.end() );
97 std::reverse( result.m_transferMatrices.begin(), result.m_transferMatrices.end() );
98
99 return result;
100 }
101
102 } // namespace gismo
103