56 USE dbcsr_api,
ONLY: dbcsr_copy,&
57 dbcsr_deallocate_matrix,&
110 #include "./base/base_uses.f90"
118 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_loc_methods'
159 zij_fm_set, para_env, cell, weights, ispin, print_loc_section, &
160 restricted, nextra, nmo, vectors_2, guess_mos)
162 INTEGER,
INTENT(IN) :: method
163 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
164 TYPE(cp_fm_type),
INTENT(IN) :: vectors
165 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: op_sm_set
166 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij_fm_set
167 TYPE(mp_para_env_type),
POINTER :: para_env
168 TYPE(cell_type),
POINTER :: cell
169 REAL(
dp),
DIMENSION(:) :: weights
170 INTEGER,
INTENT(IN) :: ispin
171 TYPE(section_vals_type),
POINTER :: print_loc_section
172 INTEGER :: restricted
173 INTEGER,
INTENT(IN),
OPTIONAL :: nextra, nmo
174 TYPE(cp_fm_type),
INTENT(IN),
OPTIONAL :: vectors_2, guess_mos
176 CHARACTER(len=*),
PARAMETER :: routinen =
'optimize_loc_berry'
178 INTEGER :: handle, max_iter, nao, nmoloc, out_each, &
180 LOGICAL :: converged, crazy_use_diag, &
181 do_jacobi_refinement, my_do_mixed
182 REAL(
dp) :: crazy_scale, eps_localization, &
183 max_crazy_angle, start_time, &
185 TYPE(cp_logger_type),
POINTER :: logger
187 CALL timeset(routinen, handle)
190 extension=
".locInfo")
195 CALL zij_matrix(vectors, op_sm_set, zij_fm_set)
197 max_iter = qs_loc_env%localized_wfn_control%max_iter
198 max_crazy_angle = qs_loc_env%localized_wfn_control%max_crazy_angle
199 crazy_use_diag = qs_loc_env%localized_wfn_control%crazy_use_diag
200 crazy_scale = qs_loc_env%localized_wfn_control%crazy_scale
201 eps_localization = qs_loc_env%localized_wfn_control%eps_localization
202 out_each = qs_loc_env%localized_wfn_control%out_each
203 target_time = qs_loc_env%target_time
204 start_time = qs_loc_env%start_time
205 do_jacobi_refinement = qs_loc_env%localized_wfn_control%jacobi_refinement
206 my_do_mixed = qs_loc_env%localized_wfn_control%do_mixed
209 ispin, print_loc_section, only_initial_out=.true.)
213 CALL jacobi_rotations(weights, zij_fm_set, vectors, para_env, max_iter=max_iter, &
214 eps_localization=eps_localization, sweeps=sweeps, &
215 out_each=out_each, target_time=target_time, start_time=start_time, &
216 restricted=restricted)
218 IF (my_do_mixed)
THEN
220 IF (
PRESENT(guess_mos))
THEN
222 eps_localization, sweeps, out_each, nextra, &
223 qs_loc_env%localized_wfn_control%do_cg_po, &
224 nmo=nmo, vectors_2=vectors_2, mos_guess=guess_mos)
227 eps_localization, sweeps, out_each, nextra, &
228 qs_loc_env%localized_wfn_control%do_cg_po, &
229 nmo=nmo, vectors_2=vectors_2)
233 eps_localization, sweeps, out_each, 0, &
234 qs_loc_env%localized_wfn_control%do_cg_po)
237 cpabort(
"GAPO works only with STATES MIXED")
243 CALL zij_matrix(vectors, op_sm_set, zij_fm_set)
244 IF (do_jacobi_refinement)
THEN
247 ispin, print_loc_section, only_initial_out=.true.)
248 CALL jacobi_rotations(weights, zij_fm_set, vectors, para_env, max_iter=max_iter, &
249 eps_localization=eps_localization, sweeps=sweeps, &
250 out_each=out_each, target_time=target_time, start_time=start_time, &
251 restricted=restricted)
254 CALL crazy_rotations(weights, zij_fm_set, vectors, max_iter=max_iter, max_crazy_angle=max_crazy_angle, &
255 crazy_scale=crazy_scale, crazy_use_diag=crazy_use_diag, &
256 eps_localization=eps_localization, iterations=sweeps, converged=converged)
258 IF (.NOT. converged)
THEN
259 IF (qs_loc_env%localized_wfn_control%jacobi_fallback)
THEN
260 IF (output_unit > 0)
WRITE (output_unit,
"(T4,A,I6,/,T4,A)") &
261 " Crazy Wannier localization not converged after ", sweeps, &
262 " iterations, switching to jacobi rotations"
263 CALL jacobi_rotations(weights, zij_fm_set, vectors, para_env, max_iter=max_iter, &
264 eps_localization=eps_localization, sweeps=sweeps, &
265 out_each=out_each, target_time=target_time, start_time=start_time, &
266 restricted=restricted)
268 IF (output_unit > 0)
WRITE (output_unit,
"(T4,A,I6,/,T4,A)") &
269 " Crazy Wannier localization not converged after ", sweeps, &
270 " iterations, and jacobi_fallback switched off "
274 CALL direct_mini(weights, zij_fm_set, vectors, max_iter=max_iter, &
275 eps_localization=eps_localization, iterations=sweeps)
277 IF (.NOT. cell%orthorhombic)
THEN
278 cpabort(
"Non-orthorhombic cell with the selected method NYI")
282 CALL zij_matrix(vectors, op_sm_set, zij_fm_set)
285 IF (output_unit > 0)
THEN
286 WRITE (output_unit,
'(T4,A,I6,A)')
" No MOS localization applied "
289 cpabort(
"Unknown localization method")
291 IF (output_unit > 0)
THEN
292 IF (sweeps <= max_iter)
WRITE (output_unit,
'(T4,A,I3,A,I6,A)')
" Localization for spin ", ispin, &
293 " converged in ", sweeps,
" iterations"
297 ispin, print_loc_section)
301 CALL timestop(handle)
319 ispin, print_loc_section)
320 TYPE(qs_environment_type),
POINTER :: qs_env
321 INTEGER,
INTENT(IN) :: method
322 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
323 TYPE(cp_fm_type),
INTENT(IN) :: vectors
324 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij_fm_set
325 INTEGER,
INTENT(IN) :: ispin
326 TYPE(section_vals_type),
POINTER :: print_loc_section
328 CHARACTER(len=*),
PARAMETER :: routinen =
'optimize_loc_pipek'
330 INTEGER :: handle, iatom, isgf, ldz, nao, natom, &
331 ncol, nmoloc, output_unit, sweeps
332 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: first_sgf, last_sgf, nsgf
333 TYPE(cp_fm_pool_p_type),
DIMENSION(:),
POINTER :: ao_ao_fm_pools
334 TYPE(cp_fm_type) :: opvec
335 TYPE(cp_fm_type),
POINTER :: ov_fm
336 TYPE(cp_logger_type),
POINTER :: logger
337 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s
338 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
339 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
341 CALL timeset(routinen, handle)
344 extension=
".locInfo")
346 NULLIFY (particle_set)
354 CALL get_qs_env(qs_env=qs_env, matrix_s=matrix_s, &
355 particle_set=particle_set, qs_kind_set=qs_kind_set)
357 natom =
SIZE(particle_set, 1)
358 ALLOCATE (first_sgf(natom))
359 ALLOCATE (last_sgf(natom))
360 ALLOCATE (nsgf(natom))
364 first_sgf=first_sgf, last_sgf=last_sgf, nsgf=nsgf)
367 CALL mpools_get(qs_env%mpools, ao_ao_fm_pools=ao_ao_fm_pools)
376 isgf = first_sgf(iatom)
380 CALL parallel_gemm(
'N',
'N', nao, nmoloc, nao, 1.0_dp, ov_fm, vectors, 0.0_dp, opvec, &
381 a_first_col=1, a_first_row=1, b_first_col=1, b_first_row=1)
383 CALL parallel_gemm(
'T',
'N', nmoloc, nmoloc, ncol, 0.5_dp, vectors, opvec, &
384 0.0_dp, zij_fm_set(iatom, 1), &
385 a_first_col=1, a_first_row=isgf, b_first_col=1, b_first_row=isgf)
387 CALL parallel_gemm(
'N',
'N', nao, nmoloc, ncol, 1.0_dp, ov_fm, vectors, 0.0_dp, opvec, &
388 a_first_col=isgf, a_first_row=1, b_first_col=1, b_first_row=isgf)
390 CALL parallel_gemm(
'T',
'N', nmoloc, nmoloc, nao, 0.5_dp, vectors, opvec, &
391 1.0_dp, zij_fm_set(iatom, 1), &
392 a_first_col=1, a_first_row=1, b_first_col=1, b_first_row=1)
401 cpabort(
"GAPO and Pipek not implemented.")
403 cpabort(
"Crazy and Pipek not implemented.")
405 cpabort(
"L1 norm and Pipek not implemented.")
407 cpabort(
"Direct and Pipek not implemented.")
409 IF (output_unit > 0)
WRITE (output_unit,
'(A,I6,A)')
" No MOS localization applied "
411 cpabort(
"Unknown localization method")
414 IF (output_unit > 0)
WRITE (output_unit,
'(A,I3,A,I6,A)')
" Localization for spin ", ispin, &
415 " converged in ", sweeps,
" iterations"
417 CALL centers_spreads_pipek(qs_loc_env, zij_fm_set, particle_set, ispin, &
420 DEALLOCATE (first_sgf, last_sgf, nsgf)
422 CALL cp_fm_release(opvec)
425 CALL timestop(handle)
441 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij_fm_set
442 TYPE(cp_fm_type),
INTENT(IN) :: vectors
445 INTEGER :: iatom, istate, jstate, natom, nstate
446 REAL(kind=
dp) :: aij, bij, ct, mii, mij, mjj, ratio, st, &
448 TYPE(cp_fm_type) :: rmat
454 tolerance = 1.0e10_dp
456 natom =
SIZE(zij_fm_set, 1)
458 DO WHILE (tolerance >= 1.0e-4_dp)
460 DO istate = 1, nstate
461 DO jstate = istate + 1, nstate
468 aij = aij + mij*(mii - mjj)
469 bij = bij + mij*mij - 0.25_dp*(mii - mjj)*(mii - mjj)
471 IF (abs(bij) > 1.e-10_dp)
THEN
473 theta = 0.25_dp*atan(ratio)
480 IF (theta >
pi*0.5_dp)
THEN
481 theta = theta -
pi*0.25_dp
482 ELSE IF (theta < -
pi*0.5_dp)
THEN
483 theta = theta +
pi*0.25_dp
497 CALL check_tolerance_real(zij_fm_set, tolerance)
501 CALL cp_fm_release(rmat)
513 SUBROUTINE check_tolerance_real(zij_fm_set, tolerance)
515 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij_fm_set
516 REAL(
dp),
INTENT(OUT) :: tolerance
518 INTEGER :: iatom, istate, jstate, natom, &
519 ncol_local, nrow_global, nrow_local
520 INTEGER,
DIMENSION(:),
POINTER :: col_indices, row_indices
521 REAL(
dp) :: grad_ij, zii, zij, zjj
522 REAL(
dp),
DIMENSION(:, :),
POINTER :: diag
523 TYPE(cp_fm_type) :: force
525 CALL cp_fm_create(force, zij_fm_set(1, 1)%matrix_struct)
528 NULLIFY (diag, col_indices, row_indices)
529 natom =
SIZE(zij_fm_set, 1)
531 ncol_local=ncol_local, nrow_global=nrow_global, &
532 row_indices=row_indices, col_indices=col_indices)
533 ALLOCATE (diag(nrow_global, natom))
536 DO istate = 1, nrow_global
537 CALL cp_fm_get_element(zij_fm_set(iatom, 1), istate, istate, diag(istate, iatom))
541 DO istate = 1, nrow_local
542 DO jstate = 1, ncol_local
545 zii = diag(row_indices(istate), iatom)
546 zjj = diag(col_indices(jstate), iatom)
547 zij = zij_fm_set(iatom, 1)%local_data(istate, jstate)
548 grad_ij = grad_ij + 4.0_dp*zij*(zjj - zii)
550 force%local_data(istate, jstate) = grad_ij
557 CALL cp_fm_release(force)
559 END SUBROUTINE check_tolerance_real
575 print_loc_section, only_initial_out)
576 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
577 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij
578 INTEGER,
INTENT(IN) :: nmoloc
579 TYPE(cell_type),
POINTER :: cell
580 REAL(
dp),
DIMENSION(:) :: weights
581 INTEGER,
INTENT(IN) :: ispin
582 TYPE(section_vals_type),
POINTER :: print_loc_section
583 LOGICAL,
INTENT(IN),
OPTIONAL :: only_initial_out
585 CHARACTER(len=default_path_length) :: file_tmp, iter
586 COMPLEX(KIND=dp) :: z
587 INTEGER :: idir, istate, jdir, nstates, &
588 output_unit, unit_out_s
589 LOGICAL :: my_only_init
590 REAL(
dp) :: avg_spread_ii, spread_i, spread_ii, &
591 sum_spread_i, sum_spread_ii
592 REAL(
dp),
DIMENSION(3) :: c, c2, cpbc
593 REAL(
dp),
DIMENSION(:, :),
POINTER :: centers
594 REAL(kind=
dp) :: imagpart, realpart
595 TYPE(cp_logger_type),
POINTER :: logger
596 TYPE(section_vals_type),
POINTER :: print_key
598 NULLIFY (centers, logger, print_key)
600 my_only_init = .false.
601 IF (
PRESENT(only_initial_out)) my_only_init = only_initial_out
603 file_tmp = trim(qs_loc_env%tag_mo)//
"_spreads_s"//trim(adjustl(cp_to_string(ispin)))
605 extension=
".locInfo")
607 middle_name=file_tmp, extension=
".data")
609 IF (unit_out_s > 0 .AND. .NOT. my_only_init)
WRITE (unit_out_s,
'(i6,/,A)') nmoloc, trim(iter)
612 cpassert(nstates >= nmoloc)
614 centers => qs_loc_env%localized_wfn_control%centers_set(ispin)%array
615 cpassert(
ASSOCIATED(centers))
616 cpassert(
SIZE(centers, 2) == nmoloc)
617 sum_spread_i = 0.0_dp
618 sum_spread_ii = 0.0_dp
619 avg_spread_ii = 0.0_dp
620 DO istate = 1, nmoloc
625 DO jdir = 1,
SIZE(zij, 2)
628 z = cmplx(realpart, imagpart,
dp)
629 spread_i = spread_i - weights(jdir)* &
630 log(realpart*realpart + imagpart*imagpart)/
twopi/
twopi
631 spread_ii = spread_ii + weights(jdir)* &
632 (1.0_dp - (realpart*realpart + imagpart*imagpart))/
twopi/
twopi
635 c(idir) = c(idir) + &
636 (cell%hmat(idir, jdir)/
twopi)*aimag(log(z))
641 centers(1:3, istate) = cpbc(1:3)
642 centers(4, istate) = spread_i
643 centers(5, istate) = spread_ii
644 sum_spread_i = sum_spread_i + centers(4, istate)
645 sum_spread_ii = sum_spread_ii + centers(5, istate)
646 IF (unit_out_s > 0 .AND. .NOT. my_only_init)
WRITE (unit_out_s,
'(I6,2F16.8)') istate, centers(4:5, istate)
648 avg_spread_ii = sum_spread_ii/real(nmoloc, kind=
dp)
652 IF (.NOT. my_only_init)
CALL print_wannier_centers(qs_loc_env, print_key, centers, logger, ispin)
654 IF (output_unit > 0)
THEN
655 WRITE (output_unit,
'(T4, A, 2x, 2A26,/,T23, A28)')
" Spread Functional ",
"sum_in -w_i ln(|z_in|^2)", &
656 "sum_in w_i(1-|z_in|^2)",
"sum_in w_i(1-|z_in|^2)/n"
657 IF (my_only_init)
THEN
658 WRITE (output_unit,
'(T4,A,T38,2F20.10,/,T38,F20.10)')
" Initial Spread (Berry) : ", &
659 sum_spread_i, sum_spread_ii, avg_spread_ii
661 WRITE (output_unit,
'(T4,A,T38,2F20.10,/,T38,F20.10)')
" Total Spread (Berry) : ", &
662 sum_spread_i, sum_spread_ii, avg_spread_ii
666 IF (unit_out_s > 0 .AND. .NOT. my_only_init)
WRITE (unit_out_s,
'(A,2F16.10)')
" Total ", sum_spread_i, sum_spread_ii
684 SUBROUTINE centers_spreads_pipek(qs_loc_env, zij_fm_set, particle_set, ispin, &
687 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
688 TYPE(cp_fm_type),
DIMENSION(:, :),
INTENT(IN) :: zij_fm_set
689 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
690 INTEGER,
INTENT(IN) :: ispin
691 TYPE(section_vals_type),
POINTER :: print_loc_section
693 CHARACTER(len=default_path_length) :: file_tmp, iter
694 INTEGER :: iatom, istate, natom, nstate, unit_out_s
695 INTEGER,
DIMENSION(:),
POINTER :: atom_of_state
697 REAL(
dp),
ALLOCATABLE,
DIMENSION(:) :: qii, ziimax
698 REAL(
dp),
ALLOCATABLE,
DIMENSION(:, :) :: diag
699 REAL(
dp),
DIMENSION(:, :),
POINTER :: centers
700 TYPE(cp_logger_type),
POINTER :: logger
701 TYPE(section_vals_type),
POINTER :: print_key
703 NULLIFY (centers, logger, print_key)
707 natom =
SIZE(zij_fm_set, 1)
709 centers => qs_loc_env%localized_wfn_control%centers_set(ispin)%array
710 cpassert(
ASSOCIATED(centers))
711 cpassert(
SIZE(centers, 2) == nstate)
713 file_tmp = trim(qs_loc_env%tag_mo)//
"_spreads_s"//trim(adjustl(cp_to_string(ispin)))
715 middle_name=file_tmp, extension=
".data", log_filename=.false.)
717 IF (unit_out_s > 0)
WRITE (unit_out_s,
'(i6,/,A)') trim(iter)
719 ALLOCATE (atom_of_state(nstate))
722 ALLOCATE (diag(nstate, natom))
726 DO istate = 1, nstate
727 CALL cp_fm_get_element(zij_fm_set(iatom, 1), istate, istate, diag(istate, iatom))
731 ALLOCATE (qii(nstate), ziimax(nstate))
736 DO istate = 1, nstate
737 qii(istate) = qii(istate) + diag(istate, iatom)*diag(istate, iatom)
738 IF (abs(diag(istate, iatom)) > ziimax(istate))
THEN
739 ziimax(istate) = abs(diag(istate, iatom))
740 atom_of_state(istate) = iatom
745 DO istate = 1, nstate
746 iatom = atom_of_state(istate)
747 r(1:3) = particle_set(iatom)%r(1:3)
748 centers(1:3, istate) = r(1:3)
749 centers(4, istate) = 1.0_dp/qii(istate)
750 IF (unit_out_s > 0)
WRITE (unit_out_s,
'(I6,F16.8)') istate,
angstrom*centers(4, istate)
755 CALL print_wannier_centers(qs_loc_env, print_key, centers, logger, ispin)
759 DEALLOCATE (qii, ziimax, atom_of_state, diag)
761 END SUBROUTINE centers_spreads_pipek
783 print_key, root, ispin, idir, state0, file_position)
785 TYPE(qs_environment_type),
POINTER :: qs_env
786 TYPE(cp_fm_type),
INTENT(IN) :: mo_coeff
787 INTEGER,
INTENT(IN) :: nstates
788 INTEGER,
DIMENSION(:),
POINTER :: state_list
789 REAL(
dp),
DIMENSION(:, :),
POINTER :: centers
790 TYPE(section_vals_type),
POINTER :: print_key
791 CHARACTER(LEN=*) :: root
792 INTEGER,
INTENT(IN),
OPTIONAL :: ispin, idir
793 INTEGER,
OPTIONAL :: state0
794 CHARACTER(LEN=default_string_length),
OPTIONAL :: file_position
796 CHARACTER(len=*),
PARAMETER :: routinen =
'qs_print_cubes'
797 CHARACTER,
DIMENSION(3),
PARAMETER :: labels = (/
'x',
'y',
'z'/)
800 CHARACTER(LEN=default_path_length) :: file_tmp, filename, my_pos
801 CHARACTER(LEN=default_string_length) :: title
802 INTEGER :: handle, ia, ie, istate, ivector, &
803 my_ispin, my_state0, unit_out_c
804 LOGICAL :: add_idir, add_spin, mpi_io
805 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
806 TYPE(cell_type),
POINTER :: cell
807 TYPE(cp_logger_type),
POINTER :: logger
808 TYPE(dft_control_type),
POINTER :: dft_control
809 TYPE(particle_list_type),
POINTER :: particles
810 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
811 TYPE(pw_c1d_gs_type) :: wf_g
812 TYPE(pw_env_type),
POINTER :: pw_env
813 TYPE(pw_pool_type),
POINTER :: auxbas_pw_pool
814 TYPE(pw_r3d_rs_type) :: wf_r
815 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
816 TYPE(qs_subsys_type),
POINTER :: subsys
818 CALL timeset(routinen, handle)
819 NULLIFY (auxbas_pw_pool, pw_env, logger)
822 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env, subsys=subsys)
825 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
827 CALL auxbas_pw_pool%create_pw(wf_r)
828 CALL auxbas_pw_pool%create_pw(wf_g)
831 IF (
PRESENT(state0)) my_state0 = state0
835 IF (
PRESENT(ispin))
THEN
840 IF (
PRESENT(idir))
THEN
846 IF (
PRESENT(file_position))
THEN
847 my_pos = file_position
850 DO istate = 1, nstates
851 ivector = state_list(istate) - my_state0
852 CALL get_qs_env(qs_env=qs_env, atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set, cell=cell, &
853 dft_control=dft_control, particle_set=particle_set, pw_env=pw_env)
856 qs_kind_set, cell, dft_control, particle_set, pw_env)
859 ivector = state_list(istate)
862 filename = root(ia:ie)//
"_"//label//
"_w"//trim(adjustl(cp_to_string(ivector)))
864 filename = root(ia:ie)//
"_w"//trim(adjustl(cp_to_string(ivector)))
869 filename = file_tmp(ia:ie)//
"_s"//trim(adjustl(cp_to_string(ispin)))
875 extension=
".cube", file_position=my_pos, log_filename=.false., &
877 IF (
SIZE(centers, 1) == 6)
THEN
878 WRITE (title,
'(A7,I5.5,A2,I1.1,A1,6F10.4)')
"WFN ", ivector,
"_s", my_ispin,
" ", &
881 WRITE (title,
'(A7,I5.5,A2,I1.1,A1,3F10.4)')
"WFN ", ivector,
"_s", my_ispin,
" ", &
885 particles=particles, &
890 CALL auxbas_pw_pool%give_back_pw(wf_r)
891 CALL auxbas_pw_pool%give_back_pw(wf_g)
892 CALL timestop(handle)
903 SUBROUTINE print_wannier_centers(qs_loc_env, print_key, center, logger, ispin)
904 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
905 TYPE(section_vals_type),
POINTER :: print_key
906 REAL(kind=
dp),
INTENT(IN) :: center(:, :)
907 TYPE(cp_logger_type),
POINTER :: logger
908 INTEGER,
INTENT(IN) :: ispin
910 CHARACTER(default_string_length) :: iter, unit_str
911 CHARACTER(LEN=default_string_length) :: my_ext, my_form
912 INTEGER :: iunit, l, nstates
913 LOGICAL :: first_time, init_traj
914 REAL(kind=
dp) :: unit_conv
916 nstates =
SIZE(center, 2)
917 my_form =
"formatted"
926 IF (first_time .OR. (.NOT. qs_loc_env%first_time))
THEN
928 middle_name=trim(qs_loc_env%tag_mo)//
"_centers_s"//trim(adjustl(cp_to_string(ispin))), &
929 log_filename=.false., on_file=.true., is_new_file=init_traj)
937 CALL print_wannier_traj(qs_loc_env, print_key, center, iunit, init_traj, unit_conv)
941 WRITE (iunit,
'(i6,/,a)') nstates, trim(iter)
943 WRITE (iunit,
'(A,4F16.8)')
"X", unit_conv*center(1:4, l)
950 END SUBROUTINE print_wannier_centers
961 SUBROUTINE print_wannier_traj(qs_loc_env, print_key, center, iunit, init_traj, unit_conv)
962 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
963 TYPE(section_vals_type),
POINTER :: print_key
964 REAL(kind=
dp),
INTENT(IN) :: center(:, :)
965 INTEGER,
INTENT(IN) :: iunit
966 LOGICAL,
INTENT(IN) :: init_traj
967 REAL(kind=
dp),
INTENT(IN) :: unit_conv
969 CHARACTER(len=default_string_length) :: iter, remark1, remark2, title
970 INTEGER :: i, iskip, natom, ntot, outformat
971 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
972 TYPE(atomic_kind_type),
POINTER :: atomic_kind
973 TYPE(cp_logger_type),
POINTER :: logger
974 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
976 NULLIFY (particle_set, atomic_kind_set, atomic_kind, logger)
978 natom =
SIZE(qs_loc_env%particle_set)
979 ntot = natom +
SIZE(center, 2)
981 ALLOCATE (atomic_kind_set(1))
982 atomic_kind => atomic_kind_set(1)
984 name=
"X", element_symbol=
"X", mass=0.0_dp)
987 particle_set(i)%atomic_kind => qs_loc_env%particle_set(i)%atomic_kind
988 particle_set(i)%r =
pbc(qs_loc_env%particle_set(i)%r, qs_loc_env%cell)
991 DO i = natom + 1, ntot
992 particle_set(i)%atomic_kind => atomic_kind
993 particle_set(i)%r =
pbc(center(1:3, i - natom), qs_loc_env%cell)
999 SELECT CASE (outformat)
1010 WRITE (iunit)
"CORD", 0, -1, iskip, &
1011 0, 0, 0, 0, 0, 0, real(0, kind=
sp), 1, 0, 0, 0, 0, 0, 0, 0, 0, 24
1012 remark1 =
"REMARK FILETYPE CORD DCD GENERATED BY CP2K"
1013 remark2 =
"REMARK Support new DCD format with cell information"
1014 WRITE (iunit) 2, remark1, remark2
1020 WRITE (unit=title, fmt=
"(A)")
" Particles+Wannier centers. Iteration:"//trim(iter)
1025 unit_conv=unit_conv)
1029 END SUBROUTINE print_wannier_traj
1044 TYPE(qs_environment_type),
POINTER :: qs_env
1045 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
1046 TYPE(section_vals_type),
POINTER :: print_loc_section
1047 INTEGER,
INTENT(IN) :: ispin
1049 INTEGER,
PARAMETER :: norder = 2
1050 LOGICAL,
PARAMETER :: debug_this_subroutine = .false.
1052 CHARACTER(len=default_path_length) :: file_tmp
1053 INTEGER :: imoment, istate, ncol_global, nm, &
1054 nmoloc, nrow_global, output_unit, &
1056 REAL(
dp),
DIMENSION(:, :),
POINTER :: centers
1057 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: moment_set
1058 REAL(kind=
dp),
DIMENSION(3) :: rcc
1059 REAL(kind=
dp),
DIMENSION(6) :: cov
1060 TYPE(cp_fm_struct_type),
POINTER :: fm_struct
1061 TYPE(cp_fm_type) :: momv, mvector, omvector
1062 TYPE(cp_fm_type),
DIMENSION(:),
POINTER :: moloc_coeff
1063 TYPE(cp_fm_type),
POINTER :: mo_localized
1064 TYPE(cp_logger_type),
POINTER :: logger
1065 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s, moments
1066 TYPE(mp_para_env_type),
POINTER :: para_env
1071 extension=
".locInfo")
1073 IF (output_unit > 0)
THEN
1074 WRITE (output_unit,
'(/,T2,A)') &
1075 '!-----------------------------------------------------------------------------!'
1076 WRITE (output_unit,
'(T12,A)')
"Computing second moments of Wannier functions..."
1077 WRITE (output_unit,
'(T2,A)') &
1078 '!-----------------------------------------------------------------------------!'
1081 file_tmp =
"_r_loc_cov_matrix_"//trim(adjustl(cp_to_string(ispin)))
1083 middle_name=file_tmp, extension=
".data")
1085 CALL get_qs_loc_env(qs_loc_env=qs_loc_env, moloc_coeff=moloc_coeff)
1086 centers => qs_loc_env%localized_wfn_control%centers_set(ispin)%array
1094 para_env => qs_loc_env%para_env
1097 ALLOCATE (moment_set(nm, nmoloc))
1100 mo_localized => moloc_coeff(ispin)
1102 DO istate = 1, nmoloc
1103 rcc = centers(1:3, istate)
1106 ALLOCATE (moments(nm))
1108 ALLOCATE (moments(imoment)%matrix)
1109 CALL dbcsr_copy(moments(imoment)%matrix, matrix_s(ispin)%matrix,
'MOM MAT')
1110 CALL dbcsr_set(moments(imoment)%matrix, 0.0_dp)
1116 CALL cp_fm_get_info(mo_localized, ncol_global=ncol_global, nrow_global=nrow_global)
1118 para_env=mo_localized%matrix_struct%para_env, &
1119 context=mo_localized%matrix_struct%context)
1121 CALL cp_fm_create(omvector, fm_struct, name=
"omvector")
1126 CALL cp_fm_to_fm(mo_localized, mvector, 1, istate, 1)
1132 moment_set(imoment, istate) = sum(momv%local_data)
1135 CALL cp_fm_release(mvector)
1136 CALL cp_fm_release(omvector)
1137 CALL cp_fm_release(momv)
1140 CALL dbcsr_deallocate_matrix(moments(imoment)%matrix)
1142 DEALLOCATE (moments)
1145 CALL para_env%sum(moment_set)
1147 IF (output_unit_sm > 0)
THEN
1148 WRITE (unit=output_unit_sm, fmt=
'(A,T13,A,I1)')
"#", &
1149 "SECOND MOMENTS OF WANNIER CENTERS FOR SPIN ", ispin
1150 WRITE (unit=output_unit_sm, fmt=
'(A,T29,A2,5(14x,A2))')
"#",
"XX",
"XY",
"XZ",
"YY",
"YZ",
"ZZ"
1151 DO istate = 1, nmoloc
1153 cov(1) = moment_set(4, istate) - moment_set(1, istate)*moment_set(1, istate)
1154 cov(2) = moment_set(5, istate) - moment_set(1, istate)*moment_set(2, istate)
1155 cov(3) = moment_set(6, istate) - moment_set(1, istate)*moment_set(3, istate)
1156 cov(4) = moment_set(7, istate) - moment_set(2, istate)*moment_set(2, istate)
1157 cov(5) = moment_set(8, istate) - moment_set(2, istate)*moment_set(3, istate)
1158 cov(6) = moment_set(9, istate) - moment_set(3, istate)*moment_set(3, istate)
1159 WRITE (unit=output_unit_sm, fmt=
'(T4,A,I6,6(T20,6F16.8))')
"Center:", istate, cov(1:6)
1160 IF (debug_this_subroutine)
THEN
1161 WRITE (unit=output_unit_sm, fmt=
'(T4,A,I5,9(T20,9F12.6))')
"Center:", istate, moment_set(1:, istate)
1168 DEALLOCATE (moment_set)
1187 TYPE(qs_environment_type),
POINTER :: qs_env
1188 TYPE(qs_loc_env_type),
POINTER :: qs_loc_env
1189 TYPE(section_vals_type),
POINTER :: print_loc_section
1190 INTEGER,
INTENT(IN) :: ispin
1192 INTEGER,
PARAMETER :: taylor_order = 1
1193 LOGICAL,
PARAMETER :: debug_this_subroutine = .false.
1195 CHARACTER(len=default_path_length) :: file_tmp
1196 COMPLEX(KIND=dp) :: z
1197 INTEGER :: icov, imoment, istate, ncol_global, nm, &
1198 nmoloc, norder, nrow_global, &
1199 output_unit, output_unit_sm
1200 REAL(
dp),
DIMENSION(:, :),
POINTER :: centers
1201 REAL(kind=
dp) :: imagpart, lx, ly, lz, realpart
1202 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: moment_set
1203 REAL(kind=
dp),
DIMENSION(3) :: rcc
1204 REAL(kind=
dp),
DIMENSION(6) :: cov
1205 TYPE(cell_type),
POINTER :: cell
1206 TYPE(cp_fm_struct_type),
POINTER :: fm_struct
1207 TYPE(cp_fm_type) :: momv, mvector, omvector
1208 TYPE(cp_fm_type),
DIMENSION(:),
POINTER :: moloc_coeff
1209 TYPE(cp_fm_type),
POINTER :: mo_localized
1210 TYPE(cp_logger_type),
POINTER :: logger
1211 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s, moments
1212 TYPE(mp_para_env_type),
POINTER :: para_env
1219 extension=
".locInfo")
1221 IF (output_unit > 0)
THEN
1222 WRITE (output_unit,
'(/,T2,A)') &
1223 '!-----------------------------------------------------------------------------!'
1224 WRITE (output_unit,
'(T12,A)')
"Computing second moments of Wannier functions..."
1225 WRITE (output_unit,
'(T2,A)') &
1226 '!-----------------------------------------------------------------------------!'
1229 file_tmp =
"_r_berry_cov_matrix_"//trim(adjustl(cp_to_string(ispin)))
1231 middle_name=file_tmp, extension=
".data")
1233 norder = 2*(2*taylor_order + 1)
1234 nm = (6 + 11*norder + 6*norder**2 + norder**3)/6 - 1
1237 CALL get_qs_loc_env(qs_loc_env=qs_loc_env, moloc_coeff=moloc_coeff, cell=cell)
1238 lx = cell%hmat(1, 1)
1239 ly = cell%hmat(2, 2)
1240 lz = cell%hmat(3, 3)
1242 centers => qs_loc_env%localized_wfn_control%centers_set(ispin)%array
1244 para_env => qs_loc_env%para_env
1250 ALLOCATE (moment_set(nm, nmoloc))
1253 mo_localized => moloc_coeff(ispin)
1255 DO istate = 1, nmoloc
1256 rcc = centers(1:3, istate)
1259 ALLOCATE (moments(nm))
1261 ALLOCATE (moments(imoment)%matrix)
1262 CALL dbcsr_copy(moments(imoment)%matrix, matrix_s(ispin)%matrix,
'MOM MAT')
1263 CALL dbcsr_set(moments(imoment)%matrix, 0.0_dp)
1269 CALL cp_fm_get_info(mo_localized, ncol_global=ncol_global, nrow_global=nrow_global)
1271 para_env=mo_localized%matrix_struct%para_env, &
1272 context=mo_localized%matrix_struct%context)
1274 CALL cp_fm_create(omvector, fm_struct, name=
"omvector")
1279 CALL cp_fm_to_fm(mo_localized, mvector, 1, istate, 1)
1284 moment_set(imoment, istate) = sum(momv%local_data)
1287 CALL cp_fm_release(mvector)
1288 CALL cp_fm_release(omvector)
1289 CALL cp_fm_release(momv)
1292 CALL dbcsr_deallocate_matrix(moments(imoment)%matrix)
1294 DEALLOCATE (moments)
1297 CALL para_env%sum(moment_set)
1299 IF (output_unit_sm > 0)
THEN
1300 WRITE (unit=output_unit_sm, fmt=
'(A,T13,A,I1)')
"#", &
1301 "SECOND MOMENTS OF WANNIER CENTERS FOR SPIN ", ispin
1302 WRITE (unit=output_unit_sm, fmt=
'(A,T29,A2,5(14x,A2))')
"#",
"XX",
"XY",
"XZ",
"YY",
"YZ",
"ZZ"
1303 DO istate = 1, nmoloc
1308 z = cmplx(realpart, imagpart,
dp)
1313 realpart = 1.0 - 0.5*
twopi*
twopi/lx/lx/lx/lx &
1314 *moment_set(20, istate)
1315 imagpart =
twopi/lx/lx*moment_set(4, istate) &
1323 z = cmplx(realpart, imagpart,
dp)
1324 cov(1) = lx*lx/
twopi*aimag(log(z)) - lx*lx/
twopi/
twopi*moment_set(1, istate)*moment_set(1, istate)
1328 realpart = 1.0 - 0.5*
twopi*
twopi/lx/ly/lx/ly &
1329 *moment_set(23, istate)
1330 imagpart =
twopi/lx/ly*moment_set(5, istate) &
1338 z = cmplx(realpart, imagpart,
dp)
1339 cov(2) = lx*ly/
twopi*aimag(log(z)) - lx*ly/
twopi/
twopi*moment_set(1, istate)*moment_set(2, istate)
1343 realpart = 1.0 - 0.5*
twopi*
twopi/lx/lz/lx/lz &
1344 *moment_set(25, istate)
1345 imagpart =
twopi/lx/lz*moment_set(6, istate) &
1353 z = cmplx(realpart, imagpart,
dp)
1354 cov(3) = lx*lz/
twopi*aimag(log(z)) - lx*lz/
twopi/
twopi*moment_set(1, istate)*moment_set(3, istate)
1358 realpart = 1.0 - 0.5*
twopi*
twopi/ly/ly/ly/ly &
1359 *moment_set(30, istate)
1360 imagpart =
twopi/ly/ly*moment_set(7, istate) &
1368 z = cmplx(realpart, imagpart,
dp)
1369 cov(4) = ly*ly/
twopi*aimag(log(z)) - ly*ly/
twopi/
twopi*moment_set(2, istate)*moment_set(2, istate)
1373 realpart = 1.0 - 0.5*
twopi*
twopi/ly/lz/ly/lz &
1374 *moment_set(32, istate)
1375 imagpart =
twopi/ly/lz*moment_set(8, istate) &
1383 z = cmplx(realpart, imagpart,
dp)
1384 cov(5) = ly*lz/
twopi*aimag(log(z)) - ly*lz/
twopi/
twopi*moment_set(2, istate)*moment_set(3, istate)
1388 realpart = 1.0 - 0.5*
twopi*
twopi/lz/lz/lz/lz &
1389 *moment_set(34, istate)
1390 imagpart =
twopi/lz/lz*moment_set(9, istate) &
1398 z = cmplx(realpart, imagpart,
dp)
1399 cov(6) = lz*lz/
twopi*aimag(log(z)) - lz*lz/
twopi/
twopi*moment_set(3, istate)*moment_set(3, istate)
1403 WRITE (unit=output_unit_sm, fmt=
'(T4,A,I6,6(T20,6F16.8))')
"Center:", istate, cov(1:6)
1404 IF (debug_this_subroutine)
THEN
1405 WRITE (unit=output_unit_sm, fmt=
'(T4,A,I5,9(T20,9F12.6))')
"Center:", istate, moment_set(1:, istate)
1412 DEALLOCATE (moment_set)
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
Define the atomic kind types and their sub types.
subroutine, public set_atomic_kind(atomic_kind, element_symbol, name, mass, kind_number, natom, atom_list, fist_potential, shell, shell_active, damping)
Set the components of an atomic kind data set.
subroutine, public deallocate_atomic_kind_set(atomic_kind_set)
Destructor routine for a set of atomic kinds.
Handles all functions related to the CELL.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
basic linear algebra operations for full matrices
subroutine, public cp_fm_rot_rows(matrix, irow, jrow, cs, sn)
Applies a planar rotation defined by cs and sn to the i'th and j'th rows.
subroutine, public cp_fm_rot_cols(matrix, icol, jcol, cs, sn)
Applies a planar rotation defined by cs and sn to the i'th and j'th columnns.
subroutine, public cp_fm_schur_product(matrix_a, matrix_b, matrix_c)
computes the schur product of two matrices c_ij = a_ij * b_ij
pool for for elements that are retained and released
subroutine, public fm_pool_create_fm(pool, element, name)
returns an element, allocating it if none is in the pool
represent the structure of a full matrix
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
subroutine, public cp_fm_get_element(matrix, irow_global, icol_global, alpha, local)
returns an element of a fm this value is valid on every cpu using this call is expensive
subroutine, public cp_fm_maxabsval(matrix, a_max, ir_max, ic_max)
find the maximum absolute value of the matrix element maxval(abs(matrix))
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
character(len=default_string_length) function, public cp_iter_string(iter_info, print_key, for_file)
returns the iteration string, a string that is useful to create unique filenames (once you trim it)
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
A wrapper around pw_to_cube() which accepts particle_list_type.
subroutine, public cp_pw_to_cube(pw, unit_nr, title, particles, stride, zero_tails, silent, mpi_io)
...
real(kind=dp) function, public cp_unit_from_cp2k(value, unit_str, defaults, power)
converts from the internal cp2k units to the given unit
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
integer, parameter, public default_path_length
integer, parameter, public sp
Machine interface based on Fortran 2003 and POSIX.
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), parameter, public twopi
Interface to the message passing library MPI.
Output Utilities for MOTION_SECTION.
subroutine, public get_output_format(section, path, my_form, my_ext)
Info on the unit to be opened to dump MD informations.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public ncoset
basic linear algebra operations for full matrixes
represent a simple array based list of the given type
Define methods related to particle_type.
subroutine, public get_particle_set(particle_set, qs_kind_set, first_sgf, last_sgf, nsgf, nmao, basis)
Get the components of a particle set.
subroutine, public write_particle_coordinates(particle_set, iunit, output_format, content, title, cell, array, unit_conv, charge_occup, charge_beta, charge_extended, print_kind)
Should be able to write a few formats e.g. xmol, and some binary format (dcd) some format can be used...
Define the data structure for the particle information.
subroutine, public deallocate_particle_set(particle_set)
Deallocate a particle set.
subroutine, public allocate_particle_set(particle_set, nparticle)
Allocate a particle set.
Definition of physical constants:
real(kind=dp), parameter, public angstrom
container for various plainwaves related things
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
subroutine, public calculate_wavefunction(mo_vectors, ivector, rho, rho_gspace, atomic_kind_set, qs_kind_set, cell, dft_control, particle_set, pw_env, basis_type, external_vector)
maps a given wavefunction on the grid
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_RI_aux_kp, matrix_s, matrix_s_RI_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, WannierCentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
Driver for the localization that should be general for all the methods available and all the definiti...
subroutine, public optimize_loc_pipek(qs_env, method, qs_loc_env, vectors, zij_fm_set, ispin, print_loc_section)
...
subroutine, public centers_second_moments_berry(qs_env, qs_loc_env, print_loc_section, ispin)
Compute the second moments of the centers using a periodic quadrupole operator.
subroutine, public centers_second_moments_loc(qs_env, qs_loc_env, print_loc_section, ispin)
Compute the second moments of the centers using the local (non-periodic) pos operators.
subroutine, public optimize_loc_berry(method, qs_loc_env, vectors, op_sm_set, zij_fm_set, para_env, cell, weights, ispin, print_loc_section, restricted, nextra, nmo, vectors_2, guess_mos)
Calculate and optimize the spread functional as calculated from the selected mos and by the definitio...
subroutine, public centers_spreads_berry(qs_loc_env, zij, nmoloc, cell, weights, ispin, print_loc_section, only_initial_out)
...
subroutine, public qs_print_cubes(qs_env, mo_coeff, nstates, state_list, centers, print_key, root, ispin, idir, state0, file_position)
write the cube files for a set of selected states
subroutine, public jacobi_rotation_pipek(zij_fm_set, vectors, sweeps)
2by2 rotation for the pipek operator in this case the z_ij numbers are reals
New version of the module for the localization of the molecular orbitals This should be able to use d...
subroutine, public get_qs_loc_env(qs_loc_env, cell, local_molecules, localized_wfn_control, moloc_coeff, op_sm_set, op_fm_set, para_env, particle_set, weights, dim_op)
...
Localization methods such as 2x2 Jacobi rotations Steepest Decents Conjugate Gradient.
subroutine, public rotate_orbitals(rmat, vectors)
...
subroutine, public jacobi_rotations(weights, zij, vectors, para_env, max_iter, eps_localization, sweeps, out_each, target_time, start_time, restricted)
wrapper for the jacobi routines, should be removed if jacobi_rot_para can deal with serial para_envs.
subroutine, public zij_matrix(vectors, op_sm_set, zij_fm_set)
...
subroutine, public direct_mini(weights, zij, vectors, max_iter, eps_localization, iterations)
use the exponential parametrization as described in to perform a direct mini Gerd Berghold et al....
subroutine, public crazy_rotations(weights, zij, vectors, max_iter, max_crazy_angle, crazy_scale, crazy_use_diag, eps_localization, iterations, converged)
yet another crazy try, computes the angles needed to rotate the orbitals first and rotates them all a...
subroutine, public approx_l1_norm_sd(C, iterations, eps, converged, sweeps)
...
subroutine, public scdm_qrfact(vectors)
...
subroutine, public jacobi_cg_edf_ls(para_env, weights, zij, vectors, max_iter, eps_localization, iter, out_each, nextra, do_cg, nmo, vectors_2, mos_guess)
combine jacobi rotations (serial) and conjugate gradient with golden section line search for partiall...
wrapper for the pools of matrixes
subroutine, public mpools_get(mpools, ao_mo_fm_pools, ao_ao_fm_pools, mo_mo_fm_pools, ao_mosub_fm_pools, mosub_mosub_fm_pools, maxao_maxmo_fm_pool, maxao_maxao_fm_pool, maxmo_maxmo_fm_pool)
returns various attributes of the mpools (notably the pools contained in it)
Calculates the moment integrals <a|r^m|b> and <a|r x d/dr|b>
subroutine, public build_local_moment_matrix(qs_env, moments, nmoments, ref_point, ref_points, basis_type)
...
types that represent a quickstep subsys
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
Utilities for string manipulations.
elemental subroutine, public xstring(string, ia, ib)
...