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 extrapolateshell(yi,yn,ipkon,inum,kon,lakon,nfield,nk,
20     &  ne,mi,ndim,orab,ielorien,co,iorienloc,cflag,
21     &  ielmat,thicke,ielprop,prop,iflag)
22!
23!     extrapolates field values at the integration points to the
24!     nodes for user-defined shell elements
25!
26!     iflag=-1: NEG-position
27!     iflag=0: MID-position
28!     iflag=1: POS-position
29!
30      implicit none
31!
32      character*1 cflag
33      character*8 lakon(*)
34!
35      integer ipkon(*),inum(*),kon(*),mi(*),ne,iorienloc,nfield,nk,i,j,
36     &     ndim,ielorien(mi(3),*),ielmat(mi(3),*),ielprop(*),iflag
37!
38      real*8 yi(ndim,mi(1),*),yn(nfield,*),orab(7,*),co(3,*),prop(*),
39     &  thicke(mi(3),*)
40!
41      do i=1,nk
42         inum(i)=0
43      enddo
44!
45      do i=1,nk
46         do j=1,nfield
47            yn(j,i)=0.d0
48         enddo
49      enddo
50!
51      do i=1,ne
52!
53         if(ipkon(i).lt.0) cycle
54!
55         if(lakon(i)(1:4).eq.'US45') then
56            call extrapolateshell_us45(yi,yn,ipkon,inum,kon,lakon,
57     &           nfield,nk,ne,mi,ndim,orab,ielorien,co,iorienloc,cflag,
58     &           ielmat,thicke,ielprop,prop,i,iflag)
59         elseif(lakon(i)(1:3).eq.'US3') then
60            call extrapolateshell_us3(yi,yn,ipkon,inum,kon,lakon,
61     &           nfield,nk,ne,mi,ndim,orab,ielorien,co,iorienloc,cflag,
62     &           ielmat,thicke,ielprop,prop,i,iflag)
63         else
64            cycle
65         endif
66!
67      enddo
68!
69!     taking the mean of nodal contributions coming from different
70!     elements having the node in common
71!
72      do i=1,nk
73         if(inum(i).gt.0) then
74            do j=1,nfield
75               yn(j,i)=yn(j,i)/inum(i)
76            enddo
77         endif
78      enddo
79!
80      return
81      end
82