FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
tri6n.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!-------------------------------------------------------------------------------
8 integer, parameter, private :: kreal = kind(0.0d0)
9
10contains
11 subroutine shapefunc_tri6n(areacoord,func)
12 real(kind=kreal), intent(in) :: areacoord(2)
13 real(kind=kreal) :: func(6)
14 real(kind=kreal) :: xi,et,st
15 xi=areacoord(1); et=areacoord(2); st=1.d0-xi-et
16
17 func(1)=st*(2.d0*st-1.d0)
18 func(2)=xi*(2.d0*xi-1.d0)
19 func(3)=et*(2.d0*et-1.d0)
20 func(4)=4.d0*xi*st
21 func(5)=4.d0*xi*et
22 func(6)=4.d0*et*st
23 end subroutine
24
25 subroutine shapederiv_tri6n(areacoord,func)
26 real(kind=kreal), intent(in) :: areacoord(2)
27 real(kind=kreal) :: func(6,2)
28 real(kind=kreal) :: xi,et,st
29 xi=areacoord(1); et=areacoord(2); st=1.d0-xi-et
30
31 func(1,1)=1.d0-4.d0*st
32 func(2,1)=4.d0*xi-1.d0
33 func(3,1)=0.d0
34 func(4,1)=4.d0*(st-xi)
35 func(5,1)=4.d0*et
36 func(6,1)=-4.d0*et
37
38 func(1,2)=1.d0-4.d0*st
39 func(2,2)=0.d0
40 func(3,2)=4.d0*et-1.d0
41 func(4,2)=-4.d0*xi
42 func(5,2)=4.d0*xi
43 func(6,2)=4.d0*(st-et)
44
45 end subroutine
46
47 subroutine shape2ndderiv_tri6n(func)
48 real(kind=kreal) :: func(6,2,2)
49 func(1,1,1) = 4.d0; func(1,1,2) = 4.d0
50 func(2,1,1) = 4.d0; func(2,1,2) = 0.d0
51 func(3,1,1) = 0.d0; func(3,1,2) = 0.d0
52 func(4,1,1) =-8.d0; func(4,1,2) = -4.d0
53 func(5,1,1) = 0.d0; func(5,1,2) = 4.d0
54 func(6,1,1) = 0.d0; func(6,1,2) = -4.d0
55
56 func(1,2,1) = 4.d0; func(1,2,2) = 4.d0
57 func(2,2,1) = 0.d0; func(2,2,2) = 0.d0
58 func(3,2,1) = 0.d0; func(3,2,2) = 4.d0
59 func(4,2,1) =-4.d0; func(4,2,2) = 0.d0
60 func(5,2,1) = 4.d0; func(5,2,2) = 0.d0
61 func(6,2,1) =-4.d0; func(6,2,2) =-8.d0
62 end subroutine
63
64end module
This module contains functions for interpolation in 6 node trianglar element (Langrange interpolation...
Definition: tri6n.f90:7
subroutine shapefunc_tri6n(areacoord, func)
Definition: tri6n.f90:12
subroutine shape2ndderiv_tri6n(func)
Definition: tri6n.f90:48
subroutine shapederiv_tri6n(areacoord, func)
Definition: tri6n.f90:26