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 adjustcontactnodes(tieset,ntie,itietri,cg,straight, 20 & co,vold,xo,yo,zo,x,y,z,nx,ny,nz,istep,iinc,iit, 21 & mi,imastop,nslavnode,islavnode,set,nset,istartset, 22 & iendset,ialset,tietol,clearini,clearslavnode,itiefac, 23 & ipkon,kon,lakon,islavsurf) 24! 25! Adjusting contact nodes if requested by the user 26! 27 implicit none 28! 29 character*8 lakon(*) 30 character*81 tieset(3,*),slavset,set(*),noset 31! 32 integer ntie, 33 & itietri(2,ntie),node,neigh(1),kneigh,i,j,k,l,m,isol,iset, 34 & idummy,itri,ll,kflag,n,nx(*),ny(*),istep,iinc,mi(*), 35 & nz(*),nstart,iit,imastop(3,*),itriangle(100),ntriangle, 36 & ntriangle_,itriold,itrinew,id,nslavnode(*),islavnode(*), 37 & ipos,nset,istartset(*),iendset(*),ialset(*),iclear, 38 & istart,ilength,nope,nopes,nelems,jj,ifaces,itiefac(2,*), 39 & jfaces,islavsurf(2,*),ifaceq(8,6),ifacet(6,4),ifacew1(4,5), 40 & ifacew2(8,5),ipkon(*),kon(*) 41! 42 real*8 cg(3,*),straight(16,*),co(3,*),vold(0:mi(2),*),p(3), 43 & dist,xo(*),yo(*),zo(*),x(*),y(*),z(*),tietol(3,*),adjust, 44 & clearini(3,9,*),clearslavnode(3,*),clear 45! 46! nodes per face for hex elements 47! 48 data ifaceq /4,3,2,1,11,10,9,12, 49 & 5,6,7,8,13,14,15,16, 50 & 1,2,6,5,9,18,13,17, 51 & 2,3,7,6,10,19,14,18, 52 & 3,4,8,7,11,20,15,19, 53 & 4,1,5,8,12,17,16,20/ 54! 55! nodes per face for tet elements 56! 57 data ifacet /1,3,2,7,6,5, 58 & 1,2,4,5,9,8, 59 & 2,3,4,6,10,9, 60 & 1,4,3,8,10,7/ 61! 62! nodes per face for linear wedge elements 63! 64 data ifacew1 /1,3,2,0, 65 & 4,5,6,0, 66 & 1,2,5,4, 67 & 2,3,6,5, 68 & 3,1,4,6/ 69! 70! nodes per face for quadratic wedge elements 71! 72 data ifacew2 /1,3,2,9,8,7,0,0, 73 & 4,5,6,10,11,12,0,0, 74 & 1,2,5,4,7,14,10,13, 75 & 2,3,6,5,8,15,11,14, 76 & 3,1,4,6,9,13,12,15/ 77! 78 do i=1,ntie 79 if(tieset(1,i)(81:81).ne.'C') cycle 80 kneigh=1 81 slavset=tieset(2,i) 82! 83! check whether a clearance has been specified by the user 84! 85c tietol(3,i)=1.2357111317d0 86 if((tietol(3,i).gt.1.2357111316d0).and. 87 & (tietol(3,i).lt.1.2357111318d0)) then 88 iclear=0 89 else 90 iclear=1 91 clear=tietol(3,i) 92 endif 93! 94! check whether an adjust node set has been defined 95! only checked in the first increment of the first step 96! 97 if(istep.eq.1) then 98 iset=0 99 if(tieset(1,i)(1:1).ne.' ') then 100 noset(1:80)=tieset(1,i)(1:80) 101 noset(81:81)=' ' 102 ipos=index(noset,' ') 103 noset(ipos:ipos)='N' 104c do iset=1,nset 105c if(set(iset).eq.noset) exit 106c enddo 107 call cident81(set,noset,nset,id) 108 iset=nset+1 109 if(id.gt.0) then 110 if(noset.eq.set(id)) then 111 iset=id 112 endif 113 endif 114 kflag=1 115 call isortii(ialset(istartset(iset)),idummy, 116 & iendset(iset)-istartset(iset)+1,kflag) 117 endif 118 endif 119! 120! search a master face for each slave node; determine the 121! distance and adjust it if requested by the user 122! 123 nstart=itietri(1,i)-1 124 n=itietri(2,i)-nstart 125 if(n.lt.kneigh) kneigh=n 126 do j=1,n 127 xo(j)=cg(1,nstart+j) 128 x(j)=xo(j) 129 nx(j)=j 130 yo(j)=cg(2,nstart+j) 131 y(j)=yo(j) 132 ny(j)=j 133 zo(j)=cg(3,nstart+j) 134 z(j)=zo(j) 135 nz(j)=j 136 enddo 137 kflag=2 138 call dsort(x,nx,n,kflag) 139 call dsort(y,ny,n,kflag) 140 call dsort(z,nz,n,kflag) 141! 142 do j=nslavnode(i)+1,nslavnode(i+1) 143 node=islavnode(j) 144! 145 do k=1,3 146 p(k)=co(k,node) 147c p(k)=co(k,node)+vold(k,node) 148 enddo 149! 150! determining the kneigh neighboring master contact 151! triangle centers of gravity 152! 153 call near3d(xo,yo,zo,x,y,z,nx,ny,nz,p(1),p(2),p(3), 154 & n,neigh,kneigh) 155! 156 isol=0 157! 158 itriold=0 159 itri=neigh(1)+itietri(1,i)-1 160 ntriangle=0 161 ntriangle_=100 162! 163 loop1: do 164 do l=1,3 165 ll=4*l-3 166 dist=straight(ll,itri)*p(1)+ 167 & straight(ll+1,itri)*p(2)+ 168 & straight(ll+2,itri)*p(3)+ 169 & straight(ll+3,itri) 170 if(dist.gt.1.d-6) then 171 itrinew=imastop(l,itri) 172 if(itrinew.eq.0) then 173c write(*,*) '**border reached' 174 exit loop1 175 elseif(itrinew.eq.itriold) then 176c write(*,*) '**solution in between triangles' 177 isol=itri 178 exit loop1 179 else 180 call nident(itriangle,itrinew,ntriangle,id) 181 if(id.gt.0) then 182 if(itriangle(id).eq.itrinew) then 183c write(*,*) '**circular path; no solution' 184 exit loop1 185 endif 186 endif 187 ntriangle=ntriangle+1 188 if(ntriangle.gt.ntriangle_) then 189c write(*,*) '**too many iterations' 190 exit loop1 191 endif 192 do k=ntriangle,id+2,-1 193 itriangle(k)=itriangle(k-1) 194 enddo 195 itriangle(id+1)=itrinew 196 itriold=itri 197 itri=itrinew 198 cycle loop1 199 endif 200 elseif(l.eq.3) then 201c write(*,*) '**regular solution' 202 isol=itri 203 exit loop1 204 endif 205 enddo 206 enddo loop1 207! 208! calculating the distance 209! 210 if(isol.ne.0) then 211 dist=straight(13,itri)*p(1)+ 212 & straight(14,itri)*p(2)+ 213 & straight(15,itri)*p(3)+ 214 & straight(16,itri) 215! 216! check for an adjust parameter (only in the first 217! increment of the first step) 218! 219 if(istep.eq.1) then 220 if(iset.ne.0) then 221! 222! check whether node belongs to the adjust node 223! set 224! 225 call nident(ialset(istartset(iset)),node, 226 & iendset(iset)-istartset(iset)+1,id) 227 if(id.gt.0) then 228 if(ialset(istartset(iset)+id-1).eq.node) then 229 do k=1,3 230 co(k,node)=co(k,node)- 231 & dist*straight(12+k,itri) 232 enddo 233 dist=0.d0 234 endif 235 endif 236 elseif(dabs(tietol(1,i)).ge.2.d0) then 237! 238! adjust parameter 239! 240 adjust=dabs(tietol(1,i))-2.d0 241 242 if(dist.le.adjust) then 243 do k=1,3 244 co(k,node)=co(k,node)- 245 & dist*straight(12+k,itri) 246 enddo 247 dist=0.d0 248 endif 249 endif 250 endif 251! 252! clearance in each slave node 253! 254 if(iclear.eq.1) then 255 do k=1,3 256 clearslavnode(k,j)=(clear-dist)*straight(12+k,itri) 257 enddo 258 endif 259! 260 endif 261 enddo 262! 263! loop over all slave faces 264! 265 if(iclear.eq.1) then 266 do jj=itiefac(1,i), itiefac(2,i) 267 ifaces=islavsurf(1,jj) 268 nelems=int(ifaces/10) 269 jfaces=ifaces - nelems*10 270! 271 if(lakon(nelems)(4:5).eq.'8R') then 272 nopes=4 273 nope=8 274 elseif(lakon(nelems)(4:4).eq.'8') then 275 nopes=4 276 nope=8 277 elseif(lakon(nelems)(4:6).eq.'20R') then 278 nopes=8 279 nope=20 280 elseif(lakon(nelems)(4:5).eq.'20') then 281 nopes=8 282 nope=20 283 elseif(lakon(nelems)(4:5).eq.'10') then 284 nopes=6 285 nope=10 286 elseif(lakon(nelems)(4:4).eq.'4') then 287 nopes=3 288 nope=4 289! 290! treatment of wedge faces 291! 292 elseif(lakon(nelems)(4:4).eq.'6') then 293 nope=6 294 if(jfaces.le.2) then 295 nopes=3 296 else 297 nopes=4 298 endif 299 elseif(lakon(nelems)(4:5).eq.'15') then 300 nope=15 301 if(jfaces.le.2) then 302 nopes=6 303 else 304 nopes=8 305 endif 306 endif 307! 308! actual position of the nodes belonging to the 309! slave surface 310! 311 istart=nslavnode(i)+1 312 ilength=nslavnode(i+1)-nslavnode(i) 313! 314 do m=1,nopes 315 if((nope.eq.20).or.(nope.eq.8))then 316 node=kon(ipkon(nelems)+ifaceq(m,jfaces)) 317 elseif((nope.eq.10).or.(nope.eq.4).or.(nope.eq.14)) 318 & then 319 node=kon(ipkon(nelems)+ifacet(m,jfaces)) 320 elseif(nope.eq.15) then 321 node=kon(ipkon(nelems)+ifacew2(m,jfaces)) 322 else 323 node=kon(ipkon(nelems)+ifacew1(m,jfaces)) 324 endif 325 call nident(islavnode(istart),node,ilength,id) 326 do k=1,3 327 clearini(k,m,jj)= 328 & clearslavnode(k,nslavnode(i)+id) 329 enddo 330 enddo 331 enddo 332 endif 333! 334! transfer the clearance from the slave nodes in islavnode 335! to the nodes per slave face in islavsurf 336! 337 enddo 338! 339 return 340 end 341 342