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 for quad-lin mortar
21!     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 dualshape6tritilde_lin(xi,et,xl,xsj,xs,shp,ns,
38     &  pslavdual,iflag)
39!
40!     iflag=2: calculate the value of the shape functions,
41!              their derivatives w.r.t. the local coordinates
42!              and the Jacobian vector (local normal to the
43!              surface)
44!     iflag=3: calculate the value of the shape functions, the
45!              value of their derivatives w.r.t. the global
46!              coordinates and the Jacobian vector (local normal
47!              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,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      shpold(4,1)=2.d0*(0.5d0-xi-et)*(1.d0-xi-et)
71      shpold(4,2)=xi*(2.d0*xi-1.d0)
72      shpold(4,3)=et*(2.d0*et-1.d0)
73      shpold(4,4)=4.d0*xi*(1.d0-xi-et)
74      shpold(4,5)=4.d0*xi*et
75      shpold(4,6)=4.d0*et*(1.d0-xi-et)
76!
77!     transformed shape functions
78!
79      shp(3,1)=1.d0-xi-et
80      shp(3,2)=xi
81      shp(3,3)=et
82      shp(3,4)=4.d0*xi*(1.d0-xi-et)
83      shp(3,5)=4.d0*xi*et
84      shp(3,6)=4.d0*et*(1.d0-xi-et)
85!
86!     Caution: derivatives and exspecially jacobian for untransformed
87!     standard basis functions are given
88!     needed for consistent integration
89!
90!     local derivatives of the shape functions: xi-derivative
91!
92      shp(1,1)=4.d0*(xi+et)-3.d0
93      shp(1,2)=4.d0*xi-1.d0
94      shp(1,3)=0.d0
95      shp(1,4)=4.d0*(1.d0-2.d0*xi-et)
96      shp(1,5)=4.d0*et
97      shp(1,6)=-4.d0*et
98!
99!     local derivatives of the shape functions: eta-derivative
100!
101      shp(2,1)=4.d0*(xi+et)-3.d0
102      shp(2,2)=0.d0
103      shp(2,3)=4.d0*et-1.d0
104      shp(2,4)=-4.d0*xi
105      shp(2,5)=4.d0*xi
106      shp(2,6)=4.d0*(1.d0-xi-2.d0*et)
107!
108!     Dual shape functions with transformed std basis
109!
110      do i=1,6
111         shp(4,i)=0.0
112         do j=1,6
113            shp(4,i)=shp(4,i)+pslavdual((i-1)*8+j,ns)*shp(3,j)
114         enddo
115      enddo
116!
117!     computation of the local derivative of the global coordinates
118!     (xs)
119!
120      do i=1,3
121        do j=1,2
122           xs(i,j)=0.d0
123           do k=1,6
124              xs(i,j)=xs(i,j)+xl(i,k)*shp(j,k)
125           enddo
126        enddo
127      enddo
128!
129!     computation of the jacobian vector
130!
131      xsj(1)=xs(2,1)*xs(3,2)-xs(3,1)*xs(2,2)
132      xsj(2)=xs(1,2)*xs(3,1)-xs(3,2)*xs(1,1)
133      xsj(3)=xs(1,1)*xs(2,2)-xs(2,1)*xs(1,2)
134!
135      if(iflag.eq.2) return
136!
137!     computation of the global derivative of the local coordinates
138!     (xsi) (inversion of xs)
139!
140      xsi(1,1)=xs(2,2)/xsj(3)
141      xsi(2,1)=-xs(2,1)/xsj(3)
142      xsi(1,2)=-xs(1,2)/xsj(3)
143      xsi(2,2)=xs(1,1)/xsj(3)
144      xsi(1,3)=-xs(2,2)/xsj(1)
145      xsi(2,3)=xs(2,1)/xsj(1)
146!
147!     computation of the global derivatives of the shape functions
148!
149      do k=1,6
150         do j=1,3
151            sh(j)=shp(1,k)*xsi(1,j)+shp(2,k)*xsi(2,j)
152         enddo
153         do j=1,3
154            shp(j,k)=sh(j)
155         enddo
156      enddo
157!
158      return
159      end
160
161