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 findslavcfd(nmpc,labmpc,ipompc,nodempc,islav,
20     &  nslav,inoslav,inomast,ics,cs,imast,nmast,co,inomat,
21     &  nr,nz,rcs,zcs,rcs0,zcs0,ncs)
22!
23!     find the slave nodes in a cyclic symmetry constraint for
24!     CFD calculations
25!
26      implicit none
27!
28      character*20 labmpc(*)
29!
30      integer i,j,nmpc,nodeslav,nodempc(3,*),ipompc(*),nslav,id,
31     &  islav(*),inoslav(*),inomast(*),nodemast,ics(*),imast(*),
32     &  nmast,inomat(*),nr(*),nz(*),noden(1),nneigh,ncs,kflag,node
33!
34      real*8 cs(17,*),co(3,*),xa(3),xn(3),rcs(*),zcs(*),rcs0(*),
35     &  zcs0(*),rp,zp,xap,yap,zap,dd
36!
37      nslav=0
38      nmast=0
39!
40!     cyclic symmetry slave nodes are the dependent nodes in a
41!     CYCLIC MPC; The fluid cyclic MPC's are assumed to be one-to-one
42!     node-MPC's (slave and master side must match)
43!
44!     they are stored and ordered in field islav
45!
46!     the fluid master nodes are stored in imast
47!
48!     inoslav contains for each fluid node the corresponding node on the
49!             slave side, if any
50!     inomast contains for each fluid node the corresponding node on the
51!             master side, if any
52!
53!     xn is an normed vector along the axis
54!     xa is a point on the axis
55!
56      do i=1,3
57         xa(i)=cs(5+i,1)
58         xn(i)=cs(8+i,1)-xa(i)
59      enddo
60      dd=dsqrt(xn(1)*xn(1)+xn(2)*xn(2)+xn(3)*xn(3))
61      do i=1,3
62         xn(i)=xn(i)/dd
63      enddo
64!
65!     catalogue all master nodes (fluid and solid) by their
66!     radial and axial coordinate
67!
68      do i=1,ncs
69         node=ics(i)
70!
71         xap=co(1,node)-xa(1)
72         yap=co(2,node)-xa(2)
73         zap=co(3,node)-xa(3)
74!
75         zcs(i)=xap*xn(1)+yap*xn(2)+zap*xn(3)
76         rcs(i)=dsqrt((xap-zcs(i)*xn(1))**2+
77     &                (yap-zcs(i)*xn(2))**2+
78     &                (zap-zcs(i)*xn(3))**2)
79      enddo
80!
81!     initialization for near2d
82!
83      do i=1,ncs
84         nr(i)=i
85         nz(i)=i
86         rcs0(i)=rcs(i)
87         zcs0(i)=zcs(i)
88      enddo
89      kflag=2
90      call dsort(rcs,nr,ncs,kflag)
91      call dsort(zcs,nz,ncs,kflag)
92!
93!     only one neighbor is looked for (one-to-one correspondence
94!     of slave and master side assumed)
95!
96      nneigh=1
97!
98      do i=1,nmpc
99         if(labmpc(i)(1:6).ne.'CYCLIC') cycle
100         nodeslav=nodempc(1,ipompc(i))
101!
102!        check whether slave node is cfd-node
103!
104         if(inomat(nodeslav).eq.0) cycle
105!
106         call nident(islav,nodeslav,nslav,id)
107         if(id.gt.0) then
108            if(islav(id).eq.nodeslav) cycle
109         endif
110!
111         xap=co(1,nodeslav)-xa(1)
112         yap=co(2,nodeslav)-xa(2)
113         zap=co(3,nodeslav)-xa(3)
114!
115         zp=xap*xn(1)+yap*xn(2)+zap*xn(3)
116         rp=dsqrt((xap-zp*xn(1))**2+(yap-zp*xn(2))**2+(zap-zp*xn(3))**2)
117!
118         call near2d(rcs0,zcs0,rcs,zcs,nr,nz,rp,zp,ncs,noden,nneigh)
119!
120         nodemast=ics(noden(1))
121!
122         if(dabs((rcs0(noden(1))-rp)**2+
123     &           (zcs0(noden(1))-zp)**2).gt.1.d-10) then
124            write(*,*) '*ERROR in findslavcfd: cyclic symmetric fluid'
125            write(*,*) '       slave and master nodes do not match'
126            stop
127         endif
128!
129         nslav=nslav+1
130         do j=nslav,id+2,-1
131            islav(j)=islav(j-1)
132         enddo
133         islav(id+1)=nodeslav
134!
135         call nident(imast,nodemast,nmast,id)
136         nmast=nmast+1
137         do j=nmast,id+2,-1
138            imast(j)=imast(j-1)
139         enddo
140         imast(id+1)=nodemast
141!
142         inomast(nodeslav)=nodemast
143         inoslav(nodemast)=nodeslav
144      enddo
145!
146c      do i=1,nslav
147c         write(*,*) i,islav(i)
148c      enddo
149c      do i=1,nmast
150c         write(*,*) i,imast(i)
151c      enddo
152c      do i=1,38
153c         write(*,*) i,inoslav(i),inomast(i)
154c      enddo
155      return
156      end
157
158