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! 20! function to evaluate transformed shape funciton \f$ shape(\xi,\eta) \f$ 21! for quad-lin mortar method, see phd-thesis Sitzmann Chapter 4.1. 22! 23! Author: Saskia Sitzmann 24! 25! [in] xi xi-coordinate 26! [in] et eta-coordinate 27! [in] xl local node coordinates 28! [out] xsj jacobian vector 29! [out] xs local derivative of the global coordinates 30! [out] shp evaluated shape functions and derivatives 31! [in] iflag flag indicating what to compute 32! 33 subroutine shape6tritilde_lin(xi,et,xl,xsj,xs,shp,iflag) 34! 35! iflag=1: calculate only the value of the shape functions 36! iflag=2: calculate the value of the shape functions, 37! their derivatives w.r.t. the local coordinates 38! and the Jacobian vector (local normal to the 39! surface) 40! iflag=3: calculate the value of the shape functions, the 41! value of their derivatives w.r.t. the global 42! coordinates and the Jacobian vector (local normal 43! to the surface) 44! iflag=4: calculate the value of the shape functions, the 45! value of their 1st and 2nd order derivatives 46! w.r.t. the local coordinates, the Jacobian vector 47! (local normal to the surface) 48! 49! shape functions and derivatives for a 6-node quadratic 50! isoparametric triangular element. 0<=xi,et<=1,xi+et<=1 51! 52 implicit none 53! 54 integer i,j,k,iflag 55! 56 real*8 shp(7,8),xs(3,7),xsi(2,3),xl(3,8),sh(3),xsj(3) 57! 58 real*8 xi,et 59! 60! 61! 62! shape functions and their glocal derivatives for an element 63! described with two local parameters and three global ones. 64! 65! shape functions 66! 67 shp(4,1)=1.d0-xi-et 68 shp(4,2)=xi 69 shp(4,3)=et 70 shp(4,4)=4.d0*xi*(1.d0-xi-et) 71 shp(4,5)=4.d0*xi*et 72 shp(4,6)=4.d0*et*(1.d0-xi-et) 73! 74! Caution: derivatives and exspecially jacobian for untransformed 75! basis functions are given 76! needed for consistent integration 77! 78 if(iflag.eq.1) return 79! 80! local derivatives of the shape functions: xi-derivative 81! 82 shp(1,1)=4.d0*(xi+et)-3.d0 83 shp(1,2)=4.d0*xi-1.d0 84 shp(1,3)=0.d0 85 shp(1,4)=4.d0*(1.d0-2.d0*xi-et) 86 shp(1,5)=4.d0*et 87 shp(1,6)=-4.d0*et 88! 89! local derivatives of the shape functions: eta-derivative 90! 91 shp(2,1)=4.d0*(xi+et)-3.d0 92 shp(2,2)=0.d0 93 shp(2,3)=4.d0*et-1.d0 94 shp(2,4)=-4.d0*xi 95 shp(2,5)=4.d0*xi 96 shp(2,6)=4.d0*(1.d0-xi-2.d0*et) 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,6 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,6 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)=0.d0 173 shp(5,2)=0.d0 174 shp(5,3)=0.d0 175 shp(5,4)=-8.d0 176 shp(5,5)=0.d0 177 shp(5,6)=0.d0 178! 179! local 2nd order derivatives of the shape functions: xi,eta-derivative 180! 181 shp(6,1)=0.d0 182 shp(6,2)=0.d0 183 shp(6,3)=0.d0 184 shp(6,4)=-4.d0 185 shp(6,5)=4.d0 186 shp(6,6)=-4.d0 187! 188! local 2nd order derivatives of the shape functions: eta,eta-derivative 189! 190 shp(7,1)=0.d0 191 shp(7,2)=0.d0 192 shp(7,3)=0.d0 193 shp(7,4)=0.d0 194 shp(7,5)=0.d0 195 shp(7,6)=-8.d0 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,6 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 213