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 interpolextnodes(ibounnod,nbounnod,co,doubleglob, 20 & integerglob,stress,ifront,nfront,ifrontrel,costruc, 21 & temp,nstep) 22! 23! interpolation of the master file stress results to the 24! bounal nodes of the crack shape(s) 25! 26 implicit none 27! 28 integer ibounnod(*),nbounnod,integerglob(*),nktet,netet,ne, 29 & nkon,nfaces,nfield,nselect,imastset,iselect(6),nterms, 30 & nelem,ialset(1),ifront(*),nfront,ifrontrel(*),kflag, 31 & iendset(1),istartset(1),konl(20),loopa,nstep, 32 & node,i,j,k,l,m,iconstant 33! 34 real*8 co(3,*),doubleglob(*),coords(3),ratio(20),costruc(3,*), 35 & stress(6,nstep,*),dist,temp(nstep,*),value(7) 36! 37 nktet=integerglob(1) 38 netet=integerglob(2) 39 ne=integerglob(3) 40 nkon=integerglob(4) 41 nfaces=integerglob(5) 42 nfield=13 43! 44! select the stresses for interpolation 45! 46 nselect=7 47 iselect(1)=1 48 do k=2,7 49 iselect(k)=k+3 50 enddo 51! 52 imastset=0 53 loopa=8 54! 55 nfront=0 56 loop:do i=1,nbounnod 57 node =ibounnod(i) 58c write(*,*) 'interpolextnodes ',node 59! 60 do j=1,3 61 coords(j)=co(j,node) 62 enddo 63! 64 call basis(doubleglob(1),doubleglob(netet+1), 65 & doubleglob(2*netet+1), 66 & doubleglob(3*netet+1),doubleglob(4*netet+1), 67 & doubleglob(5*netet+1),integerglob(6),integerglob(netet+6), 68 & integerglob(2*netet+6),doubleglob(6*netet+1), 69 & integerglob(3*netet+6),nktet,netet, 70 & doubleglob(4*nfaces+6*netet+1),nfield, 71 & doubleglob(nstep*13*nktet+4*nfaces+6*netet+1), 72 & integerglob(7*netet+6),integerglob(ne+7*netet+6), 73 & integerglob(2*ne+7*netet+6), 74 & integerglob(nkon+2*ne+7*netet+6), 75 & coords(1),coords(2),coords(3),value,ratio,iselect, 76 & nselect, 77 & istartset,iendset,ialset,imastset, 78 & integerglob(nkon+2*ne+8*netet+6),nterms,konl,nelem,loopa, 79 & dist) 80! 81! store the coordinates of the node; these may not be the same 82! als in co due to the projection of external nodes onto the 83! structure 84! 85 do j=1,3 86 costruc(j,i)=coords(j) 87 enddo 88! 89 temp(1,i)=value(1) 90 do j=1,6 91 stress(j,1,i)=value(j+1) 92 enddo 93! 94! interpolating the values of the subsequent steps 95! 96 do l=2,nstep 97 iconstant=13*(l-1)*nktet+4*nfaces+6*netet 98 do k=1,nselect 99 m=iselect(k) 100 value(k)=0.d0 101 do j=1,nterms 102 value(k)=value(k)+ratio(j)* 103 & doubleglob(iconstant+(konl(j)-1)*13+m) 104 enddo 105 enddo 106 temp(l,i)=value(1) 107 do j=1,6 108 stress(j,l,i)=value(j+1) 109 enddo 110 enddo 111! 112! check whether inside structure 113! 114 if(dist.lt.1.d-6) then 115 nfront=nfront+1 116 ifront(nfront)=node 117 ifrontrel(nfront)=i 118 endif 119! 120 enddo loop 121! 122! sorting the front nodes in ascending order 123! 124 kflag=2 125 call isortii(ifront,ifrontrel,nfront,kflag) 126! 127 return 128 end 129 130