xref: /phasta/phSolver/incompressible/BCprofile.f (revision 8746ab438bbda91291f8cdd62b94f8385f2d26f1)
1c-----------------------------------------------------------------------
2        subroutine BCprofileScale(vbc_prof,BC,yold)
3        use pvsQbi
4        include "common.h"
5        real*8 vbc_prof(nshg,3)
6        real*8 BC(nshg,ndofBC),yold(nshg,ndof)
7        real*8 offphase
8        integer factor
9
10        r_amp =rampmdot(1,1)
11        r_freq=rampmdot(2,1)
12!  Usual sinusoidal in time syn jet is given in the next line.  It assumes that you are NOT changing the time step from previous runs since it computes the
13!  current time in the sin function to be (lstep+1)* Dt where lstep is a running total of all of the runs up to now.  This will be "ok" under the following
14!  assumption:  lstep*dt_current = m * T where m is an integer and T is the period of the jet.  That is because this will evaluate to 2*m*Pi and the sin
15!  that is zero
16
17        if(abs(rampmdot(1,2)) .lt. 1.0e-5) then
18           r_time_factor = r_amp*sin(two*pi*r_freq*(lstep+1)*Delt(1)) ! BC set for next time step - all phases are computed
19        else
20           r_time_factor = rampmdot(1,2)*r_amp   ! read parameter scales Vmax
21        endif
22
23        icount = 0
24        do kk=1,nshg
25          if(ndsurf(kk).eq.18) then ! this means diaphragm for the Cube Test case
26
27            factor = idnint(rampmdot(2,2))
28
29            if(factor == 0) then
30              offphase = 0.d0
31
32            elseif(factor == 1) then
33              offphase = 1.d0
34
35            elseif(factor == 2) then
36              !Start to count from tip. If factor == 2 -> 1, 3, 5, etc
37              if(mod(ndsurf(kk)-1,factor) == 0) then
38                offphase = 1.d0
39              else
40                offphase = 0.d0
41              endif
42
43            elseif(factor == 3) then
44              !Start to count from tip. If factor == 3 -> 2, 5, 8, etc
45              if(mod(ndsurf(kk)+1,factor) == 0) then
46                offphase = 1.d0
47              else
48                offphase = 0.d0
49              endif
50
51            elseif(factor < 0) then
52              !Only one jet blowing. If factor = -5, then jet 5 only
53              !blows
54              if (ndsurf(kk) == -factor) then
55                offphase = 1.d0
56              else
57                offphase = 0.d0
58              endif
59            endif
60
61            BC(kk,3)=r_time_factor*vbc_prof(kk,1)*offphase
62            BC(kk,4)=r_time_factor*vbc_prof(kk,2)*offphase
63            BC(kk,5)=r_time_factor*vbc_prof(kk,3)*offphase
64
65            icount = icount + 1
66
67          endif
68        enddo
69
70!        if(istep.eq.0 .and. icount.ne.0)
71!     &     write(*,*) 'BCprofile count',myrank,icount
72
73        return
74        end
75
76c--------------------------------------------------------------
77        subroutine BCprofileInit(vbc_prof,x)
78
79        use pvsQbi
80        include "common.h"
81        real*8 vbc_prof(nshg,3), x(numnp,nsd)
82        real*8 rcenter(3),rnorml(3)
83
84        open(unit=789, file='bcprofile.dat',status='unknown')
85        do kk=1,nshg
86c.............Factors below are negative for desired blowing direction
87           if(ndsurf(kk).eq.18) then
88             x1=-0.025d0
89             x2=0.025d0
90             z1=-0.01d0
91             z2=0.01d0
92             vbc_prof(kk,1)=0.d0
93             vbc_prof(kk,2)=-1.6*1E7*(x(kk,1)-x1)*(x(kk,1)-x2)*
94     &                              (x(kk,3)-z1)*(x(kk,3)-z2)
95             vbc_prof(kk,3)=0.d0
96
97             write(789,987) kk,vbc_prof(kk,1),vbc_prof(kk,2),
98     &                          vbc_prof(kk,3)
99
100           else
101              vbc_prof(kk,:)=zero
102           endif
103
104        enddo
105        close(789)
106987     format(i6,3(2x,e14.7))
107
108        return
109        end
110