hex20n.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! Copyright (c) 2019 FrontISTR Commons
! This software is released under the MIT License, see LICENSE.txt
!-------------------------------------------------------------------------------
module shape_hex20n
  integer, parameter, private :: kreal = kind(0.0d0)
 
contains
  subroutine shapefunc_hex20n(localcoord,func)
    real(kind=kreal) :: localcoord(3)
    real(kind=kreal) :: func(20)
    real(kind=kreal) ri,si,ti,rp,sp,tp,rm,sm,tm
    ri=localcoord(1);  si=localcoord(2);  ti=localcoord(3)
    rp=1.0+ri; sp=1.0+si; tp=1.0+ti
    rm=1.0-ri; sm=1.0-si; tm=1.0-ti
    func(1)=-0.125*rm*sm*tm*(2.0+ri+si+ti)
    func(2)=-0.125*rp*sm*tm*(2.0-ri+si+ti)
    func(3)=-0.125*rp*sp*tm*(2.0-ri-si+ti)
    func(4)=-0.125*rm*sp*tm*(2.0+ri-si+ti)
    func(5)=-0.125*rm*sm*tp*(2.0+ri+si-ti)
    func(6)=-0.125*rp*sm*tp*(2.0-ri+si-ti)
    func(7)=-0.125*rp*sp*tp*(2.0-ri-si-ti)
    func(8)=-0.125*rm*sp*tp*(2.0+ri-si-ti)
    func(9)=0.25*(1.0-ri**2)*sm*tm
    func(10)=0.25*rp*(1.0-si**2)*tm
    func(11)=0.25*(1.0-ri**2)*sp*tm
    func(12)=0.25*rm*(1.0-si**2)*tm
    func(13)=0.25*(1.0-ri**2)*sm*tp
    func(14)=0.25*rp*(1.0-si**2)*tp
    func(15)=0.25*(1.0-ri**2)*sp*tp
    func(16)=0.25*rm*(1.0-si**2)*tp
    func(17)=0.25*rm*sm*(1.0-ti**2)
    func(18)=0.25*rp*sm*(1.0-ti**2)
    func(19)=0.25*rp*sp*(1.0-ti**2)
    func(20)=0.25*rm*sp*(1.0-ti**2)
  end subroutine
 
  subroutine shapederiv_hex20n(localcoord,func)
    real(kind=kreal) :: localcoord(3)
    real(kind=kreal) :: func(20,3)
    real(kind=kreal) ri,si,ti,rp,sp,tp,rm,sm,tm
    ri=localcoord(1);  si=localcoord(2);  ti=localcoord(3)
    rp=1.d0+ri; sp=1.d0+si; tp=1.d0+ti
    rm=1.d0-ri; sm=1.d0-si; tm=1.d0-ti
    !  FOR R-COORDINATE
    func(1,1)=-0.125*rm*sm*tm+0.125*sm*tm*(2.0+ri+si+ti)
    func(2,1)=+0.125*rp*sm*tm-0.125*sm*tm*(2.0-ri+si+ti)
    func(3,1)=+0.125*rp*sp*tm-0.125*sp*tm*(2.0-ri-si+ti)
    func(4,1)=-0.125*rm*sp*tm+0.125*sp*tm*(2.0+ri-si+ti)
    func(5,1)=-0.125*rm*sm*tp+0.125*sm*tp*(2.0+ri+si-ti)
    func(6,1)=+0.125*rp*sm*tp-0.125*sm*tp*(2.0-ri+si-ti)
    func(7,1)=+0.125*rp*sp*tp-0.125*sp*tp*(2.0-ri-si-ti)
    func(8,1)=-0.125*rm*sp*tp+0.125*sp*tp*(2.0+ri-si-ti)
    func(9,1 )=-0.50*ri*sm*tm
    func(10,1)=+0.25*(1.0-si**2)*tm
    func(11,1)=-0.50*ri*sp*tm
    func(12,1)=-0.25*(1.0-si**2)*tm
    func(13,1)=-0.50*ri*sm*tp
    func(14,1)=+0.25*(1.0-si**2)*tp
    func(15,1)=-0.50*ri*sp*tp
    func(16,1)=-0.25*(1.0-si**2)*tp
    func(17,1)=-0.25*sm*(1.0-ti**2)
    func(18,1)=+0.25*sm*(1.0-ti**2)
    func(19,1)=+0.25*sp*(1.0-ti**2)
    func(20,1)=-0.25*sp*(1.0-ti**2)
    !  FOR S-COORDINATE
    func(1,2)=-0.125*rm*sm*tm+0.125*rm*tm*(2.0+ri+si+ti)
    func(2,2)=-0.125*rp*sm*tm+0.125*rp*tm*(2.0-ri+si+ti)
    func(3,2)=+0.125*rp*sp*tm-0.125*rp*tm*(2.0-ri-si+ti)
    func(4,2)=+0.125*rm*sp*tm-0.125*rm*tm*(2.0+ri-si+ti)
    func(5,2)=-0.125*rm*sm*tp+0.125*rm*tp*(2.0+ri+si-ti)
    func(6,2)=-0.125*rp*sm*tp+0.125*rp*tp*(2.0-ri+si-ti)
    func(7,2)=+0.125*rp*sp*tp-0.125*rp*tp*(2.0-ri-si-ti)
    func(8,2)=+0.125*rm*sp*tp-0.125*rm*tp*(2.0+ri-si-ti)
    func(9,2)=-0.25*(1.0-ri**2)*tm
    func(10,2)=-0.50*rp*si*tm
    func(11,2)=+0.25*(1.0-ri**2)*tm
    func(12,2)=-0.50*rm*si*tm
    func(13,2)=-0.25*(1.0-ri**2)*tp
    func(14,2)=-0.50*rp*si*tp
    func(15,2)=+0.25*(1.0-ri**2)*tp
    func(16,2)=-0.50*rm*si*tp
    func(17,2)=-0.25*rm*(1.0-ti**2)
    func(18,2)=-0.25*rp*(1.0-ti**2)
    func(19,2)=+0.25*rp*(1.0-ti**2)
    func(20,2)=+0.25*rm*(1.0-ti**2)
    !  FOR T-COORDINATE
    func(1,3)=-0.125*rm*sm*tm+0.125*rm*sm*(2.0+ri+si+ti)
    func(2,3)=-0.125*rp*sm*tm+0.125*rp*sm*(2.0-ri+si+ti)
    func(3,3)=-0.125*rp*sp*tm+0.125*rp*sp*(2.0-ri-si+ti)
    func(4,3)=-0.125*rm*sp*tm+0.125*rm*sp*(2.0+ri-si+ti)
    func(5,3)=+0.125*rm*sm*tp-0.125*rm*sm*(2.0+ri+si-ti)
    func(6,3)=+0.125*rp*sm*tp-0.125*rp*sm*(2.0-ri+si-ti)
    func(7,3)=+0.125*rp*sp*tp-0.125*rp*sp*(2.0-ri-si-ti)
    func(8,3)=+0.125*rm*sp*tp-0.125*rm*sp*(2.0+ri-si-ti)
    func(9,3)=-0.25*(1.0-ri**2)*sm
    func(10,3)=-0.25*rp*(1.0-si**2)
    func(11,3)=-0.25*(1.0-ri**2)*sp
    func(12,3)=-0.25*rm*(1.0-si**2)
    func(13,3)=0.25*(1.0-ri**2)*sm
    func(14,3)=0.25*rp*(1.0-si**2)
    func(15,3)=0.25*(1.0-ri**2)*sp
    func(16,3)=0.25*rm*(1.0-si**2)
    func(17,3)=-0.5*rm*sm*ti
    func(18,3)=-0.5*rp*sm*ti
    func(19,3)=-0.5*rp*sp*ti
    func(20,3)=-0.5*rm*sp*ti
  end subroutine
 
end module