45 pw_integrate_function,&
62 #include "./base/base_uses.f90"
68 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'xray_diffraction'
86 TYPE(qs_environment_type),
POINTER :: qs_env
87 INTEGER,
INTENT(IN) :: unit_number
88 REAL(kind=
dp),
INTENT(IN) :: q_max
90 CHARACTER(LEN=*),
PARAMETER :: routinen =
'xray_diffraction_spectrum'
91 INTEGER,
PARAMETER :: nblock = 100
93 INTEGER :: handle, i, ig, ig_shell, ipe, ishell, &
94 jg, ng, npe, nshell, nshell_gather
95 INTEGER(KIND=int_8) :: ngpts
96 INTEGER,
DIMENSION(3) :: npts
97 INTEGER,
DIMENSION(:),
POINTER :: aux_index, ng_shell, ng_shell_gather, &
99 REAL(kind=
dp) :: cutoff, f, f2, q, rho_hard, rho_soft, &
101 REAL(kind=
dp),
DIMENSION(3) :: dg, dr
102 REAL(kind=
dp),
DIMENSION(:),
POINTER :: f2sum, f2sum_gather, f4sum, f4sum_gather, fmax, &
103 fmax_gather, fmin, fmin_gather, fsum, fsum_gather, gsq, q_shell, q_shell_gather
104 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
105 TYPE(dft_control_type),
POINTER :: dft_control
106 TYPE(mp_para_env_type),
POINTER :: para_env
107 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
108 TYPE(pw_c1d_gs_type) :: rhotot_elec_gspace
109 TYPE(pw_env_type),
POINTER :: pw_env
110 TYPE(pw_pool_type),
POINTER :: auxbas_pw_pool
111 TYPE(qs_rho_type),
POINTER :: rho
112 TYPE(rho_atom_type),
DIMENSION(:),
POINTER :: rho_atom_set
114 cpassert(
ASSOCIATED(qs_env))
116 CALL timeset(routinen, handle)
118 NULLIFY (atomic_kind_set)
120 NULLIFY (auxbas_pw_pool)
121 NULLIFY (dft_control)
123 NULLIFY (f2sum_gather)
125 NULLIFY (f4sum_gather)
127 NULLIFY (fmax_gather)
129 NULLIFY (fmin_gather)
131 NULLIFY (fsum_gather)
134 NULLIFY (ng_shell_gather)
138 NULLIFY (particle_set)
141 NULLIFY (q_shell_gather)
143 NULLIFY (rho_atom_set)
148 atomic_kind_set=atomic_kind_set, &
149 dft_control=dft_control, &
151 particle_set=particle_set, &
154 rho_atom_set=rho_atom_set)
157 auxbas_pw_pool=auxbas_pw_pool)
159 npe = para_env%num_pe
163 CALL auxbas_pw_pool%create_pw(pw=rhotot_elec_gspace)
164 CALL pw_zero(rhotot_elec_gspace)
173 dg(:) =
twopi/(npts(:)*dr(:))
178 auxbas_pw_pool=auxbas_pw_pool, &
179 rhotot_elec_gspace=rhotot_elec_gspace, &
184 rho_total = rho_hard + rho_soft
193 CALL reallocate(q_shell, 1, nblock)
194 CALL reallocate(ng_shell, 1, nblock)
200 q_shell(1) = sqrt(gsq(1))
204 cpassert(gsq(ig) >= gsq(jg))
205 IF (abs(gsq(ig) - gsq(jg)) > 1.0e-12_dp)
THEN
207 IF (nshell >
SIZE(q_shell))
THEN
208 CALL reallocate(q_shell, 1,
SIZE(q_shell) + nblock)
209 CALL reallocate(ng_shell, 1,
SIZE(ng_shell) + nblock)
220 ng_shell(nshell) = ng_shell(nshell) + 1
224 CALL reallocate(q_shell, 1, nshell)
225 CALL reallocate(ng_shell, 1, nshell)
226 CALL reallocate(fmin, 1, nshell)
227 CALL reallocate(fmax, 1, nshell)
228 CALL reallocate(fsum, 1, nshell)
229 CALL reallocate(f2sum, 1, nshell)
230 CALL reallocate(f4sum, 1, nshell)
233 DO ishell = 1, nshell
234 fmin(ishell) = huge(0.0_dp)
235 fmax(ishell) = 0.0_dp
236 fsum(ishell) = 0.0_dp
237 f2sum(ishell) = 0.0_dp
238 f4sum(ishell) = 0.0_dp
239 DO ig_shell = 1, ng_shell(ishell)
240 f = abs(rhotot_elec_gspace%array(ig + ig_shell))
241 fmin(ishell) = min(fmin(ishell), f)
242 fmax(ishell) = max(fmax(ishell), f)
243 fsum(ishell) = fsum(ishell) + f
245 f2sum(ishell) = f2sum(ishell) + f2
246 f4sum(ishell) = f4sum(ishell) + f2*f2
248 ig = ig + ng_shell(ishell)
251 CALL reallocate(nshell_pe, 0, npe - 1)
252 CALL reallocate(offset_pe, 0, npe - 1)
256 CALL para_env%gather(nshell, nshell_pe)
261 IF (unit_number > 0)
THEN
262 nshell_gather = sum(nshell_pe)
265 offset_pe(ipe) = offset_pe(ipe - 1) + nshell_pe(ipe - 1)
271 CALL reallocate(q_shell_gather, 1, nshell_gather)
272 CALL reallocate(ng_shell_gather, 1, nshell_gather)
273 CALL reallocate(fmin_gather, 1, nshell_gather)
274 CALL reallocate(fmax_gather, 1, nshell_gather)
275 CALL reallocate(fsum_gather, 1, nshell_gather)
276 CALL reallocate(f2sum_gather, 1, nshell_gather)
277 CALL reallocate(f4sum_gather, 1, nshell_gather)
279 CALL para_env%gatherv(q_shell, q_shell_gather, nshell_pe, offset_pe)
280 CALL para_env%gatherv(ng_shell, ng_shell_gather, nshell_pe, offset_pe)
281 CALL para_env%gatherv(fmax, fmax_gather, nshell_pe, offset_pe)
282 CALL para_env%gatherv(fmin, fmin_gather, nshell_pe, offset_pe)
283 CALL para_env%gatherv(fsum, fsum_gather, nshell_pe, offset_pe)
284 CALL para_env%gatherv(f2sum, f2sum_gather, nshell_pe, offset_pe)
285 CALL para_env%gatherv(f4sum, f4sum_gather, nshell_pe, offset_pe)
287 IF (
ASSOCIATED(offset_pe))
THEN
288 DEALLOCATE (offset_pe)
291 IF (
ASSOCIATED(nshell_pe))
THEN
292 DEALLOCATE (nshell_pe)
297 IF (unit_number > 0)
THEN
299 CALL reallocate(aux_index, 1, nshell_gather)
303 CALL sort(q_shell_gather, nshell_gather, aux_index)
308 CALL reallocate(q_shell, 1, nshell_gather)
309 CALL reallocate(ng_shell, 1, nshell_gather)
310 CALL reallocate(fmin, 1, nshell_gather)
311 CALL reallocate(fmax, 1, nshell_gather)
312 CALL reallocate(fsum, 1, nshell_gather)
313 CALL reallocate(f2sum, 1, nshell_gather)
314 CALL reallocate(f4sum, 1, nshell_gather)
318 q_shell(1) = q_shell_gather(1)
320 ng_shell(1) = ng_shell_gather(i)
321 fmin(1) = fmin_gather(i)
322 fmax(1) = fmax_gather(i)
323 fsum(1) = fsum_gather(i)
324 f2sum(1) = f2sum_gather(i)
325 f4sum(1) = f4sum_gather(i)
327 DO ig = 2, nshell_gather
329 IF (abs(q_shell_gather(ig) - q_shell_gather(jg)) > 1.0e-12_dp)
THEN
331 q_shell(nshell) = q_shell_gather(ig)
332 ng_shell(nshell) = ng_shell_gather(i)
333 fmin(nshell) = fmin_gather(i)
334 fmax(nshell) = fmax_gather(i)
335 fsum(nshell) = fsum_gather(i)
336 f2sum(nshell) = f2sum_gather(i)
337 f4sum(nshell) = f4sum_gather(i)
340 ng_shell(nshell) = ng_shell(nshell) + ng_shell_gather(i)
341 fmin(nshell) = min(fmin(nshell), fmin_gather(i))
342 fmax(nshell) = max(fmax(nshell), fmax_gather(i))
343 fsum(nshell) = fsum(nshell) + fsum_gather(i)
344 f2sum(nshell) = f2sum(nshell) + f2sum_gather(i)
345 f4sum(nshell) = f4sum(nshell) + f4sum_gather(i)
351 IF (
ASSOCIATED(aux_index))
THEN
352 DEALLOCATE (aux_index)
357 CALL reallocate(q_shell, 1, nshell)
358 CALL reallocate(ng_shell, 1, nshell)
359 CALL reallocate(fmin, 1, nshell)
360 CALL reallocate(fmax, 1, nshell)
361 CALL reallocate(fsum, 1, nshell)
362 CALL reallocate(f2sum, 1, nshell)
363 CALL reallocate(f4sum, 1, nshell)
367 WRITE (unit=unit_number, fmt=
"(A)") &
369 "# Coherent X-ray diffraction spectrum", &
371 WRITE (unit=unit_number, fmt=
"(A,1X,F20.10)") &
372 "# Soft electronic charge (G-space) :", rho_soft, &
373 "# Hard electronic charge (G-space) :", rho_hard, &
374 "# Total electronic charge (G-space):", rho_total, &
375 "# Density cutoff [Rydberg] :", 2.0_dp*cutoff, &
376 "# q(min) [1/Angstrom] :", q_shell(2)/
angstrom, &
377 "# q(max) [1/Angstrom] :", q_shell(nshell)/
angstrom, &
378 "# q(max) [1/Angstrom] (requested) :", q_max/
angstrom
379 WRITE (unit=unit_number, fmt=
"(A,2X,I8)") &
380 "# Number of g-vectors (grid points):", ngpts, &
381 "# Number of g-vector shells :", nshell
382 WRITE (unit=unit_number, fmt=
"(A,3(1X,I6))") &
383 "# Grid size (a,b,c) :", npts(1:3)
384 WRITE (unit=unit_number, fmt=
"(A,3F7.3)") &
385 "# dg [1/Angstrom] :", dg(1:3)/
angstrom, &
386 "# dr [Angstrom] :", dr(1:3)*
angstrom
387 WRITE (unit=unit_number, fmt=
"(A)") &
389 "# shell points q [1/A] <|F(q)|^2> Min(|F(q)|)"// &
390 " Max(|F(q)|) <|F(q)|>^2 <|F(q)|^4>"
392 DO ishell = 1, nshell
393 WRITE (unit=unit_number, fmt=
"(T2,I6,2X,I6,5(1X,F15.6),1X,ES15.6)") &
397 f2sum(ishell)/real(ng_shell(ishell), kind=
dp), &
400 (fsum(ishell)/real(ng_shell(ishell), kind=
dp))**2, &
401 f4sum(ishell)/real(ng_shell(ishell), kind=
dp)
408 IF (
ASSOCIATED(fmin))
THEN
412 IF (
ASSOCIATED(fmax))
THEN
416 IF (
ASSOCIATED(fsum))
THEN
420 IF (
ASSOCIATED(f2sum))
THEN
424 IF (
ASSOCIATED(f4sum))
THEN
428 IF (
ASSOCIATED(ng_shell))
THEN
429 DEALLOCATE (ng_shell)
432 IF (
ASSOCIATED(q_shell))
THEN
436 IF (
ASSOCIATED(fmin_gather))
THEN
437 DEALLOCATE (fmin_gather)
440 IF (
ASSOCIATED(fmax_gather))
THEN
441 DEALLOCATE (fmax_gather)
444 IF (
ASSOCIATED(fsum_gather))
THEN
445 DEALLOCATE (fsum_gather)
448 IF (
ASSOCIATED(f2sum_gather))
THEN
449 DEALLOCATE (f2sum_gather)
452 IF (
ASSOCIATED(f4sum_gather))
THEN
453 DEALLOCATE (f4sum_gather)
456 IF (
ASSOCIATED(ng_shell_gather))
THEN
457 DEALLOCATE (ng_shell_gather)
460 IF (
ASSOCIATED(q_shell_gather))
THEN
461 DEALLOCATE (q_shell_gather)
464 CALL auxbas_pw_pool%give_back_pw(rhotot_elec_gspace)
466 CALL timestop(handle)
485 rhotot_elec_gspace, q_max, rho_hard, &
488 TYPE(qs_environment_type),
POINTER :: qs_env
489 TYPE(pw_pool_type),
POINTER :: auxbas_pw_pool
490 TYPE(pw_c1d_gs_type),
INTENT(INOUT) :: rhotot_elec_gspace
491 REAL(kind=
dp),
INTENT(IN) :: q_max
492 REAL(kind=
dp),
INTENT(OUT) :: rho_hard, rho_soft
493 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: fsign
495 CHARACTER(LEN=*),
PARAMETER :: routinen =
'calculate_rhotot_elec_gspace'
497 INTEGER ::
atom, handle, iatom, ico, ico1_pgf, ico1_set, ikind, ipgf, iset, iso, iso1_pgf, &
498 iso1_set, ison, ispin, jco, jco1_pgf, jco1_set, jpgf, jset, jso, jso1_pgf, jso1_set, &
499 json, la, lb, maxco, maxso, na, natom, nb, ncoa, ncob, ncotot, nkind, nsatbas, nset, &
500 nsoa, nsob, nsotot, nspin
501 INTEGER,
DIMENSION(:),
POINTER :: atom_list, lmax, lmin, npgf, o2nindex
502 LOGICAL :: orthorhombic, paw_atom
503 REAL(kind=
dp) :: alpha, eps_rho_gspace, rho_total, scale, &
505 REAL(kind=
dp),
DIMENSION(3) :: ra
506 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: delta_cpc, pab, work, zet
507 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
508 TYPE(cell_type),
POINTER :: cell
509 TYPE(dft_control_type),
POINTER :: dft_control
510 TYPE(gto_basis_set_type),
POINTER :: basis_1c_set
511 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
512 TYPE(pw_c1d_gs_type) :: rho_elec_gspace
513 TYPE(pw_r3d_rs_type),
DIMENSION(:),
POINTER :: rho_r
514 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
515 TYPE(qs_rho_type),
POINTER :: rho
516 TYPE(rho_atom_coeff),
DIMENSION(:),
POINTER :: cpc_h, cpc_s
517 TYPE(rho_atom_type),
DIMENSION(:),
POINTER :: rho_atom_set
518 TYPE(rho_atom_type),
POINTER :: rho_atom
520 cpassert(
ASSOCIATED(qs_env))
521 cpassert(
ASSOCIATED(auxbas_pw_pool))
523 CALL timeset(routinen, handle)
526 NULLIFY (atomic_kind_set)
527 NULLIFY (qs_kind_set)
532 NULLIFY (dft_control)
536 NULLIFY (basis_1c_set)
538 NULLIFY (particle_set)
541 NULLIFY (rho_atom_set)
546 atomic_kind_set=atomic_kind_set, &
547 qs_kind_set=qs_kind_set, &
549 dft_control=dft_control, &
550 particle_set=particle_set, &
552 rho_atom_set=rho_atom_set)
555 eps_rho_gspace = dft_control%qs_control%eps_rho_gspace
556 nkind =
SIZE(atomic_kind_set)
557 nspin = dft_control%nspins
561 CALL auxbas_pw_pool%create_pw(pw=rho_elec_gspace)
563 CALL pw_zero(rhotot_elec_gspace)
566 CALL pw_zero(rho_elec_gspace)
567 CALL pw_transfer(rho_r(ispin), rho_elec_gspace)
568 IF (
PRESENT(fsign) .AND. (ispin == 2))
THEN
573 CALL pw_axpy(rho_elec_gspace, rhotot_elec_gspace, alpha=alpha)
579 CALL auxbas_pw_pool%give_back_pw(rho_elec_gspace)
581 rho_soft = pw_integrate_function(rhotot_elec_gspace, isign=-1)
585 orthorhombic=orthorhombic)
586 IF (.NOT. orthorhombic)
THEN
587 CALL cp_abort(__location__, &
588 "The calculation of XRD spectra for non-orthorhombic cells is not implemented")
591 CALL pw_scale(rhotot_elec_gspace, volume)
601 atom_list=atom_list, &
605 basis_set=basis_1c_set, &
606 basis_type=
"GAPW_1C", &
609 IF (.NOT. paw_atom) cycle
624 CALL reallocate(delta_cpc, 1, nsatbas, 1, nsatbas)
625 CALL reallocate(pab, 1, ncotot, 1, ncotot)
626 CALL reallocate(work, 1, maxso, 1, maxco)
630 atom = atom_list(iatom)
631 rho_atom => rho_atom_set(
atom)
637 ra(:) =
pbc(particle_set(iatom)%r, cell)
642 IF (
PRESENT(fsign) .AND. (ispin == 2))
THEN
647 delta_cpc = delta_cpc + alpha*(cpc_h(ispin)%r_coef - cpc_s(ispin)%r_coef)
653 ico1_set = (iset - 1)*maxco + 1
654 iso1_set = (iset - 1)*maxso + 1
658 jco1_set = (jset - 1)*maxco + 1
659 jso1_set = (jset - 1)*maxso + 1
662 DO ipgf = 1, npgf(iset)
663 ico1_pgf = ico1_set + (ipgf - 1)*ncoa
664 iso1_pgf = iso1_set + (ipgf - 1)*nsoa
665 DO jpgf = 1, npgf(jset)
666 jco1_pgf = jco1_set + (jpgf - 1)*ncob
667 jso1_pgf = jso1_set + (jpgf - 1)*nsob
668 ico = ico1_pgf +
ncoset(lmin(iset) - 1)
669 iso = iso1_pgf +
nsoset(lmin(iset) - 1)
673 DO la = lmin(iset), lmax(iset)
674 jco = jco1_pgf +
ncoset(lmin(jset) - 1)
675 jso = jso1_pgf +
nsoset(lmin(jset) - 1)
676 DO lb = lmin(jset), lmax(jset)
680 delta_cpc(ison:ison +
nso(la) - 1, json),
SIZE(delta_cpc, 1), &
685 0.0_dp, pab(ico:ico +
nco(la) - 1, jco),
SIZE(pab, 1))
698 CALL collocate_pgf_product_gspace( &
700 zeta=zet(ipgf, iset), &
703 zetb=zet(jpgf, jset), &
706 rab=(/0.0_dp, 0.0_dp, 0.0_dp/), &
712 eps_rho_gspace=eps_rho_gspace, &
713 gsq_max=q_max*q_max, &
714 pw=rhotot_elec_gspace)
721 DEALLOCATE (o2nindex)
724 rho_total = pw_integrate_function(rhotot_elec_gspace, isign=-1)/volume
726 rho_hard = rho_total - rho_soft
730 IF (
ASSOCIATED(delta_cpc))
THEN
731 DEALLOCATE (delta_cpc)
734 IF (
ASSOCIATED(work))
THEN
738 IF (
ASSOCIATED(pab))
THEN
742 CALL timestop(handle)
765 SUBROUTINE collocate_pgf_product_gspace(la_max, zeta, la_min, &
766 lb_max, zetb, lb_min, &
767 ra, rab, rab2, scale, pab, na, nb, &
768 eps_rho_gspace, gsq_max, pw)
772 INTEGER,
INTENT(IN) :: la_max
773 REAL(
dp),
INTENT(IN) :: zeta
774 INTEGER,
INTENT(IN) :: la_min, lb_max
775 REAL(
dp),
INTENT(IN) :: zetb
776 INTEGER,
INTENT(IN) :: lb_min
777 REAL(
dp),
DIMENSION(3),
INTENT(IN) :: ra, rab
778 REAL(
dp),
INTENT(IN) :: rab2, scale
779 REAL(
dp),
DIMENSION(:, :),
POINTER :: pab
780 INTEGER,
INTENT(IN) :: na, nb
781 REAL(
dp),
INTENT(IN) :: eps_rho_gspace, gsq_max
782 TYPE(pw_c1d_gs_type),
INTENT(IN) :: pw
784 CHARACTER(LEN=*),
PARAMETER :: routinen =
'collocate_pgf_product_gspace'
786 COMPLEX(dp),
DIMENSION(3) :: phasefactor
787 COMPLEX(dp),
DIMENSION(:),
POINTER :: rag, rbg
788 COMPLEX(dp),
DIMENSION(:, :, :, :),
POINTER :: cubeaxis
789 INTEGER :: ax, ay, az, bx, by, bz, handle, i, ico, &
790 ig, ig2, jco, jg, kg, la, lb, &
791 lb_cube_min, lb_grid, ub_cube_max, &
793 INTEGER,
DIMENSION(3) :: lb_cube, ub_cube
794 REAL(
dp) :: f, fa, fb, pij, prefactor, rzetp, &
796 REAL(
dp),
DIMENSION(3) :: dg, expfactor, fap, fbp, rap, rbp, rp
797 REAL(
dp),
DIMENSION(:),
POINTER :: g
799 CALL timeset(routinen, handle)
801 dg(:) =
twopi/(pw%pw_grid%npts(:)*pw%pw_grid%dr(:))
807 rbp(:) = rap(:) - rab(:)
808 rp(:) = ra(:) + rap(:)
809 twozetp = 2.0_dp*zetp
810 fap(:) = twozetp*rap(:)
811 fbp(:) = twozetp*rbp(:)
813 prefactor = scale*sqrt((
pi*rzetp)**3)*exp(-zeta*f*rab2)
814 phasefactor(:) = exp(cmplx(0.0_dp, -rp(:)*dg(:), kind=
dp))
815 expfactor(:) = exp(-0.25*rzetp*dg(:)*dg(:))
817 lb_cube(:) = pw%pw_grid%bounds(1, :)
818 ub_cube(:) = pw%pw_grid%bounds(2, :)
820 lb_cube_min = minval(lb_cube(:))
821 ub_cube_max = maxval(ub_cube(:))
823 NULLIFY (cubeaxis, g, rag, rbg)
825 CALL reallocate(cubeaxis, lb_cube_min, ub_cube_max, 1, 3, 0, la_max, 0, lb_max)
826 CALL reallocate(g, lb_cube_min, ub_cube_max)
827 CALL reallocate(rag, lb_cube_min, ub_cube_max)
828 CALL reallocate(rbg, lb_cube_min, ub_cube_max)
830 lb_grid = lbound(pw%array, 1)
831 ub_grid = ubound(pw%array, 1)
835 DO ig = lb_cube(i), ub_cube(i)
837 cubeaxis(ig, i, 0, 0) = expfactor(i)**ig2*phasefactor(i)**ig
841 DO ig = lb_cube(i), ub_cube(i)
842 g(ig) = real(ig,
dp)*dg(i)
843 rag(ig) = cmplx(fap(i), -g(ig), kind=
dp)
844 cubeaxis(ig, i, 1, 0) = rag(ig)*cubeaxis(ig, i, 0, 0)
847 fa = real(la - 1,
dp)*twozetp
848 DO ig = lb_cube(i), ub_cube(i)
849 cubeaxis(ig, i, la, 0) = rag(ig)*cubeaxis(ig, i, la - 1, 0) + &
850 fa*cubeaxis(ig, i, la - 2, 0)
855 DO ig = lb_cube(i), ub_cube(i)
856 rbg(ig) = cmplx(fbp(i), -g(ig), kind=
dp)
857 cubeaxis(ig, i, 0, 1) = rbg(ig)*cubeaxis(ig, i, 0, 0)
858 cubeaxis(ig, i, 1, 1) = rbg(ig)*cubeaxis(ig, i, 1, 0) + &
859 fa*cubeaxis(ig, i, 0, 0)
862 fb = real(lb - 1,
dp)*twozetp
863 DO ig = lb_cube(i), ub_cube(i)
864 cubeaxis(ig, i, 0, lb) = rbg(ig)*cubeaxis(ig, i, 0, lb - 1) + &
865 fb*cubeaxis(ig, i, 0, lb - 2)
866 cubeaxis(ig, i, 1, lb) = rbg(ig)*cubeaxis(ig, i, 1, lb - 1) + &
867 fb*cubeaxis(ig, i, 1, lb - 2) + &
868 fa*cubeaxis(ig, i, 0, lb - 1)
872 fa = real(la,
dp)*twozetp
873 DO ig = lb_cube(i), ub_cube(i)
874 cubeaxis(ig, i, la, 1) = rbg(ig)*cubeaxis(ig, i, la, 0) + &
875 fa*cubeaxis(ig, i, la - 1, 0)
878 fb = real(lb - 1,
dp)*twozetp
879 DO ig = lb_cube(i), ub_cube(i)
880 cubeaxis(ig, i, la, lb) = rbg(ig)*cubeaxis(ig, i, la, lb - 1) + &
881 fb*cubeaxis(ig, i, la, lb - 2) + &
882 fa*cubeaxis(ig, i, la - 1, lb - 1)
889 DO ig = lb_cube(i), ub_cube(i)
890 g(ig) = real(ig,
dp)*dg(i)
891 rbg(ig) = cmplx(fbp(i), -g(ig), kind=
dp)
892 cubeaxis(ig, i, 0, 1) = rbg(ig)*cubeaxis(ig, i, 0, 0)
895 fb = real(lb - 1,
dp)*twozetp
896 DO ig = lb_cube(i), ub_cube(i)
897 cubeaxis(ig, i, 0, lb) = rbg(ig)*cubeaxis(ig, i, 0, lb - 1) + &
898 fb*cubeaxis(ig, i, 0, lb - 2)
908 IF (la + lb == 0) cycle
909 fa = (1.0_dp/twozetp)**(la + lb)
911 DO ig = lb_cube(i), ub_cube(i)
912 cubeaxis(ig, i, la, lb) = fa*cubeaxis(ig, i, la, lb)
928 pij = prefactor*pab(na + ico, nb + jco)
930 IF (abs(pij) < eps_rho_gspace) cycle
936 DO i = lb_grid, ub_grid
937 IF (pw%pw_grid%gsq(i) > gsq_max) cycle
938 ig = pw%pw_grid%g_hat(1, i)
939 jg = pw%pw_grid%g_hat(2, i)
940 kg = pw%pw_grid%g_hat(3, i)
941 pw%array(i) = pw%array(i) + pij*cubeaxis(ig, 1, ax, bx)* &
942 cubeaxis(jg, 2, ay, by)* &
943 cubeaxis(kg, 3, az, bz)
950 DEALLOCATE (cubeaxis)
955 CALL timestop(handle)
957 END SUBROUTINE collocate_pgf_product_gspace
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
static void dgemm(const char transa, const char transb, const int m, const int n, const int k, const double alpha, const double *a, const int lda, const double *b, const int ldb, const double beta, double *c, const int ldc)
Convenient wrapper to hide Fortran nature of dgemm_, swapping a and b.
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.
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius)
...
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public krack2002
Handles all functions related to the CELL.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Defines the basic variable types.
integer, parameter, public int_8
integer, parameter, public dp
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), parameter, public twopi
Utility routines for the memory handling.
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nco
integer, dimension(:), allocatable, public nsoset
integer, dimension(:), allocatable, public ncoset
integer, dimension(:, :), allocatable, public indco
integer, dimension(:), allocatable, public nso
Define the data structure for the particle information.
subroutine, public get_paw_basis_info(basis_1c, o2nindex, n2oindex, nsatbas)
Return some info on the PAW basis derived from a GTO basis 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
This module defines the grid data type and some basic operations on it.
subroutine, public get_pw_grid_info(pw_grid, id_nr, mode, vol, dvol, npts, ngpts, ngpts_cut, dr, cutoff, orthorhombic, gvectors, gsquare)
Access to information stored in the pw_grid_type.
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
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.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_r3d_rs_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, U_of_dft_plus_u, J_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, J0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public get_rho_atom(rho_atom, cpc_h, cpc_s, rho_rad_h, rho_rad_s, drho_rad_h, drho_rad_s, vrho_rad_h, vrho_rad_s, rho_rad_h_d, rho_rad_s_d, ga_Vlocal_gb_h, ga_Vlocal_gb_s, int_scr_h, int_scr_s)
...
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...
All kind of helpful little routines.
subroutine, public xray_diffraction_spectrum(qs_env, unit_number, q_max)
Calculate the coherent X-ray diffraction spectrum using the total electronic density in reciprocal sp...
subroutine, public calculate_rhotot_elec_gspace(qs_env, auxbas_pw_pool, rhotot_elec_gspace, q_max, rho_hard, rho_soft, fsign)
The total electronic density in reciprocal space (g-space) is calculated.