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