Skip to content

Commit

Permalink
Perform check of interatomic distances
Browse files Browse the repository at this point in the history
  • Loading branch information
@mkrack
mkrack committed Oct 22, 2021
1 parent ca5348d commit eb17c48
Show file tree
Hide file tree
Showing 2 changed files with 63 additions and 19 deletions.
15 changes: 12 additions & 3 deletions src/input_cp2k_subsys.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -372,10 +372,19 @@ SUBROUTINE create_subsys_print_section(section)
CALL section_add_subsection(section, print_key)
CALL section_release(print_key)

CALL cp_print_key_section_create(print_key, __LOCATION__, "interatomic_distances", &
description="controls the output of the interatomic distances when setting up the"// &
CALL cp_print_key_section_create(print_key, __LOCATION__, "INTERATOMIC_DISTANCES", &
description="Controls the printout of the interatomic distances when setting up the "// &
"force environment", unit_str="angstrom", &
print_level=debug_print_level, filename="__STD_OUT__")
CALL keyword_create(keyword, __LOCATION__, name="CHECK_INTERATOMIC_DISTANCES", &
description="Minimum allowed distance between two atoms. "// &
"A warning is printed, if a smaller interatomic distance is encountered. "// &
"The check is disabled for the threshold value 0. "// &
"The run is aborted, if an interatomic distance is smaller than the absolute "// &
"value of a negative threshold value.", &
default_r_val=0.5_dp*bohr, unit_str="angstrom")
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)
CALL section_add_subsection(section, print_key)
CALL section_release(print_key)

Expand Down Expand Up @@ -480,7 +489,7 @@ SUBROUTINE create_subsys_print_section(section)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="EPS_GEO", &
description="Accuracy required for symmetry detection", &
default_r_val=1.e-4_dp)
default_r_val=1.0E-4_dp)
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="STANDARD_ORIENTATION", &
Expand Down
67 changes: 51 additions & 16 deletions src/particle_methods.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,7 +24,8 @@ MODULE particle_methods
real_to_scaled,&
set_cell_param
USE cp_log_handling, ONLY: cp_get_default_logger,&
cp_logger_type
cp_logger_type,&
cp_to_string
USE cp_output_handling, ONLY: cp_print_key_finished_output,&
cp_print_key_unit_nr
USE cp_units, ONLY: cp_unit_from_cp2k
Expand Down Expand Up @@ -619,7 +620,7 @@ SUBROUTINE write_particle_distances(particle_set, cell, subsys_section)
CHARACTER(LEN=default_string_length) :: unit_str
INTEGER :: handle, iatom, iw, jatom, natom
INTEGER, DIMENSION(3) :: periodic
REAL(KIND=dp) :: conv, dab
REAL(KIND=dp) :: conv, dab, dab_abort, dab_min, dab_warn
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: distance_matrix
REAL(KIND=dp), DIMENSION(3) :: rab
TYPE(cp_logger_type), POINTER :: logger
Expand All @@ -633,31 +634,63 @@ SUBROUTINE write_particle_distances(particle_set, cell, subsys_section)

CALL section_vals_val_get(subsys_section, "PRINT%INTERATOMIC_DISTANCES%UNIT", c_val=unit_str)
conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str))
IF (iw > 0) THEN
CALL get_cell(cell=cell, periodic=periodic)
natom = SIZE(particle_set)
ALLOCATE (distance_matrix(natom, natom))
distance_matrix(:, :) = 0.0_dp
DO iatom = 1, natom
DO jatom = iatom + 1, natom
rab(:) = pbc(particle_set(iatom)%r(:), &
particle_set(jatom)%r(:), cell)
dab = SQRT(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
distance_matrix(iatom, jatom) = dab*conv
distance_matrix(jatom, iatom) = distance_matrix(iatom, jatom)
END DO
CALL section_vals_val_get(subsys_section, "PRINT%INTERATOMIC_DISTANCES%CHECK_INTERATOMIC_DISTANCES", &
r_val=dab_min)
dab_abort = 0.0_dp
dab_warn = 0.0_dp
IF (dab_min > 0.0_dp) THEN
dab_warn = dab_min*conv
ELSE IF (dab_min < 0.0_dp) THEN
dab_abort = ABS(dab_min)*conv
END IF
CALL get_cell(cell=cell, periodic=periodic)
natom = SIZE(particle_set)
ALLOCATE (distance_matrix(natom, natom))
distance_matrix(:, :) = 0.0_dp
DO iatom = 1, natom
DO jatom = iatom + 1, natom
rab(:) = pbc(particle_set(iatom)%r(:), &
particle_set(jatom)%r(:), cell)
dab = SQRT(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
distance_matrix(iatom, jatom) = dab*conv
distance_matrix(jatom, iatom) = distance_matrix(iatom, jatom)
IF (dab_abort > 0.0_dp) THEN
! Stop the run for interatomic distances smaller than the requested threshold
IF (distance_matrix(iatom, jatom) < dab_abort) THEN
CALL cp_abort(__LOCATION__, "The distance between the atoms "// &
TRIM(ADJUSTL(cp_to_string(iatom, fmt="(I8)")))//" and "// &
TRIM(ADJUSTL(cp_to_string(jatom, fmt="(I8)")))//" is only "// &
TRIM(ADJUSTL(cp_to_string(distance_matrix(iatom, jatom), fmt="(F6.3)")))//" "// &
TRIM(ADJUSTL(unit_str))//" and thus smaller than the requested threshold of "// &
TRIM(ADJUSTL(cp_to_string(dab_abort, fmt="(F6.3)")))//" "// &
TRIM(ADJUSTL(unit_str)))
END IF
END IF
IF (distance_matrix(iatom, jatom) < dab_warn) THEN
! Print warning for interatomic distances smaller than 0.5 Angstrom
CALL cp_warn(__LOCATION__, "The distance between the atoms "// &
TRIM(ADJUSTL(cp_to_string(iatom, fmt="(I8)")))//" and "// &
TRIM(ADJUSTL(cp_to_string(jatom, fmt="(I8)")))//" is only "// &
TRIM(ADJUSTL(cp_to_string(distance_matrix(iatom, jatom), fmt="(F6.3)")))//" "// &
TRIM(ADJUSTL(unit_str))//" and thus smaller than the threshold of "// &
TRIM(ADJUSTL(cp_to_string(dab_warn, fmt="(F6.3)")))//" "// &
TRIM(ADJUSTL(unit_str)))
END IF
END DO
END DO

IF (iw > 0) THEN
! Print the distance matrix
WRITE (UNIT=iw, FMT="(/,/,T2,A)") &
"INTERATOMIC DISTANCES IN "//TRIM(unit_str)

CALL write_particle_matrix(distance_matrix, particle_set, iw)
END IF

CALL cp_print_key_finished_output(iw, logger, subsys_section, &
"PRINT%INTERATOMIC_DISTANCES")

IF (ALLOCATED(distance_matrix)) DEALLOCATE (distance_matrix)

CALL timestop(handle)

END SUBROUTINE write_particle_distances
Expand Down Expand Up @@ -722,7 +755,9 @@ SUBROUTINE write_particle_matrix(matrix, particle_set, iw, el_per_part, Ilist, p
(matrix(iatom, icol), icol=from, to)
END DO
END DO

DEALLOCATE (my_list)

END SUBROUTINE write_particle_matrix

! **************************************************************************************************
Expand Down

0 comments on commit eb17c48

Please sign in to comment.