BELFEM 0.9.0
Berkeley Lab Finite Element Framework
Loading...
Searching...
No Matches
fn_intpoints_gauss_tet31.hpp
Go to the documentation of this file.
1/*
2 * BELFEM -- The Berkeley Lab Finite Element Framework
3 * Copyright (c) 2026, The Regents of the University of California,
4 * through Lawrence Berkeley National Laboratory (subject to receipt of any required
5 * approvals from the U.S. Dept. of Energy). All rights reserved.
6 *
7 * Developers: Christian Messe, Gregory Giard
8 *
9 * See the top-level LICENSE file for the complete license and disclaimer.
10 */
11
12#ifndef BELFEM_FN_INTPOINTS_GAUSS_TET31_HPP
13#define BELFEM_FN_INTPOINTS_GAUSS_TET31_HPP
14
15#include "typedefs.hpp"
16#include "cl_Vector.hpp"
17#include "cl_Matrix.hpp"
18#include "fn_sum.hpp"
19
20namespace belfem
21{
22 namespace integration
23 {
24// ----------------------------------------------------------------------------
25 inline void
26 gauss_tet31(
27 Vector <real> & aWeights,
28 Matrix <real> & aPoints )
29 {
30
31
32 // source Keast: 10.1016/0045-7825(86)90059-9
33
34 aPoints.set_size( 4, 31 );
35
36 real tAlpha1 = 0.0782131923303186549;
37 real tBeta1 = 1.0 - 3.0*tAlpha1;
38
39 real tAlpha2 = 0.121843216663904411;
40 real tBeta2 = 1.0 - 3.0*tAlpha2;
41
42 real tAlpha3 = 0.332539164446420554;
43 real tBeta3 = 1.0 - 3.0*tAlpha3;
44
45 real tAlpha4 = 0.5;
46 real tBeta4 = 0.0;
47
48 real tAlpha5 = 0.1;
49 real tBeta5 = 0.2;
50 real tDelta5 = 0.6;
51
52 aPoints( 0, 0 ) = 0.25;
53 aPoints( 1, 0 ) = 0.25;
54 aPoints( 2, 0 ) = 0.25;
55
56 aPoints( 0, 1 ) = tAlpha1;
57 aPoints( 1, 1 ) = tAlpha1;
58 aPoints( 2, 1 ) = tAlpha1;
59
60 aPoints( 0, 2 ) = tAlpha1;
61 aPoints( 1, 2 ) = tAlpha1;
62 aPoints( 2, 2 ) = tBeta1;
63
64 aPoints( 0, 3 ) = tAlpha1;
65 aPoints( 1, 3 ) = tBeta1;
66 aPoints( 2, 3 ) = tAlpha1;
67
68 aPoints( 0, 4 ) = tBeta1;
69 aPoints( 1, 4 ) = tAlpha1;
70 aPoints( 2, 4 ) = tAlpha1;
71
72 aPoints( 0, 5 ) = tAlpha2;
73 aPoints( 1, 5 ) = tAlpha2;
74 aPoints( 2, 5 ) = tAlpha2;
75
76 aPoints( 0, 6 ) = tAlpha2;
77 aPoints( 1, 6 ) = tAlpha2;
78 aPoints( 2, 6 ) = tBeta2;
79
80 aPoints( 0, 7 ) = tAlpha2;
81 aPoints( 1, 7 ) = tBeta2;
82 aPoints( 2, 7 ) = tAlpha2;
83
84 aPoints( 0, 8 ) = tBeta2;
85 aPoints( 1, 8 ) = tAlpha2;
86 aPoints( 2, 8 ) = tAlpha2;
87
88 aPoints( 0, 9 ) = tAlpha3;
89 aPoints( 1, 9 ) = tAlpha3;
90 aPoints( 2, 9 ) = tAlpha3;
91
92 aPoints( 0, 10 ) = tAlpha3;
93 aPoints( 1, 10 ) = tAlpha3;
94 aPoints( 2, 10 ) = tBeta3;
95
96 aPoints( 0, 11 ) = tAlpha3;
97 aPoints( 1, 11 ) = tBeta3;
98 aPoints( 2, 11 ) = tAlpha3;
99
100 aPoints( 0, 12 ) = tBeta3;
101 aPoints( 1, 12 ) = tAlpha3;
102 aPoints( 2, 12 ) = tAlpha3;
103
104 aPoints( 0, 13 ) = tAlpha4;
105 aPoints( 1, 13 ) = tAlpha4;
106 aPoints( 2, 13 ) = tBeta4;
107
108 aPoints( 0, 14 ) = tAlpha4;
109 aPoints( 1, 14 ) = tBeta4;
110 aPoints( 2, 14 ) = tAlpha4;
111
112 aPoints( 0, 15 ) = tAlpha4;
113 aPoints( 1, 15 ) = tBeta4;
114 aPoints( 2, 15 ) = tBeta4;
115
116 aPoints( 0, 16 ) = tBeta4;
117 aPoints( 1, 16 ) = tAlpha4;
118 aPoints( 2, 16 ) = tAlpha4;
119
120 aPoints( 0, 17 ) = tBeta4;
121 aPoints( 1, 17 ) = tAlpha4;
122 aPoints( 2, 17 ) = tBeta4;
123
124 aPoints( 0, 18 ) = tBeta4;
125 aPoints( 1, 18 ) = tBeta4;
126 aPoints( 2, 18 ) = tAlpha4;
127
128 aPoints( 0, 19 ) = tAlpha5;
129 aPoints( 1, 19 ) = tAlpha5;
130 aPoints( 2, 19 ) = tBeta5;
131
132 aPoints( 0, 20 ) = tAlpha5;
133 aPoints( 1, 20 ) = tAlpha5;
134 aPoints( 2, 20 ) = tDelta5;
135
136 aPoints( 0, 21 ) = tAlpha5;
137 aPoints( 1, 21 ) = tBeta5;
138 aPoints( 2, 21 ) = tAlpha5;
139
140 aPoints( 0, 22 ) = tAlpha5;
141 aPoints( 1, 22 ) = tBeta5;
142 aPoints( 2, 22 ) = tDelta5;
143
144 aPoints( 0, 23 ) = tAlpha5;
145 aPoints( 1, 23 ) = tDelta5;
146 aPoints( 2, 23 ) = tAlpha5;
147
148 aPoints( 0, 24 ) = tAlpha5;
149 aPoints( 1, 24 ) = tDelta5;
150 aPoints( 2, 24 ) = tBeta5;
151
152 aPoints( 0, 25 ) = tBeta5;
153 aPoints( 1, 25 ) = tAlpha5;
154 aPoints( 2, 25 ) = tAlpha5;
155
156 aPoints( 0, 26 ) = tBeta5;
157 aPoints( 1, 26 ) = tAlpha5;
158 aPoints( 2, 26 ) = tDelta5;
159
160 aPoints( 0, 27 ) = tBeta5;
161 aPoints( 1, 27 ) = tDelta5;
162 aPoints( 2, 27 ) = tAlpha5;
163
164 aPoints( 0, 28 ) = tDelta5;
165 aPoints( 1, 28 ) = tAlpha5;
166 aPoints( 2, 28 ) = tAlpha5;
167
168 aPoints( 0, 29 ) = tDelta5;
169 aPoints( 1, 29 ) = tAlpha5;
170 aPoints( 2, 29 ) = tBeta5;
171
172 aPoints( 0, 30 ) = tDelta5;
173 aPoints( 1, 30 ) = tBeta5;
174 aPoints( 2, 30 ) = tAlpha5;
175
176 for( uint k=0; k<31; ++k )
177 {
178 aPoints( 3, k ) = 1.0
179 - aPoints( 0, k )
180 - aPoints( 1, k )
181 - aPoints( 2, k );
182 }
183
184 aWeights.set_size( 31 );
185
186 aWeights( 0 ) = 0.0182642234661087939;
187
188 aWeights( 1 ) = 0.0105999415244141609;
189 aWeights( 2 ) = 0.0105999415244141609;
190 aWeights( 3 ) = 0.0105999415244141609;
191 aWeights( 4 ) = 0.0105999415244141609;
192
193 aWeights( 5 ) = -0.0625177401143299494;
194 aWeights( 6 ) = -0.0625177401143299494;
195 aWeights( 7 ) = -0.0625177401143299494;
196 aWeights( 8 ) = -0.0625177401143299494;
197
198 aWeights( 9 ) = 0.00489142526307353653;
199 aWeights( 10 ) = 0.00489142526307353653;
200 aWeights( 11 ) = 0.00489142526307353653;
201 aWeights( 12 ) = 0.00489142526307353653;
202
203 aWeights( 13 ) = 0.000970017636684296702;
204 aWeights( 14 ) = 0.000970017636684296702;
205 aWeights( 15 ) = 0.000970017636684296702;
206 aWeights( 16 ) = 0.000970017636684296702;
207 aWeights( 17 ) = 0.000970017636684296702;
208 aWeights( 18 ) = 0.000970017636684296702;
209
210 aWeights( 19 ) = 0.0275573192239850917;
211 aWeights( 20 ) = 0.0275573192239850917;
212 aWeights( 21 ) = 0.0275573192239850917;
213 aWeights( 22 ) = 0.0275573192239850917;
214 aWeights( 23 ) = 0.0275573192239850917;
215 aWeights( 24 ) = 0.0275573192239850917;
216
217 aWeights( 25 ) = 0.0275573192239850917;
218 aWeights( 26 ) = 0.0275573192239850917;
219 aWeights( 27 ) = 0.0275573192239850917;
220 aWeights( 28 ) = 0.0275573192239850917;
221 aWeights( 29 ) = 0.0275573192239850917;
222 aWeights( 30 ) = 0.0275573192239850917;
223
224 aWeights( 0 ) = 0.0;
225 aWeights( 0 ) = 1./6. - sum( aWeights );
226 }
227
228// ----------------------------------------------------------------------------
229 } /* namespace integration */
230} /* namespace belfem */
231
232#endif //BELFEM_FN_INTPOINTS_GAUSS_TET31_HPP
cl_Matrix.hpp
cl_Vector.hpp
belfem::Matrix::set_size
void set_size(const size_t aNumRows, const size_t aNumCols)
Definition cl_AR_Matrix.hpp:186
belfem::Vector::set_size
void set_size(const size_t aNumRows)
change the size of the vector
Definition cl_AR_Vector.hpp:237
fn_sum.hpp
Sum of the entries of a vector.
belfem::integration
Definition fn_intpoints.cpp:576
belfem::integration::gauss_tet31
void gauss_tet31(Vector< real > &aWeights, Matrix< real > &aPoints)
Definition fn_intpoints_gauss_tet31.hpp:26
belfem
USER GUIDES:
Definition cl_Capacitor.cpp:16
belfem::uint
unsigned int uint
Definition typedefs.hpp:30
belfem::sum
auto sum(const Vector< T > &aA) -> decltype(sum(aA.vector_data()))
Definition fn_sum.hpp:49
belfem::real
double real
Definition typedefs.hpp:36