MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
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
8 
9  implicit none
10  public
11 
12 
13  double precision :: phys_gamma=5.d0/3.d0
14 
15  !> String describing the physics type of the simulation
16  character(len=name_len) :: physics_type = ""
17 
18  !> To use wider stencils in flux calculations. A value of 1 will extend it by
19  !> one cell in both directions, in any dimension
20  integer :: phys_wider_stencil = 0
21 
22  !> Whether the physics routines require diagonal ghost cells, for example for
23  !> computing a curl.
24  logical :: phys_req_diagonal = .true.
25 
26  !> Solve energy equation or not
27  logical :: phys_energy=.false.
28 
29  !> Solve total energy equation or not
30  logical :: phys_total_energy=.false.
31 
32  !> Solve internal energy instead of total energy
33  logical :: phys_internal_e=.false.
34 
35  !> Solve partially ionized one-fluid plasma
36  logical :: phys_partial_ionization=.false.
37 
38  !> Array per direction per variable, which can be used to specify that certain
39  !> fluxes have to be treated differently
40  integer, allocatable :: flux_type(:, :)
41 
42  !> Indicates a normal flux
43  integer, parameter :: flux_default = 0
44  !> Indicates the flux should be treated with tvdlf
45  integer, parameter :: flux_tvdlf = 1
46  !> Indicates dissipation should be omitted
47  integer, parameter :: flux_no_dissipation = 2
48  !> Indicates the flux should be specially treated
49  integer, parameter :: flux_special = 3
50  !> Indicates the flux should be treated with hll
51  integer, parameter :: flux_hll = 4
52 
53  procedure(sub_check_params), pointer :: phys_check_params => null()
54  procedure(sub_set_mg_bounds), pointer :: phys_set_mg_bounds => null()
55  procedure(sub_convert), pointer :: phys_to_conserved => null()
56  procedure(sub_convert), pointer :: phys_to_primitive => null()
57  procedure(sub_modify_wlr), pointer :: phys_modify_wlr => null()
58  procedure(sub_get_cmax), pointer :: phys_get_cmax => null()
59  procedure(sub_get_a2max), pointer :: phys_get_a2max => null()
60  procedure(sub_get_tcutoff), pointer :: phys_get_tcutoff => null()
61  procedure(sub_trac_after_setdt), pointer:: phys_trac_after_setdt => null()
62  procedure(sub_get_h_speed), pointer :: phys_get_h_speed => null()
63  procedure(sub_get_cbounds), pointer :: phys_get_cbounds => null()
64  procedure(sub_get_flux), pointer :: phys_get_flux => null()
65  procedure(sub_get_v), pointer :: phys_get_v => null()
66  procedure(sub_get_dt), pointer :: phys_get_dt => null()
67  procedure(sub_add_source_geom), pointer :: phys_add_source_geom => null()
68  procedure(sub_add_source), pointer :: phys_add_source => null()
69  procedure(sub_global_source), pointer :: phys_global_source_after => null()
70  procedure(sub_special_advance), pointer :: phys_special_advance => null()
71  procedure(sub_check_w), pointer :: phys_check_w => null()
72  procedure(sub_get_pthermal), pointer :: phys_get_pthermal => null()
73  procedure(sub_get_tgas), pointer :: phys_get_tgas => null()
74  procedure(sub_get_trad), pointer :: phys_get_trad => null()
75  procedure(sub_boundary_adjust), pointer :: phys_boundary_adjust => null()
76  procedure(sub_write_info), pointer :: phys_write_info => null()
77  procedure(sub_small_values), pointer :: phys_handle_small_values => null()
78  procedure(sub_get_ct_velocity), pointer :: phys_get_ct_velocity => null()
79  procedure(sub_update_faces), pointer :: phys_update_faces => null()
80  procedure(sub_face_to_center), pointer :: phys_face_to_center => null()
81  procedure(sub_implicit_update), pointer :: phys_implicit_update => null()
82  procedure(sub_evaluate_implicit),pointer:: phys_evaluate_implicit => null()
83  procedure(sub_clean_divb), pointer :: phys_clean_divb => null()
84  ! set the equilibrium variables
85  procedure(sub_set_equi_vars), pointer :: phys_set_equi_vars => null()
86  ! subroutine with no parameters which creates EUV images
87  procedure(sub_check_params), pointer :: phys_te_images => null()
88  ! to update temperature variable with partial ionization
89  procedure(sub_update_temperature), pointer :: phys_update_temperature => null()
90  procedure(sub_get_auxiliary), pointer :: phys_get_auxiliary => null()
91  procedure(sub_get_auxiliary_prim), pointer :: phys_get_auxiliary_prim => null()
92 
93  abstract interface
94 
95  subroutine sub_check_params
96  end subroutine sub_check_params
97 
98  subroutine sub_set_mg_bounds
100  use mod_usr_methods
101  end subroutine sub_set_mg_bounds
102 
103  subroutine sub_boundary_adjust(igrid,psb)
105  integer, intent(in) :: igrid
106  type(state), target :: psb(max_blocks)
107  end subroutine sub_boundary_adjust
108 
109  subroutine sub_convert(ixI^L, ixO^L, w, x)
111  integer, intent(in) :: ixI^L, ixO^L
112  double precision, intent(inout) :: w(ixI^S, nw)
113  double precision, intent(in) :: x(ixI^S, 1:^ND)
114  end subroutine sub_convert
115 
116  subroutine sub_modify_wlr(ixI^L, ixO^L, qt, wLC, wRC, wLp, wRp, s, idir)
118  integer, intent(in) :: ixI^L, ixO^L, idir
119  double precision, intent(in) :: qt
120  double precision, intent(inout) :: wLC(ixI^S,1:nw), wRC(ixI^S,1:nw)
121  double precision, intent(inout) :: wLp(ixI^S,1:nw), wRp(ixI^S,1:nw)
122  type(state) :: s
123  end subroutine sub_modify_wlr
124 
125  subroutine sub_get_cmax(w, x, ixI^L, ixO^L, idim, cmax)
127  integer, intent(in) :: ixI^L, ixO^L, idim
128  double precision, intent(in) :: w(ixI^S, nw), x(ixI^S, 1:^ND)
129  double precision, intent(inout) :: cmax(ixI^S)
130  end subroutine sub_get_cmax
131 
132  subroutine sub_get_a2max(w, x, ixI^L, ixO^L, a2max)
134  integer, intent(in) :: ixI^L, ixO^L
135  double precision, intent(in) :: w(ixI^S, nw), x(ixI^S, 1:^ND)
136  double precision, intent(inout) :: a2max(ndim)
137  end subroutine sub_get_a2max
138 
139  subroutine sub_get_tcutoff(ixI^L,ixO^L,w,x,tco_local,Tmax_local)
141  integer, intent(in) :: ixI^L, ixO^L
142  double precision, intent(inout) :: w(ixI^S, nw)
143  double precision, intent(in) :: x(ixI^S, 1:^ND)
144  double precision, intent(out) :: tco_local, Tmax_local
145  end subroutine sub_get_tcutoff
146 
147  subroutine sub_trac_after_setdt(trac_alfa,tco,T_peak,T_bott)
148  double precision, intent(in) :: trac_alfa,tco,T_peak,T_bott
149  end subroutine sub_trac_after_setdt
150 
151  subroutine sub_get_v(w,x,ixI^L,ixO^L,v)
153 
154  integer, intent(in) :: ixI^L, ixO^L
155  double precision, intent(in) :: w(ixI^S,nw), x(ixI^S,1:^ND)
156  double precision, intent(out) :: v(ixI^S,1:ndir)
157 
158  end subroutine sub_get_v
159 
160  subroutine sub_get_h_speed(wprim,x,ixI^L,ixO^L,idim,Hspeed)
162  integer, intent(in) :: ixI^L, ixO^L, idim
163  double precision, intent(in) :: wprim(ixI^S, nw)
164  double precision, intent(in) :: x(ixI^S,1:ndim)
165  double precision, intent(out) :: Hspeed(ixI^S,1:number_species)
166  end subroutine sub_get_h_speed
167 
168  subroutine sub_get_cbounds(wLC, wRC, wLp, wRp, x, ixI^L, ixO^L, idim, Hspeed, cmax, cmin)
170  use mod_variables
171  integer, intent(in) :: ixI^L, ixO^L, idim
172  double precision, intent(in) :: wLC(ixI^S, nw), wRC(ixI^S, nw)
173  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
174  double precision, intent(in) :: x(ixI^S, 1:^ND)
175  double precision, intent(inout) :: cmax(ixI^S,1:number_species)
176  double precision, intent(inout), optional :: cmin(ixI^S,1:number_species)
177  double precision, intent(in) :: Hspeed(ixI^S,1:number_species)
178  end subroutine sub_get_cbounds
179 
180  subroutine sub_get_flux(wC, w, x, ixI^L, ixO^L, idim, f)
182  integer, intent(in) :: ixI^L, ixO^L, idim
183  double precision, intent(in) :: wC(ixI^S, 1:nw)
184  double precision, intent(in) :: w(ixI^S, 1:nw)
185  double precision, intent(in) :: x(ixI^S, 1:^ND)
186  double precision, intent(out) :: f(ixI^S, nwflux)
187  end subroutine sub_get_flux
188 
189  subroutine sub_add_source_geom(qdt, dtfactor, ixI^L, ixO^L, wCT, w, x)
191  integer, intent(in) :: ixI^L, ixO^L
192  double precision, intent(in) :: qdt, dtfactor, x(ixI^S, 1:^ND)
193  double precision, intent(inout) :: wCT(ixI^S, 1:nw), w(ixI^S, 1:nw)
194  end subroutine sub_add_source_geom
195 
196  subroutine sub_add_source(qdt, dtfactor, ixI^L, ixO^L, wCT, wCTprim, w, x, &
197  qsourcesplit, active)
199  integer, intent(in) :: ixI^L, ixO^L
200  double precision, intent(in) :: qdt, dtfactor
201  double precision, intent(in) :: wCT(ixI^S, 1:nw), wCTprim(ixI^S,1:nw), x(ixI^S, 1:ndim)
202  double precision, intent(inout) :: w(ixI^S, 1:nw)
203  logical, intent(in) :: qsourcesplit
204  logical, intent(inout) :: active !< Needs to be set to true when active
205  end subroutine sub_add_source
206 
207  !> Add global source terms on complete domain (potentially implicit)
208  subroutine sub_global_source(qdt, qt, active)
210  double precision, intent(in) :: qdt !< Current time step
211  double precision, intent(in) :: qt !< Current time
212  logical, intent(inout) :: active !< Output if the source is active
213  end subroutine sub_global_source
214 
215  !> clean initial divb
216  subroutine sub_clean_divb(qdt, qt, active)
218  double precision, intent(in) :: qdt !< Current time step
219  double precision, intent(in) :: qt !< Current time
220  logical, intent(inout) :: active !< Output if the source is active
221  end subroutine sub_clean_divb
222 
223  !> set equilibrium variables other than b0 (e.g. p0 and rho0)
224  subroutine sub_set_equi_vars(igrid)
225  integer, intent(in) :: igrid
226  end subroutine sub_set_equi_vars
227 
228 
229  !> Add special advance in each advect step
230  subroutine sub_special_advance(qt, psa)
232  double precision, intent(in) :: qt !< Current time
233  type(state), target :: psa(max_blocks) !< Compute based on this state
234  end subroutine sub_special_advance
235 
236  subroutine sub_get_dt(w, ixI^L, ixO^L, dtnew, dx^D, x)
238  integer, intent(in) :: ixI^L, ixO^L
239  double precision, intent(in) :: dx^D, x(ixI^S, 1:^ND)
240  double precision, intent(in) :: w(ixI^S, 1:nw)
241  double precision, intent(inout) :: dtnew
242  end subroutine sub_get_dt
243 
244  subroutine sub_check_w(primitive, ixI^L, ixO^L, w, w_flag)
246  logical, intent(in) :: primitive
247  integer, intent(in) :: ixI^L, ixO^L
248  double precision, intent(in) :: w(ixI^S,1:nw)
249  logical, intent(inout) :: w_flag(ixI^S,1:nw)
250  end subroutine sub_check_w
251 
252  subroutine sub_get_pthermal(w,x,ixI^L,ixO^L,pth)
254  integer, intent(in) :: ixI^L, ixO^L
255  double precision, intent(in) :: w(ixI^S,nw)
256  double precision, intent(in) :: x(ixI^S,1:ndim)
257  double precision, intent(out):: pth(ixI^S)
258  end subroutine sub_get_pthermal
259 
260  subroutine sub_get_auxiliary(ixI^L,ixO^L,w,x)
262  integer, intent(in) :: ixI^L, ixO^L
263  double precision, intent(inout) :: w(ixI^S,nw)
264  double precision, intent(in) :: x(ixI^S,1:ndim)
265  end subroutine sub_get_auxiliary
266 
267  subroutine sub_get_auxiliary_prim(ixI^L,ixO^L,w)
269  integer, intent(in) :: ixI^L, ixO^L
270  double precision, intent(inout) :: w(ixI^S,nw)
271  end subroutine sub_get_auxiliary_prim
272 
273  subroutine sub_get_tgas(w,x,ixI^L,ixO^L,tgas)
275  integer, intent(in) :: ixI^L, ixO^L
276  double precision, intent(in) :: w(ixI^S,nw)
277  double precision, intent(in) :: x(ixI^S,1:ndim)
278  double precision, intent(out):: tgas(ixI^S)
279  end subroutine sub_get_tgas
280 
281  subroutine sub_get_trad(w,x,ixI^L,ixO^L,trad)
283  integer, intent(in) :: ixI^L, ixO^L
284  double precision, intent(in) :: w(ixI^S,nw)
285  double precision, intent(in) :: x(ixI^S,1:ndim)
286  double precision, intent(out):: trad(ixI^S)
287  end subroutine sub_get_trad
288 
289  subroutine sub_write_info(file_handle)
290  integer, intent(in) :: file_handle
291  end subroutine sub_write_info
292 
293  subroutine sub_small_values(primitive, w, x, ixI^L, ixO^L, subname)
295  logical, intent(in) :: primitive
296  integer, intent(in) :: ixI^L,ixO^L
297  double precision, intent(inout) :: w(ixI^S,1:nw)
298  double precision, intent(in) :: x(ixI^S,1:ndim)
299  character(len=*), intent(in) :: subname
300  end subroutine sub_small_values
301 
302  subroutine sub_get_var(ixI^L, ixO^L, w, out)
304  integer, intent(in) :: ixI^L, ixO^L
305  double precision, intent(in) :: w(ixI^S, nw)
306  double precision, intent(out) :: out(ixO^S)
307  end subroutine sub_get_var
308 
309  subroutine sub_get_ct_velocity(vcts,wLp,wRp,ixI^L,ixO^L,idim,cmax,cmin)
311 
312  integer, intent(in) :: ixI^L, ixO^L, idim
313  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
314  double precision, intent(in) :: cmax(ixI^S)
315  double precision, intent(in), optional :: cmin(ixI^S)
316  type(ct_velocity), intent(inout):: vcts
317  end subroutine sub_get_ct_velocity
318 
319  subroutine sub_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s,vcts)
321  integer, intent(in) :: ixI^L, ixO^L
322  double precision, intent(in) :: qt, qdt
323  ! cell-center primitive variables
324  double precision, intent(in) :: wprim(ixI^S,1:nw)
325  ! velocity structure
326  type(state) :: sCT, s
327  type(ct_velocity) :: vcts
328  double precision, intent(in) :: fC(ixI^S,1:nwflux,1:ndim)
329  double precision, intent(inout) :: fE(ixI^S,sdim:3)
330  end subroutine sub_update_faces
331 
332  subroutine sub_face_to_center(ixO^L,s)
334  integer, intent(in) :: ixO^L
335  type(state) :: s
336  end subroutine sub_face_to_center
337 
338  subroutine sub_evaluate_implicit(qtC,psa)
340  type(state), target :: psa(max_blocks)
341  double precision, intent(in) :: qtC
342  end subroutine sub_evaluate_implicit
343 
344  subroutine sub_implicit_update(dtfactor,qdt,qtC,psa,psb)
346  type(state), target :: psa(max_blocks)
347  type(state), target :: psb(max_blocks)
348  double precision, intent(in) :: qdt
349  double precision, intent(in) :: qtC
350  double precision, intent(in) :: dtfactor
351  end subroutine sub_implicit_update
352 
353  subroutine sub_update_temperature(ixI^L,ixO^L,wCT,w,x)
355  integer, intent(in) :: ixI^L, ixO^L
356  double precision, intent(in) :: wCT(ixI^S,nw),x(ixI^S,1:ndim)
357  double precision, intent(inout) :: w(ixI^S,nw)
358  end subroutine sub_update_temperature
359 
360  end interface
361 
362 contains
363 
364  subroutine phys_check()
365  use mod_global_parameters, only: nw, ndir
366 
369  use mod_comm_lib, only: mpistop
370 
371  if (physics_type == "") call mpistop("Error: no physics module loaded")
372 
373  call phys_hllc_check()
374  call phys_roe_check()
375 
376  ! Checks whether the required physics methods have been defined
377  if (.not. associated(phys_check_params)) &
379 
380  if (.not. associated(phys_to_conserved)) &
381  call mpistop("Error: phys_to_conserved not defined")
382 
383  if (.not. associated(phys_to_primitive)) &
384  call mpistop("Error: phys_to_primitive not defined")
385 
386  if (.not. associated(phys_modify_wlr)) &
388 
389  if (.not. associated(phys_get_cmax)) &
390  call mpistop("Error: no phys_get_cmax not defined")
391 
392  if (.not. associated(phys_get_a2max)) &
394 
395  if (.not. associated(phys_get_h_speed)) &
397 
398  if (.not. associated(phys_get_cbounds)) &
399  call mpistop("Error: no phys_get_cbounds not defined")
400 
401  if (.not. associated(phys_get_flux)) &
402  call mpistop("Error: no phys_get_flux not defined")
403 
404  if (.not. associated(phys_get_dt)) &
405  call mpistop("Error: no phys_get_dt not defined")
406 
407  if (.not. associated(phys_add_source_geom)) &
409 
410  if (.not. associated(phys_add_source)) &
412 
413  if (.not. associated(phys_check_w)) &
415 
416  if (.not. associated(phys_get_pthermal)) &
418  if (.not. associated(phys_get_auxiliary)) &
420  if (.not. associated(phys_get_auxiliary_prim)) &
422 
423  if (.not. associated(phys_boundary_adjust)) &
425 
426  if (.not. associated(phys_write_info)) &
428 
429  if (.not. associated(phys_handle_small_values)) &
431 
432  if (.not. associated(phys_get_ct_velocity)) &
434 
435  if (.not. associated(phys_update_faces)) &
437 
438  if (.not. associated(phys_face_to_center)) &
440 
441  if (.not. associated(phys_evaluate_implicit)) &
443 
444  if (.not. associated(phys_implicit_update)) &
446 
447  end subroutine phys_check
448 
449  subroutine dummy_init_params
450  end subroutine dummy_init_params
451 
453  end subroutine dummy_check_params
454 
455  subroutine dummy_modify_wlr(ixI^L, ixO^L, qt, wLC, wRC, wLp, wRp, s, idir)
457  integer, intent(in) :: ixI^L, ixO^L, idir
458  double precision, intent(in) :: qt
459  double precision, intent(inout) :: wLC(ixI^S,1:nw), wRC(ixI^S,1:nw)
460  double precision, intent(inout) :: wLp(ixI^S,1:nw), wRp(ixI^S,1:nw)
461  type(state) :: s
462  end subroutine dummy_modify_wlr
463 
464  subroutine dummy_get_h_speed(wprim,x,ixI^L,ixO^L,idim,Hspeed)
466  integer, intent(in) :: ixI^L, ixO^L, idim
467  double precision, intent(in) :: wprim(ixI^S, nw)
468  double precision, intent(in) :: x(ixI^S,1:ndim)
469  double precision, intent(out) :: Hspeed(ixI^S,1:number_species)
470  end subroutine dummy_get_h_speed
471 
472  subroutine dummy_get_a2max(w, x, ixI^L, ixO^L, a2max)
474  use mod_comm_lib, only: mpistop
475  integer, intent(in) :: ixI^L, ixO^L
476  double precision, intent(in) :: w(ixI^S, nw), x(ixI^S, 1:^ND)
477  double precision, intent(inout) :: a2max(ndim)
478  call mpistop("Error: entered dummy_get_a2max")
479  end subroutine dummy_get_a2max
480 
481  subroutine dummy_add_source_geom(qdt, dtfactor, ixI^L, ixO^L, wCT, w, x)
483  integer, intent(in) :: ixI^L, ixO^L
484  double precision, intent(in) :: qdt, dtfactor, x(ixI^S, 1:^ND)
485  double precision, intent(inout) :: wCT(ixI^S, 1:nw), w(ixI^S, 1:nw)
486  end subroutine dummy_add_source_geom
487 
488  subroutine dummy_add_source(qdt, dtfactor, ixI^L, ixO^L, wCT, wCTprim, w, x, &
489  qsourcesplit, active)
491  integer, intent(in) :: ixI^L, ixO^L
492  double precision, intent(in) :: qdt,dtfactor
493  double precision, intent(in) :: wCT(ixI^S, 1:nw), wCTprim(ixI^S,1:nw), x(ixI^S, 1:ndim)
494  double precision, intent(inout) :: w(ixI^S, 1:nw)
495  logical, intent(in) :: qsourcesplit
496  logical, intent(inout) :: active
497  ! Don't have to set active, since it starts as .false.
498  end subroutine dummy_add_source
499 
500  subroutine dummy_check_w(primitive, ixI^L, ixO^L, w, w_flag)
502  logical, intent(in) :: primitive
503  integer, intent(in) :: ixI^L, ixO^L
504  double precision, intent(in) :: w(ixI^S,1:nw)
505  logical, intent(inout) :: w_flag(ixI^S,1:nw)
506 
507  w_flag=.false. ! All okay
508  end subroutine dummy_check_w
509 
510  subroutine dummy_get_pthermal(w, x, ixI^L, ixO^L, pth)
512  use mod_comm_lib, only: mpistop
513 
514  integer, intent(in) :: ixI^L, ixO^L
515  double precision, intent(in) :: w(ixI^S, nw)
516  double precision, intent(in) :: x(ixI^S, 1:ndim)
517  double precision, intent(out):: pth(ixI^S)
518 
519  call mpistop("No get_pthermal method specified")
520  end subroutine dummy_get_pthermal
521 
522  subroutine dummy_get_auxiliary(ixI^L, ixO^L, w, x)
524  use mod_comm_lib, only: mpistop
525 
526  integer, intent(in) :: ixI^L, ixO^L
527  double precision, intent(inout) :: w(ixI^S, nw)
528  double precision, intent(in) :: x(ixI^S, 1:ndim)
529 
530  call mpistop("No get_auxiliary method specified")
531  end subroutine dummy_get_auxiliary
532 
533  subroutine dummy_get_auxiliary_prim(ixI^L, ixO^L, w)
535  use mod_comm_lib, only: mpistop
536 
537  integer, intent(in) :: ixI^L, ixO^L
538  double precision, intent(inout) :: w(ixI^S, nw)
539 
540  call mpistop("No get_auxiliary_prim method specified")
541  end subroutine dummy_get_auxiliary_prim
542 
543  subroutine dummy_boundary_adjust(igrid,psb)
545  integer, intent(in) :: igrid
546  type(state), target :: psb(max_blocks)
547  end subroutine dummy_boundary_adjust
548 
549  subroutine dummy_write_info(fh)
551  integer, intent(in) :: fh !< File handle
552  integer, dimension(MPI_STATUS_SIZE) :: st
553  integer :: er
554 
555  ! Number of physics parameters
556  integer, parameter :: n_par = 0
557 
558  call mpi_file_write(fh, n_par, 1, mpi_integer, st, er)
559  end subroutine dummy_write_info
560 
561  subroutine dummy_small_values(primitive, w, x, ixI^L, ixO^L, subname)
563  logical, intent(in) :: primitive
564  integer, intent(in) :: ixI^L,ixO^L
565  double precision, intent(inout) :: w(ixI^S,1:nw)
566  double precision, intent(in) :: x(ixI^S,1:ndim)
567  character(len=*), intent(in) :: subname
568  end subroutine dummy_small_values
569 
570  subroutine dummy_get_ct_velocity(vcts,wLp,wRp,ixI^L,ixO^L,idim,cmax,cmin)
572 
573  integer, intent(in) :: ixI^L, ixO^L, idim
574  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
575  double precision, intent(in) :: cmax(ixI^S)
576  double precision, intent(in), optional :: cmin(ixI^S)
577  type(ct_velocity), intent(inout):: vcts
578  end subroutine dummy_get_ct_velocity
579 
580  subroutine dummy_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s,vcts)
582  integer, intent(in) :: ixI^L, ixO^L
583  double precision, intent(in) :: qt, qdt
584  ! cell-center primitive variables
585  double precision, intent(in) :: wprim(ixI^S,1:nw)
586  type(state) :: sCT, s
587  type(ct_velocity) :: vcts
588  double precision, intent(in) :: fC(ixI^S,1:nwflux,1:ndim)
589  double precision, intent(inout) :: fE(ixI^S,sdim:3)
590  end subroutine dummy_update_faces
591 
592  subroutine dummy_face_to_center(ixO^L,s)
594  integer, intent(in) :: ixO^L
595  type(state) :: s
596  end subroutine dummy_face_to_center
597 
598  subroutine dummy_evaluate_implicit(qtC,psa)
600  type(state), target :: psa(max_blocks)
601  double precision, intent(in) :: qtC
602  integer :: iigrid, igrid
603 
604  ! Just copy in nul state
605  !$OMP PARALLEL DO PRIVATE(igrid)
606  do iigrid=1,igridstail; igrid=igrids(iigrid);
607  psa(igrid)%w = 0.0d0*psa(igrid)%w
608  if(stagger_grid) psa(igrid)%ws = 0.0d0*psa(igrid)%ws
609  end do
610  !$OMP END PARALLEL DO
611 
612  end subroutine dummy_evaluate_implicit
613 
614  subroutine dummy_implicit_update(dtfactor,qdt,qtC,psa,psb)
616  type(state), target :: psa(max_blocks)
617  type(state), target :: psb(max_blocks)
618  double precision, intent(in) :: qdt
619  double precision, intent(in) :: qtC
620  double precision, intent(in) :: dtfactor
621  integer :: iigrid, igrid
622 
623  ! Just copy in psb state when using the scheme without implicit part
624  !$OMP PARALLEL DO PRIVATE(igrid)
625  do iigrid=1,igridstail; igrid=igrids(iigrid);
626  psa(igrid)%w = psb(igrid)%w
627  if(stagger_grid) psa(igrid)%ws = psb(igrid)%ws
628  end do
629  !$OMP END PARALLEL DO
630 
631  end subroutine dummy_implicit_update
632 
633 end module mod_physics
clean initial divb
Definition: mod_physics.t:216
Add global source terms on complete domain (potentially implicit)
Definition: mod_physics.t:208
set equilibrium variables other than b0 (e.g. p0 and rho0)
Definition: mod_physics.t:224
Add special advance in each advect step.
Definition: mod_physics.t:230
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: mod_comm_lib.t:208
This module contains definitions of global parameters and variables and some generic functions/subrou...
logical stagger_grid
True for using stagger grid.
integer ndir
Number of spatial dimensions (components) for vector variables.
subroutine phys_hllc_check
subroutine phys_roe_check()
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_get_a2max), pointer phys_get_a2max
Definition: mod_physics.t:59
procedure(sub_get_ct_velocity), pointer phys_get_ct_velocity
Definition: mod_physics.t:78
subroutine dummy_update_faces(ixIL, ixOL, qt, qdt, wprim, fC, fE, sCT, s, vcts)
Definition: mod_physics.t:581
subroutine dummy_boundary_adjust(igrid, psb)
Definition: mod_physics.t:544
subroutine dummy_get_ct_velocity(vcts, wLp, wRp, ixIL, ixOL, idim, cmax, cmin)
Definition: mod_physics.t:571
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:56
procedure(sub_get_tgas), pointer phys_get_tgas
Definition: mod_physics.t:73
integer, parameter flux_tvdlf
Indicates the flux should be treated with tvdlf.
Definition: mod_physics.t:45
procedure(sub_small_values), pointer phys_handle_small_values
Definition: mod_physics.t:77
procedure(sub_write_info), pointer phys_write_info
Definition: mod_physics.t:76
procedure(sub_get_flux), pointer phys_get_flux
Definition: mod_physics.t:64
procedure(sub_evaluate_implicit), pointer phys_evaluate_implicit
Definition: mod_physics.t:82
logical phys_req_diagonal
Whether the physics routines require diagonal ghost cells, for example for computing a curl.
Definition: mod_physics.t:24
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition: mod_physics.t:63
procedure(sub_check_params), pointer phys_te_images
Definition: mod_physics.t:87
subroutine dummy_face_to_center(ixOL, s)
Definition: mod_physics.t:593
procedure(sub_get_dt), pointer phys_get_dt
Definition: mod_physics.t:66
logical phys_partial_ionization
Solve partially ionized one-fluid plasma.
Definition: mod_physics.t:36
logical phys_total_energy
Solve total energy equation or not.
Definition: mod_physics.t:30
integer, parameter flux_hll
Indicates the flux should be treated with hll.
Definition: mod_physics.t:51
procedure(sub_get_pthermal), pointer phys_get_pthermal
Definition: mod_physics.t:72
subroutine dummy_get_h_speed(wprim, x, ixIL, ixOL, idim, Hspeed)
Definition: mod_physics.t:465
procedure(sub_update_temperature), pointer phys_update_temperature
Definition: mod_physics.t:89
procedure(sub_get_tcutoff), pointer phys_get_tcutoff
Definition: mod_physics.t:60
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition: mod_physics.t:67
procedure(sub_check_params), pointer phys_check_params
Definition: mod_physics.t:53
integer, parameter flux_default
Indicates a normal flux.
Definition: mod_physics.t:43
procedure(sub_get_auxiliary), pointer phys_get_auxiliary
Definition: mod_physics.t:90
procedure(sub_get_auxiliary_prim), pointer phys_get_auxiliary_prim
Definition: mod_physics.t:91
integer, parameter flux_special
Indicates the flux should be specially treated.
Definition: mod_physics.t:49
procedure(sub_trac_after_setdt), pointer phys_trac_after_setdt
Definition: mod_physics.t:61
integer, parameter flux_no_dissipation
Indicates dissipation should be omitted.
Definition: mod_physics.t:47
procedure(sub_set_equi_vars), pointer phys_set_equi_vars
Definition: mod_physics.t:85
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:40
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:20
subroutine dummy_init_params
Definition: mod_physics.t:450
procedure(sub_update_faces), pointer phys_update_faces
Definition: mod_physics.t:79
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:55
character(len=name_len) physics_type
String describing the physics type of the simulation.
Definition: mod_physics.t:16
procedure(sub_check_w), pointer phys_check_w
Definition: mod_physics.t:71
procedure(sub_set_mg_bounds), pointer phys_set_mg_bounds
Definition: mod_physics.t:54
subroutine dummy_evaluate_implicit(qtC, psa)
Definition: mod_physics.t:599
procedure(sub_implicit_update), pointer phys_implicit_update
Definition: mod_physics.t:81
procedure(sub_get_trad), pointer phys_get_trad
Definition: mod_physics.t:74
subroutine dummy_add_source_geom(qdt, dtfactor, ixIL, ixOL, wCT, w, x)
Definition: mod_physics.t:482
subroutine dummy_check_params
Definition: mod_physics.t:453
subroutine dummy_get_auxiliary(ixIL, ixOL, w, x)
Definition: mod_physics.t:523
procedure(sub_clean_divb), pointer phys_clean_divb
Definition: mod_physics.t:83
subroutine dummy_implicit_update(dtfactor, qdt, qtC, psa, psb)
Definition: mod_physics.t:615
procedure(sub_boundary_adjust), pointer phys_boundary_adjust
Definition: mod_physics.t:75
double precision phys_gamma
Definition: mod_physics.t:13
procedure(sub_global_source), pointer phys_global_source_after
Definition: mod_physics.t:69
procedure(sub_add_source), pointer phys_add_source
Definition: mod_physics.t:68
procedure(sub_face_to_center), pointer phys_face_to_center
Definition: mod_physics.t:80
subroutine dummy_get_pthermal(w, x, ixIL, ixOL, pth)
Definition: mod_physics.t:511
procedure(sub_modify_wlr), pointer phys_modify_wlr
Definition: mod_physics.t:57
subroutine dummy_write_info(fh)
Definition: mod_physics.t:550
procedure(sub_get_h_speed), pointer phys_get_h_speed
Definition: mod_physics.t:62
logical phys_internal_e
Solve internal energy instead of total energy.
Definition: mod_physics.t:33
procedure(sub_get_cmax), pointer phys_get_cmax
Definition: mod_physics.t:58
procedure(sub_get_v), pointer phys_get_v
Definition: mod_physics.t:65
subroutine dummy_add_source(qdt, dtfactor, ixIL, ixOL, wCT, wCTprim, w, x, qsourcesplit, active)
Definition: mod_physics.t:490
subroutine dummy_get_auxiliary_prim(ixIL, ixOL, w)
Definition: mod_physics.t:534
subroutine dummy_get_a2max(w, x, ixIL, ixOL, a2max)
Definition: mod_physics.t:473
subroutine phys_check()
Definition: mod_physics.t:365
subroutine dummy_small_values(primitive, w, x, ixIL, ixOL, subname)
Definition: mod_physics.t:562
logical phys_energy
Solve energy equation or not.
Definition: mod_physics.t:27
subroutine dummy_check_w(primitive, ixIL, ixOL, w, w_flag)
Definition: mod_physics.t:501
procedure(sub_special_advance), pointer phys_special_advance
Definition: mod_physics.t:70
subroutine dummy_modify_wlr(ixIL, ixOL, qt, wLC, wRC, wLp, wRp, s, idir)
Definition: mod_physics.t:456
Module with all the methods that users can customize in AMRVAC.