Choral
fe_func_mod.f90
Go to the documentation of this file.
1 
31 
33 
34  use real_type
35 
36  public
37 
38 
42 
44  procedure(fescal), nopass, pointer :: u =>null()
45 
47  procedure(fescalgrad), nopass, pointer :: grad_u =>null()
48 
50  procedure(fescalgrad), nopass, pointer :: phi =>null()
51 
53  procedure(fescal), nopass, pointer :: div_phi =>null()
54 
55  end type fe_functions
56 
57 
58  ! %----------------------------------------%
59  ! | |
60  ! | INTERFACE FOR FE FUNCTIONS |
61  ! | |
62  ! %----------------------------------------%
63  abstract interface
64 
66  subroutine fescal(val, n, x, m)
67  import :: rp
68  real(RP), dimension(n), intent(out) :: val
69  real(RP), dimension(m), intent(in) :: x
70  integer , intent(in) :: n, m
71  end subroutine fescal
72 
74  subroutine fescalgrad(val, n, x, m)
75  import :: rp
76  real(RP), dimension(m,n), intent(out) :: val
77  real(RP), dimension(m) , intent(in) :: x
78  integer , intent(in) :: n, m
79  end subroutine fescalgrad
80 
81  end interface
82 
83 contains
84 
86  subroutine p1_1d_u(val, n, x, m)
87  real(RP), dimension(n), intent(out) :: val
88  real(RP), dimension(m), intent(in) :: x
89  integer , intent(in) :: n, m
90 
91  val(1) = 1._rp - x(1)
92  val(2) = x(1)
93 
94  end subroutine p1_1d_u
95 
97  subroutine p2_1d_u(val, n, x, m)
98  real(RP), dimension(n), intent(out) :: val
99  real(RP), dimension(m), intent(in) :: x
100  integer , intent(in) :: n, m
101 
102  val(1) = 2._rp * (x(1) - 1._rp )*( x(1) - 0.5_rp )
103  val(2) = 2._rp * x(1) * ( x(1) - 0.5_rp )
104  val(3) =-4._rp * x(1) * (x(1) - 1._rp )
105 
106  end subroutine p2_1d_u
107 
108 
110  subroutine p3_1d_u(val, n, x, m)
111  real(RP), dimension(n), intent(out) :: val
112  real(RP), dimension(m), intent(in) :: x
113  integer , intent(in) :: n, m
114 
115  real(RP) :: z, v2, v3, v4
116 
117  v2 = 1._rp / 3._rp
118  v3 = 2._rp / 3._rp
119  v4 = 1._rp
120 
121  z = x(1)
122  val(1) =-(z-v2)*(z-v3)*(z-v4) * 4.5_rp
123  val(2) = z*(z-v2)*(z-v3) * 4.5_rp
124  val(3) = z*(z-v3)*(z-v4) * 13.5_rp
125  val(4) =-z*(z-v2)*(z-v4) * 13.5_rp
126 
127  end subroutine p3_1d_u
128 
129 
131  subroutine p0_1d_u(val, n, x, m)
132  real(RP), dimension(n), intent(out) :: val
133  real(RP), dimension(m), intent(in) :: x
134  integer , intent(in) :: n, m
135 
136  val(1) = 1._rp
137 
138  end subroutine p0_1d_u
139 
141  subroutine p0_2d_u(val, n, x, m)
142  real(RP), dimension(n), intent(out) :: val
143  real(RP), dimension(m), intent(in) :: x
144  integer , intent(in) :: n, m
145 
146  val(1) = 1._rp
147 
148  end subroutine p0_2d_u
149 
151  subroutine p1_2d_u(val, n, x, m)
152  real(RP), dimension(n), intent(out) :: val
153  real(RP), dimension(m), intent(in) :: x
154  integer , intent(in) :: n, m
155 
156  val(1) = 1._rp - x(1) - x(2)
157  val(2) = x(1)
158  val(3) = x(2)
159 
160  end subroutine p1_2d_u
161 
162 
164  subroutine p1_1d_disc_ortho_u(val, n, x, m)
165  real(RP), dimension(n), intent(out) :: val
166  real(RP), dimension(m), intent(in) :: x
167  integer , intent(in) :: n, m
168 
169  val(1) =-1.7320508075688774_rp * ( x(1) - 0.788675134594813_rp )
170  val(2) = 1.7320508075688774_rp * ( x(1) - 0.21132486540518702_rp )
171 
172  end subroutine p1_1d_disc_ortho_u
173 
175  subroutine p1_2d_disc_ortho_u(val, n, x, m)
176  real(RP), dimension(n), intent(out) :: val
177  real(RP), dimension(m), intent(in) :: x
178  integer , intent(in) :: n, m
179 
180  val(1) = re(5,3) - 2._rp*x(1) - 2._rp*x(2)
181  val(2) =-re(1,3) + 2._rp*x(1)
182  val(3) =-re(1,3) + 2._rp*x(2)
183 
184  end subroutine p1_2d_disc_ortho_u
185 
186 
188  subroutine p2_2d_u(val, n, x, m)
189  real(RP), dimension(n), intent(out) :: val
190  real(RP), dimension(m), intent(in) :: x
191  integer , intent(in) :: n, m
192 
193  real(RP) :: a
194 
195  a = 1._rp - x(1) - x(2)
196 
197  val(1) = a * (1._rp - 2._rp*x(1) - 2._rp*x(2))
198  val(2) = x(1)* ( -1._rp + 2._rp*x(1) )
199  val(3) = x(2)* ( -1._rp + 2._rp*x(2) )
200 
201  val(4) = 4._rp * x(1) * a
202  val(5) = 4._rp * x(1) * x(2)
203  val(6) = 4._rp * x(2) * a
204 
205  end subroutine p2_2d_u
206 
208  subroutine p1_3d_u(val, n, x, m)
209  real(RP), dimension(n), intent(out) :: val
210  real(RP), dimension(m), intent(in) :: x
211  integer , intent(in) :: n, m
212 
213  val(1) = 1._rp - x(1) - x(2) - x(3)
214  val(2) = x(1)
215  val(3) = x(2)
216  val(4) = x(3)
217 
218  end subroutine p1_3d_u
219 
220 
222  subroutine p3_2d_u(val, n, x, m)
223  real(RP), dimension(n), intent(out) :: val
224  real(RP), dimension(m), intent(in) :: x
225  integer , intent(in) :: n, m
226 
227  real(RP) :: a, b, c, d, e, f, g
228 
229  a = x(1) + x(2) - 1._rp/3._rp
230  b = x(1) + x(2) - 2._rp/3._rp
231  c = x(1) + x(2)- 1._rp
232  d = x(1) - 1._rp/3._rp
233  e = x(1) - 2._rp/3._rp
234  f = x(2) - 1._rp/3._rp
235  g = x(2) - 2._rp/3._rp
236 
237  val(1) = -a*b*c*9._rp/2._rp
238  val(2) = x(1)*d*e*9._rp/2._rp
239  val(3) = x(2)*f*g*9._rp/2._rp
240  val(4) = x(1)*b*c*27._rp/2._rp
241  val(5) = -x(1)*d*c*27._rp/2._rp
242  val(6) = x(1)*x(2)*d*27._rp/2._rp
243  val(7) = x(1)*x(2)*f*27._rp/2._rp
244  val(8) = -x(2)*f*c*27._rp/2._rp
245  val(9) = x(2)*b*c*27._rp/2._rp
246  val(10) = -x(1)*x(2)*c*27._rp
247 
248  end subroutine p3_2d_u
249 
250 
252  subroutine p2_3d_u(val, n, x, m)
253  real(RP), dimension(n), intent(out) :: val
254  real(RP), dimension(m), intent(in) :: x
255  integer , intent(in) :: n, m
256 
257  val(1) = 1._rp - 3._rp*x(1) - 3._rp*x(2) - 3._rp*x(3) &
258  & + 2._rp*x(1)**2 + 2._rp*x(2)**2 + 2._rp*x(3)**2 &
259  & + 4._rp*x(2)*x(3) + 4._rp*x(1)*x(3) + 4._rp*x(1)*x(2)
260  val(2) = -1._rp*x(1)+2._rp*x(1)**2
261  val(3) = -1._rp*x(2)+2._rp*x(2)**2
262  val(4) = -1._rp*x(3)+2._rp*x(3)**2
263  val(5) = 4._rp*x(1)-4._rp*x(1)**2-4._rp*x(1)*x(3)-4._rp*x(1)*x(2)
264  val(6) = 4._rp*x(1)*x(2)
265  val(7) = 4._rp*x(2)-4._rp*x(2)**2-4._rp*x(2)*x(3)-4._rp*x(1)*x(2)
266  val(8) = 4._rp*x(3)-4._rp*x(3)**2-4._rp*x(2)*x(3)-4._rp*x(1)*x(3)
267  val(9) = 4._rp*x(1)*x(3)
268  val(10)= 4._rp*x(2)*x(3)
269 
270  end subroutine p2_3d_u
271 
272 
273 
274  ! %----------------------------------------%
275  ! | |
276  ! | FINITE ELEMENT BASE FUNCTIONS !
277  ! | SCALAR CASE : GRADIENT |
278  ! | |
279  ! %----------------------------------------%
281  subroutine p1_1d_grad_u(val, n, x, m)
282  real(RP), dimension(m,n), intent(out) :: val
283  real(RP), dimension(m) , intent(in) :: x
284  integer , intent(in) :: n, m
285 
286  val(1,1) = -1._rp
287  val(1,2) = 1._rp
288 
289  end subroutine p1_1d_grad_u
290 
291 
292 
294  subroutine p2_1d_grad_u(val, n, x, m)
295  real(RP), dimension(m,n), intent(out) :: val
296  real(RP), dimension(m) , intent(in) :: x
297  integer , intent(in) :: n, m
298 
299  val(1,1) = 4._rp*x(1) - 3._rp
300  val(1,2) = 4._rp*x(1) - 1._rp
301  val(1,3) =-8._rp*x(1) + 4._rp
302 
303  end subroutine p2_1d_grad_u
304 
305 
307  subroutine p3_1d_grad_u(val, n, x, m)
308  real(RP), dimension(m,n), intent(out) :: val
309  real(RP), dimension(m) , intent(in) :: x
310  integer , intent(in) :: n, m
311 
312  real(RP) :: z
313 
314  z = x(1)
315 
316  val(1,1) = -13.5_rp * z**2 + 18._rp * z - 5.5_rp
317  val(1,2) = 13.5_rp * z**2 - 9._rp * z + 1._rp
318  val(1,3) = 40.5_rp * z**2 - 45._rp * z + 9._rp
319  val(1,4) = -40.5_rp * z**2 + 36._rp * z - 4.5_rp
320 
321  end subroutine p3_1d_grad_u
322 
323 
324 
326  subroutine p1_2d_grad_u(val, n, x, m)
327  real(RP), dimension(m,n), intent(out) :: val
328  real(RP), dimension(m) , intent(in) :: x
329  integer , intent(in) :: n, m
330 
331  val(:,1) = (/-1._rp,-1._rp/)
332 
333  val(:,2) = (/ 1._rp, 0._rp/)
334 
335  val(:,3) = (/ 0._rp, 1._rp/)
336 
337  end subroutine p1_2d_grad_u
338 
340  subroutine p2_2d_grad_u(val, n, x, m)
341  real(RP), dimension(m,n), intent(out) :: val
342  real(RP), dimension(m) , intent(in) :: x
343  integer , intent(in) :: n, m
344 
345  real(RP) :: a
346 
347 
348  a = -3._rp + 4._rp*x(1) + 4._rp*x(2)
349  val(:,1) = (/a, a/)
350 
351  val(:,2) = (/ -1._rp + 4._rp*x(1), 0._rp/)
352 
353  val(:,3) = (/ 0._rp, -1._rp + 4._rp*x(2) /)
354 
355  val(:,4) = (/ 4._rp - 8._rp*x(1) - 4._rp*x(2), -4._rp*x(1) /)
356 
357  val(:,5) = (/ 4._rp*x(2), 4._rp*x(1) /)
358 
359  val(:,6) = (/ -4._rp*x(2), 4._rp - 4._rp*x(1) - 8._rp*x(2)/)
360 
361  end subroutine p2_2d_grad_u
362 
364  subroutine p1_3d_grad_u(val, n, x, m)
365  real(RP), dimension(m,n), intent(out) :: val
366  real(RP), dimension(m) , intent(in) :: x
367  integer , intent(in) :: n, m
368 
369  val(:,1) = (/-1._rp,-1._rp,-1._rp/)
370 
371  val(:,2) = (/ 1._rp, 0._rp, 0._rp/)
372 
373  val(:,3) = (/ 0._rp, 1._rp, 0._rp/)
374 
375  val(:,4) = (/ 0._rp, 0._rp, 1._rp/)
376 
377  end subroutine p1_3d_grad_u
378 
379 
380 
382  subroutine p2_3d_grad_u(val, n, x, m)
383  integer , intent(in) :: n, m
384  real(RP), dimension(m,n), intent(out) :: val
385  real(RP), dimension(m) , intent(in) :: x
386 
387  real(RP) :: a, b
388 
389  val = 0.0_rp
390 
391  val(:,1) = -3._rp + 4._rp * sum(x)
392 
393  val(1,2) = -1._rp + 4._rp * x(1)
394 
395  val(2,3) = -1._rp + 4._rp * x(2)
396 
397  val(3,4) = -1._rp + 4._rp * x(3)
398 
399  a = 4._rp - 8._rp*x(1) - 4._rp*x(2) - 4._rp*x(3)
400  b = -4._rp*x(1)
401  val(:,5) = (/a, b, b/)
402 
403  val(1,6) = 4._rp*x(2)
404  val(2,6) = 4._rp*x(1)
405 
406  a = - 4._rp*x(2)
407  b = 4._rp - 4._rp*x(1) - 8._rp*x(2) - 4._rp*x(3)
408  val(:,7) = (/a, b, a/)
409 
410  a = - 4._rp*x(3)
411  b = 4._rp - 4._rp*x(1) - 4._rp*x(2) - 8._rp*x(3)
412  val(:,8) = (/a, a, b/)
413 
414  val(1,9) = 4._rp*x(3)
415  val(3,9) = 4._rp*x(1)
416 
417  val(2,10) = 4._rp*x(3)
418  val(3,10) = 4._rp*x(2)
419 
420  end subroutine p2_3d_grad_u
421 
422 
424  subroutine p3_2d_grad_u(val, n, x, m)
425  real(RP), dimension(m,n), intent(out) :: val
426  real(RP), dimension(m) , intent(in) :: x
427  integer , intent(in) :: n, m
428 
429  real(RP) :: a, b
430 
431  a = -(27._rp/2._rp)*x(1)**2 - 27._rp*x(1)*x(2) + 18._rp*x(1)
432  a = a - (27._rp/2._rp)*x(2)**2 + 18._rp*x(2) - 11._rp/2._rp
433  val(:,1) = (/a, a/)
434 
435  a = (27._rp/2._rp)*x(1)**2 - 9._rp*x(1) + 1._rp
436  val(:,2) = (/a, 0._rp/)
437 
438  b = (27._rp/2._rp)*x(2)**2 - 9._rp*x(2) + 1._rp
439  val(:,3) = (/0._rp, b/)
440 
441  a = (81._rp/2._rp)*x(1)**2 + 54._rp*x(1)*x(2) - 45._rp*x(1)
442  a = a + (27._rp/2._rp)*x(2)**2 - (45._rp/2._rp)*x(2) + 9._rp
443  b = 27._rp*x(1)**2 + 27._rp*x(1)*x(2) - (45._rp/2._rp)*x(1)
444  val(:,4) = (/a, b/)
445 
446  a = -(81._rp/2._rp)*x(1)**2 - 27._rp*x(1)*x(2) + 36._rp*x(1) + (9._rp/2._rp)*x(2) - 9._rp/2._rp
447  b = -(9._rp/2._rp)*x(1)*(3._rp*x(1)-1._rp)
448  val(:,5) = (/a, b /)
449 
450  a = 27._rp*x(1)*x(2)-9._rp*x(2)/2._rp
451  b = -b
452  val(:,6) = (/a, b/)
453 
454  a = (9._rp/2._rp) * x(2) * (3._rp*x(2) - 1._rp)
455  b = 27._rp*x(1)*x(2)-(9._rp/2._rp)*x(1)
456  val(:,7) = (/a, b/)
457 
458  a = - a
459  b = -27._rp*x(1)*x(2)-(81._rp/2._rp)*x(2)**2 + 36._rp*x(2) + 4.5_rp*x(1) - 4.5_rp
460  val(:,8) = (/a, b/)
461 
462  a = 27._rp*x(1)*x(2) + 27._rp*x(2)**2 - 45._rp*x(2)/2._rp
463  b = 13.5_rp*x(1)**2 + 54._rp*x(1)*x(2) - 22.5_rp*x(1)
464  b= b + 40.5_rp*x(2)**2 - 45._rp*x(2) + 9._rp
465  val(:,9) = (/a, b/)
466 
467  a = -54._rp*x(1)*x(2) - 27._rp*x(2)**2 + 27._rp*x(2)
468  b = -27._rp*x(1)**2 - 54._rp*x(1)*x(2) + 27._rp*x(1)
469  val(:,10) = (/a, b/)
470 
471  end subroutine p3_2d_grad_u
472 
473 
474  ! %----------------------------------------%
475  ! | |
476  ! | FINITE ELEMENT BASE FUNCTIONS !
477  ! | VECTOR CASE |
478  ! | |
479  ! %----------------------------------------%
481  subroutine rt0_1d_phi(val, n, x, m)
482  real(RP), dimension(m,n), intent(out) :: val
483  real(RP), dimension(m) , intent(in) :: x
484  integer , intent(in) :: n, m
485 
486  val(1,1) = x(1) - 1._rp
487  val(1,2) = x(1)
488 
489  end subroutine rt0_1d_phi
490 
491 
494  subroutine rt0_1d_dual_phi(val, n, x, m)
495  real(RP), dimension(m,n), intent(out) :: val
496  real(RP), dimension(m) , intent(in) :: x
497  integer , intent(in) :: n, m
498 
499  val(1,1) = -20._rp*x(1)**3 + 30._rp*x(1)**2 - 9._rp*x(1) - 1._rp
500  val(1,2) = val(1,1) + 1._rp
501 
504  ! val(1,1) = -2._RP*x(1)**2 + 3._RP*x(1) - 1._RP
505  ! val(1,2) = 2._RP*x(1)**2 - x(1)
506 
507  end subroutine rt0_1d_dual_phi
508 
509 
511  subroutine rt0_2d_phi(val, n, x, m)
512  real(RP), dimension(m,n), intent(out) :: val
513  real(RP), dimension(m) , intent(in) :: x
514  integer , intent(in) :: n, m
515 
516  val(:,1) = x - (/0._rp, 1._rp/)
517  val(:,2) = x
518  val(:,3) = x - (/1._rp, 0._rp/)
519 
520  end subroutine rt0_2d_phi
521 
525  subroutine rt1_1d_phi(val, n, x, m)
526  real(RP), dimension(m,n), intent(out) :: val
527  real(RP), dimension(m) , intent(in) :: x
528  integer , intent(in) :: n, m
529 
530  val(1,1) =-2._rp*x(1)**2 + 3._rp*x(1) - 1._rp
531  val(1,2) = 2._rp*x(1)**2 - x(1)
532  val(1,3) =-4._rp*x(1)**2 + 4._rp*x(1)
533 
534  end subroutine rt1_1d_phi
535 
539  subroutine rt1_1d_2_phi(val, n, x, m)
540  real(RP), dimension(m,n), intent(out) :: val
541  real(RP), dimension(m) , intent(in) :: x
542  integer , intent(in) :: n, m
543 
544  val(1,1) =-3._rp*x(1)**2 + 4._rp*x(1) - 1._rp
545  val(1,2) = 3._rp*x(1)**2 - 2._rp*x(1)
546  val(1,3) =-6._rp*x(1)**2 + 6._rp*x(1)
547 
548  end subroutine rt1_1d_2_phi
549 
552  subroutine rt1_1d_dual_phi(val, n, x, m)
553  real(RP), dimension(m,n), intent(out) :: val
554  real(RP), dimension(m) , intent(in) :: x
555  integer , intent(in) :: n, m
556 
557  val(1,1) = 35._rp*x(1)**4 - 70._rp*x(1)**3 + 40._rp*x(1)**2 - 4._rp*x(1) - 1._rp
558  val(1,2) =-35._rp*x(1)**4 + 70._rp*x(1)**3 - 40._rp*x(1)**2 + 6._rp*x(1)
559  val(1,3) = 7._rp*x(1)**4 - 14._rp*x(1)**3 + 8._rp*x(1)**2 - 1._rp*x(1)
560 
561  !! choice 2
562  !! ref = applications/PG_star_1D_poisson/sage/PG2_1D_P2_2.sage
563  val(1,1) = 5._rp*x(1)**3 - 10._rp*x(1)**2 + 6._rp*x(1) - 1._rp
564  val(1,2) = 5._rp*x(1)**3 - 5._rp*x(1)**2 + x(1)
565  val(1,3) =-7._rp*x(1)**4 + 14._rp*x(1)**3 - 8._rp*x(1)**2 + x(1)
566 
567  end subroutine rt1_1d_dual_phi
568 
569 
572  subroutine rt1_1d_2_dual_phi(val, n, x, m)
573  real(RP), dimension(m,n), intent(out) :: val
574  real(RP), dimension(m) , intent(in) :: x
575  integer , intent(in) :: n, m
576 
577  val(1,1) = 5._rp*x(1)**3 - 10._rp*x(1)**2 + 6._rp*x(1) - 1._rp
578  val(1,2) = 5._rp*x(1)**3 - 5._rp*x(1)**2 + x(1)
579  val(1,3) = 7._rp*x(1)**4 - 14._rp*x(1)**3 + 9._rp*x(1)**2 - 2._rp* x(1)
580 
581  end subroutine rt1_1d_2_dual_phi
582 
583 
585  subroutine rt1_2d_2_phi(val, n, xx, m)
586  real(RP), dimension(m,n), intent(out) :: val
587  real(RP), dimension(m) , intent(in) :: xx
588  integer , intent(in) :: n, m
589 
590  real(RP) :: x, y, q, px, py, s3
591 
592  x = xx(1)
593  y = xx(2)
594  s3 = sqrt(3._rp)
595 
596  q = -( s3*x + ( s3 + 3._rp ) / 2._rp * y )
597  px = (1._rp + s3)/2._rp * x
598  py = -(1._rp + s3)/2._rp + s3*x + (2._rp + s3)*y
599  val(1,1) = px + q * x
600  val(2,1) = py + q * y
601 
602  q = ( s3*x + ( s3 - 3._rp ) / 2._rp * y )
603  px = (1._rp - s3)/2._rp * x
604  py = -(1._rp - s3)/2._rp - s3*x + (2._rp - s3)*y
605  val(1,2) = px + q * x
606  val(2,2) = py + q * y
607 
608  q = ( ( 3._rp + s3 )*x + &
609  & ( 3._rp - s3 )*y )/2._rp
610  px = -x
611  py = -y
612  val(1,3) = px + q * x
613  val(2,3) = py + q * y
614 
615  q = ( ( 3._rp - s3 )*x + &
616  & ( 3._rp + s3 )*y )/2._rp
617  px = -x
618  py = -y
619  val(1,4) = px + q * x
620  val(2,4) = py + q * y
621 
622  q = ( ( s3 - 3._rp)/2._rp * x + s3 * y )
623  px = (-1._rp + s3)/2._rp + (2._rp - s3)*x - s3*y
624  py = ( 1._rp - s3)/2._rp * y
625  val(1,5) = px + q * x
626  val(2,5) = py + q *y
627 
628  q = -( ( s3 + 3._rp)/2._rp * x + s3 * y )
629  px = (-1._rp - s3)/2._rp + (2._rp + s3)*x + s3*y
630  py = ( 1._rp + s3)/2._rp * y
631  val(1,6) = px + q * x
632  val(2,6) = py + q * y
633 
634  q = - ( 6._rp*x + 3._rp*y )
635  px = 6._rp*x
636  py = 3._rp*y
637  val(1,7) = px + q * x
638  val(2,7) = py + q * y
639 
640  q = - ( 3._rp*x + 6._rp*y )
641  px = 3._rp*x
642  py = 6._rp*y
643  val(1,8) = px + q * x
644  val(2,8) = py + q * y
645 
646  end subroutine rt1_2d_2_phi
647 
648 
649  ! %----------------------------------------%
650  ! | |
651  ! | FINITE ELEMENT BASE FUNCTIONS !
652  ! | VECTOR CASE : DIVERGENCE |
653  ! | |
654  ! %----------------------------------------%
656  subroutine rt0_1d_div_phi(val, n, x, m)
657  real(RP), dimension(n), intent(out) :: val
658  real(RP), dimension(m), intent(in) :: x
659  integer , intent(in) :: n, m
660 
661  val = 1._rp
662 
663  end subroutine rt0_1d_div_phi
664 
667  subroutine rt0_1d_dual_div_phi(val, n, x, m)
668  real(RP), dimension(n), intent(out) :: val
669  real(RP), dimension(m), intent(in) :: x
670  integer , intent(in) :: n, m
671 
672  val = -60._rp*x(1)**2 + 60._rp*x(1) - 9._rp
673 
676  ! val(1) = -4._RP*x(1) + 3._RP
677  ! val(2) = 4._RP*x(1) - 1._RP
678 
679  end subroutine rt0_1d_dual_div_phi
680 
682  subroutine rt0_2d_div_phi(val, n, x, m)
683  real(RP), dimension(n), intent(out) :: val
684  real(RP), dimension(m), intent(in) :: x
685  integer , intent(in) :: n, m
686 
687  val = 2._rp
688 
689  end subroutine rt0_2d_div_phi
690 
694  subroutine rt1_1d_div_phi(val, n, x, m)
695  real(RP), dimension(n), intent(out) :: val
696  real(RP), dimension(m), intent(in) :: x
697  integer , intent(in) :: n, m
698 
699  val(1) =-4._rp*x(1) + 3._rp
700  val(2) = 4._rp*x(1) - 1._rp
701  val(3) =-8._rp*x(1) + 4._rp
702 
703  end subroutine rt1_1d_div_phi
704 
708  subroutine rt1_1d_2_div_phi(val, n, x, m)
709  real(RP), dimension(n), intent(out) :: val
710  real(RP), dimension(m), intent(in) :: x
711  integer , intent(in) :: n, m
712 
713  val(1) = -6._rp*x(1) + 4._rp
714  val(2) = 6._rp*x(1) - 2._rp
715  val(3) =-12._rp*x(1) + 6._rp
716 
717  end subroutine rt1_1d_2_div_phi
718 
719 
722  subroutine rt1_1d_dual_div_phi(val, n, x, m)
723  real(RP), dimension(n), intent(out) :: val
724  real(RP), dimension(m), intent(in) :: x
725  integer , intent(in) :: n, m
726 
727  val(1) = 140._rp*x(1)**3 - 210._rp*x(1)**2 + 80._rp*x(1) - 4._rp
728  val(2) =-140._rp*x(1)**3 + 210._rp*x(1)**2 - 80._rp*x(1) + 6._rp
729  val(3) = 28._rp*x(1)**3 - 42._rp*x(1)**2 + 16._rp*x(1) - 1._rp
730 
731  !! choice 2
733  val(1) = 15._rp*x(1)**2 - 20._rp*x(1) + 6._rp
734  val(2) = 15._rp*x(1)**2 - 10._rp*x(1) + 1._rp
735  val(3) =-28._rp*x(1)**3 + 42._rp*x(1)**2 - 16._rp*x(1) + 1._rp
736 
737  end subroutine rt1_1d_dual_div_phi
738 
739 
742  subroutine rt1_1d_2_dual_div_phi(val, n, x, m)
743  real(RP), dimension(n), intent(out) :: val
744  real(RP), dimension(m), intent(in) :: x
745  integer , intent(in) :: n, m
746 
747  val(1) = 15._rp*x(1)**2 - 20._rp*x(1) + 6._rp
748  val(2) = 15._rp*x(1)**2 - 10._rp*x(1) + 1._rp
749  val(3) = 28._rp*x(1)**3 - 42._rp*x(1)**2 + 18._rp*x(1) - 2._rp
750 
751  end subroutine rt1_1d_2_dual_div_phi
752 
754  subroutine rt1_2d_2_div_phi(val, n, xx, m)
755  real(RP), dimension(n), intent(out) :: val
756  real(RP), dimension(m), intent(in) :: xx
757  integer , intent(in) :: n, m
758 
759  real(RP) :: x, y, s3
760 
761  x = xx(1)
762  y = xx(2)
763  s3 = sqrt(3._rp)
764 
765  val(1) = (5._rp+3._rp*s3)/2._rp - 3._rp*s3*x - 1.5_rp*(3._rp+s3)*y
766  val(2) = (5._rp-3._rp*s3)/2._rp + 3._rp*s3*x - 1.5_rp*(3._rp-s3)*y
767 
768  val(3) = -2._rp + 1.5_rp*( (3._rp+s3)*x + (3._rp-s3)*y )
769  val(4) = -2._rp + 1.5_rp*( (3._rp-s3)*x + (3._rp+s3)*y )
770 
771  val(5) = (5._rp-3._rp*s3)/2._rp - 1.5_rp*(3._rp-s3)*x + 3._rp*s3*y
772  val(6) = (5._rp+3._rp*s3)/2._rp - 1.5_rp*(3._rp+s3)*x - 3._rp*s3*y
773 
774  val(7) = 9._rp - 18._rp*x - 9._rp*y
775  val(8) = 9._rp - 9._rp*x - 18._rp*y
776 
777  end subroutine rt1_2d_2_div_phi
778 
779 end module fe_func_mod
780 
subroutine rt0_1d_dual_div_phi(val, n, x, m)
RT0 Dual ref choral/applications/PG_star_1D_poisson/sage/PG2_1D_P1_1.sage.
subroutine rt0_1d_phi(val, n, x, m)
RT0 1D.
DEFINITION OF FINITE ELEMENT BASIS FUNCTIONS
Definition: fe_func_mod.f90:32
subroutine p1_3d_u(val, n, x, m)
P1_3d.
subroutine p1_2d_disc_ortho_u(val, n, x, m)
P1_2d_disc_ortho.
subroutine rt1_2d_2_div_phi(val, n, xx, m)
RT1_2D.
subroutine rt1_1d_dual_phi(val, n, x, m)
RT1_1D_DUAL ref = applications/PG_star_1D_poisson/sage/PG2_1D_P2_1.sage.
subroutine rt1_1d_div_phi(val, n, x, m)
RT1_1D, considered with the usual basis for the DOFs: f(0), f(1), f(1/2) ref = applications/PG_star_1...
subroutine rt0_2d_div_phi(val, n, x, m)
RT0 2D.
subroutine p2_3d_u(val, n, x, m)
P2_3d.
subroutine p1_1d_u(val, n, x, m)
P1_1d.
Definition: fe_func_mod.f90:87
subroutine p1_1d_disc_ortho_u(val, n, x, m)
P1_1d_disc_ortho.
subroutine p3_2d_grad_u(val, n, x, m)
P3_2d.
integer, parameter, public rp
real(kind=RP) = real precision in the code REAL_TOL = epsilon to test real equality ...
Definition: real_type.F90:90
subroutine rt0_1d_dual_phi(val, n, x, m)
RT0 Dual ref choral/applications/PG_star_1D_poisson/sage/PG2_1D_P1_1.sage.
subroutine rt1_1d_2_dual_div_phi(val, n, x, m)
RT1_1D_2_DUAL ref = applications/PG_star_1D_poisson/sage/PG2_1D_P2_3.sage.
subroutine rt1_2d_2_phi(val, n, xx, m)
RT1_2D.
subroutine p1_1d_grad_u(val, n, x, m)
P1_1d.
conversion integers or rational to real
Definition: real_type.F90:153
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
Definition: real_type.F90:33
subroutine p2_2d_u(val, n, x, m)
P2_2d.
subroutine p0_2d_u(val, n, x, m)
P0_2d.
subroutine p2_3d_grad_u(val, n, x, m)
P2_3d.
subroutine rt1_1d_phi(val, n, x, m)
RT1_1D, considered with the usual basis for the DOFs: f(0), f(1), f(1/2) ref = applications/PG_star_1...
subroutine p2_1d_grad_u(val, n, x, m)
P2_1d.
subroutine p1_2d_grad_u(val, n, x, m)
P1_2d.
real(rp) function, dimension(3) grad_u(x)
gradient of the exact solution
real(rp) function u(x)
exact solution 'u'
subroutine p0_1d_u(val, n, x, m)
P0_1d.
subroutine p3_2d_u(val, n, x, m)
P3_2d.
subroutine rt1_1d_2_phi(val, n, x, m)
RT1_1D_2, considered with the alternative basis for the DOFs: f(0), f(1), ^1 f(x) dx ref = applicatio...
subroutine p1_2d_u(val, n, x, m)
P1_2d.
subroutine rt0_2d_phi(val, n, x, m)
RT0 2D.
subroutine rt1_1d_2_div_phi(val, n, x, m)
RT1 1D_2, considered with the alternative basis for the DOFs: f(0), f(1), ^1 f(x) dx ref = applicatio...
subroutine p2_2d_grad_u(val, n, x, m)
P2_2d.
subroutine rt1_1d_2_dual_phi(val, n, x, m)
RT1_1D_2_DUAL ref = applications/PG_star_1D_poisson/sage/PG2_1D_P2_3.sage.
subroutine p1_3d_grad_u(val, n, x, m)
P1_3d.
subroutine rt1_1d_dual_div_phi(val, n, x, m)
RT1_1D_DUAL ref = applications/PG_star_1D_poisson/sage/PG2_1D_P2_1.sage.
subroutine rt0_1d_div_phi(val, n, x, m)
RT0 1D.
subroutine p2_1d_u(val, n, x, m)
P2_1d.
Definition: fe_func_mod.f90:98
Functions associated with a finite element.
Definition: fe_func_mod.f90:41
subroutine p3_1d_u(val, n, x, m)
P3_1d.
subroutine p3_1d_grad_u(val, n, x, m)
P3_1d.