1!! Copyright (C) 2004 E.S. Kadantsev, M. Marques 2!! 3!! This program is free software; you can redistribute it and/or modify 4!! it under the terms of the GNU General Public License as published by 5!! the Free Software Foundation; either version 2, or (at your option) 6!! any later version. 7!! 8!! This program is distributed in the hope that it will be useful, 9!! but WITHOUT ANY WARRANTY; without even the implied warranty of 10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11!! GNU General Public License for more details. 12!! 13!! You should have received a copy of the GNU General Public License 14!! along with this program; if not, write to the Free Software 15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 16!! 02110-1301, USA. 17!! 18 19#include "global.h" 20 21module scf_tol_oct_m 22 use global_oct_m 23 use messages_oct_m 24 use namespace_oct_m 25 use parser_oct_m 26 use varinfo_oct_m 27 28 implicit none 29 30 private 31 32 integer, public, parameter :: & 33 SCF_TOL_FIXED = 0, & 34 SCF_TOL_ADAPTIVE = 1, & 35 SCF_TOL_LINEAR = 2, & 36 SCF_TOL_EXP = 3 37 38 public :: & 39 scf_tol_t, & 40 scf_tol_init, & 41 scf_tol_end, & 42 scf_tol_stop, & 43 scf_tol_step, & 44 scf_tol_final, & 45 scf_tol_obsolete_variables 46 47 type scf_tol_t 48 private 49 integer, public :: max_iter 50 integer :: scheme 51 FLOAT, public :: conv_rel_dens 52 FLOAT, public :: conv_abs_dens 53 FLOAT :: conv_threshold_use 54 FLOAT :: dynamic_tol_factor 55 FLOAT :: current_tol 56 FLOAT :: initial_tol 57 FLOAT, public :: final_tol 58 integer :: iter_window 59 end type scf_tol_t 60 61contains 62 63 !----------------------------------------------------------------- 64 subroutine scf_tol_init(this, namespace, qtot, def_maximumiter, tol_scheme) 65 type(scf_tol_t), intent(out) :: this 66 type(namespace_t), intent(in) :: namespace 67 FLOAT, intent(in) :: qtot 68 integer, optional, intent(in) :: def_maximumiter 69 integer, optional, intent(in) :: tol_scheme 70 71 integer :: def_maximumiter_ 72 73 PUSH_SUB(scf_tol_init) 74 75 !%Variable LRMaximumIter 76 !%Type integer 77 !%Default 200 78 !%Section Linear Response::SCF in LR calculations 79 !%Description 80 !% The maximum number of SCF iterations to calculate response. 81 !%End 82 def_maximumiter_ = 200 83 if(present(def_maximumiter)) def_maximumiter_ = def_maximumiter 84 85 call parse_variable(namespace, 'LRMaximumIter', def_maximumiter_, this%max_iter) 86 87 !%Variable LRConvAbsDens 88 !%Type float 89 !%Default 1e-5 90 !%Section Linear Response::SCF in LR calculations 91 !%Description 92 !% The tolerance in the absolute variation of the density response, to determine if 93 !% the SCF for linear response is converged. 94 !% <math>\varepsilon = \int {\rm d}^3r \left| \rho^{out}(\bf r) -\rho^{inp}(\bf r) \right|</math>. 95 !% A zero value means do not use this criterion. 96 !%End 97 98 call parse_variable(namespace, 'LRConvAbsDens', CNST(1e-5), this%conv_abs_dens) 99 100 !%Variable LRConvRelDens 101 !%Type float 102 !%Default 0.0 103 !%Section Linear Response::SCF in LR calculations 104 !%Description 105 !% The tolerance in the relative variation of the density response, to determine if 106 !% the SCF for linear response is converged. 107 !% <math>\varepsilon = \frac{1}{N_{\rm elec}}</math> <tt>LRConvAbsDens</tt>. 108 !% A zero value means do not use this criterion. 109 !%End 110 111 call parse_variable(namespace, 'LRConvRelDens', M_ZERO, this%conv_rel_dens) 112 113 ! value to use for adaptive tol scheme 114 if(this%conv_abs_dens <= M_ZERO) then 115 this%conv_threshold_use = this%conv_rel_dens * qtot 116 else if(this%conv_rel_dens <= M_ZERO) then 117 this%conv_threshold_use = this%conv_abs_dens 118 else 119 this%conv_threshold_use = min(this%conv_abs_dens, this%conv_rel_dens * qtot) 120 end if 121 122 if(this%conv_abs_dens <= M_ZERO .and. this%conv_rel_dens <= M_ZERO) then 123 message(1) = "Input: Not all convergence criteria can be <= 0" 124 message(2) = "Please set one of the following to a nonzero value:" 125 message(3) = "LRConvAbsDens | LRConvRelDens" 126 call messages_fatal(3, namespace=namespace) 127 end if 128 129 !%Variable LRTolScheme 130 !%Type integer 131 !%Default tol_adaptive 132 !%Section Linear Response::SCF in LR calculations 133 !%Description 134 !% The scheme used to adjust the tolerance of the solver during 135 !% the SCF iteration. For <tt>kdotp</tt> and magnetic <tt>em_resp</tt> modes, or 136 !% whenever <tt>HamiltonianVariation = V_ext_only</tt>, the 137 !% scheme is set to <tt>tol_fixed</tt>, and this variable is ignored. 138 !%Option tol_fixed 0 139 !% The solver tolerance is fixed for all the iterations; this 140 !% improves convergence but increases the computational cost 141 !%Option tol_adaptive 1 142 !% The tolerance is increased according to the level of 143 !% convergence of the SCF. 144 !%Option tol_linear 2 145 !% The tolerance decreases linearly for the first <tt>LRTolIterWindow</tt> iterations. 146 !%Option tol_exp 3 147 !% The tolerance decreases exponentially for the first <tt>LRTolIterWindow</tt> iterations. 148 !%End 149 if(present(tol_scheme)) then 150 this%scheme = tol_scheme 151 else 152 call parse_variable(namespace, 'LRTolScheme', SCF_TOL_ADAPTIVE, this%scheme) 153 end if 154 if(.not.varinfo_valid_option('LRTolScheme', this%scheme)) then 155 call messages_input_error(namespace, 'LRTolScheme') 156 end if 157 158 !%Variable LRTolInitTol 159 !%Type float 160 !%Default 1e-2 161 !%Section Linear Response::Solver 162 !%Description 163 !% This is the tolerance to determine that the linear solver has converged, 164 !% for the first SCF iteration. Ignored if <tt>LRTolScheme = fixed</tt>. 165 !%End 166 167 call parse_variable(namespace, 'LRTolInitTol', CNST(1e-2), this%initial_tol) 168 this%current_tol = this%initial_tol 169 170 !%Variable LRTolFinalTol 171 !%Type float 172 !%Default 1e-6 173 !%Section Linear Response::Solver 174 !%Description 175 !% This is the tolerance to determine that the linear solver has converged. 176 !%End 177 178 call parse_variable(namespace, 'LRTolFinalTol', CNST(1e-6), this%final_tol) 179 180 if(this%scheme == SCF_TOL_ADAPTIVE) then 181 !%Variable LRTolAdaptiveFactor 182 !%Type float 183 !%Default 0.1 184 !%Section Linear Response::SCF in LR calculations 185 !%Description 186 !% This factor controls how much the tolerance is decreased 187 !% during the self-consistency process. Smaller values mean that 188 !% tolerance is decreased faster. 189 !%End 190 191 call parse_variable(namespace, 'LRTolAdaptiveFactor', CNST(0.1), this%dynamic_tol_factor) 192 end if 193 194 if(this%scheme==SCF_TOL_LINEAR.or.this%scheme==SCF_TOL_EXP) then 195 !%Variable LRTolIterWindow 196 !%Type float 197 !%Default 10 198 !%Section Linear Response::SCF in LR calculations 199 !%Description 200 !% Number of iterations necessary to reach the final tolerance. 201 !%End 202 203 call parse_variable(namespace, 'LRTolIterWindow', 10, this%iter_window) 204 end if 205 206 POP_SUB(scf_tol_init) 207 208 end subroutine scf_tol_init 209 210 211 !----------------------------------------------------------------- 212 FLOAT function scf_tol_step(this, iter, scf_res) result(r) 213 type(scf_tol_t), intent(inout) :: this 214 integer, intent(in) :: iter 215 FLOAT, intent(in) :: scf_res 216 217 FLOAT :: logi, logf 218 219 PUSH_SUB(scf_tol_step) 220 221 if(iter == 0) this%current_tol = M_HUGE 222 223 select case(this%scheme) 224 case(SCF_TOL_FIXED) 225 r = this%final_tol 226 227 case(SCF_TOL_ADAPTIVE) 228 if(iter == 0) then 229 r = this%initial_tol 230 else 231 r = this%dynamic_tol_factor * (this%final_tol/this%conv_threshold_use)*scf_res 232 end if 233 234 case(SCF_TOL_LINEAR) 235 r = this%initial_tol + (this%final_tol - this%initial_tol) * & 236 TOFLOAT(iter) / TOFLOAT(this%iter_window) 237 238 case(SCF_TOL_EXP) 239 logi = log(this%initial_tol) 240 logf = log(this%final_tol) 241 r = logi + (logf - logi) * & 242 TOFLOAT(iter) / TOFLOAT(this%iter_window) 243 r = exp(r) 244 end select 245 246 ! tolerance can never be larger than final tolerance 247 r = max(r, this%final_tol) 248 ! tolerance always has to decrease 249 r = min(r, this%current_tol) 250 251 this%current_tol = r 252 253 POP_SUB(scf_tol_step) 254 end function scf_tol_step 255 256 257 !----------------------------------------------------------------- 258 subroutine scf_tol_stop(this) 259 type(scf_tol_t), intent(inout) :: this 260 261 PUSH_SUB(scf_tol_stop) 262 this%current_tol = M_ZERO 263 264 POP_SUB(scf_tol_stop) 265 end subroutine scf_tol_stop 266 267 268 !----------------------------------------------------------------- 269 subroutine scf_tol_end(this) 270 type(scf_tol_t), intent(inout) :: this 271 272 PUSH_SUB(scf_tol_end) 273 this%current_tol = M_ZERO 274 275 POP_SUB(scf_tol_end) 276 end subroutine scf_tol_end 277 278 279 !----------------------------------------------------------------- 280 FLOAT function scf_tol_final(this) 281 type(scf_tol_t), intent(in) :: this 282 283 scf_tol_final = this%final_tol 284 285 end function scf_tol_final 286 287 !----------------------------------------------------------------- 288 289 subroutine scf_tol_obsolete_variables(namespace, old_prefix, new_prefix) 290 type(namespace_t), intent(in) :: namespace 291 character(len=*), intent(in) :: old_prefix 292 character(len=*), intent(in) :: new_prefix 293 294 PUSH_SUB(scf_tol_obsolete_variables) 295 296 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRMaximumIter', trim(new_prefix)//'LRMaximumIter') 297 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRConvAbsDens', trim(new_prefix)//'LRConvAbsDens') 298 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRTolScheme', trim(new_prefix)//'LRTolScheme') 299 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRTolInitTol', trim(new_prefix)//'LRTolInitTol') 300 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRTolFinalTol', trim(new_prefix)//'LRTolFinalTol') 301 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRTolAdaptiveFactor', trim(new_prefix)//'LRTolAdaptiveFactor') 302 call messages_obsolete_variable(namespace, trim(old_prefix)//'LRTolIterWindow', trim(new_prefix)//'LRTolIterWindow') 303 304 POP_SUB(scf_tol_obsolete_variables) 305 end subroutine scf_tol_obsolete_variables 306 307end module scf_tol_oct_m 308 309!! Local Variables: 310!! mode: f90 311!! coding: utf-8 312!! End: 313