MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_global_parameters.t
Go to the documentation of this file.
1 !> This module contains definitions of global parameters and variables and some
2 !> generic functions/subroutines used in AMRVAC.
3 !>
4 !> \todo Move the parameters to the relevant (physics) modules
8  use mpi
9  use mod_constants
10  use mod_variables
11  use mod_basic_types
12 
13  implicit none
14  public
15 
16  ! Parameters
17 
18  character(len=*), parameter :: undefined = 'undefined'
19 
20  !> @todo Move mpi related variables to e.g. mod_comm
21 
22  !> The number of MPI tasks
23  integer :: npe
24 
25  !> The rank of the current MPI task
26  integer :: mype
27 
28  !> The MPI communicator
29  integer :: icomm
30 
31  !> A global MPI error return code
32  !> @todo Make local
33  integer :: ierrmpi
34 
35  !> MPI file handle for logfile
36  integer :: log_fh
37  !> MPI type for block including ghost cells and its size
38  integer :: type_block, size_block
39  !> MPI type for block coarsened by 2, and for its children blocks
41  !> MPI type for staggered block coarsened by 2, and for its children blocks
42  integer :: type_coarse_block_stg(^nd,2^d&), type_sub_block_stg(^ND,2^D&)
43  !> MPI type for IO: block excluding ghost cells
45  !> MPI type for IO of staggered variables
47  !> MPI type for IO: cell corner (xc) or cell center (xcc) coordinates
49  !> MPI type for IO: cell corner (wc) or cell center (wcc) variables
51 
52  !> MPI recv send variables for AMR
53  integer :: itag, irecv, isend
54  integer, dimension(:), allocatable :: recvrequest, sendrequest
55  integer, dimension(:,:), allocatable :: recvstatus, sendstatus
56  !> MPI recv send variables for staggered-variable AMR
57  integer :: itag_stg
58  integer, dimension(:), allocatable :: recvrequest_stg, sendrequest_stg
59  integer, dimension(:,:), allocatable :: recvstatus_stg, sendstatus_stg
60 
61  ! geometry and domain setups
62 
63  !> the mesh range (within a block with ghost cells)
64  integer :: ixm^ll
65 
66  !> minimum and maximum domain boundaries for each dimension
67  double precision :: xprob^l
68 
69  !> Indices for cylindrical coordinates FOR TESTS, negative value when not used:
70  integer :: r_ = -1
71  integer :: phi_ = -1
72  integer :: z_ = -1
73 
74  !> Number of spatial dimensions for grid variables
75  integer, parameter :: ndim=^nd
76 
77  !> Number of spatial dimensions (components) for vector variables
78  integer :: ndir=ndim
79 
80  !> Cartesian geometry or not
81  logical :: slab
82 
83  !> uniform Cartesian geometry or not (stretched Cartesian)
84  logical :: slab_uniform
85 
86  !> number of grid blocks in domain per dimension, in array over levels
87  integer, dimension(:), allocatable :: ng^d
88  !> extent of grid blocks in domain per dimension, in array over levels
89  double precision, dimension(:), allocatable :: dg^d
90 
91  !> number of cells for each dimension in level-one mesh
92  integer :: domain_nx^d
93 
94  !> number of cells for each dimension in grid block excluding ghostcells
95  integer :: block_nx^d
96 
97  !> Lower index of grid block arrays (always 1)
98  integer, parameter :: {ixglo^d = 1|, }
99 
100  !> Upper index of grid block arrays
101  integer :: ixghi^d
102 
103  !> Lower index of stagger grid block arrays (always 0)
104  integer, parameter :: {ixgslo^d = 0|, }
105 
106  !> Upper index of stagger grid block arrays
107  integer :: ixgshi^d
108 
109  !> Number of ghost cells surrounding a grid
110  integer :: nghostcells
111 
112  integer, parameter :: stretch_none = 0 !< No stretching
113  integer, parameter :: stretch_uni = 1 !< Unidirectional stretching from a side
114  integer, parameter :: stretch_symm = 2 !< Symmetric stretching around the center
115 
116  !> If true, adjust mod_geometry routines to account for grid stretching (but
117  !> the flux computation will not)
119  !> True if a dimension is stretched
120  logical :: stretched_dim(ndim)
121  !> What kind of stretching is used per dimension
122  integer :: stretch_type(ndim)
123  !> stretch factor between cells at AMR level 1, per dimension
124  double precision :: qstretch_baselevel(ndim)
125  !> (even) number of (symmetrically) stretched
126  !> blocks at AMR level 1, per dimension
128  !> (even) number of (symmetrically) stretched blocks per level and dimension
129  integer, allocatable :: nstretchedblocks(:,:)
130  !> physical extent of stretched border in symmetric stretching
131  double precision :: xstretch^d
132  !> Stretching factors and first cell size for each AMR level and dimension
133  double precision, allocatable :: qstretch(:,:), dxfirst(:,:), &
134  dxfirst_1mq(:,:), dxmid(:,:)
135 
136  !> grid hierarchy info (level and grid indices)
137  integer, parameter :: nodehi=^nd+1
138  integer, parameter :: plevel_=1
139  integer, parameter :: pig^d_=plevel_+^d
140 
141  integer, allocatable :: node(:,:)
142  integer, allocatable :: node_sub(:,:)
143 
144  !> grid location info (corner coordinates and grid spacing)
145  integer, parameter :: rnodehi=3*^nd
146  integer, parameter :: rpxmin0_=0
147  integer, parameter :: rpxmin^d_=rpxmin0_+^d
148  integer, parameter :: rpxmax0_=^nd
149  integer, parameter :: rpxmax^d_=rpxmax0_+^d
150  integer, parameter :: rpdx^d_=2*^nd+^d
151 
152  !> Corner coordinates
153  double precision, allocatable :: rnode(:,:)
154  double precision, allocatable :: rnode_sub(:,:)
155 
156  double precision, allocatable :: dx(:,:)
157  double precision :: dxlevel(ndim)
158 
159  ! IO related quantities
160 
161  !> Maximum number of saves that can be defined by tsave or itsave
162  integer, parameter :: nsavehi=100
163 
164  !> Number of output methods
165  integer, parameter :: nfile = 5
166 
167  !> Names of the output methods
168  character(len=40), parameter :: output_names(nfile) = &
169  ['log ', 'normal ', 'slice ', 'collapsed', 'analysis ']
170 
171  !> User parameter file
172  character(len=std_len) :: usr_filename
173 
174  !> If collapse(DIM) is true, generate output integrated over DIM
175  logical :: collapse(ndim)
176 
177  !> Save output of type N on times tsave(:, N)
178  double precision :: tsave(nsavehi,nfile)
179 
180  !> \todo Move tsavelast to amrvac.t
181  double precision :: tsavelast(nfile)
182 
183  !> Repeatedly save output of type N when dtsave(N) simulation time has passed
184  double precision :: dtsave(nfile)
185 
186  !> Save output of type N on iterations itsave(:, N)
187  integer :: itsave(nsavehi,nfile)
188 
189  !> \todo remove itsavelast?
190  integer :: itsavelast(nfile)
191 
192  !> Repeatedly save output of type N when ditsave(N) time steps have passed
193  integer :: ditsave(nfile)
194 
195  !> \todo Move to amrvac.t
196  integer :: isavet(nfile)
197 
198  !> \todo Move to amrvac.t
199  integer :: isaveit(nfile)
200 
201  !> Start of read out (not counting specified read outs)
202  double precision :: tsavestart(nfile)
203 
204  !> The level at which to produce line-integrated / collapsed output
205  integer :: collapselevel
206 
207  !> Number of saved files of each type
208  !> \todo Move to mod_input_output
209  integer :: n_saves(1:nfile)
210 
211  !> whether or not to save an output file
212  logical :: save_file(nfile)
213 
214  !> to monitor timeintegration loop at given wall-clock time intervals
215  double precision :: time_between_print
216 
217  !> accumulated wall-clock time spent on boundary conditions
218  double precision :: time_bc
219 
220  !> IO: snapshot and collapsed views output numbers/labels
222 
223  !> Constant indicating log output
224  integer, parameter :: filelog_ = 1
225 
226  !> Constant indicating regular output
227  integer, parameter :: fileout_ = 2
228 
229  !> Constant indicating slice output
230  integer, parameter :: fileslice_ = 3
231 
232  !> Constant indicating collapsed output
233  integer, parameter :: filecollapse_ = 4
234 
235  !> Constant indicating analysis output (see @ref analysis.md)
236  integer, parameter :: fileanalysis_ = 5
237 
238  !> Unit for standard input
239  integer, parameter :: unitstdin=5
240 
241  !> Unit for standard output
242  integer, parameter :: unitterm=6
243 
244  !> Unit for error messages
245  integer, parameter :: uniterr=6
246 
247  !> \todo Move to mod_input_output
248  integer, parameter :: unitpar=9
249  integer, parameter :: unitconvert=10
250  integer, parameter :: unitslice=11
251  integer, parameter :: unitsnapshot=12
252  integer, parameter :: unitcollapse=13
253  integer, parameter :: unitanalysis=14
254 
255  !> Number of auxiliary variables that are only included in output
256  integer :: nwauxio
257 
258  !> IO switches for conversion
259  logical :: nocartesian
260  logical, allocatable :: w_write(:)
261  logical, allocatable :: writelevel(:)
262  double precision :: writespshift(ndim,2)
264 
265  !> Which par files are used as input
266  character(len=std_len), allocatable :: par_files(:)
267 
268  !> Base file name for simulation output, which will be followed by a number
269  character(len=std_len) :: base_filename
270 
271  !> If not 'unavailable', resume from snapshot with this base file name
272  character(len=std_len) :: restart_from_file
273 
274  !> Which type of log to write: 'normal', 'special', 'regression_test'
275  character(len=std_len) :: typefilelog
276 
277  !> Resume from the snapshot with this index
278  integer :: snapshotini
279 
280  !> If true, restart a previous run from the latest snapshot
282 
283  !> If true and restart_from_file is given, convert snapshots to
284  !> other file formats
285  logical :: convert
286 
287  !> If true, already convert to output format during the run
288  logical :: autoconvert
289 
290  !> If true, convert from conservative to primitive variables in output
291  logical :: saveprim
292 
293  !> Which format to use when converting
294  !>
295  !> Options are: tecplot, tecplotCC, vtu, vtuCC, vtuB, vtuBCC,
296  !> tecplotmpi, tecplotCCmpi, vtumpi, vtuCCmpi, vtuBmpi, vtuBCCmpi, pvtumpi, pvtuCCmpi,
297  !> pvtuBmpi, pvtuBCCmpi, tecline, teclinempi, onegrid
298  character(len=std_len) :: convert_type
299 
300  character(len=std_len) :: collapse_type
301 
302  !> Conversion factors the primitive variables
303  double precision, allocatable :: w_convert_factor(:)
304 
305  double precision :: length_convert_factor
306 
307  !> Conversion factor for time unit
308  double precision :: time_convert_factor
309 
310  integer :: saveigrid
311 
312  !> Stores the memory and load imbalance, used in printlog
313  double precision :: xload, xmemory
314 
315  !> Save a snapshot before crash a run met unphysical values
316  logical :: crash=.false.
317 
318  ! Physics factors
319 
320  !> Physical scaling factor for length
321  double precision :: unit_length=1.d0
322 
323  !> Physical scaling factor for time
324  double precision :: unit_time=1.d0
325 
326  !> Physical scaling factor for density
327  double precision :: unit_density=1.d0
328 
329  !> Physical scaling factor for velocity
330  double precision :: unit_velocity=0.d0
331 
332  !> Physical scaling factor for temperature
333  double precision :: unit_temperature=1.d0
334 
335  !> Physical scaling factor for pressure
336  double precision :: unit_pressure=1.d0
337 
338  !> Physical scaling factor for magnetic field
339  double precision :: unit_magneticfield=1.d0
340 
341  !> Physical scaling factor for number density
342  double precision :: unit_numberdensity=1.d0
343 
344  !> error handling
346 
347  !> amplitude of background dipolar, quadrupolar, octupolar, user's field
348  double precision :: bdip=0.d0
349  double precision :: bquad=0.d0
350  double precision :: boct=0.d0
351  double precision :: busr=0.d0
352 
353  !> check and optionally fix unphysical small values (density, gas pressure)
354  logical :: check_small_values=.false.
355 
356  !> Use primitive variables when correcting small values
357  logical :: small_values_use_primitive=.false.
358 
359  !> Whether to apply small value fixes to certain variables
360  logical, allocatable :: small_values_fix_iw(:)
361 
362  !> split magnetic field as background B0 field
363  logical :: b0field=.false.
364 
365  !> Use SI units (.true.) or use cgs units (.false.)
366  logical :: si_unit=.false.
367 
368  !> Solve energy equation or not
369  logical :: phys_energy=.true.
370 
371  !> Solve polytropic process instead of solving total energy
372  logical :: solve_internal_e=.false.
373 
374  !> Enable to strictly conserve the angular momentum
375  !> (works both in cylindrical and spherical coordinates)
376  logical :: angmomfix=.false.
377 
378  !> Use particles module or not
379  logical :: use_particles=.false.
380 
381  !> Use multigrid (only available in 2D and 3D)
382  logical :: use_multigrid = .false.
383 
384  ! AMR switches
385 
386  !> The maximum number of grid blocks in a processor
387  integer :: max_blocks
388 
389  !> The maximum number of levels in the grid refinement
390  integer, parameter :: nlevelshi = 20
391 
392  !> Maximal number of AMR levels
393  integer :: refine_max_level
394 
395  !> Weights of variables used to calculate error for mesh refinement
396  double precision, allocatable :: w_refine_weight(:)
397 
398  !> Fix the AMR grid after this time
399  double precision :: tfixgrid
400 
401  !> Fix the AMR grid after this many time steps
402  integer :: itfixgrid
403 
404  !> Reconstruct the AMR grid once every ditregrid iteration(s)
405  integer :: ditregrid
406 
407  !> refinement: lohner estimate wavefilter setting
408  double precision, allocatable :: amr_wavefilter(:)
409 
410  integer :: refine_criterion
411  logical :: prolongprimitive
412  logical :: coarsenprimitive
413 
414  !> Error tolerance for refinement decision
415  double precision, allocatable :: refine_threshold(:)
416  double precision, allocatable :: derefine_ratio(:)
417 
418  !> If true, rebuild the AMR grid upon restarting
419  logical :: reset_grid
420  !> True for using stagger grid
421  logical :: stagger_grid=.false.
422 
423  !> Number of cells as buffer zone
424  !> \todo is it necessary?
425  integer :: nbufferx^d
426 
427  integer :: levmin
428  integer :: levmax
429  integer :: levmax_sub
430 
431  ! Miscellaneous
432 
433  !> problem switch allowing different setups in same usr_mod.t
434  integer :: iprob
435 
436  !> Kronecker delta tensor
437  integer :: kr(3,3)
438 
439  !> Levi-Civita tensor
440  integer :: lvc(3,3,3)
441 
442  ! Time integration aspects
443 
444  double precision :: dt
445  double precision, allocatable :: dt_grid(:)
446 
447  logical :: time_advance
448 
449  !> The Courant (CFL) number used for the simulation
450  double precision :: courantpar
451 
452  !> How to compute the CFL-limited time step.
453  !>
454  !> Options are 'maxsum': max(sum(c/dx)); 'summax': sum(max(c/dx)) and
455  !> 'minimum: max(c/dx), where the summations loop over the grid dimensions and
456  !> c is the velocity. The default 'maxsum' is the conventiontal way of
457  !> computing CFL-limited time steps.
458  character(len=std_len) :: typecourant
459 
460  !> If dtpar is positive, it sets the timestep dt, otherwise courantpar is used
461  !> to limit the time step based on the Courant condition.
462  double precision :: dtpar
463 
464  !> For resistive MHD, the time step is also limited by the diffusion time:
465  !> \f$ dt < dtdiffpar \times dx^2/eta \f$
466  double precision :: dtdiffpar
467 
468  !> The global simulation time
469  double precision :: global_time
470 
471  !> Start time for the simulation
472  double precision :: time_init
473 
474  !> End time for the simulation
475  double precision :: time_max
476 
477  !> Ending wall time (in hours) for the simulation
478  double precision :: wall_time_max
479 
480  !> Stop the simulation when the time step becomes smaller than this value
481  double precision :: dtmin
482 
483  !> If true, reset iteration count and global_time to original values, and
484  !> start writing snapshots at index 0
485  logical :: reset_time
486 
487  !> If true, reset iteration count to 0
488  logical :: reset_it
489 
490  !> If true, call initonegrid_usr upon restarting
491  logical :: firstprocess
492 
493  !> If true, wall time is up, modify snapshotnext for later overwrite
494  logical :: pass_wall_time
495 
496  !> Number of time steps taken
497  integer :: it
498 
499  !> Stop the simulation after this many time steps have been taken
500  integer :: it_max
501 
502  !> initial iteration count
503  integer :: it_init
504 
505  !> If > 1, then in the first slowsteps-1 time steps dt is reduced
506  !> by a factor \f$ 1 - (1- step/slowsteps)^2 \f$
507  integer :: slowsteps
508 
509  ! Method switches
510 
511  !> Index of the sub-step in a multi-step time integrator
512  integer :: istep
513 
514  !> How many sub-steps the time integrator takes
515  integer :: nstep
516 
517  !> Which time integrator to use
518  character(len=std_len) :: time_integrator
519 
520  !> How to apply dimensional splitting to the source terms, see
521  !> @ref disretization.md
522  character(len=std_len) :: typesourcesplit
523 
524  !> Which spatial discretization to use (per grid level)
525  character(len=std_len), allocatable :: flux_scheme(:)
526 
527  !> The spatial dicretization for the predictor step when using a two
528  !> step method
529  character(len=std_len), allocatable :: typepred1(:)
530 
531  !> Type of slope limiter used for reconstructing variables on cell edges
532  integer, allocatable :: type_limiter(:)
533 
534  !> Type of slope limiter used for computing gradients or divergences, when
535  !> typegrad or typediv are set to 'limited'
536  integer, allocatable :: type_gradient_limiter(:)
537 
538  !> \todo Remove / replace with limiter
539  integer :: typelimiter
540 
541  !> \todo Remove / replace with gradient_limiter
542  integer :: typegradlimiter
543 
544  !> Limiter used for prolongation to refined grids and ghost cells
545  character(len=std_len) :: typeprolonglimit
546 
547  !> Which type of entropy fix to use with Riemann-type solvers
548  character(len=std_len), allocatable :: typeentropy(:)
549 
550  !> Which type of TVD method to use
551  character(len=std_len) :: typetvd
552 
553  !> Which type of TVDLF method to use
554  character(len=std_len) :: typeboundspeed
555 
556  character(len=std_len) :: typeaverage
557  character(len=std_len) :: typedimsplit
558  character(len=std_len) :: geometry_name='default'
559  character(len=std_len) :: typepoly
560 
561  integer :: nxdiffusehllc
562  double precision, allocatable :: entropycoef(:)
563  double precision :: tvdlfeps
564  logical, allocatable :: loglimit(:), logflag(:)
565  logical :: flathllc,flatcd,flatsh
566  !> Use split or unsplit way to add user's source terms, default: unsplit
567  logical :: source_split_usr
568  logical :: dimsplit
569 
570  character(len=std_len) :: typediv,typegrad,typecurl
571 
572  !> global fastest wave speed needed in fd scheme and glm method
573  double precision :: cmax_global
574 
575  !> global fastest flow speed needed in glm method
576  double precision :: vmax_global
577 
578  !> need global maximal wave speed
579  logical :: need_global_cmax=.false.
580 
581  ! Boundary region parameters
582 
583  !> True for dimensions with periodic boundaries
584  logical :: periodb(ndim)
585 
586  !> Indicates whether there is a pole at a boundary
587  logical :: poleb(2,ndim)
588 
589  !> True for dimensions with aperiodic boundaries
590  logical :: aperiodb(ndim)
591 
592  !> True for save physical boundary cells in dat files
594 
595  !> True if a block has any physical boundary
596  logical, allocatable :: phyboundblock(:)
597 
598  !> Array indicating the type of boundary condition per variable and per
599  !> physical boundary
600  character(len=std_len), allocatable :: typeboundary(:, :)
601 
602  character(len=std_len) :: typeghostfill='linear',prolongation_method
603  logical :: internalboundary
604 
605  ! parameters for bc_phys
606  integer, parameter :: ismin^d=-1+2*^d
607  integer, parameter :: ismax^d=2*^d
608 
609  !> Base file name for synthetic EUV emission output
610  character(len=std_len) :: filename_euv
611  !> output image
612  logical :: image=.false.
613  !> output spectrum
614  logical :: spectrum=.false.
615  ! wavelength for output
616  integer :: wavelength
617  !> direction of light of sight
618  integer :: direction_los=3
619  !> direction of slit (for spectrum)
620  integer :: direction_slit=2
621  !> location of the slit
622  double precision :: location_slit=0.d0
623  !> resolution of the output
624  character(len=std_len) :: resolution_euv='instrument'
625 
626  !$OMP THREADPRIVATE(dxlevel)
627  !$OMP THREADPRIVATE(saveigrid)
628  !$OMP THREADPRIVATE(typelimiter,typegradlimiter)
629 
630 contains
631 
632  !> Cross product of two vectors
633  pure subroutine cross_product(ixI^L,ixO^L,a,b,axb)
634  integer, intent(in) :: ixI^L, ixO^L
635  double precision, intent(in) :: a(ixi^s,3), b(ixi^s,3)
636  double precision, intent(out) :: axb(ixi^s,3)
637  !-------------------------------------------------------------------------
638 
639  axb(ixo^s,1)=a(ixo^s,2)*b(ixo^s,3)-a(ixo^s,3)*b(ixo^s,2)
640  axb(ixo^s,2)=a(ixo^s,3)*b(ixo^s,1)-a(ixo^s,1)*b(ixo^s,3)
641  axb(ixo^s,3)=a(ixo^s,1)*b(ixo^s,2)-a(ixo^s,2)*b(ixo^s,1)
642  end subroutine cross_product
643 
644 end module mod_global_parameters
This module contains definitions of global parameters and variables and some generic functions/subrou...
character(len=std_len) prolongation_method
integer direction_los
direction of light of sight
integer, parameter rnodehi
grid location info (corner coordinates and grid spacing)
character(len=std_len) time_integrator
Which time integrator to use.
double precision vmax_global
global fastest flow speed needed in glm method
double precision, dimension(:), allocatable refine_threshold
Error tolerance for refinement decision.
character(len=std_len), dimension(:), allocatable typepred1
The spatial dicretization for the predictor step when using a two step method.
logical, dimension(:), allocatable w_write
logical small_values_use_primitive
Use primitive variables when correcting small values.
integer domain_nx
number of cells for each dimension in level-one mesh
double precision time_max
End time for the simulation.
double precision, dimension(:), allocatable w_refine_weight
Weights of variables used to calculate error for mesh refinement.
character(len=std_len) typetvd
Which type of TVD method to use.
integer, parameter unitslice
double precision unit_length
Physical scaling factor for length.
logical saveprim
If true, convert from conservative to primitive variables in output.
double precision, dimension(ndim) qstretch_baselevel
stretch factor between cells at AMR level 1, per dimension
integer max_blocks
The maximum number of grid blocks in a processor.
double precision unit_density
Physical scaling factor for density.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
integer, parameter stretch_uni
Unidirectional stretching from a side.
logical, dimension(ndim) collapse
If collapse(DIM) is true, generate output integrated over DIM.
integer type_block
MPI type for block including ghost cells and its size.
double precision, dimension(:,:), allocatable dxmid
integer type_block_wc_io
MPI type for IO: cell corner (wc) or cell center (wcc) variables.
integer, parameter filecollapse_
Constant indicating collapsed output.
integer, dimension(:), allocatable type_gradient_limiter
Type of slope limiter used for computing gradients or divergences, when typegrad or typediv are set t...
double precision xprob
minimum and maximum domain boundaries for each dimension
double precision, dimension(:), allocatable dg
extent of grid blocks in domain per dimension, in array over levels
integer, dimension(:,:), allocatable node_sub
double precision, dimension(ndim, 2) writespshift
double precision time_convert_factor
Conversion factor for time unit.
logical need_global_cmax
need global maximal wave speed
character(len=std_len) convert_type
Which format to use when converting.
double precision, dimension(:,:), allocatable rnode_sub
integer, parameter unitcollapse
integer collapselevel
The level at which to produce line-integrated / collapsed output.
integer, parameter plevel_
integer ndir
Number of spatial dimensions (components) for vector variables.
character(len=std_len) collapse_type
character(len=std_len) typedimsplit
integer nstep
How many sub-steps the time integrator takes.
integer, dimension(:), allocatable type_limiter
Type of slope limiter used for reconstructing variables on cell edges.
double precision, dimension(:), allocatable amr_wavefilter
refinement: lohner estimate wavefilter setting
logical resume_previous_run
If true, restart a previous run from the latest snapshot.
double precision, dimension(:), allocatable entropycoef
double precision wall_time_max
Ending wall time (in hours) for the simulation.
integer it_init
initial iteration count
integer it_max
Stop the simulation after this many time steps have been taken.
double precision dtmin
Stop the simulation when the time step becomes smaller than this value.
character(len=std_len) base_filename
Base file name for simulation output, which will be followed by a number.
integer itag_stg
MPI recv send variables for staggered-variable AMR.
integer direction_slit
direction of slit (for spectrum)
character(len= *), parameter undefined
integer itag
MPI recv send variables for AMR.
character(len=std_len) usr_filename
User parameter file.
integer, parameter nsavehi
Maximum number of saves that can be defined by tsave or itsave.
integer istep
Index of the sub-step in a multi-step time integrator.
character(len=std_len) typeboundspeed
Which type of TVDLF method to use.
character(len=std_len) typeprolonglimit
Limiter used for prolongation to refined grids and ghost cells.
integer npe
The number of MPI tasks.
logical b0field
split magnetic field as background B0 field
logical, dimension(:), allocatable writelevel
integer, parameter ixglo
Lower index of grid block arrays (always 1)
double precision courantpar
The Courant (CFL) number used for the simulation.
character(len=std_len), dimension(:), allocatable flux_scheme
Which spatial discretization to use (per grid level)
integer, parameter fileanalysis_
Constant indicating analysis output (see Writing a custom analysis subroutine)
integer, dimension(3, 3, 3) lvc
Levi-Civita tensor.
integer, dimension(:,:), allocatable sendstatus
integer r_
Indices for cylindrical coordinates FOR TESTS, negative value when not used:
integer, parameter unitanalysis
double precision small_density
logical, dimension(ndim) stretched_dim
True if a dimension is stretched.
logical, dimension(ndim) aperiodb
True for dimensions with aperiodic boundaries.
double precision, dimension(nfile) tsavelast
logical convert
If true and restart_from_file is given, convert snapshots to other file formats.
integer, dimension(:), allocatable ng
number of grid blocks in domain per dimension, in array over levels
integer, dimension(ndim) stretch_type
What kind of stretching is used per dimension.
Module with basic data types used in amrvac.
integer, dimension(nfile) isavet
logical pass_wall_time
If true, wall time is up, modify snapshotnext for later overwrite.
integer, parameter rpxmax0_
integer nghostcells
Number of ghost cells surrounding a grid.
integer, dimension(:,:), allocatable sendstatus_stg
double precision dtdiffpar
For resistive MHD, the time step is also limited by the diffusion time: .
logical, dimension(:), allocatable loglimit
integer, parameter nodehi
grid hierarchy info (level and grid indices)
integer ixghi
Upper index of grid block arrays.
integer type_block_io
MPI type for IO: block excluding ghost cells.
character(len=std_len) typeghostfill
double precision xload
Stores the memory and load imbalance, used in printlog.
character(len=std_len) filename_euv
Base file name for synthetic EUV emission output.
integer ditregrid
Reconstruct the AMR grid once every ditregrid iteration(s)
integer, parameter unitconvert
character(len=std_len) geometry_name
logical angmomfix
Enable to strictly conserve the angular momentum (works both in cylindrical and spherical coordinates...
logical, dimension(:), allocatable small_values_fix_iw
Whether to apply small value fixes to certain variables.
logical, dimension(2, ndim) poleb
Indicates whether there is a pole at a boundary.
double precision small_pressure
character(len=std_len) typesourcesplit
How to apply dimensional splitting to the source terms, see disretization.md.
integer, parameter fileslice_
Constant indicating slice output.
logical stagger_grid
True for using stagger grid.
integer, dimension(:), allocatable recvrequest
logical save_physical_boundary
True for save physical boundary cells in dat files.
Module for physical and numeric constants Created: 01.09.2012 Oliver Porth (Physical constants) ...
Definition: mod_constants.t:3
character(len=std_len) resolution_euv
resolution of the output
logical nocartesian
IO switches for conversion.
integer ierrmpi
A global MPI error return code.
integer ixm
the mesh range (within a block with ghost cells)
integer snapshotnext
IO: snapshot and collapsed views output numbers/labels.
character(len=std_len), dimension(:, :), allocatable typeboundary
Array indicating the type of boundary condition per variable and per physical boundary.
integer, parameter ixgslo
Lower index of stagger grid block arrays (always 0)
double precision unit_velocity
Physical scaling factor for velocity.
integer, dimension(nfile) itsavelast
logical solve_internal_e
Solve polytropic process instead of solving total energy.
character(len=std_len), dimension(:), allocatable typeentropy
Which type of entropy fix to use with Riemann-type solvers.
double precision dtpar
If dtpar is positive, it sets the timestep dt, otherwise courantpar is used to limit the time step ba...
double precision, dimension(:), allocatable dt_grid
character(len=std_len) restart_from_file
If not &#39;unavailable&#39;, resume from snapshot with this base file name.
logical image
output image
double precision, dimension(:), allocatable w_convert_factor
Conversion factors the primitive variables.
character(len=std_len) typefilelog
Which type of log to write: &#39;normal&#39;, &#39;special&#39;, &#39;regression_test&#39;.
integer, parameter unitterm
Unit for standard output.
integer, parameter filelog_
Constant indicating log output.
integer, parameter unitpar
integer nbufferx
Number of cells as buffer zone.
integer, dimension(3, 3) kr
Kronecker delta tensor.
double precision, dimension(:,:), allocatable dx
character(len=std_len) typegrad
logical check_small_values
check and optionally fix unphysical small values (density, gas pressure)
double precision time_between_print
to monitor timeintegration loop at given wall-clock time intervals
integer, dimension(:), allocatable, parameter d
logical reset_it
If true, reset iteration count to 0.
character(len=40), dimension(nfile), parameter output_names
Names of the output methods.
integer mype
The rank of the current MPI task.
double precision global_time
The global simulation time.
integer, dimension(1:nfile) n_saves
Number of saved files of each type.
integer, dimension(:,:), allocatable nstretchedblocks
(even) number of (symmetrically) stretched blocks per level and dimension
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
logical reset_time
If true, reset iteration count and global_time to original values, and start writing snapshots at ind...
integer, parameter stretch_symm
Symmetric stretching around the center.
integer, parameter unitsnapshot
logical, dimension(nfile) save_file
whether or not to save an output file
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
double precision, dimension(nfile) dtsave
Repeatedly save output of type N when dtsave(N) simulation time has passed.
double precision unit_pressure
Physical scaling factor for pressure.
double precision unit_magneticfield
Physical scaling factor for magnetic field.
logical use_particles
Use particles module or not.
double precision unit_temperature
Physical scaling factor for temperature.
This module contains variables that describe the connectivity of the mesh and also data structures fo...
double precision, dimension(nfile) tsavestart
Start of read out (not counting specified read outs)
double precision, dimension(:,:), allocatable qstretch
Stretching factors and first cell size for each AMR level and dimension.
integer, parameter uniterr
Unit for error messages.
double precision xstretch
physical extent of stretched border in symmetric stretching
double precision, dimension(:,:), allocatable dxfirst_1mq
logical firstprocess
If true, call initonegrid_usr upon restarting.
double precision tfixgrid
Fix the AMR grid after this time.
double precision length_convert_factor
integer iprob
problem switch allowing different setups in same usr_mod.t
double precision, dimension(ndim) dxlevel
double precision, dimension(nsavehi, nfile) tsave
Save output of type N on times tsave(:, N)
integer, dimension(:), allocatable sendrequest_stg
integer, dimension(:,:), allocatable recvstatus
double precision bdip
amplitude of background dipolar, quadrupolar, octupolar, user&#39;s field
integer, dimension(nfile) isaveit
integer snapshotini
Resume from the snapshot with this index.
integer, parameter rpxmin0_
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer, dimension(2^d &) type_sub_block
logical stretch_uncentered
If true, adjust mod_geometry routines to account for grid stretching (but the flux computation will n...
logical si_unit
Use SI units (.true.) or use cgs units (.false.)
logical crash
Save a snapshot before crash a run met unphysical values.
logical use_multigrid
Use multigrid (only available in 2D and 3D)
integer, dimension(:), allocatable recvrequest_stg
logical slab
Cartesian geometry or not.
integer slowsteps
If > 1, then in the first slowsteps-1 time steps dt is reduced by a factor .
character(len=std_len) typepoly
integer, parameter stretch_none
No stretching.
integer nwauxio
Number of auxiliary variables that are only included in output.
logical, dimension(:), allocatable logflag
integer, dimension(:,:), allocatable recvstatus_stg
integer icomm
The MPI communicator.
double precision small_temperature
error handling
double precision cmax_global
global fastest wave speed needed in fd scheme and glm method
character(len=std_len) typediv
integer, parameter nlevelshi
The maximum number of levels in the grid refinement.
integer type_block_io_stg
MPI type for IO of staggered variables.
integer refine_max_level
Maximal number of AMR levels.
logical, dimension(ndim) periodb
True for dimensions with periodic boundaries.
integer, parameter nfile
Number of output methods.
integer log_fh
MPI file handle for logfile.
double precision location_slit
location of the slit
logical reset_grid
If true, rebuild the AMR grid upon restarting.
logical source_split_usr
Use split or unsplit way to add user&#39;s source terms, default: unsplit.
integer, dimension(:,:), allocatable node
integer, dimension(^nd, 2^d &) type_coarse_block_stg
MPI type for staggered block coarsened by 2, and for its children blocks.
integer, parameter unitstdin
Unit for standard input.
integer itfixgrid
Fix the AMR grid after this many time steps.
integer type_coarse_block
MPI type for block coarsened by 2, and for its children blocks.
double precision unit_time
Physical scaling factor for time.
integer, dimension(:), allocatable sendrequest
character(len=std_len) typecourant
How to compute the CFL-limited time step.
double precision, dimension(:,:), allocatable dxfirst
logical phys_energy
Solve energy equation or not.
integer it
Number of time steps taken.
double precision time_bc
accumulated wall-clock time spent on boundary conditions
integer, dimension(nfile) ditsave
Repeatedly save output of type N when ditsave(N) time steps have passed.
logical, dimension(:), allocatable phyboundblock
True if a block has any physical boundary.
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
double precision, dimension(:), allocatable derefine_ratio
integer, dimension(ndim) nstretchedblocks_baselevel
(even) number of (symmetrically) stretched blocks at AMR level 1, per dimension
double precision time_init
Start time for the simulation.
character(len=std_len) typeaverage
double precision unit_numberdensity
Physical scaling factor for number density.
integer ixgshi
Upper index of stagger grid block arrays.
integer type_block_xc_io
MPI type for IO: cell corner (xc) or cell center (xcc) coordinates.
integer, parameter fileout_
Constant indicating regular output.
integer, dimension(nsavehi, nfile) itsave
Save output of type N on iterations itsave(:, N)
logical spectrum
output spectrum
logical autoconvert
If true, already convert to output format during the run.
character(len=std_len) typecurl
pure subroutine cross_product(ixIL, ixOL, a, b, axb)
Cross product of two vectors.