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 map3dto1d2d_v(yn,ipkon,inum,kon,lakon,nfield,nk, 20 & ne,nactdof) 21! 22! interpolates basic degree of freedom nodal values 23! (displacements, temperatures) to 1d/2d nodal locations 24! 25 implicit none 26! 27 logical quadratic 28! 29 character*8 lakon(*),lakonl 30! 31 integer ipkon(*),inum(*),kon(*),ne,indexe,nfield,nk,i,j,k,l, 32 & node3(8,3),node6(3,6),node8(3,8),node2d,node3d,indexe2d,ne1d2d, 33 & node3m(8,3),iflag,nactdof(nfield,*),jmax 34! 35 real*8 yn(nfield,*),ratioe(3) 36! 37! 38! 39 include "gauss.f" 40! 41 data node3 /1,4,8,5,12,20,16,17,9,11,15,13, 42 & 0,0,0,0,2,3,7,6,10,19,14,18/ 43 data node3m /1,5,8,4,17,16,20,12, 44 & 0,0,0,0,0,0,0,0, 45 & 3,7,6,2,19,14,18,10/ 46 data node6 /1,13,4,2,14,5,3,15,6,7,0,10,8,0,11,9,0,12/ 47 data node8 /1,17,5,2,18,6,3,19,7,4,20,8,9,0,13,10,0,14, 48 & 11,0,15,12,0,16/ 49 data ratioe /0.16666666666667d0,0.66666666666666d0, 50 & 0.16666666666667d0/ 51 data iflag /2/ 52! 53! removing any results in 1d/2d nodes 54! 55 ne1d2d=0 56! 57 do i=1,ne 58! 59 if(ipkon(i).lt.0) cycle 60 lakonl=lakon(i) 61 if((lakonl(7:7).eq.' ').or.(lakonl(7:7).eq.'I').or. 62 & (lakonl(1:1).ne.'C')) cycle 63 ne1d2d=1 64 indexe=ipkon(i) 65c! 66c! inactivating the 3d expansion nodes of 1d/2d elements 67c! 68c do j=1,20 69c inum(kon(indexe+j))=0 70c enddo 71! 72 if((lakonl(4:5).eq.'15').or.(lakonl(4:4).eq.'6')) then 73 if(lakonl(4:5).eq.'15') then 74 indexe2d=indexe+15 75 jmax=6 76 else 77 indexe2d=indexe+6 78 jmax=3 79 endif 80 do j=1,jmax 81 node2d=kon(indexe2d+j) 82 inum(node2d)=0 83 do k=1,nfield 84 if(nactdof(k,node2d).le.0) yn(k,node2d)=0.d0 85 enddo 86 enddo 87 elseif(lakonl(7:7).eq.'B') then 88 if(lakonl(4:5).eq.'8I') then 89 indexe2d=indexe+11 90 jmax=2 91 elseif(lakonl(4:5).eq.'8R') then 92 indexe2d=indexe+8 93 jmax=2 94 elseif(lakonl(4:5).eq.'20') then 95 indexe2d=indexe+20 96 jmax=3 97 endif 98 do j=1,jmax 99 node2d=kon(indexe2d+j) 100 inum(node2d)=0 101 do k=1,nfield 102 if(nactdof(k,node2d).le.0) yn(k,node2d)=0.d0 103 enddo 104 enddo 105 else 106 if(lakonl(4:5).eq.'8I') then 107 indexe2d=indexe+11 108 jmax=4 109 elseif((lakonl(4:5).eq.'8R').or.(lakonl(4:5).eq.'8 ')) then 110 indexe2d=indexe+8 111 jmax=4 112 elseif(lakonl(4:5).eq.'20') then 113 indexe2d=indexe+20 114 jmax=8 115 endif 116 do j=1,jmax 117 node2d=kon(indexe2d+j) 118 inum(node2d)=0 119 do k=1,nfield 120 if(nactdof(k,node2d).le.0) yn(k,node2d)=0.d0 121 enddo 122 enddo 123 endif 124! 125! inactivating the 3d expansion nodes of 1d/2d elements 126! 127 do j=1,indexe2d-indexe 128 inum(kon(indexe+j))=0 129 enddo 130! 131 enddo 132! 133! if no 1d/2d elements return 134! 135 if(ne1d2d.eq.0) return 136! 137! interpolation of 3d results on 1d/2d nodes 138! 139 do i=1,ne 140! 141 if(ipkon(i).lt.0) cycle 142 lakonl=lakon(i) 143 if((lakonl(7:7).eq.' ').or.(lakonl(7:7).eq.'I').or. 144 & (lakonl(1:1).ne.'C')) cycle 145 indexe=ipkon(i) 146! 147! check whether linear or quadratic element 148! 149 if((lakonl(4:4).eq.'6').or.(lakonl(4:4).eq.'8')) then 150 quadratic=.false. 151 else 152 quadratic=.true. 153 endif 154! 155 if((lakonl(4:5).eq.'15').or.(lakonl(4:4).eq.'6')) then 156 if(lakonl(4:5).eq.'15') then 157 indexe2d=indexe+15 158 jmax=6 159 else 160 indexe2d=indexe+6 161 jmax=3 162 endif 163 do j=1,jmax 164 node2d=kon(indexe2d+j) 165 inum(node2d)=inum(node2d)-1 166! 167! taking the mean across the thickness 168! 169 if((j.le.3).and.(quadratic)) then 170! 171! end nodes: weights 1/6,2/3 and 1/6 172! 173 do l=1,3 174 node3d=kon(indexe+node6(l,j)) 175 do k=1,nfield 176 if(nactdof(k,node2d).le.0) yn(k,node2d)= 177 & yn(k,node2d)+yn(k,node3d)*ratioe(l) 178 enddo 179 enddo 180 else 181! 182! middle nodes: weights 1/2,1/2 183! 184 do l=1,3,2 185 node3d=kon(indexe+node6(l,j)) 186 do k=1,nfield 187 if(nactdof(k,node2d).le.0) yn(k,node2d)= 188 & yn(k,node2d)+yn(k,node3d)/2.d0 189 enddo 190 enddo 191 endif 192 enddo 193 elseif(lakonl(7:7).eq.'B') then 194 if(lakonl(4:5).eq.'8I') then 195 indexe2d=indexe+11 196 jmax=2 197 elseif(lakonl(4:5).eq.'8R') then 198 indexe2d=indexe+8 199 jmax=2 200 elseif(lakonl(4:5).eq.'20') then 201 indexe2d=indexe+20 202 jmax=3 203 endif 204! 205! mean values for beam elements 206! 207 do j=1,jmax 208 node2d=kon(indexe2d+j) 209! 210! mean value of vertex values 211! 212 do l=1,4 213 inum(node2d)=inum(node2d)-1 214 if(quadratic) then 215 node3d=kon(indexe+node3(l,j)) 216 else 217 node3d=kon(indexe+node3(l,2*j-1)) 218 endif 219 do k=1,nfield 220 if(nactdof(k,node2d).le.0) yn(k,node2d)= 221 & yn(k,node2d)+yn(k,node3d) 222 enddo 223 enddo 224 enddo 225 else 226 if(lakonl(4:5).eq.'8I') then 227 indexe2d=indexe+11 228 jmax=4 229 elseif((lakonl(4:5).eq.'8R').or.(lakonl(4:5).eq.'8 ')) then 230 indexe2d=indexe+8 231 jmax=4 232 elseif(lakonl(4:5).eq.'20') then 233 indexe2d=indexe+20 234 jmax=8 235 endif 236 do j=1,jmax 237 node2d=kon(indexe2d+j) 238 inum(node2d)=inum(node2d)-1 239! 240! taking the mean across the thickness 241! 242 if((j.le.4).and.(quadratic)) then 243! 244! end nodes: weights 1/6,2/3 and 1/6 245! 246 do l=1,3 247 node3d=kon(indexe+node8(l,j)) 248 do k=1,nfield 249 if(nactdof(k,node2d).le.0) yn(k,node2d)= 250 & yn(k,node2d)+yn(k,node3d)*ratioe(l) 251 enddo 252 enddo 253 else 254! 255! middle nodes: weights 1/2,1/2 256! 257 do l=1,3,2 258 node3d=kon(indexe+node8(l,j)) 259 do k=1,nfield 260 if(nactdof(k,node2d).le.0) yn(k,node2d)= 261 & yn(k,node2d)+yn(k,node3d)/2.d0 262 enddo 263 enddo 264 endif 265 enddo 266 endif 267! 268 enddo 269! 270! taking the mean of nodal contributions coming from different 271! elements having the node in common 272! 273 do i=1,nk 274 if(inum(i).lt.0) then 275 inum(i)=-inum(i) 276 do j=1,nfield 277 if(nactdof(j,i).le.0) yn(j,i)=yn(j,i)/inum(i) 278 enddo 279 endif 280 enddo 281! 282 return 283 end 284 285