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 shape8q(xi,et,xl,xsj,xs,shp,iflag) 20! 21! shape functions and derivatives for a 8-node quadratic 22! isoparametric quadrilateral element. -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, 26! their derivatives w.r.t. the local coordinates 27! and the Jacobian vector (local normal to the 28! surface) 29! iflag=3: calculate the value of the shape functions, the 30! value of their derivatives w.r.t. the global 31! coordinates and the Jacobian vector (local normal 32! to the surface) 33! iflag=4: calculate the value of the shape functions, the 34! value of their 1st and 2nd order derivatives 35! w.r.t. the local coordinates, the Jacobian vector 36! (local normal to the surface) 37! iflag=5: calculate the value of the shape functions and 38! their derivatives w.r.t. the local coordinates 39! 40 implicit none 41! 42 integer i,j,k,iflag 43! 44 real*8 shp(7,8),xs(3,7),xsi(2,3),xl(3,8),sh(3),xsj(3),xi,et, 45 & xip,xim,xim2,etp,etm,etm2,xipet,ximet 46! 47! shape functions and their glocal derivatives for an element 48! described with two local parameters and three global ones. 49! 50 xip=1.d0+xi 51 xim=1.d0-xi 52 xim2=xip*xim 53! 54 etp=1.d0+et 55 etm=1.d0-et 56 etm2=etp*etm 57! 58 xipet=xi+et 59 ximet=xi-et 60! 61! shape functions 62! 63 shp(4,1)=xim*etm*(-xipet-1.d0)/4.d0 64 shp(4,2)=xip*etm*(ximet-1.d0)/4.d0 65 shp(4,3)=xip*etp*(xipet-1.d0)/4.d0 66 shp(4,4)=xim*etp*(-ximet-1.d0)/4.d0 67 shp(4,5)=xim2*etm/2.d0 68 shp(4,6)=xip*etm2/2.d0 69 shp(4,7)=xim2*etp/2.d0 70 shp(4,8)=xim*etm2/2.d0 71! 72 if(iflag.eq.1) return 73! 74! local derivatives of the shape functions: xi-derivative 75! 76 shp(1,1)=etm*(xi+xipet)/4.d0 77 shp(1,2)=etm*(xi+ximet)/4.d0 78 shp(1,3)=etp*(xi+xipet)/4.d0 79 shp(1,4)=etp*(xi+ximet)/4.d0 80 shp(1,5)=-xi*etm 81 shp(1,6)=etm2/2.d0 82 shp(1,7)=-xi*etp 83 shp(1,8)=-etm2/2.d0 84! 85! local derivatives of the shape functions: eta-derivative 86! 87 shp(2,1)=xim*(et+xipet)/4.d0 88 shp(2,2)=xip*(et-ximet)/4.d0 89 shp(2,3)=xip*(et+xipet)/4.d0 90 shp(2,4)=xim*(et-ximet)/4.d0 91 shp(2,5)=-xim2/2.d0 92 shp(2,6)=-et*xip 93 shp(2,7)=xim2/2.d0 94 shp(2,8)=-et*xim 95! 96 if(iflag.eq.5) return 97! 98! computation of the local derivative of the global coordinates 99! (xs) 100! 101 do i=1,3 102 do j=1,2 103 xs(i,j)=0.d0 104 do k=1,8 105 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 106 enddo 107 enddo 108 enddo 109! 110! computation of the jacobian vector 111! 112 xsj(1)=xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2) 113 xsj(2)=xs(1,2)*xs(3,1)-xs(3,2)*xs(1,1) 114 xsj(3)=xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2) 115! 116 if(iflag.eq.3) then 117! 118! computation of the global derivative of the local coordinates 119! (xsi) (inversion of xs) 120! 121 if(dabs(xsj(3)).gt.1.d-10) then 122 xsi(1,1)=xs(2,2)/xsj(3) 123 xsi(2,2)=xs(1,1)/xsj(3) 124 xsi(1,2)=-xs(1,2)/xsj(3) 125 xsi(2,1)=-xs(2,1)/xsj(3) 126 if(dabs(xsj(2)).gt.1.d-10) then 127 xsi(2,3)=xs(1,1)/(-xsj(2)) 128 xsi(1,3)=-xs(1,2)/(-xsj(2)) 129 elseif(dabs(xsj(1)).gt.1.d-10) then 130 xsi(2,3)=xs(2,1)/xsj(1) 131 xsi(1,3)=-xs(2,2)/xsj(1) 132 else 133 xsi(2,3)=0.d0 134 xsi(1,3)=0.d0 135 endif 136 elseif(dabs(xsj(2)).gt.1.d-10) then 137 xsi(1,1)=xs(3,2)/(-xsj(2)) 138 xsi(2,3)=xs(1,1)/(-xsj(2)) 139 xsi(1,3)=-xs(1,2)/(-xsj(2)) 140 xsi(2,1)=-xs(3,1)/(-xsj(2)) 141 if(dabs(xsj(1)).gt.1.d-10) then 142 xsi(1,2)=xs(3,2)/xsj(1) 143 xsi(2,2)=-xs(3,1)/xsj(1) 144 else 145 xsi(1,2)=0.d0 146 xsi(2,2)=0.d0 147 endif 148 else 149 xsi(1,2)=xs(3,2)/xsj(1) 150 xsi(2,3)=xs(2,1)/xsj(1) 151 xsi(1,3)=-xs(2,2)/xsj(1) 152 xsi(2,2)=-xs(3,1)/xsj(1) 153 xsi(1,1)=0.d0 154 xsi(2,1)=0.d0 155 endif 156! 157! computation of the global derivatives of the shape functions 158! 159 do k=1,8 160 do j=1,3 161 sh(j)=shp(1,k)*xsi(1,j)+shp(2,k)*xsi(2,j) 162 enddo 163 do j=1,3 164 shp(j,k)=sh(j) 165 enddo 166 enddo 167! 168 elseif(iflag.eq.4) then 169! 170! local 2nd order derivatives of the shape functions: xi,xi-derivative 171! 172 shp(5,1)=etm/2.d0 173 shp(5,2)=etm/2.d0 174 shp(5,3)=etp/2.d0 175 shp(5,4)=etp/2.d0 176 shp(5,5)=-etm 177 shp(5,6)=0.d0 178 shp(5,7)=-etp 179 shp(5,8)=0.d0 180! 181! local 2nd order derivatives of the shape functions: xi,eta-derivative 182! 183 shp(6,1)=(1.d0-2.d0*xipet)/4.d0 184 shp(6,2)=(-1.d0-2.d0*ximet)/4.d0 185 shp(6,3)=(1.d0+2.d0*xipet)/4.d0 186 shp(6,4)=(-1.d0+2.d0*ximet)/4.d0 187 shp(6,5)=xi 188 shp(6,6)=-et 189 shp(6,7)=-xi 190 shp(6,8)=et 191! 192! local 2nd order derivatives of the shape functions: eta,eta-derivative 193! 194 shp(7,1)=xim/2.d0 195 shp(7,2)=xip/2.d0 196 shp(7,3)=xip/2.d0 197 shp(7,4)=xim/2.d0 198 shp(7,5)=0.d0 199 shp(7,6)=-xip 200 shp(7,7)=0.d0 201 shp(7,8)=-xim 202! 203! computation of the local 2nd derivatives of the global coordinates 204! (xs) 205! 206 do i=1,3 207 do j=5,7 208 xs(i,j)=0.d0 209 do k=1,8 210 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 211 enddo 212 enddo 213 enddo 214 endif 215! 216 return 217 end 218