MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_usr_methods.t
Go to the documentation of this file.
1 !> Module with all the methods that users can customize in AMRVAC
2 !>
3 !> Each procedure pointer can be initialized in a user's mod_usr.t
5 
6  implicit none
7  public
8 
9  !> Initialize the user's settings (after initializing amrvac)
10  procedure(p_no_args), pointer :: usr_set_parameters => null()
11  !> Initialize earch grid block data
12  procedure(init_one_grid), pointer :: usr_init_one_grid => null()
13 
14  ! Boundary condition related
15  procedure(special_bc), pointer :: usr_special_bc => null()
16  procedure(internal_bc), pointer :: usr_internal_bc => null()
17 
18  ! Output related
19  procedure(p_no_args), pointer :: usr_print_log => null()
20  procedure(p_no_args), pointer :: usr_write_analysis => null()
21  procedure(transform_w), pointer :: usr_transform_w => null()
22  procedure(aux_output), pointer :: usr_aux_output => null()
23  procedure(add_aux_names), pointer :: usr_add_aux_names => null()
24  procedure(sub_modify_io), pointer :: usr_modify_output => null()
25  procedure(special_convert), pointer :: usr_special_convert => null()
26 
27  ! Called at the beginning of every time step (after determining dt)
28  procedure(process_grid), pointer :: usr_process_grid => null()
29  procedure(process_global), pointer :: usr_process_global => null()
30 
31  ! Called every time step just after advance (with w^(n+1), it^n, global_time^n)
32  procedure(process_adv_grid), pointer :: usr_process_adv_grid => null()
33  procedure(process_adv_global), pointer :: usr_process_adv_global => null()
34 
35  ! Called after initial condition before the start of the simulation
36  procedure(p_no_args), pointer :: usr_improve_initial_condition => null()
37 
38  ! Called before the start of the simulation
39  procedure(p_no_args), pointer :: usr_before_main_loop => null()
40 
41  ! Source terms
42  procedure(source), pointer :: usr_source => null()
43  procedure(get_dt), pointer :: usr_get_dt => null()
44  procedure(phys_gravity), pointer :: usr_gravity => null()
45 
46  ! Usr defined dust drag force
47  procedure(phys_dust_get_dt), pointer :: usr_dust_get_dt => null()
48  procedure(phys_dust_get_3d_dragforce), pointer :: usr_get_3d_dragforce => null()
49 
50  ! Usr defined space varying viscosity
51  procedure(phys_visco), pointer :: usr_setvisco => null()
52 
53  ! Usr defined thermal pressure for hydro & energy=.False.
54  procedure(hd_pthermal), pointer :: usr_set_pthermal => null()
55 
56  ! Refinement related procedures
57  procedure(refine_grid), pointer :: usr_refine_grid => null()
58  procedure(var_for_errest), pointer :: usr_var_for_errest => null()
59  procedure(a_refine_threshold), pointer :: usr_refine_threshold => null()
60  procedure(flag_grid), pointer :: usr_flag_grid => null()
61 
62  ! Set time-independent magnetic field for B0 splitting
63  procedure(set_b0), pointer :: usr_set_b0 => null()
64  ! Set time-independent current density for B0 splitting
65  procedure(set_j0), pointer :: usr_set_j0 => null()
66  procedure(special_resistivity), pointer :: usr_special_resistivity => null()
67 
68  ! Particles module related
69  procedure(update_payload), pointer :: usr_update_payload => null()
70  procedure(create_particles), pointer :: usr_create_particles => null()
71  procedure(particle_fields), pointer :: usr_particle_fields => null()
72  procedure(particle_analytic), pointer :: usr_particle_analytic => null()
73 
74  ! Called after the mesh has been adjuste
75  procedure(after_refine), pointer :: usr_after_refine => null()
76 
77  ! allow user to explicitly set flux at cell interfaces for finite volume scheme
78  procedure(set_flux), pointer :: usr_set_flux => null()
79 
80  ! initialize vector potential on cell edges for magnetic field
81  procedure(init_vector_potential), pointer :: usr_init_vector_potential => null()
82 
83  ! allow user to change inductive electric field, especially for boundary driven applications
84  procedure(set_electric_field), pointer :: usr_set_electric_field => null()
85 
86  abstract interface
87 
88  subroutine p_no_args()
89  end subroutine p_no_args
90 
91  !> Initialize one grid
92  subroutine init_one_grid(ixI^L,ixO^L,w,x)
94  integer, intent(in) :: ixI^L, ixO^L
95  double precision, intent(in) :: x(ixi^s,1:ndim)
96  double precision, intent(inout) :: w(ixi^s,1:nw)
97  end subroutine init_one_grid
98 
99  !> special boundary types, users must assign conservative
100  !> variables in boundaries
101  subroutine special_bc(qt,ixI^L,ixO^L,iB,w,x)
103  !> Shape of input arrays
104  integer, intent(in) :: ixI^L
105  !> Region where boundary values have to be set
106  integer, intent(in) :: ixO^L
107  !> Integer indicating direction of boundary
108  integer, intent(in) :: iB
109  double precision, intent(in) :: qt, x(ixi^s,1:ndim)
110  double precision, intent(inout) :: w(ixi^s,1:nw)
111  end subroutine special_bc
112 
113  !> internal boundary, user defined
114  !> This subroutine can be used to artificially overwrite ALL conservative
115  !> variables in a user-selected region of the mesh, and thereby act as
116  !> an internal boundary region. It is called just before external (ghost cell)
117  !> boundary regions will be set by the BC selection. Here, you could e.g.
118  !> want to introduce an extra variable (nwextra, to be distinguished from nwaux)
119  !> which can be used to identify the internal boundary region location.
120  !> Its effect should always be local as it acts on the mesh.
121  subroutine internal_bc(level,qt,ixI^L,ixO^L,w,x)
123  integer, intent(in) :: ixI^L,ixO^L,level
124  double precision, intent(in) :: qt
125  double precision, intent(inout) :: w(ixi^s,1:nw)
126  double precision, intent(in) :: x(ixi^s,1:ndim)
127  end subroutine internal_bc
128 
129  !> this subroutine is ONLY to be used for computing auxiliary variables
130  !> which happen to be non-local (like div v), and are in no way used for
131  !> flux computations. As auxiliaries, they are also not advanced
132  subroutine process_grid(igrid,level,ixI^L,ixO^L,qt,w,x)
134  integer, intent(in) :: igrid,level,ixI^L,ixO^L
135  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
136  double precision, intent(inout) :: w(ixi^s,1:nw)
137  end subroutine process_grid
138 
139  !> If defined, this routine is called before writing output, and it can
140  !> set/modify the variables in the w array.
141  subroutine sub_modify_io(ixI^L,ixO^L,qt,w,x)
143  integer, intent(in) :: ixI^L,ixO^L
144  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
145  double precision, intent(inout) :: w(ixi^s,1:nw)
146  end subroutine sub_modify_io
147 
148  !> This subroutine is called at the beginning of each time step
149  !> by each processor. No communication is specified, so the user
150  !> has to implement MPI routines if information has to be shared
151  subroutine process_global(iit,qt)
153  integer, intent(in) :: iit
154  double precision, intent(in) :: qt
155  end subroutine process_global
156 
157  !> for processing after the advance (PIC-MHD, e.g.)
158  subroutine process_adv_grid(igrid,level,ixI^L,ixO^L,qt,w,x)
160  integer, intent(in) :: igrid,level,ixI^L,ixO^L
161  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
162  double precision, intent(inout) :: w(ixi^s,1:nw)
163  end subroutine process_adv_grid
164 
165  !> for processing after the advance (PIC-MHD, e.g.)
166  subroutine process_adv_global(iit,qt)
168  integer, intent(in) :: iit
169  double precision, intent(in) :: qt
170  end subroutine process_adv_global
171 
172  !> this subroutine can be used in convert, to add auxiliary variables to the
173  !> converted output file, for further analysis using tecplot, paraview, ....
174  !> these auxiliary values need to be stored in the nw+1:nw+nwauxio slots
175  !
176  !> the array normconv can be filled in the (nw+1:nw+nwauxio) range with
177  !> corresponding normalization values (default value 1)
178  subroutine aux_output(ixI^L,ixO^L,w,x,normconv)
180  integer, intent(in) :: ixI^L,ixO^L
181  double precision, intent(in) :: x(ixi^s,1:ndim)
182  double precision :: w(ixi^s,nw+nwauxio)
183  double precision :: normconv(0:nw+nwauxio)
184  end subroutine aux_output
185 
186  !> Add names for the auxiliary variables
187  subroutine add_aux_names(varnames)
189  character(len=*) :: varnames
190  end subroutine add_aux_names
191 
192  !> Calculate w(iw)=w(iw)+qdt*SOURCE[wCT,qtC,x] within ixO for all indices
193  !> iw=iwmin...iwmax. wCT is at time qCT
194  subroutine source(qdt,ixI^L,ixO^L,iw^LIM,qtC,wCT,qt,w,x)
196  integer, intent(in) :: ixI^L, ixO^L, iw^LIM
197  double precision, intent(in) :: qdt, qtC, qt
198  double precision, intent(in) :: wCT(ixi^s,1:nw), x(ixi^s,1:ndim)
199  double precision, intent(inout) :: w(ixi^s,1:nw)
200  end subroutine source
201 
202  !> Limit "dt" further if necessary, e.g. due to the special source terms.
203  !> The getdt_courant (CFL condition) and the getdt subroutine in the AMRVACPHYS
204  !> module have already been called.
205  subroutine get_dt(w,ixI^L,ixO^L,dtnew,dx^D,x)
207  integer, intent(in) :: ixI^L, ixO^L
208  double precision, intent(in) :: dx^D, x(ixi^s,1:ndim)
209  double precision, intent(in) :: w(ixi^s,1:nw)
210  double precision, intent(inout) :: dtnew
211  end subroutine get_dt
212 
213  !> Calculate gravitational acceleration in each dimension
214  subroutine phys_gravity(ixI^L,ixO^L,wCT,x,gravity_field)
216  integer, intent(in) :: ixI^L, ixO^L
217  double precision, intent(in) :: x(ixi^s,1:ndim)
218  double precision, intent(in) :: wCT(ixi^s,1:nw)
219  double precision, intent(out) :: gravity_field(ixi^s,ndim)
220  end subroutine phys_gravity
221 
222  !> Calculate the 3d drag force of gas onto dust
223  subroutine phys_dust_get_3d_dragforce(ixI^L, ixO^L, w, x, fdrag, ptherm, vgas,dust_n_species)
225  integer, intent(in) :: ixI^L, ixO^L, dust_n_species
226  double precision, intent(in) :: x(ixi^s, 1:ndim)
227  double precision, intent(in) :: w(ixi^s, 1:nw)
228  double precision, intent(out) :: &
229  fdrag(ixi^s, 1:ndir, 1:dust_n_species)
230  double precision, intent(in) :: ptherm(ixi^s), vgas(ixi^s, ndir)
231  end subroutine phys_dust_get_3d_dragforce
232 
233  !> Calculate the time step associated with the usr drag force
234  subroutine phys_dust_get_dt(w, ixI^L, ixO^L, dtdust, dx^D, x, dust_n_species)
236  integer, intent(in) :: ixI^L, ixO^L, dust_n_species
237  double precision, intent(in) :: dx^D, x(ixi^s,1:ndim)
238  double precision, intent(in) :: w(ixi^s,1:nw)
239  double precision, intent(inout) :: dtdust(1:dust_n_species)
240  end subroutine phys_dust_get_dt
241 
242  !>Calculation anormal viscosity depending on space
243  subroutine phys_visco(ixI^L,ixO^L,x,w,mu)
245  integer, intent(in) :: ixI^L, ixO^L
246  double precision, intent(in) :: x(ixi^s,1:ndim)
247  double precision, intent(in) :: w(ixi^s,1:nw)
248  double precision, intent(out) :: mu(ixi^s)
249  end subroutine phys_visco
250 
251  !>Calculation anormal pressure for hd & energy=.False.
252  subroutine hd_pthermal(w,x,ixI^L,ixO^L,pth)
254  integer, intent(in) :: ixI^L, ixO^L
255  double precision, intent(in) :: x(ixi^s,1:ndim)
256  double precision, intent(in) :: w(ixi^s,1:nw)
257  double precision, intent(out) :: pth(ixi^s)
258  end subroutine hd_pthermal
259 
260  !> Set the "eta" array for resistive MHD based on w or the
261  !> "current" variable which has components between idirmin and 3.
262  subroutine special_resistivity(w,ixI^L,ixO^L,idirmin,x,current,eta)
264  integer, intent(in) :: ixI^L, ixO^L, idirmin
265  double precision, intent(in) :: w(ixi^s,nw), x(ixi^s,1:ndim)
266  double precision :: current(ixi^s,7-2*ndir:3), eta(ixi^s)
267  end subroutine special_resistivity
268 
269  !> Enforce additional refinement or coarsening
270  !> One can use the coordinate info in x and/or time qt=t_n and w(t_n) values w.
271  !> you must set consistent values for integers refine/coarsen:
272  !> refine = -1 enforce to not refine
273  !> refine = 0 doesn't enforce anything
274  !> refine = 1 enforce refinement
275  !> coarsen = -1 enforce to not coarsen
276  !> coarsen = 0 doesn't enforce anything
277  !> coarsen = 1 enforce coarsen
278  !> e.g. refine for negative first coordinate x < 0 as
279  !> if (any(x(ix^S,1) < zero)) refine=1
280  subroutine refine_grid(igrid,level,ixI^L,ixO^L,qt,w,x,refine,coarsen)
282  integer, intent(in) :: igrid, level, ixI^L, ixO^L
283  double precision, intent(in) :: qt, w(ixi^s,1:nw), x(ixi^s,1:ndim)
284  integer, intent(inout) :: refine, coarsen
285  end subroutine refine_grid
286 
287  !> this is the place to compute a local auxiliary variable to be used
288  !> as refinement criterion for the Lohner error estimator only
289  !> -->it is then requiring and iflag>nw
290  !> note that ixO=ixI=ixG, hence the term local (gradients need special attention!)
291  subroutine var_for_errest(ixI^L,ixO^L,iflag,w,x,var)
293  integer, intent(in) :: ixI^L,ixO^L,iflag
294  double precision, intent(in) :: w(ixi^s,1:nw), x(ixi^s,1:ndim)
295  double precision, intent(out) :: var(ixi^s)
296  end subroutine var_for_errest
297 
298  !> Here one can add a steady (time-independent) potential background field
299  subroutine set_b0(ixI^L,ixO^L,x,wB0)
301  integer, intent(in) :: ixI^L,ixO^L
302  double precision, intent(in) :: x(ixi^s,1:ndim)
303  double precision, intent(inout) :: wB0(ixi^s,1:ndir)
304  end subroutine set_b0
305 
306  !> Here one can add a time-independent background current density
307  subroutine set_j0(ixI^L,ixO^L,x,wJ0)
309  integer, intent(in) :: ixI^L,ixO^L
310  double precision, intent(in) :: x(ixi^s,1:ndim)
311  double precision, intent(inout) :: wJ0(ixi^s,7-2*ndir:ndir)
312  end subroutine set_j0
313 
314  !> regenerate w and eqpartf arrays to output into *tf.dat
315  subroutine transform_w(ixI^L,ixO^L,nw_in,w_in,x,w_out)
317  integer, intent(in) :: ixI^L, ixO^L, nw_in
318  double precision, intent(in) :: w_in(ixi^s,1:nw_in)
319  double precision, intent(in) :: x(ixi^s, 1:ndim)
320  double precision, intent(out) :: w_out(ixi^s,1:nw)
321  end subroutine transform_w
322 
323  !> use different threshold in special regions for AMR to
324  !> reduce/increase resolution there where nothing/something interesting happens.
325  subroutine a_refine_threshold(wlocal,xlocal,threshold,qt)
327  double precision, intent(in) :: wlocal(1:nw),xlocal(1:ndim),qt
328  double precision, intent(inout) :: threshold
329  end subroutine a_refine_threshold
330 
331  !> Allow user to use their own data-postprocessing procedures
332  subroutine special_convert(qunitconvert)
334  integer, intent(in) :: qunitconvert
335  character(len=20) :: userconvert_type
336  end subroutine special_convert
337 
338  !> flag=-1 : Treat all cells active, omit deactivation (onentry, default)
339  !> flag=0 : Treat as normal domain
340  !> flag=1 : Treat as passive, but reduce by safety belt
341  !> flag=2 : Always treat as passive
342  subroutine flag_grid(qt,ixI^L,ixO^L,w,x,flag)
344  integer, intent(in) :: ixI^L, ixO^L
345  integer, intent(inout) :: flag
346  double precision, intent(in) :: qt
347  double precision, intent(inout) :: w(ixi^s,1:nw)
348  double precision, intent(in) :: x(ixi^s,1:ndim)
349  end subroutine flag_grid
350 
351  !> Update payload of particles
352  subroutine update_payload(igrid,w,wold,xgrid,xpart,payload,npayload,particle_time)
354  integer, intent(in) :: igrid,npayload
355  double precision, intent(in) :: w(ixg^t,1:nw),wold(ixg^t,1:nw)
356  double precision, intent(in) :: xgrid(ixg^t,1:ndim),xpart(1:ndir),particle_time
357  double precision, intent(out) :: payload(npayload)
358  end subroutine update_payload
359 
360  subroutine create_particles(n_particles, x, v, q, m, follow)
361  integer, intent(in) :: n_particles
362  double precision, intent(out) :: x(3, n_particles)
363  double precision, intent(out) :: v(3, n_particles)
364  double precision, intent(out) :: q(n_particles)
365  double precision, intent(out) :: m(n_particles)
366  logical, intent(out) :: follow(n_particles)
367  end subroutine create_particles
368 
369  subroutine particle_fields(w, x, E, B)
371  double precision, intent(in) :: w(ixg^t,1:nw)
372  double precision, intent(in) :: x(ixg^t,1:ndim)
373  double precision, intent(out) :: E(ixg^t, ndir)
374  double precision, intent(out) :: B(ixg^t, ndir)
375  end subroutine particle_fields
376 
377  subroutine particle_analytic(ix, x, tloc, vec)
379  integer, intent(in) :: ix(ndir) !< Indices in gridvars
380  double precision, intent(in) :: x(ndir)
381  double precision, intent(in) :: tloc
382  double precision, intent(out) :: vec(ndir)
383  end subroutine particle_analytic
384 
385  subroutine after_refine(n_coarsen, n_refine)
386  integer, intent(in) :: n_coarsen
387  integer, intent(in) :: n_refine
388  end subroutine after_refine
389 
390  !> allow user to explicitly set flux at cell interfaces for finite volume scheme
391  subroutine set_flux(ixI^L,ixC^L,idim,fC)
393  integer, intent(in) :: ixI^L, ixC^L, idim
394  ! face-center flux
395  double precision,intent(inout) :: fC(ixi^s,1:nwflux,1:ndim)
396  ! For example, to set flux at bottom boundary in a 3D box for induction equation
397  ! vobs and bobs are interpolated data from original observational data for data-driven application
398  !integer :: idir
399  !if(idim==3) then
400  ! if(block%is_physical_boundary(idim*2-1)) then
401  ! do idir=1,ndir
402  ! fC(ixCmin3^%3ixC^S,mag(idir),idim)=vobs(ixCmin3+1^%3ixC^S,idim)*bobs(ixCmin3+1^%3ixC^S,idir)-vobs(ixCmin3+1^%3ixC^S,idir)*bobs(ixCmin3+1^%3ixC^S,idim)
403  ! end do
404  ! end if
405  !end if
406  end subroutine set_flux
407 
408  !> initialize vector potential on cell edges for magnetic field
409  subroutine init_vector_potential(ixI^L, ixC^L, xC, A, idir)
411 
412  integer, intent(in) :: ixI^L, ixC^L, idir
413  double precision, intent(in) :: xC(ixi^s,1:ndim)
414  double precision, intent(out) :: A(ixi^s)
415 
416  end subroutine init_vector_potential
417 
418  ! allow user to change inductive electric field, especially for boundary driven applications
419  subroutine set_electric_field(ixI^L,ixO^L,qdt,fE,s)
421  integer, intent(in) :: ixI^L, ixO^L
422  double precision, intent(in) :: qdt
423  type(state) :: s
424  double precision, intent(inout) :: fE(ixi^s,7-2*ndim:3)
425 
426  !integer :: ixC^L,ixA^L
427  ! For example, to set inductive electric field at bottom boundary in a 3D box for induction equation
428  ! v and b are from observational data for data-driven application
429 
430  !associate(w=>s%w,ws=>s%ws)
431 
432  !if(s%is_physical_boundary(5)) then
433  ! ixCmin^D=ixOmin^D-1;
434  ! ixCmax^D=ixOmax^D;
435  ! ixAmin^D=ixCmin^D;
436  ! ixAmax^D=ixCmax^D+1;
437  ! fE(nghostcells^%3ixA^S,1)=-ws(nghostcells^%3ixA^S,3)*w(nghostcells^%3ixA^S,mom(2))
438  ! fE(nghostcells^%3ixA^S,2)= ws(nghostcells^%3ixA^S,3)*w(nghostcells^%3ixA^S,mom(1))
439  ! ixAmin^D=ixCmin^D+kr(2,^D);
440  ! ixAmax^D=ixCmax^D+kr(2,^D);
441  ! fE(nghostcells^%3ixC^S,1)=0.5d0*(fE(nghostcells^%3ixC^S,1)+fE(nghostcells^%3ixA^S,1))*&
442  ! qdt*s%dsC(nghostcells^%3ixC^S,1)
443  ! ixAmin^D=ixCmin^D+kr(1,^D);
444  ! ixAmax^D=ixCmax^D+kr(1,^D);
445  ! fE(nghostcells^%3ixC^S,2)=0.5d0*(fE(nghostcells^%3ixC^S,2)+fE(nghostcells^%3ixA^S,2))*&
446  ! qdt*s%dsC(nghostcells^%3ixC^S,2)
447  !end if
448 
449  !end associate
450 
451  end subroutine set_electric_field
452 
453  end interface
454 
455 end module mod_usr_methods
procedure(special_resistivity), pointer usr_special_resistivity
This module contains definitions of global parameters and variables and some generic functions/subrou...
procedure(internal_bc), pointer usr_internal_bc
procedure(phys_dust_get_3d_dragforce), pointer usr_get_3d_dragforce
procedure(particle_fields), pointer usr_particle_fields
procedure(p_no_args), pointer usr_before_main_loop
initialize vector potential on cell edges for magnetic field
Set the "eta" array for resistive MHD based on w or the "current" variable which has components betwe...
procedure(get_dt), pointer usr_get_dt
this is the place to compute a local auxiliary variable to be used as refinement criterion for the Lo...
If defined, this routine is called before writing output, and it can set/modify the variables in the ...
procedure(process_adv_global), pointer usr_process_adv_global
Calculation anormal viscosity depending on space.
procedure(process_global), pointer usr_process_global
procedure(add_aux_names), pointer usr_add_aux_names
integer ndir
Number of spatial dimensions (components) for vector variables.
Calculation anormal pressure for hd & energy=.False.
procedure(process_adv_grid), pointer usr_process_adv_grid
Allow user to use their own data-postprocessing procedures.
procedure(sub_modify_io), pointer usr_modify_output
Add names for the auxiliary variables.
procedure(hd_pthermal), pointer usr_set_pthermal
procedure(refine_grid), pointer usr_refine_grid
regenerate w and eqpartf arrays to output into *tf.dat
allow user to explicitly set flux at cell interfaces for finite volume scheme
this subroutine is ONLY to be used for computing auxiliary variables which happen to be non-local (li...
procedure(p_no_args), pointer usr_improve_initial_condition
This subroutine is called at the beginning of each time step by each processor. No communication is s...
procedure(aux_output), pointer usr_aux_output
procedure(p_no_args), pointer usr_print_log
Module with all the methods that users can customize in AMRVAC.
Limit "dt" further if necessary, e.g. due to the special source terms. The getdt_courant (CFL conditi...
procedure(init_one_grid), pointer usr_init_one_grid
Initialize earch grid block data.
procedure(phys_gravity), pointer usr_gravity
procedure(process_grid), pointer usr_process_grid
for processing after the advance (PIC-MHD, e.g.)
procedure(flag_grid), pointer usr_flag_grid
Calculate the 3d drag force of gas onto dust.
procedure(after_refine), pointer usr_after_refine
procedure(source), pointer usr_source
Calculate gravitational acceleration in each dimension.
procedure(phys_visco), pointer usr_setvisco
procedure(p_no_args), pointer usr_set_parameters
Initialize the user&#39;s settings (after initializing amrvac)
use different threshold in special regions for AMR to reduce/increase resolution there where nothing/...
procedure(set_electric_field), pointer usr_set_electric_field
procedure(p_no_args), pointer usr_write_analysis
procedure(set_flux), pointer usr_set_flux
internal boundary, user defined This subroutine can be used to artificially overwrite ALL conservativ...
Calculate w(iw)=w(iw)+qdt*SOURCE[wCT,qtC,x] within ixO for all indices iw=iwmin...iwmax. wCT is at time qCT.
Enforce additional refinement or coarsening One can use the coordinate info in x and/or time qt=t_n a...
procedure(phys_dust_get_dt), pointer usr_dust_get_dt
procedure(init_vector_potential), pointer usr_init_vector_potential
procedure(set_b0), pointer usr_set_b0
procedure(a_refine_threshold), pointer usr_refine_threshold
Update payload of particles.
integer, parameter ndim
Number of spatial dimensions for grid variables.
procedure(particle_analytic), pointer usr_particle_analytic
procedure(update_payload), pointer usr_update_payload
procedure(transform_w), pointer usr_transform_w
integer nwauxio
Number of auxiliary variables that are only included in output.
procedure(var_for_errest), pointer usr_var_for_errest
for processing after the advance (PIC-MHD, e.g.)
flag=-1 : Treat all cells active, omit deactivation (onentry, default) flag=0 : Treat as normal domai...
procedure(create_particles), pointer usr_create_particles
special boundary types, users must assign conservative variables in boundaries
Calculate the time step associated with the usr drag force.
this subroutine can be used in convert, to add auxiliary variables to the converted output file...
procedure(set_j0), pointer usr_set_j0
procedure(special_bc), pointer usr_special_bc
procedure(special_convert), pointer usr_special_convert