MPI-AMRVAC  3.0
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 
11 
12  implicit none
13  public
14 
15 
16  double precision :: phys_gamma=5.0/3
17 
18  !> String describing the physics type of the simulation
19  character(len=name_len) :: physics_type = ""
20 
21  !> To use wider stencils in flux calculations. A value of 1 will extend it by
22  !> one cell in both directions, in any dimension
23  integer :: phys_wider_stencil = 0
24 
25  !> Whether the physics routines require diagonal ghost cells, for example for
26  !> computing a curl.
27  logical :: phys_req_diagonal = .true.
28 
29  !> Solve energy equation or not
30  logical :: phys_energy=.false.
31 
32  !> Solve total energy equation or not
33  logical :: phys_total_energy=.false.
34 
35  !> Solve internal enery instead of total energy
36  logical :: phys_internal_e=.false.
37 
38  !> Solve internal energy and total energy equations
39  logical :: phys_solve_eaux=.false.
40 
41  !> Array per direction per variable, which can be used to specify that certain
42  !> fluxes have to be treated differently
43  integer, allocatable :: flux_type(:, :)
44 
45  !> Indicates a normal flux
46  integer, parameter :: flux_default = 0
47  !> Indicates the flux should be treated with tvdlf
48  integer, parameter :: flux_tvdlf = 1
49  !> Indicates dissipation should be omitted
50  integer, parameter :: flux_no_dissipation = 2
51  !> Indicates the flux should be specially treated
52  integer, parameter :: flux_special = 3
53  !> Indicates the flux should be treated with hll
54  integer, parameter :: flux_hll = 4
55 
56  procedure(sub_check_params), pointer :: phys_check_params => null()
57  procedure(sub_set_mg_bounds), pointer :: phys_set_mg_bounds => null()
58  procedure(sub_convert), pointer :: phys_to_conserved => null()
59  procedure(sub_convert), pointer :: phys_to_primitive => null()
60  procedure(sub_modify_wlr), pointer :: phys_modify_wlr => null()
61  procedure(sub_get_cmax), pointer :: phys_get_cmax => null()
62  procedure(sub_get_a2max), pointer :: phys_get_a2max => null()
63  procedure(sub_get_tcutoff), pointer :: phys_get_tcutoff => null()
64  procedure(sub_get_h_speed), pointer :: phys_get_h_speed => null()
65  procedure(sub_get_cbounds), pointer :: phys_get_cbounds => null()
66  procedure(sub_get_flux), pointer :: phys_get_flux => null()
67  procedure(sub_energy_synchro), pointer :: phys_energy_synchro => null()
68  procedure(sub_get_v), pointer :: phys_get_v => null()
69  procedure(sub_get_dt), pointer :: phys_get_dt => null()
70  procedure(sub_add_source_geom), pointer :: phys_add_source_geom => null()
71  procedure(sub_add_source), pointer :: phys_add_source => null()
72  procedure(sub_global_source), pointer :: phys_global_source_after => null()
73  procedure(sub_special_advance), pointer :: phys_special_advance => null()
74  procedure(sub_check_w), pointer :: phys_check_w => null()
75  procedure(sub_get_pthermal), pointer :: phys_get_pthermal => null()
76  procedure(sub_get_tgas), pointer :: phys_get_tgas => null()
77  procedure(sub_get_trad), pointer :: phys_get_trad => null()
78  procedure(sub_boundary_adjust), pointer :: phys_boundary_adjust => null()
79  procedure(sub_write_info), pointer :: phys_write_info => null()
80  procedure(sub_angmomfix), pointer :: phys_angmomfix => null()
81  procedure(sub_small_values), pointer :: phys_handle_small_values => null()
82  procedure(sub_get_ct_velocity), pointer :: phys_get_ct_velocity => null()
83  procedure(sub_update_faces), pointer :: phys_update_faces => null()
84  procedure(sub_face_to_center), pointer :: phys_face_to_center => null()
85  procedure(sub_implicit_update), pointer :: phys_implicit_update => null()
86  procedure(sub_evaluate_implicit),pointer:: phys_evaluate_implicit => null()
87  procedure(sub_clean_divb), pointer :: phys_clean_divb => null()
88  ! set the equilibrium variables
89  procedure(sub_set_equi_vars), pointer :: phys_set_equi_vars => null()
90  ! subroutine with no parameters which creates EUV images
91  procedure(sub_check_params), pointer :: phys_te_images => 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_get_v(w,x,ixI^L,ixO^L,v)
149 
150  integer, intent(in) :: ixI^L, ixO^L
151  double precision, intent(in) :: w(ixI^S,nw), x(ixI^S,1:^ND)
152  double precision, intent(out) :: v(ixI^S,1:ndir)
153 
154  end subroutine sub_get_v
155 
156  subroutine sub_get_h_speed(wprim,x,ixI^L,ixO^L,idim,Hspeed)
158  integer, intent(in) :: ixI^L, ixO^L, idim
159  double precision, intent(in) :: wprim(ixI^S, nw)
160  double precision, intent(in) :: x(ixI^S,1:ndim)
161  double precision, intent(out) :: Hspeed(ixI^S,1:number_species)
162  end subroutine sub_get_h_speed
163 
164  subroutine sub_get_cbounds(wLC, wRC, wLp, wRp, x, ixI^L, ixO^L, idim, Hspeed, cmax, cmin)
166  use mod_variables
167  integer, intent(in) :: ixI^L, ixO^L, idim
168  double precision, intent(in) :: wLC(ixI^S, nw), wRC(ixI^S, nw)
169  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
170  double precision, intent(in) :: x(ixI^S, 1:^ND)
171  double precision, intent(inout) :: cmax(ixI^S,1:number_species)
172  double precision, intent(inout), optional :: cmin(ixI^S,1:number_species)
173  double precision, intent(in) :: Hspeed(ixI^S,1:number_species)
174  end subroutine sub_get_cbounds
175 
176  subroutine sub_get_flux(wC, w, x, ixI^L, ixO^L, idim, f)
178  integer, intent(in) :: ixI^L, ixO^L, idim
179  double precision, intent(in) :: wC(ixI^S, 1:nw)
180  double precision, intent(in) :: w(ixI^S, 1:nw)
181  double precision, intent(in) :: x(ixI^S, 1:^ND)
182  double precision, intent(out) :: f(ixI^S, nwflux)
183  end subroutine sub_get_flux
184 
185  subroutine sub_energy_synchro(ixI^L,ixO^L,w,x)
187  integer, intent(in) :: ixI^L,ixO^L
188  double precision, intent(in) :: x(ixI^S,1:ndim)
189  double precision, intent(inout) :: w(ixI^S,1:nw)
190  end subroutine sub_energy_synchro
191 
192  subroutine sub_add_source_geom(qdt, ixI^L, ixO^L, wCT, w, x)
194  integer, intent(in) :: ixI^L, ixO^L
195  double precision, intent(in) :: qdt, x(ixI^S, 1:^ND)
196  double precision, intent(inout) :: wCT(ixI^S, 1:nw), w(ixI^S, 1:nw)
197  end subroutine sub_add_source_geom
198 
199  subroutine sub_add_source(qdt, ixI^L, ixO^L, wCT, w, x, &
200  qsourcesplit, active, wCTprim)
202  integer, intent(in) :: ixI^L, ixO^L
203  double precision, intent(in) :: qdt
204  double precision, intent(in) :: wCT(ixI^S, 1:nw), x(ixI^S, 1:ndim)
205  double precision, intent(inout) :: w(ixI^S, 1:nw)
206  logical, intent(in) :: qsourcesplit
207  logical, intent(inout) :: active !< Needs to be set to true when active
208  double precision, intent(in), optional :: wCTprim(ixI^S,1:nw)
209  end subroutine sub_add_source
210 
211  !> Add global source terms on complete domain (potentially implicit)
212  subroutine sub_global_source(qdt, qt, active)
214  double precision, intent(in) :: qdt !< Current time step
215  double precision, intent(in) :: qt !< Current time
216  logical, intent(inout) :: active !< Output if the source is active
217  end subroutine sub_global_source
218 
219  !> clean initial divb
220  subroutine sub_clean_divb(qdt, qt, active)
222  double precision, intent(in) :: qdt !< Current time step
223  double precision, intent(in) :: qt !< Current time
224  logical, intent(inout) :: active !< Output if the source is active
225  end subroutine sub_clean_divb
226 
227  !> set equilibrium variables other than b0 (e.g. p0 and rho0)
228  subroutine sub_set_equi_vars(igrid)
229  integer, intent(in) :: igrid
230  end subroutine sub_set_equi_vars
231 
232 
233  !> Add special advance in each advect step
234  subroutine sub_special_advance(qt, psa)
236  double precision, intent(in) :: qt !< Current time
237  type(state), target :: psa(max_blocks) !< Compute based on this state
238  end subroutine sub_special_advance
239 
240  subroutine sub_get_dt(w, ixI^L, ixO^L, dtnew, dx^D, x)
242  integer, intent(in) :: ixI^L, ixO^L
243  double precision, intent(in) :: dx^D, x(ixI^S, 1:^ND)
244  double precision, intent(in) :: w(ixI^S, 1:nw)
245  double precision, intent(inout) :: dtnew
246  end subroutine sub_get_dt
247 
248  subroutine sub_check_w(primitive, ixI^L, ixO^L, w, w_flag)
250  logical, intent(in) :: primitive
251  integer, intent(in) :: ixI^L, ixO^L
252  double precision, intent(in) :: w(ixI^S,1:nw)
253  logical, intent(inout) :: w_flag(ixI^S,1:nw)
254  end subroutine sub_check_w
255 
256  subroutine sub_get_pthermal(w,x,ixI^L,ixO^L,pth)
258  integer, intent(in) :: ixI^L, ixO^L
259  double precision, intent(in) :: w(ixI^S,nw)
260  double precision, intent(in) :: x(ixI^S,1:ndim)
261  double precision, intent(out):: pth(ixI^S)
262  end subroutine sub_get_pthermal
263 
264  subroutine sub_get_tgas(w,x,ixI^L,ixO^L,tgas)
266  integer, intent(in) :: ixI^L, ixO^L
267  double precision, intent(in) :: w(ixI^S,nw)
268  double precision, intent(in) :: x(ixI^S,1:ndim)
269  double precision, intent(out):: tgas(ixI^S)
270  end subroutine sub_get_tgas
271 
272  subroutine sub_get_trad(w,x,ixI^L,ixO^L,trad)
274  integer, intent(in) :: ixI^L, ixO^L
275  double precision, intent(in) :: w(ixI^S,nw)
276  double precision, intent(in) :: x(ixI^S,1:ndim)
277  double precision, intent(out):: trad(ixI^S)
278  end subroutine sub_get_trad
279 
280  subroutine sub_write_info(file_handle)
281  integer, intent(in) :: file_handle
282  end subroutine sub_write_info
283 
284  subroutine sub_angmomfix(fC,x,wnew,ixI^L,ixO^L,idim)
286  integer, intent(in) :: ixI^L, ixO^L
287  double precision, intent(in) :: x(ixI^S,1:ndim)
288  double precision, intent(inout) :: fC(ixI^S,1:nwflux,1:ndim), wnew(ixI^S,1:nw)
289  integer, intent(in) :: idim
290  end subroutine sub_angmomfix
291 
292  subroutine sub_small_values(primitive, w, x, ixI^L, ixO^L, subname)
294  logical, intent(in) :: primitive
295  integer, intent(in) :: ixI^L,ixO^L
296  double precision, intent(inout) :: w(ixI^S,1:nw)
297  double precision, intent(in) :: x(ixI^S,1:ndim)
298  character(len=*), intent(in) :: subname
299  end subroutine sub_small_values
300 
301  subroutine sub_get_var(ixI^L, ixO^L, w, out)
303  integer, intent(in) :: ixI^L, ixO^L
304  double precision, intent(in) :: w(ixI^S, nw)
305  double precision, intent(out) :: out(ixO^S)
306  end subroutine sub_get_var
307 
308  subroutine sub_get_ct_velocity(vcts,wLp,wRp,ixI^L,ixO^L,idim,cmax,cmin)
310 
311  integer, intent(in) :: ixI^L, ixO^L, idim
312  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
313  double precision, intent(in) :: cmax(ixI^S)
314  double precision, intent(in), optional :: cmin(ixI^S)
315  type(ct_velocity), intent(inout):: vcts
316  end subroutine sub_get_ct_velocity
317 
318  subroutine sub_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s,vcts)
320  integer, intent(in) :: ixI^L, ixO^L
321  double precision, intent(in) :: qt, qdt
322  ! cell-center primitive variables
323  double precision, intent(in) :: wprim(ixI^S,1:nw)
324  ! velocity structure
325  type(state) :: sCT, s
326  type(ct_velocity) :: vcts
327  double precision, intent(in) :: fC(ixI^S,1:nwflux,1:ndim)
328  double precision, intent(inout) :: fE(ixI^S,7-2*ndim:3)
329  end subroutine sub_update_faces
330 
331  subroutine sub_face_to_center(ixO^L,s)
333  integer, intent(in) :: ixO^L
334  type(state) :: s
335  end subroutine sub_face_to_center
336 
337  subroutine sub_evaluate_implicit(qtC,psa)
339  type(state), target :: psa(max_blocks)
340  double precision, intent(in) :: qtC
341  end subroutine sub_evaluate_implicit
342 
343  subroutine sub_implicit_update(dtfactor,qdt,qtC,psa,psb)
345  type(state), target :: psa(max_blocks)
346  type(state), target :: psb(max_blocks)
347  double precision, intent(in) :: qdt
348  double precision, intent(in) :: qtC
349  double precision, intent(in) :: dtfactor
350  end subroutine sub_implicit_update
351 
352  end interface
353 
354 
355 contains
356 
357  subroutine phys_check()
358  use mod_global_parameters, only: nw, ndir
359 
363 
364  if (physics_type == "") call mpistop("Error: no physics module loaded")
365 
366  call phys_hllc_check()
367  call phys_roe_check()
368  call phys_ppm_check()
369 
370  ! Checks whether the required physics methods have been defined
371  if (.not. associated(phys_check_params)) &
373 
374  if (.not. associated(phys_to_conserved)) &
375  call mpistop("Error: phys_to_conserved not defined")
376 
377  if (.not. associated(phys_to_primitive)) &
378  call mpistop("Error: phys_to_primitive not defined")
379 
380  if (.not. associated(phys_modify_wlr)) &
382 
383  if (.not. associated(phys_get_cmax)) &
384  call mpistop("Error: no phys_get_cmax not defined")
385 
386  if (.not. associated(phys_get_a2max)) &
388 
389  if (.not. associated(phys_get_h_speed)) &
391 
392  if (.not. associated(phys_get_cbounds)) &
393  call mpistop("Error: no phys_get_cbounds not defined")
394 
395  if (.not. associated(phys_get_flux)) &
396  call mpistop("Error: no phys_get_flux not defined")
397 
398  if (.not. associated(phys_get_dt)) &
399  call mpistop("Error: no phys_get_dt not defined")
400 
401  if (.not. associated(phys_add_source_geom)) &
403 
404  if (.not. associated(phys_add_source)) &
406 
407  if (.not. associated(phys_check_w)) &
409 
410  if (.not. associated(phys_get_pthermal)) &
412 
413  if (.not. associated(phys_boundary_adjust)) &
415 
416  if (.not. associated(phys_write_info)) &
418 
419  if (.not. associated(phys_angmomfix)) &
421 
422  if (.not. associated(phys_handle_small_values)) &
424 
425  if (.not. associated(phys_get_ct_velocity)) &
427 
428  if (.not. associated(phys_update_faces)) &
430 
431  if (.not. associated(phys_face_to_center)) &
433 
434  if (.not. associated(phys_evaluate_implicit)) &
436 
437  if (.not. associated(phys_implicit_update)) &
439 
440  end subroutine phys_check
441 
442  subroutine dummy_init_params
443  end subroutine dummy_init_params
444 
446  end subroutine dummy_check_params
447 
448  subroutine dummy_modify_wlr(ixI^L, ixO^L, qt, wLC, wRC, wLp, wRp, s, idir)
450  integer, intent(in) :: ixI^L, ixO^L, idir
451  double precision, intent(in) :: qt
452  double precision, intent(inout) :: wLC(ixI^S,1:nw), wRC(ixI^S,1:nw)
453  double precision, intent(inout) :: wLp(ixI^S,1:nw), wRp(ixI^S,1:nw)
454  type(state) :: s
455  end subroutine dummy_modify_wlr
456 
457  subroutine dummy_get_h_speed(wprim,x,ixI^L,ixO^L,idim,Hspeed)
459  integer, intent(in) :: ixI^L, ixO^L, idim
460  double precision, intent(in) :: wprim(ixI^S, nw)
461  double precision, intent(in) :: x(ixI^S,1:ndim)
462  double precision, intent(out) :: Hspeed(ixI^S,1:number_species)
463  end subroutine dummy_get_h_speed
464 
465  subroutine dummy_get_a2max(w, x, ixI^L, ixO^L, a2max)
467  integer, intent(in) :: ixI^L, ixO^L
468  double precision, intent(in) :: w(ixI^S, nw), x(ixI^S, 1:^ND)
469  double precision, intent(inout) :: a2max(ndim)
470  call mpistop("Error: entered dummy_get_a2max")
471  end subroutine dummy_get_a2max
472 
473  subroutine dummy_add_source_geom(qdt, ixI^L, ixO^L, wCT, w, x)
475  integer, intent(in) :: ixI^L, ixO^L
476  double precision, intent(in) :: qdt, x(ixI^S, 1:^ND)
477  double precision, intent(inout) :: wCT(ixI^S, 1:nw), w(ixI^S, 1:nw)
478  end subroutine dummy_add_source_geom
479 
480  subroutine dummy_add_source(qdt, ixI^L, ixO^L, wCT, w, x, &
481  qsourcesplit, active, wCTprim)
483  integer, intent(in) :: ixI^L, ixO^L
484  double precision, intent(in) :: qdt
485  double precision, intent(in) :: wCT(ixI^S, 1:nw), x(ixI^S, 1:ndim)
486  double precision, intent(inout) :: w(ixI^S, 1:nw)
487  logical, intent(in) :: qsourcesplit
488  logical, intent(inout) :: active
489  double precision, intent(in), optional :: wCTprim(ixI^S,1:nw)
490  ! Don't have to set active, since it starts as .false.
491  end subroutine dummy_add_source
492 
493  subroutine dummy_check_w(primitive, ixI^L, ixO^L, w, w_flag)
495  logical, intent(in) :: primitive
496  integer, intent(in) :: ixI^L, ixO^L
497  double precision, intent(in) :: w(ixI^S,1:nw)
498  logical, intent(inout) :: w_flag(ixI^S,1:nw)
499 
500  w_flag=.false. ! All okay
501  end subroutine dummy_check_w
502 
503  subroutine dummy_get_pthermal(w, x, ixI^L, ixO^L, pth)
505 
506  integer, intent(in) :: ixI^L, ixO^L
507  double precision, intent(in) :: w(ixI^S, nw)
508  double precision, intent(in) :: x(ixI^S, 1:ndim)
509  double precision, intent(out):: pth(ixI^S)
510 
511  call mpistop("No get_pthermal method specified")
512  end subroutine dummy_get_pthermal
513 
514  subroutine dummy_boundary_adjust(igrid,psb)
516  integer, intent(in) :: igrid
517  type(state), target :: psb(max_blocks)
518  end subroutine dummy_boundary_adjust
519 
520  subroutine dummy_write_info(fh)
522  integer, intent(in) :: fh !< File handle
523  integer, dimension(MPI_STATUS_SIZE) :: st
524  integer :: er
525 
526  ! Number of physics parameters
527  integer, parameter :: n_par = 0
528 
529  call mpi_file_write(fh, n_par, 1, mpi_integer, st, er)
530  end subroutine dummy_write_info
531 
532  subroutine dummy_angmomfix(fC,x,wnew,ixI^L,ixO^L,idim)
534  double precision, intent(in) :: x(ixI^S,1:ndim)
535  double precision, intent(inout) :: fC(ixI^S,1:nwflux,1:ndim), wnew(ixI^S,1:nw)
536  integer, intent(in) :: ixI^L, ixO^L
537  integer, intent(in) :: idim
538  end subroutine dummy_angmomfix
539 
540  subroutine dummy_small_values(primitive, w, x, ixI^L, ixO^L, subname)
542  logical, intent(in) :: primitive
543  integer, intent(in) :: ixI^L,ixO^L
544  double precision, intent(inout) :: w(ixI^S,1:nw)
545  double precision, intent(in) :: x(ixI^S,1:ndim)
546  character(len=*), intent(in) :: subname
547  end subroutine dummy_small_values
548 
549  subroutine dummy_get_ct_velocity(vcts,wLp,wRp,ixI^L,ixO^L,idim,cmax,cmin)
551 
552  integer, intent(in) :: ixI^L, ixO^L, idim
553  double precision, intent(in) :: wLp(ixI^S, nw), wRp(ixI^S, nw)
554  double precision, intent(in) :: cmax(ixI^S)
555  double precision, intent(in), optional :: cmin(ixI^S)
556  type(ct_velocity), intent(inout):: vcts
557  end subroutine dummy_get_ct_velocity
558 
559  subroutine dummy_update_faces(ixI^L,ixO^L,qt,qdt,wprim,fC,fE,sCT,s,vcts)
561  integer, intent(in) :: ixI^L, ixO^L
562  double precision, intent(in) :: qt, qdt
563  ! cell-center primitive variables
564  double precision, intent(in) :: wprim(ixI^S,1:nw)
565  type(state) :: sCT, s
566  type(ct_velocity) :: vcts
567  double precision, intent(in) :: fC(ixI^S,1:nwflux,1:ndim)
568  double precision, intent(inout) :: fE(ixI^S,7-2*ndim:3)
569  end subroutine dummy_update_faces
570 
571  subroutine dummy_face_to_center(ixO^L,s)
573  integer, intent(in) :: ixO^L
574  type(state) :: s
575  end subroutine dummy_face_to_center
576 
577  subroutine dummy_evaluate_implicit(qtC,psa)
579  type(state), target :: psa(max_blocks)
580  double precision, intent(in) :: qtC
581  integer :: iigrid, igrid
582 
583  ! Just copy in nul state
584  !$OMP PARALLEL DO PRIVATE(igrid)
585  do iigrid=1,igridstail; igrid=igrids(iigrid);
586  psa(igrid)%w = 0.0d0*psa(igrid)%w
587  if(stagger_grid) psa(igrid)%ws = 0.0d0*psa(igrid)%ws
588  end do
589  !$OMP END PARALLEL DO
590 
591  end subroutine dummy_evaluate_implicit
592 
593  subroutine dummy_implicit_update(dtfactor,qdt,qtC,psa,psb)
595  type(state), target :: psa(max_blocks)
596  type(state), target :: psb(max_blocks)
597  double precision, intent(in) :: qdt
598  double precision, intent(in) :: qtC
599  double precision, intent(in) :: dtfactor
600  integer :: iigrid, igrid
601 
602  ! Just copy in psb state when using the scheme without implicit part
603  !$OMP PARALLEL DO PRIVATE(igrid)
604  do iigrid=1,igridstail; igrid=igrids(iigrid);
605  psa(igrid)%w = psb(igrid)%w
606  if(stagger_grid) psa(igrid)%ws = psb(igrid)%ws
607  end do
608  !$OMP END PARALLEL DO
609 
610  end subroutine dummy_implicit_update
611 
612 
613 
614 end module mod_physics
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:194
clean initial divb
Definition: mod_physics.t:220
Add global source terms on complete domain (potentially implicit)
Definition: mod_physics.t:212
set equilibrium variables other than b0 (e.g. p0 and rho0)
Definition: mod_physics.t:228
Add special advance in each advect step.
Definition: mod_physics.t:234
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_ppm_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:62
procedure(sub_get_ct_velocity), pointer phys_get_ct_velocity
Definition: mod_physics.t:82
subroutine dummy_update_faces(ixIL, ixOL, qt, qdt, wprim, fC, fE, sCT, s, vcts)
Definition: mod_physics.t:560
subroutine dummy_boundary_adjust(igrid, psb)
Definition: mod_physics.t:515
subroutine dummy_get_ct_velocity(vcts, wLp, wRp, ixIL, ixOL, idim, cmax, cmin)
Definition: mod_physics.t:550
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:59
procedure(sub_get_tgas), pointer phys_get_tgas
Definition: mod_physics.t:76
integer, parameter flux_tvdlf
Indicates the flux should be treated with tvdlf.
Definition: mod_physics.t:48
procedure(sub_small_values), pointer phys_handle_small_values
Definition: mod_physics.t:81
procedure(sub_write_info), pointer phys_write_info
Definition: mod_physics.t:79
procedure(sub_get_flux), pointer phys_get_flux
Definition: mod_physics.t:66
procedure(sub_evaluate_implicit), pointer phys_evaluate_implicit
Definition: mod_physics.t:86
subroutine dummy_angmomfix(fC, x, wnew, ixIL, ixOL, idim)
Definition: mod_physics.t:533
logical phys_req_diagonal
Whether the physics routines require diagonal ghost cells, for example for computing a curl.
Definition: mod_physics.t:27
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition: mod_physics.t:65
subroutine dummy_add_source_geom(qdt, ixIL, ixOL, wCT, w, x)
Definition: mod_physics.t:474
procedure(sub_check_params), pointer phys_te_images
Definition: mod_physics.t:91
logical phys_solve_eaux
Solve internal energy and total energy equations.
Definition: mod_physics.t:39
subroutine dummy_add_source(qdt, ixIL, ixOL, wCT, w, x, qsourcesplit, active, wCTprim)
Definition: mod_physics.t:482
subroutine dummy_face_to_center(ixOL, s)
Definition: mod_physics.t:572
procedure(sub_get_dt), pointer phys_get_dt
Definition: mod_physics.t:69
logical phys_total_energy
Solve total energy equation or not.
Definition: mod_physics.t:33
integer, parameter flux_hll
Indicates the flux should be treated with hll.
Definition: mod_physics.t:54
procedure(sub_get_pthermal), pointer phys_get_pthermal
Definition: mod_physics.t:75
subroutine dummy_get_h_speed(wprim, x, ixIL, ixOL, idim, Hspeed)
Definition: mod_physics.t:458
procedure(sub_get_tcutoff), pointer phys_get_tcutoff
Definition: mod_physics.t:63
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition: mod_physics.t:70
procedure(sub_check_params), pointer phys_check_params
Definition: mod_physics.t:56
integer, parameter flux_default
Indicates a normal flux.
Definition: mod_physics.t:46
integer, parameter flux_special
Indicates the flux should be specially treated.
Definition: mod_physics.t:52
integer, parameter flux_no_dissipation
Indicates dissipation should be omitted.
Definition: mod_physics.t:50
procedure(sub_set_equi_vars), pointer phys_set_equi_vars
Definition: mod_physics.t:89
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:43
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:23
subroutine dummy_init_params
Definition: mod_physics.t:443
procedure(sub_update_faces), pointer phys_update_faces
Definition: mod_physics.t:83
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:58
character(len=name_len) physics_type
String describing the physics type of the simulation.
Definition: mod_physics.t:19
procedure(sub_check_w), pointer phys_check_w
Definition: mod_physics.t:74
procedure(sub_set_mg_bounds), pointer phys_set_mg_bounds
Definition: mod_physics.t:57
subroutine dummy_evaluate_implicit(qtC, psa)
Definition: mod_physics.t:578
procedure(sub_energy_synchro), pointer phys_energy_synchro
Definition: mod_physics.t:67
procedure(sub_implicit_update), pointer phys_implicit_update
Definition: mod_physics.t:85
procedure(sub_get_trad), pointer phys_get_trad
Definition: mod_physics.t:77
subroutine dummy_check_params
Definition: mod_physics.t:446
procedure(sub_clean_divb), pointer phys_clean_divb
Definition: mod_physics.t:87
subroutine dummy_implicit_update(dtfactor, qdt, qtC, psa, psb)
Definition: mod_physics.t:594
procedure(sub_boundary_adjust), pointer phys_boundary_adjust
Definition: mod_physics.t:78
double precision phys_gamma
Definition: mod_physics.t:16
procedure(sub_global_source), pointer phys_global_source_after
Definition: mod_physics.t:72
procedure(sub_add_source), pointer phys_add_source
Definition: mod_physics.t:71
procedure(sub_face_to_center), pointer phys_face_to_center
Definition: mod_physics.t:84
subroutine dummy_get_pthermal(w, x, ixIL, ixOL, pth)
Definition: mod_physics.t:504
procedure(sub_modify_wlr), pointer phys_modify_wlr
Definition: mod_physics.t:60
subroutine dummy_write_info(fh)
Definition: mod_physics.t:521
procedure(sub_get_h_speed), pointer phys_get_h_speed
Definition: mod_physics.t:64
logical phys_internal_e
Solve internal enery instead of total energy.
Definition: mod_physics.t:36
procedure(sub_angmomfix), pointer phys_angmomfix
Definition: mod_physics.t:80
procedure(sub_get_cmax), pointer phys_get_cmax
Definition: mod_physics.t:61
procedure(sub_get_v), pointer phys_get_v
Definition: mod_physics.t:68
subroutine dummy_get_a2max(w, x, ixIL, ixOL, a2max)
Definition: mod_physics.t:466
subroutine phys_check()
Definition: mod_physics.t:358
subroutine dummy_small_values(primitive, w, x, ixIL, ixOL, subname)
Definition: mod_physics.t:541
logical phys_energy
Solve energy equation or not.
Definition: mod_physics.t:30
subroutine dummy_check_w(primitive, ixIL, ixOL, w, w_flag)
Definition: mod_physics.t:494
procedure(sub_special_advance), pointer phys_special_advance
Definition: mod_physics.t:73
subroutine dummy_modify_wlr(ixIL, ixOL, qt, wLC, wRC, wLp, wRp, s, idir)
Definition: mod_physics.t:449
Module with all the methods that users can customize in AMRVAC.