Octopus
xc_oep_photon.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!! Copyright (C) 2012-2013 M. Gruning, P. Melo, M. Oliveira
3!!
4!! This program is free software; you can redistribute it and/or modify
5!! it under the terms of the GNU General Public License as published by
6!! the Free Software Foundation; either version 2, or (at your option)
7!! any later version.
8!!
9!! This program is distributed in the hope that it will be useful,
10!! but WITHOUT ANY WARRANTY; without even the implied warranty of
11!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12!! GNU General Public License for more details.
13!!
14!! You should have received a copy of the GNU General Public License
15!! along with this program; if not, write to the Free Software
16!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
17!! 02110-1301, USA.
18!!
19
20#include "global.h"
21
23 use comm_oct_m
24 use debug_oct_m
28 use global_oct_m
29 use grid_oct_m
32 use, intrinsic :: iso_fortran_env
38 use math_oct_m
40 use mesh_oct_m
42 use mpi_oct_m
46 use parser_oct_m
50 use space_oct_m
57 use xc_oct_m
58 use xc_oep_oct_m
59 use xc_f03_lib_m
61
62 implicit none
63
64 private
65 public :: &
71
72 type, extends(xc_oep_t) :: xc_oep_photon_t
73 private
74 logical :: coc_translation
75 type(photon_mode_t), public :: pt
76 type(lr_t) :: photon_lr
77 logical, public :: rm_ee_interaction
78 end type xc_oep_photon_t
79
80contains
81
82 ! ---------------------------------------------------------
83 subroutine xc_oep_photon_init(oep, namespace, family, gr, st, mc, space)
84 type(xc_oep_photon_t), intent(inout) :: oep
85 type(namespace_t), intent(in) :: namespace
86 integer, intent(in) :: family
87 type(grid_t), intent(inout) :: gr
88 type(states_elec_t), intent(in) :: st
89 type(multicomm_t), intent(in) :: mc
90 class(space_t), intent(in) :: space
91
92 push_sub(xc_oep_photon_init)
93
94 if(bitand(family, xc_family_oep) == 0) then
95 oep%level = oep_level_none
96 pop_sub(xc_oep_photon_init)
97 return
98 end if
99
100 oep%type = oep_type_photons
101
102 call messages_obsolete_variable(namespace, 'OEP_Level', 'OEPLevel')
103 call parse_variable(namespace, 'OEPLevel', oep_level_kli, oep%level)
104 if(.not. varinfo_valid_option('OEPLevel', oep%level)) call messages_input_error(namespace, 'OEPLevel')
105
106 if(oep%level /= oep_level_none) then
107 call messages_experimental("EnablePhotons = yes")
108 call photon_mode_init(oep%pt, namespace, space%dim)
109 call photon_mode_set_n_electrons(oep%pt, st%qtot)
110 if (oep%pt%nmodes > 1) then
111 call messages_not_implemented('Photon OEP for more than one photon mode.')
112 end if
113 if (oep%level == oep_level_full .and. st%d%nspin /= unpolarized) then
114 call messages_not_implemented('Spin-polarized calculations with photon OEP.')
115 end if
116
117 if (states_are_complex(st)) then
118 call messages_not_implemented('Photon OEP with complex wavefunctions.')
119 end if
120
121 if (st%nik > st%d%ispin .and. oep%level == oep_level_full) then
122 call messages_not_implemented("Full OEP for periodic systems", namespace=namespace)
123 end if
124 if (st%nik > st%d%ispin .and. st%d%ispin==spinors) then
125 call messages_not_implemented("OEP for periodic systems with spinors", namespace=namespace)
126 end if
127
128 safe_allocate(oep%pt%correlator(1:gr%np, 1:oep%pt%nmodes))
129 oep%pt%correlator = m_zero
130
131 call photon_mode_compute_dipoles(oep%pt, gr)
132 end if
133
134 if(oep%level == oep_level_full) then
135
136 call messages_experimental("Full OEP")
137 !%Variable OEPMixing
138 !%Type float
139 !%Default 1.0
140 !%Section Hamiltonian::XC
141 !%Description
142 !% The linear mixing factor used to solve the Sternheimer
143 !% equation in the full OEP procedure.
144 !%End
145 call messages_obsolete_variable(namespace, 'OEP_Mixing', 'OEPMixing')
146 call parse_variable(namespace, 'OEPMixing', m_one, oep%mixing)
147
148 !%Variable OEPMixingScheme
149 !%Type integer
150 !%Default 1.0
151 !%Section Hamiltonian::XC
152 !%Description
153 !%Different Mixing Schemes are possible
154 !%Option OEP_MIXING_SCHEME_CONST 1
155 !%Use a constant
156 !%Reference: S. Kuemmel and J. Perdew, <i>Phys. Rev. Lett.</i> <b>90</b>, 4, 043004 (2003)
157 !%Option OEP_MIXING_SCHEME_BB 2
158 !%Use the Barzilai-Borwein (BB) Method
159 !%Reference: T. W. Hollins, S. J. Clark, K. Refson, and N. I. Gidopoulos,
160 !%<i>Phys. Rev. B</i> <b>85<\b>, 235126 (2012)
161 !%Option OEP_MIXING_SCHEME_DENS 3
162 !%Use the inverse of the electron density
163 !%Reference: S. Kuemmel and J. Perdew, <i>Phys. Rev. B</i> <b>68</b>, 035103 (2003)
164 !%End
165 call parse_variable(namespace, 'OEPMixingScheme', oep_mixing_scheme_const, oep%mixing_scheme)
166
167 if (oep%mixing_scheme == oep_mixing_scheme_bb) then
168 safe_allocate(oep%vxc_old(1:gr%np,st%d%ispin))
169 safe_allocate(oep%ss_old(1:gr%np,st%d%ispin))
170 oep%vxc_old = m_zero
171 oep%ss_old = m_zero
172 end if
173
174 oep%norm2ss = m_zero
175 end if
177 ! obtain the spin factors
178 call xc_oep_spinfactor(oep, st%d%nspin)
179
180 ! This variable will keep vxc across iterations
181 if ((st%d%ispin==3) .or. oep%level == oep_level_full) then
182 safe_allocate(oep%vxc(1:gr%np,st%d%nspin))
183 else
184 safe_allocate(oep%vxc(1:gr%np,1:min(st%d%nspin, 2)))
185 end if
186 oep%vxc = m_zero
187
188 !%Variable KLIPhotonCOC
189 !%Type logical
190 !%Default .false.
191 !%Section Hamiltonian::XC
192 !%Description
193 !% Activate the center of charge translation of the electric dipole operator which should
194 !% avoid the dependence of the photon KLI on an permanent dipole.
195 !%End
196
197 ! when performing full OEP, we need to solve a linear equation
198 call scf_tol_init(oep%scftol, namespace, st%qtot, def_maximumiter=10)
199 call linear_solver_init(oep%solver, namespace, gr, states_are_real(st), mc, space)
200 call lr_init(oep%lr)
201 call lr_init(oep%photon_lr)
202 call parse_variable(namespace, 'KLIPhotonCOC', .false., oep%coc_translation)
203
204 !%Variable OEPRemoveElectron
205 !%Type logical
206 !%Default .false.
207 !%Section Hamiltonian::XC
208 !%Description
209 !% Remove electron-electron interaction in OEP-Photon calculations
210 !%End
211
212 call parse_variable(namespace, 'OEPRemoveElectron', .false., oep%rm_ee_interaction)
213
214 ! the linear equation has to be more converged if we are to attain the required precision
215 !oep%lr%conv_abs_dens = oep%lr%conv_abs_dens / (oep%mixing)
216
217 if(st%d%kpt%parallel) then
218 call messages_not_implemented("OEP parallel in spin/k-points", namespace=namespace)
219 end if
220
221 pop_sub(xc_oep_photon_init)
222 end subroutine xc_oep_photon_init
223
224
225 ! ---------------------------------------------------------
226 subroutine xc_oep_photon_end(oep)
227 type(xc_oep_photon_t), intent(inout) :: oep
228
229 push_sub(xc_oep_photon_end)
230
231 call xc_oep_end(oep)
232
233 if(oep%level /= oep_level_none) then
234 call lr_dealloc(oep%photon_lr)
235 call photon_mode_end(oep%pt)
236 end if
237
238 pop_sub(xc_oep_photon_end)
239 end subroutine xc_oep_photon_end
240
241#include "xc_oep_qed_inc.F90"
242
243#include "undef.F90"
244#include "real.F90"
245#include "xc_kli_photon_inc.F90"
246#include "xc_oep_photon_inc.F90"
247
248#include "undef.F90"
249#include "complex.F90"
250#include "xc_kli_photon_inc.F90"
251#include "xc_oep_photon_inc.F90"
252
253end module xc_oep_photon_oct_m
254
255!! Local Variables:
256!! mode: f90
257!! coding: utf-8
258!! End:
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
integer, parameter, public unpolarized
Parameters...
integer, parameter, public spinors
real(real64), parameter, public m_zero
Definition: global.F90:188
real(real64), parameter, public m_one
Definition: global.F90:189
This module implements the underlying real-space grid.
Definition: grid.F90:117
subroutine, public lr_init(lr)
subroutine, public lr_dealloc(lr)
subroutine, public linear_solver_init(this, namespace, gr, states_are_real, mc, space)
This module is intended to contain "only mathematical" functions and procedures.
Definition: math.F90:115
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
subroutine, public messages_not_implemented(feature, namespace)
Definition: messages.F90:1113
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1045
subroutine, public messages_input_error(namespace, var, details, row, column)
Definition: messages.F90:713
subroutine, public messages_experimental(name, namespace)
Definition: messages.F90:1085
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:145
subroutine, public photon_mode_compute_dipoles(this, mesh)
Computes the polarization dipole.
subroutine, public photon_mode_end(this)
subroutine, public photon_mode_set_n_electrons(this, qtot)
subroutine, public photon_mode_init(this, namespace, dim, photon_free)
subroutine, public scf_tol_init(this, namespace, qtot, def_maximumiter, tol_scheme)
Definition: scf_tol.F90:160
pure logical function, public states_are_complex(st)
pure logical function, public states_are_real(st)
This module handles spin dimensions of the states and the k-point distribution.
Definition: xc.F90:114
integer, parameter, public oep_level_none
the OEP levels
Definition: xc_oep.F90:172
subroutine, public xc_oep_spinfactor(oep, nspin)
A couple of auxiliary functions for oep.
Definition: xc_oep.F90:396
subroutine, public xc_oep_end(oep)
Definition: xc_oep.F90:356
integer, parameter, public oep_type_photons
Definition: xc_oep.F90:184
integer, parameter, public oep_mixing_scheme_bb
Definition: xc_oep.F90:178
integer, parameter, public oep_level_full
Definition: xc_oep.F90:172
integer, parameter, public oep_mixing_scheme_const
Mixing schemes.
Definition: xc_oep.F90:178
integer, parameter, public oep_level_kli
Definition: xc_oep.F90:172
subroutine, public zxc_oep_photon_calc(oep, namespace, xcs, gr, hm, st, space, ex, ec, vxc)
This file handles the evaluation of the OEP potential, in the KLI or full OEP as described in S....
subroutine, public xc_oep_photon_init(oep, namespace, family, gr, st, mc, space)
subroutine, public xc_oep_photon_end(oep)
subroutine, public dxc_oep_photon_calc(oep, namespace, xcs, gr, hm, st, space, ex, ec, vxc)
This file handles the evaluation of the OEP potential, in the KLI or full OEP as described in S....