forked from truongdangqe/Unstruct2D
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathdependentVars.f90
142 lines (120 loc) · 4.23 KB
/
dependentVars.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
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
!> @file dependentVars.f90
!!
!! Computation of dependent variables under the assumption of ideal gas
!! with constant properties.
!
! *****************************************************************************
!
! (c) J. Blazek, CFD Consulting & Analysis, www.cfd-ca.de
! Created February 25, 2014
! Last modification: May 24, 2014
!
! *****************************************************************************
!
! This program is free software; you can redistribute it and/or
! modify it under the terms of the GNU General Public License
! as published by the Free Software Foundation; either version 2
! of the License, or (at your option) any later version.
!
! This program is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
! GNU General Public License for more details.
!
! You should have received a copy of the GNU General Public License
! along with this program; if not, write to the Free Software
! Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
!
! *****************************************************************************
!> Computes values of dependent variables (pressure, temperature, speed
!! of sound, specific heat ratio, specific heat coeff. at const. pressure)
!! from conservative variables at all grid points. Additionally, laminar
!! viscosity and heat conductivity coefficients are computed in the case
!! of viscous flow.
!!
subroutine DependentVarsAll
use ModDataTypes
use ModGeometry
use ModPhysics
implicit none
! local variables
integer :: i
real(rtype) :: gam1, rgas, g1cp, rhoq, s1, s2, s12, rat, cppr
! *****************************************************************************
gam1 = gamma - 1.D0
rgas = gam1*cpgas/gamma
g1cp = gam1*cpgas
! Euler equations
if (kequs == "E") then
do i=1,nnodes
rhoq = cv(2,i)*cv(2,i) + cv(3,i)*cv(3,i)
dv(1,i) = gam1*(cv(4,i)-0.5D0*rhoq/cv(1,i))
dv(2,i) = dv(1,i)/(rgas*cv(1,i))
dv(3,i) = Sqrt(g1cp*dv(2,i))
dv(4,i) = gamma
dv(5,i) = cpgas
enddo
! Navier-Stokes equations
else
s1 = 110.D0
s2 = 288.16D0
s12 = 1.D0 + s1/s2
cppr = cpgas/prlam
do i=1,nnodes
rhoq = cv(2,i)*cv(2,i) + cv(3,i)*cv(3,i)
dv(1,i) = gam1*(cv(4,i)-0.5D0*rhoq/cv(1,i))
dv(2,i) = dv(1,i)/(rgas*cv(1,i))
dv(3,i) = Sqrt(g1cp*dv(2,i))
dv(4,i) = gamma
dv(5,i) = cpgas
rat = Sqrt(dv(2,i)/s2)*s12/(1.D0+s1/dv(2,i))
dv(6,i) = refvisc*rat
dv(7,i) = dv(6,i)*cppr
enddo
endif
end subroutine DependentVarsAll
! =============================================================================
!> Computes values of dependent variables (pressure, temperature, speed
!! of sound, specific heat ratio, specific heat coeff. at const. pressure)
!! from conservative variables at the node i. Additionally, laminar
!! viscosity and heat conductivity coefficients are computed in the case
!! of viscous flow.
!!
!! @param i node index
!!
subroutine DependentVarsOne( i )
use ModDataTypes
use ModPhysics
implicit none
! parameters
integer, intent(in) :: i
! local variables
real(rtype) :: gam1, rgas, g1cp, rhoq, s1, s2, s12, rat
! *****************************************************************************
gam1 = gamma - 1.D0
rgas = gam1*cpgas/gamma
g1cp = gam1*cpgas
! Euler equations
if (kequs == "E") then
rhoq = cv(2,i)*cv(2,i) + cv(3,i)*cv(3,i)
dv(1,i) = gam1*(cv(4,i)-0.5D0*rhoq/cv(1,i))
dv(2,i) = dv(1,i)/(rgas*cv(1,i))
dv(3,i) = Sqrt(g1cp*dv(2,i))
dv(4,i) = gamma
dv(5,i) = cpgas
! Navier-Stokes equations
else
s1 = 110.D0
s2 = 288.16D0
s12 = 1.D0 + s1/s2
rhoq = cv(2,i)*cv(2,i) + cv(3,i)*cv(3,i)
dv(1,i) = gam1*(cv(4,i)-0.5D0*rhoq/cv(1,i))
dv(2,i) = dv(1,i)/(rgas*cv(1,i))
dv(3,i) = Sqrt(g1cp*dv(2,i))
dv(4,i) = gamma
dv(5,i) = cpgas
rat = Sqrt(dv(2,i)/s2)*s12/(1.D0+s1/dv(2,i))
dv(6,i) = refvisc*rat
dv(7,i) = dv(6,i)*(cpgas/prlam)
endif
end subroutine DependentVarsOne