MPI-AMRVAC  2.0
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, t^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 before the start of the simulation
36  procedure(p_no_args), pointer :: usr_before_main_loop => null()
37 
38  ! Source terms
39  procedure(source), pointer :: usr_source => null()
40  procedure(get_dt), pointer :: usr_get_dt => null()
41  procedure(phys_gravity), pointer :: usr_gravity => null()
42 
43  ! Usr defined space varying viscosity
44  procedure(phys_visco), pointer :: usr_setvisco => null()
45 
46  ! Refinement related procedures
47  procedure(refine_grid), pointer :: usr_refine_grid => null()
48  procedure(var_for_errest), pointer :: usr_var_for_errest => null()
49  procedure(a_refine_threshold), pointer :: usr_refine_threshold => null()
50  procedure(flag_grid), pointer :: usr_flag_grid => null()
51 
52  ! Set time-independent magnetic field for B0 splitting
53  procedure(set_b0), pointer :: usr_set_b0 => null()
54  ! Set time-independent current density for B0 splitting
55  procedure(set_j0), pointer :: usr_set_j0 => null()
56  procedure(special_resistivity), pointer :: usr_special_resistivity => null()
57 
58  ! Particles module related
59  procedure(update_payload), pointer :: usr_update_payload => null()
60  procedure(create_particles), pointer :: usr_create_particles => null()
61  procedure(particle_fields), pointer :: usr_particle_fields => null()
62  procedure(particle_analytic), pointer :: usr_particle_analytic => null()
63 
64  ! Called after the mesh has been adjuste
65  procedure(after_refine), pointer :: usr_after_refine => null()
66 
67  ! allow user to explicitly set flux at cell interfaces for finite volume scheme
68  procedure(set_flux), pointer :: usr_set_flux => null()
69 
70  abstract interface
71 
72  subroutine p_no_args()
73  end subroutine p_no_args
74 
75  !> Initialize one grid
76  subroutine init_one_grid(ixI^L,ixO^L,w,x)
78  integer, intent(in) :: ixI^L, ixO^L
79  double precision, intent(in) :: x(ixi^s,1:ndim)
80  double precision, intent(inout) :: w(ixi^s,1:nw)
81  end subroutine init_one_grid
82 
83  !> special boundary types, users must assign conservative
84  !> variables in boundaries
85  subroutine special_bc(qt,ixI^L,ixO^L,iB,w,x)
87  !> Shape of input arrays
88  integer, intent(in) :: ixI^L
89  !> Region where boundary values have to be set
90  integer, intent(in) :: ixO^L
91  !> Integer indicating direction of boundary
92  integer, intent(in) :: iB
93  double precision, intent(in) :: qt, x(ixi^s,1:ndim)
94  double precision, intent(inout) :: w(ixi^s,1:nw)
95  end subroutine special_bc
96 
97  !> internal boundary, user defined
98  !
99  !> This subroutine can be used to artificially overwrite ALL conservative
100  !> variables in a user-selected region of the mesh, and thereby act as
101  !> an internal boundary region. It is called just before external (ghost cell)
102  !> boundary regions will be set by the BC selection. Here, you could e.g.
103  !> want to introduce an extra variable (nwextra, to be distinguished from nwaux)
104  !> which can be used to identify the internal boundary region location.
105  !> Its effect should always be local as it acts on the mesh.
106  subroutine internal_bc(level,qt,ixI^L,ixO^L,w,x)
108  integer, intent(in) :: ixI^L,ixO^L,level
109  double precision, intent(in) :: qt
110  double precision, intent(inout) :: w(ixi^s,1:nw)
111  double precision, intent(in) :: x(ixi^s,1:ndim)
112  end subroutine internal_bc
113 
114  !> this subroutine is ONLY to be used for computing auxiliary variables
115  !> which happen to be non-local (like div v), and are in no way used for
116  !> flux computations. As auxiliaries, they are also not advanced
117  subroutine process_grid(igrid,level,ixI^L,ixO^L,qt,w,x)
119  integer, intent(in) :: igrid,level,ixI^L,ixO^L
120  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
121  double precision, intent(inout) :: w(ixi^s,1:nw)
122  end subroutine process_grid
123 
124  !> If defined, this routine is called before writing output, and it can
125  !> set/modify the variables in the w array.
126  subroutine sub_modify_io(ixI^L,ixO^L,qt,w,x)
128  integer, intent(in) :: ixI^L,ixO^L
129  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
130  double precision, intent(inout) :: w(ixi^s,1:nw)
131  end subroutine sub_modify_io
132 
133  !> This subroutine is called at the beginning of each time step
134  !> by each processor. No communication is specified, so the user
135  !> has to implement MPI routines if information has to be shared
136  subroutine process_global(iit,qt)
138  integer, intent(in) :: iit
139  double precision, intent(in) :: qt
140  end subroutine process_global
141 
142  !> for processing after the advance (PIC-MHD, e.g.)
143  subroutine process_adv_grid(igrid,level,ixI^L,ixO^L,qt,w,x)
145  integer, intent(in) :: igrid,level,ixI^L,ixO^L
146  double precision, intent(in) :: qt,x(ixi^s,1:ndim)
147  double precision, intent(inout) :: w(ixi^s,1:nw)
148  end subroutine process_adv_grid
149 
150  !> for processing after the advance (PIC-MHD, e.g.)
151  subroutine process_adv_global(iit,qt)
153  integer, intent(in) :: iit
154  double precision, intent(in) :: qt
155  end subroutine process_adv_global
156 
157  !> this subroutine can be used in convert, to add auxiliary variables to the
158  !> converted output file, for further analysis using tecplot, paraview, ....
159  !> these auxiliary values need to be stored in the nw+1:nw+nwauxio slots
160  !
161  !> the array normconv can be filled in the (nw+1:nw+nwauxio) range with
162  !> corresponding normalization values (default value 1)
163  subroutine aux_output(ixI^L,ixO^L,w,x,normconv)
165  integer, intent(in) :: ixI^L,ixO^L
166  double precision, intent(in) :: x(ixi^s,1:ndim)
167  double precision :: w(ixi^s,nw+nwauxio)
168  double precision :: normconv(0:nw+nwauxio)
169  end subroutine aux_output
170 
171  !> Add names for the auxiliary variables
172  subroutine add_aux_names(varnames)
174  character(len=*) :: varnames
175  end subroutine add_aux_names
176 
177  !> Calculate w(iw)=w(iw)+qdt*SOURCE[wCT,qtC,x] within ixO for all indices
178  !> iw=iwmin...iwmax. wCT is at time qCT
179  subroutine source(qdt,ixI^L,ixO^L,iw^LIM,qtC,wCT,qt,w,x)
181  integer, intent(in) :: ixI^L, ixO^L, iw^LIM
182  double precision, intent(in) :: qdt, qtC, qt
183  double precision, intent(in) :: wCT(ixi^s,1:nw), x(ixi^s,1:ndim)
184  double precision, intent(inout) :: w(ixi^s,1:nw)
185  end subroutine source
186 
187  !> Limit "dt" further if necessary, e.g. due to the special source terms.
188  !> The getdt_courant (CFL condition) and the getdt subroutine in the AMRVACPHYS
189  !> module have already been called.
190  subroutine get_dt(w,ixI^L,ixO^L,dtnew,dx^D,x)
192  integer, intent(in) :: ixI^L, ixO^L
193  double precision, intent(in) :: dx^D, x(ixi^s,1:ndim)
194  double precision, intent(in) :: w(ixi^s,1:nw)
195  double precision, intent(inout) :: dtnew
196  end subroutine get_dt
197 
198  !> Calculate gravitational acceleration in each dimension
199  subroutine phys_gravity(ixI^L,ixO^L,wCT,x,gravity_field)
201  integer, intent(in) :: ixI^L, ixO^L
202  double precision, intent(in) :: x(ixi^s,1:ndim)
203  double precision, intent(in) :: wCT(ixi^s,1:nw)
204  double precision, intent(out) :: gravity_field(ixi^s,ndim)
205  end subroutine phys_gravity
206 
207  !>Calculation anormal viscosity depending on space
208  subroutine phys_visco(ixI^L,ixO^L,x,w,mu)
210  integer, intent(in) :: ixI^L, ixO^L
211  double precision, intent(in) :: x(ixi^s,1:ndim)
212  double precision, intent(in) :: w(ixi^s,1:nw)
213  double precision, intent(out) :: mu(ixi^s)
214  end subroutine phys_visco
215 
216  !> Set the "eta" array for resistive MHD based on w or the
217  !> "current" variable which has components between idirmin and 3.
218  subroutine special_resistivity(w,ixI^L,ixO^L,idirmin,x,current,eta)
220  integer, intent(in) :: ixI^L, ixO^L, idirmin
221  double precision, intent(in) :: w(ixi^s,nw), x(ixi^s,1:ndim)
222  double precision :: current(ixi^s,7-2*ndir:3), eta(ixi^s)
223  end subroutine special_resistivity
224 
225  !> Enforce additional refinement or coarsening
226  !> One can use the coordinate info in x and/or time qt=t_n and w(t_n) values w.
227  !> you must set consistent values for integers refine/coarsen:
228  !> refine = -1 enforce to not refine
229  !> refine = 0 doesn't enforce anything
230  !> refine = 1 enforce refinement
231  !> coarsen = -1 enforce to not coarsen
232  !> coarsen = 0 doesn't enforce anything
233  !> coarsen = 1 enforce coarsen
234  !> e.g. refine for negative first coordinate x < 0 as
235  !> if (any(x(ix^S,1) < zero)) refine=1
236  subroutine refine_grid(igrid,level,ixI^L,ixO^L,qt,w,x,refine,coarsen)
238  integer, intent(in) :: igrid, level, ixI^L, ixO^L
239  double precision, intent(in) :: qt, w(ixi^s,1:nw), x(ixi^s,1:ndim)
240  integer, intent(inout) :: refine, coarsen
241  end subroutine refine_grid
242 
243  !> this is the place to compute a local auxiliary variable to be used
244  !> as refinement criterion for the Lohner error estimator only
245  !> -->it is then requiring and iflag>nw
246  !> note that ixO=ixI=ixG, hence the term local (gradients need special attention!)
247  subroutine var_for_errest(ixI^L,ixO^L,iflag,w,x,var)
249  integer, intent(in) :: ixI^L,ixO^L,iflag
250  double precision, intent(in) :: w(ixi^s,1:nw), x(ixi^s,1:ndim)
251  double precision, intent(out) :: var(ixi^s)
252  end subroutine var_for_errest
253 
254  !> Here one can add a steady (time-independent) potential background field
255  subroutine set_b0(ixI^L,ixO^L,x,wB0)
257  integer, intent(in) :: ixI^L,ixO^L
258  double precision, intent(in) :: x(ixi^s,1:ndim)
259  double precision, intent(inout) :: wB0(ixi^s,1:ndir)
260  end subroutine set_b0
261 
262  !> Here one can add a time-independent background current density
263  subroutine set_j0(ixI^L,ixO^L,x,wJ0)
265  integer, intent(in) :: ixI^L,ixO^L
266  double precision, intent(in) :: x(ixi^s,1:ndim)
267  double precision, intent(inout) :: wJ0(ixi^s,7-2*ndir:ndir)
268  end subroutine set_j0
269 
270  !> regenerate w and eqpartf arrays to output into *tf.dat
271  subroutine transform_w(ixI^L,ixO^L,nw_in,w_in,x,w_out)
273  integer, intent(in) :: ixI^L, ixO^L, nw_in
274  double precision, intent(in) :: w_in(ixi^s,1:nw_in)
275  double precision, intent(in) :: x(ixi^s, 1:ndim)
276  double precision, intent(out) :: w_out(ixi^s,1:nw)
277  end subroutine transform_w
278 
279  !> use different threshold in special regions for AMR to
280  !> reduce/increase resolution there where nothing/something interesting happens.
281  subroutine a_refine_threshold(wlocal,xlocal,threshold,qt)
283  double precision, intent(in) :: wlocal(1:nw),xlocal(1:ndim),qt
284  double precision, intent(inout) :: threshold
285  end subroutine a_refine_threshold
286 
287  !> Allow user to use their own data-postprocessing procedures
288  subroutine special_convert(qunitconvert)
290  integer, intent(in) :: qunitconvert
291  character(len=20) :: userconvert_type
292  end subroutine special_convert
293 
294  !> flag=-1 : Treat all cells active, omit deactivation (onentry, default)
295  !> flag=0 : Treat as normal domain
296  !> flag=1 : Treat as passive, but reduce by safety belt
297  !> flag=2 : Always treat as passive
298  subroutine flag_grid(qt,ixI^L,ixO^L,w,x,flag)
300  integer, intent(in) :: ixI^L, ixO^L
301  integer, intent(inout) :: flag
302  double precision, intent(in) :: qt
303  double precision, intent(inout) :: w(ixi^s,1:nw)
304  double precision, intent(in) :: x(ixi^s,1:ndim)
305  end subroutine flag_grid
306 
307  !> Update payload of particles
308  subroutine update_payload(igrid,w,wold,xgrid,xpart,payload,npayload,particle_time)
310  integer, intent(in) :: igrid,npayload
311  double precision, intent(in) :: w(ixg^t,1:nw),wold(ixg^t,1:nw)
312  double precision, intent(in) :: xgrid(ixg^t,1:ndim),xpart(1:ndir),particle_time
313  double precision, intent(out) :: payload(npayload)
314  end subroutine update_payload
315 
316  subroutine create_particles(n_particles, x, v, q, m, follow)
317  integer, intent(in) :: n_particles
318  double precision, intent(out) :: x(3, n_particles)
319  double precision, intent(out) :: v(3, n_particles)
320  double precision, intent(out) :: q(n_particles)
321  double precision, intent(out) :: m(n_particles)
322  logical, intent(out) :: follow(n_particles)
323  end subroutine create_particles
324 
325  subroutine particle_fields(w, x, E, B)
327  double precision, intent(in) :: w(ixg^t,1:nw)
328  double precision, intent(in) :: x(ixg^t,1:ndim)
329  double precision, intent(out) :: E(ixg^t, ndir)
330  double precision, intent(out) :: B(ixg^t, ndir)
331  end subroutine particle_fields
332 
333  subroutine particle_analytic(ix, x, tloc, vec)
335  integer, intent(in) :: ix(ndir) !< Indices in gridvars
336  double precision, intent(in) :: x(ndir)
337  double precision, intent(in) :: tloc
338  double precision, intent(out) :: vec(ndir)
339  end subroutine particle_analytic
340 
341  subroutine after_refine(n_coarsen, n_refine)
342  integer, intent(in) :: n_coarsen
343  integer, intent(in) :: n_refine
344  end subroutine after_refine
345 
346  !> allow user to explicitly set flux at cell interfaces for finite volume scheme
347  subroutine set_flux(ixI^L,ixC^L,idim,fC)
349  integer, intent(in) :: ixI^L, ixC^L, idim
350  ! face-center flux
351  double precision,intent(inout) :: fC(ixi^s,1:nwflux,1:ndim)
352  ! For example, to set flux at bottom boundary in a 3D box for induction equation
353  ! vobs and bobs are interpolated data from original observational data for data-driven application
354  !integer :: idir
355  !if(idim==3) then
356  ! if(block%is_physical_boundary(idim*2-1)) then
357  ! do idir=1,ndir
358  ! 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)
359  ! end do
360  ! end if
361  !end if
362  end subroutine set_flux
363 
364  end interface
365 
366 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(particle_fields), pointer usr_particle_fields
procedure(p_no_args), pointer usr_before_main_loop
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.
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(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...
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
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(p_no_args), pointer usr_write_analysis
procedure(set_flux), pointer usr_set_flux
internal boundary, user defined
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(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
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