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, gr, ions, ext_partners, st, ks, hm, restart_load)
477 type(
scf_t),
intent(inout) :: scf
480 type(
grid_t),
intent(inout) :: gr
481 type(
ions_t),
intent(in) :: ions
484 type(
v_ks_t),
intent(inout) :: ks
486 type(
restart_t),
intent(in) :: restart_load
496 message(1) =
'Unable to read density. Density will be calculated from states.'
499 if (
bitand(ks%xc_family, xc_family_oep) == 0)
then
500 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
503 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
510 call hm%ks_pot%load(restart_load, gr, ierr)
512 message(1) =
'Unable to read Vhxc. Vhxc will be calculated from states.'
515 call hm%update(gr, namespace, space, ext_partners)
516 if (
bitand(ks%xc_family, xc_family_oep) /= 0)
then
518 do is = 1, st%d%nspin
519 ks%oep%vxc(1:gr%np, is) = hm%ks_pot%vhxc(1:gr%np, is) - hm%ks_pot%vhartree(1:gr%np)
521 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners)
528 if (scf%mix_field == option__mixfield__density .or. scf%mix_field == option__mixfield__potential)
then
529 call mix_load(namespace, restart_load, scf%smix, gr, ierr)
532 message(1) =
"Unable to read mixing information. Mixing will start from scratch."
538 call lda_u_load(restart_load, hm%lda_u, st, hm%energy%dft_u, ierr)
540 message(1) =
"Unable to read DFT+U information. DFT+U data will be calculated from states."
560 subroutine scf_start(scf, namespace, gr, ions, st, ks, hm, outp, verbosity)
561 type(
scf_t),
intent(inout) :: scf
563 type(
grid_t),
intent(inout) :: gr
564 type(
ions_t),
intent(inout) :: ions
566 type(
v_ks_t),
intent(inout) :: ks
568 type(
output_t),
optional,
intent(in) :: outp
569 integer,
optional,
intent(in) :: verbosity
575 if(scf%forced_finish)
then
576 message(1) =
"Previous clean stop, not doing SCF and quitting."
580 if (.not. hm%is_hermitian())
then
581 message(1) =
"Trying to run a SCF calculation for a non-hermitian Hamiltonian. This is not supported."
587 scf%output_during_scf = .false.
588 scf%output_forces = .false.
589 scf%calc_current = .false.
591 if (
present(outp))
then
594 if (outp%what(option__output__stress))
then
595 scf%calc_stress = .
true.
598 scf%output_during_scf = outp%duringscf
601 if (outp%duringscf .and. outp%what(option__output__forces))
then
602 scf%output_forces = .
true.
606 safe_allocate(scf%rhoout(1:gr%np, 1:st%d%nspin))
607 safe_allocate(scf%rhoin (1:gr%np, 1:st%d%nspin))
609 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoin)
612 if (scf%calc_force .or. scf%output_forces)
then
614 safe_allocate(scf%vhxc_old(1:gr%np, 1:st%d%nspin))
615 call lalg_copy(gr%np, st%d%nspin, hm%ks_pot%vhxc, scf%vhxc_old)
619 select case(scf%mix_field)
620 case(option__mixfield__potential)
623 case(option__mixfield__density)
626 case(option__mixfield__states)
629 allocate(
wfs_elec_t::scf%psioutb (st%group%block_start:st%group%block_end, st%d%kpt%start:st%d%kpt%end))
631 do iqn = st%d%kpt%start, st%d%kpt%end
632 do ib = st%group%block_start, st%group%block_end
633 call st%group%psib(ib, iqn)%copy_to(scf%psioutb(ib, iqn))
641 if (scf%mix_field /= option__mixfield__states)
call lda_u_mixer_set_vin(hm%lda_u, scf%lda_u_mix)
645 if ( scf%verbosity_ /= verb_no )
then
646 if(scf%max_iter > 0)
then
647 write(
message(1),
'(a)')
'Info: Starting SCF iteration.'
649 write(
message(1),
'(a)')
'Info: No SCF iterations will be done.'
656 scf%converged_current = .false.
666 character(len=*),
intent(in) :: dir
667 character(len=*),
intent(in) :: fname
670 character(len=12) :: label
673 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write')
674 write(iunit,
'(a)', advance =
'no')
'#iter energy '
675 label =
'energy_diff'
676 write(iunit,
'(1x,a)', advance =
'no') label
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
688 label =
'OEP norm2ss'
689 write(iunit,
'(1x,a)', advance =
'no') label
692 write(iunit,
'(a)')
''
702 subroutine scf_run(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, &
703 verbosity, iters_done, restart_dump)
704 type(
scf_t),
intent(inout) :: scf
708 type(
grid_t),
intent(inout) :: gr
709 type(
ions_t),
intent(inout) :: ions
712 type(
v_ks_t),
intent(inout) :: ks
714 type(
output_t),
optional,
intent(in) :: outp
715 integer,
optional,
intent(in) :: verbosity
716 integer,
optional,
intent(out) :: iters_done
717 type(
restart_t),
optional,
intent(in) :: restart_dump
724 call scf_start(scf, namespace, gr, ions, st, ks, hm, outp, verbosity)
727 do iter = 1, scf%max_iter
729 call scf_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
732 completed =
scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, iters_done)
734 if(scf%forced_finish .or. completed)
then
739 call scf_finish(scf, namespace, space, gr, ions, ext_partners, st, ks, hm, iter, outp)
745 subroutine scf_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, &
747 type(
scf_t),
intent(inout) :: scf
751 type(
grid_t),
intent(inout) :: gr
752 type(
ions_t),
intent(inout) :: ions
755 type(
v_ks_t),
intent(inout) :: ks
757 integer,
intent(in) :: iter
758 type(
output_t),
optional,
intent(in) :: outp
759 type(
restart_t),
optional,
intent(in) :: restart_dump
761 integer :: iqn, ib, ierr
764 logical :: is_crit_conv
765 real(real64) :: etime, itime
774 scf%eigens%converged = 0
777 call hm%update_span(gr%spacing(1:space%dim), minval(st%eigenval(:, :)), namespace)
783 call iterator%start(scf%criterion_list)
784 do while (iterator%has_next())
785 crit => iterator%get_next()
789 if (scf%calc_force .or. scf%output_forces)
then
791 scf%vhxc_old(1:gr%np, 1:st%d%nspin) = hm%ks_pot%vhxc(1:gr%np, 1:st%d%nspin)
796 if (
allocated(hm%vberry))
then
799 ks%frozen_hxc = .
true.
801 call berry_perform_internal_scf(scf%berry, namespace, space, scf%eigens, gr, st, hm, iter, ks, ions, ext_partners)
803 ks%frozen_hxc = .false.
805 scf%eigens%converged = 0
806 call scf%eigens%run(namespace, gr, st, hm, space, ext_partners, iter)
809 scf%matvec = scf%matvec + scf%eigens%matvec
818 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoout)
820 select case (scf%mix_field)
821 case (option__mixfield__potential)
822 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
825 case (option__mixfield__density)
827 case(option__mixfield__states)
828 do iqn = st%d%kpt%start, st%d%kpt%end
829 do ib = st%group%block_start, st%group%block_end
830 call st%group%psib(ib, iqn)%copy_data_to(gr%np, scf%psioutb(ib, iqn))
835 if (scf%mix_field /= option__mixfield__states .and. scf%mix_field /= option__mixfield__none)
then
842 if (
present(outp))
then
844 if (outp%duringscf .and. outp%what_now(option__output__forces, iter))
then
845 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old=scf%vhxc_old)
850 call iterator%start(scf%criterion_list)
851 do while (iterator%has_next())
852 crit => iterator%get_next()
857 scf%converged_last = scf%converged_current
859 scf%converged_current = scf%check_conv .and. &
860 (.not. scf%conv_eigen_error .or. all(scf%eigens%converged >= st%nst_conv))
862 call iterator%start(scf%criterion_list)
863 do while (iterator%has_next())
864 crit => iterator%get_next()
865 call crit%is_converged(is_crit_conv)
866 scf%converged_current = scf%converged_current .and. is_crit_conv
871 scf%finish = scf%converged_last .and. scf%converged_current
877 select case (scf%mix_field)
878 case (option__mixfield__density)
880 call mixing(namespace, scf%smix)
883 if (minval(st%rho(1:gr%np, 1:st%d%spin_channels)) < -1e-6_real64)
then
884 write(
message(1),*)
'Negative density after mixing. Minimum value = ', &
885 minval(st%rho(1:gr%np, 1:st%d%spin_channels))
889 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
891 case (option__mixfield__potential)
893 call mixing(namespace, scf%smix)
899 case(option__mixfield__states)
900 do iqn = st%d%kpt%start, st%d%kpt%end
901 do ib = st%group%block_start, st%group%block_end
907 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
909 case (option__mixfield__none)
910 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_current=scf%output_during_scf)
916 if (scf%finish .and. st%modelmbparticles%nparticle > 0)
then
920 if (
present(outp) .and.
present(restart_dump))
then
923 if ( (scf%finish .or. (modulo(iter, outp%restart_write_interval) == 0) &
924 .or. iter == scf%max_iter .or. scf%forced_finish) )
then
926 call states_elec_dump(scf%restart_dump, space, st, gr, hm%kpoints, ierr, iter=iter)
928 message(1) =
'Unable to write states wavefunctions.'
934 message(1) =
'Unable to write density.'
939 call lda_u_dump(scf%restart_dump, namespace, hm%lda_u, st, gr, ierr)
941 message(1) =
'Unable to write DFT+U information.'
946 select case (scf%mix_field)
947 case (option__mixfield__density)
948 call mix_dump(namespace, scf%restart_dump, scf%smix, gr, ierr)
950 message(1) =
'Unable to write mixing information.'
953 case (option__mixfield__potential)
954 call hm%ks_pot%dump(scf%restart_dump, gr, ierr)
956 message(1) =
'Unable to write Vhxc.'
960 call mix_dump(namespace, scf%restart_dump, scf%smix, gr, ierr)
962 message(1) =
'Unable to write mixing information.'
981 character(len=50) :: str
982 real(real64) :: dipole(1:space%dim)
988 write(str,
'(a,i5)')
'SCF CYCLE ITER #' ,iter
992 ' rel_ev = ', scf%evsum_diff/(abs(scf%evsum_out)+1e-20)
993 write(
message(2),
'(a,es15.2,2(a,es9.2))') &
994 ' ediff = ', scf%energy_diff,
' abs_dens = ', scf%abs_dens_diff, &
995 ' rel_dens = ', scf%abs_dens_diff/st%qtot
998 write(
message(1),
'(a,i6)')
'Matrix vector products: ', scf%eigens%matvec
999 write(
message(2),
'(a,i6)')
'Converged eigenvectors: ', sum(scf%eigens%converged(1:st%nik))
1003 if (
allocated(hm%vberry))
then
1005 call write_dipole(st, hm, space, dipole, namespace=namespace)
1018 write(
message(2),
'(a,i5,a,f14.2)')
'Elapsed time for SCF step ', iter,
':', etime
1028 write(
message(1),
'(a,i4,a,es15.8, a,es9.2, a, f7.1, a)') &
1031 ' : abs_dens', scf%abs_dens_diff, &
1032 ' : etime ', etime,
's'
1042 character(len=*),
intent(in) :: dir
1043 character(len=*),
intent(in) :: fname
1049 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write', position=
'append')
1051 call iterator%start(scf%criterion_list)
1052 do while (iterator%has_next())
1053 crit => iterator%get_next()
1058 write(iunit,
'(2es13.5)', advance =
'no') crit%val_abs, crit%val_rel
1061 write(iunit,
'(es13.5)', advance =
'no') crit%val_rel
1069 write(iunit,
'(es13.5)', advance =
'no') ks%oep%norm2ss
1072 write(iunit,
'(a)')
''
1079 logical function scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, &
1080 iters_done)
result(completed)
1081 type(
scf_t),
intent(inout) :: scf
1084 type(
grid_t),
intent(inout) :: gr
1085 type(
ions_t),
intent(inout) :: ions
1087 type(
v_ks_t),
intent(inout) :: ks
1089 integer,
intent(in) :: iter
1090 type(
output_t),
optional,
intent(in) :: outp
1091 integer,
optional,
intent(out) :: iters_done
1093 character(len=MAX_PATH_LEN) :: dirname
1094 integer(int64) :: what_i
1101 if(
present(iters_done)) iters_done = iter
1103 write(
message(1),
'(a, i4, a)')
'Info: SCF converged in ', iter,
' iterations'
1111 if (
present(outp))
then
1112 if (any(outp%what) .and. outp%duringscf)
then
1113 do what_i = lbound(outp%what, 1), ubound(outp%what, 1)
1114 if (outp%what_now(what_i, iter))
then
1115 write(dirname,
'(a,a,i4.4)') trim(outp%iter_dir),
"scf.", iter
1116 call output_all(outp, namespace, space, dirname, gr, ions, iter, st, hm, ks)
1117 call output_modelmb(outp, namespace, space, dirname, gr, ions, iter, st)
1125 call lalg_copy(gr%np, st%d%nspin, st%rho, scf%rhoin)
1128 if (scf%mix_field /= option__mixfield__none)
then
1129 if (scf%smix%ns_restart > 0)
then
1130 if (mod(iter, scf%smix%ns_restart) == 0)
then
1131 message(1) =
"Info: restarting mixing."
1138 select case(scf%mix_field)
1139 case(option__mixfield__potential)
1142 case (option__mixfield__density)
1147 if (scf%mix_field /= option__mixfield__states)
then
1158 subroutine scf_finish(scf, namespace, space, gr, ions, ext_partners, st, ks, hm, iter, outp)
1159 type(
scf_t),
intent(inout) :: scf
1162 type(
grid_t),
intent(inout) :: gr
1163 type(
ions_t),
intent(inout) :: ions
1166 type(
v_ks_t),
intent(inout) :: ks
1168 integer,
intent(in) :: iter
1169 type(
output_t),
optional,
intent(in) :: outp
1180 if ((scf%max_iter > 0 .and. scf%mix_field == option__mixfield__potential) .or. scf%calc_current)
then
1181 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, &
1182 calc_current=scf%calc_current)
1185 select case(scf%mix_field)
1186 case(option__mixfield__states)
1188 do iqn = st%d%kpt%start, st%d%kpt%end
1189 do ib = st%group%block_start, st%group%block_end
1190 call scf%psioutb(ib, iqn)%end()
1195 deallocate(scf%psioutb)
1198 safe_deallocate_a(scf%rhoout)
1199 safe_deallocate_a(scf%rhoin)
1201 if (scf%max_iter > 0 .and. any(scf%eigens%converged < st%nst))
then
1202 write(
message(1),
'(a)')
'Some of the states are not fully converged!'
1203 if (all(scf%eigens%converged >= st%nst_conv))
then
1204 write(
message(2),
'(a)')
'But all requested states to converge are converged.'
1211 if (.not.scf%finish)
then
1212 write(
message(1),
'(a,i4,a)')
'SCF *not* converged after ', iter - 1,
' iterations.'
1216 write(
message(1),
'(a,i10)')
'Info: Number of matrix-vector products: ', scf%matvec
1219 if (scf%calc_force)
then
1220 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old=scf%vhxc_old)
1223 if (scf%calc_stress)
call stress_calculate(namespace, gr, hm, st, ions, ks, ext_partners)
1226 if (scf%mix_field == option__mixfield__potential)
then
1231 if(
present(outp))
then
1238 if (space%is_periodic() .and. st%nik > st%d%nspin)
then
1241 ions, gr, hm%kpoints, hm%phase, vec_pot = hm%hm_base%uniform_vector_potential, &
1242 vec_pot_var = hm%hm_base%vector_potential)
1246 if (ks%vdw%vdw_correction == option__vdwcorrection__vdw_ts)
then
1250 safe_deallocate_a(scf%vhxc_old)
1258 character(len=*),
intent(in) :: dir, fname
1261 real(real64) :: dipole(1:space%dim)
1262 real(real64) :: ex_virial
1268 iunit =
io_open(trim(dir) //
"/" // trim(fname), namespace, action=
'write')
1274 if (space%is_periodic())
then
1275 call hm%kpoints%write_info(iunit=iunit)
1283 write(iunit,
'(a, i4, a)')
'SCF converged in ', iter,
' iterations'
1285 write(iunit,
'(a)')
'SCF *not* converged!'
1287 write(iunit,
'(1x)')
1289 if(any(scf%eigens%converged < st%nst))
then
1290 write(iunit,
'(a)')
'Some of the states are not fully converged!'
1291 if (all(scf%eigens%converged >= st%nst_conv))
then
1292 write(iunit,
'(a)')
'But all requested states to converge are converged.'
1297 write(iunit,
'(1x)')
1299 if (space%is_periodic())
then
1301 write(iunit,
'(1x)')
1311 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1314 if (
mpi_world%is_root())
write(iunit,
'(1x)')
1317 if(st%d%ispin ==
spinors .and. space%dim == 3 .and. &
1320 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1326 if(
mpi_world%is_root())
write(iunit,
'(1x)')
1329 if(scf%calc_dipole)
then
1336 hm%xc%functional(
func_c,1)%family == xc_family_none .and. st%d%ispin /=
spinors)
then
1343 write(iunit,
'(1x)')
1348 if(scf%max_iter > 0)
then
1349 write(iunit,
'(a)')
'Convergence:'
1350 call iterator%start(scf%criterion_list)
1351 do while (iterator%has_next())
1352 crit => iterator%get_next()
1353 call crit%write_info(iunit)
1362 write(iunit,
'(a)')
'Photon observables:'
1363 write(iunit,
'(6x, a, es15.8,a,es15.8,a)')
'Photon number = ', ks%oep_photon%pt%number(1)
1364 write(iunit,
'(6x, a, es15.8,a,es15.8,a)')
'Photon ex. = ', ks%oep_photon%pt%ex
1371 if (scf%calc_stress)
then
1372 call output_stress(iunit, space%periodic_dim, st%stress_tensors, all_terms=.false.)
1373 call output_pressure(iunit, space%periodic_dim, st%stress_tensors%total)
1378 if(scf%calc_partial_charges)
then
1413 real(real64) :: mem_tmp
1417 if(
conf%report_memory)
then
1419 call mpi_world%allreduce(mem, mem_tmp, 1, mpi_double_precision, mpi_sum)
1421 write(
message(1),
'(a,f14.2)')
'Memory usage [Mbytes] :', mem
1431 type(
scf_t),
intent(inout) :: scf
1437 select type (criterion)
1439 scf%energy_in = hm%energy%total
1454 type(
scf_t),
intent(inout) :: scf
1457 type(
grid_t),
intent(in) :: gr
1458 real(real64),
intent(in) :: rhoout(:,:), rhoin(:,:)
1462 real(real64),
allocatable :: tmp(:)
1466 select type (criterion)
1468 scf%energy_diff = abs(hm%energy%total - scf%energy_in)
1471 scf%abs_dens_diff =
m_zero
1472 safe_allocate(tmp(1:gr%np))
1473 do is = 1, st%d%nspin
1474 tmp(:) = abs(rhoin(1:gr%np, is) - rhoout(1:gr%np, is))
1475 scf%abs_dens_diff = scf%abs_dens_diff +
dmf_integrate(gr, tmp)
1477 safe_deallocate_a(tmp)
1481 scf%evsum_diff = abs(scf%evsum_out - scf%evsum_in)
1482 scf%evsum_in = scf%evsum_out
1492 subroutine write_dipole(st, hm, space, dipole, iunit, namespace)
1496 real(real64),
intent(in) :: dipole(:)
1497 integer,
optional,
intent(in) :: iunit
1498 type(
namespace_t),
optional,
intent(in) :: namespace
1503 call output_dipole(dipole, space%dim, iunit=iunit, namespace=namespace)
1505 if (space%is_periodic())
then
1506 message(1) =
"Defined only up to quantum of polarization (e * lattice vector)."
1507 message(2) =
"Single-point Berry's phase method only accurate for large supercells."
1510 if (hm%kpoints%full%npoints > 1)
then
1512 "WARNING: Single-point Berry's phase method for dipole should not be used when there is more than one k-point."
1513 message(2) =
"Instead, finite differences on k-points (not yet implemented) are needed."
1518 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.
Define which routines can be seen from the outside.
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
subroutine, public lda_u_update_occ_matrices(this, namespace, mesh, st, phase, energy)
integer, parameter, public dft_u_acbn0
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, mesh, ierr)
subroutine, public mix_init(smix, namespace, space, der, d1, d2, def_, func_type_, prefix_)
Initialise mix_t instance.
subroutine, public mix_load(namespace, restart, smix, mesh, ierr)
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_finish(scf, namespace, space, gr, ions, ext_partners, st, ks, hm, iter, outp)
subroutine, public scf_load(scf, namespace, space, gr, ions, ext_partners, st, ks, hm, restart_load)
Loading of restarting data of the SCF cycle.
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)
subroutine, public scf_start(scf, namespace, gr, ions, st, ks, hm, outp, verbosity)
Preparation of the SCF cycle.
integer, parameter, public verb_full
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_iter(scf, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, iter, outp, restart_dump)
logical function, public scf_iter_finish(scf, namespace, space, gr, ions, st, ks, hm, iter, outp, iters_done)
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, st, mesh, ierr)
subroutine, public states_elec_dump_rho(restart, 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