LCOV - code coverage report
Current view: top level - src - qs_scf_csr_write.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:15a58fb) Lines: 209 216 96.8 %
Date: 2025-02-18 08:24:35 Functions: 5 5 100.0 %

          Line data    Source code
       1             : !--------------------------------------------------------------------------------------------------!
       2             : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3             : !   Copyright 2000-2025 CP2K developers group <https://cp2k.org>                                   !
       4             : !                                                                                                  !
       5             : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6             : !--------------------------------------------------------------------------------------------------!
       7             : 
       8             : ! **************************************************************************************************
       9             : !> \brief Functions to print the KS and S matrix in the CSR format to file
      10             : !> \par History
      11             : !>      Started as a copy from the relevant part of qs_scf_post_gpw
      12             : !> \author Fabian Ducry (05.2020)
      13             : ! **************************************************************************************************
      14             : MODULE qs_scf_csr_write
      15             :    USE cp_dbcsr_api,                    ONLY: &
      16             :         dbcsr_convert_dbcsr_to_csr, dbcsr_copy, dbcsr_create, &
      17             :         dbcsr_csr_create_and_convert_complex, dbcsr_csr_create_from_dbcsr, &
      18             :         dbcsr_csr_dbcsr_blkrow_dist, dbcsr_csr_destroy, dbcsr_csr_type, dbcsr_csr_write, &
      19             :         dbcsr_desymmetrize, dbcsr_finalize, dbcsr_get_block_p, dbcsr_has_symmetry, dbcsr_p_type, &
      20             :         dbcsr_put_block, dbcsr_release, dbcsr_set, dbcsr_type, dbcsr_type_antisymmetric, &
      21             :         dbcsr_type_no_symmetry, dbcsr_type_real_8, dbcsr_type_symmetric
      22             :    USE cp_dbcsr_cp2k_link,              ONLY: cp_dbcsr_alloc_block_from_nbl
      23             :    USE cp_dbcsr_operations,             ONLY: dbcsr_allocate_matrix_set,&
      24             :                                               dbcsr_deallocate_matrix_set
      25             :    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
      26             :                                               cp_logger_get_default_io_unit,&
      27             :                                               cp_logger_type
      28             :    USE cp_output_handling,              ONLY: cp_p_file,&
      29             :                                               cp_print_key_finished_output,&
      30             :                                               cp_print_key_should_output,&
      31             :                                               cp_print_key_unit_nr
      32             :    USE input_section_types,             ONLY: section_vals_get_subs_vals,&
      33             :                                               section_vals_type,&
      34             :                                               section_vals_val_get
      35             :    USE kinds,                           ONLY: default_path_length,&
      36             :                                               dp
      37             :    USE kpoint_methods,                  ONLY: rskp_transform
      38             :    USE kpoint_types,                    ONLY: get_kpoint_info,&
      39             :                                               kpoint_type
      40             :    USE qs_environment_types,            ONLY: get_qs_env,&
      41             :                                               qs_environment_type
      42             :    USE qs_neighbor_list_types,          ONLY: get_iterator_info,&
      43             :                                               get_neighbor_list_set_p,&
      44             :                                               neighbor_list_iterate,&
      45             :                                               neighbor_list_iterator_create,&
      46             :                                               neighbor_list_iterator_p_type,&
      47             :                                               neighbor_list_iterator_release,&
      48             :                                               neighbor_list_set_p_type
      49             : #include "./base/base_uses.f90"
      50             : 
      51             :    IMPLICIT NONE
      52             :    PRIVATE
      53             : 
      54             :    ! Global parameters
      55             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_scf_csr_write'
      56             :    PUBLIC :: write_ks_matrix_csr, &
      57             :              write_s_matrix_csr
      58             : 
      59             : ! **************************************************************************************************
      60             : 
      61             : CONTAINS
      62             : 
      63             : !**************************************************************************************************
      64             : !> \brief writing the KS matrix in csr format into a file
      65             : !> \param qs_env qs environment
      66             : !> \param input the input
      67             : !> \par History
      68             : !>       Moved to module qs_scf_csr_write (05.2020)
      69             : !> \author Mohammad Hossein Bani-Hashemian
      70             : ! **************************************************************************************************
      71       17161 :    SUBROUTINE write_ks_matrix_csr(qs_env, input)
      72             :       TYPE(qs_environment_type), POINTER                 :: qs_env
      73             :       TYPE(section_vals_type), POINTER                   :: input
      74             : 
      75             :       CHARACTER(len=*), PARAMETER :: routineN = 'write_ks_matrix_csr'
      76             : 
      77             :       INTEGER                                            :: handle, output_unit
      78             :       LOGICAL                                            :: do_kpoints, do_ks_csr_write, real_space
      79             :       TYPE(cp_logger_type), POINTER                      :: logger
      80       17161 :       TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER       :: matrix_ks
      81             :       TYPE(kpoint_type), POINTER                         :: kpoints
      82             :       TYPE(section_vals_type), POINTER                   :: dft_section
      83             : 
      84       17161 :       CALL timeset(routineN, handle)
      85             : 
      86             :       NULLIFY (dft_section)
      87             : 
      88       17161 :       logger => cp_get_default_logger()
      89       17161 :       output_unit = cp_logger_get_default_io_unit(logger)
      90             : 
      91       17161 :       dft_section => section_vals_get_subs_vals(input, "DFT")
      92             :       do_ks_csr_write = BTEST(cp_print_key_should_output(logger%iter_info, dft_section, &
      93       17161 :                                                          "PRINT%KS_CSR_WRITE"), cp_p_file)
      94             : 
      95       17161 :       IF (do_ks_csr_write) THEN
      96          54 :          CALL get_qs_env(qs_env=qs_env, kpoints=kpoints, matrix_ks_kp=matrix_ks, do_kpoints=do_kpoints)
      97             :          CALL section_vals_val_get(dft_section, "PRINT%KS_CSR_WRITE%REAL_SPACE", &
      98          54 :                                    l_val=real_space)
      99             : 
     100          54 :          IF (do_kpoints .AND. .NOT. real_space) THEN
     101             :             CALL write_matrix_kp_csr(mat=matrix_ks, dft_section=dft_section, &
     102          10 :                                      kpoints=kpoints, prefix="KS")
     103             :          ELSE
     104             :             CALL write_matrix_csr(dft_section, mat=matrix_ks, kpoints=kpoints, do_kpoints=do_kpoints, &
     105          44 :                                   prefix="KS")
     106             :          END IF
     107             :       END IF
     108             : 
     109       17161 :       CALL timestop(handle)
     110             : 
     111       17161 :    END SUBROUTINE write_ks_matrix_csr
     112             : 
     113             : !**************************************************************************************************
     114             : !> \brief writing the overlap matrix in csr format into a file
     115             : !> \param qs_env qs environment
     116             : !> \param input the input
     117             : !> \par History
     118             : !>      Moved to module qs_scf_csr_write
     119             : !> \author Mohammad Hossein Bani-Hashemian
     120             : ! **************************************************************************************************
     121       17161 :    SUBROUTINE write_s_matrix_csr(qs_env, input)
     122             :       TYPE(qs_environment_type), POINTER                 :: qs_env
     123             :       TYPE(section_vals_type), POINTER                   :: input
     124             : 
     125             :       CHARACTER(len=*), PARAMETER :: routineN = 'write_s_matrix_csr'
     126             : 
     127             :       INTEGER                                            :: handle, output_unit
     128             :       LOGICAL                                            :: do_kpoints, do_s_csr_write, real_space
     129             :       TYPE(cp_logger_type), POINTER                      :: logger
     130       17161 :       TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER       :: matrix_s
     131             :       TYPE(kpoint_type), POINTER                         :: kpoints
     132             :       TYPE(section_vals_type), POINTER                   :: dft_section
     133             : 
     134       17161 :       CALL timeset(routineN, handle)
     135             : 
     136             :       NULLIFY (dft_section)
     137             : 
     138       17161 :       logger => cp_get_default_logger()
     139       17161 :       output_unit = cp_logger_get_default_io_unit(logger)
     140             : 
     141       17161 :       dft_section => section_vals_get_subs_vals(input, "DFT")
     142             :       do_s_csr_write = BTEST(cp_print_key_should_output(logger%iter_info, dft_section, &
     143       17161 :                                                         "PRINT%S_CSR_WRITE"), cp_p_file)
     144             : 
     145       17161 :       IF (do_s_csr_write) THEN
     146          52 :          CALL get_qs_env(qs_env=qs_env, kpoints=kpoints, matrix_s_kp=matrix_s, do_kpoints=do_kpoints)
     147             :          CALL section_vals_val_get(dft_section, "PRINT%S_CSR_WRITE%REAL_SPACE", &
     148          52 :                                    l_val=real_space)
     149             : 
     150          52 :          IF (do_kpoints .AND. .NOT. real_space) THEN
     151             :             CALL write_matrix_kp_csr(mat=matrix_s, dft_section=dft_section, &
     152          10 :                                      kpoints=kpoints, prefix="S")
     153             :          ELSE
     154             :             CALL write_matrix_csr(dft_section, mat=matrix_s, kpoints=kpoints, do_kpoints=do_kpoints, &
     155          42 :                                   prefix="S")
     156             :          END IF
     157             :       END IF
     158             : 
     159       17161 :       CALL timestop(handle)
     160             : 
     161       17161 :    END SUBROUTINE write_s_matrix_csr
     162             : 
     163             : ! **************************************************************************************************
     164             : !> \brief helper function to print the real space representation of KS or S matrix to file
     165             : !> \param dft_section the dft_section
     166             : !> \param mat Hamiltonian or overlap matrix
     167             : !> \param kpoints Kpoint environment
     168             : !> \param prefix string to distinguish between KS and S matrix
     169             : !> \param do_kpoints Whether it is a gamma-point or k-point calculation
     170             : !> \par History
     171             : !>       Moved most of the code from write_ks_matrix_csr and write_s_matrix_csr
     172             : !>       Removed the code for printing k-point dependent matrices and added
     173             : !>              printing of the real space representation
     174             : ! **************************************************************************************************
     175          86 :    SUBROUTINE write_matrix_csr(dft_section, mat, kpoints, prefix, do_kpoints)
     176             :       TYPE(section_vals_type), INTENT(IN), POINTER       :: dft_section
     177             :       TYPE(dbcsr_p_type), DIMENSION(:, :), INTENT(IN), &
     178             :          POINTER                                         :: mat
     179             :       TYPE(kpoint_type), INTENT(IN), POINTER             :: kpoints
     180             :       CHARACTER(*), INTENT(in)                           :: prefix
     181             :       LOGICAL, INTENT(IN)                                :: do_kpoints
     182             : 
     183             :       CHARACTER(len=*), PARAMETER                        :: routineN = 'write_matrix_csr'
     184             : 
     185             :       CHARACTER(LEN=default_path_length)                 :: file_name, fileformat, subs_string
     186             :       INTEGER                                            :: handle, ic, ispin, ncell, nspin, &
     187             :                                                             output_unit, unit_nr
     188          86 :       INTEGER, ALLOCATABLE, DIMENSION(:, :)              :: index_to_cell
     189          86 :       INTEGER, ALLOCATABLE, DIMENSION(:, :, :)           :: cell_to_index
     190          86 :       INTEGER, DIMENSION(:, :), POINTER                  :: i2c
     191          86 :       INTEGER, DIMENSION(:, :, :), POINTER               :: c2i
     192             :       LOGICAL                                            :: bin, do_symmetric, uptr
     193             :       REAL(KIND=dp)                                      :: thld
     194             :       TYPE(cp_logger_type), POINTER                      :: logger
     195             :       TYPE(dbcsr_csr_type)                               :: mat_csr
     196          86 :       TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER       :: mat_nosym
     197             :       TYPE(dbcsr_type)                                   :: matrix_nosym
     198             :       TYPE(neighbor_list_set_p_type), DIMENSION(:), &
     199          86 :          POINTER                                         :: sab_nl
     200             : 
     201          86 :       CALL timeset(routineN, handle)
     202             : 
     203          86 :       logger => cp_get_default_logger()
     204          86 :       output_unit = cp_logger_get_default_io_unit(logger)
     205             : 
     206          86 :       subs_string = "PRINT%"//prefix//"_CSR_WRITE"
     207             : 
     208          86 :       CALL section_vals_val_get(dft_section, subs_string//"%THRESHOLD", r_val=thld)
     209          86 :       CALL section_vals_val_get(dft_section, subs_string//"%UPPER_TRIANGULAR", l_val=uptr)
     210          86 :       CALL section_vals_val_get(dft_section, subs_string//"%BINARY", l_val=bin)
     211             : 
     212          86 :       IF (bin) THEN
     213           2 :          fileformat = "UNFORMATTED"
     214             :       ELSE
     215          84 :          fileformat = "FORMATTED"
     216             :       END IF
     217             : 
     218          86 :       nspin = SIZE(mat, 1)
     219          86 :       ncell = SIZE(mat, 2)
     220             : 
     221          86 :       IF (do_kpoints) THEN
     222             : 
     223           2 :          i2c => kpoints%index_to_cell
     224           2 :          c2i => kpoints%cell_to_index
     225             : 
     226           2 :          NULLIFY (sab_nl)
     227           2 :          CALL get_kpoint_info(kpoints, sab_nl=sab_nl)
     228           2 :          CALL get_neighbor_list_set_p(neighbor_list_sets=sab_nl, symmetric=do_symmetric)
     229             : 
     230             :          ! desymmetrize the KS or S matrices if necessary
     231           2 :          IF (do_symmetric) THEN
     232             :             CALL desymmetrize_rs_matrix(mat, mat_nosym, cell_to_index, index_to_cell, kpoints)
     233           2 :             ncell = SIZE(index_to_cell, 2) ! update the number of cells
     234             :          ELSE
     235             :             ALLOCATE (cell_to_index(LBOUND(c2i, 1):UBOUND(c2i, 1), &
     236             :                                     LBOUND(c2i, 2):UBOUND(c2i, 2), &
     237           0 :                                     LBOUND(c2i, 3):UBOUND(c2i, 3)))
     238             :             cell_to_index(LBOUND(c2i, 1):UBOUND(c2i, 1), &
     239             :                           LBOUND(c2i, 2):UBOUND(c2i, 2), &
     240           0 :                           LBOUND(c2i, 3):UBOUND(c2i, 3)) = c2i
     241             : 
     242           0 :             ALLOCATE (index_to_cell(3, ncell))
     243           0 :             index_to_cell(1:3, 1:ncell) = i2c
     244             : 
     245           0 :             mat_nosym => mat
     246             :          END IF
     247             : 
     248             :          ! print the index to cell mapping to the output
     249           2 :          IF (output_unit > 0) THEN
     250           1 :             WRITE (output_unit, "(/,A15,T15,I4,A)") prefix//" CSR write| ", &
     251           2 :                ncell, " periodic images"
     252           1 :             WRITE (output_unit, "(T7,A,T17,A,T24,A,T31,A)") "Number", "X", "Y", "Z"
     253          28 :             DO ic = 1, ncell
     254          28 :                WRITE (output_unit, "(T8,I3,T15,I3,T22,I3,T29,I3)") ic, index_to_cell(:, ic)
     255             :             END DO
     256             :          END IF
     257             :       END IF
     258             : 
     259             :       ! write the csr file(s)
     260         172 :       DO ispin = 1, nspin
     261         310 :          DO ic = 1, ncell
     262         138 :             IF (do_kpoints) THEN
     263          54 :                CALL dbcsr_copy(matrix_nosym, mat_nosym(ispin, ic)%matrix)
     264          54 :                WRITE (file_name, '(2(A,I0))') prefix//"_SPIN_", ispin, "_R_", ic
     265             :             ELSE
     266          84 :                IF (dbcsr_has_symmetry(mat(ispin, ic)%matrix)) THEN
     267          84 :                   CALL dbcsr_desymmetrize(mat(ispin, ic)%matrix, matrix_nosym)
     268             :                ELSE
     269           0 :                   CALL dbcsr_copy(matrix_nosym, mat(ispin, ic)%matrix)
     270             :                END IF
     271          84 :                WRITE (file_name, '(A,I0)') prefix//"_SPIN_", ispin
     272             :             END IF
     273             :             ! Convert dbcsr to csr
     274             :             CALL dbcsr_csr_create_from_dbcsr(matrix_nosym, &
     275         138 :                                              mat_csr, dbcsr_csr_dbcsr_blkrow_dist)
     276         138 :             CALL dbcsr_convert_dbcsr_to_csr(matrix_nosym, mat_csr)
     277             :             ! Write to file
     278             :             unit_nr = cp_print_key_unit_nr(logger, dft_section, subs_string, &
     279             :                                            extension=".csr", middle_name=TRIM(file_name), &
     280         138 :                                            file_status="REPLACE", file_form=fileformat)
     281         138 :             CALL dbcsr_csr_write(mat_csr, unit_nr, upper_triangle=uptr, threshold=thld, binary=bin)
     282             : 
     283         138 :             CALL cp_print_key_finished_output(unit_nr, logger, dft_section, subs_string)
     284         138 :             CALL dbcsr_csr_destroy(mat_csr)
     285         224 :             CALL dbcsr_release(matrix_nosym)
     286             :          END DO
     287             :       END DO
     288             : 
     289             :       ! clean up
     290          86 :       IF (do_kpoints) THEN
     291           2 :          DEALLOCATE (cell_to_index, index_to_cell)
     292           2 :          IF (do_symmetric) THEN
     293           4 :             DO ispin = 1, nspin
     294          58 :                DO ic = 1, ncell
     295          56 :                   CALL dbcsr_release(mat_nosym(ispin, ic)%matrix)
     296             :                END DO
     297             :             END DO
     298           2 :             CALL dbcsr_deallocate_matrix_set(mat_nosym)
     299             :          END IF
     300             :       END IF
     301          86 :       CALL timestop(handle)
     302             : 
     303         172 :    END SUBROUTINE write_matrix_csr
     304             : 
     305             : ! **************************************************************************************************
     306             : !> \brief helper function to print the k-dependent KS or S matrix to file
     307             : !> \param mat Hamiltonian or overlap matrix for k-point calculations
     308             : !> \param dft_section the dft_section
     309             : !> \param kpoints Kpoint environment
     310             : !> \param prefix string to distinguish between KS and S matrix
     311             : !> \par History
     312             : !>       Moved the code from write_matrix_csr to write_matrix_kp_csr
     313             : !> \author Fabian Ducry
     314             : ! **************************************************************************************************
     315          20 :    SUBROUTINE write_matrix_kp_csr(mat, dft_section, kpoints, prefix)
     316             :       TYPE(dbcsr_p_type), DIMENSION(:, :), INTENT(IN), &
     317             :          POINTER                                         :: mat
     318             :       TYPE(section_vals_type), INTENT(IN), POINTER       :: dft_section
     319             :       TYPE(kpoint_type), INTENT(IN), POINTER             :: kpoints
     320             :       CHARACTER(*), INTENT(in)                           :: prefix
     321             : 
     322             :       CHARACTER(len=*), PARAMETER :: routineN = 'write_matrix_kp_csr'
     323             : 
     324             :       CHARACTER(LEN=default_path_length)                 :: file_name, fileformat, subs_string
     325             :       INTEGER                                            :: handle, igroup, ik, ikp, ispin, kplocal, &
     326             :                                                             nkp_groups, output_unit, unit_nr
     327             :       INTEGER, DIMENSION(2)                              :: kp_range
     328          20 :       INTEGER, DIMENSION(:, :), POINTER                  :: kp_dist
     329             :       LOGICAL                                            :: bin, uptr, use_real_wfn
     330             :       REAL(KIND=dp)                                      :: thld
     331          20 :       REAL(KIND=dp), DIMENSION(:, :), POINTER            :: xkp
     332             :       TYPE(cp_logger_type), POINTER                      :: logger
     333             :       TYPE(dbcsr_csr_type)                               :: mat_csr
     334             :       TYPE(dbcsr_type)                                   :: matrix_nosym
     335             :       TYPE(dbcsr_type), POINTER                          :: imatrix, imatrix_nosym, rmatrix, &
     336             :                                                             rmatrix_nosym
     337             :       TYPE(neighbor_list_set_p_type), DIMENSION(:), &
     338          20 :          POINTER                                         :: sab_nl
     339             : 
     340          20 :       CALL timeset(routineN, handle)
     341             : 
     342          20 :       logger => cp_get_default_logger()
     343          20 :       output_unit = cp_logger_get_default_io_unit(logger)
     344             : 
     345          20 :       subs_string = "PRINT%"//prefix//"_CSR_WRITE"
     346             : 
     347          20 :       CALL section_vals_val_get(dft_section, subs_string//"%THRESHOLD", r_val=thld)
     348          20 :       CALL section_vals_val_get(dft_section, subs_string//"%UPPER_TRIANGULAR", l_val=uptr)
     349          20 :       CALL section_vals_val_get(dft_section, subs_string//"%BINARY", l_val=bin)
     350             : 
     351          20 :       IF (bin) THEN
     352          20 :          fileformat = "UNFORMATTED"
     353             :       ELSE
     354           0 :          fileformat = "FORMATTED"
     355             :       END IF
     356             : 
     357          20 :       NULLIFY (sab_nl)
     358             : 
     359             :       !  Calculate the Hamiltonian at the k-points
     360             :       CALL get_kpoint_info(kpoints, xkp=xkp, use_real_wfn=use_real_wfn, kp_range=kp_range, &
     361          20 :                            nkp_groups=nkp_groups, kp_dist=kp_dist, sab_nl=sab_nl)
     362             : 
     363          20 :       ALLOCATE (rmatrix)
     364             :       CALL dbcsr_create(rmatrix, template=mat(1, 1)%matrix, &
     365          20 :                         matrix_type=dbcsr_type_symmetric)
     366          20 :       CALL cp_dbcsr_alloc_block_from_nbl(rmatrix, sab_nl)
     367             : 
     368          20 :       IF (.NOT. use_real_wfn) THEN
     369             :          ! Allocate temporary variables
     370          16 :          ALLOCATE (rmatrix_nosym, imatrix, imatrix_nosym)
     371             :          CALL dbcsr_create(rmatrix_nosym, template=mat(1, 1)%matrix, &
     372          16 :                            matrix_type=dbcsr_type_no_symmetry)
     373             :          CALL dbcsr_create(imatrix, template=mat(1, 1)%matrix, &
     374          16 :                            matrix_type=dbcsr_type_antisymmetric)
     375             :          CALL dbcsr_create(imatrix_nosym, template=mat(1, 1)%matrix, &
     376          16 :                            matrix_type=dbcsr_type_no_symmetry)
     377          16 :          CALL cp_dbcsr_alloc_block_from_nbl(rmatrix_nosym, sab_nl)
     378          16 :          CALL cp_dbcsr_alloc_block_from_nbl(imatrix, sab_nl)
     379          16 :          CALL cp_dbcsr_alloc_block_from_nbl(imatrix_nosym, sab_nl)
     380             :       END IF
     381             : 
     382          20 :       kplocal = kp_range(2) - kp_range(1) + 1
     383          56 :       DO ikp = 1, kplocal
     384          92 :          DO ispin = 1, SIZE(mat, 1)
     385         140 :             DO igroup = 1, nkp_groups
     386             :                ! number of current kpoint
     387          68 :                ik = kp_dist(1, igroup) + ikp - 1
     388          68 :                CALL dbcsr_set(rmatrix, 0.0_dp)
     389          68 :                IF (use_real_wfn) THEN
     390             :                   ! FT of the matrix
     391             :                   CALL rskp_transform(rmatrix=rmatrix, rsmat=mat, ispin=ispin, &
     392           4 :                                       xkp=xkp(1:3, ik), cell_to_index=kpoints%cell_to_index, sab_nl=sab_nl)
     393             :                   ! Convert to desymmetrized csr matrix
     394           4 :                   CALL dbcsr_desymmetrize(rmatrix, matrix_nosym)
     395           4 :                   CALL dbcsr_csr_create_from_dbcsr(matrix_nosym, mat_csr, dbcsr_csr_dbcsr_blkrow_dist)
     396           4 :                   CALL dbcsr_convert_dbcsr_to_csr(matrix_nosym, mat_csr)
     397           4 :                   CALL dbcsr_release(matrix_nosym)
     398             :                ELSE
     399             :                   ! FT of the matrix
     400          64 :                   CALL dbcsr_set(imatrix, 0.0_dp)
     401             :                   CALL rskp_transform(rmatrix=rmatrix, cmatrix=imatrix, rsmat=mat, ispin=ispin, &
     402          64 :                                       xkp=xkp(1:3, ik), cell_to_index=kpoints%cell_to_index, sab_nl=sab_nl)
     403             : 
     404             :                   ! Desymmetrize and sum the real and imaginary part into
     405             :                   ! cmatrix
     406          64 :                   CALL dbcsr_desymmetrize(rmatrix, rmatrix_nosym)
     407          64 :                   CALL dbcsr_desymmetrize(imatrix, imatrix_nosym)
     408             :                   ! Convert to csr
     409             :                   CALL dbcsr_csr_create_and_convert_complex(rmatrix=rmatrix_nosym, &
     410             :                                                             imatrix=imatrix_nosym, &
     411             :                                                             csr_mat=mat_csr, &
     412          64 :                                                             dist_format=dbcsr_csr_dbcsr_blkrow_dist)
     413             :                END IF
     414             :                ! Write to file
     415          68 :                WRITE (file_name, '(2(A,I0))') prefix//"_SPIN_", ispin, "_K_", ik
     416             :                unit_nr = cp_print_key_unit_nr(logger, dft_section, subs_string, &
     417             :                                               extension=".csr", middle_name=TRIM(file_name), &
     418          68 :                                               file_status="REPLACE", file_form=fileformat)
     419          68 :                CALL dbcsr_csr_write(mat_csr, unit_nr, upper_triangle=uptr, threshold=thld, binary=bin)
     420             : 
     421          68 :                CALL cp_print_key_finished_output(unit_nr, logger, dft_section, subs_string)
     422             : 
     423         104 :                CALL dbcsr_csr_destroy(mat_csr)
     424             :             END DO
     425             :          END DO
     426             :       END DO
     427          20 :       CALL dbcsr_release(rmatrix)
     428          20 :       DEALLOCATE (rmatrix)
     429          20 :       IF (.NOT. use_real_wfn) THEN
     430          16 :          CALL dbcsr_release(rmatrix_nosym)
     431          16 :          CALL dbcsr_release(imatrix)
     432          16 :          CALL dbcsr_release(imatrix_nosym)
     433          16 :          DEALLOCATE (rmatrix_nosym, imatrix, imatrix_nosym)
     434             :       END IF
     435          20 :       CALL timestop(handle)
     436             : 
     437          20 :    END SUBROUTINE write_matrix_kp_csr
     438             : 
     439             : ! **************************************************************************************************
     440             : !> \brief Desymmetrizes the KS or S matrices which are stored in symmetric !matrices
     441             : !> \param mat Hamiltonian or overlap matrices
     442             : !> \param mat_nosym Desymmetrized Hamiltonian or overlap matrices
     443             : !> \param cell_to_index Mapping of cell indices to linear RS indices
     444             : !> \param index_to_cell Mapping of linear RS indices to cell indices
     445             : !> \param kpoints Kpoint environment
     446             : !> \author Fabian Ducry
     447             : ! **************************************************************************************************
     448           2 :    SUBROUTINE desymmetrize_rs_matrix(mat, mat_nosym, cell_to_index, index_to_cell, kpoints)
     449             :       TYPE(dbcsr_p_type), DIMENSION(:, :), INTENT(IN), &
     450             :          POINTER                                         :: mat
     451             :       TYPE(dbcsr_p_type), DIMENSION(:, :), &
     452             :          INTENT(INOUT), POINTER                          :: mat_nosym
     453             :       INTEGER, ALLOCATABLE, DIMENSION(:, :, :), &
     454             :          INTENT(OUT)                                     :: cell_to_index
     455             :       INTEGER, ALLOCATABLE, DIMENSION(:, :), INTENT(OUT) :: index_to_cell
     456             :       TYPE(kpoint_type), INTENT(IN), POINTER             :: kpoints
     457             : 
     458             :       CHARACTER(len=*), PARAMETER :: routineN = 'desymmetrize_rs_matrix'
     459             : 
     460             :       INTEGER                                            :: handle, iatom, ic, icn, icol, irow, &
     461             :                                                             ispin, jatom, ncell, nomirror, nspin, &
     462             :                                                             nx, ny, nz
     463             :       INTEGER, DIMENSION(3)                              :: cell
     464           2 :       INTEGER, DIMENSION(:, :), POINTER                  :: i2c
     465           2 :       INTEGER, DIMENSION(:, :, :), POINTER               :: c2i
     466             :       LOGICAL                                            :: found, lwtr
     467           2 :       REAL(KIND=dp), DIMENSION(:, :), POINTER            :: block
     468             :       TYPE(neighbor_list_iterator_p_type), &
     469           2 :          DIMENSION(:), POINTER                           :: nl_iterator
     470             :       TYPE(neighbor_list_set_p_type), DIMENSION(:), &
     471           2 :          POINTER                                         :: sab_nl
     472             : 
     473           2 :       CALL timeset(routineN, handle)
     474             : 
     475           2 :       i2c => kpoints%index_to_cell
     476           2 :       c2i => kpoints%cell_to_index
     477             : 
     478           2 :       ncell = SIZE(i2c, 2)
     479           2 :       nspin = SIZE(mat, 1)
     480             : 
     481           2 :       nx = MAX(ABS(LBOUND(c2i, 1)), ABS(UBOUND(c2i, 1)))
     482           2 :       ny = MAX(ABS(LBOUND(c2i, 2)), ABS(UBOUND(c2i, 3)))
     483           2 :       nz = MAX(ABS(LBOUND(c2i, 3)), ABS(UBOUND(c2i, 3)))
     484          10 :       ALLOCATE (cell_to_index(-nx:nx, -ny:ny, -nz:nz))
     485             :       cell_to_index(LBOUND(c2i, 1):UBOUND(c2i, 1), &
     486             :                     LBOUND(c2i, 2):UBOUND(c2i, 2), &
     487          84 :                     LBOUND(c2i, 3):UBOUND(c2i, 3)) = c2i
     488             : 
     489             :       ! identify cells with no mirror img
     490             :       nomirror = 0
     491          52 :       DO ic = 1, ncell
     492         200 :          cell = i2c(:, ic)
     493          50 :          IF (cell_to_index(-cell(1), -cell(2), -cell(3)) == 0) &
     494           6 :             nomirror = nomirror + 1
     495             :       END DO
     496             : 
     497             :       ! create the mirror imgs
     498           6 :       ALLOCATE (index_to_cell(3, ncell + nomirror))
     499         202 :       index_to_cell(:, 1:ncell) = i2c
     500             : 
     501             :       nomirror = 0 ! count the imgs without mirror
     502          52 :       DO ic = 1, ncell
     503         200 :          cell = index_to_cell(:, ic)
     504          52 :          IF (cell_to_index(-cell(1), -cell(2), -cell(3)) == 0) THEN
     505           4 :             nomirror = nomirror + 1
     506          16 :             index_to_cell(:, ncell + nomirror) = -cell
     507           4 :             cell_to_index(-cell(1), -cell(2), -cell(3)) = ncell + nomirror
     508             :          END IF
     509             :       END DO
     510           2 :       ncell = ncell + nomirror
     511             : 
     512           2 :       CALL get_kpoint_info(kpoints, sab_nl=sab_nl)
     513             :       ! allocate the nonsymmetric matrices
     514           2 :       NULLIFY (mat_nosym)
     515           2 :       CALL dbcsr_allocate_matrix_set(mat_nosym, nspin, ncell)
     516           4 :       DO ispin = 1, nspin
     517          58 :          DO ic = 1, ncell
     518          54 :             ALLOCATE (mat_nosym(ispin, ic)%matrix)
     519             :             CALL dbcsr_create(matrix=mat_nosym(ispin, ic)%matrix, &
     520             :                               template=mat(1, 1)%matrix, &
     521             :                               matrix_type=dbcsr_type_no_symmetry, &
     522          54 :                               data_type=dbcsr_type_real_8)
     523             :             CALL cp_dbcsr_alloc_block_from_nbl(mat_nosym(ispin, ic)%matrix, &
     524          54 :                                                sab_nl, desymmetrize=.TRUE.)
     525          56 :             CALL dbcsr_set(mat_nosym(ispin, ic)%matrix, 0.0_dp)
     526             :          END DO
     527             :       END DO
     528             : 
     529           4 :       DO ispin = 1, nspin
     530             :          ! desymmetrize the matrix for real space printing
     531           2 :          CALL neighbor_list_iterator_create(nl_iterator, sab_nl)
     532         506 :          DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
     533         504 :             CALL get_iterator_info(nl_iterator, iatom=iatom, jatom=jatom, cell=cell)
     534             : 
     535         504 :             ic = cell_to_index(cell(1), cell(2), cell(3))
     536         504 :             icn = cell_to_index(-cell(1), -cell(2), -cell(3))
     537         504 :             CPASSERT(icn > 0)
     538             : 
     539         504 :             irow = iatom
     540         504 :             icol = jatom
     541         504 :             lwtr = .FALSE.
     542             :             ! always copy from the top
     543         504 :             IF (iatom > jatom) THEN
     544         216 :                irow = jatom
     545         216 :                icol = iatom
     546         216 :                lwtr = .TRUE.
     547             :             END IF
     548             : 
     549             :             CALL dbcsr_get_block_p(matrix=mat(ispin, ic)%matrix, &
     550         504 :                                    row=irow, col=icol, block=block, found=found)
     551         504 :             CPASSERT(found)
     552             : 
     553             :             ! copy to M(R) at (iatom,jatom)
     554             :             ! copy to M(-R) at (jatom,iatom)
     555         506 :             IF (lwtr) THEN
     556             :                CALL dbcsr_put_block(matrix=mat_nosym(ispin, ic)%matrix, &
     557        4536 :                                     row=iatom, col=jatom, block=TRANSPOSE(block))
     558             :                CALL dbcsr_put_block(matrix=mat_nosym(ispin, icn)%matrix, &
     559         216 :                                     row=jatom, col=iatom, block=block)
     560             :             ELSE
     561             :                CALL dbcsr_put_block(matrix=mat_nosym(ispin, ic)%matrix, &
     562         288 :                                     row=iatom, col=jatom, block=block)
     563             :                CALL dbcsr_put_block(matrix=mat_nosym(ispin, icn)%matrix, &
     564        6048 :                                     row=jatom, col=iatom, block=TRANSPOSE(block))
     565             :             END IF
     566             :          END DO
     567           4 :          CALL neighbor_list_iterator_release(nl_iterator)
     568             :       END DO
     569             : 
     570           4 :       DO ispin = 1, nspin
     571          58 :          DO ic = 1, ncell
     572          56 :             CALL dbcsr_finalize(mat_nosym(ispin, ic)%matrix)
     573             :          END DO
     574             :       END DO
     575             : 
     576           2 :       CALL timestop(handle)
     577             : 
     578           2 :    END SUBROUTINE desymmetrize_rs_matrix
     579             : 
     580             : END MODULE qs_scf_csr_write

Generated by: LCOV version 1.15