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 attach_2d(pneigh,pnode,nterms,ratio,dist,xil,etl)
20!
21!     ataches node with coordinates in "pnode" to the face containing
22!     "nterms" nodes with coordinates in field "pneigh" (nterms < 9).
23!     cave: the coordinates are stored in pneigh(1..3,*)
24!
25      implicit none
26!
27      integer nterms,i,j,k,imin,jmin,im,jm
28!
29      real*8 ratio(9),pneigh(3,9),pnode(3),a,xi(-1:1,-1:1),
30     &  et(-1:1,-1:1),p(3),distmin,d1,dist,xil,etl
31!
32      d1=1.d0
33!
34      xi(0,0)=0.d0
35      et(0,0)=0.d0
36      call distattach_2d(xi(0,0),et(0,0),pneigh,pnode,a,p,
37     &     ratio,nterms)
38      distmin=a
39      imin=0
40      jmin=0
41!
42      do k=1,8
43!
44!     initialisation
45!
46         d1=d1/10.d0
47!
48         do i=-1,1
49            do j=-1,1
50               if((i.eq.0).and.(j.eq.0)) cycle
51!
52               xi(i,j)=xi(0,0)+i*d1
53               et(i,j)=et(0,0)+j*d1
54!
55!              check whether inside the (-1,1)x(-1,1) domain
56!
57               if((xi(i,j).le.1.d0).and.
58     &              (xi(i,j).ge.-1.d0).and.
59     &              (et(i,j).le.1.d0).and.
60     &              (et(i,j).ge.-1.d0)) then
61                  call distattach_2d(xi(i,j),et(i,j),pneigh,pnode,a,p,
62     &                 ratio,nterms)
63!
64!                 checking for smallest initial distance
65!
66                  if(a.lt.distmin) then
67                     distmin=a
68                     imin=i
69                     jmin=j
70                  endif
71               endif
72!
73            enddo
74         enddo
75!
76!     minimizing the distance from the face to the node
77!
78         do
79!
80!     exit if minimum found
81!
82            if((imin.eq.0).and.(jmin.eq.0)) exit
83!
84!           new center of 3x3 matrix
85!
86            xi(0,0)=xi(imin,jmin)
87            et(0,0)=et(imin,jmin)
88!
89            im=imin
90            jm=jmin
91!
92            imin=0
93            jmin=0
94!
95            do i=-1,1
96               do j=-1,1
97                  if((i+im.lt.-1).or.(i+im.gt.1).or.
98     &                 (j+jm.lt.-1).or.(j+jm.gt.1)) then
99!
100                     xi(i,j)=xi(0,0)+i*d1
101                     et(i,j)=et(0,0)+j*d1
102!
103!              check whether inside the (-1,1)x(-1,1) domain
104!
105                     if((xi(i,j).le.1.d0).and.
106     &                  (xi(i,j).ge.-1.d0).and.
107     &                  (et(i,j).le.1.d0).and.
108     &                  (et(i,j).ge.-1.d0)) then
109                        call distattach_2d(xi(i,j),et(i,j),pneigh,
110     &                       pnode,a,p,ratio,nterms)
111!
112!                       check for new minimum
113!
114                        if(a.lt.distmin) then
115                           distmin=a
116                           imin=i
117                           jmin=j
118                        endif
119                     endif
120!
121                  endif
122               enddo
123            enddo
124         enddo
125      enddo
126!
127      call distattach_2d(xi(0,0),et(0,0),pneigh,pnode,a,p,
128     &     ratio,nterms)
129!
130      do i=1,3
131        pnode(i)=p(i)
132      enddo
133!
134      dist=dsqrt(a)
135!
136      if(nterms.eq.3) then
137         if(xi(0,0)+et(0,0).le.0.d0) then
138            xil=(xi(0,0)+1.d0)/2.d0
139            etl=(et(0,0)+1.d0)/2.d0
140         else
141            xil=(1.d0-et(0,0))/2.d0
142            etl=(1.d0-xi(0,0))/2.d0
143         endif
144      elseif(nterms.eq.4) then
145         xil=xi(0,0)
146         etl=et(0,0)
147      elseif(nterms.eq.6) then
148         if(xi(0,0)+et(0,0).le.0.d0) then
149            xil=(xi(0,0)+1.d0)/2.d0
150            etl=(et(0,0)+1.d0)/2.d0
151         else
152            xil=(1.d0-et(0,0))/2.d0
153            etl=(1.d0-xi(0,0))/2.d0
154         endif
155      elseif(nterms.eq.8) then
156         xil=xi(0,0)
157         etl=et(0,0)
158      endif
159!
160c      write(*,*) 'end attach_2d.f'
161      return
162      end
163
164