-
Notifications
You must be signed in to change notification settings - Fork 1
/
qs_energy.F
200 lines (160 loc) · 8.63 KB
/
qs_energy.F
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
!--------------------------------------------------------------------------------------------------!
! CP2K: A general program to perform molecular dynamics simulations !
! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
! !
! SPDX-License-Identifier: GPL-2.0-or-later !
!--------------------------------------------------------------------------------------------------!
! **************************************************************************************************
!> \brief Perform a QUICKSTEP wavefunction optimization (single point)
!> \par History
!> none
!> \author MK (29.10.2002)
! **************************************************************************************************
MODULE qs_energy
USE almo_scf, ONLY: almo_entry_scf
USE cp_control_types, ONLY: dft_control_type
USE cp_external_control, ONLY: external_control
USE dm_ls_scf, ONLY: ls_scf
USE energy_corrections, ONLY: energy_correction
USE excited_states, ONLY: excited_state_energy
USE input_constants, ONLY: smeagol_runtype_emtransport
USE input_section_types, ONLY: section_vals_get,&
section_vals_get_subs_vals,&
section_vals_type,&
section_vals_val_get
USE lri_environment_methods, ONLY: lri_print_stat
USE mp2, ONLY: mp2_main
USE qs_active_space_methods, ONLY: active_space_main
USE qs_energy_init, ONLY: qs_energies_init
USE qs_energy_types, ONLY: qs_energy_type
USE qs_energy_utils, ONLY: qs_energies_properties
USE qs_environment_methods, ONLY: qs_env_rebuild_pw_env
USE qs_environment_types, ONLY: get_qs_env,&
qs_environment_type
USE qs_harris_utils, ONLY: harris_energy_correction
USE qs_ks_methods, ONLY: qs_ks_update_qs_env
USE qs_matrix_w, ONLY: compute_matrix_w
USE qs_nonscf, ONLY: nonscf
USE qs_scf, ONLY: scf
USE scf_control_types, ONLY: scf_control_type
#include "./base/base_uses.f90"
IMPLICIT NONE
PRIVATE
! *** Global parameters ***
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_energy'
PUBLIC :: qs_energies
CONTAINS
! **************************************************************************************************
!> \brief Driver routine for QUICKSTEP single point wavefunction optimization.
!> \param qs_env ...
!> \param consistent_energies ...
!> \param calc_forces ...
!> \date 29.10.2002
!> \par History
!> - consistent_energies option added (25.08.2005, TdK)
!> - introduced driver for energy in order to properly decide between
!> SCF or RTP (fschiff 02.09)
!> \author MK
!> \version 1.0
! **************************************************************************************************
SUBROUTINE qs_energies(qs_env, consistent_energies, calc_forces)
TYPE(qs_environment_type), POINTER :: qs_env
LOGICAL, INTENT(IN), OPTIONAL :: consistent_energies, calc_forces
CHARACTER(len=*), PARAMETER :: routineN = 'qs_energies'
INTEGER :: handle
LOGICAL :: do_consistent_energies, &
do_excited_state, loverlap_deltat, &
my_calc_forces, run_rtp
TYPE(dft_control_type), POINTER :: dft_control
TYPE(qs_energy_type), POINTER :: energy
TYPE(scf_control_type), POINTER :: scf_control
TYPE(section_vals_type), POINTER :: excited_state_section
CALL timeset(routineN, handle)
my_calc_forces = .FALSE.
IF (PRESENT(calc_forces)) my_calc_forces = calc_forces
do_consistent_energies = .FALSE.
IF (PRESENT(consistent_energies)) do_consistent_energies = consistent_energies
CALL qs_env_rebuild_pw_env(qs_env)
CALL get_qs_env(qs_env=qs_env, run_rtp=run_rtp)
IF (.NOT. run_rtp) THEN
NULLIFY (dft_control, energy)
CALL qs_energies_init(qs_env, my_calc_forces)
CALL get_qs_env(qs_env=qs_env, dft_control=dft_control, scf_control=scf_control, energy=energy)
! *** check if only overlap matrix is needed for couplings
loverlap_deltat = .FALSE.
NULLIFY (excited_state_section)
excited_state_section => section_vals_get_subs_vals(qs_env%input, "DFT%EXCITED_STATES")
CALL section_vals_get(excited_state_section, explicit=do_excited_state)
IF (do_excited_state) THEN
CALL section_vals_val_get(excited_state_section, "OVERLAP_DELTAT", &
l_val=loverlap_deltat)
END IF
! *** Perform a SCF run ***
IF (.NOT. loverlap_deltat) THEN
IF (scf_control%non_selfconsistent .AND. .NOT. scf_control%force_scf_calculation) THEN
CALL nonscf(qs_env)
ELSE IF (dft_control%qs_control%do_ls_scf) THEN
CALL ls_scf(qs_env)
ELSE IF (dft_control%qs_control%do_almo_scf) THEN
CALL almo_entry_scf(qs_env, calc_forces=my_calc_forces)
ELSE
! current-induced forces
IF (dft_control%smeagol_control%smeagol_enabled .AND. &
dft_control%smeagol_control%run_type == smeagol_runtype_emtransport) THEN
dft_control%smeagol_control%emforces = my_calc_forces
END IF
CALL scf(qs_env)
END IF
END IF
IF (do_consistent_energies) THEN
CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., just_energy=.FALSE.)
END IF
IF (.NOT. (dft_control%qs_control%do_ls_scf .OR. dft_control%qs_control%do_almo_scf)) THEN
! Compute MP2 energy
CALL qs_energies_mp2(qs_env, my_calc_forces)
IF (.NOT. ASSOCIATED(qs_env%mp2_env)) THEN
! do not overwrite w matrix computed by SMEAGOL (current-induced forces)
IF (.NOT. (dft_control%smeagol_control%smeagol_enabled .AND. &
dft_control%smeagol_control%run_type == smeagol_runtype_emtransport)) THEN
! if calculate forces, time to compute the w matrix
CALL compute_matrix_w(qs_env, my_calc_forces)
END IF
END IF
END IF
! Check for energy correction
IF (qs_env%harris_method) THEN
CALL harris_energy_correction(qs_env, my_calc_forces)
END IF
! Do active space calculation
CALL active_space_main(qs_env)
! Check for energy correction
IF (qs_env%energy_correction) THEN
CALL energy_correction(qs_env, ec_init=.TRUE., calculate_forces=.FALSE.)
END IF
IF (.NOT. loverlap_deltat) THEN
CALL qs_energies_properties(qs_env, calc_forces)
CALL excited_state_energy(qs_env, calculate_forces=.FALSE.)
END IF
IF (dft_control%qs_control%lrigpw) THEN
CALL lri_print_stat(qs_env)
END IF
END IF
CALL timestop(handle)
END SUBROUTINE qs_energies
! **************************************************************************************************
!> \brief Enters the mp2 part of cp2k
!> \param qs_env ...
!> \param calc_forces ...
! **************************************************************************************************
SUBROUTINE qs_energies_mp2(qs_env, calc_forces)
TYPE(qs_environment_type), POINTER :: qs_env
LOGICAL, INTENT(IN) :: calc_forces
LOGICAL :: should_stop
! Compute MP2 energy
IF (ASSOCIATED(qs_env%mp2_env)) THEN
CALL external_control(should_stop, "MP2", target_time=qs_env%target_time, &
start_time=qs_env%start_time)
CALL mp2_main(qs_env=qs_env, calc_forces=calc_forces)
END IF
END SUBROUTINE qs_energies_mp2
END MODULE qs_energy