forked from johannesgerer/jburkardt-f
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathfem_basis.html
323 lines (280 loc) · 9.34 KB
/
fem_basis.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
<html>
<head>
<title>
FEM_BASIS - Finite Element Basis Functions on an M-dimensional Simplex
</title>
</head>
<body bgcolor="#EEEEEE" link="#CC0000" alink="#FF3300" vlink="#000055">
<h1 align = "center">
FEM_BASIS <br> Finite Element Basis Functions <br> on an M-dimensional Simplex
</h1>
<hr>
<p>
<b>FEM_BASIS</b>
is a FORTRAN90 library which
can evaluate basis functions associated with an M-dimensional simplex
(a 1D interval, a 2D triangle, a 3D tetrahedron, and the higher-dimensional
generalizations).
</p>
<p>
A complete polynomial space is generated, up to the user-specified degree.
Thus, for instance, in 2D, simply by specifying the degree, the user can
request a space of:
<pre>
degree 0, 1 function, 1;
degree 1, 3 functions, 1, x, y;
degree 2, 6 functions, 1, x, y, x^2, xy, y^2;
degree 3, 10 functions, 1, x, y, x^2, xy, y^2, x^3, x^2y, xy^2, y^3
</pre>
and so on. Corresponding families of basis functions in higher dimensions
are generated by specifying the degree.
</p>
<p>
The feature of this library is to generate a Lagrange form of the basis
set of the given degree, and to evaluate those basis functions at any
point.
</p>
<h3 align = "center">
The Triangular Case:
</h3>
<p>
We will discuss the triangular case in detail; it is general enough
to show all the details of the procedure, while still being simple
enough to describe in pictures and tables.
</p>
<p>
Given the maximum degree D for the polynomial basis defined
on a reference triangle, we have ( ( D + 1 ) * ( D + 2 ) ) / 2 monomials
of degree at most D. Normally, the reference triangle is supplied
with barycentric coordinates (xsi1,xsi2,xsi3) which are nonnegative
and add to 1. For our purposes, we will use barycentric coordinates
that have been multiplied by D. This means that all the coordinates
we are interested in will be integers, which we will now identify by
(I,J,K) = D * (xsi1,xsi2,xsi3). We can then define a grid of integer
coordinates so that 0 <= I, J, K <= D and I+J+K = D, with
(I,J,K) corresponding to
<ul>
<li>
the basis point (X,Y)(I,J,K) = ( I/D, J/D );
</li>
<li>
the basis monomial P(I,J,K)(X,Y) = X^I Y^J.
</li>
</ul>
</p>
<p>
For example, with D = 2, we have simply:
<pre>
F
|\
C-E
|\|\
A-B-D
</pre>
with
<pre>
I J K X Y P(I,J,K)(X,Y)
A (0 0 2) (0.0, 0.0) 1
B (1 0 1) (0.5, 0.0) x
C (0 1 1) (0.0, 0.5) y
D (2 0 0) (1.0, 0.0) x^2
E (1 1 0) (0.5, 0.5) x y
F (0 2 0) (0.0, 1.0) y^2
</pre>
Note that in this arrangement, there is an obvious relationship between
a particular choice of (I,J,K) and the corresponding point (X,Y), as well
as a relationship between (I,J,K) and the polynomial P(I,J,K)(X,Y). However,
there is no obvious connection between the point (X,Y) and the polynomial
P(I,J,K)(X,Y). This will change when we move to the Lagrange basis.
</p>
<p>
Instead of the monomials P(I,J,K)(X,Y), we want a set of
polynomials L(I,J,K)(X,Y) which span the same space, but have
the Lagrange property, namely L(I,J,K)(X,Y) is 1 if (X,Y) is
equal to (X,Y)(I,J,K), and 0 if (X,Y) is equal to any other
of the basis points. One effect of this change is to tie each basis point
(X,Y) to a particular basis polynomial L(I,J,K)(X,Y), namely, the polynomial
which is 1 there.
</p>
<p>
This is easily arranged. Given an index (I,J,K), we compute a polynomial
that has:
<ol>
<li>
I factors of the form (X-0) * (X-1/D) * ... * (X-(I-1)/D);
</li>
<li>
J factors of the form (Y-0) * (Y-1/D) * ... * (Y-(J-1)/D);
</li>
<li>
K factors of the form (X+Y-(D-0)/D) * (X+Y-(D-1)/D) * ... * (X+Y-(D-(K-1))/D).
</li>
</ol>
</p>
<p>
This polynomial is the product of I+J+K linear factors, in other words,
it is a polynomial of degree D. This polynomial is 0 at all basis points
except (X,Y)(I,J,K). If we divide this polynomial by its value at
the basis point, we arrive at the desired Lagrange polynomial
L(I,J,K)(X,Y).
</p>
<p>
This method is applicable for a basis set of any polynomial degree;
moreover, it is a straightforward matter to extend it to higher dimensions.
</p>
<h3 align = "center">
The Triangular Prism:
</h3>
<p>
A (right) <i>triangular prism</i> is a 3D shape which can be formed
by drawing a triangle in the XY plane and then "lifting" the triangle
along a line in the Z direction a fixed distance.
</p>
<p>
A finite element basis family can be defined for this shape, using
the product of basis functions in the XY triangle and basis functions
for the Z line. The basis functions are then defined by scaled
barycentric coordinates I(1), I(2), I(3) for the triangle, and
an independent set of scaled barycentric coordinates J(1) and J(2)
for the line. A function is provided to evaluate basis functions
for such an element.
</p>
<h3 align = "center">
Licensing:
</h3>
<p>
The computer code and data files described and made available on this
web page are distributed under
<a href = "../../txt/gnu_lgpl.txt">the GNU LGPL license.</a>
</p>
<h3 align = "center">
Languages:
</h3>
<p>
<b>FEM_BASIS</b> is available in
<a href = "../../cpp_src/fem_basis/fem_basis.html">a C++ version</a> and
<a href = "../../f_src/fem_basis/fem_basis.html">a FORTRAN90 version</a> and
<a href = "../../m_src/fem_basis/fem_basis.html">a MATLAB version</a> and
</p>
<h3 align = "center">
Related Data and Programs:
</h3>
<p>
<a href = "../../m_src/fem_basis_t3_display/fem_basis_t3_display.html">
FEM_BASIS_T3_DISPLAY</a>,
a MATLAB program which
displays a basis function associated with a 3-node triangle "T3" mesh.
</p>
<p>
<a href = "../../m_src/fem_basis_t6_display/fem_basis_t6_display.html">
FEM_BASIS_T6_DISPLAY</a>,
a MATLAB program which
displays a basis function associated with a 6-node triangle "T6" mesh.
</p>
<p>
<a href = "../../f_src/fem2d_pack/fem2d_pack.html">
FEM2D_PACK</a>,
a FORTRAN90 library which
contains utilities for finite element calculations.
</p>
<p>
<a href = "../../f_src/fem3d_pack/fem3d_pack.html">
FEM3D_PACK</a>,
a FORTRAN90 library which
contains utilities for 3D finite element calculations.
</p>
<h3 align = "center">
Reference:
</h3>
<p>
<ol>
<li>
Gilbert Strang, George Fix,<br>
An Analysis of the Finite Element Method,<br>
Cambridge, 1973,<br>
ISBN: 096140888X,<br>
LC: TA335.S77.
</li>
<li>
Olgierd Zienkiewicz,<br>
The Finite Element Method,<br>
Sixth Edition,<br>
Butterworth-Heinemann, 2005,<br>
ISBN: 0750663200,<br>
TA640.2.Z54.
</li>
</ol>
</p>
<h3 align = "center">
Source Code:
</h3>
<p>
<ul>
<li>
<a href = "fem_basis.f90">fem_basis.f90</a>, the source code.
</li>
<li>
<a href = "fem_basis.sh">fem_basis.sh</a>,
BASH commands to compile the source code.
</li>
</ul>
</p>
<h3 align = "center">
Examples and Tests:
</h3>
<p>
<ul>
<li>
<a href = "fem_basis_prb.f90">fem_basis_prb.f90</a>,
a sample calling program.
</li>
<li>
<a href = "fem_basis_prb.sh">fem_basis_prb.sh</a>,
BASH commands to compile and run the sample program.
</li>
<li>
<a href = "fem_basis_prb_output.txt">fem_basis_prb_output.txt</a>,
the output file.
</li>
</ul>
</p>
<h3 align = "center">
List of Routines:
</h3>
<p>
<ul>
<li>
<b>FEM_BASIS_1D</b> evaluates an arbitrary 1D basis function.
</li>
<li>
<b>FEM_BASIS_2D</b> evaluates an arbitrary triangular basis function.
</li>
<li>
<b>FEM_BASIS_3D</b> evaluates an arbitrary tetrahedral basis function.
</li>
<li>
<b>FEM_BASIS_MD</b> evaluates an arbitrary M-dimensional basis function.
</li>
<li>
<b>FEM_BASIS_PRISM_TRIANGLE</b> evaluates a triangular prism basis function.
</li>
<li>
<b>R8_FRACTION</b> uses real arithmetic on an integer ratio.
</li>
<li>
<b>TIMESTAMP</b> prints the current YMDHMS date as a time stamp.
</li>
</ul>
</p>
<p>
You can go up one level to <a href = "../f_src.html">
the FORTRAN90 source codes</a>.
</p>
<hr>
<i>
Last revised on 09 October 2010.
</i>
<!-- John Burkardt -->
</body>
<!-- Initial HTML skeleton created by HTMLINDEX. -->
</html>