22 call mpistop(
"Hyperdiffusivity only implemented for Cartesian uniform grid")
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)
37 double precision :: tmp(ixI^S), tmp2(ixI^S)
38 integer :: hx^L, hx1f^L, hx1b^L, hx2b^L
40 double precision,
parameter :: eps=1e-8
43 hxmin^
d=iximin^
d+2*
kr(idimm,^
d);
44 hxmax^
d=iximax^
d-
kr(idimm,^
d);
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);
51 tmp(hx^s)= abs(3d0 * (var(hx^s) - var(hx1b^s)) - (var(hx1f^s)-var(hx2b^s)))
53 tmp2(hx^s)= abs((var(hx^s) - var(hx1b^s)))
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);
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)
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
81 ixomin^
d=ixomin^
d+
kr(idimm,^
d);
87 hx1b^l=ixo^l-
kr(ii,^
d);
89 nu_vel(ixo^s) = nu_vel(ixo^s) + (vel(ixo^s,ii) - vel(hx1b^s,ii))/
dxlevel(ii)
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))
95 where(nu_vel(ixo^s) < 0d0)
96 nu_vel(ixo^s) = -nu_vel(ixo^s)
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)
118 integer :: hxf^L, hxb^L
123 hxf^l=ixo^l+
kr(idimm,^
d);
124 hxb^l=ixo^l-
kr(idimm,^
d);
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)))
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)
143 integer :: hxf^L, hxb^L
148 hxf^l=ixo^l+
kr(idimm,^
d);
149 hxb^l=ixo^l-
kr(idimm,^
d);
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)))
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)
169 integer :: hxfi^L, hxbi^L, hxfifj^L, hxbibj^L, hxfibj^L, hxbifj^L
174 hxfi^l=ixo^l+
kr(idimm2,^
d);
175 hxbi^l=ixo^l-
kr(idimm2,^
d);
177 hxfifj^l=hxfi^l+
kr(idimm,^
d);
178 hxfibj^l=hxfi^l-
kr(idimm,^
d);
180 hxbifj^l=hxbi^l+
kr(idimm,^
d);
181 hxbibj^l=hxbi^l-
kr(idimm,^
d);
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)))
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)
204 integer :: hxfi^L, hxbi^L, hxfifj^L, hxbibj^L, hxfibj^L, hxbifj^L
209 hxfi^l=ixo^l+
kr(idimm2,^
d);
210 hxbi^l=ixo^l-
kr(idimm2,^
d);
212 hxfifj^l=hxfi^l+
kr(idimm,^
d);
213 hxfibj^l=hxfi^l-
kr(idimm,^
d);
215 hxbifj^l=hxbi^l+
kr(idimm,^
d);
216 hxbibj^l=hxbi^l-
kr(idimm,^
d);
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)))
Module with geometry-related routines (e.g., divergence, curl)
integer, parameter cartesian
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 hyperdiffusivity_init()
subroutine second_same_deriv2(ixil, ixol, nu_hyper, var2, var, idimm, res)