FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_Jacob_341.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_341 ( hecMESH, iElem, DET, W, N, NX, NY, NZ)
11 use hecmw_util
12 implicit none
13
14 type(hecmwst_local_mesh):: hecMESH
15 integer(kind=kint):: iElem
16 real(kind=kreal):: det
17 real(kind=kreal):: w(4),n(4,4),nx(4,4),ny(4,4),nz(4,4)
18
19 integer(kind=kint):: i, ii, j, iLocal
20 real(kind=kreal):: dum
21 real(kind=kreal):: xx(4), yy(4), zz(4)
22 real(kind=kreal):: nr(4), ns(4), nt(4)
23 real(kind=kreal):: jacob(3,3), jinv(3,3)
24 real(kind=kreal),parameter:: alpha = 0.58541020, beta = 0.13819660
25
26 !!
27 !! ******** Set values of coordinates********
28 !!
29 do i = 1, 4
30 ii = 4 - i
31 ilocal = hecmesh%elem_node_item( 4*ielem -ii )
32 xx(i) = hecmesh%node( ilocal*3 -2 )
33 yy(i) = hecmesh%node( ilocal*3 -1 )
34 zz(i) = hecmesh%node( ilocal*3 )
35 w(i) = 0.25d0
36 end do
37
38 !!
39 !! ******** Set values to shape functions ********
40 !!
41 n(1,1) = alpha
42 n(2,1) = beta
43 n(3,1) = beta
44 n(4,1) = beta
45
46 n(1,2) = beta
47 n(2,2) = alpha
48 n(3,2) = beta
49 n(4,2) = beta
50
51 n(1,3) = beta
52 n(2,3) = beta
53 n(3,3) = alpha
54 n(4,3) = beta
55
56 n(1,4) = beta
57 n(2,4) = beta
58 n(3,4) = beta
59 n(4,4) = alpha
60 !!
61 !! ******** Derivative of shape functions ********
62 !!
63 !! ----------- For R-Coordinate -------------
64 nr(1) = 1.0
65 nr(2) = 0.0
66 nr(3) = 0.0
67 nr(4) =-1.0
68 !! ----------- For S-Coordinate -------------
69 ns(1) = 0.0
70 ns(2) = 1.0
71 ns(3) = 0.0
72 ns(4) =-1.0
73 !! ----------- For T-Coordinate -------------
74 nt(1) = 0.0
75 nt(2) = 0.0
76 nt(3) = 1.0
77 nt(4) =-1.0
78 !!
79 !! ******** Jacobi matrix calculation********
80 !!
81 jacob(1,1) = nr(1)*xx(1)+nr(2)*xx(2)+nr(3)*xx(3)+nr(4)*xx(4)
82 jacob(2,1) = ns(1)*xx(1)+ns(2)*xx(2)+ns(3)*xx(3)+ns(4)*xx(4)
83 jacob(3,1) = nt(1)*xx(1)+nt(2)*xx(2)+nt(3)*xx(3)+nt(4)*xx(4)
84 jacob(1,2) = nr(1)*yy(1)+nr(2)*yy(2)+nr(3)*yy(3)+nr(4)*yy(4)
85 jacob(2,2) = ns(1)*yy(1)+ns(2)*yy(2)+ns(3)*yy(3)+ns(4)*yy(4)
86 jacob(3,2) = nt(1)*yy(1)+nt(2)*yy(2)+nt(3)*yy(3)+nt(4)*yy(4)
87 jacob(1,3) = nr(1)*zz(1)+nr(2)*zz(2)+nr(3)*zz(3)+nr(4)*zz(4)
88 jacob(2,3) = ns(1)*zz(1)+ns(2)*zz(2)+ns(3)*zz(3)+ns(4)*zz(4)
89 jacob(3,3) = nt(1)*zz(1)+nt(2)*zz(2)+nt(3)*zz(3)+nt(4)*zz(4)
90
91 det = jacob(1,1) * jacob(2,2) * jacob(3,3) &
92 & + jacob(1,2) * jacob(2,3) * jacob(3,1) &
93 & + jacob(1,3) * jacob(2,1) * jacob(3,2) &
94 & - jacob(1,3) * jacob(2,2) * jacob(3,1) &
95 & - jacob(1,2) * jacob(2,1) * jacob(3,3) &
96 & - jacob(1,1) * jacob(2,3) * jacob(3,2)
97 !!
98 !! ******** Inverse Jacobi matrix calculation ********
99 !!
100 dum = 1.0 / det
101 jinv(1,1) = dum*( jacob(2,2)*jacob(3,3) -jacob(2,3)*jacob(3,2) )
102 jinv(2,1) = dum*( -jacob(2,1)*jacob(3,3) +jacob(2,3)*jacob(3,1) )
103 jinv(3,1) = dum*( jacob(2,1)*jacob(3,2) -jacob(2,2)*jacob(3,1) )
104 jinv(1,2) = dum*( -jacob(1,2)*jacob(3,3) +jacob(1,3)*jacob(3,2) )
105 jinv(2,2) = dum*( jacob(1,1)*jacob(3,3) -jacob(1,3)*jacob(3,1) )
106 jinv(3,2) = dum*( -jacob(1,1)*jacob(3,2) +jacob(1,2)*jacob(3,1) )
107 jinv(1,3) = dum*( jacob(1,2)*jacob(2,3) -jacob(1,3)*jacob(2,2) )
108 jinv(2,3) = dum*( -jacob(1,1)*jacob(2,3) +jacob(1,3)*jacob(2,1) )
109 jinv(3,3) = dum*( jacob(1,1)*jacob(2,2) -jacob(1,2)*jacob(2,1) )
110
111 do i = 1, 4
112 do j = 1, 4
113
114 nx(i,j) = jinv(1,1)*nr(i) + jinv(1,2)*ns(i) +jinv(1,3)*nt(i)
115 ny(i,j) = jinv(2,1)*nr(i) + jinv(2,2)*ns(i) +jinv(2,3)*nt(i)
116 nz(i,j) = jinv(3,1)*nr(i) + jinv(3,2)*ns(i) +jinv(3,3)*nt(i)
117
118 end do
119 end do
120
121 end subroutine hecmw_jacob_341
122
123end module hecmw_jacob341
Jacobian calculation.
subroutine hecmw_jacob_341(hecmesh, ielem, det, w, n, nx, ny, nz)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal