openPMD output: Implement missing scientific metadata (#5096)

Co-authored-by: AI Agent <ai-agent@localhost>
Co-authored-by: Blablador <noreply@hzdr.de>
This commit is contained in:
Franz Pöschel 2026-04-27 16:37:16 +02:00 committed by GitHub
parent 90d8299cb6
commit e7dba46940
No known key found for this signature in database
GPG key ID: B5690EEEBB952194
5 changed files with 514 additions and 62 deletions

View file

@ -41,13 +41,16 @@ MODULE cp_output_handling_openpmd
section_vals_type, &
section_vals_val_get
USE kinds, ONLY: default_path_length, &
default_string_length
default_string_length, &
dp
USE machine, ONLY: m_mov
USE memory_utilities, ONLY: reallocate
USE message_passing, ONLY: mp_file_delete, &
mp_file_get_amode, &
mp_file_type
#ifdef __OPENPMD
USE physcon, ONLY: seconds
USE cp2k_info, ONLY: cp2k_version
USE openpmd_api, ONLY: &
openpmd_access_create, &
openpmd_attributable_type, openpmd_iteration_type, openpmd_mesh_type, &
@ -73,11 +76,13 @@ MODULE cp_output_handling_openpmd
#ifdef __OPENPMD
TYPE :: cp_openpmd_per_call_value_type
TYPE(openpmd_series_type) :: series = openpmd_series_type()
TYPE(openpmd_iteration_type) :: iteration = openpmd_iteration_type()
! TYPE(openpmd_mesh_type) :: mesh = openpmd_mesh_type()
TYPE(openpmd_series_type) :: series = openpmd_series_type()
TYPE(openpmd_iteration_type) :: iteration = openpmd_iteration_type()
! TYPE(openpmd_mesh_type) :: mesh = openpmd_mesh_type()
! TYPE(openpmd_particle_species_type) :: particle_species = openpmd_particle_species_type()
CHARACTER(len=default_string_length) :: name_prefix = "" ! e.g. 'WFN_00008_1'
CHARACTER(len=default_string_length) :: name_prefix = "" ! e.g. 'WFN_00008_1'
REAL(kind=dp), DIMENSION(7) :: unit_dimension = [0, 0, 0, 0, 0, 0, 0]
REAL(kind=dp) :: unit_si = 1
END TYPE cp_openpmd_per_call_value_type
TYPE :: cp_openpmd_per_call_type
@ -401,18 +406,30 @@ CONTAINS
!> \param series ...
!> \param middle_name ...
!> \param logger ...
!> \param sim_time Current simulation time in atomic units (used for timeOffset)
!> \return ...
! **************************************************************************************************
FUNCTION cp_openpmd_create_unit_nr_entry(openpmd_file_index, middle_name, logger) RESULT(res)
FUNCTION cp_openpmd_create_unit_nr_entry( &
openpmd_file_index, &
middle_name, &
logger, &
openpmd_unit_dimension, &
openpmd_unit_si, &
sim_time) RESULT(res)
INTEGER :: openpmd_file_index
CHARACTER(len=*), INTENT(IN) :: middle_name
TYPE(cp_logger_type), POINTER :: logger
REAL(kind=dp), DIMENSION(7), OPTIONAL, INTENT(IN) :: openpmd_unit_dimension
REAL(kind=dp), OPTIONAL, INTENT(IN) :: openpmd_unit_si
REAL(kind=dp), OPTIONAL, INTENT(IN) :: sim_time
TYPE(cp_openpmd_per_call_value_type) :: res
LOGICAL, SAVE :: opened_new_iteration = .FALSE.
LOGICAL :: opened_new_iteration
TYPE(openpmd_attributable_type) :: attr
TYPE(cp_openpmd_per_callsite_value_type), POINTER :: opmd
opened_new_iteration = .FALSE.
opmd => cp_openpmd_per_callsite(openpmd_file_index)%value
res%series = opmd%output_series
@ -420,9 +437,23 @@ CONTAINS
opened_new_iteration = cp_advance_iteration_number(logger, opmd)
res%iteration = opmd%output_series%write_iteration(opmd%iteration_counter%flat_iteration)
CALL res%iteration%set_time_unit_SI(seconds)
res%name_prefix = TRIM(middle_name)
IF (opened_new_iteration) THEN
IF (PRESENT(sim_time)) THEN
CALL res%iteration%set_time(sim_time)
END IF
IF (PRESENT(openpmd_unit_dimension)) THEN
res%unit_dimension = openpmd_unit_dimension
END IF
IF (PRESENT(openpmd_unit_si)) THEN
res%unit_si = openpmd_unit_si
END IF
attr = res%iteration%as_attributable()
CALL attr%set_attribute_vec_int( &
"ndim_iteration_index", &
@ -441,7 +472,7 @@ CONTAINS
INTEGER :: file_index
CHARACTER(:), ALLOCATABLE :: merged_config
CHARACTER(len=default_string_length), SAVE :: basename_copied = ' '
CHARACTER(len=default_string_length) :: basename_copied
TYPE(cp_openpmd_per_callsite_value_type) :: emplace_new
INTEGER :: handle
@ -449,6 +480,7 @@ CONTAINS
TYPE(openpmd_iteration_type) :: iteration
INTEGER :: i
basename_copied = ' '
basename_copied(1:LEN_TRIM(openpmd_basename)) = TRIM(openpmd_basename)
file_index = cp_openpmd_get_index_filedata(basename_copied)
@ -482,6 +514,9 @@ CONTAINS
emplace_new%output_series = openpmd_series_create( &
filename, openpmd_access_create, config=merged_config)
END IF
CALL emplace_new%output_series%set_software("CP2K", cp2k_version)
DEALLOCATE (merged_config)
file_index = cp_openpmd_add_filedata(basename_copied, emplace_new)
END FUNCTION cp_openpmd_get_openpmd_file_entry
@ -538,12 +573,22 @@ CONTAINS
!> all MPI libraries.
!> \param openpmd_basename Used to associate an identifier to each callsite of this module
!> \param use_openpmd ...
!> \param sim_time Current simulation time in atomic units (used for timeOffset)
!> \return ...
! **************************************************************************************************
FUNCTION cp_openpmd_print_key_unit_nr(logger, basis_section, print_key_path, &
middle_name, ignore_should_output, &
mpi_io, &
fout, openpmd_basename) RESULT(res)
FUNCTION cp_openpmd_print_key_unit_nr( &
logger, &
basis_section, &
print_key_path, &
middle_name, &
ignore_should_output, &
mpi_io, &
fout, &
openpmd_basename, &
openpmd_unit_dimension, &
openpmd_unit_si, &
sim_time) RESULT(res)
TYPE(cp_logger_type), POINTER :: logger
TYPE(section_vals_type), INTENT(IN) :: basis_section
CHARACTER(len=*), INTENT(IN), OPTIONAL :: print_key_path
@ -553,6 +598,9 @@ CONTAINS
CHARACTER(len=default_path_length), INTENT(OUT), &
OPTIONAL :: fout
CHARACTER(len=*), INTENT(IN), OPTIONAL :: openpmd_basename
REAL(kind=dp), DIMENSION(7), OPTIONAL, INTENT(IN) :: openpmd_unit_dimension
REAL(kind=dp), OPTIONAL, INTENT(IN) :: openpmd_unit_si
REAL(kind=dp), OPTIONAL, INTENT(IN) :: sim_time
INTEGER :: res
#ifdef __OPENPMD
@ -612,10 +660,19 @@ CONTAINS
openpmd_basename, filename, openpmd_config, logger, my_mpi_io)
OPEN (newunit=res, status='scratch', action='write')
openpmd_call_index = cp_openpmd_add_unit_nr( &
res, &
cp_openpmd_create_unit_nr_entry( &
openpmd_file_index, middle_name, logger))
BLOCK
TYPE(cp_openpmd_per_call_value_type) :: output_info
output_info = cp_openpmd_create_unit_nr_entry( &
openpmd_file_index, &
middle_name, &
logger, &
openpmd_unit_dimension, &
openpmd_unit_si, &
sim_time)
openpmd_call_index = cp_openpmd_add_unit_nr( &
res, output_info)
END BLOCK
ELSE
res = -1
@ -629,6 +686,10 @@ CONTAINS
MARK_USED(mpi_io)
MARK_USED(fout)
MARK_USED(openpmd_basename)
MARK_USED(openpmd_basename)
MARK_USED(openpmd_unit_dimension)
MARK_USED(openpmd_unit_si)
MARK_USED(sim_time)
res = 0
CPABORT("CP2K compiled without the openPMD-api")
#endif

View file

@ -93,6 +93,7 @@ MODULE openpmd_api
PROCEDURE, PUBLIC :: as_attributable => openpmd_series_as_attributable
PROCEDURE, PUBLIC :: close => openpmd_series_close
PROCEDURE, PUBLIC :: present => openpmd_series_present
PROCEDURE, PUBLIC :: set_software => openpmd_series_set_software
PROCEDURE, PUBLIC :: write_iteration => openpmd_series_write_iteration
PROCEDURE, PUBLIC :: get_iteration => openpmd_series_get_iteration
END TYPE openpmd_series_type
@ -107,6 +108,8 @@ MODULE openpmd_api
PROCEDURE, PUBLIC :: close => openpmd_iteration_close
PROCEDURE, PUBLIC :: open => openpmd_iteration_open
PROCEDURE, PUBLIC :: closed => openpmd_iteration_closed
PROCEDURE, PUBLIC :: set_time_unit_SI => openpmd_iteration_set_time_unit_SI
PROCEDURE, PUBLIC :: set_time => openpmd_iteration_set_time
END TYPE openpmd_iteration_type
TYPE openpmd_mesh_type
@ -117,7 +120,9 @@ MODULE openpmd_api
PROCEDURE, PUBLIC :: set_axis_labels => openpmd_mesh_set_axis_labels
PROCEDURE, PUBLIC :: set_grid_global_offset => openpmd_mesh_set_grid_global_offset
PROCEDURE, PUBLIC :: set_grid_spacing => openpmd_mesh_set_grid_spacing
PROCEDURE, PUBLIC :: set_grid_unit_SI => openpmd_mesh_set_grid_unit_SI
PROCEDURE, PUBLIC :: set_position => openpmd_mesh_set_position
PROCEDURE, PUBLIC :: set_unit_dimension => openpmd_mesh_set_unit_dimension
END TYPE openpmd_mesh_type
TYPE openpmd_particle_species_type
@ -134,6 +139,7 @@ MODULE openpmd_api
PROCEDURE, PUBLIC :: make_empty => openpmd_record_component_make_empty
PROCEDURE, PUBLIC :: make_constant_zero => openpmd_record_component_make_constant_zero
PROCEDURE, PUBLIC :: reset_dataset => openpmd_record_component_reset_dataset
PROCEDURE, PUBLIC :: set_unit_SI => openpmd_record_component_set_unit_SI
#:for type_f, type_enum in zip(dataset_types_f, dataset_types_enum)
#:for dim in dimensions
PROCEDURE, PRIVATE :: store_chunk_${dim}$d_${type_enum}$ &
@ -173,6 +179,7 @@ MODULE openpmd_api
CONTAINS
PROCEDURE, PUBLIC :: as_record_component => openpmd_record_as_record_component
PROCEDURE, PUBLIC :: get_component => openpmd_record_get_component
PROCEDURE, PUBLIC :: set_unit_dimension => openpmd_record_set_unit_dimension
END TYPE openpmd_record_type
#:for dim in dimensions
@ -422,6 +429,35 @@ MODULE openpmd_api
series%c_ptr = C_NULL_PTR
END SUBROUTINE openpmd_series_close
! **************************************************************************************************
!> \brief ...
!> \param this ...
!> \param name ...
!> \param version ...
! **************************************************************************************************
SUBROUTINE openpmd_series_set_software(this, name, version)
CLASS(openpmd_series_type), INTENT(IN) :: this
CHARACTER(len=*), INTENT(IN) :: name
CHARACTER(len=*), INTENT(IN), OPTIONAL :: version
INTERFACE
SUBROUTINE openpmd_c_series_set_software(series, name, version) &
BIND(C, name="openPMD_Series_setSoftware")
IMPORT :: C_PTR, C_CHAR
TYPE(C_PTR), VALUE :: series
CHARACTER(kind=C_CHAR), DIMENSION(*) :: name
CHARACTER(kind=C_CHAR), DIMENSION(*) :: version
END SUBROUTINE openpmd_c_series_set_software
END INTERFACE
CPASSERT(C_ASSOCIATED(this%c_ptr))
IF (PRESENT(version)) THEN
CALL openpmd_c_series_set_software(this%c_ptr, name//C_NULL_CHAR, version//C_NULL_CHAR)
ELSE
CALL openpmd_c_series_set_software(this%c_ptr, name//C_NULL_CHAR, C_NULL_CHAR)
END IF
END SUBROUTINE openpmd_series_set_software
! **************************************************************************************************
!> \brief ...
!> \param series ...
@ -572,6 +608,38 @@ MODULE openpmd_api
CALL openpmd_c_attributable_series_flush(this%c_ptr, my_backendconfig//C_NULL_CHAR)
END SUBROUTINE openpmd_attributable_series_flush
SUBROUTINE openpmd_iteration_set_time_unit_SI(this, timeUnitSI)
CLASS(openpmd_iteration_type) :: this
REAL(kind=C_DOUBLE), INTENT(IN), VALUE :: timeUnitSI
INTERFACE
SUBROUTINE openpmd_c_iteration_set_time_unit_SI(iteration, timeUnitSI) &
BIND(C, NAME="openPMD_Iteration_setTimeUnitSI")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: iteration
REAL(kind=C_DOUBLE), VALUE :: timeUnitSI
END SUBROUTINE openpmd_c_iteration_set_time_unit_SI
END INTERFACE
CALL openpmd_c_iteration_set_time_unit_SI(this%c_ptr, timeUnitSI)
END SUBROUTINE openpmd_iteration_set_time_unit_SI
SUBROUTINE openpmd_iteration_set_time(this, time)
CLASS(openpmd_iteration_type) :: this
REAL(kind=C_DOUBLE), INTENT(IN), VALUE :: time
INTERFACE
SUBROUTINE openpmd_c_iteration_set_time(iteration, time) &
BIND(C, NAME="openPMD_Iteration_setTime")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: iteration
REAL(kind=C_DOUBLE), VALUE :: time
END SUBROUTINE openpmd_c_iteration_set_time
END INTERFACE
CALL openpmd_c_iteration_set_time(this%c_ptr, time)
END SUBROUTINE openpmd_iteration_set_time
FUNCTION openpmd_json_merge(into, from, mpi_comm) RESULT(merged)
CHARACTER(len=*), INTENT(IN) :: into, from
CLASS(mp_comm_type), INTENT(in), OPTIONAL :: mpi_comm
@ -873,6 +941,22 @@ MODULE openpmd_api
this%c_ptr, dtype, SIZE(extent), C_LOC(extent(1)), 1, my_config//C_NULL_CHAR)
END SUBROUTINE openpmd_record_component_reset_dataset
SUBROUTINE openpmd_record_component_set_unit_SI(this, unitSI)
CLASS(openpmd_record_component_type) :: this
REAL(kind=C_DOUBLE), INTENT(IN), VALUE :: unitSI
INTERFACE
SUBROUTINE openpmd_c_record_component_set_unit_SI(record_component, unitSI) &
BIND(C, NAME="openPMD_RecordComponent_setUnitSI")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: record_component
REAL(kind=C_DOUBLE), VALUE :: unitSI
END SUBROUTINE openpmd_c_record_component_set_unit_SI
END INTERFACE
CALL openpmd_c_record_component_set_unit_SI(this%c_ptr, unitSI)
END SUBROUTINE openpmd_record_component_set_unit_SI
#:for type_f, type_enum in zip(dataset_types_f, dataset_types_enum)
#:for dim in [1, 2, 3]
! e.g. ':,:,:'
@ -1032,6 +1116,29 @@ MODULE openpmd_api
CPASSERT(C_ASSOCIATED(record_component%c_ptr))
END FUNCTION openpmd_record_get_component
! L = 0, //!< length
! M, //!< mass
! T, //!< time
! I, //!< electric current
! theta, //!< thermodynamic temperature
! N, //!< amount of substance
! J //!< luminous intensity
SUBROUTINE openpmd_record_set_unit_dimension(this, unitDimension)
CLASS(openpmd_record_type) :: this
REAL(kind=C_DOUBLE), DIMENSION(7), TARGET :: unitDimension
INTERFACE
SUBROUTINE openpmd_c_record_set_unit_dimension(record, unitDimension) &
BIND(C, NAME="openPMD_Record_setUnitDimension")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: record
REAL(kind=C_DOUBLE), DIMENSION(*) :: unitDimension
END SUBROUTINE openpmd_c_record_set_unit_dimension
END INTERFACE
CALL openpmd_c_record_set_unit_dimension(this%c_ptr, unitDimension)
END SUBROUTINE openpmd_record_set_unit_dimension
#:set f_routines = ["openpmd_mesh_set_grid_global_offset", "openpmd_mesh_set_grid_spacing", "openpmd_mesh_set_position"]
#:set c_functions = ["openPMD_Mesh_setGridGlobalOffset", "openPMD_Mesh_setGridSpacing", "openPMD_Mesh_setPosition"]
#:set f_vector_types = ["REAL(KIND=dp)" for _ in range(3)]
@ -1058,6 +1165,45 @@ MODULE openpmd_api
END SUBROUTINE ${f_routine}$
#:endfor
SUBROUTINE openpmd_mesh_set_grid_unit_SI(this, gridUnitSI)
CLASS(openpmd_mesh_type) :: this
REAL(kind=C_DOUBLE), INTENT(IN), VALUE :: gridUnitSI
INTERFACE
SUBROUTINE openpmd_c_mesh_set_grid_unit_SI(mesh, gridUnitSI) &
BIND(C, NAME="openPMD_Mesh_setGridUnitSI")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: mesh
REAL(kind=C_DOUBLE), VALUE :: gridUnitSI
END SUBROUTINE openpmd_c_mesh_set_grid_unit_SI
END INTERFACE
CALL openpmd_c_mesh_set_grid_unit_SI(this%c_ptr, gridUnitSI)
END SUBROUTINE openpmd_mesh_set_grid_unit_SI
! L = 0, //!< length
! M, //!< mass
! T, //!< time
! I, //!< electric current
! theta, //!< thermodynamic temperature
! N, //!< amount of substance
! J //!< luminous intensity
SUBROUTINE openpmd_mesh_set_unit_dimension(this, unitDimension)
CLASS(openpmd_mesh_type) :: this
REAL(kind=C_DOUBLE), DIMENSION(7), TARGET :: unitDimension
INTERFACE
SUBROUTINE openpmd_c_mesh_set_unit_dimension(mesh, unitDimension) &
BIND(C, NAME="openPMD_Mesh_setUnitDimension")
IMPORT :: C_PTR, C_DOUBLE
TYPE(C_PTR), VALUE :: mesh
REAL(kind=C_DOUBLE), DIMENSION(*) :: unitDimension
END SUBROUTINE openpmd_c_mesh_set_unit_dimension
END INTERFACE
CALL openpmd_c_mesh_set_unit_dimension(this%c_ptr, unitDimension)
END SUBROUTINE openpmd_mesh_set_unit_dimension
#endif
! In two lines of this file, `make pretty` applies a wrong indentation and keeps it for the rest of this file.
! Hence, we end up a bit further to the right than we should.

View file

@ -132,6 +132,10 @@ int openPMD_Series_upcast_to_Attributable(
int openPMD_Series_present(openPMD_Series series);
int openPMD_Series_setSoftware(
// in
openPMD_Series series, char const *name, char const *version);
char const *openPMD_get_default_extension();
/*************************
@ -183,6 +187,14 @@ int openPMD_Iteration_closed(
// in
openPMD_Iteration iteration);
int openPMD_Iteration_setTimeUnitSI(
// in
openPMD_Iteration iteration, double const timeUnitSI);
int openPMD_Iteration_setTime(
// in
openPMD_Iteration iteration, double const time);
/****************
* Mesh members *
****************/
@ -215,6 +227,14 @@ int openPMD_Mesh_setPosition(
// in
openPMD_Mesh mesh, double const *labels, int len_labels, int invert);
int openPMD_Mesh_setGridUnitSI(
// in
openPMD_Mesh mesh, double const gridUnitSI);
int openPMD_Mesh_setUnitDimension(
// in
openPMD_Mesh mesh, double const *unitDimension);
/***************************
* RecordComponent members *
***************************/
@ -239,6 +259,10 @@ int openPMD_RecordComponent_makeConstant(
openPMD_RecordComponent rc, openPMD_Datatype, int dimensions,
int const *extent, int invert, void const *value);
int openPMD_RecordComponent_setUnitSI(
// in
openPMD_RecordComponent rc, double const unitSI);
int openPMD_RecordComponent_storeChunk(
// in
openPMD_RecordComponent rc, openPMD_Datatype, int dimensions,
@ -293,6 +317,10 @@ int openPMD_Record_get_Component(
// out
openPMD_RecordComponent *rc);
int openPMD_Record_setUnitDimension(
// in
openPMD_Record record, double const *unitDimension);
/***********
* Helpers *
***********/
@ -316,6 +344,7 @@ char *openPMD_json_merge(char const *into, char const *from,
#endif
#include <any>
#include <array>
#include <cctype>
#include <cstdio>
#include <fstream>
@ -650,6 +679,14 @@ int openPMD_Series_present(openPMD_Series series_param) {
return res ? 1 : 0;
}
int openPMD_Series_setSoftware(openPMD_Series series_param, char const *name,
char const *version) {
auto series = reinterpret_cast<openPMD::Series *>(series_param);
std::string version_str = version ? std::string(version) : std::string();
series->setSoftware(std::string(name), version_str);
return 0;
}
char const *openPMD_get_default_extension() {
constexpr static char const *preferred_order[] = {"bp5", "bp4", "bp", "h5",
"json"};
@ -722,7 +759,6 @@ int openPMD_Iteration_get_particle_species(
int openPMD_Iteration_open(openPMD_Iteration iteration_param) {
auto iteration = reinterpret_cast<openPMD::Iteration *>(iteration_param);
std::cout << "OPENING ITERATION" << std::endl;
iteration->open();
return 0;
}
@ -738,6 +774,22 @@ int openPMD_Iteration_closed(openPMD_Iteration iteration_param) {
return iteration->closed();
}
int openPMD_Iteration_setTimeUnitSI(
// in
openPMD_Iteration iteration_param, double const timeUnitSI) {
auto iteration = reinterpret_cast<openPMD::Iteration *>(iteration_param);
iteration->setTimeUnitSI(timeUnitSI);
return 0;
}
int openPMD_Iteration_setTime(
// in
openPMD_Iteration iteration_param, double const time) {
auto iteration = reinterpret_cast<openPMD::Iteration *>(iteration_param);
iteration->setTime(time);
return 0;
}
int openPMD_Mesh_upcast_to_RecordComponent(
// in
openPMD_Mesh mesh_param,
@ -792,6 +844,24 @@ int openPMD_Mesh_setPosition(
return 0;
}
int openPMD_Mesh_setGridUnitSI(
// in
openPMD_Mesh mesh_param, double const gridUnitSI) {
auto mesh = reinterpret_cast<openPMD::Mesh *>(mesh_param);
mesh->setGridUnitSI(gridUnitSI);
return 0;
}
int openPMD_Mesh_setUnitDimension(
// in
openPMD_Mesh mesh_param, double const *unitDimension) {
auto mesh = reinterpret_cast<openPMD::Mesh *>(mesh_param);
std::array<double, 7> unitDim{};
std::copy_n(unitDimension, 7, unitDim.begin());
mesh->setUnitDimension(unitDim);
return 0;
}
int openPMD_MeshRecordComponent_upcast_to_RecordComponent(
// in
openPMD_MeshRecordComponent mrc,
@ -846,6 +916,14 @@ int openPMD_RecordComponent_makeConstant(
return 0;
}
int openPMD_RecordComponent_setUnitSI(
// in
openPMD_RecordComponent rc_param, double const unitSI) {
auto rc = reinterpret_cast<openPMD::RecordComponent *>(rc_param);
rc->setUnitSI(unitSI);
return 0;
}
int openPMD_RecordComponent_storeChunk(
// in
openPMD_RecordComponent rc_param, openPMD_Datatype dt, int dimensions,
@ -928,6 +1006,16 @@ int openPMD_Record_get_Component(
return 0;
}
int openPMD_Record_setUnitDimension(
// in
openPMD_Record record_param, double const *unitDimension) {
auto record = reinterpret_cast<openPMD::Record *>(record_param);
std::array<double, 7> unitDim{};
std::copy_n(unitDimension, 7, unitDim.begin());
record->setUnitDimension(unitDim);
return 0;
}
char *openPMD_json_merge(char const *into, char const *from,
MPI_Comm maybe_comm) {
#if OPENPMDAPI_VERSION_GE(0, 17, 0)

View file

@ -22,6 +22,7 @@ MODULE realspace_grid_openpmd
file_amode_rdonly, file_offset, mp_comm_type, mp_file_descriptor_type, mp_file_type, &
mp_file_type_free, mp_file_type_hindexed_make_chv, mp_file_type_set_view_chv, &
mpi_character_size
USE physcon, ONLY: a_bohr, e_charge, seconds
#if defined(__parallel)
#if defined(__MPI_F08)
USE mpi_f08, ONLY: mpi_allreduce, mpi_integer, mpi_bxor, mpi_allgather
@ -172,11 +173,15 @@ CONTAINS
position_offset = species%get_record("positionOffset")
position = species%get_record("position")
! length
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])
DO k = 1, SIZE(dims)
position_offset_component = position_offset%get_component(dims(k))
CALL position_offset_component%make_constant_zero(openpmd_type_int, global_extent)
CALL position_offset_component%set_unit_SI(a_bohr) ! doesnt really matter as it is zero
position_component = position%get_component(dims(k))
CALL position_component%reset_dataset(openpmd_type_double, global_extent)
CALL position_component%set_unit_SI(a_bohr)
unresolved_write_buffers(k) = &
position_component%store_chunk_span_1d_double(global_offset, local_extent)
write_buffers(k)%buffer => unresolved_write_buffers(k)%resolve_double(DEALLOCATE=.FALSE.)
@ -185,7 +190,10 @@ CONTAINS
IF (PRESENT(particles_zeff)) THEN
charge = species%get_record("charge")
charge_component = charge%as_record_component()
! charge ~ time x current
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])
CALL charge_component%reset_dataset(openpmd_type_double, global_extent)
CALL charge_component%set_unit_SI(e_charge)
unresolved_charge_write_buffer = charge_component%store_chunk_span_1d_double(global_offset, local_extent)
charge_write_buffer%buffer => unresolved_charge_write_buffer%resolve_double(DEALLOCATE=.FALSE.)
END IF
@ -405,7 +413,10 @@ CONTAINS
CALL mesh%set_position([0.5_dp, 0.5_dp, 0.5_dp])
CALL mesh%set_grid_global_offset([0._dp, 0._dp, 0._dp])
CALL mesh%set_grid_spacing(grid_spacing)
CALL mesh%set_grid_unit_SI(a_bohr)
CALL mesh%set_unit_dimension(openpmd_data%unit_dimension)
scalar_mesh = mesh%as_record_component()
CALL scalar_mesh%set_unit_SI(openpmd_data%unit_si)
CALL scalar_mesh%reset_dataset(openpmd_type_double, global_extent)
! shortcut

View file

@ -100,7 +100,8 @@ MODULE qs_scf_post_gpw
USE orbital_pointers, ONLY: indso
USE particle_list_types, ONLY: particle_list_type
USE particle_types, ONLY: particle_type
USE physcon, ONLY: angstrom,&
USE physcon, ONLY: a_bohr,&
angstrom,&
evolt
USE population_analyses, ONLY: lowdin_population_analysis,&
mulliken_population_analysis
@ -230,6 +231,16 @@ MODULE qs_scf_post_gpw
INTEGER, PARAMETER :: grid_output_cubes = 1, grid_output_openpmd = 2
REAL(kind=dp), DIMENSION(7), PARAMETER :: openpmd_unit_dimension_density = &
[-3, 0, 0, 0, 0, 0, 0]
REAL(kind=dp), DIMENSION(7), PARAMETER :: openpmd_unit_dimension_dimensionless = &
[0, 0, 0, 0, 0, 0, 0]
REAL(kind=dp), DIMENSION(7), PARAMETER :: openpmd_unit_dimension_wavefunction = &
[-1.5_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp, 0.0_dp]
REAL(kind=dp), PARAMETER :: openpmd_unit_si_density = a_bohr**(-3)
REAL(kind=dp), PARAMETER :: openpmd_unit_si_dimensionless = 1.0_dp
REAL(kind=dp), PARAMETER :: openpmd_unit_si_wavefunction = a_bohr**(-1.5_dp)
! Generic information on whether a certain output section has been activated
! or not, and on whether it has been activated in the Cube or openPMD variant.
! Create with function cube_or_openpmd(), see there for further details.
@ -338,12 +349,35 @@ CONTAINS
!> \param mpi_io ...
!> \param fout ...
!> \param openpmd_basename ...
!> \param openpmd_unit_dimension ...
!> \param openpmd_unit_si ...
!> \param sim_time ...
!> \return ...
! **************************************************************************************************
FUNCTION cp_forward_print_key_unit_nr(self, logger, basis_section, print_key_path, extension, &
middle_name, local, log_filename, ignore_should_output, file_form, file_position, &
file_action, file_status, do_backup, on_file, is_new_file, mpi_io, &
fout, openpmd_basename) RESULT(res)
FUNCTION cp_forward_print_key_unit_nr( &
self, &
logger, &
basis_section, &
print_key_path, &
extension, &
middle_name, &
local, &
log_filename, &
ignore_should_output, &
file_form, &
file_position, &
file_action, &
file_status, &
do_backup, &
on_file, &
is_new_file, &
mpi_io, &
fout, &
openpmd_basename, &
openpmd_unit_dimension, &
openpmd_unit_si, &
sim_time) RESULT(res)
CLASS(cp_section_key), INTENT(IN) :: self
TYPE(cp_logger_type), POINTER :: logger
TYPE(section_vals_type), INTENT(IN) :: basis_section
@ -359,6 +393,9 @@ CONTAINS
CHARACTER(len=default_path_length), INTENT(OUT), &
OPTIONAL :: fout
CHARACTER(len=*), INTENT(IN), OPTIONAL :: openpmd_basename
REAL(kind=dp), DIMENSION(7), OPTIONAL, INTENT(IN) :: openpmd_unit_dimension
REAL(kind=dp), OPTIONAL, INTENT(IN) :: openpmd_unit_si
REAL(kind=dp), OPTIONAL, INTENT(IN) :: sim_time
INTEGER :: res
IF (self%grid_output == grid_output_cubes) THEN
@ -370,9 +407,18 @@ CONTAINS
file_status=file_status, do_backup=do_backup, on_file=on_file, &
is_new_file=is_new_file, mpi_io=mpi_io, fout=fout)
ELSE
res = cp_openpmd_print_key_unit_nr(logger, basis_section, print_key_path, &
middle_name=middle_name, ignore_should_output=ignore_should_output, &
mpi_io=mpi_io, fout=fout, openpmd_basename=openpmd_basename)
res = cp_openpmd_print_key_unit_nr( &
logger, &
basis_section, &
print_key_path, &
middle_name=middle_name, &
ignore_should_output=ignore_should_output, &
mpi_io=mpi_io, &
fout=fout, &
openpmd_basename=openpmd_basename, &
openpmd_unit_dimension=openpmd_unit_dimension, &
openpmd_unit_si=openpmd_unit_si, &
sim_time=sim_time)
END IF
END FUNCTION cp_forward_print_key_unit_nr
@ -1442,9 +1488,19 @@ CONTAINS
WRITE (filename, '(a4,I5.5,a1,I1.1)') "WFN_", ivector, "_", ispin
mpi_io = .TRUE.
unit_nr = mo_section%print_key_unit_nr(logger, input, mo_section%absolute_section_key, extension=".cube", &
middle_name=TRIM(filename), file_position=my_pos_cube, log_filename=.FALSE., &
mpi_io=mpi_io, openpmd_basename="dft-mo")
unit_nr = mo_section%print_key_unit_nr( &
logger, &
input, &
mo_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
openpmd_basename="dft-mo", &
openpmd_unit_dimension=openpmd_unit_dimension_wavefunction, &
openpmd_unit_si=openpmd_unit_si_wavefunction, &
sim_time=qs_env%sim_time)
WRITE (title, *) "WAVEFUNCTION ", ivector, " spin ", ispin, " i.e. HOMO - ", ivector - homo
CALL mo_section%write_pw(wf_r, unit_nr, title, particles=particles, &
stride=section_get_ivals(dft_section, mo_section%concat_to_relative("%STRIDE")), &
@ -1542,9 +1598,19 @@ CONTAINS
WRITE (filename, '(a4,I5.5,a1,I1.1)') "WFN_", index_mo, "_", ispin
mpi_io = .TRUE.
unit_nr = mo_section%print_key_unit_nr(logger, input, mo_section%absolute_section_key, extension=".cube", &
middle_name=TRIM(filename), file_position=my_pos_cube, log_filename=.FALSE., &
mpi_io=mpi_io, openpmd_basename="dft-mo")
unit_nr = mo_section%print_key_unit_nr( &
logger, &
input, &
mo_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
openpmd_basename="dft-mo", &
openpmd_unit_dimension=openpmd_unit_dimension_wavefunction, &
openpmd_unit_si=openpmd_unit_si_wavefunction, &
sim_time=qs_env%sim_time)
WRITE (title, *) "WAVEFUNCTION ", index_mo, " spin ", ispin, " i.e. LUMO + ", ifirst + ivector - 2
CALL mo_section%write_pw(wf_r, unit_nr, title, particles=particles, &
stride=section_get_ivals(dft_section, mo_section%concat_to_relative("%STRIDE")), &
@ -1886,9 +1952,19 @@ CONTAINS
WRITE (title, *) "ELF spin ", ispin
mpi_io = .TRUE.
unit_nr = elf_section_key%print_key_unit_nr( &
logger, input, elf_section_key%absolute_section_key, extension=".cube", &
middle_name=TRIM(filename), file_position=my_pos_cube, log_filename=.FALSE., &
mpi_io=mpi_io, fout=mpi_filename, openpmd_basename="dft-elf")
logger, &
input, &
elf_section_key%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-elf", &
openpmd_unit_dimension=openpmd_unit_dimension_dimensionless, &
openpmd_unit_si=openpmd_unit_si_dimensionless, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2617,10 +2693,20 @@ CONTAINS
rho_total_rspace = pw_integrate_function(rho_elec_rspace, isign=-1)
filename = "TOTAL_ELECTRON_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-total-electron-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-total-electron-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2662,10 +2748,20 @@ CONTAINS
rho_total_rspace = pw_integrate_function(rho_elec_rspace, isign=-1)
filename = "TOTAL_SPIN_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-total-spin-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-total-spin-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2702,10 +2798,20 @@ CONTAINS
CALL pw_axpy(rho_r(2), rho_elec_rspace)
filename = "ELECTRON_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-electron-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-electron-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2724,10 +2830,20 @@ CONTAINS
CALL pw_axpy(rho_r(2), rho_elec_rspace, alpha=-1.0_dp)
filename = "SPIN_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-spin-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-spin-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2746,10 +2862,20 @@ CONTAINS
ELSE
filename = "ELECTRON_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-electron-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-electron-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2797,10 +2923,20 @@ CONTAINS
filename = "ELECTRON_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-electron-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-electron-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)
@ -2838,10 +2974,20 @@ CONTAINS
filename = "SPIN_DENSITY"
mpi_io = .TRUE.
unit_nr = e_density_section%print_key_unit_nr(logger, input, e_density_section%absolute_section_key, &
extension=".cube", middle_name=TRIM(filename), &
file_position=my_pos_cube, log_filename=.FALSE., mpi_io=mpi_io, &
fout=mpi_filename, openpmd_basename="dft-spin-density")
unit_nr = e_density_section%print_key_unit_nr( &
logger, &
input, &
e_density_section%absolute_section_key, &
extension=".cube", &
middle_name=TRIM(filename), &
file_position=my_pos_cube, &
log_filename=.FALSE., &
mpi_io=mpi_io, &
fout=mpi_filename, &
openpmd_basename="dft-spin-density", &
openpmd_unit_dimension=openpmd_unit_dimension_density, &
openpmd_unit_si=openpmd_unit_si_density, &
sim_time=qs_env%sim_time)
IF (output_unit > 0) THEN
IF (.NOT. mpi_io .AND. e_density_section%grid_output == grid_output_cubes) THEN
INQUIRE (UNIT=unit_nr, NAME=filename)