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