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