MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_finite_difference.t
Go to the documentation of this file.
1 !> Module with finite difference methods for fluxes
3 
4  implicit none
5  private
6 
7  public :: fd
8  public :: centdiff
9 
10 contains
11 
12  subroutine fd(qdt,dtfactor,ixI^L,ixO^L,idims^LIM,qtC,sCT,qt,snew,fC,fE,dxs,x)
13  use mod_physics
14  use mod_source, only: addsource2
17  use mod_usr_methods
18 
19  double precision, intent(in) :: qdt, dtfactor, 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
22 
23  type(state) :: sct, snew
24  double precision, dimension(ixI^S,1:nwflux,1:ndim), intent(out) :: fc
25  double precision, dimension(ixI^S,sdim:3) :: fe
26 
27  double precision, dimension(ixI^S,1:nwflux) :: fct
28  double precision, dimension(ixI^S,1:nw) :: fm, fp, fmr, fpl, wprim
29  ! left and right constructed status in conservative form
30  double precision, dimension(ixI^S,1:nw) :: wlc, wrc
31  ! left and right constructed status in primitive form, needed for better performance
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
40 
41  associate(wct=>sct%w,wnew=>snew%w)
42 
43  fc=0.d0
44  wprim=wct
45  call phys_to_primitive(ixi^l,ixi^l,wprim,x)
46 
47  b0i = 0 ! fd uses centered values in phys_get_flux
48  do idims= idims^lim
49 
50  !b0i=idims
51 
52  ! Get fluxes for the whole grid (mesh+nghostcells)
53  {^d& ixmin^d = ixomin^d - nghostcells * kr(idims,^d)\}
54  {^d& ixmax^d = ixomax^d + nghostcells * kr(idims,^d)\}
55 
56  hxo^l=ixo^l-kr(idims,^d);
57 
58  if(stagger_grid) then
59  ! ct needs all transverse cells
60  ixcmax^d=ixomax^d+nghostcells-nghostcells*kr(idims,^d); ixcmin^d=hxomin^d-nghostcells+nghostcells*kr(idims,^d);
61  ixmax^d=ixmax^d+nghostcells-nghostcells*kr(idims,^d); ixmin^d=ixmin^d-nghostcells+nghostcells*kr(idims,^d);
62  kxcmin^d=iximin^d; kxcmax^d=iximax^d-kr(idims,^d);
63  kxr^l=kxc^l+kr(idims,^d);
64  ! wRp and wLp are defined at the same locations, and will correspond to
65  ! the left and right reconstructed values at a cell face. Their indexing
66  ! is similar to cell-centered values, but in direction idims they are
67  ! shifted half a cell towards the 'lower' direction.
68  wrp(kxc^s,1:nw)=wprim(kxr^s,1:nw)
69  wlp(kxc^s,1:nw)=wprim(kxc^s,1:nw)
70  else
71  ! ixC is centered index in the idims direction from ixOmin-1/2 to ixOmax+1/2
72  ixcmax^d=ixomax^d; ixcmin^d=hxomin^d;
73  end if
74 
75  call phys_get_flux(wct,wprim,x,ixg^ll,ix^l,idims,fct)
76 
77  do iw=iwstart,nwflux
78  ! Lax-Friedrich splitting:
79  fp(ix^s,iw) = half * (fct(ix^s,iw) + tvdlfeps * cmax_global * wct(ix^s,iw))
80  fm(ix^s,iw) = half * (fct(ix^s,iw) - tvdlfeps * cmax_global * wct(ix^s,iw))
81  end do ! iw loop
82 
83  ! now do the reconstruction of fp and fm:
84  call reconstructl(ixi^l,ixc^l,idims,fp,fpl)
85  call reconstructr(ixi^l,ixc^l,idims,fm,fmr)
86 
87  fc(ixc^s,iwstart:nwflux,idims) = fpl(ixc^s,iwstart:nwflux) + fmr(ixc^s,iwstart:nwflux)
88 
89  if(stagger_grid) then
90  ! apply limited reconstruction for left and right status at cell interfaces
91  call reconstruct_lr(ixi^l,ixc^l,ixc^l,idims,wprim,wlc,wrc,wlp,wrp,x,dxs(idims))
92  if(h_correction) then
93  call phys_get_h_speed(wprim,x,ixi^l,ixo^l,idims,hspeed)
94  end if
95  call phys_get_cbounds(wlc,wrc,wlp,wrp,x,ixi^l,ixc^l,idims,hspeed,cmaxc,cminc)
96  call phys_get_ct_velocity(vcts,wlp,wrp,ixi^l,ixc^l,idims,cmaxc,cminc)
97  end if
98 
99  end do !idims loop
100  !b0i=0
101 
102  if(stagger_grid) call phys_update_faces(ixi^l,ixo^l,qt,qdt,wprim,fc,fe,sct,snew,vcts)
103 
104  if(slab_uniform) then
105  dxinv=-qdt/dxs
106  do idims= idims^lim
107  hxo^l=ixo^l-kr(idims,^d);
108  if(local_timestep) then
109  do iw=iwstart,nwflux
110  fc(ixi^s,iw,idims) = -block%dt(ixi^s)*dtfactor/dxs(idims) * fc(ixi^s,iw,idims)
111  wnew(ixo^s,iw)=wnew(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
112  end do ! iw loop
113  else
114  do iw=iwstart,nwflux
115  fc(ixi^s,iw,idims) = dxinv(idims) * fc(ixi^s,iw,idims)
116  wnew(ixo^s,iw)=wnew(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
117  end do ! iw loop
118  endif
119  end do ! Next idims
120  else
121  inv_volume=1.d0/block%dvolume(ixo^s)
122  do idims= idims^lim
123  hxo^l=ixo^l-kr(idims,^d);
124  if(local_timestep) then
125  do iw=iwstart,nwflux
126  fc(ixi^s,iw,idims)=-block%dt(ixi^s)*dtfactor*fc(ixi^s,iw,idims)*block%surfaceC(ixi^s,idims)
127  wnew(ixo^s,iw)=wnew(ixo^s,iw)+ &
128  (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
129  end do ! iw loop
130  else
131  do iw=iwstart,nwflux
132  fc(ixi^s,iw,idims)=-qdt*fc(ixi^s,iw,idims)*block%surfaceC(ixi^s,idims)
133  wnew(ixo^s,iw)=wnew(ixo^s,iw)+ &
134  (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
135  end do ! iw loop
136  endif ! local_timestep
137  end do ! Next idims
138  end if
139 
140  if (.not.slab.and.idimsmin==1) call phys_add_source_geom(qdt,dtfactor,ixi^l,ixo^l,wct,wnew,x)
141 
142  if(stagger_grid) call phys_face_to_center(ixo^l,snew)
143 
144  ! check and optionally correct unphysical values
145  if(fix_small_values) then
146  call phys_handle_small_values(.false.,wnew,x,ixi^l,ixo^l,'fd')
147  if(crash) then
148  ! replace erroneous values with values at previous step
149  wnew=wct
150  if(stagger_grid) snew%ws=sct%ws
151  end if
152  endif
153 
154  call addsource2(qdt*dble(idimsmax-idimsmin+1)/dble(ndim), &
155  dtfactor*dble(idimsmax-idimsmin+1)/dble(ndim), &
156  ixi^l,ixo^l,1,nw,qtc,wct,wprim,qt,wnew,x,.false.,active)
157 
158  end associate
159 
160  end subroutine fd
161 
162  subroutine reconstructl(ixI^L,iL^L,idims,w,wLC)
164  use mod_mp5
165  use mod_limiter
166  use mod_comm_lib, only: mpistop
167 
168  integer, intent(in) :: ixI^L, iL^L, idims
169  double precision, intent(in) :: w(ixI^S,1:nw)
170 
171  double precision, intent(out) :: wLC(ixI^S,1:nw)
172 
173  double precision :: ldw(ixI^S), dwC(ixI^S)
174  integer :: jxR^L, ixC^L, jxC^L, kxC^L, iw
175  double precision :: a2max
176 
177  select case (type_limiter(block%level))
178  case (limiter_mp5)
179  call mp5limiterl(ixi^l,il^l,idims,w,wlc)
180  case (limiter_weno5)
181  call weno5limiterl(ixi^l,il^l,idims,w,wlc,1)
182  case (limiter_weno5nm)
183  call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,1)
184  case (limiter_wenoz5)
185  call weno5limiterl(ixi^l,il^l,idims,w,wlc,2)
186  case (limiter_wenoz5nm)
187  call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,2)
188  case (limiter_wenozp5)
189  call weno5limiterl(ixi^l,il^l,idims,w,wlc,3)
190  case (limiter_wenozp5nm)
191  call weno5nmlimiterl(ixi^l,il^l,idims,w,wlc,3)
192  case default
193 
194  kxcmin^d=iximin^d; kxcmax^d=iximax^d-kr(idims,^d);
195 
196  wlc(kxc^s,iwstart:nwflux) = w(kxc^s,iwstart:nwflux)
197 
198  jxr^l=il^l+kr(idims,^d);
199 
200  ixcmax^d=jxrmax^d; ixcmin^d=ilmin^d-kr(idims,^d);
201  jxc^l=ixc^l+kr(idims,^d);
202 
203  do iw=iwstart,nwflux
204  dwc(ixc^s)=w(jxc^s,iw)-w(ixc^s,iw)
205 
206  if(need_global_a2max) then
207  a2max=a2max_global(idims)
208  else
209  select case(idims)
210  case(1)
211  a2max=schmid_rad1
212  {^iftwod
213  case(2)
214  a2max=schmid_rad2}
215  {^ifthreed
216  case(2)
217  a2max=schmid_rad2
218  case(3)
219  a2max=schmid_rad3}
220  case default
221  call mpistop("idims is wrong in mod_limiter")
222  end select
223  end if
224 
225  call dwlimiter2(dwc,ixi^l,ixc^l,idims,type_limiter(block%level),ldw=ldw,a2max=a2max)
226 
227  wlc(il^s,iw)=wlc(il^s,iw)+half*ldw(il^s)
228  end do
229  end select
230 
231  end subroutine reconstructl
232 
233  subroutine reconstructr(ixI^L,iL^L,idims,w,wRC)
235  use mod_mp5
236  use mod_limiter
237  use mod_comm_lib, only: mpistop
238 
239  integer, intent(in) :: ixI^L, iL^L, idims
240  double precision, intent(in) :: w(ixI^S,1:nw)
241 
242  double precision, intent(out) :: wRC(ixI^S,1:nw)
243 
244  double precision :: rdw(ixI^S), dwC(ixI^S)
245  integer :: jxR^L, ixC^L, jxC^L, kxC^L, kxR^L, iw
246  double precision :: a2max
247 
248  select case (type_limiter(block%level))
249  case (limiter_mp5)
250  call mp5limiterr(ixi^l,il^l,idims,w,wrc)
251  case (limiter_weno5)
252  call weno5limiterr(ixi^l,il^l,idims,w,wrc,1)
253  case (limiter_weno5nm)
254  call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,1)
255  case (limiter_wenoz5)
256  call weno5limiterr(ixi^l,il^l,idims,w,wrc,2)
257  case (limiter_wenoz5nm)
258  call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,2)
259  case (limiter_wenozp5)
260  call weno5limiterr(ixi^l,il^l,idims,w,wrc,3)
261  case (limiter_wenozp5nm)
262  call weno5nmlimiterr(ixi^l,il^l,idims,w,wrc,3)
263  case default
264 
265  kxcmin^d=iximin^d; kxcmax^d=iximax^d-kr(idims,^d);
266  kxr^l=kxc^l+kr(idims,^d);
267 
268  wrc(kxc^s,iwstart:nwflux)=w(kxr^s,iwstart:nwflux)
269 
270  jxr^l=il^l+kr(idims,^d);
271  ixcmax^d=jxrmax^d; ixcmin^d=ilmin^d-kr(idims,^d);
272  jxc^l=ixc^l+kr(idims,^d);
273 
274  do iw=iwstart,nwflux
275  dwc(ixc^s)=w(jxc^s,iw)-w(ixc^s,iw)
276 
277  if(need_global_a2max) then
278  a2max=a2max_global(idims)
279  else
280  select case(idims)
281  case(1)
282  a2max=schmid_rad1
283  {^iftwod
284  case(2)
285  a2max=schmid_rad2}
286  {^ifthreed
287  case(2)
288  a2max=schmid_rad2
289  case(3)
290  a2max=schmid_rad3}
291  case default
292  call mpistop("idims is wrong in mod_limiter")
293  end select
294  end if
295 
296  call dwlimiter2(dwc,ixi^l,ixc^l,idims,type_limiter(block%level),rdw=rdw,a2max=a2max)
297 
298  wrc(il^s,iw)=wrc(il^s,iw)-half*rdw(jxr^s)
299  end do
300  end select
301 
302  end subroutine reconstructr
303 
304  subroutine centdiff(method,qdt,dtfactor,ixI^L,ixO^L,idims^LIM,qtC,sCT,qt,s,fC,fE,dxs,x)
305 
306  ! Advance the flow variables from global_time to global_time+qdt within ixO^L by
307  ! fourth order centered differencing in space
308  ! for the dw/dt+dF_i(w)/dx_i=S type equation.
309  ! wCT contains the time centered variables at time qtC for flux and source.
310  ! w is the old value at qt on input and the new value at qt+qdt on output.
311  use mod_physics
314  use mod_source, only: addsource2
315  use mod_usr_methods
316  use mod_variables
317  use mod_comm_lib, only: mpistop
318 
319  integer, intent(in) :: method
320  integer, intent(in) :: ixi^l, ixo^l, idims^lim
321  double precision, intent(in) :: qdt, dtfactor, qtc, qt, dxs(ndim)
322  type(state) :: sct, s
323  double precision, intent(in) :: x(ixi^s,1:ndim)
324  double precision :: fc(ixi^s,1:nwflux,1:ndim)
325  double precision :: fe(ixi^s,sdim:3)
326 
327  double precision :: v(ixi^s,ndim), f(ixi^s, nwflux)
328 
329  double precision, dimension(ixI^S,1:nw) :: wprim, wlc, wrc
330  ! left and right constructed status in primitive form, needed for better performance
331  double precision, dimension(ixI^S,1:nw) :: wlp, wrp
332  double precision, dimension(ixI^S) :: vlc, cmaxlc, cmaxrc
333  double precision, dimension(ixI^S,1:number_species) :: cmaxc
334  double precision, dimension(ixI^S,1:number_species) :: cminc
335  double precision, dimension(ixI^S) :: hspeed
336  double precision, dimension(ixO^S) :: inv_volume
337 
338  double precision :: dxinv(1:ndim)
339  integer :: idims, iw, ix^l, hxo^l, ixc^l, jxc^l, hxc^l, kxc^l, kkxc^l, kkxr^l
340  type(ct_velocity) :: vcts
341  logical :: transport, new_cmax, patchw(ixi^s), active
342 
343  associate(wct=>sct%w,w=>s%w)
344  ! two extra layers are needed in each direction for which fluxes are added.
345  ix^l=ixo^l;
346  do idims= idims^lim
347  ix^l=ix^l^ladd2*kr(idims,^d);
348  end do
349 
350  if (ixi^l^ltix^l|.or.|.or.) then
351  call mpistop("Error in centdiff: Non-conforming input limits")
352  end if
353 
354  fc=0.d0
355  wprim=wct
356  call phys_to_primitive(ixi^l,ixi^l,wprim,x)
357 
358  ! get fluxes
359  do idims= idims^lim
360  b0i=idims
361 
362  ix^l=ixo^l^ladd2*kr(idims,^d);
363  hxo^l=ixo^l-kr(idims,^d);
364 
365  if(stagger_grid) then
366  ! ct needs all transverse cells
367  ixcmax^d=ixomax^d+nghostcells-nghostcells*kr(idims,^d);
368  ixcmin^d=hxomin^d-nghostcells+nghostcells*kr(idims,^d);
369  ixmax^d=ixmax^d+nghostcells-nghostcells*kr(idims,^d);
370  ixmin^d=ixmin^d-nghostcells+nghostcells*kr(idims,^d);
371  else
372  ! ixC is centered index in the idims direction from ixOmin-1/2 to ixOmax+1/2
373  ixcmax^d=ixomax^d; ixcmin^d=hxomin^d;
374  end if
375  hxc^l=ixc^l-kr(idims,^d);
376  jxc^l=ixc^l+kr(idims,^d);
377  kxc^l=ixc^l+2*kr(idims,^d);
378 
379  kkxcmin^d=iximin^d; kkxcmax^d=iximax^d-kr(idims,^d);
380  kkxr^l=kkxc^l+kr(idims,^d);
381  wrp(kkxc^s,1:nwflux)=wprim(kkxr^s,1:nwflux)
382  wlp(kkxc^s,1:nwflux)=wprim(kkxc^s,1:nwflux)
383 
384  ! apply limited reconstruction for left and right status at cell interfaces
385  call reconstruct_lr(ixi^l,ixc^l,ixc^l,idims,wprim,wlc,wrc,wlp,wrp,x,dxs(idims))
386 
387  if(stagger_grid) then
388  if(h_correction) then
389  call phys_get_h_speed(wprim,x,ixi^l,ixo^l,idims,hspeed)
390  end if
391  call phys_get_cbounds(wlc,wrc,wlp,wrp,x,ixi^l,ixc^l,idims,hspeed,cmaxc,cminc)
392  call phys_get_ct_velocity(vcts,wlp,wrp,ixi^l,ixc^l,idims,cmaxc,cminc)
393  end if
394 
395  ! Calculate velocities from upwinded values
396  call phys_get_cmax(wlc,x,ixg^ll,ixc^l,idims,cmaxlc)
397  call phys_get_cmax(wrc,x,ixg^ll,ixc^l,idims,cmaxrc)
398  ! now take the maximum of left and right states
399  vlc(ixc^s)=max(cmaxrc(ixc^s),cmaxlc(ixc^s))
400 
401  call phys_get_flux(wct,wprim,x,ixi^l,ix^l,idims,f)
402 
403  ! Center flux to interface
404  if(method==fs_cd) then
405  fc(ixc^s,iwstart:nwflux,idims)=half*(f(ixc^s,iwstart:nwflux)+f(jxc^s,iwstart:nwflux))
406  else
407  ! f_i+1/2= (-f_(i+2) +7 f_(i+1) + 7 f_i - f_(i-1))/12
408  fc(ixc^s,iwstart:nwflux,idims)=(-f(kxc^s,iwstart:nwflux)+7.0d0*(f(jxc^s,iwstart:nwflux) + &
409  f(ixc^s,iwstart:nwflux))-f(hxc^s,iwstart:nwflux))/12.0d0
410 
411  do iw=iwstart,nwflux
412  ! add rempel dissipative flux, only second order version for now
413  fc(ixc^s,iw,idims)=fc(ixc^s,iw,idims)-tvdlfeps*half*vlc(ixc^s) &
414  *(wrc(ixc^s,iw)-wlc(ixc^s,iw))
415  end do
416  end if
417 
418  end do !next idims
419  b0i=0
420 
421  if(stagger_grid) call phys_update_faces(ixi^l,ixo^l,qt,qdt,wprim,fc,fe,sct,s,vcts)
422 
423  if(slab_uniform) then
424  dxinv=-qdt/dxs
425  do idims= idims^lim
426  hxo^l=ixo^l-kr(idims,^d);
427  if(local_timestep) then
428  do iw=iwstart,nwflux
429  fc(ixi^s,iw,idims)=-block%dt(ixi^s)*dtfactor/dxs(idims)*fc(ixi^s,iw,idims)
430  ! result: f_(i+1/2)-f_(i-1/2) = [-f_(i+2)+8(f_(i+1)-f_(i-1))+f_(i-2)]/12
431  w(ixo^s,iw)=w(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
432  end do !next iw
433  else
434  do iw=iwstart,nwflux
435  fc(ixi^s,iw,idims)=dxinv(idims)*fc(ixi^s,iw,idims)
436  ! result: f_(i+1/2)-f_(i-1/2) = [-f_(i+2)+8(f_(i+1)-f_(i-1))+f_(i-2)]/12
437  w(ixo^s,iw)=w(ixo^s,iw)+(fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))
438  end do !next iw
439  endif
440  end do ! Next idims
441  else
442  inv_volume=1.d0/block%dvolume
443  do idims= idims^lim
444  hxo^l=ixo^l-kr(idims,^d);
445  if(local_timestep) then
446  do iw=iwstart,nwflux
447  fc(ixi^s,iw,idims)=-block%dt(ixi^s)*dtfactor*block%surfaceC(ixi^s,idims)*fc(ixi^s,iw,idims)
448  w(ixo^s,iw)=w(ixo^s,iw)+ &
449  (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
450  end do !next iw
451  else
452  do iw=iwstart,nwflux
453  fc(ixi^s,iw,idims)=-qdt*block%surfaceC(ixi^s,idims)*fc(ixi^s,iw,idims)
454  w(ixo^s,iw)=w(ixo^s,iw)+ &
455  (fc(ixo^s,iw,idims)-fc(hxo^s,iw,idims))*inv_volume(ixo^s)
456  end do !next iw
457  endif
458  end do ! Next idims
459  end if
460 
461  if (.not.slab.and.idimsmin==1) call phys_add_source_geom(qdt,dtfactor,ixi^l,ixo^l,wct,w,x)
462 
463  if(stagger_grid) call phys_face_to_center(ixo^l,s)
464 
465  ! check and optionally correct unphysical values
466  if(fix_small_values) then
467  call phys_handle_small_values(.false.,w,x,ixi^l,ixo^l,'centdiff')
468  endif
469 
470  call addsource2(qdt*dble(idimsmax-idimsmin+1)/dble(ndim), &
471  dtfactor*dble(idimsmax-idimsmin+1)/dble(ndim), &
472  ixi^l,ixo^l,1,nw,qtc,wct,wprim,qt,w,x,.false.,active)
473 
474  end associate
475  end subroutine centdiff
476 
477 end module mod_finite_difference
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: mod_comm_lib.t:208
Module with finite difference methods for fluxes.
subroutine, public centdiff(method, qdt, dtfactor, ixIL, ixOL, idimsLIM, qtC, sCT, qt, s, fC, fE, dxs, x)
subroutine reconstructr(ixIL, iLL, idims, w, wRC)
subroutine reconstructl(ixIL, iLL, idims, w, wLC)
subroutine, public fd(qdt, dtfactor, ixIL, ixOL, idimsLIM, qtC, sCT, qt, snew, 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.
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 local_timestep
each cell has its own timestep or not
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.
integer, parameter sdim
starting dimension for electric field
logical crash
Save a snapshot before crash a run met unphysical values.
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.
Definition: mod_limiter.t:2
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...
Definition: mod_limiter.t:129
double precision d
Definition: mod_limiter.t:14
integer, parameter limiter_wenozp5
Definition: mod_limiter.t:35
integer, parameter limiter_weno5
Definition: mod_limiter.t:31
integer, parameter limiter_wenoz5
Definition: mod_limiter.t:33
integer, parameter limiter_wenoz5nm
Definition: mod_limiter.t:34
integer, parameter limiter_wenozp5nm
Definition: mod_limiter.t:36
integer, parameter limiter_mp5
Definition: mod_limiter.t:28
integer, parameter limiter_weno5nm
Definition: mod_limiter.t:32
Module containing the MP5 (fifth order) flux scheme.
Definition: mod_mp5.t:2
subroutine, public mp5limiterr(ixIL, iLL, idims, w, wRC)
Definition: mod_mp5.t:301
subroutine, public mp5limiterl(ixIL, iLL, idims, w, wLC)
Definition: mod_mp5.t:204
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_get_ct_velocity), pointer phys_get_ct_velocity
Definition: mod_physics.t:83
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:59
procedure(sub_small_values), pointer phys_handle_small_values
Definition: mod_physics.t:82
procedure(sub_get_flux), pointer phys_get_flux
Definition: mod_physics.t:68
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition: mod_physics.t:67
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition: mod_physics.t:72
procedure(sub_update_faces), pointer phys_update_faces
Definition: mod_physics.t:84
procedure(sub_face_to_center), pointer phys_face_to_center
Definition: mod_physics.t:85
procedure(sub_get_h_speed), pointer phys_get_h_speed
Definition: mod_physics.t:66
procedure(sub_get_cmax), pointer phys_get_cmax
Definition: mod_physics.t:61
Module for handling split source terms (split from the fluxes)
Definition: mod_source.t:2
subroutine, public addsource2(qdt, dtfactor, ixIL, ixOL, iwLIM, qtC, wCT, wCTprim, qt, w, x, qsourcesplit, src_active)
Add source within ixO for iws: w=w+qdt*S[wCT].
Definition: mod_source.t:132
Module with all the methods that users can customize in AMRVAC.
integer nw
Total number of variables.
Definition: mod_variables.t:23
integer iwstart
Definition: mod_variables.t:38
integer nwflux
Number of flux variables.
Definition: mod_variables.t:8