Line data    Source code

       1             : !--------------------------------------------------------------------------------------------------!
       2             : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3             : !   Copyright 2000-2024 CP2K developers group <https://cp2k.org>                                   !
       4             : !                                                                                                  !
       5             : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6             : !--------------------------------------------------------------------------------------------------!
       7             : MODULE optbas_fenv_manipulation
       8             :    USE atomic_kind_types,               ONLY: atomic_kind_type,&
       9             :                                               get_atomic_kind
      10             :    USE basis_set_container_types,       ONLY: get_basis_from_container
      11             :    USE basis_set_types,                 ONLY: gto_basis_set_type,&
      12             :                                               init_orb_basis_set
      13             :    USE cp_blacs_env,                    ONLY: cp_blacs_env_type
      14             :    USE cp_control_types,                ONLY: dft_control_type
      15             :    USE cp_dbcsr_operations,             ONLY: copy_dbcsr_to_fm
      16             :    USE cp_fm_basic_linalg,              ONLY: cp_fm_upper_to_full
      17             :    USE cp_fm_cholesky,                  ONLY: cp_fm_cholesky_decompose,&
      18             :                                               cp_fm_cholesky_invert
      19             :    USE cp_fm_pool_types,                ONLY: cp_fm_pool_p_type
      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_type
      26             :    USE cp_log_handling,                 ONLY: cp_to_string
      27             :    USE cp_output_handling,              ONLY: debug_print_level
      28             :    USE dbcsr_api,                       ONLY: dbcsr_get_info,&
      29             :                                               dbcsr_p_type,&
      30             :                                               dbcsr_type
      31             :    USE input_section_types,             ONLY: section_vals_get,&
      32             :                                               section_vals_get_subs_vals,&
      33             :                                               section_vals_type,&
      34             :                                               section_vals_val_get,&
      35             :                                               section_vals_val_set
      36             :    USE kinds,                           ONLY: default_string_length
      37             :    USE message_passing,                 ONLY: mp_para_env_type
      38             :    USE optimize_basis_types,            ONLY: basis_optimization_type,&
      39             :                                               flex_basis_type
      40             :    USE particle_types,                  ONLY: particle_type
      41             :    USE qs_density_matrices,             ONLY: calculate_density_matrix
      42             :    USE qs_energy_init,                  ONLY: qs_energies_init
      43             :    USE qs_environment_types,            ONLY: get_qs_env,&
      44             :                                               qs_environment_type
      45             :    USE qs_interactions,                 ONLY: init_interaction_radii
      46             :    USE qs_kind_types,                   ONLY: qs_kind_type
      47             :    USE qs_ks_methods,                   ONLY: qs_ks_update_qs_env
      48             :    USE qs_ks_types,                     ONLY: qs_ks_did_change
      49             :    USE qs_matrix_pools,                 ONLY: mpools_get
      50             :    USE qs_mo_io,                        ONLY: read_mo_set_from_restart
      51             :    USE qs_mo_types,                     ONLY: init_mo_set,&
      52             :                                               mo_set_type
      53             :    USE qs_rho_methods,                  ONLY: qs_rho_update_rho
      54             :    USE qs_rho_types,                    ONLY: qs_rho_get,&
      55             :                                               qs_rho_type
      56             :    USE string_utilities,                ONLY: uppercase
      57             : #include "./base/base_uses.f90"
      58             : 
      59             :    IMPLICIT NONE
      60             :    PRIVATE
      61             : 
      62             :    PUBLIC :: modify_input_settings, &
      63             :              allocate_mo_sets, &
      64             :              update_basis_set, &
      65             :              calculate_ks_matrix, &
      66             :              calculate_overlap_inverse
      67             : 
      68             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'optbas_fenv_manipulation'
      69             : 
      70             : CONTAINS
      71             : 
      72             : ! **************************************************************************************************
      73             : !> \brief change settings in the training input files to initialize
      74             : !>        all needed structures and adjust settings to basis optimization
      75             : !> \param basis_optimization ...
      76             : !> \param bas_id ...
      77             : !> \param input_file ...
      78             : !> \author Florian Schiffmann
      79             : ! **************************************************************************************************
      80          27 :    SUBROUTINE modify_input_settings(basis_optimization, bas_id, input_file)
      81             :       TYPE(basis_optimization_type)                      :: basis_optimization
      82             :       INTEGER                                            :: bas_id
      83             :       TYPE(section_vals_type), POINTER                   :: input_file
      84             : 
      85             :       CHARACTER(LEN=default_string_length)               :: atom
      86             :       CHARACTER(LEN=default_string_length), &
      87           9 :          DIMENSION(:), POINTER                           :: abasinfo, obasinfo
      88             :       INTEGER                                            :: ibasis, ikind, jkind, nbasis, nkind
      89             :       TYPE(section_vals_type), POINTER                   :: dft_section, feval_section, &
      90             :                                                             kind_section, subsys_section
      91             : 
      92          18 :       feval_section => section_vals_get_subs_vals(input_file, "FORCE_EVAL")
      93           9 :       dft_section => section_vals_get_subs_vals(feval_section, "DFT")
      94           9 :       subsys_section => section_vals_get_subs_vals(feval_section, "SUBSYS")
      95           9 :       kind_section => section_vals_get_subs_vals(subsys_section, "KIND")
      96             : 
      97             :       CALL section_vals_val_set(feval_section, "PRINT%DISTRIBUTION%_SECTION_PARAMETERS_", &
      98           9 :                                 i_val=debug_print_level)
      99             :       CALL section_vals_val_set(dft_section, "SCF%PRINT%TOTAL_DENSITIES%_SECTION_PARAMETERS_", &
     100           9 :                                 i_val=debug_print_level)
     101             :       CALL section_vals_val_set(dft_section, "SCF%PRINT%DETAILED_ENERGY%_SECTION_PARAMETERS_", &
     102           9 :                                 i_val=debug_print_level)
     103             : 
     104             :       ! add the new basis file containing the templates to the basis file list
     105           9 :       CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", n_rep_val=nbasis)
     106             :       CALL section_vals_val_set(dft_section, "BASIS_SET_FILE_NAME", i_rep_val=nbasis + 1, &
     107           9 :                                 c_val=basis_optimization%work_basis_file)
     108             : 
     109             :       ! Set the auxilarry basis in the kind sections
     110           9 :       CALL section_vals_get(kind_section, n_repetition=nkind)
     111          24 :       DO ikind = 1, nkind
     112             :          CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", &
     113          15 :                                    c_val=atom, i_rep_section=ikind)
     114          15 :          CALL uppercase(atom)
     115          15 :          CALL section_vals_val_get(kind_section, "BASIS_SET", n_rep_val=nbasis, i_rep_section=ikind)
     116          15 :          IF (nbasis > 1) THEN
     117             :             CALL cp_abort(__LOCATION__, &
     118           0 :                           "Basis set optimization: Only one single BASIS_SET allowed per KIND in the reference input")
     119             :          END IF
     120             :          CALL section_vals_val_get(kind_section, "BASIS_SET", &
     121          15 :                                    c_vals=obasinfo, i_rep_val=1, i_rep_section=ikind)
     122          15 :          ALLOCATE (abasinfo(2))
     123          15 :          abasinfo(1) = "AUX_OPT"
     124          15 :          IF (SIZE(obasinfo) == 1) THEN
     125          15 :             abasinfo(2) = obasinfo(1)
     126             :          ELSE
     127           0 :             abasinfo(2) = obasinfo(2)
     128             :          END IF
     129             :          CALL section_vals_val_set(kind_section, "BASIS_SET", &
     130          15 :                                    c_vals_ptr=abasinfo, i_rep_val=2, i_rep_section=ikind)
     131          15 :          CALL section_vals_val_get(kind_section, "BASIS_SET", n_rep_val=nbasis, i_rep_section=ikind)
     132          15 :          CPASSERT(nbasis == 2)
     133             : 
     134          60 :          DO jkind = 1, basis_optimization%nkind
     135          21 :             IF (atom == basis_optimization%kind_basis(jkind)%element) THEN
     136             : 
     137          15 :                NULLIFY (abasinfo)
     138             :                CALL section_vals_val_get(kind_section, "BASIS_SET", &
     139          15 :                                          c_vals=abasinfo, i_rep_val=2, i_rep_section=ikind)
     140          15 :                ibasis = basis_optimization%combination(bas_id, jkind)
     141          15 :                CPASSERT(SIZE(abasinfo) == 2)
     142          15 :                CPASSERT(abasinfo(1) == "AUX_OPT")
     143          15 :                abasinfo(2) = TRIM(ADJUSTL(basis_optimization%kind_basis(jkind)%flex_basis(ibasis)%basis_name))
     144          15 :                EXIT
     145             :             END IF
     146             :          END DO
     147             :       END DO
     148             : 
     149           9 :    END SUBROUTINE modify_input_settings
     150             : 
     151             : ! **************************************************************************************************
     152             : !> \brief ...
     153             : !> \param qs_env ...
     154             : ! **************************************************************************************************
     155           9 :    SUBROUTINE allocate_mo_sets(qs_env)
     156             :       TYPE(qs_environment_type), POINTER                 :: qs_env
     157             : 
     158             :       INTEGER                                            :: ispin
     159             :       INTEGER, DIMENSION(2)                              :: nelectron_spin
     160             :       LOGICAL                                            :: natom_mismatch
     161           9 :       TYPE(atomic_kind_type), DIMENSION(:), POINTER      :: atomic_kind_set
     162           9 :       TYPE(cp_fm_pool_p_type), DIMENSION(:), POINTER     :: ao_mo_fm_pools
     163             :       TYPE(dft_control_type), POINTER                    :: dft_control
     164           9 :       TYPE(mo_set_type), DIMENSION(:), POINTER           :: mos
     165             :       TYPE(mp_para_env_type), POINTER                    :: para_env
     166           9 :       TYPE(particle_type), DIMENSION(:), POINTER         :: particle_set
     167           9 :       TYPE(qs_kind_type), DIMENSION(:), POINTER          :: qs_kind_set
     168             :       TYPE(section_vals_type), POINTER                   :: dft_section
     169             : 
     170           9 :       NULLIFY (para_env)
     171             :       CALL get_qs_env(qs_env=qs_env, &
     172             :                       dft_control=dft_control, &
     173             :                       mos=mos, nelectron_spin=nelectron_spin, &
     174             :                       atomic_kind_set=atomic_kind_set, &
     175             :                       qs_kind_set=qs_kind_set, &
     176             :                       particle_set=particle_set, &
     177           9 :                       para_env=para_env)
     178           9 :       dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
     179             : 
     180           9 :       CALL mpools_get(qs_env%mpools, ao_mo_fm_pools=ao_mo_fm_pools)
     181          18 :       DO ispin = 1, dft_control%nspins
     182          18 :          IF (.NOT. ASSOCIATED(mos(ispin)%mo_coeff)) THEN
     183             :             CALL init_mo_set(mos(ispin), &
     184             :                              fm_pool=ao_mo_fm_pools(ispin)%pool, &
     185           9 :                              name="qs_env%mo"//TRIM(ADJUSTL(cp_to_string(ispin))))
     186             :          END IF
     187             :       END DO
     188             : 
     189             :       CALL read_mo_set_from_restart(mos, atomic_kind_set, qs_kind_set, particle_set, para_env, &
     190             :                                     id_nr=0, multiplicity=dft_control%multiplicity, dft_section=dft_section, &
     191           9 :                                     natom_mismatch=natom_mismatch)
     192             : 
     193           9 :    END SUBROUTINE allocate_mo_sets
     194             : 
     195             : ! **************************************************************************************************
     196             : !> \brief ...
     197             : !> \param qs_env ...
     198             : ! **************************************************************************************************
     199           9 :    SUBROUTINE calculate_ks_matrix(qs_env)
     200             :       TYPE(qs_environment_type), POINTER                 :: qs_env
     201             : 
     202             :       INTEGER                                            :: ispin
     203           9 :       TYPE(dbcsr_p_type), DIMENSION(:), POINTER          :: rho_ao
     204             :       TYPE(dft_control_type), POINTER                    :: dft_control
     205             :       TYPE(qs_rho_type), POINTER                         :: rho
     206             : 
     207           9 :       NULLIFY (rho, dft_control, rho_ao)
     208             : 
     209           9 :       CALL qs_energies_init(qs_env, .FALSE.)
     210           9 :       CALL get_qs_env(qs_env, rho=rho, dft_control=dft_control)
     211           9 :       CALL qs_rho_get(rho, rho_ao=rho_ao)
     212          18 :       DO ispin = 1, dft_control%nspins
     213          18 :          CALL calculate_density_matrix(qs_env%mos(ispin), rho_ao(ispin)%matrix)
     214             :       END DO
     215           9 :       CALL qs_rho_update_rho(rho, qs_env)
     216           9 :       CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.TRUE.)
     217           9 :       qs_env%requires_mo_derivs = .FALSE.
     218           9 :       CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE.)
     219             : 
     220           9 :    END SUBROUTINE calculate_ks_matrix
     221             : 
     222             : ! **************************************************************************************************
     223             : !> \brief ...
     224             : !> \param matrix_s ...
     225             : !> \param matrix_s_inv ...
     226             : !> \param para_env ...
     227             : !> \param context ...
     228             : ! **************************************************************************************************
     229           9 :    SUBROUTINE calculate_overlap_inverse(matrix_s, matrix_s_inv, para_env, context)
     230             :       TYPE(dbcsr_type), POINTER                          :: matrix_s
     231             :       TYPE(cp_fm_type), INTENT(OUT)                      :: matrix_s_inv
     232             :       TYPE(mp_para_env_type), POINTER                    :: para_env
     233             :       TYPE(cp_blacs_env_type), POINTER                   :: context
     234             : 
     235             :       INTEGER                                            :: nao
     236             :       TYPE(cp_fm_struct_type), POINTER                   :: fm_struct_tmp
     237             :       TYPE(cp_fm_type)                                   :: work1
     238             : 
     239           9 :       CALL dbcsr_get_info(matrix_s, nfullrows_total=nao)
     240             :       CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=nao, ncol_global=nao, &
     241           9 :                                para_env=para_env, context=context)
     242             : 
     243           9 :       CALL cp_fm_create(matrix_s_inv, matrix_struct=fm_struct_tmp)
     244           9 :       CALL cp_fm_create(work1, matrix_struct=fm_struct_tmp)
     245           9 :       CALL copy_dbcsr_to_fm(matrix_s, matrix_s_inv)
     246           9 :       CALL cp_fm_upper_to_full(matrix_s_inv, work1)
     247           9 :       CALL cp_fm_cholesky_decompose(matrix_s_inv)
     248           9 :       CALL cp_fm_cholesky_invert(matrix_s_inv)
     249           9 :       CALL cp_fm_upper_to_full(matrix_s_inv, work1)
     250           9 :       CALL cp_fm_struct_release(fm_struct_tmp)
     251           9 :       CALL cp_fm_release(work1)
     252             : 
     253           9 :    END SUBROUTINE calculate_overlap_inverse
     254             : 
     255             : ! **************************************************************************************************
     256             : !> \brief ...
     257             : !> \param opt_bas ...
     258             : !> \param bas_id ...
     259             : !> \param basis_type ...
     260             : !> \param qs_env ...
     261             : ! **************************************************************************************************
     262         234 :    SUBROUTINE update_basis_set(opt_bas, bas_id, basis_type, qs_env)
     263             :       TYPE(basis_optimization_type)                      :: opt_bas
     264             :       INTEGER                                            :: bas_id
     265             :       CHARACTER(*)                                       :: basis_type
     266             :       TYPE(qs_environment_type), POINTER                 :: qs_env
     267             : 
     268             :       CHARACTER(default_string_length)                   :: elem
     269             :       INTEGER                                            :: ibasis, ikind, jkind
     270         234 :       TYPE(atomic_kind_type), DIMENSION(:), POINTER      :: atomic_kind_set
     271             :       TYPE(dft_control_type), POINTER                    :: dft_control
     272             :       TYPE(gto_basis_set_type), POINTER                  :: gto_basis
     273         234 :       TYPE(qs_kind_type), DIMENSION(:), POINTER          :: qs_kind_set
     274             : 
     275             :       CALL get_qs_env(qs_env=qs_env, dft_control=dft_control, &
     276         234 :                       atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set)
     277         639 :       DO ikind = 1, SIZE(qs_kind_set)
     278        1449 :          DO jkind = 1, opt_bas%nkind
     279         810 :             CALL get_atomic_kind(atomic_kind_set(ikind), name=elem)
     280         810 :             CALL uppercase(elem)
     281        1215 :             IF (elem == opt_bas%kind_basis(jkind)%element) THEN
     282         405 :                ibasis = opt_bas%combination(bas_id, jkind)
     283             :                CALL get_basis_from_container(qs_kind_set(ikind)%basis_sets, basis_set=gto_basis, &
     284         405 :                                              basis_type=basis_type)
     285         405 :                CALL transfer_data_to_gto(gto_basis, opt_bas%kind_basis(jkind)%flex_basis(ibasis))
     286         405 :                CALL init_orb_basis_set(gto_basis)
     287             :             END IF
     288             :          END DO
     289             :       END DO
     290             : 
     291         234 :       CALL init_interaction_radii(dft_control%qs_control, qs_kind_set)
     292             : 
     293         234 :    END SUBROUTINE update_basis_set
     294             : 
     295             : ! **************************************************************************************************
     296             : !> \brief ...
     297             : !> \param gto_basis ...
     298             : !> \param basis ...
     299             : ! **************************************************************************************************
     300         405 :    SUBROUTINE transfer_data_to_gto(gto_basis, basis)
     301             :       TYPE(gto_basis_set_type), POINTER                  :: gto_basis
     302             :       TYPE(flex_basis_type)                              :: basis
     303             : 
     304             :       INTEGER                                            :: ipgf, iset, ishell
     305             : 
     306         810 :       DO iset = 1, basis%nsets
     307        1827 :          DO ishell = 1, basis%subset(iset)%ncon_tot
     308       11781 :             DO ipgf = 1, basis%subset(iset)%nexp
     309       11376 :                gto_basis%gcc(ipgf, ishell, iset) = basis%subset(iset)%coeff(ipgf, ishell)
     310             :             END DO
     311             :          END DO
     312        3645 :          DO ipgf = 1, basis%subset(iset)%nexp
     313        3240 :             gto_basis%zet(ipgf, iset) = basis%subset(iset)%exps(ipgf)
     314             :          END DO
     315             :       END DO
     316             : 
     317         405 :    END SUBROUTINE transfer_data_to_gto
     318             : 
     319             : END MODULE optbas_fenv_manipulation