41 use,
intrinsic :: iso_fortran_env
110 integer,
public,
parameter :: &
118 integer,
public :: max_iter
120 real(real64),
public :: lmm_r
123 logical :: conv_eigen_error
124 logical :: check_conv
127 logical :: calc_force
128 logical,
public :: calc_stress
129 logical :: calc_dipole
130 logical :: calc_partial_charges
132 type(mixfield_t),
pointer :: mixfield
133 type(eigensolver_t) :: eigens
135 logical :: forced_finish
136 type(lda_u_mixer_t) :: lda_u_mix
137 type(vtau_mixer_t) :: vtau_mix
138 type(berry_t) :: berry
141 type(restart_t),
public :: restart_load, restart_dump
143 type(criterion_list_t),
public :: criterion_list
144 real(real64) :: energy_in, energy_diff, abs_dens_diff, evsum_in, evsum_out, evsum_diff
147 logical :: converged_current, converged_last
148 integer :: verbosity_
150 real(real64),
allocatable :: rhoout(:,:), rhoin(:,:)
151 real(real64),
allocatable :: vhxc_old(:,:)
152 class(wfs_elec_t),
allocatable :: psioutb(:, :)
153 logical :: output_forces, calc_current, output_during_scf, finish
159 subroutine scf_init(scf, namespace, gr, ions, st, mc, hm, space)
160 type(scf_t),
intent(inout) :: scf
161 type(grid_t),
intent(in) :: gr
162 type(namespace_t),
intent(in) :: namespace
163 type(ions_t),
intent(in) :: ions
164 type(states_elec_t),
intent(in) :: st
165 type(multicomm_t),
intent(in) :: mc
166 type(hamiltonian_elec_t),
intent(inout) :: hm
167 class(space_t),
intent(in) :: space
170 integer :: mixdefault
171 type(type_t) :: mix_type
172 class(convergence_criterion_t),
pointer :: crit
173 type(criterion_iterator_t) :: iter
174 logical :: deactivate_oracle
197 if (
allocated(hm%vberry))
then
204 call iter%start(scf%criterion_list)
205 do while (iter%has_next())
206 crit => iter%get_next()
209 call crit%set_pointers(scf%energy_diff, scf%energy_in)
211 call crit%set_pointers(scf%abs_dens_diff, st%qtot)
213 call crit%set_pointers(scf%evsum_diff, scf%evsum_out)
220 if(.not. scf%check_conv .and. scf%max_iter < 0)
then
221 call messages_write(
"All convergence criteria are disabled. Octopus is cowardly refusing")
226 call messages_write(
"Please set one of the following variables to a positive value:")
249 call parse_variable(namespace,
'ConvEigenError', .false., scf%conv_eigen_error)
251 if(scf%max_iter < 0) scf%max_iter = huge(scf%max_iter)
257 call eigensolver_init(scf%eigens, namespace, gr, st, hm, mc, space, deactivate_oracle)
259 if(scf%eigens%es_type /=
rs_evo)
then
283 mixdefault = option__mixfield__potential
286 call parse_variable(namespace,
'MixField', mixdefault, scf%mix_field)
290 if (scf%mix_field == option__mixfield__potential .and. hm%theory_level ==
independent_particles)
then
291 call messages_write(
'Input: Cannot mix the potential for non-interacting particles.')
295 if (scf%mix_field == option__mixfield__potential .and. hm%pcm%run_pcm)
then
296 call messages_write(
'Input: You have selected to mix the potential.', new_line = .
true.)
297 call messages_write(
' This might produce convergence problems for solvated systems.', new_line = .
true.)
302 if(scf%mix_field == option__mixfield__density &
305 call messages_write(
'Input: You have selected to mix the density with OEP or MGGA XC functionals.', new_line = .
true.)
306 call messages_write(
' This might produce convergence problems. Mix the potential instead.')
310 if(scf%mix_field == option__mixfield__states)
then
315 select case(scf%mix_field)
316 case (option__mixfield__potential, option__mixfield__density)
318 case(option__mixfield__states)
325 if (scf%mix_field /= option__mixfield__none)
then
326 call mix_init(scf%smix, namespace, space, gr%der, scf%mixdim1, st%d%nspin, func_type_ = mix_type)
330 if (scf%mix_field /= option__mixfield__states .and. scf%mix_field /= option__mixfield__none )
then
336 if(scf%mix_field == option__mixfield__potential)
then
342 scf%mix_field = option__mixfield__none
354 call parse_variable(namespace,
'SCFCalculateForces', .not. ions%only_user_def, scf%calc_force)
356 if(scf%calc_force .and. gr%der%boundaries%spiralBC)
then
357 message(1) =
'Forces cannot be calculated when using spiral boundary conditions.'
358 write(
message(2),
'(a)')
'Please use SCFCalculateForces = no.'
361 if(scf%calc_force)
then
362 if (
allocated(hm%ep%b_field) .or.
allocated(hm%ep%a_static))
then
363 write(
message(1),
'(a)')
'The forces are currently not properly calculated if static'
364 write(
message(2),
'(a)')
'magnetic fields or static vector potentials are present.'
365 write(
message(3),
'(a)')
'Please use SCFCalculateForces = no.'
378 call parse_variable(namespace,
'SCFCalculateStress', .false. , scf%calc_stress)
392 call parse_variable(namespace,
'SCFCalculateDipole', .not. space%is_periodic(), scf%calc_dipole)
393 if (
allocated(hm%vberry)) scf%calc_dipole = .
true.
403 call parse_variable(namespace,
'SCFCalculatePartialCharges', .false., scf%calc_partial_charges)
404 if (scf%calc_partial_charges)
call messages_experimental(
'SCFCalculatePartialCharges', namespace=namespace)
406 rmin = ions%min_distance()
422 scf%forced_finish = .false.
430 type(
scf_t),
intent(inout) :: scf
439 if(scf%mix_field /= option__mixfield__none)
call mix_end(scf%smix)
441 nullify(scf%mixfield)
443 if (scf%mix_field /= option__mixfield__states)
then
448 call iter%start(scf%criterion_list)
449 do while (iter%has_next())
450 crit => iter%get_next()
451 safe_deallocate_p(crit)
460 type(
scf_t),
intent(inout) :: scf
466 if (scf%mix_field /= option__mixfield__states)
then
476 subroutine scf_load(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, restart_load)
477 type(
scf_t),
intent(inout) :: scf
481 type(
grid_t),
intent(inout) :: gr
482 type(
ions_t),
intent(in) :: ions
485 type(
v_ks_t),
intent(inout) :: ks
487 type(
restart_t),
intent(in) :: restart_load
497 message(1) =
'Unable to read density. Density will be calculated from states.'
500 if (
bitand(ks%xc_family, xc_family_oep) == 0)
then
501 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
504 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
511 call hm%ks_pot%load(restart_load, space, gr, ierr)
513 message(1) =
'Unable to read Vhxc. Vhxc will be calculated from states.'
516 call hm%update(gr, namespace, space, ext_partners)
517 if (
bitand(ks%xc_family, xc_family_oep) /= 0)
then
519 do is = 1, st%d%nspin
520 ks%oep%vxc(1:gr%np, is) = hm%ks_pot%vhxc(1:gr%np, is) - hm%ks_pot%vhartree(1:gr%np)
522 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
529 if (scf%mix_field == option__mixfield__density .or. scf%mix_field == option__mixfield__potential)
then
530 call mix_load(namespace, restart_load, scf%smix, space, gr, ierr)
533 message(1) =
"Unable to read mixing information. Mixing will start from scratch."
539 call lda_u_load(restart_load, hm%lda_u, st, hm%energy%dft_u, ierr)
541 message(1) =
"Unable to read DFT+U information. DFT+U data will be calculated from states."
561 subroutine scf_start(scf, namespace, space, gr, ions, st, ks, hm, outp, verbosity)
562 type(
scf_t),
intent(inout) :: scf
565 type(
grid_t),
intent(inout) :: gr
566 type(
ions_t),
intent(inout) :: ions
568 type(
v_ks_t),
intent(inout) :: ks
570 type(
output_t),
optional,
intent(in) :: outp
571 integer,
optional,
intent(in) :: verbosity
577 if(scf%forced_finish)
then
578 message(1) =
"Previous clean stop, not doing SCF and quitting."
582 if (.not. hm%is_hermitian())
then
583 message(1) =
"Trying to run a SCF calculation for a non-hermitian Hamiltonian. This is not supported."
589 scf%output_during_scf = .false.
590 scf%output_forces = .false.
591 scf%calc_current = .false.
593 if (
present(outp))
then
596 if (outp%what(option__output__stress))
then
597 scf%calc_stress = .
true.
600 scf%output_during_scf = outp%duringscf
603 if (outp%duringscf .and. outp%what(option__output__forces))
then
604 scf%output_forces = .
true.
608 safe_allocate(scf%rhoout(1:gr%np, 1:st%d%nspin))
609 safe_allocate(scf%rhoin (1:gr%np, 1:st%d%nspin))
611 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoin)
614 if (scf%calc_force .or. scf%output_forces)
then
616 safe_allocate(scf%vhxc_old(1:gr%np, 1:st%d%nspin))
617 call lalg_copy(gr%np, st%d%nspin, hm%ks_pot%vhxc, scf%vhxc_old)
621 select case(scf%mix_field)
622 case(option__mixfield__potential)
625 case(option__mixfield__density)
628 case(option__mixfield__states)
631 allocate(
wfs_elec_t::scf%psioutb (st%group%block_start:st%group%block_end, st%d%kpt%start:st%d%kpt%end))
633 do iqn = st%d%kpt%start, st%d%kpt%end
634 do ib = st%group%block_start, st%group%block_end
635 call st%group%psib(ib, iqn)%copy_to(scf%psioutb(ib, iqn))
643 if (scf%mix_field /= option__mixfield__states)
call lda_u_mixer_set_vin(hm%lda_u, scf%lda_u_mix)
647 if ( scf%verbosity_ /= verb_no )
then
648 if(scf%max_iter > 0)
then
649 write(
message(1),
'(a)')
'Info: Starting SCF iteration.'
651 write(
message(1),
'(a)')
'Info: No SCF iterations will be done.'
658 scf%converged_current = .false.
668 character(len=*),
intent(in) :: dir
669 character(len=*),
intent(in) :: fname
672 character(len=12) :: label
675 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write')
676 write(iunit,
'(a)', advance =
'no')
'#iter energy '
677 label =
'energy_diff'
678 write(iunit,
'(1x,a)', advance =
'no') label
680 write(iunit,
'(1x,a)', advance =
'no') label
682 write(iunit,
'(1x,a)', advance =
'no') label
684 write(iunit,
'(1x,a)', advance =
'no') label
686 write(iunit,
'(1x,a)', advance =
'no') label
690 label =
'OEP norm2ss'
691 write(iunit,
'(1x,a)', advance =
'no') label
694 write(iunit,
'(a)')
''
704 subroutine scf_run(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, &
705 verbosity, iters_done, restart_dump)
706 type(
scf_t),
intent(inout) :: scf
710 type(
grid_t),
intent(inout) :: gr
711 type(
ions_t),
intent(inout) :: ions
714 type(
v_ks_t),
intent(inout) :: ks
716 type(
output_t),
optional,
intent(in) :: outp
717 integer,
optional,
intent(in) :: verbosity
718 integer,
optional,
intent(out) :: iters_done
719 type(
restart_t),
optional,
intent(in) :: restart_dump
726 call scf_start(scf, namespace, space, gr, ions, st, ks, hm, outp, verbosity)
729 do iter = 1, scf%max_iter
731 call scf_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
732 verbosity, restart_dump)
734 completed =
scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, verbosity, iters_done)
736 if(scf%forced_finish .or. completed)
then
741 call scf_finish(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
742 verbosity, iters_done, restart_dump)
748 subroutine scf_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
749 verbosity, restart_dump)
750 type(
scf_t),
intent(inout) :: scf
754 type(
grid_t),
intent(inout) :: gr
755 type(
ions_t),
intent(inout) :: ions
758 type(
v_ks_t),
intent(inout) :: ks
760 integer,
intent(in) :: iter
761 type(
output_t),
optional,
intent(in) :: outp
762 integer,
optional,
intent(in) :: verbosity
763 type(
restart_t),
optional,
intent(in) :: restart_dump
765 integer :: iqn, ib, ierr
768 logical :: is_crit_conv
769 real(real64) :: etime, itime
778 scf%eigens%converged = 0
781 call hm%update_span(gr%spacing(1:space%dim), minval(st%eigenval(:, :)), namespace)
787 call iterator%start(scf%criterion_list)
788 do while (iterator%has_next())
789 crit => iterator%get_next()
793 if (scf%calc_force .or. scf%output_forces)
then
795 scf%vhxc_old(1:gr%np, 1:st%d%nspin) = hm%ks_pot%vhxc(1:gr%np, 1:st%d%nspin)
800 if (
allocated(hm%vberry))
then
803 ks%frozen_hxc = .
true.
805 call berry_perform_internal_scf(scf%berry, namespace, space, scf%eigens, gr, st, hm, iter, ks, ions, ext_partners)
807 ks%frozen_hxc = .false.
809 scf%eigens%converged = 0
810 call scf%eigens%run(namespace, gr, st, hm, space, ext_partners, iter)
813 scf%matvec = scf%matvec + scf%eigens%matvec
822 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoout)
824 select case (scf%mix_field)
825 case (option__mixfield__potential)
826 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
829 case (option__mixfield__density)
831 case(option__mixfield__states)
832 do iqn = st%d%kpt%start, st%d%kpt%end
833 do ib = st%group%block_start, st%group%block_end
834 call st%group%psib(ib, iqn)%copy_data_to(gr%np, scf%psioutb(ib, iqn))
839 if (scf%mix_field /= option__mixfield__states .and. scf%mix_field /= option__mixfield__none)
then
846 if (
present(outp))
then
848 if (outp%duringscf .and. outp%what_now(option__output__forces, iter))
then
849 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old=scf%vhxc_old)
854 call iterator%start(scf%criterion_list)
855 do while (iterator%has_next())
856 crit => iterator%get_next()
861 scf%converged_last = scf%converged_current
863 scf%converged_current = scf%check_conv .and. &
864 (.not. scf%conv_eigen_error .or. all(scf%eigens%converged >= st%nst_conv))
866 call iterator%start(scf%criterion_list)
867 do while (iterator%has_next())
868 crit => iterator%get_next()
869 call crit%is_converged(is_crit_conv)
870 scf%converged_current = scf%converged_current .and. is_crit_conv
875 scf%finish = scf%converged_last .and. scf%converged_current
881 select case (scf%mix_field)
882 case (option__mixfield__density)
884 call mixing(namespace, scf%smix)
887 if (minval(st%rho(1:gr%np, 1:st%d%spin_channels)) < -1e-6_real64)
then
888 write(
message(1),*)
'Negative density after mixing. Minimum value = ', &
889 minval(st%rho(1:gr%np, 1:st%d%spin_channels))
893 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
895 case (option__mixfield__potential)
897 call mixing(namespace, scf%smix)
903 case(option__mixfield__states)
904 do iqn = st%d%kpt%start, st%d%kpt%end
905 do ib = st%group%block_start, st%group%block_end
911 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
913 case (option__mixfield__none)
914 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
920 if (scf%finish .and. st%modelmbparticles%nparticle > 0)
then
924 if (
present(outp) .and.
present(restart_dump))
then
927 if ( (scf%finish .or. (modulo(iter, outp%restart_write_interval) == 0) &
928 .or. iter == scf%max_iter .or. scf%forced_finish) )
then
930 call states_elec_dump(scf%restart_dump, space, st, gr, hm%kpoints, ierr, iter=iter)
932 message(1) =
'Unable to write states wavefunctions.'
938 message(1) =
'Unable to write density.'
943 call lda_u_dump(scf%restart_dump, namespace, hm%lda_u, st, gr, ierr)
945 message(1) =
'Unable to write DFT+U information.'
950 select case (scf%mix_field)
951 case (option__mixfield__density)
952 call mix_dump(namespace, scf%restart_dump, scf%smix, space, gr, ierr)
954 message(1) =
'Unable to write mixing information.'
957 case (option__mixfield__potential)
958 call hm%ks_pot%dump(scf%restart_dump, space, gr, ierr)
960 message(1) =
'Unable to write Vhxc.'
964 call mix_dump(namespace, scf%restart_dump, scf%smix, space, gr, ierr)
966 message(1) =
'Unable to write mixing information.'
985 character(len=50) :: str
986 real(real64) :: dipole(1:space%dim)
992 write(str,
'(a,i5)')
'SCF CYCLE ITER #' ,iter
996 ' rel_ev = ', scf%evsum_diff/(abs(scf%evsum_out)+1e-20)
997 write(
message(2),
'(a,es15.2,2(a,es9.2))') &
998 ' ediff = ', scf%energy_diff,
' abs_dens = ', scf%abs_dens_diff, &
999 ' rel_dens = ', scf%abs_dens_diff/st%qtot
1002 write(
message(1),
'(a,i6)')
'Matrix vector products: ', scf%eigens%matvec
1003 write(
message(2),
'(a,i6)')
'Converged eigenvectors: ', sum(scf%eigens%converged(1:st%nik))
1007 if (
allocated(hm%vberry))
then
1009 call write_dipole(st, hm, space, dipole, namespace=namespace)
1022 write(
message(2),
'(a,i5,a,f14.2)')
'Elapsed time for SCF step ', iter,
':', etime
1032 write(
message(1),
'(a,i4,a,es15.8, a,es9.2, a, f7.1, a)') &
1035 ' : abs_dens', scf%abs_dens_diff, &
1036 ' : etime ', etime,
's'
1046 character(len=*),
intent(in) :: dir
1047 character(len=*),
intent(in) :: fname
1053 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write', position=
'append')
1055 call iterator%start(scf%criterion_list)
1056 do while (iterator%has_next())
1057 crit => iterator%get_next()
1062 write(iunit,
'(2es13.5)', advance =
'no') crit%val_abs, crit%val_rel
1065 write(iunit,
'(es13.5)', advance =
'no') crit%val_rel
1073 write(iunit,
'(es13.5)', advance =
'no') ks%oep%norm2ss
1076 write(iunit,
'(a)')
''
1083 logical function scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, &
1084 verbosity, iters_done)
result(completed)
1085 type(
scf_t),
intent(inout) :: scf
1088 type(
grid_t),
intent(inout) :: gr
1089 type(
ions_t),
intent(inout) :: ions
1091 type(
v_ks_t),
intent(inout) :: ks
1093 integer,
intent(in) :: iter
1094 type(
output_t),
optional,
intent(in) :: outp
1095 integer,
optional,
intent(in) :: verbosity
1096 integer,
optional,
intent(out) :: iters_done
1098 character(len=MAX_PATH_LEN) :: dirname
1099 integer(int64) :: what_i
1106 if(
present(iters_done)) iters_done = iter
1108 write(
message(1),
'(a, i4, a)')
'Info: SCF converged in ', iter,
' iterations'
1116 if (
present(outp))
then
1117 if (any(outp%what) .and. outp%duringscf)
then
1118 do what_i = lbound(outp%what, 1), ubound(outp%what, 1)
1119 if (outp%what_now(what_i, iter))
then
1120 write(dirname,
'(a,a,i4.4)') trim(outp%iter_dir),
"scf.", iter
1121 call output_all(outp, namespace, space, dirname, gr, ions, iter, st, hm, ks)
1122 call output_modelmb(outp, namespace, space, dirname, gr, ions, iter, st)
1130 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoin)
1133 if (scf%mix_field /= option__mixfield__none)
then
1134 if (scf%smix%ns_restart > 0)
then
1135 if (mod(iter, scf%smix%ns_restart) == 0)
then
1136 message(1) =
"Info: restarting mixing."
1143 select case(scf%mix_field)
1144 case(option__mixfield__potential)
1147 case (option__mixfield__density)
1152 if (scf%mix_field /= option__mixfield__states)
then
1163 subroutine scf_finish(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
1164 verbosity, iters_done, restart_dump)
1165 type(
scf_t),
intent(inout) :: scf
1169 type(
grid_t),
intent(inout) :: gr
1170 type(
ions_t),
intent(inout) :: ions
1173 type(
v_ks_t),
intent(inout) :: ks
1175 integer,
intent(in) :: iter
1176 type(
output_t),
optional,
intent(in) :: outp
1177 integer,
optional,
intent(in) :: verbosity
1178 integer,
optional,
intent(out) :: iters_done
1179 type(
restart_t),
optional,
intent(in) :: restart_dump
1190 if ((scf%max_iter > 0 .and. scf%mix_field == option__mixfield__potential) .or. scf%calc_current)
then
1191 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, &
1192 calc_current=scf%calc_current)
1195 select case(scf%mix_field)
1196 case(option__mixfield__states)
1198 do iqn = st%d%kpt%start, st%d%kpt%end
1199 do ib = st%group%block_start, st%group%block_end
1200 call scf%psioutb(ib, iqn)%end()
1205 deallocate(scf%psioutb)
1208 safe_deallocate_a(scf%rhoout)
1209 safe_deallocate_a(scf%rhoin)
1211 if (scf%max_iter > 0 .and. any(scf%eigens%converged < st%nst))
then
1212 write(
message(1),
'(a)')
'Some of the states are not fully converged!'
1213 if (all(scf%eigens%converged >= st%nst_conv))
then
1214 write(
message(2),
'(a)')
'But all requested states to converge are converged.'
1221 if (.not.scf%finish)
then
1222 write(
message(1),
'(a,i4,a)')
'SCF *not* converged after ', iter - 1,
' iterations.'
1226 write(
message(1),
'(a,i10)')
'Info: Number of matrix-vector products: ', scf%matvec
1229 if (scf%calc_force)
then
1230 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old=scf%vhxc_old)
1233 if (scf%calc_stress)
call stress_calculate(namespace, gr, hm, st, ions, ks, ext_partners)
1236 if (scf%mix_field == option__mixfield__potential)
then
1241 if(
present(outp))
then
1248 if (space%is_periodic() .and. st%nik > st%d%nspin)
then
1251 ions, gr, hm%kpoints, hm%phase, vec_pot = hm%hm_base%uniform_vector_potential, &
1252 vec_pot_var = hm%hm_base%vector_potential)
1256 if (ks%vdw%vdw_correction == option__vdwcorrection__vdw_ts)
then
1260 safe_deallocate_a(scf%vhxc_old)
1268 character(len=*),
intent(in) :: dir, fname
1271 real(real64) :: dipole(1:space%dim)
1272 real(real64) :: ex_virial
1278 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write')
1284 if (space%is_periodic())
then
1285 call hm%kpoints%write_info(iunit=iunit)
1293 write(iunit,
'(a, i4, a)')
'SCF converged in ', iter,
' iterations'
1295 write(iunit,
'(a)')
'SCF *not* converged!'
1297 write(iunit,
'(1x)')
1299 if(any(scf%eigens%converged < st%nst))
then
1300 write(iunit,
'(a)')
'Some of the states are not fully converged!'
1301 if (all(scf%eigens%converged >= st%nst_conv))
then
1302 write(iunit,
'(a)')
'But all requested states to converge are converged.'
1307 write(iunit,
'(1x)')
1309 if (space%is_periodic())
then
1311 write(iunit,
'(1x)')
1321 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1324 if (
mpi_world%is_root())
write(iunit,
'(1x)')
1327 if(st%d%ispin ==
spinors .and. space%dim == 3 .and. &
1330 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1336 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1339 if(scf%calc_dipole)
then
1346 hm%xc%functional(
func_c,1)%family == xc_family_none .and. st%d%ispin /=
spinors)
then
1353 write(iunit,
'(1x)')
1358 if(scf%max_iter > 0)
then
1359 write(iunit,
'(a)')
'Convergence:'
1360 call iterator%start(scf%criterion_list)
1361 do while (iterator%has_next())
1362 crit => iterator%get_next()
1363 call crit%write_info(iunit)
1372 write(iunit,
'(a)')
'Photon observables:'
1373 write(iunit,
'(6x, a, es15.8,a,es15.8,a)')
'Photon number = ', ks%oep_photon%pt%number(1)
1374 write(iunit,
'(6x, a, es15.8,a,es15.8,a)')
'Photon ex. = ', ks%oep_photon%pt%ex
1381 if (scf%calc_stress)
then
1382 call output_stress(iunit, space%periodic_dim, st%stress_tensors, all_terms=.false.)
1383 call output_pressure(iunit, space%periodic_dim, st%stress_tensors%total)
1388 if(scf%calc_partial_charges)
then
1423 real(real64) :: mem_tmp
1427 if(
conf%report_memory)
then
1429 call mpi_world%allreduce(mem, mem_tmp, 1, mpi_double_precision, mpi_sum)
1431 write(
message(1),
'(a,f14.2)')
'Memory usage [Mbytes] :', mem
1441 type(
scf_t),
intent(inout) :: scf
1447 select type (criterion)
1449 scf%energy_in = hm%energy%total
1464 type(
scf_t),
intent(inout) :: scf
1467 type(
grid_t),
intent(in) :: gr
1468 real(real64),
intent(in) :: rhoout(:,:), rhoin(:,:)
1472 real(real64),
allocatable :: tmp(:)
1476 select type (criterion)
1478 scf%energy_diff = abs(hm%energy%total - scf%energy_in)
1481 scf%abs_dens_diff =
m_zero
1482 safe_allocate(tmp(1:gr%np))
1483 do is = 1, st%d%nspin
1484 tmp(:) = abs(rhoin(1:gr%np, is) - rhoout(1:gr%np, is))
1485 scf%abs_dens_diff = scf%abs_dens_diff +
dmf_integrate(gr, tmp)
1487 safe_deallocate_a(tmp)
1491 scf%evsum_diff = abs(scf%evsum_out - scf%evsum_in)
1492 scf%evsum_in = scf%evsum_out
1502 subroutine write_dipole(st, hm, space, dipole, iunit, namespace)
1506 real(real64),
intent(in) :: dipole(:)
1507 integer,
optional,
intent(in) :: iunit
1508 type(
namespace_t),
optional,
intent(in) :: namespace
1513 call output_dipole(dipole, space%dim, iunit=iunit, namespace=namespace)
1515 if (space%is_periodic())
then
1516 message(1) =
"Defined only up to quantum of polarization (e * lattice vector)."
1517 message(2) =
"Single-point Berry's phase method only accurate for large supercells."
1520 if (hm%kpoints%full%npoints > 1)
then
1522 "WARNING: Single-point Berry's phase method for dipole should not be used when there is more than one k-point."
1523 message(2) =
"Instead, finite differences on k-points (not yet implemented) are needed."
1528 message(1) =
"Single-point Berry's phase dipole calculation not correct without integer occupations."
batchified version of the BLAS axpy routine:
scale a batch by a constant or vector
Copies a vector x, to a vector y.
This module implements common operations on batches of mesh functions.
subroutine, public berry_perform_internal_scf(this, namespace, space, eigensolver, gr, st, hm, iter, ks, ions, ext_partners)
subroutine, public berry_init(this, namespace)
subroutine, public calc_dipole(dipole, space, mesh, st, ions)
subroutine, public criteria_factory_init(list, namespace, check_conv)
This module implements a calculator for the density and defines related functions.
subroutine, public density_calc(st, gr, density, istin)
Computes the density from the orbitals in st.
integer, parameter, public rs_evo
subroutine, public eigensolver_init(eigens, namespace, gr, st, hm, mc, space, deactivate_oracle)
subroutine, public eigensolver_end(eigens)
integer, parameter, public unpolarized
Parameters...
integer, parameter, public spinors
subroutine, public energy_calc_total(namespace, space, hm, gr, st, ext_partners, iunit, full)
This subroutine calculates the total energy of the system. Basically, it adds up the KS eigenvalues,...
subroutine, public energy_calc_virial_ex(der, vxc, st, ex)
subroutine, public energy_calc_eigenvalues(namespace, hm, der, st)
subroutine, public forces_write_info(iunit, ions, dir, namespace)
subroutine, public forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old, t, dt)
real(real64), parameter, public m_zero
integer, parameter, public hartree_fock
integer, parameter, public independent_particles
Theory level.
integer, parameter, public generalized_kohn_sham_dft
real(real64), parameter, public lmm_r_single_atom
Default local magnetic moments sphere radius for an isolated system.
integer, parameter, public kohn_sham_dft
type(conf_t), public conf
Global instance of Octopus configuration.
character(len= *), parameter, public static_dir
real(real64), parameter, public m_half
real(real64), parameter, public m_one
This module implements the underlying real-space grid.
subroutine, public grid_write_info(gr, iunit, namespace)
subroutine, public hamiltonian_elec_update_pot(this, mesh, accumulate)
Update the KS potential of the electronic Hamiltonian.
This module defines classes and functions for interaction partners.
subroutine, public io_close(iunit, grp)
subroutine, public io_debug_on_the_fly(namespace)
check if debug mode should be enabled or disabled on the fly
subroutine, public io_mkdir(fname, namespace, parents)
integer function, public io_open(file, namespace, action, status, form, position, die, recl, grp)
integer, parameter, public kpoints_path
A module to handle KS potential, without the external potential.
subroutine, public lda_u_dump(restart, namespace, this, st, mesh, ierr)
subroutine, public lda_u_write_u(this, iunit, namespace)
subroutine, public lda_u_load(restart, this, st, dftu_energy, ierr, occ_only, u_only)
subroutine, public lda_u_write_v(this, iunit, namespace)
subroutine, public lda_u_mixer_set_vin(this, mixer)
subroutine, public lda_u_mixer_init(this, mixer, st)
subroutine, public lda_u_mixer_clear(mixer, smix)
subroutine, public lda_u_mixer_init_auxmixer(this, namespace, mixer, smix, st)
subroutine, public lda_u_mixer_get_vnew(this, mixer, st)
subroutine, public lda_u_mixer_set_vout(this, mixer)
subroutine, public lda_u_mixer_end(mixer, smix)
integer, parameter, public dft_u_none
integer, parameter, public dft_u_acbn0
subroutine, public lda_u_update_occ_matrices(this, namespace, mesh, st, hm_base, phase, energy)
subroutine, public write_magnetic_moments(mesh, st, ions, boundaries, lmm_r, iunit, namespace)
subroutine, public write_total_xc_torque(iunit, mesh, vxc, st)
This module is intended to contain "only mathematical" functions and procedures.
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
character(len=512), private msg
subroutine, public messages_warning(no_lines, all_nodes, namespace)
subroutine, public messages_obsolete_variable(namespace, name, rep)
subroutine, public messages_new_line()
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
subroutine, public messages_input_error(namespace, var, details, row, column)
subroutine, public messages_experimental(name, namespace)
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
real(real64) pure function, public mix_coefficient(this)
subroutine, public mixing(namespace, smix)
Main entry-point to SCF mixer.
subroutine, public mix_get_field(this, mixfield)
subroutine, public mix_dump(namespace, restart, smix, space, mesh, ierr)
subroutine, public mix_load(namespace, restart, smix, space, mesh, ierr)
subroutine, public mix_init(smix, namespace, space, der, d1, d2, def_, func_type_, prefix_)
Initialise mix_t instance.
subroutine, public mix_end(smix)
subroutine, public mix_clear(smix)
subroutine, public modelmb_sym_all_states(space, mesh, st)
type(mpi_grp_t), public mpi_world
This module handles the communicators for the various parallelization strategies.
this module contains the low-level part of the output system
subroutine, public output_modelmb(outp, namespace, space, dir, gr, ions, iter, st)
this module contains the output system
logical function, public output_needs_current(outp, states_are_real)
subroutine, public output_all(outp, namespace, space, dir, gr, ions, iter, st, hm, ks)
subroutine, public partial_charges_compute_and_print_charges(mesh, st, ions, iunit)
Computes and write partial charges to a file.
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
logical function, public clean_stop(comm)
returns true if a file named stop exists
integer, parameter, public restart_flag_mix
integer, parameter, public restart_flag_rho
integer, parameter, public restart_flag_vhxc
subroutine, public scf_load(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, restart_load)
Loading of restarting data of the SCF cycle.
logical function, public scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, verbosity, iters_done)
subroutine write_dipole(st, hm, space, dipole, iunit, namespace)
subroutine scf_update_initial_quantity(scf, hm, criterion)
Update the quantity at the begining of a SCF cycle.
subroutine scf_update_diff_quantity(scf, hm, st, gr, rhoout, rhoin, criterion)
Update the quantity at the begining of a SCF cycle.
subroutine, public scf_state_info(namespace, st)
subroutine, public scf_print_mem_use(namespace)
subroutine, public scf_mix_clear(scf)
integer, parameter, public verb_full
subroutine, public scf_finish(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, verbosity, iters_done, restart_dump)
subroutine, public scf_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, verbosity, restart_dump)
integer, parameter, public verb_compact
subroutine, public scf_init(scf, namespace, gr, ions, st, mc, hm, space)
subroutine, public scf_end(scf)
subroutine, public scf_run(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, verbosity, iters_done, restart_dump)
Legacy version of the SCF code.
subroutine, public scf_start(scf, namespace, space, gr, ions, st, ks, hm, outp, verbosity)
Preparation of the SCF cycle.
logical pure function, public smear_is_semiconducting(this)
pure logical function, public states_are_real(st)
This module defines routines to write information about states.
subroutine, public states_elec_write_eigenvalues(nst, st, space, kpoints, error, st_start, compact, iunit, namespace)
write the eigenvalues for some states to a file.
subroutine, public states_elec_write_gaps(iunit, st, space)
calculate gaps and write to a file.
subroutine, public states_elec_write_bandstructure(dir, namespace, nst, st, ions, mesh, kpoints, phase, vec_pot, vec_pot_var)
calculate and write the bandstructure
subroutine, public states_elec_fermi(st, namespace, mesh, compute_spin)
calculate the Fermi level for the states in this object
real(real64) function, public states_elec_eigenvalues_sum(st, alt_eig)
function to calculate the eigenvalues sum using occupations as weights
This module handles reading and writing restart information for the states_elec_t.
subroutine, public states_elec_dump(restart, space, st, mesh, kpoints, ierr, iter, lr, verbose)
subroutine, public states_elec_load_rho(restart, space, st, mesh, ierr)
subroutine, public states_elec_dump_rho(restart, space, st, mesh, ierr, iter)
This module implements the calculation of the stress tensor.
subroutine, public output_pressure(iunit, space_dim, total_stress_tensor)
subroutine, public stress_calculate(namespace, gr, hm, st, ions, ks, ext_partners)
This computes the total stress on the lattice.
subroutine, public output_stress(iunit, space_dim, stress_tensors, all_terms)
subroutine, public symmetries_write_info(this, space, iunit, namespace)
type(type_t), public type_float
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
character(len=20) pure function, public units_abbrev(this)
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
type(unit_system_t), public units_inp
the units systems for reading and writing
This module is intended to contain simple general-purpose utility functions and procedures.
subroutine, public output_dipole(dipole, ndim, iunit, namespace)
subroutine, public v_ks_write_info(ks, iunit, namespace)
subroutine, public v_ks_update_dftu_energy(ks, namespace, hm, st, int_dft_u)
Update the value of <\psi | V_U | \psi>, where V_U is the DFT+U potential.
subroutine, public v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval, time, calc_energy, calc_current, force_semilocal)
Tkatchenko-Scheffler pairwise method for van der Waals (vdW, dispersion) interactions.
subroutine, public vdw_ts_write_c6ab(this, ions, dir, fname, namespace)
subroutine, public vtau_mixer_end(mixer, smix)
subroutine, public vtau_mixer_init_auxmixer(namespace, mixer, smix, hm, np, nspin)
subroutine, public vtau_mixer_set_vout(mixer, hm)
subroutine, public vtau_mixer_get_vnew(mixer, hm)
subroutine, public vtau_mixer_clear(mixer, smix)
subroutine, public vtau_mixer_set_vin(mixer, hm)
This module provices a simple timer class which can be used to trigger the writing of a restart file ...
logical function, public walltimer_alarm(comm, print)
indicate whether time is up
integer, parameter, public xc_family_nc_mgga
integer, parameter, public func_c
integer, parameter, public oep_level_full
integer, parameter, public oep_level_kli
subroutine scf_write_static(dir, fname)
subroutine create_convergence_file(dir, fname)
subroutine scf_write_iter(namespace)
subroutine write_convergence_file(dir, fname)
Extension of space that contains the knowledge of the spin dimension.
Description of the grid, containing information on derivatives, stencil, and symmetries.
Stores all communicators and groups.
some variables used for the SCF cycle
abstract class for states
The states_elec_t class contains all electronic wave functions.
batches of electronic states