SFEMaNS  version 5.3
Reference documentation for SFEMaNS
GP_2D_p1.f90
Go to the documentation of this file.
1 !===
2 !Authors: Jean-Luc Guermond, Lugi Quartapelle, Copyright 1994
3 !Author: Jean-Luc Guermond, Copyright December 24th 2018
4 !===
5 MODULE gp_2d_p1
6  PRIVATE
8 CONTAINS
9  SUBROUTINE element_2d_p1 (w, d, p, n_w, l_G)
10  !===Triangular element with linear interpolation
11  !===and three Gauss integration points
12  !===w(n_w, l_G) : values of shape functions at Gauss points
13  !===d(2, n_w, l_G) : derivatives values of shape functions at Gauss points
14  !===p(l_G) : weight for Gaussian quadrature at Gauss points
15  IMPLICIT NONE
16  INTEGER, INTENT(IN) :: n_w, l_G
17  REAL(KIND=8), DIMENSION( n_w, l_G), INTENT(OUT) :: w
18  REAL(KIND=8), DIMENSION(2, n_w, l_G), INTENT(OUT) :: d
19  REAL(KIND=8), DIMENSION(l_G), INTENT(OUT) :: p
20  REAL(KIND=8), DIMENSION(l_G) :: xx, yy
21  INTEGER :: j
22  REAL(KIND=8) :: zero = 0, one = 1, two = 2, three = 3, six = 6
23  REAL(KIND=8) :: f1, f2, f3, x, y
24  f1(x, y) = one - x - y
25  f2(x, y) = x
26  f3(x, y) = y
27  xx(1) = one/six; xx(2) = two/three; xx(3) = one/six
28  yy(1) = one/six; yy(2) = one/six; yy(3) = two/three
29  DO j = 1, l_g
30  w(1, j) = f1(xx(j), yy(j))
31  d(1, 1, j) = - one
32  d(2, 1, j) = - one
33  w(2, j) = f2(xx(j), yy(j))
34  d(1, 2, j) = one
35  d(2, 2, j) = zero
36  w(3, j) = f3(xx(j), yy(j))
37  d(1, 3, j) = zero
38  d(2, 3, j) = one
39  p(j) = one/six
40  ENDDO
41  END SUBROUTINE element_2d_p1
42 
43  SUBROUTINE element_2d_p1_boundary (face, d, w, n_w, l_Gs)
44  !===Triangular element with linear interpolation
45  !===and three Gauss integration points
46  !===w(n_w, l_Gs) : values of shape functions at Gauss points
47  !===d(2, n_w, l_Gs) : derivatives values of shape functions at Gauss points
48  !===p(l_Gs) : weight for Gaussian quadrature at Gauss points
49  ! 3
50  ! 1 2 with orientation 1->2, 1->3, 2->3 (lowest to highest index convention)
51  IMPLICIT NONE
52  INTEGER, INTENT(IN) :: n_w, l_Gs
53  INTEGER, INTENT(IN) :: face
54  REAL(KIND=8), DIMENSION(2, n_w, l_Gs), INTENT(OUT) :: d
55  REAL(KIND=8), DIMENSION( n_w, l_Gs), INTENT(OUT) :: w
56  REAL(KIND=8), DIMENSION(l_Gs) :: xx, yy
57  INTEGER :: j
58  REAL(KIND=8) :: zero = 0, one = 1, three=3, half = 0.5d0
59  REAL(KIND=8) :: f1, f2, f3, x, y
60  f1(x, y) = 1-x-y
61  f2(x, y) = x
62  f3(x, y) = y
63  IF (face==1) THEN
64  xx(1) = half + half/sqrt(three)
65  xx(2) = half - half/sqrt(three)
66  yy(1) = half - half/sqrt(three)
67  yy(2) = half + half/sqrt(three)
68  ELSE IF (face==2) THEN
69  xx(1) = zero
70  xx(2) = zero
71  yy(1) = half - half/sqrt(three)
72  yy(2) = half + half/sqrt(three)
73  ELSE
74  xx(1) = half - half/sqrt(three)
75  xx(2) = half + half/sqrt(three)
76  yy(1) = zero
77  yy(2) = zero
78  END IF
79  DO j = 1, l_gs
80  w(1,j) = f1(xx(j), yy(j))
81  d(1, 1, j) = -one
82  d(2, 1, j) = -one
83  w(2,j) = f2(xx(j), yy(j))
84  d(1, 2, j) = one
85  d(2, 2, j) = zero
86  w(3,j) = f3(xx(j), yy(j))
87  d(1, 3, j) = zero
88  d(2, 3, j) = one
89  ENDDO
90  END SUBROUTINE element_2d_p1_boundary
91 
92  SUBROUTINE element_1d_p1 (w, d, p, n_ws, l_Gs)
93  !===One-dimensional element with linear interpolation
94  !===and two Gauss integration points
95  !===w(n_w, l_G) : values of shape functions at Gauss points
96  !===d(1, n_w, l_G) : derivatives values of shape functions at Gauss points
97  !===p(l_G) : weight for Gaussian quadrature at Gauss points
98  IMPLICIT NONE
99  INTEGER, INTENT(IN) :: n_ws, l_Gs
100  REAL(KIND=8), DIMENSION( n_ws, l_Gs), INTENT(OUT) :: w
101  REAL(KIND=8), DIMENSION(1, n_ws, l_Gs), INTENT(OUT) :: d
102  REAL(KIND=8), DIMENSION(l_Gs), INTENT(OUT) :: p
103  REAL(KIND=8), DIMENSION(l_Gs) :: xx
104  INTEGER :: j
105  REAL(KIND=8) :: one = 1.d0, two = 2.d0, three = 3.d0
106  REAL(KIND=8) :: f1, f2, x
107  f1(x) = (one - x)/two
108  f2(x) = (x + one)/two
109  xx(1) = - one/sqrt(three)
110  xx(2) = + one/sqrt(three)
111  DO j = 1, l_gs
112  w(1, j) = f1(xx(j))
113  d(1, 1, j) = - one/two
114  w(2, j) = f2(xx(j))
115  d(1, 2, j) = + one/two
116  p(j) = one
117  ENDDO
118  END SUBROUTINE element_1d_p1
119 
120  SUBROUTINE element_1d_p1_at_nodes (d, n_ws)
121  IMPLICIT NONE
122  INTEGER, INTENT(IN) :: n_ws
123  REAL(KIND=8), DIMENSION(n_ws, n_ws), INTENT(OUT) :: d
124  INTEGER :: j
125  REAL(KIND=8) :: one = 1.d0, two = 2.d0
126  DO j = 1, n_ws
127  d(1, j) = - one/two
128  d(2, j) = + one/two
129  ENDDO
130  END SUBROUTINE element_1d_p1_at_nodes
131 
132 END MODULE gp_2d_p1
133 
subroutine, public element_1d_p1(w, d, p, n_ws, l_Gs)
Definition: GP_2D_p1.f90:93
subroutine, public element_2d_p1(w, d, p, n_w, l_G)
Definition: GP_2D_p1.f90:10
subroutine, public element_1d_p1_at_nodes(d, n_ws)
Definition: GP_2D_p1.f90:121
subroutine, public element_2d_p1_boundary(face, d, w, n_w, l_Gs)
Definition: GP_2D_p1.f90:44