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  character(len=*), parameter :: undefined = 'undefined'
18 
19  !> @todo Move mpi related variables to e.g. mod_comm
20 
21  !> The number of MPI tasks
22  integer :: npe
23 
24  !> The rank of the current MPI task
25  integer :: mype
26 
27  !> The MPI communicator
28  integer :: icomm
29 
30  !> A global MPI error return code
31  !> @todo Make local
32  integer :: ierrmpi
33 
34  !> MPI file handle for logfile
35  integer :: log_fh
36  !> MPI type for block including ghost cells and its size
37  integer :: type_block, size_block
38  !> MPI type for block coarsened by 2, and for its children blocks
40  !> MPI type for staggered block coarsened by 2, and for its children blocks
41  integer :: type_coarse_block_stg(^nd,2^d&), type_sub_block_stg(^ND,2^D&)
42  !> MPI type for IO: block excluding ghost cells
44  !> MPI type for IO of staggered variables
46  !> MPI type for IO: cell corner (xc) or cell center (xcc) coordinates
48  !> MPI type for IO: cell corner (wc) or cell center (wcc) variables
50 
51 
52  ! geometry and domain setups
53 
54  !> the mesh range (within a block with ghost cells)
55  integer :: ixm^ll
56 
57  !> minimum and maximum domain boundaries for each dimension
58  double precision :: xprob^l
59 
60  !> Indices for cylindrical coordinates FOR TESTS, negative value when not used:
61  integer :: r_ = -1
62  integer :: phi_ = -1
63  integer :: z_ = -1
64 
65  !> Number of spatial dimensions for grid variables
66  integer, parameter :: ndim=^nd
67 
68  !> Number of spatial dimensions (components) for vector variables
69  integer :: ndir=ndim
70 
71  !> Cartesian geometry or not
72  logical :: slab
73 
74  !> uniform Cartesian geometry or not (stretched Cartesian)
75  logical :: slab_uniform
76 
77  !> number of grid blocks in domain per dimension, in array over levels
78  integer, dimension(:), allocatable :: ng^d
79  !> extent of grid blocks in domain per dimension, in array over levels
80  double precision, dimension(:), allocatable :: dg^d
81 
82  !> number of cells for each dimension in level-one mesh
83  integer :: domain_nx^d
84 
85  !> number of cells for each dimension in grid block excluding ghostcells
86  integer :: block_nx^d
87 
88  !> Lower index of grid block arrays (always 1)
89  integer, parameter :: {ixglo^d = 1|, }
90 
91  !> Upper index of grid block arrays
92  integer :: ixghi^d
93 
94  !> Lower index of stagger grid block arrays (always 0)
95  integer, parameter :: {ixgslo^d = 0|, }
96 
97  !> Upper index of stagger grid block arrays
98  integer :: ixgshi^d
99 
100  !> Number of ghost cells surrounding a grid
101  integer :: nghostcells
102 
103  integer, parameter :: stretch_none = 0 !< No stretching
104  integer, parameter :: stretch_uni = 1 !< Unidirectional stretching from a side
105  integer, parameter :: stretch_symm = 2 !< Symmetric stretching around the center
106 
107  !> If true, adjust mod_geometry routines to account for grid stretching (but
108  !> the flux computation will not)
110  !> True if a dimension is stretched
111  logical :: stretched_dim(ndim)
112  !> What kind of stretching is used per dimension
113  integer :: stretch_type(ndim)
114  !> stretch factor between cells at AMR level 1, per dimension
115  double precision :: qstretch_baselevel(ndim)
116  !> (even) number of (symmetrically) stretched
117  !> blocks at AMR level 1, per dimension
119  !> (even) number of (symmetrically) stretched blocks per level and dimension
120  integer, allocatable :: nstretchedblocks(:,:)
121  !> physical extent of stretched border in symmetric stretching
122  double precision :: xstretch^d
123  !> Stretching factors and first cell size for each AMR level and dimension
124  double precision, allocatable :: qstretch(:,:), dxfirst(:,:), &
125  dxfirst_1mq(:,:), dxmid(:,:)
126 
127  !> grid hierarchy info (level and grid indices)
128  integer, parameter :: nodehi=^nd+1
129  integer, parameter :: plevel_=1
130  integer, parameter :: pig^d_=plevel_+^d
131 
132  integer, allocatable :: node(:,:)
133  integer, allocatable :: node_sub(:,:)
134 
135  !> grid location info (corner coordinates and grid spacing)
136  integer, parameter :: rnodehi=3*^nd
137  integer, parameter :: rpxmin0_=0
138  integer, parameter :: rpxmin^d_=rpxmin0_+^d
139  integer, parameter :: rpxmax0_=^nd
140  integer, parameter :: rpxmax^d_=rpxmax0_+^d
141  integer, parameter :: rpdx^d_=2*^nd+^d
142 
143  !> Corner coordinates
144  double precision, allocatable :: rnode(:,:)
145  double precision, allocatable :: rnode_sub(:,:)
146 
147  double precision, allocatable :: dx(:,:)
148  double precision :: dxlevel(ndim)
149 
150  ! IO related quantities
151 
152  !> Maximum number of saves that can be defined by tsave or itsave
153  integer, parameter :: nsavehi=100
154 
155  !> Number of output methods
156  integer, parameter :: nfile = 5
157 
158  !> index number of the latest existing data file
159  integer :: index_latest_data
160 
161  !> Names of the output methods
162  character(len=40), parameter :: output_names(nfile) = &
163  ['log ', 'normal ', 'slice ', 'collapsed', 'analysis ']
164 
165  !> User parameter file
166  character(len=std_len) :: usr_filename
167 
168  !> If collapse(DIM) is true, generate output integrated over DIM
169  logical :: collapse(ndim)
170 
171  !> Save output of type N on times tsave(:, N)
172  double precision :: tsave(nsavehi,nfile)
173 
174  !> \todo Move tsavelast to amrvac.t
175  double precision :: tsavelast(nfile)
176 
177  !> Repeatedly save output of type N when dtsave(N) simulation time has passed
178  double precision :: dtsave(nfile)
179 
180  !> Save output of type N on iterations itsave(:, N)
181  integer :: itsave(nsavehi,nfile)
182 
183  !> \todo remove itsavelast?
184  integer :: itsavelast(nfile)
185 
186  !> Repeatedly save output of type N when ditsave(N) time steps have passed
187  integer :: ditsave(nfile)
188 
189  !> \todo Move to amrvac.t
190  integer :: isavet(nfile)
191 
192  !> \todo Move to amrvac.t
193  integer :: isaveit(nfile)
194 
195  !> Start of read out (not counting specified read outs)
196  double precision :: tsavestart(nfile)
197 
198  !> The level at which to produce line-integrated / collapsed output
199  integer :: collapselevel
200 
201  !> Number of saved files of each type
202  !> \todo Move to mod_input_output
203  integer :: n_saves(1:nfile)
204 
205  !> whether or not to save an output file
206  logical :: save_file(nfile)
207 
208  !> to monitor timeintegration loop at given wall-clock time intervals
209  double precision :: time_between_print
210 
211  !> accumulated wall-clock time spent on boundary conditions
212  double precision :: time_bc
213 
214  !> IO: snapshot and collapsed views output numbers/labels
216 
217  !> Constant indicating log output
218  integer, parameter :: filelog_ = 1
219 
220  !> Constant indicating regular output
221  integer, parameter :: fileout_ = 2
222 
223  !> Constant indicating slice output
224  integer, parameter :: fileslice_ = 3
225 
226  !> Constant indicating collapsed output
227  integer, parameter :: filecollapse_ = 4
228 
229  !> Constant indicating analysis output (see @ref analysis.md)
230  integer, parameter :: fileanalysis_ = 5
231 
232  !> Unit for standard input
233  integer, parameter :: unitstdin=5
234 
235  !> Unit for standard output
236  integer, parameter :: unitterm=6
237 
238  !> Unit for error messages
239  integer, parameter :: uniterr=6
240 
241  !> file handle for IO
242  integer, parameter :: unitpar=9
243  integer, parameter :: unitconvert=10
244  integer, parameter :: unitslice=11
245  integer, parameter :: unitsnapshot=12
246  integer, parameter :: unitcollapse=13
247  integer, parameter :: unitanalysis=14
248 
249  !> Number of auxiliary variables that are only included in output
250  integer :: nwauxio
251 
252  !> IO switches for conversion
253  logical :: nocartesian
254  logical, allocatable :: w_write(:)
255  logical, allocatable :: writelevel(:)
256  double precision :: writespshift(ndim,2)
258 
259  !> Which par files are used as input
260  character(len=std_len), allocatable :: par_files(:)
261 
262  !> Base file name for simulation output, which will be followed by a number
263  character(len=std_len) :: base_filename
264 
265  !> If not 'unavailable', resume from snapshot with this base file name
266  character(len=std_len) :: restart_from_file
267 
268  !> Which type of log to write: 'normal', 'special', 'regression_test'
269  character(len=std_len) :: typefilelog
270 
271  !> Resume from the snapshot with this index
272  integer :: snapshotini
273 
274  !> If true, restart a previous run from the latest snapshot
276 
277  !> If true and restart_from_file is given, convert snapshots to
278  !> other file formats
279  logical :: convert
280 
281  !> If true, already convert to output format during the run
282  logical :: autoconvert
283 
284  !> If true, convert from conservative to primitive variables in output
285  logical :: saveprim
286 
287  !> Which format to use when converting
288  !>
289  !> Options are: tecplot, tecplotCC, vtu, vtuCC, vtuB, vtuBCC,
290  !> tecplotmpi, tecplotCCmpi, vtumpi, vtuCCmpi, vtuBmpi, vtuBCCmpi, pvtumpi, pvtuCCmpi,
291  !> pvtuBmpi, pvtuBCCmpi, tecline, teclinempi, onegrid
292  character(len=std_len) :: convert_type
293 
294  character(len=std_len) :: collapse_type
295 
296  !> Conversion factors the primitive variables
297  double precision, allocatable :: w_convert_factor(:)
298 
299  double precision :: length_convert_factor
300 
301  !> Conversion factor for time unit
302  double precision :: time_convert_factor
303 
304  !> Stores the memory and load imbalance, used in printlog
305  double precision :: xload, xmemory
306 
307  !> Save a snapshot before crash a run met unphysical values
308  logical :: crash=.false.
309 
310  ! Physics factors
311 
312  !> Physical scaling factor for length
313  double precision :: unit_length=1.d0
314 
315  !> Physical scaling factor for time
316  double precision :: unit_time=1.d0
317 
318  !> Physical scaling factor for density
319  double precision :: unit_density=1.d0
320 
321  !> Physical scaling factor for velocity
322  double precision :: unit_velocity=1.d0
323 
324  !> Physical scaling factor for temperature
325  double precision :: unit_temperature=1.d0
326 
327  !> Physical scaling factor for pressure
328  double precision :: unit_pressure=1.d0
329 
330  !> Physical scaling factor for magnetic field
331  double precision :: unit_magneticfield=1.d0
332 
333  !> Physical scaling factor for number density
334  double precision :: unit_numberdensity=1.d0
335 
336  !> Physical scaling factor for charge
337  double precision :: unit_charge=1.d0
338 
339  !> Physical scaling factor for mass
340  double precision :: unit_mass=1.d0
341 
342  !> Normalised speed of light
343  double precision :: c_norm=1.d0
344 
345  !> error handling
347 
348  !> amplitude of background dipolar, quadrupolar, octupolar, user's field
349  double precision :: bdip=0.d0
350  double precision :: bquad=0.d0
351  double precision :: boct=0.d0
352  double precision :: busr=0.d0
353 
354  !> check and optionally fix unphysical small values (density, gas pressure)
355  logical :: check_small_values=.true.
356  logical :: fix_small_values=.false.
357 
358  !> split magnetic field as background B0 field
359  logical :: b0field=.false.
360 
361  !> Use SI units (.true.) or use cgs units (.false.)
362  logical :: si_unit=.false.
363 
364  !> Use TRAC (Johnston 2019 ApJL, 873, L22) for MHD or 1D HD
365  logical :: phys_trac=.false.
366  integer :: phys_trac_type=1
367  double precision :: phys_trac_mask
368 
369  !> Enable to strictly conserve the angular momentum
370  !> (works both in cylindrical and spherical coordinates)
371  logical :: angmomfix=.false.
372 
373  !> Use particles module or not
374  logical :: use_particles=.false.
375 
376  !> Use multigrid (only available in 2D and 3D)
377  logical :: use_multigrid = .false.
378 
379  ! AMR switches
380 
381  !> The maximum number of grid blocks in a processor
382  integer :: max_blocks
383 
384  !> The maximum number of levels in the grid refinement
385  integer, parameter :: nlevelshi = 20
386 
387  !> Maximal number of AMR levels
388  integer :: refine_max_level
389 
390  !> Weights of variables used to calculate error for mesh refinement
391  double precision, allocatable :: w_refine_weight(:)
392 
393  !> Fix the AMR grid after this time
394  double precision :: tfixgrid
395 
396  !> Whether to apply flux conservation at refinement boundaries
397  logical :: fix_conserve_global = .true.
398 
399  !> Fix the AMR grid after this many time steps
400  integer :: itfixgrid
401 
402  !> Reconstruct the AMR grid once every ditregrid iteration(s)
403  integer :: ditregrid
404 
405  !> refinement: lohner estimate wavefilter setting
406  double precision, allocatable :: amr_wavefilter(:)
407 
408  integer :: refine_criterion
409  logical :: prolongprimitive=.false.
410  logical :: coarsenprimitive=.false.
411 
412  !> Error tolerance for refinement decision
413  double precision, allocatable :: refine_threshold(:)
414  double precision, allocatable :: derefine_ratio(:)
415 
416  !> If true, rebuild the AMR grid upon restarting
417  logical :: reset_grid
418  !> True for using stagger grid
419  logical :: stagger_grid=.false.
420  !> True for record electric field
421  logical :: record_electric_field=.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 
446  logical :: time_advance
447 
448  !> The Courant (CFL) number used for the simulation
449  double precision :: courantpar
450 
451  !> How to compute the CFL-limited time step.
452  !>
453  !> Options are 'maxsum': max(sum(c/dx)); 'summax': sum(max(c/dx)) and
454  !> 'minimum: max(c/dx), where the summations loop over the grid dimensions and
455  !> c is the velocity. The default 'maxsum' is the conventiontal way of
456  !> computing CFL-limited time steps.
457  character(len=std_len) :: typecourant
458 
459  !> If dtpar is positive, it sets the timestep dt, otherwise courantpar is used
460  !> to limit the time step based on the Courant condition.
461  double precision :: dtpar
462 
463  !> For resistive MHD, the time step is also limited by the diffusion time:
464  !> \f$ dt < dtdiffpar \times dx^2/eta \f$
465  double precision :: dtdiffpar
466 
467  !> The global simulation time
468  double precision :: global_time
469 
470  !> Start time for the simulation
471  double precision :: time_init
472 
473  !> End time for the simulation
474  double precision :: time_max
475 
476  !> Ending wall time (in hours) for the simulation
477  double precision :: wall_time_max
478 
479  !> Stop the simulation when the time step becomes smaller than this value
480  double precision :: dtmin
481 
482  !> Force timeloop exit when final dt < dtmin
483  logical :: final_dt_exit
484 
485  !> If true, reset iteration count and global_time to original values, and
486  !> start writing snapshots at index 0
487  logical :: reset_time
488 
489  !> If true, reset iteration count to 0
490  logical :: reset_it
491 
492  !> If true, allow final dt reduction for matching time_max on output
494 
495  !> If true, call initonegrid_usr upon restarting
496  logical :: firstprocess
497 
498  !> If true, wall time is up, modify snapshotnext for later overwrite
499  logical :: pass_wall_time
500 
501  !> If true, do H-correction to fix the carbuncle problem at grid-aligned shocks
502  logical :: h_correction=.false.
503 
504  !> Number of time steps taken
505  integer :: it
506 
507  !> Stop the simulation after this many time steps have been taken
508  integer :: it_max
509 
510  !> initial iteration count
511  integer :: it_init
512 
513  !> If > 1, then in the first slowsteps-1 time steps dt is reduced
514  !> by a factor \f$ 1 - (1- step/slowsteps)^2 \f$
515  integer :: slowsteps
516 
517  ! Method switches
518 
519  !> Index of the sub-step in a multi-step time integrator
520  integer :: istep
521 
522  !> How many sub-steps the time integrator takes
523  integer :: nstep
524 
525  !> Which flux scheme of spatial discretization to use (per grid level)
526  integer, allocatable :: flux_method(:)
527 
528  !> The spatial discretization for the predictor step when using a two
529  !> step PC method
530  integer, allocatable :: typepred1(:)
531 
532  !> flux schemes
533  integer, parameter :: fs_hll=1
534  integer, parameter :: fs_hllc=2
535  integer, parameter :: fs_hlld=3
536  integer, parameter :: fs_hllcd=4
537  integer, parameter :: fs_tvdlf=5
538  integer, parameter :: fs_tvdmu=6
539  integer, parameter :: fs_tvd=7
540  integer, parameter :: fs_hancock=8
541  integer, parameter :: fs_cd=9
542  integer, parameter :: fs_cd4=10
543  integer, parameter :: fs_fd=11
544  integer, parameter :: fs_source=12
545  integer, parameter :: fs_nul=13
546 
547  !> time stepper type
548  integer :: t_stepper=0
549  integer, parameter :: onestep=1
550  integer, parameter :: twostep=2
551  integer, parameter :: threestep=3
552  integer, parameter :: fourstep=4
553  integer, parameter :: fivestep=5
554 
555  !> time integrator method
556  integer :: t_integrator=0
557  integer, parameter :: forward_euler=1
558  integer, parameter :: predictor_corrector=2
559  integer, parameter :: ssprk3=3
560  integer, parameter :: ssprk4=4
561  integer, parameter :: ssprk5=5
562 
563  integer, parameter :: imex_euler=6
564  integer, parameter :: imex_sp=7
565  integer, parameter :: rk2_alf=8
566  integer, parameter :: ssprk2=9
567  integer, parameter :: imex_midpoint=10
568  integer, parameter :: imex_trapezoidal=11
569  integer, parameter :: imex_222=12
570 
571  integer, parameter :: rk3_bt=13
572  integer, parameter :: imex_ars3=14
573  integer, parameter :: imex_232=15
574  integer, parameter :: imex_cb3a=16
575 
576  integer, parameter :: rk4=17
577  integer, parameter :: jameson=18
578 
579  !> Type of slope limiter used for reconstructing variables on cell edges
580  integer, allocatable :: type_limiter(:)
581 
582  !> Type of slope limiter used for computing gradients or divergences, when
583  !> typegrad or typediv are set to 'limited'
584  integer, allocatable :: type_gradient_limiter(:)
585 
586  !> background magnetic field location indicator
587  integer :: b0i=0
588 
589  !> Limiter used for prolongation to refined grids and ghost cells
590  integer :: prolong_limiter=0
591 
592  !> Which type of entropy fix to use with Riemann-type solvers
593  character(len=std_len), allocatable :: typeentropy(:)
594 
595  !> Which type of TVD method to use
596  character(len=std_len) :: typetvd
597 
598  !> bound (left/min and right.max) speed of Riemann fan
599  integer :: boundspeed
600 
601  character(len=std_len) :: typeaverage
602  character(len=std_len) :: typedimsplit
603  character(len=std_len) :: geometry_name='default'
604  character(len=std_len) :: typepoly
605 
606  integer :: nxdiffusehllc
607  double precision, allocatable :: entropycoef(:)
608  double precision :: tvdlfeps
609  logical, allocatable :: loglimit(:), logflag(:)
610  logical :: flathllc,flatcd,flatsh
611  !> Use split or unsplit way to add user's source terms, default: unsplit
612  logical :: source_split_usr
613  logical :: dimsplit
614 
615  !> RK2(alfa) method parameters from Butcher tableau
616  double precision :: rk_a21,rk_b1,rk_b2
617  !> IMEX-222(lambda) one-parameter family of schemes
618  double precision :: imex222_lambda
619  !> SSPRK choice of methods (both threestep and fourstep, Shu-Osher 2N* implementation)
620  !> also fivestep SSPRK54
621  integer :: ssprk_order
625  !> RK3 Butcher table
626  integer :: rk3_switch
628  !> IMEX_ARS3 parameter ars_gamma
629  double precision :: ars_gamma
630  !> IMEX_232 choice and parameters
631  integer :: imex_switch
633  double precision :: imex_b3,imex_c2,imex_c3
634  !> IMEX_CB3a extra parameters
635  double precision :: imex_a22, imex_a33, imex_ha32
636  !> whether IMEX in use or not
637  logical :: use_imex_scheme
638 
639  character(len=std_len) :: typediv,typegrad
640 
641  !> global fastest wave speed needed in fd scheme and glm method
642  double precision :: cmax_global
643 
644  !> global fastest flow speed needed in glm method
645  double precision :: vmax_global
646 
647  !> global largest a2 for schmid scheme
648  double precision :: a2max_global(ndim)
649 
650  !> need global maximal wave speed
651  logical :: need_global_cmax=.false.
652 
653  !> global value for schmid scheme
654  logical :: need_global_a2max=.false.
655 
656  ! Boundary region parameters
657 
658  !> True for dimensions with periodic boundaries
659  logical :: periodb(ndim)
660 
661  !> Indicates whether there is a pole at a boundary
662  logical :: poleb(2,ndim)
663 
664  !> True for dimensions with aperiodic boundaries
665  logical :: aperiodb(ndim)
666 
667  !> True for save physical boundary cells in dat files
669 
670  !> True if a block has any physical boundary
671  logical, allocatable :: phyboundblock(:)
672 
673  !> Array indicating the type of boundary condition per variable and per
674  !> physical boundary
675  integer, allocatable :: typeboundary(:, :)
676  !> boundary condition types
677  integer, parameter :: bc_special=1
678  integer, parameter :: bc_cont=2
679  integer, parameter :: bc_symm=3
680  integer, parameter :: bc_asymm=4
681  integer, parameter :: bc_periodic=5
682  integer, parameter :: bc_aperiodic=6
683  integer, parameter :: bc_noinflow=7
684  integer, parameter :: bc_data=8
685  integer, parameter :: bc_character=9
686 
687  !> whether copy values instead of interpolation in ghost cells of finer blocks
688  logical :: ghost_copy=.false.
689 
690  !> if there is an internal boundary
691  logical :: internalboundary
692 
693  !> Base file name for synthetic EUV emission output
694  character(len=std_len) :: filename_euv
695  !> output EUV image
696  logical :: image_euv=.false.
697  !> wavelength for output
698  integer :: wavelength
699  !> resolution of the EUV image
700  character(len=std_len) :: resolution_euv='instrument'
701  !> Base file name for synthetic SXR emission output
702  character(len=std_len) :: filename_sxr
703  !> output SXR image
704  logical :: image_sxr=.false.
705  ! minimum and maximum energy of SXR (keV)
706  integer :: emin_sxr,emax_sxr
707  !> resolution of the SXR image
708  character(len=std_len) :: resolution_sxr='instrument'
709  !> direction of the line of sight (LOS)
710  integer :: los_theta,los_phi
711  !> rotation of image
712  integer :: image_rotate=0
713  !> Base file name for synthetic EUV spectrum output
714  character(len=std_len) :: filename_spectrum
715  !> output EUV spectral profile
716  logical :: spectrum_euv=.false.
717  !> wave length for spectrum
718  integer :: spectrum_wl
719  !> spectral window
721  !> location of the slit
722  double precision :: location_slit=0.d0
723  !> direction of the slit
724  integer :: direction_slit
725  !> resolution of the spectrum
726  character(len=std_len) :: resolution_spectrum='instrument'
727 
728  !> Block pointer for using one block and its previous state
729  type(state), pointer :: block
730 
731  !$OMP THREADPRIVATE(block,dxlevel,b0i)
732 
733 contains
734 
735  !> Cross product of two vectors
736  pure subroutine cross_product(ixI^L,ixO^L,a,b,axb)
737  integer, intent(in) :: ixi^l, ixo^l
738  double precision, intent(in) :: a(ixi^s,3), b(ixi^s,3)
739  double precision, intent(out) :: axb(ixi^s,3)
740 
741  axb(ixo^s,1)=a(ixo^s,2)*b(ixo^s,3)-a(ixo^s,3)*b(ixo^s,2)
742  axb(ixo^s,2)=a(ixo^s,3)*b(ixo^s,1)-a(ixo^s,1)*b(ixo^s,3)
743  axb(ixo^s,3)=a(ixo^s,1)*b(ixo^s,2)-a(ixo^s,2)*b(ixo^s,1)
744  end subroutine cross_product
745 
746 end module mod_global_parameters
Module with basic data types used in amrvac.
This module contains variables that describe the connectivity of the mesh and also data structures fo...
Module for physical and numeric constants Created: 01.09.2012 Oliver Porth (Physical constants)
Definition: mod_constants.t:3
This module contains definitions of global parameters and variables and some generic functions/subrou...
character(len=std_len), dimension(:), allocatable typeentropy
Which type of entropy fix to use with Riemann-type solvers.
double precision, dimension(nfile) tsavelast
double precision, dimension(:), allocatable w_convert_factor
Conversion factors the primitive variables.
type(state), pointer block
Block pointer for using one block and its previous state.
double precision xload
Stores the memory and load imbalance, used in printlog.
integer, parameter unitslice
integer nstep
How many sub-steps the time integrator takes.
logical h_correction
If true, do H-correction to fix the carbuncle problem at grid-aligned shocks.
integer it_max
Stop the simulation after this many time steps have been taken.
logical, dimension(ndim) aperiodb
True for dimensions with aperiodic boundaries.
logical internalboundary
if there is an internal boundary
character(len=std_len) filename_sxr
Base file name for synthetic SXR emission output.
integer spectrum_wl
wave length for spectrum
logical nocartesian
IO switches for conversion.
integer, dimension(:), allocatable typepred1
The spatial discretization for the predictor step when using a two step PC method.
double precision dtdiffpar
For resistive MHD, the time step is also limited by the diffusion time: .
character(len=std_len) typegrad
integer ixgshi
Upper index of stagger grid block arrays.
logical reset_it
If true, reset iteration count to 0.
double precision unit_charge
Physical scaling factor for charge.
integer, parameter bc_noinflow
integer type_coarse_block
MPI type for block coarsened by 2, and for its children blocks.
integer, parameter fs_tvdlf
double precision small_pressure
integer ixghi
Upper index of grid block arrays.
integer, parameter stretch_uni
Unidirectional stretching from a side.
character(len=std_len) geometry_name
integer, dimension(3, 3, 3) lvc
Levi-Civita tensor.
double precision unit_time
Physical scaling factor for time.
double precision unit_density
Physical scaling factor for density.
character(len=std_len) typepoly
integer, parameter fs_hlld
logical source_split_usr
Use split or unsplit way to add user's source terms, default: unsplit.
integer, parameter imex_euler
integer domain_nx
number of cells for each dimension in level-one mesh
integer, parameter unitpar
file handle for IO
character(len=std_len) filename_spectrum
Base file name for synthetic EUV spectrum output.
logical, dimension(nfile) save_file
whether or not to save an output file
integer, parameter bc_asymm
logical resume_previous_run
If true, restart a previous run from the latest snapshot.
double precision global_time
The global simulation time.
integer type_block_xc_io
MPI type for IO: cell corner (xc) or cell center (xcc) coordinates.
integer, dimension(nsavehi, nfile) itsave
Save output of type N on iterations itsave(:, N)
double precision unit_mass
Physical scaling factor for mass.
logical use_imex_scheme
whether IMEX in use or not
integer istep
Index of the sub-step in a multi-step time integrator.
logical spectrum_euv
output EUV spectral profile
double precision time_max
End time for the simulation.
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, parameter threestep
double precision xstretch
physical extent of stretched border in symmetric stretching
logical, dimension(:), allocatable logflag
double precision time_init
Start time for the simulation.
logical stretch_uncentered
If true, adjust mod_geometry routines to account for grid stretching (but the flux computation will n...
double precision phys_trac_mask
logical firstprocess
If true, call initonegrid_usr upon restarting.
integer, parameter imex_232
integer snapshotini
Resume from the snapshot with this index.
double precision small_temperature
error handling
double precision xprob
minimum and maximum domain boundaries for each dimension
integer it
Number of time steps taken.
character(len=std_len) filename_euv
Base file name for synthetic EUV emission output.
logical, dimension(:), allocatable loglimit
double precision, dimension(:), allocatable dg
extent of grid blocks in domain per dimension, in array over levels
integer, parameter fivestep
integer, parameter bc_character
integer it_init
initial iteration count
integer, dimension(:, :), allocatable typeboundary
Array indicating the type of boundary condition per variable and per physical boundary.
integer ditregrid
Reconstruct the AMR grid once every ditregrid iteration(s)
logical saveprim
If true, convert from conservative to primitive variables in output.
double precision ars_gamma
IMEX_ARS3 parameter ars_gamma.
double precision unit_numberdensity
Physical scaling factor for number density.
integer, parameter fs_hllcd
character(len=std_len) convert_type
Which format to use when converting.
double precision unit_pressure
Physical scaling factor for pressure.
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer itfixgrid
Fix the AMR grid after this many time steps.
integer, parameter filecollapse_
Constant indicating collapsed output.
logical angmomfix
Enable to strictly conserve the angular momentum (works both in cylindrical and spherical coordinates...
integer prolong_limiter
Limiter used for prolongation to refined grids and ghost cells.
double precision, dimension(:), allocatable amr_wavefilter
refinement: lohner estimate wavefilter setting
double precision unit_length
Physical scaling factor for length.
integer, parameter nlevelshi
The maximum number of levels in the grid refinement.
double precision location_slit
location of the slit
logical save_physical_boundary
True for save physical boundary cells in dat files.
double precision vmax_global
global fastest flow speed needed in glm method
logical stagger_grid
True for using stagger grid.
double precision cmax_global
global fastest wave speed needed in fd scheme and glm method
double precision time_convert_factor
Conversion factor for time unit.
logical, dimension(:), allocatable phyboundblock
True if a block has any physical boundary.
integer, dimension(:,:), allocatable nstretchedblocks
(even) number of (symmetrically) stretched blocks per level and dimension
integer, parameter imex_ars3
integer, dimension(^nd, 2^d &) type_coarse_block_stg
MPI type for staggered block coarsened by 2, and for its children blocks.
integer, parameter bc_data
logical use_particles
Use particles module or not.
integer, parameter fs_tvdmu
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
integer, parameter jameson
integer icomm
The MPI communicator.
integer, dimension(:), allocatable ng
number of grid blocks in domain per dimension, in array over levels
logical reset_time
If true, reset iteration count and global_time to original values, and start writing snapshots at ind...
double precision bdip
amplitude of background dipolar, quadrupolar, octupolar, user's field
integer b0i
background magnetic field location indicator
integer, parameter imex_trapezoidal
integer, parameter nsavehi
Maximum number of saves that can be defined by tsave or itsave.
logical image_sxr
output SXR image
integer mype
The rank of the current MPI task.
double precision dtpar
If dtpar is positive, it sets the timestep dt, otherwise courantpar is used to limit the time step ba...
integer, dimension(1:nfile) n_saves
Number of saved files of each type.
character(len=std_len) typediv
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
integer image_rotate
rotation of image
integer type_block_io
MPI type for IO: block excluding ghost cells.
double precision, dimension(nfile) tsavestart
Start of read out (not counting specified read outs)
integer, dimension(nfile) ditsave
Repeatedly save output of type N when ditsave(N) time steps have passed.
integer, parameter plevel_
integer, dimension(2^d &) type_sub_block
logical, dimension(ndim) collapse
If collapse(DIM) is true, generate output integrated over DIM.
integer, parameter unitstdin
Unit for standard input.
integer, dimension(:), allocatable, parameter d
character(len=std_len) usr_filename
User parameter file.
logical need_global_a2max
global value for schmid scheme
logical ghost_copy
whether copy values instead of interpolation in ghost cells of finer blocks
double precision length_convert_factor
integer, parameter nodehi
grid hierarchy info (level and grid indices)
integer ndir
Number of spatial dimensions (components) for vector variables.
integer, parameter uniterr
Unit for error messages.
double precision imex222_lambda
IMEX-222(lambda) one-parameter family of schemes.
double precision courantpar
The Courant (CFL) number used for the simulation.
double precision wall_time_max
Ending wall time (in hours) for the simulation.
integer ixm
the mesh range (within a block with ghost cells)
integer ierrmpi
A global MPI error return code.
logical autoconvert
If true, already convert to output format during the run.
integer, dimension(:), allocatable flux_method
Which flux scheme of spatial discretization to use (per grid level)
character(len=std_len) collapse_type
logical slab
Cartesian geometry or not.
integer slowsteps
If > 1, then in the first slowsteps-1 time steps dt is reduced by a factor .
integer ssprk_order
SSPRK choice of methods (both threestep and fourstep, Shu-Osher 2N* implementation) also fivestep SSP...
double precision, dimension(:,:), allocatable qstretch
Stretching factors and first cell size for each AMR level and dimension.
integer, parameter bc_periodic
integer type_block_wc_io
MPI type for IO: cell corner (wc) or cell center (wcc) variables.
character(len=std_len) resolution_sxr
resolution of the SXR image
integer, parameter bc_special
boundary condition types
integer snapshotnext
IO: snapshot and collapsed views output numbers/labels.
double precision unit_magneticfield
Physical scaling factor for magnetic field.
integer, parameter unitanalysis
integer, parameter rk2_alf
logical, dimension(ndim) stretched_dim
True if a dimension is stretched.
integer, dimension(ndim) nstretchedblocks_baselevel
(even) number of (symmetrically) stretched blocks at AMR level 1, per dimension
integer npe
The number of MPI tasks.
integer nwauxio
Number of auxiliary variables that are only included in output.
double precision unit_velocity
Physical scaling factor for velocity.
integer index_latest_data
index number of the latest existing data file
integer, dimension(nfile) itsavelast
integer, parameter stretch_none
No stretching.
double precision, dimension(ndim) qstretch_baselevel
stretch factor between cells at AMR level 1, per dimension
double precision, dimension(ndim, 2) writespshift
integer imex_switch
IMEX_232 choice and parameters.
integer, parameter fs_hll
flux schemes
double precision time_between_print
to monitor timeintegration loop at given wall-clock time intervals
integer, parameter unitterm
Unit for standard output.
double precision, dimension(nfile) dtsave
Repeatedly save output of type N when dtsave(N) simulation time has passed.
logical, dimension(:), allocatable w_write
character(len=std_len) typecourant
How to compute the CFL-limited time step.
integer iprob
problem switch allowing different setups in same usr_mod.t
character(len=std_len) restart_from_file
If not 'unavailable', resume from snapshot with this base file name.
logical, dimension(ndim) periodb
True for dimensions with periodic boundaries.
double precision c_norm
Normalised speed of light.
logical b0field
split magnetic field as background B0 field
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer, parameter filelog_
Constant indicating log output.
character(len=40), dimension(nfile), parameter output_names
Names of the output methods.
double precision unit_temperature
Physical scaling factor for temperature.
integer, dimension(:), allocatable type_limiter
Type of slope limiter used for reconstructing variables on cell edges.
integer, parameter imex_sp
logical final_dt_reduction
If true, allow final dt reduction for matching time_max on output.
integer, parameter bc_cont
logical si_unit
Use SI units (.true.) or use cgs units (.false.)
logical, dimension(:), allocatable writelevel
integer, parameter unitcollapse
integer, parameter fileout_
Constant indicating regular output.
character(len=std_len) typetvd
Which type of TVD method to use.
integer nbufferx
Number of cells as buffer zone.
double precision, dimension(:), allocatable entropycoef
integer, parameter twostep
double precision time_bc
accumulated wall-clock time spent on boundary conditions
double precision, dimension(:,:), allocatable dx
integer, parameter rpxmax0_
integer type_block
MPI type for block including ghost cells and its size.
integer, parameter onestep
double precision rk_a21
RK2(alfa) method parameters from Butcher tableau.
integer, parameter predictor_corrector
double precision tfixgrid
Fix the AMR grid after this time.
integer, parameter unitsnapshot
pure subroutine cross_product(ixIL, ixOL, a, b, axb)
Cross product of two vectors.
integer nghostcells
Number of ghost cells surrounding a grid.
double precision, dimension(:), allocatable w_refine_weight
Weights of variables used to calculate error for mesh refinement.
integer, parameter bc_symm
integer, parameter forward_euler
logical phys_trac
Use TRAC (Johnston 2019 ApJL, 873, L22) for MHD or 1D HD.
logical fix_conserve_global
Whether to apply flux conservation at refinement boundaries.
character(len=std_len) typedimsplit
character(len=std_len) typeaverage
character(len= *), parameter undefined
logical image_euv
output EUV image
double precision, dimension(nsavehi, nfile) tsave
Save output of type N on times tsave(:, N)
double precision spectrum_window_max
logical need_global_cmax
need global maximal wave speed
logical convert
If true and restart_from_file is given, convert snapshots to other file formats.
character(len=std_len) resolution_euv
resolution of the EUV image
integer wavelength
wavelength for output
integer collapselevel
The level at which to produce line-integrated / collapsed output.
logical reset_grid
If true, rebuild the AMR grid upon restarting.
integer, dimension(ndim) stretch_type
What kind of stretching is used per dimension.
logical crash
Save a snapshot before crash a run met unphysical values.
integer t_stepper
time stepper type
double precision, dimension(:,:), allocatable rnode_sub
logical use_multigrid
Use multigrid (only available in 2D and 3D)
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
integer, parameter fourstep
character(len=std_len) resolution_spectrum
resolution of the spectrum
double precision, dimension(:), allocatable refine_threshold
Error tolerance for refinement decision.
character(len=std_len) base_filename
Base file name for simulation output, which will be followed by a number.
integer, parameter rnodehi
grid location info (corner coordinates and grid spacing)
double precision small_density
double precision spectrum_window_min
spectral window
double precision dtmin
Stop the simulation when the time step becomes smaller than this value.
integer refine_max_level
Maximal number of AMR levels.
integer, parameter fs_hllc
integer, parameter fileslice_
Constant indicating slice output.
double precision, dimension(:), allocatable derefine_ratio
integer, parameter nfile
Number of output methods.
integer max_blocks
The maximum number of grid blocks in a processor.
integer, parameter stretch_symm
Symmetric stretching around the center.
integer, parameter fileanalysis_
Constant indicating analysis output (see Writing a custom analysis subroutine)
integer r_
Indices for cylindrical coordinates FOR TESTS, negative value when not used:
integer rk3_switch
RK3 Butcher table.
integer, parameter rpxmin0_
integer, parameter bc_aperiodic
character(len=std_len) typefilelog
Which type of log to write: 'normal', 'special', 'regression_test'.
double precision imex_a22
IMEX_CB3a extra parameters.
integer direction_slit
direction of the slit
integer type_block_io_stg
MPI type for IO of staggered variables.
logical, dimension(2, ndim) poleb
Indicates whether there is a pole at a boundary.
integer boundspeed
bound (left/min and right.max) speed of Riemann fan
logical pass_wall_time
If true, wall time is up, modify snapshotnext for later overwrite.
integer, parameter imex_cb3a
integer, parameter ixgslo
Lower index of stagger grid block arrays (always 0)
logical record_electric_field
True for record electric field.
integer, parameter imex_222
integer, parameter unitconvert
double precision, dimension(:,:), allocatable dxfirst
integer t_integrator
time integrator method
integer, parameter fs_hancock
logical final_dt_exit
Force timeloop exit when final dt < dtmin.
integer, dimension(nfile) isaveit
integer, dimension(:), allocatable type_gradient_limiter
Type of slope limiter used for computing gradients or divergences, when typegrad or typediv are set t...
integer, dimension(:,:), allocatable node
integer, dimension(:,:), allocatable node_sub
integer, parameter fs_source
double precision, dimension(ndim) dxlevel
integer, dimension(nfile) isavet
integer, parameter imex_midpoint
double precision, dimension(:,:), allocatable dxfirst_1mq
double precision, dimension(ndim) a2max_global
global largest a2 for schmid scheme
integer los_theta
direction of the line of sight (LOS)
logical check_small_values
check and optionally fix unphysical small values (density, gas pressure)
double precision, dimension(:,:), allocatable dxmid
integer, parameter ixglo
Lower index of grid block arrays (always 1)
integer log_fh
MPI file handle for logfile.