80#include "./base/base_uses.f90"
86 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'almo_scf_qs'
115 name_new, size_keys, symmetry_new, &
116 spin_key, init_domains)
120 CHARACTER(len=*),
INTENT(IN) :: name_new
121 INTEGER,
DIMENSION(2),
INTENT(IN) :: size_keys
122 CHARACTER,
INTENT(IN) :: symmetry_new
123 INTEGER,
INTENT(IN) :: spin_key
124 LOGICAL,
INTENT(IN) :: init_domains
126 CHARACTER(len=*),
PARAMETER :: routinen =
'matrix_almo_create'
128 INTEGER :: dimen, handle, hold, iatom, iblock_col, &
129 iblock_row, imol, mynode, natoms, &
130 nblkrows_tot, nlength, nmols, row
131 INTEGER,
DIMENSION(:),
POINTER :: blk_distr, blk_sizes, block_sizes_new, col_blk_size, &
132 col_distr_new, col_sizes_new, distr_new_array, row_blk_size, row_distr_new, row_sizes_new
133 LOGICAL :: active, one_dim_is_mo, tr
134 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: new_block
149 CALL timeset(routinen, handle)
164 natoms = almo_scf_env%natoms
165 nmols = almo_scf_env%nmolecules
183 ALLOCATE (block_sizes_new(natoms), distr_new_array(natoms))
184 block_sizes_new(:) = blk_sizes(:)
185 distr_new_array(:) = blk_distr(:)
188 ALLOCATE (block_sizes_new(nmols), distr_new_array(nmols))
189 block_sizes_new(:) = 0
191 block_sizes_new(almo_scf_env%domain_index_of_atom(iatom)) = &
192 block_sizes_new(almo_scf_env%domain_index_of_atom(iatom)) + &
196 distr_new_array(imol) = &
197 blk_distr(almo_scf_env%first_atom_of_domain(imol))
200 cpabort(
"Illegal distribution")
213 ALLOCATE (block_sizes_new(nlength))
214 block_sizes_new(:) = 0
225 ALLOCATE (block_sizes_new(nlength))
227 block_sizes_new(:) = almo_scf_env%nocc_of_domain(:, spin_key)
230 WHERE (block_sizes_new == 0) block_sizes_new = 1
232 block_sizes_new(:) = almo_scf_env%nvirt_disc_of_domain(:, spin_key)
234 block_sizes_new(:) = almo_scf_env%nvirt_full_of_domain(:, spin_key)
236 block_sizes_new(:) = almo_scf_env%nvirt_of_domain(:, spin_key)
239 cpabort(
"Illegal distribution")
244 cpabort(
"Illegal dimension")
249 ALLOCATE (distr_new_array(nlength))
252 distr_new_array(:) = blk_distr(:)
256 distr_new_array(imol) = &
257 blk_distr(almo_scf_env%first_atom_of_domain(imol))
264 row_sizes_new => block_sizes_new
265 row_distr_new => distr_new_array
267 col_sizes_new => block_sizes_new
268 col_distr_new => distr_new_array
274 row_dist=row_distr_new, col_dist=col_distr_new, &
279 dist_new, symmetry_new, &
280 row_sizes_new, col_sizes_new, reuse_arrays=.true.)
285 IF (init_domains)
THEN
290 row_blk_size=row_blk_size, col_blk_size=col_blk_size)
350 DO row = 1, nblkrows_tot
356 IF (hold .EQ. mynode)
THEN
360 one_dim_is_mo = .false.
364 IF (one_dim_is_mo)
THEN
365 IF (almo_scf_env%nocc_of_domain(row, spin_key) == 0) active = .false.
368 one_dim_is_mo = .false.
372 IF (one_dim_is_mo)
THEN
373 IF (almo_scf_env%nvirt_of_domain(row, spin_key) == 0) active = .false.
376 one_dim_is_mo = .false.
380 IF (one_dim_is_mo)
THEN
381 IF (almo_scf_env%nvirt_disc_of_domain(row, spin_key) == 0) active = .false.
384 one_dim_is_mo = .false.
388 IF (one_dim_is_mo)
THEN
389 IF (almo_scf_env%nvirt_full_of_domain(row, spin_key) == 0) active = .false.
393 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
394 new_block(:, :) = 1.0_dp
396 DEALLOCATE (new_block)
407 CALL timestop(handle)
423 INTEGER :: mat_distr_aos
425 CHARACTER(len=*),
PARAMETER :: routinen =
'matrix_qs_to_almo'
430 CALL timeset(routinen, handle)
433 SELECT CASE (mat_distr_aos)
439 CALL dbcsr_create(matrix_qs_nosym, template=matrix_qs, matrix_type=dbcsr_type_no_symmetry)
457 CALL timestop(handle)
470 SUBROUTINE matrix_almo_to_qs(matrix_almo, matrix_qs, mat_distr_aos)
472 INTEGER,
INTENT(IN) :: mat_distr_aos
474 CHARACTER(len=*),
PARAMETER :: routinen =
'matrix_almo_to_qs'
479 CALL timeset(routinen, handle)
482 SELECT CASE (mat_distr_aos)
484 CALL dbcsr_copy(matrix_qs, matrix_almo, keep_sparsity=.true.)
486 CALL dbcsr_create(matrix_almo_redist, template=matrix_qs)
489 CALL dbcsr_copy(matrix_qs, matrix_almo_redist, keep_sparsity=.true.)
495 CALL timestop(handle)
497 END SUBROUTINE matrix_almo_to_qs
513 INTEGER :: mat_distr_aos
514 REAL(kind=
dp) :: eps_filter
516 CHARACTER(len=*),
PARAMETER :: routinen =
'init_almo_ks_matrix_via_qs'
518 INTEGER :: handle, ispin, nspin
519 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_qs_ks, matrix_qs_s
525 CALL timeset(routinen, handle)
531 dft_control=dft_control, &
532 matrix_s=matrix_qs_s, &
533 matrix_ks=matrix_qs_ks, &
537 nspin = dft_control%nspins
540 IF (.NOT.
ASSOCIATED(matrix_qs_ks))
THEN
543 ALLOCATE (matrix_qs_ks(ispin)%matrix)
545 template=matrix_qs_s(1)%matrix)
547 CALL dbcsr_set(matrix_qs_ks(ispin)%matrix, 0.0_dp)
549 CALL set_ks_env(ks_env, matrix_ks=matrix_qs_ks)
554 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
558 CALL timestop(handle)
575 CHARACTER(len=*),
PARAMETER :: routinen =
'construct_qs_mos'
577 INTEGER :: handle, ispin, ncol_fm, nrow_fm
584 CALL timeset(routinen, handle)
587 NULLIFY (mos, fm_struct_tmp, scf_env)
593 CALL get_qs_env(qs_env, dft_control=dft_control, mos=mos)
595 CALL dbcsr_get_info(almo_scf_env%matrix_t(1), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
598 DO ispin = 1, almo_scf_env%nspins
607 nelectron=almo_scf_env%nelectrons_total, &
608 n_el_f=real(almo_scf_env%nelectrons_total,
dp), &
610 flexible_electron_count=dft_control%relax_multiplicity)
613 context=almo_scf_env%blacs_env, &
614 para_env=almo_scf_env%para_env)
616 CALL cp_fm_create(mo_fm_copy, fm_struct_tmp, name=
"t_orthogonal_converted_to_fm")
620 CALL init_mo_set(mos(ispin), fm_ref=mo_fm_copy, name=
'fm_mo')
626 CALL timestop(handle)
642 INTEGER,
INTENT(IN) :: mat_distr_aos
644 CHARACTER(len=*),
PARAMETER :: routinen =
'almo_dm_to_qs_env'
646 INTEGER :: handle, ispin, nspins
650 CALL timeset(routinen, handle)
652 NULLIFY (rho, rho_ao)
653 nspins =
SIZE(matrix_p)
659 CALL matrix_almo_to_qs(matrix_p(ispin), &
660 rho_ao(ispin)%matrix, &
666 CALL timestop(handle)
684 SUBROUTINE almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, smear, kTS_sum)
687 REAL(kind=
dp) :: energy_total
688 INTEGER,
INTENT(IN) :: mat_distr_aos
689 LOGICAL,
INTENT(IN),
OPTIONAL :: smear
690 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: kts_sum
692 CHARACTER(len=*),
PARAMETER :: routinen =
'almo_dm_to_qs_ks'
696 REAL(kind=
dp) :: entropic_term
699 CALL timeset(routinen, handle)
701 IF (
PRESENT(smear))
THEN
707 IF (
PRESENT(kts_sum))
THEN
708 entropic_term = kts_sum
710 entropic_term = 0.0_dp
722 energy%total = energy%total - energy%kTS + entropic_term
725 energy_total = energy%total
727 CALL timestop(handle)
729 END SUBROUTINE almo_dm_to_qs_ks
748 mat_distr_aos, smear, kTS_sum)
751 TYPE(
dbcsr_type),
DIMENSION(:) :: matrix_p, matrix_ks
752 REAL(kind=
dp) :: energy_total, eps_filter
753 INTEGER,
INTENT(IN) :: mat_distr_aos
754 LOGICAL,
INTENT(IN),
OPTIONAL :: smear
755 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: kts_sum
757 CHARACTER(len=*),
PARAMETER :: routinen =
'almo_dm_to_almo_ks'
759 INTEGER :: handle, ispin, nspins
761 REAL(kind=
dp) :: entropic_term
762 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_qs_ks
764 CALL timeset(routinen, handle)
766 IF (
PRESENT(smear))
THEN
772 IF (
PRESENT(kts_sum))
THEN
773 entropic_term = kts_sum
775 entropic_term = 0.0_dp
779 CALL almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, &
781 kts_sum=entropic_term)
783 nspins =
SIZE(matrix_ks)
786 CALL get_qs_env(qs_env, matrix_ks=matrix_qs_ks)
788 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
792 CALL timestop(handle)
808 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: energy, energy_singles_corr
814 IF (
PRESENT(energy_singles_corr))
THEN
815 qs_energy%singles_corr = energy_singles_corr
820 IF (
PRESENT(energy))
THEN
841 CHARACTER(len=*),
PARAMETER :: routinen =
'almo_scf_construct_quencher'
844 INTEGER :: col, contact_atom_1, contact_atom_2, domain_col, domain_map_local_entries, &
845 domain_row, global_entries, global_list_length, grid1, groupid, handle, hold, iatom, &
846 iatom2, iblock_col, iblock_row, idomain, idomain2, ientry, igrid, ineig, ineighbor, &
847 inode, inode2, ipair, ispin, jatom, jatom2, jdomain2, local_list_length, &
848 max_domain_neighbors, max_neig, mynode, nblkcols_tot, nblkrows_tot, nblks, ndomains, &
849 neig_temp, nnode2, nnodes, row, unit_nr
850 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: current_number_neighbors, domain_entries_cpu, &
851 domain_map_global, domain_map_local, first_atom_of_molecule, global_list, &
852 last_atom_of_molecule, list_length_cpu, list_offset_cpu, local_list, offset_for_cpu
853 INTEGER,
ALLOCATABLE,
DIMENSION(:, :) :: domain_grid, domain_neighbor_list, &
854 domain_neighbor_list_excessive
855 INTEGER,
DIMENSION(:),
POINTER :: col_blk_size, row_blk_size
856 LOGICAL :: already_listed, block_active, &
857 delayed_increment, found, &
859 REAL(kind=
dp) :: contact1_radius, contact2_radius, &
861 r0, r1, s0, s1, trial_distance_squared
862 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: new_block
863 REAL(kind=
dp),
DIMENSION(3) :: rab
864 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: p_old_block
873 DIMENSION(:),
POINTER :: nl_iterator, nl_iterator2
878 CALL timeset(routinen, handle)
882 IF (logger%para_env%is_source())
THEN
888 ndomains = almo_scf_env%ndomains
891 particle_set=particle_set, &
892 molecule_set=molecule_set, &
900 ALLOCATE (first_atom_of_molecule(almo_scf_env%nmolecules))
901 ALLOCATE (last_atom_of_molecule(almo_scf_env%nmolecules))
903 mol_to_first_atom=first_atom_of_molecule, &
904 mol_to_last_atom=last_atom_of_molecule)
909 template=almo_scf_env%matrix_s(1), &
910 matrix_type=dbcsr_type_no_symmetry)
913 IF (sym .EQ. dbcsr_type_no_symmetry)
THEN
914 CALL dbcsr_copy(matrix_s_sym, almo_scf_env%matrix_s(1))
920 ALLOCATE (almo_scf_env%quench_t(almo_scf_env%nspins))
921 ALLOCATE (almo_scf_env%domain_map(almo_scf_env%nspins))
928 matrix_qs=matrix_s(1)%matrix, &
929 almo_scf_env=almo_scf_env, &
930 name_new=
"T_QUENCHER", &
932 symmetry_new=dbcsr_type_no_symmetry, &
934 init_domains=.false.)
938 CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), distribution=dist, &
939 nblkrows_total=nblkrows_tot, nblkcols_total=nblkcols_tot)
941 CALL group%set_handle(groupid)
945 local_list_length = 0
951 iatom=iatom2, jatom=jatom2, inode=inode2, nnode=nnode2)
953 IF (inode2 == 1)
THEN
954 local_list_length = local_list_length + nnode2
960 ALLOCATE (local_list(2*local_list_length))
962 local_list_length = 0
966 iatom=iatom2, jatom=jatom2)
967 local_list(2*local_list_length + 1) = iatom2
968 local_list(2*local_list_length + 2) = jatom2
969 local_list_length = local_list_length + 1
974 ALLOCATE (list_length_cpu(nnodes), list_offset_cpu(nnodes))
975 CALL group%allgather(2*local_list_length, list_length_cpu)
978 list_offset_cpu(1) = 0
980 list_offset_cpu(inode) = list_offset_cpu(inode - 1) + &
981 list_length_cpu(inode - 1)
983 global_list_length = list_offset_cpu(nnodes) + list_length_cpu(nnodes)
986 ALLOCATE (global_list(global_list_length))
987 CALL group%allgatherv(local_list, global_list, &
988 list_length_cpu, list_offset_cpu)
989 DEALLOCATE (list_length_cpu, list_offset_cpu)
990 DEALLOCATE (local_list)
993 ALLOCATE (current_number_neighbors(almo_scf_env%ndomains))
994 current_number_neighbors(:) = 0
995 global_list_length = global_list_length/2
996 DO ipair = 1, global_list_length
997 iatom2 = global_list(2*(ipair - 1) + 1)
998 jatom2 = global_list(2*(ipair - 1) + 2)
999 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
1000 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
1002 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
1004 IF (idomain2 .NE. jdomain2)
THEN
1005 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
1008 max_domain_neighbors = maxval(current_number_neighbors)
1011 ALLOCATE (domain_neighbor_list_excessive(ndomains, max_domain_neighbors))
1012 current_number_neighbors(:) = 1
1013 DO ipair = 1, ndomains
1014 domain_neighbor_list_excessive(ipair, 1) = ipair
1016 DO ipair = 1, global_list_length
1017 iatom2 = global_list(2*(ipair - 1) + 1)
1018 jatom2 = global_list(2*(ipair - 1) + 2)
1019 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
1020 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
1021 already_listed = .false.
1022 DO ineighbor = 1, current_number_neighbors(idomain2)
1023 IF (domain_neighbor_list_excessive(idomain2, ineighbor) .EQ. jdomain2)
THEN
1024 already_listed = .true.
1028 IF (.NOT. already_listed)
THEN
1030 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
1031 domain_neighbor_list_excessive(idomain2, current_number_neighbors(idomain2)) = jdomain2
1033 IF (idomain2 .NE. jdomain2)
THEN
1034 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
1035 domain_neighbor_list_excessive(jdomain2, current_number_neighbors(jdomain2)) = idomain2
1039 DEALLOCATE (global_list)
1041 max_domain_neighbors = maxval(current_number_neighbors)
1042 ALLOCATE (domain_neighbor_list(ndomains, max_domain_neighbors))
1043 domain_neighbor_list(:, :) = 0
1044 domain_neighbor_list(:, :) = domain_neighbor_list_excessive(:, 1:max_domain_neighbors)
1045 DEALLOCATE (domain_neighbor_list_excessive)
1047 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1048 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(max_domain_neighbors*ndomains, 2))
1049 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1050 almo_scf_env%domain_map(ispin)%index1(:) = 0
1051 domain_map_local_entries = 0
1057 row_blk_size=row_blk_size, col_blk_size=col_blk_size)
1059 DO row = 1, nblkrows_tot
1060 DO col = 1, current_number_neighbors(row)
1063 iblock_col = domain_neighbor_list(row, col)
1065 iblock_row, iblock_col, hold)
1067 IF (hold .EQ. mynode)
THEN
1071 domain_row = almo_scf_env%domain_index_of_ao_block(iblock_row)
1073 domain_col = almo_scf_env%domain_index_of_mo_block(iblock_col)
1075 SELECT CASE (almo_scf_env%constraint_type)
1078 block_active = .false.
1086 IF (domain_row == domain_col)
THEN
1087 block_active = .true.
1093 cpabort(
"Illegal: atomic domains and molecular groups")
1103 cpabort(
"Illegal: molecular domains and atomic groups")
1108 IF (domain_row == domain_col)
THEN
1109 block_active = .true.
1116 IF (block_active)
THEN
1118 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1119 new_block(:, :) = 1.0_dp
1120 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1121 DEALLOCATE (new_block)
1123 IF (domain_map_local_entries .GE. max_domain_neighbors*almo_scf_env%ndomains)
THEN
1124 cpabort(
"weird... max_domain_neighbors is exceeded")
1126 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1127 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1128 domain_map_local_entries = domain_map_local_entries + 1
1143 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1147 overlap = maxval(abs(p_old_block))
1156 cpabort(
"atomic domains and molecular groups - NYI")
1167 cpabort(
"molecular domains and atomic groups - NYI")
1173 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1177 overlap = maxval(abs(p_old_block))
1186 s0 = -log10(abs(almo_scf_env%quencher_s0))
1187 s1 = -log10(abs(almo_scf_env%quencher_s1))
1188 IF (overlap .EQ. 0.0_dp)
THEN
1189 overlap = -log10(abs(almo_scf_env%eps_filter)) + 100.0_dp
1191 overlap = -log10(overlap)
1193 IF (s0 .LT. 0.0_dp)
THEN
1194 cpabort(
"S0 is less than zero")
1196 IF (s1 .LE. 0.0_dp)
THEN
1197 cpabort(
"S1 is less than or equal to zero")
1199 IF (s0 .GE. s1)
THEN
1200 cpabort(
"S0 is greater than or equal to S1")
1204 IF (overlap .LT. s1)
THEN
1205 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1206 IF (overlap .LE. s0)
THEN
1207 new_block(:, :) = 1.0_dp
1211 new_block(:, :) = 1.0_dp/(1.0_dp + exp(-(s0 - s1)/(s0 - overlap) - (s0 - s1)/(overlap - s1)))
1216 IF (abs(new_block(1, 1)) .GT. abs(almo_scf_env%eps_filter))
THEN
1217 IF (domain_map_local_entries .GE. max_domain_neighbors*almo_scf_env%ndomains)
THEN
1218 cpabort(
"weird... max_domain_neighbors is exceeded")
1220 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1221 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1222 domain_map_local_entries = domain_map_local_entries + 1
1225 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1226 DEALLOCATE (new_block)
1243 IF (domain_row == domain_col)
THEN
1245 contact_atom_1 = first_atom_of_molecule(domain_row)
1246 contact_atom_2 = first_atom_of_molecule(domain_col)
1251 DO iatom = first_atom_of_molecule(domain_row), last_atom_of_molecule(domain_row)
1252 DO jatom = first_atom_of_molecule(domain_col), last_atom_of_molecule(domain_col)
1253 rab(:) =
pbc(particle_set(iatom)%r(:), particle_set(jatom)%r(:), cell)
1254 trial_distance_squared = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
1257 contact_atom_1 = iatom
1258 contact_atom_2 = jatom
1262 cpassert(contact_atom_1 .GT. 0)
1270 cpabort(
"atomic domains and molecular groups - NYI")
1281 cpabort(
"molecular domains and atomic groups - NYI")
1287 rab(:) =
pbc(particle_set(domain_row)%r(:), particle_set(domain_col)%r(:), cell)
1288 distance = sqrt(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
1289 contact_atom_1 = domain_row
1290 contact_atom_2 = domain_col
1298 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1299 rcov=contact1_radius)
1300 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1301 rcov=contact2_radius)
1303 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1304 rvdw=contact1_radius)
1305 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1306 rvdw=contact2_radius)
1308 cpabort(
"Illegal quencher_radius_type")
1318 r0 = almo_scf_env%quencher_r0_factor*(contact1_radius + contact2_radius)
1320 r1 = almo_scf_env%quencher_r1_factor*(contact1_radius + contact2_radius)
1323 IF (r0 .LT. 0.0_dp)
THEN
1324 cpabort(
"R0 is less than zero")
1326 IF (r1 .LE. 0.0_dp)
THEN
1327 cpabort(
"R1 is less than or equal to zero")
1329 IF (r0 .GT. r1)
THEN
1330 cpabort(
"R0 is greater than or equal to R1")
1334 IF (distance .LT. r1)
THEN
1335 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1336 IF (distance .LE. r0)
THEN
1337 new_block(:, :) = 1.0_dp
1344 new_block(:, :) = 1.0_dp/(1.0_dp + exp((r1 - r0)/(r0 - distance) + (r1 - r0)/(r1 - distance)))
1350 IF (abs(new_block(1, 1)) .GT. abs(almo_scf_env%eps_filter))
THEN
1351 IF (domain_map_local_entries .GE. max_domain_neighbors*almo_scf_env%ndomains)
THEN
1352 cpabort(
"weird... max_domain_neighbors is exceeded")
1354 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1355 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1356 domain_map_local_entries = domain_map_local_entries + 1
1359 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1360 DEALLOCATE (new_block)
1364 cpabort(
"Illegal constraint type")
1372 DEALLOCATE (domain_neighbor_list)
1373 DEALLOCATE (current_number_neighbors)
1379 almo_scf_env%eps_filter)
1383 IF (nblks .NE. domain_map_local_entries)
THEN
1384 cpabort(
"number of blocks is wrong")
1388 ALLOCATE (domain_entries_cpu(nnodes), offset_for_cpu(nnodes))
1389 CALL group%allgather(2*domain_map_local_entries, domain_entries_cpu)
1392 offset_for_cpu(1) = 0
1393 DO inode = 2, nnodes
1394 offset_for_cpu(inode) = offset_for_cpu(inode - 1) + &
1395 domain_entries_cpu(inode - 1)
1397 global_entries = offset_for_cpu(nnodes) + domain_entries_cpu(nnodes)
1400 ALLOCATE (domain_map_global(global_entries))
1401 ALLOCATE (domain_map_local(2*domain_map_local_entries))
1402 DO ientry = 1, domain_map_local_entries
1403 domain_map_local(2*(ientry - 1) + 1) = almo_scf_env%domain_map(ispin)%pairs(ientry, 1)
1404 domain_map_local(2*ientry) = almo_scf_env%domain_map(ispin)%pairs(ientry, 2)
1406 CALL group%allgatherv(domain_map_local, domain_map_global, &
1407 domain_entries_cpu, offset_for_cpu)
1408 DEALLOCATE (domain_entries_cpu, offset_for_cpu)
1409 DEALLOCATE (domain_map_local)
1411 DEALLOCATE (almo_scf_env%domain_map(ispin)%index1)
1412 DEALLOCATE (almo_scf_env%domain_map(ispin)%pairs)
1413 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1414 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(global_entries/2, 2))
1415 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1416 almo_scf_env%domain_map(ispin)%index1(:) = 0
1422 max_neig = max_domain_neighbors
1423 max_neig_fails = .true.
1424 max_neig_loop:
DO WHILE (max_neig_fails)
1425 ALLOCATE (domain_grid(almo_scf_env%ndomains, 0:max_neig))
1426 domain_grid(:, :) = 0
1428 domain_grid(:, 0) = 1
1430 global_entries = global_entries/2
1431 DO ientry = 1, global_entries
1433 grid1 = domain_map_global(2*ientry)
1435 ineig = domain_map_global(2*(ientry - 1) + 1)
1437 IF (domain_grid(grid1, 0) .GT. max_neig)
THEN
1440 DEALLOCATE (domain_grid)
1441 max_neig = max_neig*2
1446 delayed_increment = .false.
1447 DO igrid = 1, domain_grid(grid1, 0)
1449 IF (ineig .LT. domain_grid(grid1, igrid))
THEN
1453 ineig = domain_grid(grid1, igrid)
1454 domain_grid(grid1, igrid) = neig_temp
1456 IF (domain_grid(grid1, igrid) .EQ. 0)
THEN
1458 domain_grid(grid1, igrid) = ineig
1460 delayed_increment = .true.
1464 IF (delayed_increment)
THEN
1465 domain_grid(grid1, 0) = domain_grid(grid1, 0) + 1
1468 cpabort(
"all records must be inserted")
1471 max_neig_fails = .false.
1472 END DO max_neig_loop
1473 DEALLOCATE (domain_map_global)
1476 DO idomain = 1, almo_scf_env%ndomains
1477 DO ineig = 1, domain_grid(idomain, 0) - 1
1478 almo_scf_env%domain_map(ispin)%pairs(ientry, 1) = domain_grid(idomain, ineig)
1479 almo_scf_env%domain_map(ispin)%pairs(ientry, 2) = idomain
1482 almo_scf_env%domain_map(ispin)%index1(idomain) = ientry
1484 DEALLOCATE (domain_grid)
1487 IF (almo_scf_env%nspins .EQ. 2)
THEN
1489 almo_scf_env%quench_t(1))
1490 almo_scf_env%domain_map(2)%pairs(:, :) = &
1491 almo_scf_env%domain_map(1)%pairs(:, :)
1492 almo_scf_env%domain_map(2)%index1(:) = &
1493 almo_scf_env%domain_map(1)%index1(:)
1500 DEALLOCATE (first_atom_of_molecule)
1501 DEALLOCATE (last_atom_of_molecule)
1521 CALL timestop(handle)
1538 CHARACTER(len=*),
PARAMETER :: routinen =
'calculate_w_matrix_almo'
1540 INTEGER :: handle, ispin
1541 REAL(kind=
dp) :: scaling
1542 TYPE(
dbcsr_type) :: tmp_nn1, tmp_no1, tmp_oo1, tmp_oo2
1544 CALL timeset(routinen, handle)
1546 IF (almo_scf_env%nspins == 1)
THEN
1552 DO ispin = 1, almo_scf_env%nspins
1554 CALL dbcsr_create(tmp_nn1, template=almo_scf_env%matrix_s(1), &
1555 matrix_type=dbcsr_type_no_symmetry)
1556 CALL dbcsr_create(tmp_no1, template=almo_scf_env%matrix_t(ispin), &
1557 matrix_type=dbcsr_type_no_symmetry)
1558 CALL dbcsr_create(tmp_oo1, template=almo_scf_env%matrix_sigma_inv(ispin), &
1559 matrix_type=dbcsr_type_no_symmetry)
1560 CALL dbcsr_create(tmp_oo2, template=almo_scf_env%matrix_sigma_inv(ispin), &
1561 matrix_type=dbcsr_type_no_symmetry)
1563 CALL dbcsr_copy(tmp_nn1, almo_scf_env%matrix_ks(ispin))
1565 CALL dbcsr_multiply(
"N",
"N", scaling, tmp_nn1, almo_scf_env%matrix_t(ispin), &
1566 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1568 CALL dbcsr_multiply(
"T",
"N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_no1, &
1569 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1571 CALL dbcsr_multiply(
"N",
"N", 1.0_dp, tmp_oo1, almo_scf_env%matrix_sigma_inv(ispin), &
1572 0.0_dp, tmp_oo2, filter_eps=almo_scf_env%eps_filter)
1574 CALL dbcsr_multiply(
"N",
"N", 1.0_dp, almo_scf_env%matrix_sigma_inv(ispin), tmp_oo2, &
1575 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1577 CALL dbcsr_multiply(
"N",
"N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_oo1, &
1578 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1580 CALL dbcsr_multiply(
"N",
"T", 1.0_dp, tmp_no1, almo_scf_env%matrix_t(ispin), &
1581 0.0_dp, tmp_nn1, filter_eps=almo_scf_env%eps_filter)
1582 CALL matrix_almo_to_qs(tmp_nn1, matrix_w(ispin)%matrix, almo_scf_env%mat_distr_aos)
1591 CALL timestop(handle)
double distance_squared(double *A, double *B)
Interface between ALMO SCF and QS.
subroutine, public almo_scf_update_ks_energy(qs_env, energy, energy_singles_corr)
update qs_env total energy
subroutine, public construct_qs_mos(qs_env, almo_scf_env)
Create MOs in the QS env to be able to return ALMOs to QS.
subroutine, public almo_dm_to_almo_ks(qs_env, matrix_p, matrix_ks, energy_total, eps_filter, mat_distr_aos, smear, kts_sum)
uses the ALMO density matrix to compute ALMO KS matrix and the new energy
subroutine, public calculate_w_matrix_almo(matrix_w, almo_scf_env)
Compute matrix W (energy-weighted density matrix) that is needed for the evaluation of forces.
subroutine, public matrix_almo_create(matrix_new, matrix_qs, almo_scf_env, name_new, size_keys, symmetry_new, spin_key, init_domains)
create the ALMO matrix templates
subroutine, public init_almo_ks_matrix_via_qs(qs_env, matrix_ks, mat_distr_aos, eps_filter)
Initialization of the QS and ALMO KS matrix.
subroutine, public almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
return density matrix to the qs_env
subroutine, public matrix_qs_to_almo(matrix_qs, matrix_almo, mat_distr_aos)
convert between two types of matrices: QS style to ALMO style
subroutine, public almo_scf_construct_quencher(qs_env, almo_scf_env)
Creates the matrix that imposes absolute locality on MOs.
Types for all ALMO-based methods.
integer, parameter, public almo_mat_dim_occ
integer, parameter, public almo_mat_dim_virt_full
integer, parameter, public almo_mat_dim_aobasis
integer, parameter, public almo_mat_dim_virt
integer, parameter, public almo_mat_dim_virt_disc
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Handles all functions related to the CELL.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_distribution_release(dist)
...
subroutine, public dbcsr_distribution_new(dist, template, group, pgrid, row_dist, col_dist, reuse_arrays)
...
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_get_block_p(matrix, row, col, block, found, row_size, col_size)
...
subroutine, public dbcsr_multiply(transa, transb, alpha, matrix_a, matrix_b, beta, matrix_c, first_row, last_row, first_column, last_column, first_k, last_k, retain_sparsity, filter_eps, flop)
...
subroutine, public dbcsr_get_info(matrix, nblkrows_total, nblkcols_total, nfullrows_total, nfullcols_total, nblkrows_local, nblkcols_local, nfullrows_local, nfullcols_local, my_prow, my_pcol, local_rows, local_cols, proc_row_dist, proc_col_dist, row_blk_size, col_blk_size, row_blk_offset, col_blk_offset, distribution, name, matrix_type, group)
...
subroutine, public dbcsr_get_stored_coordinates(matrix, row, column, processor)
...
subroutine, public dbcsr_work_create(matrix, nblks_guess, sizedata_guess, n, work_mutable)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_finalize(matrix)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_release(matrix)
...
subroutine, public dbcsr_complete_redistribute(matrix, redist)
...
integer function, public dbcsr_get_num_blocks(matrix)
...
subroutine, public dbcsr_put_block(matrix, row, col, block, summation)
...
subroutine, public dbcsr_distribution_get(dist, row_dist, col_dist, nrows, ncols, has_threads, group, mynode, numnodes, nprows, npcols, myprow, mypcol, pgrid, subgroups_defined, prow_group, pcol_group)
...
Routines that link DBCSR and CP2K concepts together.
subroutine, public cp_dbcsr_alloc_block_from_nbl(matrix, sab_orb, desymmetrize)
allocate the blocks of a dbcsr based on the neighbor list
DBCSR operations in CP2K.
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_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 ...
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
Defines the basic variable types.
integer, parameter, public dp
Interface to the message passing library MPI.
Define the data structure for the molecule information.
subroutine, public get_molecule_set_info(molecule_set, atom_to_mol, mol_to_first_atom, mol_to_last_atom, mol_to_nelectrons, mol_to_nbasis, mol_to_charge, mol_to_multiplicity)
returns information about molecules in the set.
Define the data structure for the particle information.
Perform a QUICKSTEP wavefunction optimization (single point)
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_pp, 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, harris_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, eeq, rhs)
Get the QUICKSTEP environment.
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, harris_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, wanniercentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs)
Set the QUICKSTEP environment.
routines that build the Kohn-Sham matrix (i.e calculate the coulomb and xc parts
subroutine, public qs_ks_update_qs_env(qs_env, calculate_forces, just_energy, print_active)
updates the Kohn Sham matrix of the given qs_env (facility method)
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
subroutine, public qs_ks_did_change(ks_env, s_mstruct_changed, rho_changed, potential_changed, full_reset)
tells that some of the things relevant to the ks calculation did change. has to be called when change...
Definition and initialisation of the mo data type.
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
subroutine, public deallocate_mo_set(mo_set)
Deallocate a wavefunction data structure.
subroutine, public init_mo_set(mo_set, fm_pool, fm_ref, fm_struct, name)
initializes an allocated mo_set. eigenvalues, mo_coeff, occupation_numbers are valid only after this ...
Define the neighbor list data types and the corresponding functionality.
subroutine, public neighbor_list_iterator_create(iterator_set, nl, search, nthread)
Neighbor list iterator functions.
subroutine, public neighbor_list_iterator_release(iterator_set)
...
integer function, public neighbor_list_iterate(iterator_set, mepos)
...
subroutine, public get_iterator_info(iterator_set, mepos, ikind, jkind, nkind, ilist, nlist, inode, nnode, iatom, jatom, r, cell)
...
methods of the rho structure (defined in qs_rho_types)
subroutine, public qs_rho_update_rho(rho_struct, qs_env, rho_xc_external, local_rho_set, task_list_external, task_list_external_soft, pw_env_external, para_env_external)
updates rho_r and rho_g to the rhorho_ao. if use_kinetic_energy_density also computes tau_r and tau_g...
superstucture that hold various representations of the density and keeps track of which ones are vali...
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
module that contains the definitions of the scf types
subroutine, public scf_env_create(scf_env)
allocates and initialize an scf_env
Type defining parameters related to the simulation cell.
keeps the information about the structure of a full matrix
type of a logger, at the moment it contains just a print level starting at which level it should be l...
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
keeps the density in various representations, keeping track of which ones are valid.