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 shape8hu(xi,et,ze,xl,xsj,shp,iflag) 20! 21! shape functions and derivatives for a 8-node linear isoparametric 22! solid element 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! author: Otto-Ernst Bernhardi 32! 33 implicit none 34! 35 integer i,j,k,iflag 36! 37 real*8 shp(4,23),xs(3,3),xsi(3,3),xl(3,23),sh(3),xsi0(3,3) 38! 39 real*8 xi,et,ze,xsj,omg,omh,omr,opg,oph,opr 40! 41! 42! 43 if(iflag.gt.2) then 44! 45! local derivatives at center point: xi-derivative 46! 47 shp(1,1)=-1.0d0/8.d0 48 shp(1,2)=1.0d0/8.d0 49 shp(1,3)=1.0d0/8.d0 50 shp(1,4)=-1.0d0/8.d0 51 shp(1,5)=-1.0d0/8.d0 52 shp(1,6)=1.0d0/8.d0 53 shp(1,7)=1.0d0/8.d0 54 shp(1,8)=-1.0d0/8.d0 55! 56! local derivatives at center point: eta-derivative 57! 58 shp(2,1)=-1.0d0/8.d0 59 shp(2,2)=-1.0d0/8.d0 60 shp(2,3)=1.0d0/8.d0 61 shp(2,4)=1.0d0/8.d0 62 shp(2,5)=-1.0d0/8.d0 63 shp(2,6)=-1.0d0/8.d0 64 shp(2,7)=1.0d0/8.d0 65 shp(2,8)=1.0d0/8.d0 66! 67! local derivatives at center point: zeta-derivative 68! 69 shp(3,1)=-1.0d0/8.d0 70 shp(3,2)=-1.0d0/8.d0 71 shp(3,3)=-1.0d0/8.d0 72 shp(3,4)=-1.0d0/8.d0 73 shp(3,5)=1.0d0/8.d0 74 shp(3,6)=1.0d0/8.d0 75 shp(3,7)=1.0d0/8.d0 76 shp(3,8)=1.0d0/8.d0 77! 78! computation of the local derivative of the global coordinates 79! (xs) at center point of element. 80! 81 do i=1,3 82 do j=1,3 83 xs(i,j)=0.d0 84 do k=1,8 85 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 86 enddo 87 enddo 88 enddo 89! 90! computation of the jacobian determinant at center point 91! 92 xsj=xs(1,1)*(xs(2,2)*xs(3,3)-xs(2,3)*xs(3,2)) 93 & -xs(1,2)*(xs(2,1)*xs(3,3)-xs(2,3)*xs(3,1)) 94 & +xs(1,3)*(xs(2,1)*xs(3,2)-xs(2,2)*xs(3,1)) 95! 96! computation of the global derivative of the local coordinates 97! at center point of element. 98! 99 xsi0(1,1)=(xs(2,2)*xs(3,3)-xs(3,2)*xs(2,3))/xsj 100 xsi0(1,2)=(xs(1,3)*xs(3,2)-xs(1,2)*xs(3,3))/xsj 101 xsi0(1,3)=(xs(1,2)*xs(2,3)-xs(2,2)*xs(1,3))/xsj 102 xsi0(2,1)=(xs(2,3)*xs(3,1)-xs(2,1)*xs(3,3))/xsj 103 xsi0(2,2)=(xs(1,1)*xs(3,3)-xs(3,1)*xs(1,3))/xsj 104 xsi0(2,3)=(xs(1,3)*xs(2,1)-xs(1,1)*xs(2,3))/xsj 105 xsi0(3,1)=(xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2))/xsj 106 xsi0(3,2)=(xs(1,2)*xs(3,1)-xs(1,1)*xs(3,2))/xsj 107 xsi0(3,3)=(xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2))/xsj 108! 109 endif 110! 111! shape functions and their global derivatives 112! 113 omg=1.d0-xi 114 omh=1.d0-et 115 omr=1.d0-ze 116 opg=1.d0+xi 117 oph=1.d0+et 118 opr=1.d0+ze 119! 120! shape functions 121! 122 shp(4, 1)=omg*omh*omr/8.d0 123 shp(4, 2)=opg*omh*omr/8.d0 124 shp(4, 3)=opg*oph*omr/8.d0 125 shp(4, 4)=omg*oph*omr/8.d0 126 shp(4, 5)=omg*omh*opr/8.d0 127 shp(4, 6)=opg*omh*opr/8.d0 128 shp(4, 7)=opg*oph*opr/8.d0 129 shp(4, 8)=omg*oph*opr/8.d0 130c 131c change on 190315: set shape functions to 132c zero in order to obtain convergence in 133c contact calculations with c3d8i 134c 135 shp(4, 9)=0.d0 136 shp(4,10)=0.d0 137 shp(4,11)=0.d0 138c shp(4, 9)=1.0d0-xi*xi 139c shp(4,10)=1.0d0-et*et 140c shp(4,11)=1.0d0-ze*ze 141! 142 if(iflag.eq.1) return 143! 144! local derivatives of the shape functions: xi-derivative 145! 146 shp(1, 1)=-omh*omr/8.d0 147 shp(1, 2)=omh*omr/8.d0 148 shp(1, 3)=oph*omr/8.d0 149 shp(1, 4)=-oph*omr/8.d0 150 shp(1, 5)=-omh*opr/8.d0 151 shp(1, 6)=omh*opr/8.d0 152 shp(1, 7)=oph*opr/8.d0 153 shp(1, 8)=-oph*opr/8.d0 154 shp(1, 9)=-2.d0*xi 155 shp(1,10)=0.d0 156 shp(1,11)=0.d0 157! 158! local derivatives of the shape functions: eta-derivative 159! 160 shp(2, 1)=-omg*omr/8.d0 161 shp(2, 2)=-opg*omr/8.d0 162 shp(2, 3)=opg*omr/8.d0 163 shp(2, 4)=omg*omr/8.d0 164 shp(2, 5)=-omg*opr/8.d0 165 shp(2, 6)=-opg*opr/8.d0 166 shp(2, 7)=opg*opr/8.d0 167 shp(2, 8)=omg*opr/8.d0 168 shp(2, 9)=0.d0 169 shp(2,10)=-2.d0*et 170 shp(2,11)=0.d0 171! 172! local derivatives of the shape functions: zeta-derivative 173! 174 shp(3, 1)=-omg*omh/8.d0 175 shp(3, 2)=-opg*omh/8.d0 176 shp(3, 3)=-opg*oph/8.d0 177 shp(3, 4)=-omg*oph/8.d0 178 shp(3, 5)=omg*omh/8.d0 179 shp(3, 6)=opg*omh/8.d0 180 shp(3, 7)=opg*oph/8.d0 181 shp(3, 8)=omg*oph/8.d0 182 shp(3, 9)=0.d0 183 shp(3,10)=0.d0 184 shp(3,11)=-2.d0*ze 185! 186! computation of the local derivative of the global coordinates 187! (xs). Incompatible modes are not included here. 188! 189 do i=1,3 190 do j=1,3 191 xs(i,j)=0.d0 192 do k=1,8 193 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 194 enddo 195 enddo 196 enddo 197! 198! computation of the jacobian determinant 199! 200 xsj=xs(1,1)*(xs(2,2)*xs(3,3)-xs(2,3)*xs(3,2)) 201 & -xs(1,2)*(xs(2,1)*xs(3,3)-xs(2,3)*xs(3,1)) 202 & +xs(1,3)*(xs(2,1)*xs(3,2)-xs(2,2)*xs(3,1)) 203! 204 if(iflag.eq.2) return 205! 206! computation of the global derivative of the local coordinates 207! (xsi) (inversion of xs) 208! 209 xsi(1,1)=(xs(2,2)*xs(3,3)-xs(3,2)*xs(2,3))/xsj 210 xsi(1,2)=(xs(1,3)*xs(3,2)-xs(1,2)*xs(3,3))/xsj 211 xsi(1,3)=(xs(1,2)*xs(2,3)-xs(2,2)*xs(1,3))/xsj 212 xsi(2,1)=(xs(2,3)*xs(3,1)-xs(2,1)*xs(3,3))/xsj 213 xsi(2,2)=(xs(1,1)*xs(3,3)-xs(3,1)*xs(1,3))/xsj 214 xsi(2,3)=(xs(1,3)*xs(2,1)-xs(1,1)*xs(2,3))/xsj 215 xsi(3,1)=(xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2))/xsj 216 xsi(3,2)=(xs(1,2)*xs(3,1)-xs(1,1)*xs(3,2))/xsj 217 xsi(3,3)=(xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2))/xsj 218! 219! computation of the global derivatives of the shape functions 220! 221 do k=1,8 222 do j=1,3 223 sh(j)=shp(1,k)*xsi(1,j)+shp(2,k)*xsi(2,j)+shp(3,k)*xsi(3,j) 224 enddo 225 do j=1,3 226 shp(j,k)=sh(j) 227 enddo 228 enddo 229 do k=9,11 230 do j=1,3 231 sh(j)=shp(1,k)*xsi0(1,j)+shp(2,k)*xsi0(2,j)+shp(3,k)*xsi0(3,j) 232 enddo 233 do j=1,3 234 shp(j,k)=sh(j) 235 enddo 236 enddo 237! 238 return 239 end 240