MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_hyperdiffusivity.t
Go to the documentation of this file.
1!> Subroutines for Roe-type Riemann solver for HD
3
4
5 implicit none
6
7 !all public
8 !private
9 !public :: hyperdiffusivity_init
10
11
12contains
13
16 use mod_geometry
17
18 !print*, slab_uniform !!THIS IS FALSE, why??
19
20 !if(coordinate .ne. Cartesian .or. .not. slab_uniform) then
21 if(coordinate .ne. cartesian) then
22 call mpistop("Hyperdiffusivity only implemented for Cartesian uniform grid")
23 endif
24
26
27 end subroutine hyperdiffusivity_init
28
29
30 subroutine hyp_coeff(ixI^L, ixO^L, var, idimm, nu_hyp)
32 integer, intent(in) :: ixI^L, idimm
33 integer, intent(out) :: ixO^L
34 double precision, intent(in) :: var(ixI^S)
35 double precision, intent(out) :: nu_hyp(ixI^S)
36
37 double precision :: tmp(ixI^S), tmp2(ixI^S)
38 integer :: hx^L, hx1f^L, hx1b^L, hx2b^L
39
40 double precision, parameter :: eps=1e-8
41
42
43 hxmin^d=iximin^d+2*kr(idimm,^d);
44 hxmax^d=iximax^d-kr(idimm,^d);
45
46 hx1f^l=hx^l+kr(idimm,^d);
47 hx1b^l=hx^l-kr(idimm,^d);
48 hx2b^l=hx^l-2*kr(idimm,^d);
49
50 !store d3 in tmp
51 tmp(hx^s)= abs(3d0 * (var(hx^s) - var(hx1b^s)) - (var(hx1f^s)-var(hx2b^s)))
52 !store d1 in tmp2
53 tmp2(hx^s)= abs((var(hx^s) - var(hx1b^s)))
54
55 ixomin^d=hxmin^d+kr(idimm,^d);
56 ixomax^d=hxmax^d-kr(idimm,^d);
57 hx1f^l=ixo^l+kr(idimm,^d);
58 hx1b^l=ixo^l-kr(idimm,^d);
59
60 nu_hyp(ixo^s) = max(tmp(hx1b^s),tmp(ixo^s),tmp(hx1f^s))/max(tmp2(hx1b^s),tmp2(ixo^s),tmp2(hx1f^s),eps)
61
62 !print*, "HYP IXO ", ixO^L
63
64 end subroutine hyp_coeff
65
66
67
68 subroutine div_vel_coeff(ixI^L, ixO^L, vel, idimm, nu_vel)
70 integer, intent(in) :: ixI^L, idimm
71 integer, intent(out) :: ixO^L
72 double precision, intent(in) :: vel(ixI^S,1:ndir)
73 double precision, intent(out) :: nu_vel(ixI^S)
74 integer :: hx1f^L, hx1b^L
75
76 integer :: ii
77
78 ixomin^d=iximin^d+1;
79 ixomax^d=iximax^d-1;
80
81 ixomin^d=ixomin^d+kr(idimm,^d);
82
83 nu_vel(ixo^s) = 0d0
84
85 ! ndim should always be smaller or equal to ndir !
86 do ii = 1, ndim
87 hx1b^l=ixo^l-kr(ii,^d);
88 if(ii==idimm) then
89 nu_vel(ixo^s) = nu_vel(ixo^s) + (vel(ixo^s,ii) - vel(hx1b^s,ii))/dxlevel(ii)
90 else
91 hx1f^l=ixo^l+kr(ii,^d);
92 nu_vel(ixo^s) = nu_vel(ixo^s) + (vel(hx1f^s,ii) - vel(hx1b^s,ii))/(2d0*dxlevel(ii))
93 endif
94
95 where(nu_vel(ixo^s) < 0d0)
96 nu_vel(ixo^s) = -nu_vel(ixo^s)
97 elsewhere
98 nu_vel(ixo^s) = 0d0
99 endwhere
100
101 enddo
102
103 !print*, "DIV_VEL IXO ", ixO^L
104
105 end subroutine div_vel_coeff
106
107
108
109 !var has cell centered values
110 !nu_hyper is defined at the interfaces
111 subroutine second_same_deriv(ixI^L, ixO^L, nu_hyper, var, idimm, res)
113 integer, intent(in) :: ixI^L, idimm
114 integer, intent(out) :: ixO^L
115 double precision, intent(in) :: nu_hyper(ixI^S),var(ixI^S)
116 double precision, intent(out) :: res(ixI^S)
117
118 integer :: hxf^L, hxb^L
119
120 ixomin^d=iximin^d+3;
121 ixomax^d=iximax^d-3;
122
123 hxf^l=ixo^l+kr(idimm,^d);
124 hxb^l=ixo^l-kr(idimm,^d);
125
126 res(ixo^s) = 1d0/(dxlevel(idimm)**2)*&
127 (nu_hyper(hxf^s) * (var(hxf^s)-var(ixo^s))-&
128 nu_hyper(ixo^s) * (var(ixo^s)-var(hxb^s)))
129 !print*, "SECOND SAME DERIV IXO ", ixO^L
130
131 end subroutine second_same_deriv
132
133 !var has cell centered values
134 !var2 has cell centered values
135 !nu_hyper is defined at the interfaces
136 subroutine second_same_deriv2(ixI^L, ixO^L, nu_hyper, var2, var, idimm, res)
138 integer, intent(in) :: ixI^L, idimm
139 integer, intent(out) :: ixO^L
140 double precision, intent(in) :: nu_hyper(ixI^S),var(ixI^S), var2(ixI^S)
141 double precision, intent(out) :: res(ixI^S)
142
143 integer :: hxf^L, hxb^L
144
145 ixomin^d=iximin^d+3;
146 ixomax^d=iximax^d-3;
147
148 hxf^l=ixo^l+kr(idimm,^d);
149 hxb^l=ixo^l-kr(idimm,^d);
150
151 res(ixo^s) = 1d0/(2d0*dxlevel(idimm)**2)*&
152 (nu_hyper(hxf^s) * (var2(hxf^s)+var2(ixo^s)) * (var(hxf^s)-var(ixo^s))-&
153 nu_hyper(ixo^s) * (var2(hxb^s)+var2(ixo^s)) * (var(ixo^s)-var(hxb^s)))
154 !print*, "SECOND SAME DERIV2 IXO ", ixO^L
155
156 end subroutine second_same_deriv2
157
158 !idimm inner derivative, idimm2 outer
159 ! deriv_idimm2(nu * deriv_idimm (u) )
160 !var has cell centered values
161 !nu_hyper is defined at the interfaces
162 subroutine second_cross_deriv(ixI^L, ixO^L, nu_hyper, var, idimm, idimm2, res)
164 integer, intent(in) :: ixI^L, idimm,idimm2
165 integer, intent(out) :: ixO^L
166 double precision, intent(in) :: nu_hyper(ixI^S),var(ixI^S)
167 double precision, intent(out) :: res(ixI^S)
168
169 integer :: hxfi^L, hxbi^L, hxfifj^L, hxbibj^L, hxfibj^L, hxbifj^L
170
171 ixomin^d=iximin^d+3;
172 ixomax^d=iximax^d-3;
173
174 hxfi^l=ixo^l+kr(idimm2,^d);
175 hxbi^l=ixo^l-kr(idimm2,^d);
176
177 hxfifj^l=hxfi^l+kr(idimm,^d);
178 hxfibj^l=hxfi^l-kr(idimm,^d);
179
180 hxbifj^l=hxbi^l+kr(idimm,^d);
181 hxbibj^l=hxbi^l-kr(idimm,^d);
182
183 res(ixo^s) = 1d0/(8d0*dxlevel(idimm) * dxlevel(idimm2))*&
184 ((nu_hyper(hxfifj^s) + nu_hyper(hxfi^s)) * (var(hxfifj^s)-var(hxfibj^s))-&
185 (nu_hyper(hxbifj^s) + nu_hyper(hxbi^s)) * (var(hxbifj^s)-var(hxbibj^s)))
186
187 !print*, "SECOND CROSS DERIV IXO ", ixO^L
188
189 end subroutine second_cross_deriv
190
191
192 !idimm inner derivative, idimm2 outer
193 ! deriv_idimm2(nu * var2 * deriv_idimm (u) )
194 !var has cell centered values
195 !var2 has cell centered values
196 !nu_hyper is defined at the interfaces (with numbering index left center): center_i-1 interface_i center_i
197 subroutine second_cross_deriv2(ixI^L, ixO^L, nu_hyper, var2, var, idimm, idimm2, res)
199 integer, intent(in) :: ixI^L, idimm,idimm2
200 integer, intent(out) :: ixO^L
201 double precision, intent(in) :: nu_hyper(ixI^S),var(ixI^S),var2(ixI^S)
202 double precision, intent(out) :: res(ixI^S)
203
204 integer :: hxfi^L, hxbi^L, hxfifj^L, hxbibj^L, hxfibj^L, hxbifj^L
205
206 ixomin^d=iximin^d+3;
207 ixomax^d=iximax^d-3;
208
209 hxfi^l=ixo^l+kr(idimm2,^d);
210 hxbi^l=ixo^l-kr(idimm2,^d);
211
212 hxfifj^l=hxfi^l+kr(idimm,^d);
213 hxfibj^l=hxfi^l-kr(idimm,^d);
214
215 hxbifj^l=hxbi^l+kr(idimm,^d);
216 hxbibj^l=hxbi^l-kr(idimm,^d);
217
218 res(ixo^s) = 1d0/(8d0*dxlevel(idimm) * dxlevel(idimm2))*&
219 ((nu_hyper(hxfifj^s) + nu_hyper(hxfi^s)) * var2(hxfi^s) * (var(hxfifj^s)-var(hxfibj^s))-&
220 (nu_hyper(hxbifj^s) + nu_hyper(hxbi^s)) * var2(hxbi^s) * (var(hxbifj^s)-var(hxbibj^s)))
221
222 !print*, "SECOND CROSS DERIV2 IXO ", ixO^L
223
224 end subroutine second_cross_deriv2
225
226end module mod_hyperdiffusivity
Module with geometry-related routines (e.g., divergence, curl)
Definition mod_geometry.t:2
integer coordinate
Definition mod_geometry.t:7
integer, parameter cartesian
Definition mod_geometry.t:8
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, parameter ndim
Number of spatial dimensions for grid variables.
double precision, dimension(:), allocatable, parameter d
integer nghostcells
Number of ghost cells surrounding a grid.
double precision, dimension(^nd) dxlevel
store unstretched cell size of current level
Subroutines for Roe-type Riemann solver for HD.
subroutine second_same_deriv(ixil, ixol, nu_hyper, var, idimm, res)
subroutine div_vel_coeff(ixil, ixol, vel, idimm, nu_vel)
subroutine second_cross_deriv2(ixil, ixol, nu_hyper, var2, var, idimm, idimm2, res)
subroutine hyp_coeff(ixil, ixol, var, idimm, nu_hyp)
subroutine second_cross_deriv(ixil, ixol, nu_hyper, var, idimm, idimm2, res)
subroutine second_same_deriv2(ixil, ixol, nu_hyper, var2, var, idimm, res)