LCOV - code coverage report
Current view: top level - src - qs_scf_wfn_mix.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:42dac4a) Lines: 90.7 % 107 97
Test Date: 2025-07-25 12:55:17 Functions: 100.0 % 1 1

            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 Does all kind of post scf calculations for GPW/GAPW
      10              : !> \par History
      11              : !>      Started as a copy from the relevant part of qs_scf
      12              : !> \author Joost VandeVondele (10.2003)
      13              : ! **************************************************************************************************
      14              : MODULE qs_scf_wfn_mix
      15              :    USE cp_dbcsr_api,                    ONLY: dbcsr_p_type
      16              :    USE cp_dbcsr_operations,             ONLY: copy_fm_to_dbcsr
      17              :    USE cp_files,                        ONLY: close_file,&
      18              :                                               open_file
      19              :    USE cp_fm_basic_linalg,              ONLY: cp_fm_scale_and_add
      20              :    USE cp_fm_struct,                    ONLY: cp_fm_struct_create,&
      21              :                                               cp_fm_struct_release,&
      22              :                                               cp_fm_struct_type
      23              :    USE cp_fm_types,                     ONLY: cp_fm_create,&
      24              :                                               cp_fm_release,&
      25              :                                               cp_fm_to_fm,&
      26              :                                               cp_fm_type
      27              :    USE input_constants,                 ONLY: wfn_mix_orig_external,&
      28              :                                               wfn_mix_orig_occ,&
      29              :                                               wfn_mix_orig_virtual
      30              :    USE input_section_types,             ONLY: section_vals_get,&
      31              :                                               section_vals_get_subs_vals,&
      32              :                                               section_vals_type,&
      33              :                                               section_vals_val_get
      34              :    USE kinds,                           ONLY: default_path_length,&
      35              :                                               dp
      36              :    USE message_passing,                 ONLY: mp_para_env_type
      37              :    USE particle_types,                  ONLY: particle_type
      38              :    USE qs_kind_types,                   ONLY: qs_kind_type
      39              :    USE qs_mo_io,                        ONLY: read_mos_restart_low,&
      40              :                                               write_mo_set_to_restart
      41              :    USE qs_mo_methods,                   ONLY: calculate_orthonormality
      42              :    USE qs_mo_types,                     ONLY: deallocate_mo_set,&
      43              :                                               duplicate_mo_set,&
      44              :                                               get_mo_set,&
      45              :                                               mo_set_type
      46              :    USE qs_scf_types,                    ONLY: qs_scf_env_type,&
      47              :                                               special_diag_method_nr
      48              : #include "./base/base_uses.f90"
      49              : 
      50              :    IMPLICIT NONE
      51              :    PRIVATE
      52              : 
      53              :    ! Global parameters
      54              :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_scf_wfn_mix'
      55              :    PUBLIC :: wfn_mix
      56              : 
      57              : CONTAINS
      58              : 
      59              : ! **************************************************************************************************
      60              : !> \brief writes a new 'mixed' set of mos to restart file, without touching the current MOs
      61              : !> \param mos ...
      62              : !> \param particle_set ...
      63              : !> \param dft_section ...
      64              : !> \param qs_kind_set ...
      65              : !> \param para_env ...
      66              : !> \param output_unit ...
      67              : !> \param unoccupied_orbs ...
      68              : !> \param scf_env ...
      69              : !> \param matrix_s ...
      70              : !> \param marked_states ...
      71              : !> \param for_rtp ...
      72              : ! **************************************************************************************************
      73         9687 :    SUBROUTINE wfn_mix(mos, particle_set, dft_section, qs_kind_set, para_env, output_unit, &
      74              :                       unoccupied_orbs, scf_env, matrix_s, marked_states, for_rtp)
      75              : 
      76              :       TYPE(mo_set_type), DIMENSION(:), POINTER           :: mos
      77              :       TYPE(particle_type), DIMENSION(:), POINTER         :: particle_set
      78              :       TYPE(section_vals_type), POINTER                   :: dft_section
      79              :       TYPE(qs_kind_type), DIMENSION(:), POINTER          :: qs_kind_set
      80              :       TYPE(mp_para_env_type), POINTER                    :: para_env
      81              :       INTEGER                                            :: output_unit
      82              :       TYPE(cp_fm_type), DIMENSION(:), INTENT(IN), &
      83              :          OPTIONAL, POINTER                               :: unoccupied_orbs
      84              :       TYPE(qs_scf_env_type), OPTIONAL, POINTER           :: scf_env
      85              :       TYPE(dbcsr_p_type), DIMENSION(:), OPTIONAL, &
      86              :          POINTER                                         :: matrix_s
      87              :       INTEGER, DIMENSION(:, :, :), OPTIONAL, POINTER     :: marked_states
      88              :       LOGICAL, OPTIONAL                                  :: for_rtp
      89              : 
      90              :       CHARACTER(len=*), PARAMETER                        :: routineN = 'wfn_mix'
      91              : 
      92              :       CHARACTER(LEN=default_path_length)                 :: read_file_name
      93              :       INTEGER :: handle, homo, i_rep, ispin, mark_ind, mark_number, n_rep, nmo, orig_mo_index, &
      94              :          orig_spin_index, orig_type, restart_unit, result_mo_index, result_spin_index
      95              :       LOGICAL                                            :: explicit, is_file, my_for_rtp, &
      96              :                                                             overwrite_mos, reverse_mo_index
      97              :       REAL(KIND=dp)                                      :: orig_scale, orthonormality, result_scale
      98              :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct_vector
      99              :       TYPE(cp_fm_type)                                   :: matrix_x, matrix_y
     100         9687 :       TYPE(mo_set_type), ALLOCATABLE, DIMENSION(:)       :: mos_new, mos_orig_ext
     101              :       TYPE(section_vals_type), POINTER                   :: update_section, wfn_mix_section
     102              : 
     103         9687 :       CALL timeset(routineN, handle)
     104         9687 :       wfn_mix_section => section_vals_get_subs_vals(dft_section, "PRINT%WFN_MIX")
     105         9687 :       CALL section_vals_get(wfn_mix_section, explicit=explicit)
     106              : 
     107              :       ! only perform action if explicitly required
     108         9687 :       IF (explicit) THEN
     109              : 
     110           24 :          IF (PRESENT(for_rtp)) THEN
     111            4 :             my_for_rtp = for_rtp
     112              :          ELSE
     113              :             my_for_rtp = .FALSE.
     114              :          END IF
     115              : 
     116           24 :          IF (output_unit > 0) THEN
     117           12 :             WRITE (output_unit, '()')
     118           12 :             WRITE (output_unit, '(T2,A)') "Performing wfn mixing"
     119           12 :             WRITE (output_unit, '(T2,A)') "====================="
     120              :          END IF
     121              : 
     122          112 :          ALLOCATE (mos_new(SIZE(mos)))
     123           64 :          DO ispin = 1, SIZE(mos)
     124           64 :             CALL duplicate_mo_set(mos_new(ispin), mos(ispin))
     125              :          END DO
     126              : 
     127              :          ! a single vector matrix structure
     128           24 :          NULLIFY (fm_struct_vector)
     129              :          CALL cp_fm_struct_create(fm_struct_vector, template_fmstruct=mos(1)%mo_coeff%matrix_struct, &
     130           24 :                                   ncol_global=1)
     131           24 :          CALL cp_fm_create(matrix_x, fm_struct_vector, name="x")
     132           24 :          CALL cp_fm_create(matrix_y, fm_struct_vector, name="y")
     133           24 :          CALL cp_fm_struct_release(fm_struct_vector)
     134              : 
     135           24 :          update_section => section_vals_get_subs_vals(wfn_mix_section, "UPDATE")
     136           24 :          CALL section_vals_get(update_section, n_repetition=n_rep)
     137           24 :          CALL section_vals_get(update_section, explicit=explicit)
     138           24 :          IF (.NOT. explicit) n_rep = 0
     139              : 
     140              :          ! Mix the MOs as : y = ay + bx
     141           54 :          DO i_rep = 1, n_rep
     142              :             ! The occupied MO that will be modified or saved, 'y'
     143           30 :             CALL section_vals_val_get(update_section, "RESULT_MO_INDEX", i_rep_section=i_rep, i_val=result_mo_index)
     144           30 :             CALL section_vals_val_get(update_section, "RESULT_MARKED_STATE", i_rep_section=i_rep, i_val=mark_number)
     145           30 :             CALL section_vals_val_get(update_section, "RESULT_SPIN_INDEX", i_rep_section=i_rep, i_val=result_spin_index)
     146              :             ! result_scale is the 'a' coefficient
     147           30 :             CALL section_vals_val_get(update_section, "RESULT_SCALE", i_rep_section=i_rep, r_val=result_scale)
     148              : 
     149           30 :             mark_ind = 1
     150           30 :             IF (mark_number .GT. 0) result_mo_index = marked_states(mark_number, result_spin_index, mark_ind)
     151              : 
     152              :             ! The MO that will be added to the previous one, 'x'
     153              :             CALL section_vals_val_get(update_section, "ORIG_TYPE", i_rep_section=i_rep, &
     154           30 :                                       i_val=orig_type)
     155           30 :             CALL section_vals_val_get(update_section, "ORIG_MO_INDEX", i_rep_section=i_rep, i_val=orig_mo_index)
     156           30 :             CALL section_vals_val_get(update_section, "ORIG_MARKED_STATE", i_rep_section=i_rep, i_val=mark_number)
     157           30 :             CALL section_vals_val_get(update_section, "ORIG_SPIN_INDEX", i_rep_section=i_rep, i_val=orig_spin_index)
     158              :             ! orig_scale is the 'b' coefficient
     159           30 :             CALL section_vals_val_get(update_section, "ORIG_SCALE", i_rep_section=i_rep, r_val=orig_scale)
     160              : 
     161           30 :             IF (orig_type == wfn_mix_orig_virtual) mark_ind = 2
     162           30 :             IF (mark_number .GT. 0) orig_mo_index = marked_states(mark_number, orig_spin_index, mark_ind)
     163              : 
     164           30 :             CALL section_vals_val_get(wfn_mix_section, "OVERWRITE_MOS", l_val=overwrite_mos)
     165              : 
     166           30 :             CALL section_vals_val_get(update_section, "REVERSE_MO_INDEX", i_rep_section=i_rep, l_val=reverse_mo_index)
     167              : 
     168              :             ! First get a copy of the proper orig
     169              :             ! Origin is in an occupied MO
     170           30 :             IF (orig_type == wfn_mix_orig_occ) THEN
     171            4 :                IF (reverse_mo_index) THEN
     172              :                   CALL cp_fm_to_fm(mos(orig_spin_index)%mo_coeff, matrix_x, 1, &
     173            0 :                                    orig_mo_index, 1)
     174              :                ELSE
     175              :                   CALL cp_fm_to_fm(mos(orig_spin_index)%mo_coeff, matrix_x, 1, &
     176            4 :                                    mos(orig_spin_index)%nmo - orig_mo_index + 1, 1)
     177              :                END IF
     178              :                ! Origin is an unoccupied MO
     179           26 :             ELSE IF (orig_type == wfn_mix_orig_virtual) THEN
     180           22 :                CALL get_mo_set(mos(orig_spin_index), homo=homo, nmo=nmo)
     181              :                ! check whether the MO is already available in the original MO set
     182           22 :                IF (orig_mo_index + homo <= nmo) THEN
     183            0 :                   CALL cp_fm_to_fm(mos(orig_spin_index)%mo_coeff, matrix_x, 1, orig_mo_index + homo, 1)
     184           22 :                ELSE IF (.NOT. ASSOCIATED(unoccupied_orbs)) THEN
     185              :                   CALL cp_abort(__LOCATION__, &
     186              :                                 "If ORIG_TYPE is set to VIRTUAL, the array unoccupied_orbs must be associated! "// &
     187            0 :                                 "For instance, ask in the SCF section to compute virtual orbitals after the GS optimization.")
     188              :                ELSE
     189           22 :                   CALL cp_fm_to_fm(unoccupied_orbs(orig_spin_index), matrix_x, 1, orig_mo_index, 1)
     190              :                END IF
     191              : 
     192              :                ! Origin is to be read from an external .wfn file
     193            4 :             ELSE IF (orig_type == wfn_mix_orig_external) THEN
     194              :                CALL section_vals_val_get(update_section, "ORIG_EXT_FILE_NAME", i_rep_section=i_rep, &
     195            4 :                                          c_val=read_file_name)
     196            4 :                IF (read_file_name == "EMPTY") &
     197              :                   CALL cp_abort(__LOCATION__, &
     198              :                                 "If ORIG_TYPE is set to EXTERNAL, a file name should be set in ORIG_EXT_FILE_NAME "// &
     199            0 :                                 "so that it can be used as the orginal MO.")
     200              : 
     201           18 :                ALLOCATE (mos_orig_ext(SIZE(mos)))
     202           10 :                DO ispin = 1, SIZE(mos)
     203           10 :                   CALL duplicate_mo_set(mos_orig_ext(ispin), mos(ispin))
     204              :                END DO
     205              : 
     206            4 :                IF (para_env%is_source()) THEN
     207            2 :                   INQUIRE (FILE=TRIM(read_file_name), exist=is_file)
     208            2 :                   IF (.NOT. is_file) &
     209              :                      CALL cp_abort(__LOCATION__, &
     210            0 :                                    "Reference file not found! Name of the file CP2K looked for: "//TRIM(read_file_name))
     211              : 
     212              :                   CALL open_file(file_name=read_file_name, &
     213              :                                  file_action="READ", &
     214              :                                  file_form="UNFORMATTED", &
     215              :                                  file_status="OLD", &
     216            2 :                                  unit_number=restart_unit)
     217              :                END IF
     218              :                CALL read_mos_restart_low(mos_orig_ext, para_env=para_env, qs_kind_set=qs_kind_set, &
     219              :                                          particle_set=particle_set, natom=SIZE(particle_set, 1), &
     220            4 :                                          rst_unit=restart_unit)
     221            4 :                IF (para_env%is_source()) CALL close_file(unit_number=restart_unit)
     222              : 
     223            4 :                IF (reverse_mo_index) THEN
     224              :                   CALL cp_fm_to_fm(mos_orig_ext(orig_spin_index)%mo_coeff, matrix_x, 1, &
     225            0 :                                    orig_mo_index, 1)
     226              :                ELSE
     227              :                   CALL cp_fm_to_fm(mos_orig_ext(orig_spin_index)%mo_coeff, matrix_x, 1, &
     228            4 :                                    mos_orig_ext(orig_spin_index)%nmo - orig_mo_index + 1, 1)
     229              :                END IF
     230           10 :                DO ispin = 1, SIZE(mos_orig_ext)
     231           10 :                   CALL deallocate_mo_set(mos_orig_ext(ispin))
     232              :                END DO
     233            4 :                DEALLOCATE (mos_orig_ext)
     234              :             END IF
     235              : 
     236              :             ! Second, get a copy of the target
     237           30 :             IF (reverse_mo_index) THEN
     238              :                CALL cp_fm_to_fm(mos_new(result_spin_index)%mo_coeff, matrix_y, &
     239            0 :                                 1, result_mo_index, 1)
     240              :             ELSE
     241              :                CALL cp_fm_to_fm(mos_new(result_spin_index)%mo_coeff, matrix_y, &
     242           30 :                                 1, mos_new(result_spin_index)%nmo - result_mo_index + 1, 1)
     243              :             END IF
     244              : 
     245              :             ! Third, perform the mix
     246           30 :             CALL cp_fm_scale_and_add(result_scale, matrix_y, orig_scale, matrix_x)
     247              : 
     248              :             ! and copy back in the result mos
     249          144 :             IF (reverse_mo_index) THEN
     250              :                CALL cp_fm_to_fm(matrix_y, mos_new(result_spin_index)%mo_coeff, &
     251            0 :                                 1, 1, result_mo_index)
     252              :             ELSE
     253              :                CALL cp_fm_to_fm(matrix_y, mos_new(result_spin_index)%mo_coeff, &
     254           30 :                                 1, 1, mos_new(result_spin_index)%nmo - result_mo_index + 1)
     255              :             END IF
     256              :          END DO
     257              : 
     258           24 :          CALL cp_fm_release(matrix_x)
     259           24 :          CALL cp_fm_release(matrix_y)
     260              : 
     261           24 :          IF (my_for_rtp) THEN
     262           12 :             DO ispin = 1, SIZE(mos_new)
     263            8 :                CALL cp_fm_to_fm(mos_new(ispin)%mo_coeff, mos(ispin)%mo_coeff)
     264            8 :                IF (mos_new(1)%use_mo_coeff_b) &
     265            0 :                   CALL copy_fm_to_dbcsr(mos_new(ispin)%mo_coeff, mos_new(ispin)%mo_coeff_b)
     266            8 :                IF (mos(1)%use_mo_coeff_b) &
     267            4 :                   CALL copy_fm_to_dbcsr(mos_new(ispin)%mo_coeff, mos(ispin)%mo_coeff_b)
     268              :             END DO
     269              :          ELSE
     270           20 :             IF (scf_env%method == special_diag_method_nr) THEN
     271            0 :                CALL calculate_orthonormality(orthonormality, mos)
     272              :             ELSE
     273           20 :                CALL calculate_orthonormality(orthonormality, mos, matrix_s(1)%matrix)
     274              :             END IF
     275              : 
     276           20 :             IF (output_unit > 0) THEN
     277           10 :                WRITE (output_unit, '()')
     278              :                WRITE (output_unit, '(T2,A,T61,E20.4)') &
     279           10 :                   "Maximum deviation from MO S-orthonormality", orthonormality
     280           10 :                WRITE (output_unit, '(T2,A)') "Writing new MOs to file"
     281              :             END IF
     282              : 
     283              :             ! *** Write WaveFunction restart file ***
     284              : 
     285           52 :             DO ispin = 1, SIZE(mos_new)
     286           32 :                IF (overwrite_mos) THEN
     287           12 :                   CALL cp_fm_to_fm(mos_new(ispin)%mo_coeff, mos(ispin)%mo_coeff)
     288           12 :                   IF (mos_new(1)%use_mo_coeff_b) &
     289            4 :                      CALL copy_fm_to_dbcsr(mos_new(ispin)%mo_coeff, mos_new(ispin)%mo_coeff_b)
     290              :                END IF
     291           32 :                IF (mos(1)%use_mo_coeff_b) &
     292           24 :                   CALL copy_fm_to_dbcsr(mos_new(ispin)%mo_coeff, mos(ispin)%mo_coeff_b)
     293              :             END DO
     294           20 :             CALL write_mo_set_to_restart(mos_new, particle_set, dft_section=dft_section, qs_kind_set=qs_kind_set)
     295              :          END IF
     296              : 
     297           64 :          DO ispin = 1, SIZE(mos_new)
     298           64 :             CALL deallocate_mo_set(mos_new(ispin))
     299              :          END DO
     300           72 :          DEALLOCATE (mos_new)
     301              : 
     302              :       END IF
     303              : 
     304         9687 :       CALL timestop(handle)
     305              : 
     306        19374 :    END SUBROUTINE wfn_mix
     307              : 
     308              : END MODULE qs_scf_wfn_mix
        

Generated by: LCOV version 2.0-1