FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_Jacob_241.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
6
8contains
9
10 subroutine hecmw_jacob_241 ( hecMESH, iElem, DET, W, N, NX, NY )
11 use hecmw_util
12
13 implicit none
14
15 type(hecmwst_local_mesh):: hecMESH
16 integer(kind=kint):: iElem
17 real(kind=kreal):: det
18 real(kind=kreal):: w(4), n(4,4), nx(4,4), ny(4,4)
19
20 integer(kind=kint):: i, j, jj, iLocal
21 integer(kind=kint):: LX, LY
22 real(kind=kreal):: dum
23 real(kind=kreal):: xx(4), yy(4)
24 real(kind=kreal):: xg(2), hr(2,2,4), hs(2,2,4), ht(2,2,4)
25 real(kind=kreal):: h(2,2,4),bx(2,2,4),by(2,2,4),bz(2,2,4)
26 real(kind=kreal):: ri, si, rp, sp, rm, sm
27 real(kind=kreal):: xj11,xj12,xj21,xj22
28 real(kind=kreal):: xji11,xji12,xji21,xji22
29
30 data xg/-0.5773502691896d0, 0.5773502691896d0/
31
32 do j = 1, 4
33 jj = 4 - j
34 ilocal = hecmesh%elem_node_item ( 4*ielem -jj )
35 xx(j) = hecmesh%node( ilocal*2 -1 )
36 yy(j) = hecmesh%node( ilocal*2 )
37 w(j) = 1.0d0
38 end do
39 !C
40 !C*LOOP OVER ALL INTEGRATION POINTS
41 do lx=1,2
42 ri=xg(lx)
43 do ly=1,2
44 si=xg(ly)
45 rp=1.0+ri
46 sp=1.0+si
47 rm=1.0-ri
48 sm=1.0-si
49 !!
50 !! ******** Shape functions ********
51 !!
52 h(lx,ly,1)=0.25*rm*sm
53 h(lx,ly,2)=0.25*rp*sm
54 h(lx,ly,3)=0.25*rp*sp
55 h(lx,ly,4)=0.25*rm*sp
56 !!
57 !! ******** Derivative of shape functions ********
58 !!
59 !! ----------- For R-Coordinate -------------
60 hr(lx,ly,1)=-.25*sm
61 hr(lx,ly,2)= .25*sm
62 hr(lx,ly,3)= .25*sp
63 hr(lx,ly,4)=-.25*sp
64 !! ----------- For S-Coordinate -------------
65 hs(lx,ly,1)=-.25*rm
66 hs(lx,ly,2)=-.25*rp
67 hs(lx,ly,3)= .25*rp
68 hs(lx,ly,4)= .25*rm
69 !!
70 !! ******** Jacobi matrix calculation********
71 !!
72 xj11=0.0
73 xj21=0.0
74 xj12=0.0
75 xj22=0.0
76 do i=1,4
77 xj11=xj11+hr(lx,ly,i)*xx(i)
78 xj21=xj21+hs(lx,ly,i)*xx(i)
79 xj12=xj12+hr(lx,ly,i)*yy(i)
80 xj22=xj22+hs(lx,ly,i)*yy(i)
81 end do
82 det=xj11*xj22-xj21*xj12
83 !!
84 !! ******** Inverse Jacobi matrix calculation ********
85 !!
86 dum=1.0/det
87 xji11= xj22*dum
88 xji12=-xj12*dum
89 xji21=-xj21*dum
90 xji22= xj11*dum
91 do j=1, 4
92 bx(lx,ly,j)=xji11*hr(lx,ly,j)+xji12*hs(lx,ly,j)
93 by(lx,ly,j)=xji21*hr(lx,ly,j)+xji22*hs(lx,ly,j)
94 end do
95 !C
96 end do
97 end do
98 !C
99
100 j = 1
101 do lx = 1, 2
102 do ly = 1, 2
103
104 n(j,1) = h(lx,ly,1)
105 n(j,2) = h(lx,ly,2)
106 n(j,3) = h(lx,ly,3)
107 n(j,4) = h(lx,ly,4)
108
109 nx(j,1) = bx(lx,ly,1)
110 nx(j,2) = bx(lx,ly,2)
111 nx(j,3) = bx(lx,ly,3)
112 nx(j,4) = bx(lx,ly,4)
113
114 ny(j,1) = by(lx,ly,1)
115 ny(j,2) = by(lx,ly,2)
116 ny(j,3) = by(lx,ly,3)
117 ny(j,4) = by(lx,ly,4)
118
119 j = j + 1
120 end do
121 end do
122
123 !C
124 !C
125 end subroutine hecmw_jacob_241
126
127end module hecmw_jacob241
Jacobian calculation.
subroutine hecmw_jacob_241(hecmesh, ielem, det, w, n, nx, ny)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal