1!$Id:$ 2 subroutine ptiend(ie,ix,ib,ip,x,ma,nie,nen,nen1,ndm,numel) 3 4! * * F E A P * * A Finite Element Analysis Program 5 6!.... Copyright (c) 1984-2017: Regents of the University of California 7! All rights reserved 8 9!-----[--.----+----.----+----.-----------------------------------------] 10! Purpose: Connect line elements to boundaries of 2-d solids 11 12! Inputs: 13! ie(nie,*) - Assembly information for material sets 14! ib(*) - Node indicators?? 15! x(ndm,*) - Nodal coordinates of mesh 16! ma - Material set number of line element to tie 17! nie - Dimension of ie array 18! nen - Number of nodes on element 19! nen1 - Dimension of ix array 20! ndm - Spatial dimension of problem 21! numel - Number of elements in mesh 22 23! Outputs: 24! ix(nen1,*)- Element nodal connection list 25! ip(*) - Nodal number list 26!-----[--.----+----.----+----.-----------------------------------------] 27 28 implicit none 29 30 include 'iofile.h' 31 32 integer ma,nie,nen,nen1,ndm,numel 33 34 integer n,m, ii,i1,i2,ij1, j1,j2,ij2, n1,n2,nn, m1,m2 35 real*8 dx,dy, dx1,dx2, dy1,dy2, tol 36 37 integer ie(nie,*),ix(nen1,*),ib(*),ip(*),iord1(50),iord2(50) 38 real*8 x(ndm,*) 39 40 save 41 42 data tol/ 1.d-08 / 43 44! Loop over elements: search for line type 'ma' 45 46 do n = 1,numel 47 if(ix(nen1,n).eq.ma) then 48 49! Find edge on boundary of line: must have non-zero length 50 51 call pltord(ix(1,n),ie(nie-1,ix(nen1,n)),ij1,iord1) 52 do i1 = 1,ij1 53 i2 = mod(i1,ij1) + 1 54 n1 = ix(iord1(i1),n) 55 n2 = ix(iord1(i2),n) 56 if(ib(n1).ne.0 .and. ib(n2).ne.0) then 57 dx1 = x(1,n2) - x(1,n1) 58 dy1 = x(2,n2) - x(2,n1) 59 if(abs(dx1) + abs(dy1) .gt. tol) then 60 61! Search other types of elements for contiguous edge 62 63 do m = 1,numel 64 if(ix(nen1,m).ne.ma) then 65 call pltord(ix(1,m),ie(nie-1,ix(nen1,m)),ij2,iord2) 66 do j1 = 1,ij2 67 j2 = mod(j1,ij2) + 1 68 m1 = ix(iord2(j1),m) 69 m2 = ix(iord2(j2),m) 70 if(ib(m1).ne.0 .and. ib(m2).ne.0) then 71 72! Check that node pairs have same coordinates 73 74 dx = abs(x(1,m1)-x(1,n2))+abs(x(1,m2)-x(1,n1)) 75 dy = abs(x(2,m1)-x(2,n2))+abs(x(2,m2)-x(2,n1)) 76 if(dx + dy .lt. tol) then 77 78! Check side has non-zero length 79 80 dx2 = x(1,m2) - x(1,m1) 81 dy2 = x(2,m2) - x(2,m1) 82 if(abs(dx2) + abs(dy2) .gt. tol) then 83 84! Contiguous edge has opposite normals 85 86 if(abs(dx1+dx2)+abs(dy1+dy2).lt.tol) then 87 ix(iord2(i1),n) = m2 88 ix(iord2(i2),n) = m1 89 ib(m1) = -1 90 ib(m2) = -1 91 ib(n1) = -1 92 ib(n2) = -1 93 ip(n1) = m2 94 ip(n2) = m1 95 do nn = 1,numel 96 do ii = 1,nen 97 if(ix(ii,nn).eq.n1) ix(ii,nn) = m2 98 if(ix(ii,nn).eq.n2) ix(ii,nn) = m1 99 end do 100 end do 101 endif 102 endif 103 endif 104 endif 105 end do 106 endif 107 end do 108 endif 109 endif 110 end do 111 endif 112 end do 113 114 end 115