26#if defined(__parallel)
28 USE mpi_f08,
ONLY: mpi_allreduce, mpi_integer, mpi_bxor, mpi_allgather
30 USE mpi,
ONLY: mpi_allreduce, mpi_integer, mpi_bxor, mpi_allgather
35 openpmd_attributable_type, &
36 openpmd_dynamic_memory_view_type_1d, &
37 openpmd_dynamic_memory_view_type_3d, &
39 openpmd_particle_species_type, &
40 openpmd_record_component_type, &
41 openpmd_record_type, openpmd_type_double, openpmd_type_int
53#include "../base/base_uses.f90"
62 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'realspace_grid_openpmd'
63 LOGICAL,
PARAMETER,
PRIVATE :: debug_this_module = .false.
64 LOGICAL,
PRIVATE :: parses_linebreaks = .false., &
67 TYPE cp_openpmd_write_buffer_1d
68 REAL(KIND=
dp),
POINTER :: buffer(:)
69 END TYPE cp_openpmd_write_buffer_1d
82 SUBROUTINE pw_get_atom_types(particles_z, res_atom_types, res_atom_counts, res_len)
83 INTEGER,
DIMENSION(:),
INTENT(IN) :: particles_z
84 INTEGER,
ALLOCATABLE,
DIMENSION(:),
INTENT(OUT) :: res_atom_types, res_atom_counts
85 INTEGER,
INTENT(OUT) :: res_len
87 INTEGER :: current_atom_number, i
88 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: particles_z_sorted
91 ALLOCATE (particles_z_sorted(
SIZE(particles_z)))
92 particles_z_sorted(:) = particles_z(:)
95 ALLOCATE (res_atom_types(min(118,
SIZE(particles_z))))
96 ALLOCATE (res_atom_counts(min(118,
SIZE(particles_z))))
97 current_atom_number = -1
99 DO i = 1,
SIZE(particles_z_sorted)
100 IF (particles_z_sorted(i) /= current_atom_number)
THEN
101 res_len = res_len + 1
102 current_atom_number = particles_z_sorted(i)
103 res_atom_types(res_len) = current_atom_number
104 res_atom_counts(res_len) = 1
106 res_atom_counts(res_len) = res_atom_counts(res_len) + 1
110 END SUBROUTINE pw_get_atom_types
122 SUBROUTINE pw_write_particle_species( &
131 INTEGER,
DIMENSION(:),
INTENT(IN) :: particles_z
132 REAL(KIND=
dp),
DIMENSION(:, :),
INTENT(IN) :: particles_r
133 REAL(KIND=
dp),
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_zeff
134 INTEGER,
INTENT(IN) :: atom_type, atom_count
135 TYPE(cp_openpmd_per_call_value_type) :: openpmd_data
136 LOGICAL :: do_write_data
138 CHARACTER(len=1),
DIMENSION(3),
PARAMETER :: dims = [
"x",
"y",
"z"]
140 CHARACTER(len=3) :: atom_type_as_string
141 CHARACTER(len=default_string_length) :: species_name
143 INTEGER,
DIMENSION(1) :: global_extent, global_offset, &
145 TYPE(cp_openpmd_write_buffer_1d) :: charge_write_buffer
146 TYPE(cp_openpmd_write_buffer_1d),
DIMENSION(3) :: write_buffers
147 TYPE(openpmd_attributable_type) :: attr
148 TYPE(openpmd_dynamic_memory_view_type_1d) :: unresolved_charge_write_buffer
149 TYPE(openpmd_dynamic_memory_view_type_1d), &
150 DIMENSION(3) :: unresolved_write_buffers
151 TYPE(openpmd_particle_species_type) :: species
152 TYPE(openpmd_record_component_type) :: charge_component, position_component, &
153 position_offset_component
154 TYPE(openpmd_record_type) :: charge, position, position_offset
159 global_extent(1) = atom_count
160 IF (do_write_data)
THEN
162 local_extent(1) = atom_count
168 WRITE (atom_type_as_string,
'(I3)') atom_type
169 species_name = trim(openpmd_data%name_prefix)//
"-"//adjustl(atom_type_as_string)
171 CALL openpmd_data%iteration%open()
172 species = openpmd_data%iteration%get_particle_species(trim(species_name))
174 position_offset = species%get_record(
"positionOffset")
175 position = species%get_record(
"position")
177 CALL position%set_unit_dimension([1.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp])
179 position_offset_component = position_offset%get_component(dims(k))
180 CALL position_offset_component%make_constant_zero(openpmd_type_int, global_extent)
181 CALL position_offset_component%set_unit_SI(
a_bohr)
182 position_component = position%get_component(dims(k))
183 CALL position_component%reset_dataset(openpmd_type_double, global_extent)
184 CALL position_component%set_unit_SI(
a_bohr)
185 unresolved_write_buffers(k) = &
186 position_component%store_chunk_span_1d_double(global_offset, local_extent)
187 write_buffers(k)%buffer => unresolved_write_buffers(k)%resolve_double(deallocate=.false.)
190 IF (
PRESENT(particles_zeff))
THEN
191 charge = species%get_record(
"charge")
192 charge_component = charge%as_record_component()
194 CALL charge%set_unit_dimension([0.0_dp, 0.0_dp, 1.0_dp, 1.0_dp, 0.0_dp, 0.0_dp, 0.0_dp])
195 CALL charge_component%reset_dataset(openpmd_type_double, global_extent)
196 CALL charge_component%set_unit_SI(
e_charge)
197 unresolved_charge_write_buffer = charge_component%store_chunk_span_1d_double(global_offset, local_extent)
198 charge_write_buffer%buffer => unresolved_charge_write_buffer%resolve_double(deallocate=.false.)
203 write_buffers(k)%buffer = unresolved_write_buffers(k)%resolve_double(deallocate=.true.)
205 IF (
PRESENT(particles_zeff))
THEN
206 charge_write_buffer%buffer = unresolved_charge_write_buffer%resolve_double(deallocate=.true.)
208 IF (do_write_data)
THEN
210 DO i = 1,
SIZE(particles_z)
211 IF (particles_z(i) == atom_type)
THEN
213 write_buffers(k)%buffer(j) = particles_r(k, i)
215 IF (
PRESENT(particles_zeff))
THEN
216 charge_write_buffer%buffer(j) = particles_zeff(i)
222 attr = openpmd_data%iteration%as_attributable()
223 CALL attr%series_flush(
"hdf5.independent_stores = true")
224 END SUBROUTINE pw_write_particle_species
237 SUBROUTINE pw_write_particles( &
247 INTEGER,
DIMENSION(:),
INTENT(IN) :: particles_z
248 REAL(KIND=
dp),
DIMENSION(:, :),
INTENT(IN) :: particles_r
249 REAL(KIND=
dp),
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_zeff
250 INTEGER,
DIMENSION(:),
INTENT(IN) :: atom_types, atom_counts
251 INTEGER,
INTENT(IN),
TARGET :: num_atom_types
252 TYPE(cp_openpmd_per_call_value_type) :: openpmd_data
253 TYPE(mp_comm_type),
OPTIONAL :: gid
255 INTEGER :: i, mpi_rank
256 LOGICAL :: do_write_data
258 IF (
PRESENT(gid))
THEN
259 CALL gid%get_rank(mpi_rank)
260 do_write_data = mpi_rank == 0
262 do_write_data = .true.
264 DO i = 1, num_atom_types
265 CALL pw_write_particle_species( &
275 END SUBROUTINE pw_write_particles
302 TYPE(pw_r3d_rs_type),
INTENT(IN) :: pw
304 CHARACTER(*),
INTENT(IN),
OPTIONAL :: title
305 REAL(KIND=
dp),
DIMENSION(:, :),
INTENT(IN), &
306 OPTIONAL :: particles_r
307 INTEGER,
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_z
308 REAL(KIND=
dp),
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_zeff
309 INTEGER,
DIMENSION(:),
OPTIONAL,
POINTER :: stride
310 LOGICAL,
INTENT(IN),
OPTIONAL :: zero_tails, silent, mpi_io
312 CHARACTER(len=*),
PARAMETER :: routineN =
'pw_to_openpmd'
314 CHARACTER(LEN=default_string_length) :: my_title
315 INTEGER :: count1, count2, count3, handle, i, I1, &
316 I2, I3, iat, L1, L2, L3, my_rank, &
317 my_stride(3), np, num_atom_types, &
319 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: atom_counts, atom_types
320 INTEGER,
DIMENSION(3) :: global_extent, local_extent, offset
321 LOGICAL :: be_silent, my_zero_tails, parallel_write
322 REAL(KIND=
dp),
DIMENSION(3) :: grid_spacing
323 REAL(KIND=
dp),
POINTER :: write_buffer(:, :, :)
324 TYPE(cp_openpmd_per_call_value_type) :: openpmd_data
325 TYPE(mp_comm_type) :: gid
326 TYPE(openpmd_attributable_type) :: attr
327 TYPE(openpmd_dynamic_memory_view_type_3d) :: unresolved_write_buffer
328 TYPE(openpmd_mesh_type) :: mesh
329 TYPE(openpmd_record_component_type) :: scalar_mesh
331 CALL timeset(routinen, handle)
333 my_zero_tails = .false.
335 parallel_write = .false.
336 gid = pw%pw_grid%para%group
337 IF (
PRESENT(title)) my_title = trim(title)
338 IF (
PRESENT(zero_tails)) my_zero_tails = zero_tails
339 IF (
PRESENT(silent)) be_silent = silent
340 IF (
PRESENT(mpi_io)) parallel_write = mpi_io
342 IF (
PRESENT(stride))
THEN
343 IF (
SIZE(stride) /= 1 .AND.
SIZE(stride) /= 3)
THEN
344 CALL cp_abort(__location__,
"STRIDE keyword can accept only 1 "// &
345 "(the same for X,Y,Z) or 3 values. Correct your input file.")
347 IF (
SIZE(stride) == 1)
THEN
349 my_stride(i) = stride(1)
352 my_stride = stride(1:3)
354 cpassert(my_stride(1) > 0)
355 cpassert(my_stride(2) > 0)
356 cpassert(my_stride(3) > 0)
361 cpassert(
PRESENT(particles_z) .EQV.
PRESENT(particles_r))
363 IF (
PRESENT(particles_z))
THEN
364 CALL pw_get_atom_types(particles_z,
atom_types, atom_counts, num_atom_types)
365 cpassert(
SIZE(particles_z) ==
SIZE(particles_r, dim=2))
366 np =
SIZE(particles_z)
374 grid_spacing(i) = sqrt(sum(pw%pw_grid%dh(:, i)**2))*real(my_stride(i),
dp)
377 IF (
PRESENT(particles_z))
THEN
378 IF (parallel_write)
THEN
379 CALL pw_write_particles( &
390 CALL pw_write_particles( &
403 global_extent(iat) = (pw%pw_grid%npts(iat) + my_stride(iat) - 1)/my_stride(iat)
405 offset(iat) = ((pw%pw_grid%bounds_local(1, iat) - pw%pw_grid%bounds(1, iat) + my_stride(iat) - 1)/my_stride(iat))
408 local_extent(iat) = ((pw%pw_grid%bounds_local(2, iat) + 1 - pw%pw_grid%bounds(1, iat) + my_stride(iat) - 1)/my_stride(iat))
410 local_extent = local_extent - offset
412 mesh = openpmd_data%iteration%get_mesh(trim(openpmd_data%name_prefix))
413 CALL mesh%set_axis_labels([
"x",
"y",
"z"])
414 CALL mesh%set_position([0.5_dp, 0.5_dp, 0.5_dp])
415 CALL mesh%set_grid_global_offset([ &
416 pw%pw_grid%bounds(1, 1)*grid_spacing(1), &
417 pw%pw_grid%bounds(1, 2)*grid_spacing(2), &
418 pw%pw_grid%bounds(1, 3)*grid_spacing(3)])
419 CALL mesh%set_grid_spacing(grid_spacing)
420 CALL mesh%set_grid_unit_SI(
a_bohr)
421 CALL mesh%set_unit_dimension(openpmd_data%unit_dimension)
422 scalar_mesh = mesh%as_record_component()
423 CALL scalar_mesh%set_unit_SI(openpmd_data%unit_si)
424 CALL scalar_mesh%reset_dataset(openpmd_type_double, global_extent)
431 l1 = pw%pw_grid%bounds(1, 1) + offset(1)*my_stride(1)
432 l2 = pw%pw_grid%bounds_local(1, 2)
433 l3 = pw%pw_grid%bounds_local(1, 3)
436 u1 = pw%pw_grid%bounds(1, 1) + (offset(1) + local_extent(1) - 1)*my_stride(1)
437 u2 = pw%pw_grid%bounds_local(2, 2)
438 u3 = pw%pw_grid%bounds_local(2, 3)
440 my_rank = pw%pw_grid%para%group%mepos
441 num_pe = pw%pw_grid%para%group%num_pe
443 IF (all(my_stride == 1))
THEN
444 CALL scalar_mesh%store_chunk(pw%array(l1:u1, l2:u2, l3:u3), offset)
446 attr = openpmd_data%iteration%as_attributable()
447 CALL attr%series_flush(
"hdf5.independent_stores = false")
450 DO i3 = l3, u3, my_stride(3)
454 unresolved_write_buffer = scalar_mesh%store_chunk_span_3d_double( &
455 [offset(1), offset(2), offset(3) + count3], &
456 [local_extent(1), local_extent(2), 1])
457 write_buffer => unresolved_write_buffer%resolve_double(deallocate=.true.)
460 cpassert(
ASSOCIATED(write_buffer))
461 cpassert(
SIZE(write_buffer, 1) == local_extent(1))
462 cpassert(
SIZE(write_buffer, 2) == local_extent(2))
463 cpassert(
SIZE(write_buffer, 3) == 1)
466 DO i2 = l2, u2, my_stride(2)
470 DO i1 = l1, u1, my_stride(1)
471 write_buffer(count1 + 1, count2 + 1, 1) = pw%array(i1, i2, i3)
482 CALL timestop(handle)
515 CHARACTER(*),
INTENT(IN),
OPTIONAL :: title
516 REAL(kind=
dp),
DIMENSION(:, :),
INTENT(IN), &
517 OPTIONAL :: particles_r
518 INTEGER,
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_z
519 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN),
OPTIONAL :: particles_zeff
520 INTEGER,
DIMENSION(:),
OPTIONAL,
POINTER :: stride
521 LOGICAL,
INTENT(IN),
OPTIONAL :: zero_tails, silent, mpi_io
526 mark_used(particles_r)
527 mark_used(particles_z)
528 mark_used(particles_zeff)
530 mark_used(zero_tails)
533 cpabort(
"CP2K compiled without the openPMD-api")
Define the atom type and its sub types.
Utility routines to open and close files. Tracking of preconnections.
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
type(cp_openpmd_per_call_value_type) function, public cp_openpmd_get_value_unit_nr(key)
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
Interface to the message passing library MPI.
type(mp_file_descriptor_type) function, public mp_file_type_hindexed_make_chv(count, lengths, displs)
Creates an indexed MPI type for arrays of strings using bytes for spacing (hindexed type)
subroutine, public mp_file_type_free(type_descriptor)
Releases the type used for MPI I/O.
integer, parameter, public mpi_character_size
integer, parameter, public file_offset
integer, parameter, public file_amode_rdonly
subroutine, public mp_file_type_set_view_chv(fh, offset, type_descriptor)
Uses a previously created indexed MPI character type to tell the MPI processes how to partition (set_...
Definition of physical constants:
real(kind=dp), parameter, public a_bohr
real(kind=dp), parameter, public e_charge
real(kind=dp), parameter, public seconds
integer, parameter, public pw_mode_local
Generate Gaussian cube files.
subroutine, public pw_to_openpmd(pw, unit_nr, title, particles_r, particles_z, particles_zeff, stride, zero_tails, silent, mpi_io)
...
All kind of helpful little routines.