MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_physics.t
Go to the documentation of this file.
1 !> This module defines the procedures of a physics module. It contains function
2 !> pointers for the various supported routines. An actual physics module has to
3 !> set these pointers to its implementation of these routines.
4 module mod_physics
5  use mod_global_parameters, only: name_len, max_nw
9 
10  implicit none
11  public
12 
13  !> String describing the physics type of the simulation
14  character(len=name_len) :: physics_type = ""
15 
16  !> To use wider stencils in flux calculations. A value of 1 will extend it by
17  !> one cell in both directions, in any dimension
18  integer :: phys_wider_stencil = 0
19 
20  !> Whether the physics routines require diagonal ghost cells, for example for
21  !> computing a curl.
22  logical :: phys_req_diagonal = .true.
23 
24  !> Solve energy equation or not
25  logical :: phys_energy=.false.
26 
27  !> Solve total energy equation or not
28  logical :: phys_total_energy=.false.
29 
30  !> Solve internal enery instead of total energy
31  logical :: phys_internal_e=.false.
32 
33  !> Solve internal energy and total energy equations
34  logical :: phys_solve_eaux=.false.
35 
36  !> Array per direction per variable, which can be used to specify that certain
37  !> fluxes have to be treated differently
38  integer, allocatable :: flux_type(:, :)
39 
40  !> Indicates a normal flux
41  integer, parameter :: flux_default = 0
42  !> Indicates the flux should be treated with tvdlf
43  integer, parameter :: flux_tvdlf = 1
44  !> Indicates dissipation should be omitted
45  integer, parameter :: flux_no_dissipation = 2
46  !> Indicates the flux should be specially treated
47  integer, parameter :: flux_special = 3
48 
49  !> Type for special methods defined per variable
51  integer :: test
52  !> If this is set, use the routine as a capacity function when adding fluxes
53  procedure(sub_get_var), pointer, nopass :: inv_capacity => null()
54  end type iw_methods
55 
56  !> Special methods defined per variable
57  type(iw_methods) :: phys_iw_methods(max_nw)
58 
59  procedure(sub_check_params), pointer :: phys_check_params => null()
60  procedure(sub_convert), pointer :: phys_to_conserved => null()
61  procedure(sub_convert), pointer :: phys_to_primitive => null()
62  procedure(sub_convert), pointer :: phys_ei_to_e => null()
63  procedure(sub_convert), pointer :: phys_e_to_ei => null()
64  procedure(sub_modify_wlr), pointer :: phys_modify_wlr => null()
65  procedure(sub_get_cmax), pointer :: phys_get_cmax => null()
66  procedure(sub_get_a2max), pointer :: phys_get_a2max => null()
67  procedure(sub_get_tcutoff), pointer :: phys_get_tcutoff => null()
68  procedure(sub_get_cbounds), pointer :: phys_get_cbounds => null()
69  procedure(sub_get_flux), pointer :: phys_get_flux => null()
70  procedure(sub_energy_synchro), pointer :: phys_energy_synchro => null()
71  procedure(sub_get_v_idim), pointer :: phys_get_v_idim => null()
72  procedure(sub_get_dt), pointer :: phys_get_dt => null()
73  procedure(sub_add_source_geom), pointer :: phys_add_source_geom => null()
74  procedure(sub_add_source), pointer :: phys_add_source => null()
75  procedure(sub_global_source), pointer :: phys_global_source_before => null()
76  procedure(sub_global_source), pointer :: phys_global_source_after => null()
77  procedure(sub_get_aux), pointer :: phys_get_aux => null()
78  procedure(sub_check_w), pointer :: phys_check_w => null()
79  procedure(sub_get_pthermal), pointer :: phys_get_pthermal => null()
80  procedure(sub_boundary_adjust), pointer :: phys_boundary_adjust => null()
81  procedure(sub_write_info), pointer :: phys_write_info => null()
82  procedure(sub_angmomfix), pointer :: phys_angmomfix => null()
83  procedure(sub_small_values), pointer :: phys_handle_small_values => null()
84  procedure(sub_update_faces), pointer :: phys_update_faces => null()
85  procedure(sub_face_to_center), pointer :: phys_face_to_center => null()
86  procedure(sub_implicit_update), pointer :: phys_implicit_update => null()
87  procedure(sub_evaluate_implicit),pointer:: phys_evaluate_implicit => null()
88 
89  abstract interface
90 
91  subroutine sub_check_params
92  end subroutine sub_check_params
93 
96  end subroutine sub_boundary_adjust
97 
98  subroutine sub_convert(ixI^L, ixO^L, w, x)
100  integer, intent(in) :: ixI^L, ixO^L
101  double precision, intent(inout) :: w(ixi^s, nw)
102  double precision, intent(in) :: x(ixi^s, 1:^nd)
103  end subroutine sub_convert
104 
105  subroutine sub_modify_wlr(ixI^L, ixO^L, qt, wLC, wRC, wLp, wRp, s, idir)
107  integer, intent(in) :: ixI^L, ixO^L, idir
108  double precision, intent(in) :: qt
109  double precision, intent(inout) :: wLC(ixi^s,1:nw), wRC(ixi^s,1:nw)
110  double precision, intent(inout) :: wLp(ixi^s,1:nw), wRp(ixi^s,1:nw)
111  type(state) :: s
112  end subroutine sub_modify_wlr
113 
114  subroutine sub_get_cmax(w, x, ixI^L, ixO^L, idim, cmax)
116  integer, intent(in) :: ixI^L, ixO^L, idim
117  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
118  double precision, intent(inout) :: cmax(ixi^s)
119  end subroutine sub_get_cmax
120 
121  subroutine sub_get_a2max(w, x, ixI^L, ixO^L, a2max)
123  integer, intent(in) :: ixI^L, ixO^L
124  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
125  double precision, intent(inout) :: a2max(ndim)
126  end subroutine sub_get_a2max
127 
128  subroutine sub_get_tcutoff(ixI^L,ixO^L,w,x,tco_local,Tmax_local)
130  integer, intent(in) :: ixI^L, ixO^L
131  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
132  double precision, intent(out) :: tco_local, Tmax_local
133  end subroutine sub_get_tcutoff
134 
135  subroutine sub_get_v_idim(w,x,ixI^L,ixO^L,idim,v)
137 
138  integer, intent(in) :: ixI^L, ixO^L, idim
139  double precision, intent(in) :: w(ixi^s,nw), x(ixi^s,1:^nd)
140  double precision, intent(out) :: v(ixi^s)
141 
142  end subroutine sub_get_v_idim
143 
144  subroutine sub_get_cbounds(wLC, wRC, wLp, wRp, x, ixI^L, ixO^L, idim, cmax, cmin)
146  integer, intent(in) :: ixI^L, ixO^L, idim
147  double precision, intent(in) :: wLC(ixi^s, nw), wRC(ixi^s, nw)
148  double precision, intent(in) :: wLp(ixi^s, nw), wRp(ixi^s, nw)
149  double precision, intent(in) :: x(ixi^s, 1:^nd)
150  double precision, intent(inout) :: cmax(ixi^s)
151  double precision, intent(inout), optional :: cmin(ixi^s)
152  end subroutine sub_get_cbounds
153 
154  subroutine sub_get_flux(wC, w, x, ixI^L, ixO^L, idim, f)
156  integer, intent(in) :: ixI^L, ixO^L, idim
157  double precision, intent(in) :: wC(ixi^s, 1:nw)
158  double precision, intent(in) :: w(ixi^s, 1:nw)
159  double precision, intent(in) :: x(ixi^s, 1:^nd)
160  double precision, intent(out) :: f(ixi^s, nwflux)
161  end subroutine sub_get_flux
162 
163  subroutine sub_energy_synchro(qdt,ixI^L,ixO^L,wCT,w,x)
165  integer, intent(in) :: ixI^L,ixO^L
166  double precision, intent(in) :: qdt,wCT(ixi^s,1:nw),x(ixi^s,1:ndim)
167  double precision, intent(inout) :: w(ixi^s,1:nw)
168  end subroutine sub_energy_synchro
169 
170  subroutine sub_add_source_geom(qdt, ixI^L, ixO^L, wCT, w, x)
172  integer, intent(in) :: ixI^L, ixO^L
173  double precision, intent(in) :: qdt, x(ixi^s, 1:^nd)
174  double precision, intent(inout) :: wCT(ixi^s, 1:nw), w(ixi^s, 1:nw)
175  end subroutine sub_add_source_geom
176 
177  subroutine sub_add_source(qdt, ixI^L, ixO^L, wCT, w, x, &
178  qsourcesplit, active)
180  integer, intent(in) :: ixI^L, ixO^L
181  double precision, intent(in) :: qdt
182  double precision, intent(in) :: wCT(ixi^s, 1:nw), x(ixi^s, 1:ndim)
183  double precision, intent(inout) :: w(ixi^s, 1:nw)
184  logical, intent(in) :: qsourcesplit
185  logical, intent(inout) :: active !< Needs to be set to true when active
186  end subroutine sub_add_source
187 
188  !> Add global source terms on complete domain (potentially implicit)
189  subroutine sub_global_source(qdt, qt, active)
191  double precision, intent(in) :: qdt !< Current time step
192  double precision, intent(in) :: qt !< Current time
193  logical, intent(inout) :: active !< Output if the source is active
194  end subroutine sub_global_source
195 
196  subroutine sub_get_dt(w, ixI^L, ixO^L, dtnew, dx^D, x)
198  integer, intent(in) :: ixI^L, ixO^L
199  double precision, intent(in) :: dx^D, x(ixi^s, 1:^nd)
200  double precision, intent(in) :: w(ixi^s, 1:nw)
201  double precision, intent(inout) :: dtnew
202  end subroutine sub_get_dt
203 
204  subroutine sub_get_aux(clipping,w,x,ixI^L,ixO^L,subname)
206  integer, intent(in) :: ixI^L, ixO^L
207  double precision, intent(in) :: x(ixi^s,1:ndim)
208  double precision, intent(inout) :: w(ixi^s,nw)
209  logical, intent(in) :: clipping
210  character(len=*) :: subname
211  end subroutine sub_get_aux
212 
213  subroutine sub_check_w(primitive, ixI^L, ixO^L, w, w_flag)
215  logical, intent(in) :: primitive
216  integer, intent(in) :: ixI^L, ixO^L
217  double precision, intent(in) :: w(ixi^s,1:nw)
218  logical, intent(inout) :: w_flag(ixi^s,1:nw)
219  end subroutine sub_check_w
220 
221  subroutine sub_get_pthermal(w,x,ixI^L,ixO^L,pth)
223  integer, intent(in) :: ixI^L, ixO^L
224  double precision, intent(in) :: w(ixi^s,nw)
225  double precision, intent(in) :: x(ixi^s,1:ndim)
226  double precision, intent(out):: pth(ixi^s)
227  end subroutine sub_get_pthermal
228 
229  subroutine sub_write_info(file_handle)
230  integer, intent(in) :: file_handle
231  end subroutine sub_write_info
232 
233  subroutine sub_angmomfix(fC,x,wnew,ixI^L,ixO^L,idim)
235  integer, intent(in) :: ixI^L, ixO^L
236  double precision, intent(in) :: x(ixi^s,1:ndim)
237  double precision, intent(inout) :: fC(ixi^s,1:nwflux,1:ndim), wnew(ixi^s,1:nw)
238  integer, intent(in) :: idim
239  end subroutine sub_angmomfix
240 
241  subroutine sub_small_values(primitive, w, x, ixI^L, ixO^L, subname)
243  logical, intent(in) :: primitive
244  integer, intent(in) :: ixI^L,ixO^L
245  double precision, intent(inout) :: w(ixi^s,1:nw)
246  double precision, intent(in) :: x(ixi^s,1:ndim)
247  character(len=*), intent(in) :: subname
248  end subroutine sub_small_values
249 
250  subroutine sub_get_var(ixI^L, ixO^L, w, out)
252  integer, intent(in) :: ixI^L, ixO^L
253  double precision, intent(in) :: w(ixi^s, nw)
254  double precision, intent(out) :: out(ixo^s)
255  end subroutine sub_get_var
256 
257  subroutine sub_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s)
259  integer, intent(in) :: ixI^L, ixO^L
260  double precision, intent(in) :: qt, qdt
261  ! cell-center primitive variables
262  double precision, intent(in) :: wprim(ixi^s,1:nw)
263  ! velocity structure
264  type(state) :: sCT, s
265  double precision, intent(in) :: fC(ixi^s,1:nwflux,1:ndim)
266  double precision, intent(inout) :: fE(ixi^s,7-2*ndim:3)
267  end subroutine sub_update_faces
268 
269  subroutine sub_face_to_center(ixO^L,s)
271  integer, intent(in) :: ixO^L
272  type(state) :: s
273  end subroutine sub_face_to_center
274 
275  subroutine sub_evaluate_implicit(qtC,psa)
277  type(state), target :: psa(max_blocks)
278  double precision, intent(in) :: qtC
279  end subroutine sub_evaluate_implicit
280 
281  subroutine sub_implicit_update(dtfactor,qdt,qtC,psa,psb)
283  type(state), target :: psa(max_blocks)
284  type(state), target :: psb(max_blocks)
285  double precision, intent(in) :: qdt
286  double precision, intent(in) :: qtC
287  double precision, intent(in) :: dtfactor
288  end subroutine sub_implicit_update
289 
290  end interface
291 
292 contains
293 
294  subroutine phys_check()
296 
300 
301  if (physics_type == "") call mpistop("Error: no physics module loaded")
302 
303  call phys_hllc_check()
304  call phys_roe_check()
305  call phys_ppm_check()
306 
307  ! Checks whether the required physics methods have been defined
308  if (.not. associated(phys_check_params)) &
310 
311  if (.not. associated(phys_to_conserved)) &
312  call mpistop("Error: phys_to_conserved not defined")
313 
314  if (.not. associated(phys_to_primitive)) &
315  call mpistop("Error: phys_to_primitive not defined")
316 
317  if (.not. associated(phys_modify_wlr)) &
319 
320  if (.not. associated(phys_get_cmax)) &
321  call mpistop("Error: no phys_get_cmax not defined")
322 
323  if (.not. associated(phys_get_a2max)) &
325 
326  if (.not. associated(phys_get_cbounds)) &
327  call mpistop("Error: no phys_get_cbounds not defined")
328 
329  if (.not. associated(phys_get_flux)) &
330  call mpistop("Error: no phys_get_flux not defined")
331 
332  if (.not. associated(phys_get_dt)) &
333  call mpistop("Error: no phys_get_dt not defined")
334 
335  if (.not. associated(phys_add_source_geom)) &
337 
338  if (.not. associated(phys_add_source)) &
340 
341  if (.not. associated(phys_get_aux)) &
343 
344  if (.not. associated(phys_check_w)) &
346 
347  if (.not. associated(phys_get_pthermal)) &
349 
350  if (.not. associated(phys_boundary_adjust)) &
352 
353  if (.not. associated(phys_write_info)) &
355 
356  if (.not. associated(phys_angmomfix)) &
358 
359  if (.not. associated(phys_handle_small_values)) &
361 
362  if (.not. associated(phys_update_faces)) &
364 
365  if (.not. associated(phys_face_to_center)) &
367 
368  if (.not. associated(phys_evaluate_implicit)) &
370 
371  if (.not. associated(phys_implicit_update)) &
373 
374  end subroutine phys_check
375 
376  subroutine dummy_init_params
377  end subroutine dummy_init_params
378 
379  subroutine dummy_check_params
380  end subroutine dummy_check_params
381 
382  subroutine dummy_modify_wlr(ixI^L, ixO^L, qt, wLC, wRC, wLp, wRp, s, idir)
384  integer, intent(in) :: ixI^L, ixO^L, idir
385  double precision, intent(in) :: qt
386  double precision, intent(inout) :: wLC(ixi^s,1:nw), wRC(ixi^s,1:nw)
387  double precision, intent(inout) :: wLp(ixi^s,1:nw), wRp(ixi^s,1:nw)
388  type(state) :: s
389  end subroutine dummy_modify_wlr
390 
391  subroutine dummy_get_a2max(w, x, ixI^L, ixO^L, a2max)
393  integer, intent(in) :: ixI^L, ixO^L
394  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
395  double precision, intent(inout) :: a2max(ndim)
396  call mpistop("Error: entered dummy_get_a2max")
397  end subroutine dummy_get_a2max
398 
399  subroutine dummy_add_source_geom(qdt, ixI^L, ixO^L, wCT, w, x)
401  integer, intent(in) :: ixI^L, ixO^L
402  double precision, intent(in) :: qdt, x(ixi^s, 1:^nd)
403  double precision, intent(inout) :: wCT(ixi^s, 1:nw), w(ixi^s, 1:nw)
404  end subroutine dummy_add_source_geom
405 
406  subroutine dummy_add_source(qdt, ixI^L, ixO^L, wCT, w, x, &
407  qsourcesplit, active)
409  integer, intent(in) :: ixI^L, ixO^L
410  double precision, intent(in) :: qdt
411  double precision, intent(in) :: wCT(ixi^s, 1:nw), x(ixi^s, 1:ndim)
412  double precision, intent(inout) :: w(ixi^s, 1:nw)
413  logical, intent(in) :: qsourcesplit
414  logical, intent(inout) :: active
415  ! Don't have to set active, since it starts as .false.
416  end subroutine dummy_add_source
417 
418  subroutine dummy_get_aux(clipping,w,x,ixI^L,ixO^L,subname)
420  integer, intent(in) :: ixI^L, ixO^L
421  double precision, intent(in) :: x(ixi^s,1:ndim)
422  double precision, intent(inout) :: w(ixi^s,nw)
423  logical, intent(in) :: clipping
424  character(len=*) :: subname
425  end subroutine dummy_get_aux
426 
427  subroutine dummy_check_w(primitive, ixI^L, ixO^L, w, w_flag)
429  logical, intent(in) :: primitive
430  integer, intent(in) :: ixI^L, ixO^L
431  double precision, intent(in) :: w(ixi^s,1:nw)
432  logical, intent(inout) :: w_flag(ixi^s,1:nw)
433 
434  w_flag=.false. ! All okay
435  end subroutine dummy_check_w
436 
437  subroutine dummy_get_pthermal(w, x, ixI^L, ixO^L, pth)
439 
440  integer, intent(in) :: ixI^L, ixO^L
441  double precision, intent(in) :: w(ixi^s, nw)
442  double precision, intent(in) :: x(ixi^s, 1:ndim)
443  double precision, intent(out):: pth(ixi^s)
444 
445  call mpistop("No get_pthermal method specified")
446  end subroutine dummy_get_pthermal
447 
448  subroutine dummy_boundary_adjust
449  end subroutine dummy_boundary_adjust
450 
451  subroutine dummy_write_info(fh)
453  integer, intent(in) :: fh !< File handle
454  integer, dimension(MPI_STATUS_SIZE) :: st
455  integer :: er
456 
457  ! Number of physics parameters
458  integer, parameter :: n_par = 0
459 
460  call mpi_file_write(fh, n_par, 1, mpi_integer, st, er)
461  end subroutine dummy_write_info
462 
463  subroutine dummy_angmomfix(fC,x,wnew,ixI^L,ixO^L,idim)
465  double precision, intent(in) :: x(ixi^s,1:ndim)
466  double precision, intent(inout) :: fC(ixi^s,1:nwflux,1:ndim), wnew(ixi^s,1:nw)
467  integer, intent(in) :: ixI^L, ixO^L
468  integer, intent(in) :: idim
469  end subroutine dummy_angmomfix
470 
471  subroutine dummy_small_values(primitive, w, x, ixI^L, ixO^L, subname)
473  logical, intent(in) :: primitive
474  integer, intent(in) :: ixI^L,ixO^L
475  double precision, intent(inout) :: w(ixi^s,1:nw)
476  double precision, intent(in) :: x(ixi^s,1:ndim)
477  character(len=*), intent(in) :: subname
478  end subroutine dummy_small_values
479 
480  subroutine dummy_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s)
482  integer, intent(in) :: ixI^L, ixO^L
483  double precision, intent(in) :: qt, qdt
484  ! cell-center primitive variables
485  double precision, intent(in) :: wprim(ixi^s,1:nw)
486  type(state) :: sCT, s
487  double precision, intent(in) :: fC(ixi^s,1:nwflux,1:ndim)
488  double precision, intent(inout) :: fE(ixi^s,7-2*ndim:3)
489  end subroutine dummy_update_faces
490 
491  subroutine dummy_face_to_center(ixO^L,s)
493  integer, intent(in) :: ixO^L
494  type(state) :: s
495  end subroutine dummy_face_to_center
496 
497  subroutine dummy_evaluate_implicit(qtC,psa)
499  type(state), target :: psa(max_blocks)
500  double precision, intent(in) :: qtC
501  integer :: iigrid, igrid
502 
503  ! Just copy in nul state
504  !$OMP PARALLEL DO PRIVATE(igrid)
505  do iigrid=1,igridstail; igrid=igrids(iigrid);
506  psa(igrid)%w = 0.0d0*psa(igrid)%w
507  if(stagger_grid) psa(igrid)%ws = 0.0d0*psa(igrid)%ws
508  end do
509  !$OMP END PARALLEL DO
510 
511  end subroutine dummy_evaluate_implicit
512 
513  subroutine dummy_implicit_update(dtfactor,qdt,qtC,psa,psb)
515  type(state), target :: psa(max_blocks)
516  type(state), target :: psb(max_blocks)
517  double precision, intent(in) :: qdt
518  double precision, intent(in) :: qtC
519  double precision, intent(in) :: dtfactor
520  integer :: iigrid, igrid
521 
522  ! Just copy in psb state when using the scheme without implicit part
523  !$OMP PARALLEL DO PRIVATE(igrid)
524  do iigrid=1,igridstail; igrid=igrids(iigrid);
525  psa(igrid)%w = psb(igrid)%w
526  if(stagger_grid) psa(igrid)%ws = psb(igrid)%ws
527  end do
528  !$OMP END PARALLEL DO
529 
530  end subroutine dummy_implicit_update
531 
532 end module mod_physics
type(iw_methods), dimension(max_nw) phys_iw_methods
Special methods defined per variable.
Definition: mod_physics.t:57
subroutine dummy_modify_wlr(ixIL, ixOL, qt, wLC, wRC, wLp, wRp, s, idir)
Definition: mod_physics.t:383
This module contains definitions of global parameters and variables and some generic functions/subrou...
subroutine dummy_angmomfix(fC, x, wnew, ixIL, ixOL, idim)
Definition: mod_physics.t:464
subroutine phys_hllc_check
procedure(sub_global_source), pointer phys_global_source_after
Definition: mod_physics.t:76
procedure(sub_energy_synchro), pointer phys_energy_synchro
Definition: mod_physics.t:70
subroutine dummy_boundary_adjust
Definition: mod_physics.t:449
integer max_blocks
The maximum number of grid blocks in a processor.
character(len=name_len) physics_type
String describing the physics type of the simulation.
Definition: mod_physics.t:14
subroutine dummy_init_params
Definition: mod_physics.t:377
logical phys_total_energy
Solve total energy equation or not.
Definition: mod_physics.t:28
Add global source terms on complete domain (potentially implicit)
Definition: mod_physics.t:189
subroutine dummy_write_info(fh)
Definition: mod_physics.t:452
Type for special methods defined per variable.
Definition: mod_physics.t:50
integer ndir
Number of spatial dimensions (components) for vector variables.
procedure(sub_check_w), pointer phys_check_w
Definition: mod_physics.t:78
procedure(sub_face_to_center), pointer phys_face_to_center
Definition: mod_physics.t:85
procedure(sub_small_values), pointer phys_handle_small_values
Definition: mod_physics.t:83
procedure(sub_get_a2max), pointer phys_get_a2max
Definition: mod_physics.t:66
procedure(sub_boundary_adjust), pointer phys_boundary_adjust
Definition: mod_physics.t:80
logical phys_req_diagonal
Whether the physics routines require diagonal ghost cells, for example for computing a curl...
Definition: mod_physics.t:22
subroutine dummy_implicit_update(dtfactor, qdt, qtC, psa, psb)
Definition: mod_physics.t:514
procedure(sub_update_faces), pointer phys_update_faces
Definition: mod_physics.t:84
procedure(sub_convert), pointer phys_e_to_ei
Definition: mod_physics.t:63
procedure(sub_get_tcutoff), pointer phys_get_tcutoff
Definition: mod_physics.t:67
logical stagger_grid
True for using stagger grid.
procedure(sub_angmomfix), pointer phys_angmomfix
Definition: mod_physics.t:82
subroutine dummy_add_source_geom(qdt, ixIL, ixOL, wCT, w, x)
Definition: mod_physics.t:400
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:60
subroutine phys_check()
Definition: mod_physics.t:295
subroutine dummy_add_source(qdt, ixIL, ixOL, wCT, w, x, qsourcesplit, active)
Definition: mod_physics.t:408
integer, parameter flux_default
Indicates a normal flux.
Definition: mod_physics.t:41
procedure(sub_get_flux), pointer phys_get_flux
Definition: mod_physics.t:69
subroutine dummy_get_a2max(w, x, ixIL, ixOL, a2max)
Definition: mod_physics.t:392
procedure(sub_evaluate_implicit), pointer phys_evaluate_implicit
Definition: mod_physics.t:87
subroutine dummy_evaluate_implicit(qtC, psa)
Definition: mod_physics.t:498
integer, parameter flux_no_dissipation
Indicates dissipation should be omitted.
Definition: mod_physics.t:45
procedure(sub_global_source), pointer phys_global_source_before
Definition: mod_physics.t:75
procedure(sub_get_dt), pointer phys_get_dt
Definition: mod_physics.t:72
logical phys_internal_e
Solve internal enery instead of total energy.
Definition: mod_physics.t:31
subroutine phys_roe_check()
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition: mod_physics.t:73
subroutine dummy_check_w(primitive, ixIL, ixOL, w, w_flag)
Definition: mod_physics.t:428
subroutine phys_ppm_check
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:198
subroutine dummy_update_faces(ixIL, ixOL, qt, qdt, wprim, fC, fE, sCT, s)
Definition: mod_physics.t:481
procedure(sub_check_params), pointer phys_check_params
Definition: mod_physics.t:59
integer, parameter flux_special
Indicates the flux should be specially treated.
Definition: mod_physics.t:47
subroutine dummy_face_to_center(ixOL, s)
Definition: mod_physics.t:492
procedure(sub_implicit_update), pointer phys_implicit_update
Definition: mod_physics.t:86
subroutine dummy_check_params
Definition: mod_physics.t:380
integer, parameter ndim
Number of spatial dimensions for grid variables.
subroutine dummy_get_pthermal(w, x, ixIL, ixOL, pth)
Definition: mod_physics.t:438
integer phys_wider_stencil
To use wider stencils in flux calculations. A value of 1 will extend it by one cell in both direction...
Definition: mod_physics.t:18
logical phys_solve_eaux
Solve internal energy and total energy equations.
Definition: mod_physics.t:34
procedure(sub_get_v_idim), pointer phys_get_v_idim
Definition: mod_physics.t:71
procedure(sub_get_cmax), pointer phys_get_cmax
Definition: mod_physics.t:65
logical phys_energy
Solve energy equation or not.
Definition: mod_physics.t:25
procedure(sub_convert), pointer phys_ei_to_e
Definition: mod_physics.t:62
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_add_source), pointer phys_add_source
Definition: mod_physics.t:74
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:61
subroutine dummy_small_values(primitive, w, x, ixIL, ixOL, subname)
Definition: mod_physics.t:472
procedure(sub_get_aux), pointer phys_get_aux
Definition: mod_physics.t:77
integer, dimension(:, :), allocatable flux_type
Array per direction per variable, which can be used to specify that certain fluxes have to be treated...
Definition: mod_physics.t:38
procedure(sub_get_pthermal), pointer phys_get_pthermal
Definition: mod_physics.t:79
subroutine dummy_get_aux(clipping, w, x, ixIL, ixOL, subname)
Definition: mod_physics.t:419
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition: mod_physics.t:68
procedure(sub_modify_wlr), pointer phys_modify_wlr
Definition: mod_physics.t:64
procedure(sub_write_info), pointer phys_write_info
Definition: mod_physics.t:81
integer, parameter flux_tvdlf
Indicates the flux should be treated with tvdlf.
Definition: mod_physics.t:43