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