13 integer,
intent(in) :: ixi^
l,ixo^
l, nwc
14 double precision,
intent(in) :: w(ixi^s, 1:nw)
15 double precision,
intent(in) :: x(ixi^s,1:
ndim)
16 double precision :: wnew(ixo^s, 1:nwc)
24 character(len=40) :: file_suffix
25 character(len=40) :: dataset_names(
max_nw)
41 character(len=*),
intent(in):: aux_variable_names
42 integer,
intent(in) :: nwc
43 character(len=name_len) :: names(1:nwc)
45 integer:: space_position,iw
46 character(len=name_len):: wname
47 character(len=std_len):: scanstring
51 scanstring=trim(adjustl(aux_variable_names))
52 space_position=index(scanstring,
' ')
54 do while (space_position==1)
55 scanstring=scanstring(2:)
56 space_position=index(scanstring,
' ')
58 wname=scanstring(:space_position-1)
59 scanstring=scanstring(space_position+1:)
60 space_position=index(scanstring,
' ')
62 names(iw)=trim(adjustl(wname))
70 integer,
intent(in) :: nwc
71 character(len=*),
intent(in):: aux_variable_names
72 character(len=*),
intent(in) :: file_suffix
75 function phys_convert_vars(ixI^L, ixO^L, w, x, nwc)
result(wnew)
77 integer,
intent(in) :: ixI^L, ixO^L, nwc
78 double precision,
intent(in) :: w(ixI^S, 1:nw)
79 double precision,
intent(in) :: x(ixI^S,1:ndim)
80 double precision :: wnew(ixO^S, 1:nwc)
81 end function phys_convert_vars
89 integer,
intent(in) :: nwc
90 character(len=*),
intent(in) :: dataset_names(:)
91 character(len=*),
intent(in) :: file_suffix
94 function phys_convert_vars(ixI^L, ixO^L, w, x, nwc)
result(wnew)
96 integer,
intent(in) :: ixI^L, ixO^L, nwc
97 double precision,
intent(in) :: w(ixI^S, 1:nw)
98 double precision,
intent(in) :: x(ixI^S,1:ndim)
99 double precision :: wnew(ixO^S, 1:nwc)
100 end function phys_convert_vars
106 call mpistop(
"INCREASE max_nw ")
110 temp%phys_convert_vars => phys_convert_vars
112 temp%file_suffix = file_suffix
113 temp%dataset_names(1:nwc) = dataset_names(1:nwc)
122 do while(
associated(temp))
123 call convert_dat_generic(temp%nwc, temp%dataset_names, temp%file_suffix, temp%phys_convert_vars)
136 integer :: file_handle, igrid, Morton_no, iwrite
137 integer :: ipe, ix_buffer(2*ndim+1), n_values
138 integer :: ixI^L, ixO^L, n_ghost(2*ndim)
139 integer :: ixOs^L,n_values_stagger
140 integer :: iorecvstatus(MPI_STATUS_SIZE)
141 integer :: ioastatus(MPI_STATUS_SIZE)
142 integer :: igrecvstatus(MPI_STATUS_SIZE)
143 integer :: istatus(MPI_STATUS_SIZE)
145 integer(kind=MPI_OFFSET_KIND) :: offset_tree_info
146 integer(kind=MPI_OFFSET_KIND) :: offset_block_data
147 integer(kind=MPI_OFFSET_KIND) :: offset_offsets
148 double precision,
allocatable :: w_buffer(:)
149 double precision,
allocatable:: converted_vars(:^D&,:)
151 integer :: idim, itag
152 integer,
allocatable :: block_ig(:, :)
153 integer,
allocatable :: block_lvl(:)
154 integer(kind=MPI_OFFSET_KIND),
allocatable :: block_offset(:)
156 integer,
intent(in) :: nwc
157 character(len=*),
intent(in) :: dataset_names(:)
158 character(len=*),
intent(in) :: file_suffix
161 function convert_vars(ixI^L, ixO^L,w,x,nwc)
result(wres)
163 integer,
intent(in) :: ixI^L, ixO^L, nwc
164 double precision,
intent(in) :: x(ixI^S,1:ndim)
165 double precision,
intent(in) :: w(ixI^S,1:nw)
166 double precision :: wres(ixO^S,1:nwc)
167 end function convert_vars
176 allocate(w_buffer(n_values))
179 allocate(block_ig(ndim,
nleafs))
180 allocate(block_lvl(
nleafs))
181 allocate(block_offset(
nleafs+1))
188 offset_tree_info = -1
189 offset_block_data = -1
191 offset_block_data, dataset_names, nwc)
193 call mpi_file_get_position(file_handle, offset_tree_info,
ierrmpi)
200 igrid =
sfc(1, morton_no)
201 ipe =
sfc(2, morton_no)
204 block_ig(:, morton_no) = [ pnode%ig^d ]
205 block_lvl(morton_no) = pnode%level
206 block_offset(morton_no) = 0
209 call mpi_file_write(file_handle, block_lvl,
size(block_lvl), &
212 call mpi_file_write(file_handle, block_ig,
size(block_ig), &
216 call mpi_file_get_position(file_handle, offset_offsets,
ierrmpi)
217 call mpi_file_write(file_handle, block_offset(1:
nleafs),
nleafs, &
220 call mpi_file_get_position(file_handle, offset_block_data,
ierrmpi)
223 if (offset_block_data - offset_tree_info /= &
225 nleafs * ((1+ndim) * size_int + 2 * size_int))
then
227 print *,
"Warning: MPI_OFFSET type /= 8 bytes"
228 print *,
"This *could* cause problems when reading .dat files"
232 block_offset(1) = offset_block_data
252 if(ps(igrid)%is_physical_boundary(2*idim-1)) n_ghost(idim)=
nghostcells
254 if(ps(igrid)%is_physical_boundary(2*idim)) n_ghost(ndim+idim)=
nghostcells
258 {ixomin^d = ixmlo^d - n_ghost(^d)\}
259 {ixomax^d = ixmhi^d + n_ghost(ndim+^d)\}
264 {iximin^d =
ixglo^d\}
265 {iximax^d =
ixghi^d\}
267 w_buffer(1:n_values) = pack(convert_vars(ixi^l, ixo^l,ps(igrid)%w(ixi^s, 1:nw), ps(igrid)%x(ixi^s, 1:ndim),nwc), .true.)
269 ix_buffer(1) = n_values
270 ix_buffer(2:) = n_ghost
273 call mpi_send(ix_buffer, 2*ndim+1, &
275 call mpi_send(w_buffer, n_values, &
279 call mpi_file_write(file_handle, ix_buffer(2:), &
280 2*ndim, mpi_integer, istatus,
ierrmpi)
281 call mpi_file_write(file_handle, w_buffer, &
282 n_values, mpi_double_precision, istatus,
ierrmpi)
285 block_offset(iwrite+1) = block_offset(iwrite) + &
286 int(n_values, mpi_offset_kind) * size_double + &
298 call mpi_recv(ix_buffer, 2*ndim+1, mpi_integer, ipe, itag,
icomm,&
300 n_values = ix_buffer(1)
302 call mpi_recv(w_buffer, n_values, mpi_double_precision,&
305 call mpi_file_write(file_handle, ix_buffer(2:), &
306 2*ndim, mpi_integer, istatus,
ierrmpi)
307 call mpi_file_write(file_handle, w_buffer, &
308 n_values, mpi_double_precision, istatus,
ierrmpi)
311 block_offset(iwrite+1) = block_offset(iwrite) + &
312 int(n_values, mpi_offset_kind) * size_double + &
318 call mpi_file_seek(file_handle, offset_offsets, mpi_seek_set,
ierrmpi)
319 call mpi_file_write(file_handle, block_offset(1:
nleafs),
nleafs, &
323 call mpi_file_seek(file_handle, 0_mpi_offset_kind, mpi_seek_set,
ierrmpi)
325 offset_block_data, dataset_names, nwc)
327 call mpi_file_close(file_handle,
ierrmpi)
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
subroutine write_forest(file_handle)
subroutine convert_dat_generic(nwc, dataset_names, file_suffix, convert_vars)
subroutine add_convert_method2(phys_convert_vars, nwc, aux_variable_names, file_suffix)
type(convert_vars_method), pointer head_convert_vars_methods
character(len=name_len) function, dimension(1:nwc) get_names_from_string(aux_variable_names, nwc)
subroutine init_convert()
subroutine add_convert_method(phys_convert_vars, nwc, dataset_names, file_suffix)
Module with basic grid data structures.
integer, dimension(:), allocatable, save sfc_to_igrid
Go from a Morton number to an igrid index (for a single processor)
integer, dimension(:), allocatable, save morton_start
First Morton number per processor.
integer, save nleafs
Number of leaf block.
integer, dimension(:), allocatable, save morton_stop
Last Morton number per processor.
integer, dimension(:,:), allocatable, save sfc
Array to go from a Morton number to an igrid and processor index. Sfc(1:3, MN) contains [igrid,...
integer, save nparents
Number of parent blocks.
type(tree_node_ptr), dimension(:,:), allocatable, save igrid_to_node
Array to go from an [igrid, ipe] index to a node pointer.
This module contains definitions of global parameters and variables and some generic functions/subrou...
type(state), pointer block
Block pointer for using one block and its previous state.
integer ixghi
Upper index of grid block arrays.
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical save_physical_boundary
True for save physical boundary cells in dat files.
integer icomm
The MPI communicator.
integer mype
The rank of the current MPI task.
integer ierrmpi
A global MPI error return code.
integer snapshotnext
IO: snapshot and collapsed views output numbers/labels.
integer npe
The number of MPI tasks.
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer nghostcells
Number of ghost cells surrounding a grid.
double precision, dimension(ndim) dxlevel
integer, parameter ixglo
Lower index of grid block arrays (always 1)
integer, parameter max_nw
Maximum number of variables.
The data structure that contains information about a tree node/grid block.