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