1! 2! CalculiX - A 3-dimensional finite element program 3! Copyright (C) 1998-2021 Guido Dhondt 4! 5! This program is free software; you can redistribute it and/or 6! modify it under the terms of the GNU General Public License as 7! published by the Free Software Foundation(version 2); 8! 9! 10! This program is distributed in the hope that it will be useful, 11! but WITHOUT ANY WARRANTY; without even the implied warranty of 12! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13! GNU General Public License for more details. 14! 15! You should have received a copy of the GNU General Public License 16! along with this program; if not, write to the Free Software 17! Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 18! 19 subroutine shape15w(xi,et,ze,xl,xsj,shp,iflag) 20! 21! shape functions and derivatives for a 15-node quadratic 22! isoparametric wedge element. 0<=xi,et<=1,-1<=ze<=1,xi+et<=1. 23! 24! iflag=1: calculate only the value of the shape functions 25! iflag=2: calculate the value of the shape functions and 26! the Jacobian determinant 27! iflag=3: calculate the value of the shape functions, the 28! value of their derivatives w.r.t. the global 29! coordinates and the Jacobian determinant 30! 31! 32! Copyright (c) 2003 WB 33! 34! Written February 2003 on the basis of the Guido's shape function files 35! 36 implicit none 37! 38 integer i,j,k,iflag 39! 40 real*8 shp(4,15),xs(3,3),xsi(3,3),xl(3,15),sh(3) 41! 42 real*8 xi,et,ze,xsj,a 43! 44! 45! 46! shape functions and their glocal derivatives 47! 48 a=1.d0-xi-et 49! 50! shape functions 51! 52 shp(4,1)=-0.5d0*a*(1.d0-ze)*(2.d0*xi+2.d0*et+ze) 53 shp(4,2)=0.5d0*xi*(1.d0-ze)*(2.d0*xi-2.d0-ze) 54 shp(4,3)=0.5d0*et*(1.d0-ze)*(2.d0*et-2.d0-ze) 55 shp(4,4)=-0.5d0*a*(1.d0+ze)*(2.d0*xi+2.d0*et-ze) 56 shp(4,5)=0.5d0*xi*(1.d0+ze)*(2.d0*xi-2.d0+ze) 57 shp(4,6)=0.5d0*et*(1.d0+ze)*(2.d0*et-2.d0+ze) 58 shp(4,7)=2.d0*xi*a*(1.d0-ze) 59 shp(4,8)=2.d0*xi*et*(1.d0-ze) 60 shp(4,9)=2.d0*et*a*(1.d0-ze) 61 shp(4,10)=2.d0*xi*a*(1.d0+ze) 62 shp(4,11)=2.d0*xi*et*(1.d0+ze) 63 shp(4,12)=2.d0*et*a*(1.d0+ze) 64 shp(4,13)= a*(1.d0-ze*ze) 65 shp(4,14)=xi*(1.d0-ze*ze) 66 shp(4,15)=et*(1.d0-ze*ze) 67! 68 if(iflag.eq.1) return 69! 70! local derivatives of the shape functions: xi-derivative 71! 72 shp(1,1)= 0.5d0*(1.d0-ze)*(4.0d0*xi+4.0d0*et+ze-2.d0) 73 shp(1,2)= 0.5d0*(1.d0-ze)*(4.0d0*xi-ze-2.d0) 74 shp(1,3)= 0.d0 75 shp(1,4)= 0.5d0*(1.d0+ze)*(4.0d0*xi+4.0d0*et-ze-2.d0) 76 shp(1,5)= 0.5d0*(1.d0+ze)*(4.0d0*xi+ze-2.d0) 77 shp(1,6)= 0.d0 78 shp(1,7)= 2.d0*(1.d0-ze)*(1.d0-2.d0*xi-et) 79 shp(1,8)= 2.d0*et*(1.d0-ze) 80 shp(1,9)= -2.d0*et*(1.d0-ze) 81 shp(1,10)= 2.d0*(1.d0+ze)*(1.d0-2.d0*xi-et) 82 shp(1,11)= 2.d0*et*(1.d0+ze) 83 shp(1,12)= -2.d0*et*(1.d0+ze) 84 shp(1,13)= -(1.d0-ze*ze) 85 shp(1,14)= (1.d0-ze*ze) 86 shp(1,15)= 0.d0 87! 88! local derivatives of the shape functions: eta-derivative 89! 90 shp(2,1)= 0.5d0*(1.d0-ze)*(4.0d0*xi+4.0d0*et+ze-2.d0) 91 shp(2,2)= 0.d0 92 shp(2,3)= 0.5d0*(1.d0-ze)*(4.0d0*et-ze-2.d0) 93 shp(2,4)= 0.5d0*(1.d0+ze)*(4.0d0*xi+4.0d0*et-ze-2.d0) 94 shp(2,5)= 0.d0 95 shp(2,6)= 0.5d0*(1.d0+ze)*(4.0d0*et+ze-2.d0) 96 shp(2,7)=-2.d0*xi*(1.d0-ze) 97 shp(2,8)= 2.d0*xi*(1.d0-ze) 98 shp(2,9)= 2.d0*(1.d0-ze)*(1.d0-xi-2.d0*et) 99 shp(2,10)=-2.d0*xi*(1.d0+ze) 100 shp(2,11)= 2.d0*xi*(1.d0+ze) 101 shp(2,12)= 2.d0*(1.d0+ze)*(1.d0-xi-2.d0*et) 102 shp(2,13)=-(1.d0-ze*ze) 103 shp(2,14)= 0.0d0 104 shp(2,15)= (1.d0-ze*ze) 105! 106! local derivatives of the shape functions: zeta-derivative 107! 108 shp(3,1)= a*(xi+et+ze-0.5d0) 109 shp(3,2)= xi*(-xi+ze+0.5d0) 110 shp(3,3)= et*(-et+ze+0.5d0) 111 shp(3,4)= a*(-xi-et+ze+0.5d0) 112 shp(3,5)= xi*(xi+ze-0.5d0) 113 shp(3,6)= et*(et+ze-0.5d0) 114 shp(3,7)= -2*xi*a 115 shp(3,8)= -2*xi*et 116 shp(3,9)= -2*et*a 117 shp(3,10)= 2*xi*a 118 shp(3,11)= 2*xi*et 119 shp(3,12)= 2*et*a 120 shp(3,13)=-2*a*ze 121 shp(3,14)=-2*xi*ze 122 shp(3,15)=-2*et*ze 123! 124! computation of the local derivative of the global coordinates 125! (xs) 126! 127 do i=1,3 128 do j=1,3 129 xs(i,j)=0.d0 130 do k=1,15 131 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 132 enddo 133 enddo 134 enddo 135! 136! computation of the jacobian determinant 137! 138 xsj=xs(1,1)*(xs(2,2)*xs(3,3)-xs(2,3)*xs(3,2)) 139 & -xs(1,2)*(xs(2,1)*xs(3,3)-xs(2,3)*xs(3,1)) 140 & +xs(1,3)*(xs(2,1)*xs(3,2)-xs(2,2)*xs(3,1)) 141! 142 if(iflag.eq.2) return 143! 144! computation of the global derivative of the local coordinates 145! (xsi) (inversion of xs) 146! 147 xsi(1,1)=(xs(2,2)*xs(3,3)-xs(3,2)*xs(2,3))/xsj 148 xsi(1,2)=(xs(1,3)*xs(3,2)-xs(1,2)*xs(3,3))/xsj 149 xsi(1,3)=(xs(1,2)*xs(2,3)-xs(2,2)*xs(1,3))/xsj 150 xsi(2,1)=(xs(2,3)*xs(3,1)-xs(2,1)*xs(3,3))/xsj 151 xsi(2,2)=(xs(1,1)*xs(3,3)-xs(3,1)*xs(1,3))/xsj 152 xsi(2,3)=(xs(1,3)*xs(2,1)-xs(1,1)*xs(2,3))/xsj 153 xsi(3,1)=(xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2))/xsj 154 xsi(3,2)=(xs(1,2)*xs(3,1)-xs(1,1)*xs(3,2))/xsj 155 xsi(3,3)=(xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2))/xsj 156! 157! computation of the global derivatives of the shape functions 158! 159 do k=1,15 160 do j=1,3 161 sh(j)=shp(1,k)*xsi(1,j)+shp(2,k)*xsi(2,j)+shp(3,k)*xsi(3,j) 162 enddo 163 do j=1,3 164 shp(j,k)=sh(j) 165 enddo 166 enddo 167! 168 return 169 end 170