1 // GetDP - Copyright (C) 1997-2021 P. Dular and C. Geuzaine, University of Liege
2 //
3 // See the LICENSE.txt file for license information. Please report all
4 // issues on https://gitlab.onelab.info/getdp/getdp/issues.
5
6 #include "GetDPConfig.h"
7 #include "ProData.h"
8 #include "BF.h"
9 #include "MallocUtils.h"
10 #include "Message.h"
11 #include "Cal_Quantity.h"
12
13 #if defined(HAVE_KERNEL)
14 #include "Get_DofOfElement.h"
15 #include "Pos_FemInterpolation.h"
16 #endif
17
18 extern struct Problem Problem_S ;
19 extern struct CurrentData Current ;
20
21 /* ------------------------------------------------------------------------ */
22 /* B F _ S u b F u n c t i o n */
23 /* ------------------------------------------------------------------------ */
24
BF_SubFunction(struct Element * Element,int NumExpression,int Dim,double s[])25 void BF_SubFunction(struct Element * Element, int NumExpression,
26 int Dim, double s[])
27 {
28 struct Value Value ;
29
30 Get_ValueOfExpressionByIndex(NumExpression, NULL, 0., 0., 0., &Value) ;
31
32 switch (Dim) {
33 case 1 :
34 *s *= Value.Val[0] ;
35 break ;
36 case 3 :
37 s[0] *= Value.Val[0] ;
38 s[1] *= Value.Val[0] ;
39 s[2] *= Value.Val[0] ;
40 break ;
41 }
42 }
43
44 /* ------------------------------------------------------------------------ */
45 /* B F _ R e g i o n */
46 /* ------------------------------------------------------------------------ */
47
BF_Region(struct Element * Element,int NumRegion,double u,double v,double w,double * s)48 void BF_Region(struct Element * Element, int NumRegion,
49 double u, double v, double w, double *s)
50 {
51 *s = 1. ;
52
53 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
54 BF_SubFunction(Element, Element->NumSubFunction[0][NumRegion-1], 1, s) ;
55 }
56
BF_dRegion(struct Element * Element,int NumRegion,double u,double v,double w,double * s)57 void BF_dRegion(struct Element * Element, int NumRegion,
58 double u, double v, double w, double *s)
59 {
60 *s = 1. ;
61
62 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
63 BF_SubFunction(Element, Element->NumSubFunction[2][NumRegion-1], 1, s) ;
64 else
65 *s = 0. ;
66 }
67
68 /* ------------------------------------------------------------------------ */
69 /* B F _ R e g i o n X , Y , Z */
70 /* ------------------------------------------------------------------------ */
71
BF_RegionX(struct Element * Element,int NumRegion,double u,double v,double w,double s[])72 void BF_RegionX(struct Element * Element, int NumRegion,
73 double u, double v, double w, double s[])
74 {
75 s[1] = s[2] = 0. ; s[0] = 1. ;
76
77 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
78 BF_SubFunction(Element, Element->NumSubFunction[0][NumRegion-1], 3, s) ;
79 }
80
BF_RegionY(struct Element * Element,int NumRegion,double u,double v,double w,double s[])81 void BF_RegionY(struct Element * Element, int NumRegion,
82 double u, double v, double w, double s[])
83 {
84
85 s[0] = s[2] = 0. ; s[1] = 1. ;
86
87 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
88 BF_SubFunction(Element, Element->NumSubFunction[0][NumRegion-1], 3, s) ;
89 }
90
BF_RegionZ(struct Element * Element,int NumRegion,double u,double v,double w,double s[])91 void BF_RegionZ(struct Element * Element, int NumRegion,
92 double u, double v, double w, double s[])
93 {
94 s[0] = s[1] = 0. ; s[2] = 1. ;
95
96 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
97 BF_SubFunction(Element, Element->NumSubFunction[0][NumRegion-1], 3, s) ;
98 }
99
BF_dRegionX(struct Element * Element,int NumRegion,double u,double v,double w,double s[])100 void BF_dRegionX(struct Element * Element, int NumRegion,
101 double u, double v, double w, double s[])
102 {
103 s[1] = s[2] = 0. ; s[0] = 1. ; /* Patrick (a finaliser) */
104
105 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
106 BF_SubFunction(Element, Element->NumSubFunction[2][NumRegion-1], 3, s) ;
107 else
108 s[0] = 0. ;
109 }
110
BF_dRegionY(struct Element * Element,int NumRegion,double u,double v,double w,double s[])111 void BF_dRegionY(struct Element * Element, int NumRegion,
112 double u, double v, double w, double s[])
113 {
114 s[0] = s[2] = 0. ; s[1] = 1. ; /* Patrick (a finaliser) */
115
116 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
117 BF_SubFunction(Element, Element->NumSubFunction[2][NumRegion-1], 3, s) ;
118 else
119 s[1] = 0. ;
120 }
121
BF_dRegionZ(struct Element * Element,int NumRegion,double u,double v,double w,double s[])122 void BF_dRegionZ(struct Element * Element, int NumRegion,
123 double u, double v, double w, double s[])
124 {
125 /* Patrick (a finaliser)
126 s[0] = s[1] = 0. ; s[2] = 1. ;
127 */
128 s[0] = s[2] = 0. ; s[1] = -1. ;
129
130 if (Element->NumSubFunction[0][NumRegion-1] >= 0)
131 BF_SubFunction(Element, Element->NumSubFunction[2][NumRegion-1], 3, s) ;
132 else
133 s[1] = 0. ;
134 }
135
136 /* ------------------------------------------------------------------------ */
137 /* B F _ Z e r o */
138 /* ------------------------------------------------------------------------ */
139
BF_Zero(struct Element * Element,int Num,double u,double v,double w,double * s)140 void BF_Zero(struct Element * Element, int Num,
141 double u, double v, double w, double *s)
142 {
143 s[0] = s[1] = s[2] = 0. ;
144 }
145
BF_One(struct Element * Element,int Num,double u,double v,double w,double * s)146 void BF_One(struct Element * Element, int Num,
147 double u, double v, double w, double *s)
148 {
149 s[0] = 1. ; s[1] = s[2] = 0. ;
150 }
151
BF_OneZ(struct Element * Element,int Num,double u,double v,double w,double * s)152 void BF_OneZ(struct Element * Element, int Num,
153 double u, double v, double w, double *s)
154 {
155 s[0] = s[1] = 0. ; s[2] = 1. ;
156 }
157
158 /* ------------------------------------------------------------------------ */
159 /* B F _ I n i t G l o b a l */
160 /* ------------------------------------------------------------------------ */
161
BF_InitGlobal(struct GlobalBasisFunction * GlobalBasisFunction_P)162 void BF_InitGlobal(struct GlobalBasisFunction * GlobalBasisFunction_P)
163 {
164 struct QuantityStorage * QuantityStorage_P ;
165 struct Formulation * Formulation_P ;
166
167 QuantityStorage_P =
168 GlobalBasisFunction_P->QuantityStorage =
169 (struct QuantityStorage *)Malloc(sizeof(struct QuantityStorage)) ;
170
171 QuantityStorage_P->NumLastElementForFunctionSpace = 0 ;
172
173 Formulation_P = (struct Formulation*)
174 List_Pointer(Problem_S.Formulation,
175 GlobalBasisFunction_P->FormulationIndex) ;
176 QuantityStorage_P->DefineQuantity = (struct DefineQuantity*)
177 List_Pointer(Formulation_P->DefineQuantity,
178 GlobalBasisFunction_P->DefineQuantityIndex) ;
179 QuantityStorage_P->FunctionSpace = (struct FunctionSpace*)
180 List_Pointer(Problem_S.FunctionSpace,
181 QuantityStorage_P->DefineQuantity->FunctionSpaceIndex) ;
182 QuantityStorage_P->TypeQuantity = QuantityStorage_P->FunctionSpace->Type ;
183 }
184
185 /* ------------------------------------------------------------------------ */
186 /* B F _ G l o b a l , B F _ d G l o b a l , B F _ d I n v G l o b a l */
187 /* ------------------------------------------------------------------------ */
188
BF_Global_OP(struct Element * Element,int NumGlobal,double u,double v,double w,double * s,int type_OP)189 void BF_Global_OP(struct Element * Element, int NumGlobal,
190 double u, double v, double w, double *s, int type_OP)
191 {
192 #if !defined(HAVE_KERNEL)
193 Message::Error("BF_Global_OP requires Kernel");
194 #else
195 struct Value Value ;
196 struct GlobalBasisFunction * GlobalBasisFunction_P ;
197 struct QuantityStorage * QuantityStorage_P ;
198 int Save_NbrHar;
199
200 GlobalBasisFunction_P = Element->GlobalBasisFunction[NumGlobal-1] ;
201
202 if (!GlobalBasisFunction_P->QuantityStorage)
203 BF_InitGlobal(GlobalBasisFunction_P) ; /* Init QuantityStorage */
204
205 QuantityStorage_P = GlobalBasisFunction_P->QuantityStorage ;
206
207 if (QuantityStorage_P->NumLastElementForFunctionSpace != Element->Num) {
208 QuantityStorage_P->NumLastElementForFunctionSpace = Element->Num ;
209
210 Get_DofOfElement
211 (Element, QuantityStorage_P->FunctionSpace, QuantityStorage_P,
212 QuantityStorage_P->DefineQuantity->IndexInFunctionSpace) ;
213 }
214
215 Save_NbrHar = Current.NbrHar;
216 Current.NbrHar = 1; /* for real basis function */
217 Pos_FemInterpolation
218 (Element,
219 NULL,
220 GlobalBasisFunction_P->QuantityStorage,
221 QUANTITY_SIMPLE, type_OP, -1, 0,
222 u, v, w, 0., 0., 0., Value.Val, &Value.Type, 0) ;
223 Current.NbrHar = Save_NbrHar;
224
225 switch (Value.Type) {
226 case SCALAR :
227 s[0] = Value.Val[0] ;
228 break ;
229 case VECTOR :
230 s[0] = Value.Val[0] ; s[1] = Value.Val[1] ; s[2] = Value.Val[2] ;
231 break ;
232 default :
233 Message::Error("Bad type of value for Global BasisFunction") ;
234 }
235 #endif
236 }
237
238 /* ------------------------------------------------------------------------ */
239
BF_Global(struct Element * Element,int NumGlobal,double u,double v,double w,double * s)240 void BF_Global(struct Element * Element, int NumGlobal,
241 double u, double v, double w, double *s)
242 {
243 BF_Global_OP(Element, NumGlobal, u, v, w, s, NOOP);
244 }
245
BF_dGlobal(struct Element * Element,int NumGlobal,double u,double v,double w,double * s)246 void BF_dGlobal(struct Element * Element, int NumGlobal,
247 double u, double v, double w, double *s )
248 {
249 BF_Global_OP(Element, NumGlobal, u, v, w, s, EXTDER);
250 }
251
BF_dInvGlobal(struct Element * Element,int NumGlobal,double u,double v,double w,double * s)252 void BF_dInvGlobal(struct Element * Element, int NumGlobal,
253 double u, double v, double w, double *s )
254 {
255 BF_Global_OP(Element, NumGlobal, u, v, w, s, EXTDERINV);
256 }
257