-
Notifications
You must be signed in to change notification settings - Fork 4
/
initialisation.f90
114 lines (111 loc) · 3.45 KB
/
initialisation.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
! This Source Code Form is subject to the terms of the Mozilla Public License, v. 2.0. If a copy of the MPL was not distributed with this file, You can obtain one at https://mozilla.org/MPL/2.0/.
subroutine initblfv(muinf,xc,yc,w,field,im,gh,jm)
!
implicit none
! variables for dimension -----------------------------------------
integer :: im,jm,gh
! required arguments ----------------------------------------------
real(8),intent(in) :: muinf
real(8),dimension(1-gh:im+gh ,1-gh:jm+gh ),intent(in) :: xc
real(8),dimension(1-gh:im+gh ,1-gh:jm+gh ),intent(in) :: yc
! Returned objects ------------------------------------------------
real(8),dimension(1-gh:im+gh ,1-gh:jm+gh , 5 ),intent(inout) :: w
real(8),dimension( jm , gh , 5 ),intent(inout) :: field
! Non-required arguments -------------------------------------------
real(8) :: eta,ubl,ue ! BL parameter
real(8) :: ro, roE, rou, xdeltam1,rou2
integer :: i,j
!
! w contains state inf adim
ro = w(3,3,1)
rou = w(3,3,2)
roE = w(3,3,5)
ue = rou/ro
rou2 = 0.5d0*rou*ue
xdeltam1 = sqrt(rou/muinf)
!init w
do j=1,jm+gh
do i=1-gh, im+gh
eta = yc(i,j)*xdeltam1/sqrt(xc(i,j))
if (eta.lt.1.d0) then
ubl = 2.d0*eta -2.d0*eta*eta*eta + eta*eta*eta*eta
else
ubl = 1.d0
endif
ubl = ubl * ue
w(i,j,2) = ro*ubl
w(i,j,5) = roE - rou2 + 0.5d0*ro*ubl*ubl
enddo
enddo
w(:,0,2) = - w(:,1,2)
!set field
do j=1,jm
do i=1-gh,0
eta = yc(i,j)*xdeltam1*sqrt(xc(i,j))
if (eta.lt.1.d0) then
ubl = 2.d0*eta -2.d0*eta*eta*eta + eta*eta*eta*eta
else
ubl = 1.d0
endif
ubl = ubl * ue
field(j,1-i,2) = ro*ubl
field(j,1-i,5) = roE - rou2 + 0.5d0*ro*ubl*ubl
! write(209876,*) "i,j, field(j,i)", i,j,
enddo
enddo
end subroutine initblfv
subroutine initbl(muinf,xc,yc,w,field,im,gh,jm)
!
implicit none
! variables for dimension -----------------------------------------
integer :: im,jm,gh
! required arguments ----------------------------------------------
real(8),intent(in) :: muinf
real(8),dimension(1-gh:im+gh ),intent(in) :: xc
real(8),dimension( 1-gh:jm+gh ),intent(in) :: yc
! Returned objects ------------------------------------------------
real(8),dimension(1-gh:im+gh ,1-gh:jm+gh , 5 ),intent(inout) :: w
real(8),dimension( jm , gh , 5 ),intent(inout) :: field
! Non-required arguments -------------------------------------------
real(8) :: eta,ubl,ue ! BL parameter
real(8) :: ro, roE, rou, xdeltam1,rou2
integer :: i,j
!
! w contains state inf adim
ro = w(3,3,1)
rou = w(3,3,2)
roE = w(3,3,5)
ue = rou/ro
rou2 = 0.5d0*rou*ue
xdeltam1 = sqrt(rou/muinf)
!init w
do j=1,jm+gh
do i=1-gh, im+gh
eta = yc(j)*xdeltam1/sqrt(xc(i))
if (eta.lt.1.d0) then
ubl = 2.d0*eta -2.d0*eta*eta*eta + eta*eta*eta*eta
else
ubl = 1.d0
endif
ubl = ubl * ue
w(i,j,2) = ro*ubl
w(i,j,5) = roE - rou2 + 0.5d0*ro*ubl*ubl
enddo
enddo
w(:,0,2) = - w(:,1,2)
!set field
do j=1,jm
do i=1-gh,0
eta = yc(j)*xdeltam1*sqrt(xc(i))
if (eta.lt.1.d0) then
ubl = 2.d0*eta -2.d0*eta*eta*eta + eta*eta*eta*eta
else
ubl = 1.d0
endif
ubl = ubl * ue
field(j,1-i,2) = ro*ubl
field(j,1-i,5) = roE - rou2 + 0.5d0*ro*ubl*ubl
! write(209876,*) "i,j, field(j,i)", i,j,
enddo
enddo
end subroutine initbl