Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Local FCI solver interface for active-space embedding.
10 : ! **************************************************************************************************
11 : MODULE qs_active_space_fci
12 : #if defined(__LIBFCI)
13 : USE ISO_C_BINDING, ONLY: C_CHAR, &
14 : C_DOUBLE, &
15 : C_INT, &
16 : C_NULL_CHAR
17 : #endif
18 : USE kinds, ONLY: dp
19 : USE message_passing, ONLY: mp_para_env_type
20 : USE qs_active_space_types, ONLY: active_space_type
21 : #if defined(__LIBFCI)
22 : USE qs_active_space_utils, ONLY: eri_to_array, &
23 : subspace_matrix_to_array
24 : #endif
25 : #include "./base/base_uses.f90"
26 :
27 : IMPLICIT NONE
28 : PRIVATE
29 :
30 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_active_space_fci'
31 :
32 : PUBLIC :: solve_active_space_fci
33 :
34 : #if defined(__LIBFCI)
35 : INTERFACE
36 : FUNCTION libfci_solve(nspins, norb, nelec, ms2, h1_a, h1_b, eri_aa, eri_ab, eri_bb, &
37 : max_iter, threshold, pspace, project_spin, verbose, energy, &
38 : rdm1_a, rdm1_b, err_msg, err_msg_len) RESULT(status) &
39 : BIND(C, name="libfci_solve")
40 : IMPORT :: C_CHAR, C_DOUBLE, C_INT
41 : INTEGER(C_INT), VALUE :: nspins, norb, nelec, ms2
42 : REAL(C_DOUBLE), DIMENSION(*), INTENT(IN) :: h1_a, h1_b, eri_aa, eri_ab, eri_bb
43 : INTEGER(C_INT), VALUE :: max_iter, pspace, project_spin, verbose
44 : REAL(C_DOUBLE), VALUE :: threshold
45 : REAL(C_DOUBLE), INTENT(OUT) :: energy
46 : REAL(C_DOUBLE), DIMENSION(*), INTENT(OUT) :: rdm1_a, rdm1_b
47 : CHARACTER(KIND=C_CHAR), DIMENSION(*), INTENT(OUT) :: err_msg
48 : INTEGER(C_INT), VALUE :: err_msg_len
49 : INTEGER(C_INT) :: status
50 : END FUNCTION libfci_solve
51 : END INTERFACE
52 : #endif
53 :
54 : CONTAINS
55 :
56 : ! **************************************************************************************************
57 : !> \brief Solve the current active-space Hamiltonian with the local FCI kernel.
58 : !> \param active_space_env active-space environment
59 : !> \param para_env parallel environment
60 : !> \param p_act_mo_a alpha or spin-summed active-space 1-RDM
61 : !> \param p_act_mo_b beta active-space 1-RDM for unrestricted calculations
62 : ! **************************************************************************************************
63 22 : SUBROUTINE solve_active_space_fci(active_space_env, para_env, p_act_mo_a, p_act_mo_b)
64 : TYPE(active_space_type), INTENT(INOUT), POINTER :: active_space_env
65 : TYPE(mp_para_env_type), INTENT(IN), POINTER :: para_env
66 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
67 : INTENT(OUT) :: p_act_mo_a
68 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:), &
69 : INTENT(OUT), OPTIONAL :: p_act_mo_b
70 :
71 : #if defined(__LIBFCI)
72 : CHARACTER(KIND=C_CHAR), DIMENSION(512) :: c_error
73 : CHARACTER(LEN=512) :: error_text
74 : INTEGER :: max_iter, ms2, n2, n4, nmo_active, &
75 : nspins, status
76 : LOGICAL :: ionode, restricted_orbitals
77 : REAL(KIND=dp) :: energy_active, threshold
78 22 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eri_aa, eri_ab, eri_bb, fock_a, fock_b, &
79 22 : p_beta
80 :
81 22 : nmo_active = active_space_env%nmo_active
82 22 : nspins = active_space_env%nspins
83 22 : restricted_orbitals = active_space_env%restricted_orbitals
84 22 : n2 = nmo_active*nmo_active
85 22 : n4 = n2*n2
86 22 : ionode = para_env%is_source()
87 :
88 132 : ALLOCATE (fock_a(n2), eri_aa(n4), p_act_mo_a(n2))
89 154 : ALLOCATE (fock_b(MAX(1, n2)), eri_ab(MAX(1, n4)), eri_bb(MAX(1, n4)), p_beta(MAX(1, n2)))
90 22 : fock_b(:) = 0.0_dp
91 22 : eri_ab(:) = 0.0_dp
92 22 : eri_bb(:) = 0.0_dp
93 22 : p_beta(:) = 0.0_dp
94 :
95 : ASSOCIATE (act_indices => active_space_env%active_orbitals(:, 1))
96 22 : CALL subspace_matrix_to_array(active_space_env%fock_sub(1), fock_a, act_indices, act_indices)
97 : END ASSOCIATE
98 22 : CALL eri_to_array(active_space_env%eri, eri_aa, active_space_env%active_orbitals, 1, 1)
99 :
100 22 : IF (nspins == 2) THEN
101 14 : IF (.NOT. PRESENT(p_act_mo_b)) THEN
102 0 : CPABORT("Missing beta output buffer for local active-space FCI solver.")
103 : END IF
104 28 : ALLOCATE (p_act_mo_b(n2))
105 : ASSOCIATE (act_indices => active_space_env%active_orbitals(:, 2))
106 14 : CALL subspace_matrix_to_array(active_space_env%fock_sub(2), fock_b(1:n2), act_indices, act_indices)
107 : END ASSOCIATE
108 14 : IF (restricted_orbitals) THEN
109 2056 : eri_ab(1:n4) = eri_aa
110 2056 : eri_bb(1:n4) = eri_aa
111 : ELSE
112 6 : CALL eri_to_array(active_space_env%eri, eri_ab(1:n4), active_space_env%active_orbitals, 1, 2)
113 6 : CALL eri_to_array(active_space_env%eri, eri_bb(1:n4), active_space_env%active_orbitals, 2, 2)
114 : END IF
115 : END IF
116 :
117 22 : status = 0
118 22 : energy_active = 0.0_dp
119 22 : error_text = ""
120 11286 : c_error(:) = C_NULL_CHAR
121 22 : ms2 = active_space_env%multiplicity - 1
122 22 : max_iter = 2048
123 : ! LibFCI tests the squared residual norm, not the energy error.
124 22 : threshold = active_space_env%eps_fci
125 :
126 22 : IF (ionode) THEN
127 : status = libfci_solve(INT(nspins, C_INT), INT(nmo_active, C_INT), &
128 : INT(active_space_env%nelec_active, C_INT), INT(ms2, C_INT), &
129 : fock_a, fock_b, eri_aa, eri_ab, eri_bb, INT(max_iter, C_INT), &
130 : REAL(threshold, C_DOUBLE), INT(200, C_INT), INT(1, C_INT), INT(0, C_INT), &
131 11 : energy_active, p_act_mo_a, p_beta, c_error, INT(SIZE(c_error), C_INT))
132 11 : error_text = c_string(c_error)
133 : END IF
134 :
135 22 : CALL para_env%bcast(status, para_env%source)
136 22 : CALL para_env%bcast(error_text, para_env%source)
137 22 : IF (status /= 0) THEN
138 0 : IF (LEN_TRIM(error_text) > 0) THEN
139 0 : CPABORT("Local active-space FCI solver failed: "//TRIM(error_text))
140 : ELSE
141 0 : CPABORT("Local active-space FCI solver failed.")
142 : END IF
143 : END IF
144 :
145 22 : CALL para_env%bcast(energy_active, para_env%source)
146 22 : CALL para_env%bcast(p_act_mo_a, para_env%source)
147 22 : active_space_env%energy_active = energy_active
148 22 : IF (nspins == 2) THEN
149 238 : p_act_mo_b(:) = p_beta(1:n2)
150 14 : CALL para_env%bcast(p_act_mo_b, para_env%source)
151 : END IF
152 :
153 22 : DEALLOCATE (fock_a, fock_b, eri_aa, eri_ab, eri_bb, p_beta)
154 : #else
155 : MARK_USED(active_space_env)
156 : MARK_USED(para_env)
157 : MARK_USED(p_act_mo_a)
158 : IF (PRESENT(p_act_mo_b)) THEN
159 : MARK_USED(p_act_mo_b)
160 : END IF
161 : CPABORT("AS_SOLVER FCI requires CP2K to be built with LibFCI support.")
162 : #endif
163 :
164 22 : END SUBROUTINE solve_active_space_fci
165 :
166 : #if defined(__LIBFCI)
167 : ! **************************************************************************************************
168 : !> \brief Convert a C null-terminated string buffer to a Fortran string.
169 : !> \param buffer ...
170 : !> \return ...
171 : ! **************************************************************************************************
172 11 : FUNCTION c_string(buffer) RESULT(text)
173 : CHARACTER(KIND=C_CHAR), DIMENSION(:), INTENT(IN) :: buffer
174 : CHARACTER(LEN=512) :: text
175 :
176 : INTEGER :: i
177 :
178 11 : text = ""
179 11 : DO i = 1, MIN(SIZE(buffer), LEN(text))
180 11 : IF (buffer(i) == C_NULL_CHAR) EXIT
181 11 : text(i:i) = buffer(i)
182 : END DO
183 11 : END FUNCTION c_string
184 : #endif
185 :
186 : END MODULE qs_active_space_fci
|
