Octopus
eigensolver.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
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
22 use batch_oct_m
24 use debug_oct_m
30 use global_oct_m
31 use grid_oct_m
34 use, intrinsic :: iso_fortran_env
38 use loct_oct_m
39 use mesh_oct_m
43 use mpi_oct_m
47 use parser_oct_m
50 use smear_oct_m
51 use space_oct_m
58 use unit_oct_m
61 use xc_oct_m
62
63 implicit none
64
65 private
66 public :: &
70
71 type eigensolver_t
72 private
73 integer, public :: es_type
74
75 real(real64), public :: tolerance
76 integer, public :: es_maxiter
77
78 real(real64) :: imag_time
79
81 real(real64), allocatable, public :: diff(:, :)
82 integer, public :: matvec
83 integer, allocatable, public :: converged(:)
84
86 type(preconditioner_t), public :: pre
87
89 type(subspace_t) :: sdiag
90
91 integer :: rmmdiis_minimization_iter
92
93 logical, public :: folded_spectrum
94
95 ! cg options
96 logical, public :: orthogonalize_to_all
97 integer, public :: conjugate_direction
98 logical, public :: additional_terms
99 real(real64), public :: energy_change_threshold
100
101 ! Chebyshev filtering options
102 type(eigen_chebyshev_t), public :: cheby_params
103
104 type(exponential_t) :: exponential_operator
105 contains
106 procedure :: run => eigensolver_run
107 end type eigensolver_t
108
109
110 integer, public, parameter :: &
111 RS_PLAN = 11, &
112 rs_cg = 5, &
113 rs_evo = 9, &
115 rs_chebyshev = 12
116
117contains
118
119 ! ---------------------------------------------------------
120 subroutine eigensolver_init(eigens, namespace, gr, st, mc, space)
121 type(eigensolver_t), intent(out) :: eigens
122 type(namespace_t), intent(in) :: namespace
123 type(grid_t), intent(in) :: gr
124 type(states_elec_t), intent(in) :: st
125 type(multicomm_t), intent(in) :: mc
126 class(space_t), intent(in) :: space
127
128 integer :: default_iter, default_es
129 real(real64) :: default_tol
130 real(real64) :: mem
131
132 push_sub(eigensolver_init)
133
134 !%Variable Eigensolver
135 !%Type integer
136 !%Section SCF::Eigensolver
137 !%Description
138 !% Which eigensolver to use to obtain the lowest eigenvalues and
139 !% eigenfunctions of the Kohn-Sham Hamiltonian. The default is
140 !% conjugate gradients (<tt>cg</tt>), except that when parallelization in states is
141 !% enabled, the default is <tt>chebyshev_filter</tt>.
142 !%Option cg 5
143 !% Conjugate-gradients algorithm.
144 !%Option plan 11
145 !% Preconditioned Lanczos scheme. Ref: Y. Saad, A. Stathopoulos, J. Chelikowsky, K. Wu and S. Ogut,
146 !% "Solution of Large Eigenvalue Problems in Electronic Structure Calculations", <i>BIT</i> <b>36</b>, 1 (1996).
147 !%Option evolution 9
148 !% (Experimental) Propagation in imaginary time.
149 !%Option rmmdiis 10
150 !% Residual minimization scheme, direct inversion in the
151 !% iterative subspace eigensolver, based on the implementation of
152 !% Kresse and Furthm&uuml;ller [<i>Phys. Rev. B</i> <b>54</b>, 11169
153 !% (1996)]. This eigensolver requires almost no orthogonalization
154 !% so it can be considerably faster than the other options for
155 !% large systems. To improve its performance a large number of <tt>ExtraStates</tt>
156 !% are required (around 10-20% of the number of occupied states).
157 !% Note: with <tt>unocc</tt>, you will need to stop the calculation
158 !% by hand, since the highest states will probably never converge.
159 !% Usage with more than one block of states per node is experimental, unfortunately.
160 !%Option chebyshev_filter 12
161 !% A Chebyshev-filtered subspace iteration method, which avoids most of the explicit computation of
162 !% eigenvectors and results in a significant speedup over iterative diagonalization methods.
163 !% This method may be viewed as an approach to solve the original nonlinear Kohn-Sham equation by a
164 !% nonlinear subspace iteration technique, without emphasizing the intermediate linearized Kohn-Sham
165 !% eigenvalue problems. For further details, see [Zhou et. al.](http://dx.doi.org/10.1016/j.jcp.2014.06.056)
166 !%End
167
168 if(st%parallel_in_states) then
169 default_es = rs_chebyshev
170 else
171 default_es = rs_cg
172 end if
173
174 call parse_variable(namespace, 'Eigensolver', default_es, eigens%es_type)
176 if(st%parallel_in_states .and. .not. eigensolver_parallel_in_states(eigens)) then
177 message(1) = "The selected eigensolver is not parallel in states."
178 message(2) = "Please use the rmmdiis or Chebyshev filtering eigensolvers."
179 call messages_fatal(2, namespace=namespace)
180 end if
181
182 call messages_obsolete_variable(namespace, 'EigensolverVerbose')
183 call messages_obsolete_variable(namespace, 'EigensolverSubspaceDiag', 'SubspaceDiagonalization')
185 default_iter = 25
187 select case(eigens%es_type)
188 case(rs_cg)
189 !%Variable CGOrthogonalizeAll
190 !%Type logical
191 !%Default no
192 !%Section SCF::Eigensolver
193 !%Description
194 !% Used by the cg solver only.
195 !% During the cg iterations, the current band can be orthogonalized
196 !% against all other bands or only against the lower bands. Orthogonalizing
197 !% against all other bands can improve convergence properties, whereas
198 !% orthogonalizing against lower bands needs less operations.
199 !% Moreover, orthogonalizing against all bands can make converging
200 !% the highest band or unoccupied bands more difficult.
201 !%End
202 call parse_variable(namespace, 'CGOrthogonalizeAll', .false., eigens%orthogonalize_to_all)
204 !%Variable CGDirection
205 !%Type integer
206 !%Section SCF::Eigensolver
207 !%Description
208 !% Used by the cg solver only.
209 !% The conjugate direction is updated using a certain coefficient to the previous
210 !% direction. This coeffiction can be computed in different ways. The default is
211 !% to use Fletcher-Reeves (FR), an alternative is Polak-Ribiere (PR).
212 !%Option fletcher 1
213 !% The coefficient for Fletcher-Reeves consists of the current norm of the
214 !% steepest descent vector divided by that of the previous iteration.
215 !%Option polak 2
216 !% For the Polak-Ribiere scheme, a product of the current with the previous
217 !% steepest descent vector is subtracted in the nominator.
218 !%End
219 call parse_variable(namespace, 'CGDirection', option__cgdirection__fletcher, eigens%conjugate_direction)
220
221 !%Variable CGAdditionalTerms
222 !%Type logical
223 !%Section SCF::Eigensolver
224 !%Default no
225 !%Description
226 !% Used by the cg solver only.
227 !% Add additional terms during the line minimization, see PTA92, eq. 5.31ff.
228 !% These terms can improve convergence for some systems, but they are quite costly.
229 !% If you experience convergence problems, you might try out this option.
230 !% This feature is still experimental.
231 !%End
232 call parse_variable(namespace, 'CGAdditionalTerms', .false., eigens%additional_terms)
233 if(eigens%additional_terms) then
234 call messages_experimental("The additional terms for the CG eigensolver are not tested for all cases.")
235 end if
236
237 !%Variable CGEnergyChangeThreshold
238 !%Type float
239 !%Section SCF::Eigensolver
240 !%Default 0.1
241 !%Description
242 !% Used by the cg solver only.
243 !% For each band, the CG iterations are stopped when the change in energy is smaller than the
244 !% change in the first iteration multiplied by this factor. This limits the number of CG
245 !% iterations for each band, while still showing good convergence for the SCF cycle. The criterion
246 !% is discussed in Sec. V.B.6 of Payne et al. (1992), Rev. Mod. Phys. 64, 4.
247 !% The default value is 0.1, which is usually a good choice for LDA and GGA potentials. If you
248 !% are solving the OEP equation, you might want to set this value to 1e-3 or smaller. In general,
249 !% smaller values might help if you experience convergence problems.
250 !% For very small convergence tolerances, choose 0 to disable this criterion.
251 !%End
252 call parse_variable(namespace, 'CGEnergyChangeThreshold', 0.1_real64, eigens%energy_change_threshold)
253
254 case(rs_plan)
255 case(rs_evo)
256 call messages_experimental("imaginary-time evolution eigensolver")
257
258 default_iter = 1
259
260 !%Variable EigensolverImaginaryTime
261 !%Type float
262 !%Default 0.1
263 !%Section SCF::Eigensolver
264 !%Description
265 !% The imaginary-time step that is used in the imaginary-time evolution
266 !% method (<tt>Eigensolver = evolution</tt>) to obtain the lowest eigenvalues/eigenvectors.
267 !% It must satisfy <tt>EigensolverImaginaryTime > 0</tt>.
268 !% Increasing this value can make the propagation faster, but could lead to unstable propagations.
269 !%End
270 call parse_variable(namespace, 'EigensolverImaginaryTime', 0.1_real64, eigens%imag_time)
271 if(eigens%imag_time <= m_zero) call messages_input_error(namespace, 'EigensolverImaginaryTime')
272
273 call exponential_init(eigens%exponential_operator, namespace)
274
275 if(st%smear%method /= smear_semiconductor .and. st%smear%method /= smear_fixed_occ) then
276 message(1) = "Smearing of occupations is incompatible with imaginary time evolution."
277 call messages_fatal(1, namespace=namespace)
278 end if
279
280 case(rs_rmmdiis)
281 default_iter = 5
282
283 !%Variable EigensolverMinimizationIter
284 !%Type integer
285 !%Default 0
286 !%Section SCF::Eigensolver
287 !%Description
288 !% During the first iterations, the RMMDIIS eigensolver requires
289 !% some steepest-descent minimizations to improve
290 !% convergence. This variable determines the number of those
291 !% minimizations.
292 !%End
293
294 call parse_variable(namespace, 'EigensolverMinimizationIter', 0, eigens%rmmdiis_minimization_iter)
295
296 if(gr%use_curvilinear) call messages_experimental("RMMDIIS eigensolver for curvilinear coordinates")
297
298 case(rs_chebyshev)
299 !%Variable ChebyshevFilterLanczosOrder
300 !%Type integer
301 !%Default 5
302 !%Section SCF::Eigensolver
303 !%Description
304 !% Used by the Chebyshev filter only.
305 !% The number of Lanczos iterations used to construct the tridiagonal matrix,
306 !% from which the upper bound of H is estimated.
307 !% A value in the range 4 <= ChebyshevFilterLanczosOrder <= 10 is reasonable.
308 !% Values greater than 10 will raise an assertion.
309 !%End
310 call parse_variable(namespace, 'ChebyshevFilterLanczosOrder', default_chebyshev_params%n_lanczos, &
311 eigens%cheby_params%n_lanczos)
312
313 ! TODO Rewrite this, having changed the behaviour
314 !%Variable ChebyshevFilterDegree
315 !%Type integer
316 !%Default 25
317 !%Section SCF::Eigensolver
318 !%Description
319 !% Used by the Chebyshev filter only.
320 !% The degree of the Chebyshev polynomial used to dampen the interval of eigenstates one is not interested in.
321 !% If used in conjunction with "OptimizeChebyshevFilterDegree", which is the default, "ChebyshevFilterDegree" defines the
322 !% the maximum Chebyshev polynomial degree Octopus will consider when determining an approximate, optimal degree.
323 !%
324 !% If not used in conjunction with "OptimizeChebyshevFilterDegree", this defines the polynomial degree used for all SCF steps.
325 !% The larger the degree, the larger the separation between the subspaces, which will reduce the number of SCF iterations
326 !% required to reach convergence. However, ChebyshevFilterDegree also directly correlates with the number of matrix-vector
327 !% products performed on the Hamiltonian per SCF step, and so larger values result in slower SCF cycles.
328 !% A value in the range 8 <= ChebyshevFilterDegree <= 20 is reasonable.
329 !%End
330 call parse_variable(namespace, 'ChebyshevFilterDegree', default_chebyshev_params%degree, eigens%cheby_params%degree)
331
332 !%Variable OptimizeChebyshevFilterDegree
333 !%Type logical
334 !%Default yes
335 !%Section SCF::Eigensolver
336 !%Description
337 !% Used by the Chebyshev filter only.
338 !% Octopus generates a best-estimate for the degree of the Chebyshev polynomial used to filter the subspace.
339 !% A separate estimate is generated for each state block, per SCF iteration. Note that if running parallelism
340 !% over states, the block/batch containing the largest eigenstates will converge more slowly and will
341 !% therefore use a larger degree relative to all other batches. One should therefore avoid setting "ChebyshevFilterDegree"
342 !% to an excessive value (for example > 50). For more details regarding how the degree is estimated, one can refer to Section 5.4
343 !% in "Computer Physics Communications 187 (2015) 98–105" (http://dx.doi.org/10.1016/j.cpc.2014.10.015).
344 !%End
345 call parse_variable(namespace, 'OptimizeChebyshevFilterDegree', default_chebyshev_params%optimize_degree, &
346 eigens%cheby_params%optimize_degree)
347
348 !%Variable ChebyshevFilterBoundMixing
349 !%Type float
350 !%Default 0.5
351 !%Section SCF::Eigensolver
352 !%Description
353 !% In the first application of the filter, for the first SCF step, the initial lower bound estimate
354 !% is defined as a linear combination of the smallest and largest eigenvalues of the Hamiltonian.
355 !% The bound mixing determines the proportion of linear mixing, beta:
356 !% $b_{lower} = beta min{\lambda} + (\beta - 1) max{\lambda}$
357 !% of the smallest and largest eigenvalues.
358 !%End
359 call parse_variable(namespace, 'ChebyshevFilterBoundMixing', default_chebyshev_params%bound_mixing, &
360 eigens%cheby_params%bound_mixing)
361
362 !%Variable ChebyshevFilterNIter
363 !%Type integer
364 !%Default 5
365 !%Section SCF::Eigensolver
366 !%Description
367 !% The max number of iterations in the first SCF step of the Chebyshev solver. In practice this
368 !% does not need to exceed 10, as the eigenstates will vary alot between the first and second
369 !% SCF steps.
370 !%End
371 call parse_variable(namespace, 'ChebyshevFilterNIter', default_chebyshev_params%n_iter, &
372 eigens%cheby_params%n_iter)
373
374 case default
375 call messages_input_error(namespace, 'Eigensolver')
376 end select
377
378 call messages_print_with_emphasis(msg='Eigensolver', namespace=namespace)
379
380 call messages_print_var_option("Eigensolver", eigens%es_type, namespace=namespace)
381
382 call messages_obsolete_variable(namespace, 'EigensolverInitTolerance', 'EigensolverTolerance')
383 call messages_obsolete_variable(namespace, 'EigensolverFinalTolerance', 'EigensolverTolerance')
384 call messages_obsolete_variable(namespace, 'EigensolverFinalToleranceIteration')
385
386 ! this is an internal option that makes the solver use the
387 ! folded operator (H-shift)^2 to converge first eigenvalues around
388 ! the values of shift
389 ! c.f. L. W. Wang and A. Zunger
390 ! JCP 100, 2394 (1994); doi: http://dx.doi.org/10.1063/1.466486
391 eigens%folded_spectrum = .false.
392
393 !%Variable EigensolverTolerance
394 !%Type float
395 !%Section SCF::Eigensolver
396 !%Description
397 !% This is the tolerance for the eigenvectors. The default is 1e-7.
398 !% However, if <tt>ConvRelDens</tt> is defined, the default is set to a tenth of its value.
399 !%End
400 default_tol = 1e-7_real64
401 if (parse_is_defined(namespace, 'ConvRelDens')) then
402 call parse_variable(namespace, 'ConvRelDens', 1e-6_real64, default_tol)
403 default_tol = default_tol / 10.0_real64
404 end if
405 call parse_variable(namespace, 'EigensolverTolerance', default_tol, eigens%tolerance)
406
407 !%Variable EigensolverMaxIter
408 !%Type integer
409 !%Section SCF::Eigensolver
410 !%Description
411 !% Determines the maximum number of iterations that the
412 !% eigensolver will perform if the desired tolerance is not
413 !% achieved. The default is 25 iterations for all eigensolvers
414 !% except for <tt>rmdiis</tt>, which performs only 5 iterations.
415 !% Increasing this value for <tt>rmdiis</tt> increases the convergence speed,
416 !% at the cost of an increased memory footprint.
417 !%
418 !% In the case of imaginary time propatation, this variable controls the number of iterations
419 !% for which the Hxc potential is frozen. Default is 1 for the imaginary time evolution.
420 !%End
421 call parse_variable(namespace, 'EigensolverMaxIter', default_iter, eigens%es_maxiter)
422 if(eigens%es_maxiter < 1) call messages_input_error(namespace, 'EigensolverMaxIter')
423
424 if(eigens%es_maxiter > default_iter) then
425 call messages_write('You have specified a large number of eigensolver iterations (')
426 call messages_write(eigens%es_maxiter)
427 call messages_write(').', new_line = .true.)
428 call messages_write('This is not a good idea as it might slow down convergence, even for', new_line = .true.)
429 call messages_write('independent particles, as subspace diagonalization will not be used', new_line = .true.)
430 call messages_write('often enough.')
431 call messages_warning(namespace=namespace)
432 end if
433
434 if (any(eigens%es_type == (/rs_plan, rs_cg, rs_rmmdiis/))) then
435 call preconditioner_init(eigens%pre, namespace, gr, mc, space)
436 end if
437
438 safe_allocate(eigens%diff(1:st%nst, 1:st%nik))
439 eigens%diff(1:st%nst, 1:st%nik) = 0
440
441 safe_allocate(eigens%converged(1:st%nik))
442 eigens%converged(1:st%nik) = 0
443 eigens%matvec = 0
444
445 call eigens%sdiag%init(namespace, st)
446
447 ! print memory requirements
448 select case(eigens%es_type)
449 case(rs_rmmdiis)
450 call messages_write('Info: The rmmdiis eigensolver requires ')
451 mem = (m_two*eigens%es_maxiter - m_one)*st%block_size*real(gr%np_part, real64)
452 if(states_are_real(st)) then
453 mem = mem*8.0_real64
454 else
455 mem = mem*16.0_real64
456 end if
457 call messages_write(mem, units = unit_megabytes, fmt = '(f9.1)')
458 call messages_write(' of additional')
459 call messages_new_line()
460 call messages_write(' memory. This amount can be reduced by decreasing the value')
461 call messages_new_line()
462 call messages_write(' of the variable StatesBlockSize (currently set to ')
463 call messages_write(st%block_size)
464 call messages_write(').')
465 call messages_info()
466 end select
467
468 call messages_print_with_emphasis(namespace=namespace)
469
470 pop_sub(eigensolver_init)
471 end subroutine eigensolver_init
472
473
474 ! ---------------------------------------------------------
475 subroutine eigensolver_end(eigens)
476 type(eigensolver_t), intent(inout) :: eigens
477
478 push_sub(eigensolver_end)
479
480 select case(eigens%es_type)
481 case(rs_plan, rs_cg, rs_rmmdiis)
482 call preconditioner_end(eigens%pre)
483 end select
484
485 safe_deallocate_a(eigens%converged)
486 safe_deallocate_a(eigens%diff)
487
488 pop_sub(eigensolver_end)
489 end subroutine eigensolver_end
490
491
492 ! ---------------------------------------------------------
493 subroutine eigensolver_run(eigens, namespace, gr, st, hm, iter, conv, nstconv)
494 class(eigensolver_t), intent(inout) :: eigens
495 type(namespace_t), intent(in) :: namespace
496 type(grid_t), intent(in) :: gr
497 type(states_elec_t), intent(inout) :: st
498 type(hamiltonian_elec_t), intent(inout) :: hm
499 integer, intent(in) :: iter
500 logical, optional, intent(out) :: conv
501 integer, optional, intent(in) :: nstconv
502 ! !< the convergence criteria
503
504 integer :: ik, ist, nstconv_
505#ifdef HAVE_MPI
506 logical :: conv_reduced
507 integer :: outcount, lmatvec
508 real(real64), allocatable :: ldiff(:), leigenval(:)
509 real(real64), allocatable :: ldiff_out(:), leigenval_out(:)
510 integer, allocatable :: lconv(:)
511#endif
512
513 call profiling_in("EIGEN_SOLVER")
514 push_sub(eigensolver_run)
515
516 if(present(conv)) conv = .false.
517 if(present(nstconv)) then
518 nstconv_ = nstconv
519 else
520 nstconv_ = st%nst
521 end if
522
523 eigens%matvec = 0
524
525 if(mpi_grp_is_root(mpi_world) .and. eigensolver_has_progress_bar(eigens) .and. .not. debug%info) then
526 call loct_progress_bar(-1, st%lnst*st%d%kpt%nlocal)
527 end if
528
529 ik_loop: do ik = st%d%kpt%start, st%d%kpt%end
530 if (states_are_real(st)) then
531 call deigensolver_run(eigens, namespace, gr, st, hm, iter, ik)
532 else
533 call zeigensolver_run(eigens, namespace, gr, st, hm, iter, ik)
534 end if
535
536 if (.not. eigens%folded_spectrum) then
537 ! recheck convergence after subspace diagonalization, since states may have reordered
538 eigens%converged(ik) = 0
539 do ist = 1, st%nst
540 if(eigens%diff(ist, ik) < eigens%tolerance) then
541 eigens%converged(ik) = ist
542 else
543 exit
544 end if
545 end do
546 end if
547 end do ik_loop
548
549 if(mpi_grp_is_root(mpi_world) .and. eigensolver_has_progress_bar(eigens) .and. .not. debug%info) then
550 write(stdout, '(1x)')
551 end if
552
553 if(present(conv)) conv = all(eigens%converged(st%d%kpt%start:st%d%kpt%end) >= nstconv_)
554
555#ifdef HAVE_MPI
556 if (st%d%kpt%parallel) then
557 if (present(conv)) then
558 call st%d%kpt%mpi_grp%allreduce(conv, conv_reduced, 1, mpi_logical, mpi_land)
559 conv = conv_reduced
560 end if
561
562 lmatvec = eigens%matvec
563 call st%d%kpt%mpi_grp%allreduce(lmatvec, eigens%matvec, 1, mpi_integer, mpi_sum)
564
565 safe_allocate(lconv(1:st%d%kpt%nlocal))
566 lconv(1:st%d%kpt%nlocal) = eigens%converged(st%d%kpt%start:st%d%kpt%end)
567 call lmpi_gen_allgatherv(st%d%kpt%nlocal, lconv, outcount, eigens%converged, st%d%kpt%mpi_grp)
568 assert(outcount == st%nik)
569 safe_deallocate_a(lconv)
570
571 ! every node needs to know all eigenvalues (and diff)
572 safe_allocate(ldiff(1:st%d%kpt%nlocal))
573 safe_allocate(leigenval(1:st%d%kpt%nlocal))
574 safe_allocate(ldiff_out(1:st%nik))
575 safe_allocate(leigenval_out(1:st%nik))
576 do ist = st%st_start, st%st_end
577 ldiff(1:st%d%kpt%nlocal) = eigens%diff(ist, st%d%kpt%start:st%d%kpt%end)
578 leigenval(1:st%d%kpt%nlocal) = st%eigenval(ist, st%d%kpt%start:st%d%kpt%end)
579 call lmpi_gen_allgatherv(st%d%kpt%nlocal, ldiff, outcount, ldiff_out, st%d%kpt%mpi_grp)
580 eigens%diff(ist, :) = ldiff_out
581 assert(outcount == st%nik)
582 call lmpi_gen_allgatherv(st%d%kpt%nlocal, leigenval, outcount, leigenval_out, st%d%kpt%mpi_grp)
583 st%eigenval(ist, :) = leigenval_out
584 assert(outcount == st%nik)
585 end do
586 safe_deallocate_a(ldiff)
587 safe_deallocate_a(ldiff_out)
588 safe_deallocate_a(leigenval)
589 safe_deallocate_a(leigenval_out)
590 end if
591#endif
592
593 pop_sub(eigensolver_run)
594 call profiling_out("EIGEN_SOLVER")
595 end subroutine eigensolver_run
596
597
598 ! ---------------------------------------------------------
599 logical function eigensolver_parallel_in_states(this) result(par_stat)
600 type(eigensolver_t), intent(in) :: this
601
603
604 par_stat = .false.
605
606 select case(this%es_type)
608 par_stat = .true.
609 end select
610
613
614
615 ! ---------------------------------------------------------
616 logical function eigensolver_has_progress_bar(this) result(has)
617 type(eigensolver_t), intent(in) :: this
618
620
621 has = .false.
622
623 select case(this%es_type)
624 case(rs_rmmdiis, rs_cg)
625 has = .true.
626 end select
627
630
631
632#include "undef.F90"
633#include "real.F90"
634#include "eigensolver_inc.F90"
635#include "eigen_plan_inc.F90"
636#include "eigen_evolution_inc.F90"
637
638#include "undef.F90"
639#include "complex.F90"
640#include "eigensolver_inc.F90"
641#include "eigen_plan_inc.F90"
642#include "eigen_evolution_inc.F90"
643
644end module eigensolver_oct_m
645
646
647!! Local Variables:
648!! mode: f90
649!! coding: utf-8
650!! End:
This module implements batches of mesh functions.
Definition: batch.F90:133
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:116
type(debug_t), save, public debug
Definition: debug.F90:156
This module calculates the derivatives (gradients, Laplacians, etc.) of a function.
type(eigen_chebyshev_t), public, protected default_chebyshev_params
Default Chebyshev input parameters Arguments 1 and 2 taken from 10.1016/j.jcp.2006....
subroutine eigensolver_run(eigens, namespace, gr, st, hm, iter, conv, nstconv)
subroutine zeigensolver_run(eigens, namespace, mesh, st, hm, iter, ik)
subroutine, public eigensolver_init(eigens, namespace, gr, st, mc, space)
integer, parameter, public rs_evo
integer, parameter, public rs_cg
subroutine deigensolver_run(eigens, namespace, mesh, st, hm, iter, ik)
logical function eigensolver_parallel_in_states(this)
integer, parameter, public rs_chebyshev
logical function eigensolver_has_progress_bar(this)
subroutine, public eigensolver_end(eigens)
integer, parameter, public rs_rmmdiis
subroutine, public exponential_init(te, namespace, full_batch)
real(real64), parameter, public m_two
Definition: global.F90:190
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
A module to handle KS potential, without the external potential.
This module defines functions over batches of mesh functions.
Definition: mesh_batch.F90:116
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_print_with_emphasis(msg, iunit, namespace)
Definition: messages.F90:920
character(len=512), private msg
Definition: messages.F90:165
subroutine, public messages_warning(no_lines, all_nodes, namespace)
Definition: messages.F90:537
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1045
subroutine, public messages_new_line()
Definition: messages.F90:1134
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:160
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:414
subroutine, public messages_input_error(namespace, var, details, row, column)
Definition: messages.F90:713
subroutine, public messages_experimental(name, namespace)
Definition: messages.F90:1085
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:616
This module contains some common usage patterns of MPI routines.
Definition: mpi_lib.F90:115
logical function mpi_grp_is_root(grp)
Is the current MPI process of grpcomm, root.
Definition: mpi.F90:430
type(mpi_grp_t), public mpi_world
Definition: mpi.F90:266
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:145
logical function, public parse_is_defined(namespace, name)
Definition: parser.F90:502
subroutine, public preconditioner_end(this)
subroutine, public preconditioner_init(this, namespace, gr, mc, space)
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
Definition: profiling.F90:623
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
Definition: profiling.F90:552
integer, parameter, public smear_semiconductor
Definition: smear.F90:171
integer, parameter, public smear_fixed_occ
Definition: smear.F90:171
pure logical function, public states_are_real(st)
This module handles spin dimensions of the states and the k-point distribution.
This module provides routines for communicating states when using states parallelization.
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:132
This module defines the unit system, used for input and output.
type(unit_t), public unit_megabytes
For large amounts of data (natural code units are bytes)
Description of the grid, containing information on derivatives, stencil, and symmetries.
Definition: grid.F90:168
The states_elec_t class contains all electronic wave functions.
int true(void)