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 shape funciton dual shape funciton for 21! quad-quad mortar method \f$ shape(\xi,\eta) \f$ 22! see phd-thesis Sitzmann Chapter 4.1. 23! 24! Author: Saskia Sitzmann 25! 26! [in] xi xi-coordinate 27! [in] et eta-coordinate 28! [in] xl local node coordinates 29! [out] xsj jacobian vector 30! [out] xs local derivative of the global coordinates 31! [out] shp evaluated shape functions and derivatives 32! [in] ns current slave face 33! [in] pslavdual (:,i)coefficients \f$ \alpha_{ij}\f$, 34! \f$ 1,j=1,..8\f$ for dual shape functions for face i 35! [in] iflag flag indicating what to compute 36! 37 subroutine dualshape8qtilde(xi,et,xl,xsj,xs,shp,ns, 38 & pslavdual,iflag) 39! 40! shape functions and derivatives for a 8-node quadratic 41! isoparametric quadrilateral element. -1<=xi,et<=1 42! 43! iflag=2: calculate the value of the shape functions, 44! their derivatives w.r.t. the local coordinates 45! and the Jacobian vector (local normal to the 46! surface) 47! iflag=3: calculate the value of the shape functions, the 48! value of their derivatives w.r.t. the global 49! coordinates and the Jacobian vector (local normal 50! to the surface) 51! 52 implicit none 53! 54 integer i,j,k,iflag,ns 55! 56 real*8 shp(7,8),xs(3,7),xsi(2,3),xl(3,8),sh(3),xsj(3), 57 & shpold(7,8),alpha,pslavdual(64,*) 58! 59 real*8 xi,et 60! 61! 62! 63! shape functions and their glocal derivatives for an element 64! described with two local parameters and three global ones. 65! 66 alpha=1.0/5.0 67! 68! standard shape functions 69! 70! shape functions 71! 72 shpold(4,1)=(1.d0-xi)*(1.d0-et)*(-xi-et-1.d0)/4.d0 73 shpold(4,2)=(1.d0+xi)*(1.d0-et)*(xi-et-1.d0)/4.d0 74 shpold(4,3)=(1.d0+xi)*(1.d0+et)*(xi+et-1.d0)/4.d0 75 shpold(4,4)=(1.d0-xi)*(1.d0+et)*(-xi+et-1.d0)/4.d0 76 shpold(4,5)=(1.d0-xi*xi)*(1.d0-et)/2.d0 77 shpold(4,6)=(1.d0+xi)*(1.d0-et*et)/2.d0 78 shpold(4,7)=(1.d0-xi*xi)*(1.d0+et)/2.d0 79 shpold(4,8)=(1.d0-xi)*(1.d0-et*et)/2.d0 80! 81! transformed standard shape functions 82! 83 shp(3,1)=1.0*shpold(4,1)+alpha*shpold(4,5)+alpha*shpold(4,8) 84 shp(3,2)=1.0*shpold(4,2)+alpha*shpold(4,5)+alpha*shpold(4,6) 85 shp(3,3)=1.0*shpold(4,3)+alpha*shpold(4,6)+alpha*shpold(4,7) 86 shp(3,4)=1.0*shpold(4,4)+alpha*shpold(4,7)+alpha*shpold(4,8) 87 shp(3,5)=(1.0-2.0*alpha)*shpold(4,5) 88 shp(3,6)=(1.0-2.0*alpha)*shpold(4,6) 89 shp(3,7)=(1.0-2.0*alpha)*shpold(4,7) 90 shp(3,8)=(1.0-2.0*alpha)*shpold(4,8) 91! 92! Caution: derivatives and exspecially jacobian for untransformed 93! basis functions are given 94! needed for consistent integration 95! 96! local derivatives of the shape functions: xi-derivative 97! 98 shp(1,1)=(1.d0-et)*(2.d0*xi+et)/4.d0 99 shp(1,2)=(1.d0-et)*(2.d0*xi-et)/4.d0 100 shp(1,3)=(1.d0+et)*(2.d0*xi+et)/4.d0 101 shp(1,4)=(1.d0+et)*(2.d0*xi-et)/4.d0 102 shp(1,5)=-xi*(1.d0-et) 103 shp(1,6)=(1.d0-et*et)/2.d0 104 shp(1,7)=-xi*(1.d0+et) 105 shp(1,8)=-(1.d0-et*et)/2.d0 106! 107! local derivatives of the shape functions: eta-derivative 108! 109 shp(2,1)=(1.d0-xi)*(2.d0*et+xi)/4.d0 110 shp(2,2)=(1.d0+xi)*(2.d0*et-xi)/4.d0 111 shp(2,3)=(1.d0+xi)*(2.d0*et+xi)/4.d0 112 shp(2,4)=(1.d0-xi)*(2.d0*et-xi)/4.d0 113 shp(2,5)=-(1.d0-xi*xi)/2.d0 114 shp(2,6)=-et*(1.d0+xi) 115 shp(2,7)=(1.d0-xi*xi)/2.d0 116 shp(2,8)=-et*(1.d0-xi) 117! 118! Dual shape functions with transformed std basis 119! 120 do i=1,8 121 shp(4,i)=0.0 122 do j=1,8 123 shp(4,i)=shp(4,i)+pslavdual((i-1)*8+j,ns)*shp(3,j) 124 enddo 125 enddo 126! 127! computation of the local derivative of the global coordinates 128! (xs) 129! 130 do i=1,3 131 do j=1,2 132 xs(i,j)=0.d0 133 do k=1,8 134 xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k) 135 enddo 136 enddo 137 enddo 138! 139! computation of the jacobian vector 140! 141 xsj(1)=xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2) 142 xsj(2)=xs(1,2)*xs(3,1)-xs(3,2)*xs(1,1) 143 xsj(3)=xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2) 144! 145 if(iflag.eq.2) return 146! 147! computation of the global derivative of the local coordinates 148! (xsi) (inversion of xs) 149! 150 if(dabs(xsj(3)).gt.1.d-10) then 151 xsi(1,1)=xs(2,2)/xsj(3) 152 xsi(2,2)=xs(1,1)/xsj(3) 153 xsi(1,2)=-xs(1,2)/xsj(3) 154 xsi(2,1)=-xs(2,1)/xsj(3) 155 if(dabs(xsj(2)).gt.1.d-10) then 156 xsi(2,3)=xs(1,1)/xsj(2) 157 xsi(1,3)=-xs(1,2)/xsj(2) 158 elseif(dabs(xsj(1)).gt.1.d-10) then 159 xsi(2,3)=xs(2,1)/xsj(1) 160 xsi(1,3)=-xs(2,2)/xsj(1) 161 else 162 xsi(2,3)=0.d0 163 xsi(1,3)=0.d0 164 endif 165 elseif(dabs(xsj(2)).gt.1.d-10) then 166 xsi(1,1)=xs(3,2)/xsj(2) 167 xsi(2,3)=xs(1,1)/xsj(2) 168 xsi(1,3)=-xs(1,2)/xsj(2) 169 xsi(2,1)=-xs(3,1)/xsj(2) 170 if(dabs(xsj(1)).gt.1.d-10) then 171 xsi(1,2)=xs(3,2)/xsj(1) 172 xsi(2,2)=-xs(3,1)/xsj(1) 173 else 174 xsi(1,2)=0.d0 175 xsi(2,2)=0.d0 176 endif 177 else 178 xsi(1,2)=xs(3,2)/xsj(1) 179 xsi(2,3)=xs(2,1)/xsj(1) 180 xsi(1,3)=-xs(2,2)/xsj(1) 181 xsi(2,2)=-xs(3,1)/xsj(1) 182 xsi(1,1)=0.d0 183 xsi(2,1)=0.d0 184 endif 185! 186! computation of the global derivatives of the shape functions 187! 188 do k=1,8 189 do j=1,3 190 sh(j)=shp(1,k)*xsi(1,j)+shp(2,k)*xsi(2,j) 191 enddo 192 do j=1,3 193 shp(j,k)=sh(j) 194 enddo 195 enddo 196! 197 return 198 end 199 200