MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_physicaldata.t
Go to the documentation of this file.
2 implicit none
3 save
4
5 type state
6 !> Variables, normally cell center conservative values
7 double precision, dimension(:^D&,:), allocatable :: w
8 !> Variables, cell face values
9 double precision, dimension(:^D&,:), allocatable :: ws
10 !> Variables, cell edge values
11 double precision, dimension(:^D&,:), allocatable :: we
12 !> Variables, cell corner values
13 double precision, dimension(:^D&,:), allocatable :: wc
14 !> extra variables do not need ghost cell and equation flux
15 double precision, dimension(:^D&,:), pointer :: wextra=>null()
16 !> Time-independent magnetic field at cell center and cell interface
17 double precision, dimension(:^D&,:,:), pointer :: b0=>null()
18 !> Time-independent electric current density at cell center
19 double precision, dimension(:^D&,:), pointer :: j0=>null()
20 !> Time-independent equi vars (B0 is not into this array)
21 double precision, dimension(:^D&,:,:), pointer :: equi_vars=>null()
22 !> Cell-center positions
23 double precision, dimension(:^D&,:), pointer :: x=>null()
24 !> Cell sizes in coordinate units
25 double precision, dimension(:^D&,:), pointer :: dx=>null()
26 !> Cell local timesteps
27 double precision, dimension(:^D&), pointer :: dt=>null()
28 !> Cell sizes at cell center in length unit
29 double precision, dimension(:^D&,:), pointer :: ds=>null()
30 !> Cell sizes at cell face in length unit
31 double precision, dimension(:^D&,:), pointer :: dsc=>null()
32 !> Volumes of a cell
33 double precision, dimension(:^D&), pointer :: dvolume=>null()
34 !> Areas of cell-center surfaces
35 double precision, dimension(:^D&,:), pointer :: surface=>null()
36 !> Areas of cell-face surfaces
37 double precision, dimension(:^D&,:), pointer :: surfacec=>null()
38 !> special values for a block
39 double precision, dimension(:), pointer :: special_values=>null()
40 !> ID of a grid block
41 integer :: igrid=-1
42 !> index range of block array in cell centers
43 integer :: ixg^l
44 !> index range of block array in cell faces
45 integer :: ixgs^l
46 !> level of AMR
47 integer :: level
48 !> If it face a physical boundary
49 logical, dimension(:), pointer :: is_physical_boundary(:) =>null()
50 end type state
51
52{^nooned
54 !> Variables, normally center
55 double precision, dimension(:^DE&,:), allocatable :: w
56 !> Variables for the cornerpositions on the slice
57 double precision, dimension(:^DE&,:), allocatable :: wc
58 !> Variables, normally center, one level coarser representative
59 double precision, dimension(:^DE&,:), allocatable :: wcoarse
60 !> Cell-center positions
61 double precision, dimension(:^DE&,:), allocatable :: x
62 !> Corner positions on the slice
63 double precision, dimension(:^DE&,:), allocatable :: xc
64 !> Cell-center positions, one level coarser representative
65 double precision, dimension(:^DE&,:), allocatable :: xcoarse
66 !> Cell sizes
67 double precision, dimension(:^DE&,:), allocatable :: dx
68 !> Cell sizes, one level coarser
69 double precision, dimension(:^D&,:), allocatable :: dxcoarse
70 !> Cell sizes in length unit
71 double precision, dimension(:^D&,:), allocatable :: ds
72 !> Volumes of a cell
73 double precision, dimension(:^DE&), allocatable :: dvolume
74 !> Volumes of a cell, one level coarser representative
75 double precision, dimension(:^DE&), allocatable :: dvolumecoarse
76 !> Areas of cell-center surfaces
77 double precision, dimension(:^DE&,:), allocatable :: surface
78 !> Areas of cell-face surfaces
79 double precision, dimension(:^DE&,:), allocatable :: surfacec
80 !> ID of a grid block
81 integer :: igrid=-1
82 end type state_sub
83}
84{^ifoned
85 type state_sub
86 !> Variables, normally center
87 double precision, dimension(:), allocatable :: w
88 !> Variables for the cornerpositions on the slice
89 double precision, dimension(:), allocatable :: wc
90 !> Variables, normally center, one level coarser representative
91 double precision, dimension(:), allocatable :: wcoarse
92 !> Cell-center positions
93 double precision, dimension(:), allocatable :: x
94 !> Corner positions on the slice
95 double precision, dimension(:), allocatable :: xc
96 !> Cell-center positions, one level coarser representative
97 double precision, dimension(:), allocatable :: xcoarse
98 !> ID of a grid block
99 integer :: igrid=-1
100 end type state_sub
101}
103 !> Variables new state
104 double precision, dimension(:^D&,:), allocatable :: w
105 !> Variables old state
106 double precision, dimension(:^D&,:), allocatable :: wold
107 end type grid_field
108 !> buffer for pole boundary
110
111 !> array of physical states for all blocks on my processor
112 type(state), dimension(:), allocatable, target :: ps
113 !> array of physical states, temp 1 for multi-step time integrator
114 type(state), dimension(:), allocatable, target :: ps1
115 !> array of physical states, temp 2 for multi-step time integrator
116 type(state), dimension(:), allocatable, target :: ps2
117 !> array of physical states, temp 3 for multi-step time integrator
118 type(state), dimension(:), allocatable, target :: ps3
119 !> array of physical states, temp 4 for multi-step time integrator
120 type(state), dimension(:), allocatable, target :: ps4
121 !> array of physical blocks, one level coarser representative
122 type(state), dimension(:), allocatable, target :: psc
123
124 !> array of physical blocks in reduced dimension
125 type(state_sub), dimension(:), allocatable, target :: ps_sub
126
127{^ifoned
128 double precision, dimension(:), allocatable :: collapseddata
129}
130{^nooned
131 double precision, dimension(:^DE&,:), allocatable :: collapseddata
132}
133 !> array of physical blocks of meshed fields for particles
134 type(grid_field), dimension(:), allocatable, target :: gridvars
135
136 !> velocities store for constrained transport
138 double precision, dimension(:^D&,:), allocatable :: vnorm,cbarmin,cbarmax
139 double precision, dimension(:^D&,:,:), allocatable :: vbarc,vbarlc,vbarrc
140 end type ct_velocity
141
142end module mod_physicaldata
type(state), dimension(:), allocatable, target ps4
array of physical states, temp 4 for multi-step time integrator
type(state_sub), dimension(:), allocatable, target ps_sub
array of physical blocks in reduced dimension
type(grid_field), dimension(:), allocatable, target gridvars
array of physical blocks of meshed fields for particles
type(state), dimension(:), allocatable, target ps2
array of physical states, temp 2 for multi-step time integrator
type(state), dimension(:), allocatable, target ps
array of physical states for all blocks on my processor
type(state), dimension(:), allocatable, target psc
array of physical blocks, one level coarser representative
type(state), dimension(:), allocatable, target ps1
array of physical states, temp 1 for multi-step time integrator
type(state) pole_buf
buffer for pole boundary
double precision, dimension(:^de &,:), allocatable collapseddata
type(state), dimension(:), allocatable, target ps3
array of physical states, temp 3 for multi-step time integrator
velocities store for constrained transport