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