xref: /phasta/phSolver/common/setSuction_Duct3.f90 (revision 1016729149754f57cd03fe576ba6fd0f1723ab31)
1c========================================================================
2c Dynamically controll the suction level based on main contraction mdot,
3c called by itrdrv.f
4c========================================================================
5        subroutine setSuction_Duct3(x,BC,y, ilwork)
6
7        use wallData ! wnorm
8        use timedata
9        include "common.h"
10        include "mpif.h"
11        include "auxmpi.h"
12
13        integer i, j, k, nn, nL
14        integer ilwork(nlwork)
15        integer nSuctionNode
16        integer, allocatable :: suctionNodeMap(:)
17!        logical inPatch
18        real*8 xcoor,ycoor,zcoor
19        real*8 BC(nshg,ndofBC)
20        real*8 BC3(numnp, 5)                      !used for hack to get suction right on part boundaries
21        real*8 x(nshg,nsd), y(nshg,ndof)
22!       real*8 xmin_t, xmax_t, xmin_b, xmax_b     !ramp geometry controls
23!       real*8 xmin_su, xmax_su, xmin_sl, xmax_sl !side wall patch limits
24!       real*8 t_upper                            !upper thickness of suction area
25!       real*8 ymin_su, ymax_sl
26        real*8 nvel_s, nvel_t, nvel_b             !normal velocity on the side walls, top, and bottom
27
28!        real*8 delX, delY, X1, Y1     !Variables used to map upper surface or trumpet
29!        real*8 Xline(30,2)
30!        integer suctionIDs(4)
31
32!        real*8 wallBCout(nshg, 6)
33
34        if(myrank.eq.0)write(*,*)'Setting side suctions'
35
36        nvel_b     = suctionVbottom      !normal velocity at the bottom patch
37        nvel_sl    = suctionVside_lower  !normal velocity of lower side wall patches
38        nvel_su    = suctionVside_upper  !normal velocity of upper side wall patches
39        nvel_t     = suctionVtop         !normal velocity at the top patch
40
41        allocate(suctionNodeMap(nshg))
42        call sfID2np(isetSuctionID_Duct, nSuctionNode, suctionNodeMap)  !to figure out whether it's the top, bottom, or side surface.
43
44        smallNum=1.0e-5
45        xsideBotMax=0.0428625+smallNum
46        xsideTopMax=0.130175+smallNum
47
48        do i = 1, nSuctionNode
49          nn = suctionNodeMap(i) ! map face node index to global node
50          xcoor = x(nn,1)
51          ycoor = x(nn,2)
52          !zcoor = x(nn,3)
53
54          !wnormx = wnorm(nn, 1)
55          wnormy = wnorm(nn, 2)
56          wnormz = wnorm(nn, 3)
57
58          !Test if the point lies in one of the suction patches
59
60          if((abs(nvel_sl) >  0  .or. abs(nvel_su) > 0) .and.  !if either upper or lower side walls are on...
61     &       (abs(wnormz) >=  0.999)) then                     ! and the node normal is in the z direction
62            if((ycoor > 0) .and. (xcoor .le. xsideTopMax)) then !hack to prevent edges from turning on.
63              BC(nn, 3) = nvel_su*wnorm(nn, 1) ! set xVel
64              BC(nn, 4) = nvel_su*wnorm(nn, 2) ! set yVel
65              BC(nn, 5) = nvel_su*wnorm(nn, 3) ! set zVel
66            elseif (xcoor .le. xsideBotMax) then
67              BC(nn, 3) = nvel_sl*wnorm(nn, 1) ! set xVel
68              BC(nn, 4) = nvel_sl*wnorm(nn, 2) ! set yVel
69              BC(nn, 5) = nvel_sl*wnorm(nn, 3) ! set zVel
70            endif
71          else if(abs(wnormz) > 0.1) then      ! side walls should have a larger znormal than this - it's ambiguous what to do so do nothing
72            cycle
73          else if(abs(nvel_t) > 0      .and. !top suction patch
74     &                 wnormy > 0.001) then  !face is pointing up (in the +y direction)
75            BC(nn, 3) = nvel_t*wnorm(nn, 1) ! set xVel using wall normal
76            BC(nn, 4) = nvel_t*wnorm(nn, 2) ! set yVel
77            BC(nn, 5) = nvel_t*wnorm(nn, 3) ! set zVel
78
79          else if(abs(nvel_b) > 0       .and. !bottom suction patches
80     &                wnormy  < -0.001) then  !face is pointing down (in the -y direction)
81            BC(nn, 3) = nvel_b*wnorm(nn, 1) ! set xVel using wall normal
82            BC(nn, 4) = nvel_b*wnorm(nn, 2) ! set yVel
83            BC(nn, 5) = nvel_b*wnorm(nn, 3) ! set zVel
84          endif
85
86        enddo
87
88        if(numpe.gt.1) then
89            BC3(:,1:3)=BC(:,3:5)
90            do i=1,nshg
91                bc3mag=BC3(i,1)**2+BC3(i,2)**2+BC3(i,3)**2
92                if(bc3mag.gt.0) BC3(i,4)=one
93            enddo
94            call commu(BC3,ilwork,5,'in ')
95! This accumulated and failed                 call commu(BC(:,3:5),ilwork,3,'in ')
96            do i=1,nshg
97                bcmag=BC(i,3)**2+BC(i,4)**2+BC(i,5)**2
98                if((bcmag.eq.0).and.(BC3(i,4).ne.0)) then !node is slave and has not been correctly set.
99                     BC(i,3:5)=BC3(i,1:3)/BC3(i,4)  !division by BC3 is necessary to account for more than two blocks sharing the same node
100                endif
101            enddo
102       endif
103
104        !Debugging loop
105!        if(nSuctionNode.gt.0) then
106!          icount=0
107!          do i=1,nSuctionNode
108!             nn=suctionNodeMap(i)
109!             dz=dabs(dabs(x(nn,3))-0.0762)
110!             if(dz.lt.0.0001) then
111!               icount=icount+1
112!               if(icount.eq.1)
113!     & write(myrank+1000,*) 'writing wnorm and BC setSuctionNG3'
114!               write(myrank+1000,1234) nn,i,
115!     &           x(nn,1), x(nn,2), x(nn,3),
116!     &           wnorm(nn,1),wnorm(nn,2),wnorm(nn,3),
117!     &           BC(nn,3),BC(nn,4),BC(nn,5)
118!             endif
119!          enddo
120!         if(icount.gt.0) close(myrank+1000)
121!        endif
122!1234           format(i5,2x,i5,9(2x,e14.7))
123
124      end
125
126
127
128
129
130
131