-
Notifications
You must be signed in to change notification settings - Fork 2
/
VFCartFE.h
205 lines (180 loc) · 9.6 KB
/
VFCartFE.h
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
#ifndef VFVFCartFEH
#define VFVFCartFEH
/*
Implement a simple P1 Lagrange element on a rectangle aligned with axis directions.
(c) 2010-2012 B. Bourdin, [email protected]
Uses 3 points / 5th order quadrature from
A. Ern and J.-L. Guermond "Theory and Practice of Finite Elements"
Table 8.1 p. 359
*/
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt ng; /* number of integration points */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal lx; /* length of the element */
PetscReal ly; /* length of the element */
PetscReal lz; /* length of the element */
/*PetscReal *weight;*/ /* integration weight */
PetscReal weight[3]; /* integration weight */
PetscReal phi[1][1][2][3]; /* phi[i][g] = value of i^th basis function at g^th integration point */
PetscReal dphi[1][1][2][3]; /* dphi[i][g] = value of the derivative of the i^th basis function at g^th integration point */
} VFCartFEElement1D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt ng; /* number of integration points */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal lx; /* length of the element */
PetscReal ly; /* length of the element */
PetscReal lz; /* length of the element */
/*PetscReal *weight;*/ /* integration weight */
PetscReal weight[9]; /* integration weight */
PetscReal phi[1][2][2][9]; /* phi[j][i][g] = value of (i,j)^th basis function at g^th integration point */
PetscReal dphi[1][2][2][2][9]; /* phi[j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j)^th basis function at g^th integration point */
} VFCartFEElement2D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt ng; /* number of integration points */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal lx; /* length of the element */
PetscReal ly; /* length of the element */
PetscReal lz; /* length of the element */
/*PetscReal *weight;*/ /* integration weight */
PetscReal weight[9]; /* integration weight */
PetscReal phi[1][2][2][9]; /* phi[j][i][g] = value of (i,j)^th basis function at g^th integration point */
PetscReal dphi[1][2][2][2][9]; /* phi[j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j)^th basis function at g^th integration point */
} VFCartFEElement22D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt ng; /* number of integration points */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal lx; /* length of the element */
PetscReal ly; /* length of the element */
PetscReal lz; /* length of the element */
/*PetscReal *weight;*/ /* integration weight */
PetscReal weight[27]; /* integration weight */
PetscReal phi[2][2][2][27]; /* phi[k][j][i][g] = value of (i,j,k)^th basis function at g^th integration point */
PetscReal dphi[2][2][2][3][27]; /* phi[k][j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j,k)^th basis function at g^th integration point */
} VFCartFEElement3D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt ng; /* number of integration points */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal lx; /* length of the element */
PetscReal ly; /* length of the element */
PetscReal lz; /* length of the element */
PetscReal *weight; /* integration weight */
PetscReal ****phi; /* phi[k][j][i][g] = value of (i,j,k)^th basis function at g^th integration point */
PetscReal *****dphi; /* phi[k][j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j,k)^th basis function at g^th integration point */
} CartesianElement;
typedef struct {
PetscReal x;
PetscReal y;
PetscReal z;
} coord3d;
typedef enum {
NONE,
ZERO,
ONE,
FIXED,
} BCTYPE;
typedef enum {
X0,
X1,
Y0,
Y1,
Z0,
Z1
} FACE;
typedef enum {
X0Z0,
X1Z0,
Y0Z0,
Y1Z0,
X0Z1,
X1Z1,
Y0Z1,
Y1Z1,
X0Y0,
X0Y1,
X1Y0,
X1Y1
} EDGE;
typedef enum {
X0Y0Z0,
X1Y0Z0,
X0Y1Z0,
X1Y1Z0,
X0Y0Z1,
X1Y0Z1,
X0Y1Z1,
X1Y1Z1
} VERTEX;
typedef struct {
BCTYPE face[6];
PetscReal faceValue[6];
BCTYPE edge[12];
PetscReal edgeValue[12];
BCTYPE vertex[8];
PetscReal vertexValue[8];
} VFBC;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal phi[1][1][2]; /* phi[i][g] = value of i^th basis function at g^th integration point */
PetscReal dphi[1][1][2]; /* dphi[i][g] = value of the derivative of the i^th basis function at g^th integration point */
} CartFEElement1D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal phi[1][2][2]; /* phi[j][i][g] = value of (i,j)^th basis function at g^th integration point */
PetscReal dphi[1][2][2][2]; /* phi[j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j)^th basis function at g^th integration point */
} CartFEElement2D;
typedef struct {
PetscInt dim; /* dimension of the space */
PetscInt nphix; /* number of basis functions along the x axis */
PetscInt nphiy; /* number of basis functions along the y axis */
PetscInt nphiz; /* number of basis functions along the z axis */
PetscReal phi[2][2][2]; /* phi[k][j][i][g] = value of (i,j,k)^th basis function at g^th integration point */
PetscReal dphi[2][2][2][3]; /* phi[k][j][i][l][g] = value of the derivative w.r.t. x_l of the (i,j,k)^th basis function at g^th integration point */
} CartFEElement3D;
#ifndef VFCartFEC
extern PetscErrorCode VFCartFEInit();
extern PetscErrorCode VFCartFEElement1DCreate1D(VFCartFEElement1D *e);
extern PetscErrorCode VFCartFEElement1DInit(VFCartFEElement1D *e,PetscReal lx);
extern PetscErrorCode VFCartFEElement2DCreate(VFCartFEElement2D *e);
extern PetscErrorCode VFCartFEElement2DInit(VFCartFEElement2D *e,PetscReal lx,PetscReal ly);
extern PetscErrorCode VFCartFEElement3DCreate(VFCartFEElement3D *e);
extern PetscErrorCode VFCartFEElement3DInit(VFCartFEElement3D *e,PetscReal lx,PetscReal ly,PetscReal lz);
extern PetscErrorCode DAReadCoordinatesHDF5(DM da,const char filename[]);
extern PetscErrorCode VFBCCreate(VFBC *bc,PetscInt dof);
extern PetscErrorCode VFBCSetFromOptions(VFBC *bc,const char prefix[],PetscInt dof);
extern PetscErrorCode VFBCView(VFBC *bc,PetscViewer viewer,PetscInt dof);
extern PetscErrorCode VecSetFromBC(Vec BCVec,VFBC *BC);
extern PetscErrorCode VecApplyDirichletBC(Vec RHS,Vec BCU,VFBC *BC);
extern PetscErrorCode ResidualApplyDirichletBC(Vec Residual,Vec U,Vec BCU,VFBC *BC);
extern PetscErrorCode GradientApplyDirichletBC(Vec gradient,VFBC *BC);
extern PetscErrorCode MatApplyDirichletBC(Mat K,VFBC *BC);
extern PetscErrorCode MatApplyDirichletBCRowCol(Mat K,VFBC *BC);
extern PetscErrorCode VecApplyDirichletFlowBC(Vec RHS,Vec BCU,VFBC *BC,PetscReal *BCpres);
extern PetscErrorCode CartFEElement1DCreate(CartFEElement1D *e);
extern PetscErrorCode CartFEElement1DInit(CartFEElement1D *e,PetscReal x,PetscReal lx);
extern PetscErrorCode CartFEElement2DCreate(CartFEElement2D *e);
extern PetscErrorCode CartFEElement2DInit(CartFEElement2D *e,PetscReal x,PetscReal y,PetscReal lx,PetscReal ly);
extern PetscErrorCode CartFEElement3DCreate(CartFEElement3D *e);
extern PetscErrorCode CartFEElement3DInit(CartFEElement3D *e,PetscReal x,PetscReal y,PetscReal z,PetscReal lx,PetscReal ly,PetscReal lz);
#endif
#endif /* VFCartFEH */