12 subroutine fd(qdt,ixI^L,ixO^L,idims^LIM,qtC,sCT,qt,snew,fC,fE,dxs,x)
19 double precision,
intent(in) :: qdt, qtc, qt, dxs(
ndim)
20 integer,
intent(in) :: ixi^
l, ixo^
l, idims^lim
21 double precision,
dimension(ixI^S,1:ndim),
intent(in) :: x
23 type(state) :: sct, snew
24 double precision,
dimension(ixI^S,1:nwflux,1:ndim),
intent(out) :: fc
25 double precision,
dimension(ixI^S,7-2*ndim:3) :: fe
27 double precision,
dimension(ixI^S,1:nwflux) :: fct
28 double precision,
dimension(ixI^S,1:nw) :: fm, fp, fmr, fpl, wprim
30 double precision,
dimension(ixI^S,1:nw) :: wlc, wrc
32 double precision,
dimension(ixI^S,1:nw) :: wlp, wrp
33 double precision,
dimension(ixI^S) :: cmaxc, cminc
34 double precision,
dimension(ixI^S) :: hspeed
35 double precision,
dimension(ixO^S) :: inv_volume
36 double precision,
dimension(1:ndim) :: dxinv
37 logical :: transport, active
38 integer :: idims, iw, ixc^
l, ix^
l, hxo^
l, kxc^
l, kxr^
l
39 type(ct_velocity) :: vcts
41 associate(wct=>sct%w,wnew=>snew%w)
56 hxo^
l=ixo^
l-
kr(idims,^
d);
62 kxcmin^
d=iximin^
d; kxcmax^
d=iximax^
d-
kr(idims,^
d);
63 kxr^
l=kxc^
l+
kr(idims,^
d);
68 wrp(kxc^s,1:nw)=wprim(kxr^s,1:nw)
69 wlp(kxc^s,1:nw)=wprim(kxc^s,1:nw)
72 ixcmax^
d=ixomax^
d; ixcmin^
d=hxomin^
d;
87 fc(ixc^s,iwstart:nwflux,idims) = fpl(ixc^s,iwstart:nwflux) + fmr(ixc^s,iwstart:nwflux)
91 call reconstruct_lr(ixi^
l,ixc^
l,ixc^
l,idims,wprim,wlc,wrc,wlp,wrp,x,dxs(idims))
95 call phys_get_cbounds(wlc,wrc,wlp,wrp,x,ixi^
l,ixc^
l,idims,hspeed,cmaxc,cminc)
107 hxo^
l=ixo^
l-
kr(idims,^
d);
109 fc(ixi^s,iw,idims) = dxinv(idims) * fc(ixi^s,iw,idims)
110 wnew(ixo^s,iw)=wnew(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
114 inv_volume=1.d0/
block%dvolume(ixo^s)
116 hxo^
l=ixo^
l-
kr(idims,^
d);
118 fc(ixi^s,iw,idims)=-qdt*fc(ixi^s,iw,idims)*
block%surfaceC(ixi^s,idims)
119 wnew(ixo^s,iw)=wnew(ixo^s,iw)+ &
120 (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
135 ixi^
l,ixo^
l,1,nw,qtc,wct,qt,wnew,x,.false.,active,wprim)
150 integer,
intent(in) :: ixI^L, iL^L, idims
151 double precision,
intent(in) :: w(ixI^S,1:nw)
153 double precision,
intent(out) :: wLC(ixI^S,1:nw)
155 double precision :: ldw(ixI^S), dwC(ixI^S)
156 integer :: jxR^L, ixC^L, jxC^L, kxC^L, iw
157 double precision :: a2max
163 call weno5limiterl(ixi^l,il^l,idims,w,wlc,1)
165 call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,1)
167 call weno5limiterl(ixi^l,il^l,idims,w,wlc,2)
169 call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,2)
171 call weno5limiterl(ixi^l,il^l,idims,w,wlc,3)
173 call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,3)
176 kxcmin^
d=iximin^
d; kxcmax^
d=iximax^
d-
kr(idims,^
d);
178 wlc(kxc^s,iwstart:nwflux) = w(kxc^s,iwstart:nwflux)
180 jxr^l=il^l+
kr(idims,^
d);
182 ixcmax^
d=jxrmax^
d; ixcmin^
d=ilmin^
d-
kr(idims,^
d);
183 jxc^l=ixc^l+
kr(idims,^
d);
186 dwc(ixc^s)=w(jxc^s,iw)-w(ixc^s,iw)
203 call mpistop(
"idims is wrong in mod_limiter")
209 wlc(il^s,iw)=wlc(il^s,iw)+half*ldw(il^s)
220 integer,
intent(in) :: ixI^L, iL^L, idims
221 double precision,
intent(in) :: w(ixI^S,1:nw)
223 double precision,
intent(out) :: wRC(ixI^S,1:nw)
225 double precision :: rdw(ixI^S), dwC(ixI^S)
226 integer :: jxR^L, ixC^L, jxC^L, kxC^L, kxR^L, iw
227 double precision :: a2max
233 call weno5limiterr(ixi^l,il^l,idims,w,wrc,1)
235 call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,1)
237 call weno5limiterr(ixi^l,il^l,idims,w,wrc,2)
239 call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,2)
241 call weno5limiterr(ixi^l,il^l,idims,w,wrc,3)
243 call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,3)
246 kxcmin^
d=iximin^
d; kxcmax^
d=iximax^
d-
kr(idims,^
d);
247 kxr^l=kxc^l+
kr(idims,^
d);
249 wrc(kxc^s,iwstart:nwflux)=w(kxr^s,iwstart:nwflux)
251 jxr^l=il^l+
kr(idims,^
d);
252 ixcmax^
d=jxrmax^
d; ixcmin^
d=ilmin^
d-
kr(idims,^
d);
253 jxc^l=ixc^l+
kr(idims,^
d);
256 dwc(ixc^s)=w(jxc^s,iw)-w(ixc^s,iw)
273 call mpistop(
"idims is wrong in mod_limiter")
279 wrc(il^s,iw)=wrc(il^s,iw)-half*rdw(jxr^s)
285 subroutine centdiff(method,qdt,ixI^L,ixO^L,idims^LIM,qtC,sCT,qt,s,fC,fE,dxs,x)
299 integer,
intent(in) :: method
300 integer,
intent(in) :: ixi^
l, ixo^
l, idims^lim
301 double precision,
intent(in) :: qdt, qtc, qt, dxs(
ndim)
302 type(state) :: sct, s
303 double precision,
intent(in) :: x(ixi^s,1:
ndim)
305 double precision :: fe(ixi^s,7-2*
ndim:3)
307 double precision :: v(ixi^s,
ndim), f(ixi^s,
nwflux)
309 double precision,
dimension(ixI^S,1:nw) :: wprim, wlc, wrc
311 double precision,
dimension(ixI^S,1:nw) :: wlp, wrp
312 double precision,
dimension(ixI^S) :: vlc, cmaxlc, cmaxrc
313 double precision,
dimension(ixI^S,1:number_species) :: cmaxc
314 double precision,
dimension(ixI^S,1:number_species) :: cminc
315 double precision,
dimension(ixI^S) :: hspeed
316 double precision,
dimension(ixO^S) :: inv_volume
318 double precision :: dxinv(1:
ndim)
319 integer :: idims, iw, ix^
l, hxo^
l, ixc^
l, jxc^
l, hxc^
l, kxc^
l, kkxc^
l, kkxr^
l
320 type(ct_velocity) :: vcts
321 logical :: transport, new_cmax, patchw(ixi^s), active
323 associate(wct=>sct%w,w=>s%w)
327 ix^
l=ix^
l^ladd2*
kr(idims,^
d);
330 if (ixi^
l^ltix^
l|.or.|.or.)
then
331 call mpistop(
"Error in centdiff: Non-conforming input limits")
342 ix^
l=ixo^
l^ladd2*
kr(idims,^
d);
343 hxo^
l=ixo^
l-
kr(idims,^
d);
353 ixcmax^
d=ixomax^
d; ixcmin^
d=hxomin^
d;
355 hxc^
l=ixc^
l-
kr(idims,^
d);
356 jxc^
l=ixc^
l+
kr(idims,^
d);
357 kxc^
l=ixc^
l+2*
kr(idims,^
d);
359 kkxcmin^
d=iximin^
d; kkxcmax^
d=iximax^
d-
kr(idims,^
d);
360 kkxr^
l=kkxc^
l+
kr(idims,^
d);
365 call reconstruct_lr(ixi^
l,ixc^
l,ixc^
l,idims,wprim,wlc,wrc,wlp,wrp,x,dxs(idims))
371 call phys_get_cbounds(wlc,wrc,wlp,wrp,x,ixi^
l,ixc^
l,idims,hspeed,cmaxc,cminc)
379 vlc(ixc^s)=max(cmaxrc(ixc^s),cmaxlc(ixc^s))
384 if(method==
fs_cd)
then
393 fc(ixc^s,iw,idims)=fc(ixc^s,iw,idims)-
tvdlfeps*half*vlc(ixc^s) &
394 *(wrc(ixc^s,iw)-wlc(ixc^s,iw))
406 hxo^
l=ixo^
l-
kr(idims,^
d);
408 fc(ixi^s,iw,idims)=dxinv(idims)*fc(ixi^s,iw,idims)
410 w(ixo^s,iw)=w(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
414 inv_volume=1.d0/
block%dvolume
416 hxo^
l=ixo^
l-
kr(idims,^
d);
418 fc(ixi^s,iw,idims)=-qdt*
block%surfaceC(ixi^s,idims)*fc(ixi^s,iw,idims)
419 w(ixo^s,iw)=w(ixo^s,iw)+ &
420 (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
435 ixi^
l,ixo^
l,1,
nw,qtc,wct,qt,w,x,.false.,active,wprim)
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Module with finite difference methods for fluxes.
subroutine, public fd(qdt, ixIL, ixOL, idimsLIM, qtC, sCT, qt, snew, fC, fE, dxs, x)
subroutine reconstructr(ixIL, iLL, idims, w, wRC)
subroutine reconstructl(ixIL, iLL, idims, w, wLC)
subroutine, public centdiff(method, qdt, ixIL, ixOL, idimsLIM, qtC, sCT, qt, s, fC, fE, dxs, x)
Module with finite volume methods for fluxes.
subroutine, public reconstruct_lr(ixIL, ixLL, ixRL, idims, w, wLC, wRC, wLp, wRp, x, dxdim)
Determine the upwinded wLC(ixL) and wRC(ixR) from w. the wCT is only used when PPM is exploited.
This module contains definitions of global parameters and variables and some generic functions/subrou...
type(state), pointer block
Block pointer for using one block and its previous state.
logical h_correction
If true, do H-correction to fix the carbuncle problem at grid-aligned shocks.
double precision tvdlfeps
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical stagger_grid
True for using stagger grid.
double precision cmax_global
global fastest wave speed needed in fd scheme and glm method
integer b0i
background magnetic field location indicator
integer, dimension(:), allocatable, parameter d
logical need_global_a2max
global value for schmid scheme
logical slab
Cartesian geometry or not.
integer, dimension(:), allocatable type_limiter
Type of slope limiter used for reconstructing variables on cell edges.
integer nghostcells
Number of ghost cells surrounding a grid.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
double precision, dimension(ndim) a2max_global
global largest a2 for schmid scheme
Module with slope/flux limiters.
subroutine dwlimiter2(dwC, ixIL, ixCL, idims, typelim, ldw, rdw, a2max)
Limit the centered dwC differences within ixC for iw in direction idim. The limiter is chosen accordi...
integer, parameter limiter_wenozp5
integer, parameter limiter_weno5
integer, parameter limiter_wenoz5
integer, parameter limiter_wenoz5nm
integer, parameter limiter_wenozp5nm
integer, parameter limiter_mp5
integer, parameter limiter_weno5nm
Module containing the MP5 (fifth order) flux scheme.
subroutine, public mp5limiterr(ixIL, iLL, idims, w, wRC)
subroutine, public mp5limiterl(ixIL, iLL, idims, w, wLC)
This module defines the procedures of a physics module. It contains function pointers for the various...
procedure(sub_get_ct_velocity), pointer phys_get_ct_velocity
procedure(sub_convert), pointer phys_to_primitive
procedure(sub_small_values), pointer phys_handle_small_values
procedure(sub_get_flux), pointer phys_get_flux
procedure(sub_get_cbounds), pointer phys_get_cbounds
logical phys_solve_eaux
Solve internal energy and total energy equations.
procedure(sub_add_source_geom), pointer phys_add_source_geom
procedure(sub_update_faces), pointer phys_update_faces
procedure(sub_energy_synchro), pointer phys_energy_synchro
procedure(sub_face_to_center), pointer phys_face_to_center
procedure(sub_get_h_speed), pointer phys_get_h_speed
procedure(sub_get_cmax), pointer phys_get_cmax
Module for handling split source terms (split from the fluxes)
subroutine, public addsource2(qdt, ixIL, ixOL, iwLIM, qtC, wCT, qt, w, x, qsourcesplit, src_active, wCTprim)
Add source within ixO for iws: w=w+qdt*S[wCT].
Module with all the methods that users can customize in AMRVAC.
integer nw
Total number of variables.
integer nwflux
Number of flux variables.