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 map3dtolayer(yn,ipkon,kon,lakon,nfield, 20 & ne,co,ielmat,mi) 21! 22! interpolates 3d field nodal values to nodal values in the 23! layers of composite materials 24! 25 implicit none 26! 27 character*8 lakon(*) 28! 29 integer ipkon(*),kon(*),ne,nfield,i,j,k,l,mi(*),nelem, 30 & ielmat(mi(3),*),nlayer,iflag,nbot20(8),nmid20(8),ntop20(8), 31 & nbot15(6),nmid15(6),ntop15(6),ibot,itop,nope,nopes,nopeexp, 32 & indexe,konl(20),imid,node,nodebot,nodetop 33! 34 real*8 yn(nfield,*),shp(4,20),xsj(3),co(3,*),xl(3,20), 35 & xi,et,ze,dd,dt,xi20(8),et20(8),xi15(6),et15(6) 36! 37! 38! 39 data iflag /1/ 40! 41 data nbot20 /1,2,3,4,9,10,11,12/ 42 data nmid20 /17,18,19,20,0,0,0,0/ 43 data ntop20 /5,6,7,8,13,14,15,16/ 44! 45 data nbot15 /1,2,3,7,8,9/ 46 data nmid15 /13,14,15,0,0,0/ 47 data ntop15 /4,5,6,10,11,12/ 48! 49 data xi20 /-1.d0,1.d0,1.d0,-1.d0,0.d0,1.d0,0.d0,-1.d0/ 50 data et20 /-1.d0,-1.d0,1.d0,1.d0,-1.d0,0.d0,1.d0,0.d0/ 51! 52 data xi15 /0.d0,1.d0,0.d0,0.5d0,0.5d0,0.d0/ 53 data et15 /0.d0,0.d0,1.d0,0.d0,0.5d0,0.5d0/ 54! 55 include "gauss.f" 56! 57 do nelem=1,ne 58 if(lakon(nelem)(7:8).eq.'LC') then 59! 60! composite materials 61! 62! determining the number of layers 63! 64 nlayer=0 65 do k=1,mi(3) 66 if(ielmat(k,nelem).ne.0) then 67 nlayer=nlayer+1 68 endif 69 enddo 70! 71 indexe=ipkon(nelem) 72! 73 if(lakon(nelem)(4:5).eq.'20') then 74 nope=20 75 nopes=8 76 nopeexp=28 77 elseif(lakon(nelem)(4:5).eq.'15') then 78 nope=15 79 nopes=6 80 nopeexp=21 81 endif 82! 83 do i=1,nope 84 konl(i)=kon(indexe+i) 85 do j=1,3 86 xl(j,i)=co(j,konl(i)) 87 enddo 88 enddo 89! 90 do i=1,nopes 91 if(lakon(nelem)(4:5).eq.'20') then 92 ibot=nbot20(i) 93 imid=nmid20(i) 94 itop=ntop20(i) 95 xi=xi20(i) 96 et=et20(i) 97 else 98 ibot=nbot15(i) 99 imid=nmid15(i) 100 itop=ntop15(i) 101 xi=xi15(i) 102 et=et15(i) 103 endif 104! 105 nodebot=konl(ibot) 106 nodetop=konl(itop) 107 dd=sqrt((co(1,nodebot)-co(1,nodetop))**2+ 108 & (co(2,nodebot)-co(2,nodetop))**2+ 109 & (co(3,nodebot)-co(3,nodetop))**2) 110 do j=0,nlayer-1 111! 112! bottom node 113! 114 node=kon(indexe+nopeexp+j*nope+ibot) 115 dt=sqrt((co(1,nodebot)-co(1,node))**2+ 116 & (co(2,nodebot)-co(2,node))**2+ 117 & (co(3,nodebot)-co(3,node))**2) 118 ze=2.d0*dt/dd-1.d0 119c write(*,*) 'map3dtolayer',node,xi,et,ze 120! 121! determining the value of the shape functions 122! 123 if(lakon(nelem)(4:5).eq.'20') then 124 call shape20h(xi,et,ze,xl,xsj,shp,iflag) 125 elseif(lakon(nelem)(4:5).eq.'15') then 126 call shape15w(xi,et,ze,xl,xsj,shp,iflag) 127 endif 128! 129 do k=1,nfield 130 yn(k,node)=0.d0 131 do l=1,nope 132 yn(k,node)=yn(k,node)+ 133 & shp(4,l)*yn(k,konl(l)) 134c write(*,*) 'xi',xi,et,ze,node 135c write(*,*) l,shp(4,l),yn(k,konl(l)) 136 enddo 137 enddo 138! 139! top node 140! 141 node=kon(indexe+nopeexp+j*nope+itop) 142 dt=sqrt((co(1,nodebot)-co(1,node))**2+ 143 & (co(2,nodebot)-co(2,node))**2+ 144 & (co(3,nodebot)-co(3,node))**2) 145 ze=2.d0*dt/dd-1.d0 146c write(*,*) 'map3dtolayer',node,xi,et,ze 147! 148! determining the value of the shape functions 149! 150 if(lakon(nelem)(4:5).eq.'20') then 151 call shape20h(xi,et,ze,xl,xsj,shp,iflag) 152 elseif(lakon(nelem)(4:5).eq.'15') then 153 call shape15w(xi,et,ze,xl,xsj,shp,iflag) 154 endif 155! 156 do k=1,nfield 157 yn(k,node)=0.d0 158 do l=1,nope 159 yn(k,node)=yn(k,node)+ 160 & shp(4,l)*yn(k,konl(l)) 161 enddo 162 enddo 163! 164! middle node, if any 165! 166 if(i.gt.nopes/2) cycle 167! 168 node=kon(indexe+nopeexp+j*nope+imid) 169 dt=sqrt((co(1,nodebot)-co(1,node))**2+ 170 & (co(2,nodebot)-co(2,node))**2+ 171 & (co(3,nodebot)-co(3,node))**2) 172 ze=2.d0*dt/dd-1.d0 173c write(*,*) 'map3dtolayer',node,xi,et,ze 174! 175! determining the value of the shape functions 176! 177 if(lakon(nelem)(4:5).eq.'20') then 178 call shape20h(xi,et,ze,xl,xsj,shp,iflag) 179 elseif(lakon(nelem)(4:5).eq.'15') then 180 call shape15w(xi,et,ze,xl,xsj,shp,iflag) 181 endif 182! 183 do k=1,nfield 184 yn(k,node)=0.d0 185 do l=1,nope 186 yn(k,node)=yn(k,node)+ 187 & shp(4,l)*yn(k,konl(l)) 188 enddo 189 enddo 190 enddo 191 enddo 192 endif 193 enddo 194! 195 return 196 end 197