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