58#include "./base/base_uses.f90"
64 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'atom_types'
67 INTEGER,
PARAMETER ::
lmat = 5
74 INTEGER,
PARAMETER :: nmax = 25
80 INTEGER,
DIMENSION(0:lmat) :: nbas = 0
81 INTEGER,
DIMENSION(0:lmat) :: nprim = 0
82 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: am => null()
83 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: cm => null()
84 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: as => null()
85 INTEGER,
DIMENSION(:, :),
POINTER :: ns => null()
86 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: bf => null()
87 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: dbf => null()
88 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: ddbf => null()
89 REAL(kind=
dp) :: eps_eig = 0.0_dp
91 LOGICAL :: geometrical = .false.
92 REAL(kind=
dp) :: aval = 0.0_dp, cval = 0.0_dp
93 INTEGER,
DIMENSION(0:lmat) :: start = 0
99 CHARACTER(LEN=2) :: symbol =
""
100 CHARACTER(LEN=default_string_length) :: pname =
""
101 INTEGER,
DIMENSION(0:lmat) :: econf = 0
102 REAL(
dp) :: zion = 0.0_dp
103 REAL(
dp) :: rc = 0.0_dp
105 REAL(
dp),
DIMENSION(5) :: cl = 0.0_dp
106 INTEGER,
DIMENSION(0:lmat) :: nl = 0
107 REAL(
dp),
DIMENSION(0:lmat) :: rcnl = 0.0_dp
108 REAL(
dp),
DIMENSION(4, 4, 0:lmat) :: hnl = 0.0_dp
110 LOGICAL :: soc = .false.
111 REAL(
dp),
DIMENSION(4, 4, 0:lmat) :: knl = 0.0_dp
114 LOGICAL :: nlcc = .false.
115 INTEGER :: nexp_nlcc = 0
116 REAL(kind=
dp),
DIMENSION(10) :: alpha_nlcc = 0.0_dp
117 INTEGER,
DIMENSION(10) :: nct_nlcc = 0
118 REAL(kind=
dp),
DIMENSION(4, 10) :: cval_nlcc = 0.0_dp
120 LOGICAL :: lsdpot = .false.
121 INTEGER :: nexp_lsd = 0
122 REAL(kind=
dp),
DIMENSION(10) :: alpha_lsd = 0.0_dp
123 INTEGER,
DIMENSION(10) :: nct_lsd = 0
124 REAL(kind=
dp),
DIMENSION(4, 10) :: cval_lsd = 0.0_dp
126 LOGICAL :: lpotextended = .false.
127 INTEGER :: nexp_lpot = 0
128 REAL(kind=
dp),
DIMENSION(10) :: alpha_lpot = 0.0_dp
129 INTEGER,
DIMENSION(10) :: nct_lpot = 0
130 REAL(kind=
dp),
DIMENSION(4, 10) :: cval_lpot = 0.0_dp
134 CHARACTER(LEN=2) :: symbol =
""
135 CHARACTER(LEN=default_string_length) :: pname =
""
136 INTEGER,
DIMENSION(0:lmat) :: econf = 0
137 REAL(
dp) :: zion = 0.0_dp
140 INTEGER,
DIMENSION(1:15) :: nrloc = 0
141 REAL(
dp),
DIMENSION(1:15) :: aloc = 0.0_dp
142 REAL(
dp),
DIMENSION(1:15) :: bloc = 0.0_dp
143 INTEGER,
DIMENSION(0:10) :: npot = 0
144 INTEGER,
DIMENSION(1:15, 0:10) :: nrpot = 0
145 REAL(
dp),
DIMENSION(1:15, 0:10) :: apot = 0.0_dp
146 REAL(
dp),
DIMENSION(1:15, 0:10) :: bpot = 0.0_dp
150 CHARACTER(LEN=2) :: symbol =
""
151 CHARACTER(LEN=default_string_length) :: pname =
""
152 INTEGER,
DIMENSION(0:lmat) :: econf = 0
153 REAL(
dp) :: zion = 0.0_dp
155 LOGICAL :: has_nonlocal = .false.
156 INTEGER :: n_nonlocal = 0
157 LOGICAL,
DIMENSION(0:5) :: is_nonlocal = .false.
158 REAL(kind=
dp),
DIMENSION(nmax) :: a_nonlocal = 0.0_dp
159 REAL(kind=
dp),
DIMENSION(nmax, 0:lmat) :: h_nonlocal = 0.0_dp
160 REAL(kind=
dp),
DIMENSION(nmax, nmax, 0:lmat) :: c_nonlocal = 0.0_dp
161 INTEGER :: n_local = 0
162 REAL(kind=
dp) :: ac_local = 0.0_dp
163 REAL(kind=
dp),
DIMENSION(nmax) :: a_local = 0.0_dp
164 REAL(kind=
dp),
DIMENSION(nmax) :: c_local = 0.0_dp
165 LOGICAL :: has_nlcc = .false.
166 INTEGER :: n_nlcc = 0
167 REAL(kind=
dp),
DIMENSION(nmax) :: a_nlcc = 0.0_dp
168 REAL(kind=
dp),
DIMENSION(nmax) :: c_nlcc = 0.0_dp
172 INTEGER :: ppot_type = 0
173 LOGICAL :: confinement = .false.
174 INTEGER :: conf_type = 0
175 REAL(
dp) :: acon = 0.0_dp
176 REAL(
dp) :: rcon = 0.0_dp
177 REAL(
dp) :: scon = 0.0_dp
188 REAL(kind=
dp) :: scale_coulomb = 0.0_dp
189 REAL(kind=
dp) :: scale_longrange = 0.0_dp
190 REAL(kind=
dp) :: omega = 0.0_dp
191 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: kernel
192 LOGICAL :: do_gh = .false.
199 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: occ = 0.0_dp
200 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: core = 0.0_dp
201 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: occupation = 0.0_dp
202 INTEGER :: maxl_occ = 0
203 INTEGER,
DIMENSION(0:lmat) :: maxn_occ = 0
204 INTEGER :: maxl_calc = 0
205 INTEGER,
DIMENSION(0:lmat) :: maxn_calc = 0
206 INTEGER :: multiplicity = 0
207 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: occa = 0.0_dp, occb = 0.0_dp
213 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: int => null()
217 INTEGER :: status = 0
218 INTEGER :: ppstat = 0
219 LOGICAL :: eri_coulomb = .false.
220 LOGICAL :: eri_exchange = .false.
221 LOGICAL :: all_nu = .false.
222 INTEGER,
DIMENSION(0:lmat) :: n = 0, nne = 0
223 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: ovlp => null(), kin => null(), core => null(), clsd => null()
224 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: utrans => null(), uptrans => null()
225 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: hnl => null()
226 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: conf => null()
227 TYPE(
eri),
DIMENSION(100) :: ceri =
eri()
228 TYPE(
eri),
DIMENSION(100) :: eeri =
eri()
229 INTEGER :: dkhstat = 0
230 INTEGER :: zorastat = 0
231 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: tzora => null()
232 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: hdkh => null()
238 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: wfn => null(), wfna => null(), wfnb => null()
239 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: pmat => null(), pmata => null(), pmatb => null()
240 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: ener => null(), enera => null(), enerb => null()
241 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: refene => null(), refchg => null(), refnod => null()
242 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: wrefene => null(), wrefchg => null(), wrefnod => null()
243 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: crefene => null(), crefchg => null(), crefnod => null()
244 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: wpsir0 => null(), tpsir0 => null()
245 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: rcmax => null()
246 CHARACTER(LEN=2),
DIMENSION(:, :, :),
POINTER :: reftype => null()
252 INTEGER,
DIMENSION(0:lmat) :: n = 0
253 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: op => null()
259 REAL(kind=
dp),
DIMENSION(:),
POINTER :: op => null()
265 TYPE atom_energy_type
266 REAL(kind=
dp) :: etot = 0.0_dp
267 REAL(kind=
dp) :: eband = 0.0_dp
268 REAL(kind=
dp) :: ekin = 0.0_dp
269 REAL(kind=
dp) :: epot = 0.0_dp
270 REAL(kind=
dp) :: ecore = 0.0_dp
271 REAL(kind=
dp) :: elsd = 0.0_dp
272 REAL(kind=
dp) :: epseudo = 0.0_dp
273 REAL(kind=
dp) :: eploc = 0.0_dp
274 REAL(kind=
dp) :: epnl = 0.0_dp
275 REAL(kind=
dp) :: exc = 0.0_dp
276 REAL(kind=
dp) :: ecoulomb = 0.0_dp
277 REAL(kind=
dp) :: eexchange = 0.0_dp
278 REAL(kind=
dp) :: econfinement = 0.0_dp
279 END TYPE atom_energy_type
284 REAL(kind=
dp) :: damping = 0.0_dp
285 REAL(kind=
dp) :: eps_scf = 0.0_dp
286 REAL(kind=
dp) :: eps_diis = 0.0_dp
287 INTEGER :: max_iter = 0
288 INTEGER :: n_diis = 0
296 LOGICAL :: pp_calc = .false.
298 LOGICAL :: do_zmp = .false., doread = .false., read_vxc = .false., dm = .false.
299 CHARACTER(LEN=default_string_length) :: ext_file =
"", ext_vxc_file =
"", &
300 zmp_restart_file =
""
307 REAL(kind=
dp) :: lambda = 0.0_dp
308 REAL(kind=
dp) :: rho_diff_integral = 0.0_dp
309 REAL(kind=
dp) :: weight = 0.0_dp, zmpgrid_tol = 0.0_dp, zmpvxcgrid_tol = 0.0_dp
316 TYPE(atom_energy_type) :: energy = atom_energy_type()
346 MODULE PROCEDURE read_ecp_potential_file, &
347 read_ecp_potential_files
385 INTEGER,
INTENT(IN) :: zval
386 CHARACTER(LEN=2) :: btyp
388 CHARACTER(LEN=*),
PARAMETER :: routinen =
'init_atom_basis'
389 INTEGER,
PARAMETER :: nua = 40, nup = 16
390 REAL(kind=
dp),
DIMENSION(nua),
PARAMETER :: ugbs = [0.007299_dp, 0.013705_dp, 0.025733_dp, &
391 0.048316_dp, 0.090718_dp, 0.170333_dp, 0.319819_dp, 0.600496_dp, 1.127497_dp, 2.117000_dp,&
392 3.974902_dp, 7.463317_dp, 14.013204_dp, 26.311339_dp, 49.402449_dp, 92.758561_dp, &
393 174.164456_dp, 327.013024_dp, 614.003114_dp, 1152.858743_dp, 2164.619772_dp, &
394 4064.312984_dp, 7631.197056_dp, 14328.416324_dp, 26903.186074_dp, 50513.706789_dp, &
395 94845.070265_dp, 178082.107320_dp, 334368.848683_dp, 627814.487663_dp, 1178791.123851_dp, &
396 2213310.684886_dp, 4155735.557141_dp, 7802853.046713_dp, 14650719.428954_dp, &
397 27508345.793637_dp, 51649961.080194_dp, 96978513.342764_dp, 182087882.613702_dp, &
400 CHARACTER(LEN=default_string_length) :: basis_fn, basis_name
401 INTEGER :: basistype, handle, i, j, k, l, ll, m, &
402 ngp, nl, nr, nu, quadtype
403 INTEGER,
DIMENSION(0:lmat) :: starti
404 INTEGER,
DIMENSION(:),
POINTER :: nqm, num_gto, num_slater, sindex
405 REAL(kind=
dp) :: al, amax, aval, cval, ear, pf, rk
406 REAL(kind=
dp),
DIMENSION(:),
POINTER :: expo
409 CALL timeset(routinen, handle)
416 NULLIFY (basis%am, basis%cm, basis%as, basis%ns, basis%bf, basis%dbf, basis%ddbf)
423 cpabort(
"The number of radial grid points must be greater than zero.")
427 basis%geometrical = .false.
434 SELECT CASE (basistype)
439 IF (num_gto(1) < 1)
THEN
441 IF (btyp ==
"AE")
THEN
443 ELSE IF (btyp ==
"PP")
THEN
450 ALLOCATE (basis%am(nu, 0:
lmat))
452 basis%am(1:nu, i) = ugbs(1:nu)
456 DO i = 1,
SIZE(num_gto)
457 basis%nbas(i - 1) = num_gto(i)
459 basis%nprim = basis%nbas
460 m = maxval(basis%nbas)
461 ALLOCATE (basis%am(m, 0:
lmat))
464 IF (basis%nbas(l) > 0)
THEN
476 cpabort(
"Invalid angular quantum number l found for Gaussian basis set")
478 cpassert(
SIZE(expo) >= basis%nbas(l))
479 DO i = 1, basis%nbas(l)
480 basis%am(i, l) = expo(i)
487 m = maxval(basis%nbas)
488 ALLOCATE (basis%bf(nr, m, 0:
lmat))
489 ALLOCATE (basis%dbf(nr, m, 0:
lmat))
490 ALLOCATE (basis%ddbf(nr, m, 0:
lmat))
495 DO i = 1, basis%nbas(l)
498 rk = basis%grid%rad(k)
499 ear = exp(-al*basis%grid%rad(k)**2)
500 basis%bf(k, i, l) = rk**l*ear
501 basis%dbf(k, i, l) = (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear
502 basis%ddbf(k, i, l) = (real(l*(l - 1),
dp)*rk**(l - 2) - &
503 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))*ear
511 IF (num_gto(1) < 1)
THEN
512 IF (btyp ==
"AE")
THEN
515 ELSE IF (btyp ==
"PP")
THEN
518 ELSE IF (btyp ==
"AA")
THEN
520 amax = cval**(basis%nbas(0) - 1)
521 basis%nbas(0) = nint((log(amax)/log(1.6_dp)))
524 basis%nbas(1) = basis%nbas(0) - 4
525 basis%nbas(2) = basis%nbas(0) - 8
526 basis%nbas(3) = basis%nbas(0) - 12
527 IF (
lmat > 3) basis%nbas(4:
lmat) = 0
528 ELSE IF (btyp ==
"AP")
THEN
531 basis%nbas = nint((log(amax)/log(1.6_dp)))
538 basis%nprim = basis%nbas
541 DO i = 1,
SIZE(num_gto)
542 basis%nbas(i - 1) = num_gto(i)
544 basis%nprim = basis%nbas
548 DO i = 1,
SIZE(sindex)
549 starti(i - 1) = sindex(i)
550 cpassert(sindex(i) >= 0)
555 m = maxval(basis%nbas)
556 ALLOCATE (basis%am(m, 0:
lmat))
559 DO i = 1, basis%nbas(l)
560 ll = i - 1 + starti(l)
561 basis%am(i, l) = aval*cval**(ll)
565 basis%geometrical = .true.
572 m = maxval(basis%nbas)
573 ALLOCATE (basis%bf(nr, m, 0:
lmat))
574 ALLOCATE (basis%dbf(nr, m, 0:
lmat))
575 ALLOCATE (basis%ddbf(nr, m, 0:
lmat))
580 DO i = 1, basis%nbas(l)
583 rk = basis%grid%rad(k)
584 ear = exp(-al*basis%grid%rad(k)**2)
585 basis%bf(k, i, l) = rk**l*ear
586 basis%dbf(k, i, l) = (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear
587 basis%ddbf(k, i, l) = (real(l*(l - 1),
dp)*rk**(l - 2) - &
588 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))*ear
597 CALL read_basis_set(
ptable(zval)%symbol, basis, basis_name, basis_fn, &
602 m = maxval(basis%nbas)
603 ALLOCATE (basis%bf(nr, m, 0:
lmat))
604 ALLOCATE (basis%dbf(nr, m, 0:
lmat))
605 ALLOCATE (basis%ddbf(nr, m, 0:
lmat))
610 DO i = 1, basis%nprim(l)
613 rk = basis%grid%rad(k)
614 ear = exp(-al*basis%grid%rad(k)**2)
615 DO j = 1, basis%nbas(l)
616 basis%bf(k, j, l) = basis%bf(k, j, l) + rk**l*ear*basis%cm(i, j, l)
617 basis%dbf(k, j, l) = basis%dbf(k, j, l) &
618 + (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear*basis%cm(i, j, l)
619 basis%ddbf(k, j, l) = basis%ddbf(k, j, l) + &
620 (real(l*(l - 1),
dp)*rk**(l - 2) - 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))* &
621 ear*basis%cm(i, j, l)
630 IF (num_slater(1) < 1)
THEN
631 cpabort(
"Invalid number (less than 1) Slater-type functions found.")
634 DO i = 1,
SIZE(num_slater)
635 basis%nbas(i - 1) = num_slater(i)
637 basis%nprim = basis%nbas
638 m = maxval(basis%nbas)
639 ALLOCATE (basis%as(m, 0:
lmat), basis%ns(m, 0:
lmat))
643 IF (basis%nbas(l) > 0)
THEN
655 cpabort(
"Invalid angular quantum number l found for Slater basis set")
657 cpassert(
SIZE(expo) >= basis%nbas(l))
658 DO i = 1, basis%nbas(l)
659 basis%as(i, l) = expo(i)
672 cpabort(
"Invalid angular quantum number l found for Slater basis set")
674 cpassert(
SIZE(nqm) >= basis%nbas(l))
675 DO i = 1, basis%nbas(l)
676 basis%ns(i, l) = nqm(i)
683 m = maxval(basis%nbas)
684 ALLOCATE (basis%bf(nr, m, 0:
lmat))
685 ALLOCATE (basis%dbf(nr, m, 0:
lmat))
686 ALLOCATE (basis%ddbf(nr, m, 0:
lmat))
691 DO i = 1, basis%nbas(l)
694 pf = (2._dp*al)**nl*sqrt(2._dp*al/
fac(2*nl))
696 rk = basis%grid%rad(k)
697 ear = rk**(nl - 1)*exp(-al*rk)
698 basis%bf(k, i, l) = pf*ear
699 basis%dbf(k, i, l) = pf*(real(nl - 1,
dp)/rk - al)*ear
700 basis%ddbf(k, i, l) = pf*(real((nl - 2)*(nl - 1),
dp)/rk/rk &
701 - al*real(2*(nl - 1),
dp)/rk + al*al)*ear
707 cpabort(
"Numerical basis set type not yet implemented.")
709 cpabort(
"Unknown basis set type specified. Check the code!")
712 CALL timestop(handle)
723 CHARACTER(LEN=*),
PARAMETER :: routinen =
'init_atom_basis_default_pp'
724 INTEGER,
PARAMETER :: nua = 40, nup = 20
725 REAL(kind=
dp),
DIMENSION(nua),
PARAMETER :: ugbs = [0.007299_dp, 0.013705_dp, 0.025733_dp, &
726 0.048316_dp, 0.090718_dp, 0.170333_dp, 0.319819_dp, 0.600496_dp, 1.127497_dp, 2.117000_dp,&
727 3.974902_dp, 7.463317_dp, 14.013204_dp, 26.311339_dp, 49.402449_dp, 92.758561_dp, &
728 174.164456_dp, 327.013024_dp, 614.003114_dp, 1152.858743_dp, 2164.619772_dp, &
729 4064.312984_dp, 7631.197056_dp, 14328.416324_dp, 26903.186074_dp, 50513.706789_dp, &
730 94845.070265_dp, 178082.107320_dp, 334368.848683_dp, 627814.487663_dp, 1178791.123851_dp, &
731 2213310.684886_dp, 4155735.557141_dp, 7802853.046713_dp, 14650719.428954_dp, &
732 27508345.793637_dp, 51649961.080194_dp, 96978513.342764_dp, 182087882.613702_dp, &
735 INTEGER :: handle, i, k, l, m, ngp, nr, nu, quadtype
736 REAL(kind=
dp) :: al, ear, rk
738 CALL timeset(routinen, handle)
740 NULLIFY (basis%am, basis%cm, basis%as, basis%ns, basis%bf, basis%dbf, basis%ddbf)
749 basis%geometrical = .false.
753 basis%eps_eig = 1.e-12_dp
759 ALLOCATE (basis%am(nu, 0:
lmat))
761 basis%am(1:nu, i) = ugbs(1:nu)
765 m = maxval(basis%nbas)
766 ALLOCATE (basis%bf(nr, m, 0:
lmat))
767 ALLOCATE (basis%dbf(nr, m, 0:
lmat))
768 ALLOCATE (basis%ddbf(nr, m, 0:
lmat))
773 DO i = 1, basis%nbas(l)
776 rk = basis%grid%rad(k)
777 ear = exp(-al*basis%grid%rad(k)**2)
778 basis%bf(k, i, l) = rk**l*ear
779 basis%dbf(k, i, l) = (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear
780 basis%ddbf(k, i, l) = (real(l*(l - 1),
dp)*rk**(l - 2) - &
781 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))*ear
786 CALL timestop(handle)
800 REAL(kind=
dp),
DIMENSION(:),
INTENT(IN) :: r, rab
802 INTEGER :: i, j, k, l, m, n1, n2, n3, ngp, nl, nr, &
804 REAL(kind=
dp) :: al, ear, pf, rk
806 NULLIFY (gbasis%am, gbasis%cm, gbasis%as, gbasis%ns, gbasis%bf, gbasis%dbf, gbasis%ddbf)
809 gbasis%basis_type = basis%basis_type
810 gbasis%nbas(0:
lmat) = basis%nbas(0:
lmat)
811 gbasis%nprim(0:
lmat) = basis%nprim(0:
lmat)
812 IF (
ASSOCIATED(basis%am))
THEN
813 n1 =
SIZE(basis%am, 1)
814 n2 =
SIZE(basis%am, 2)
815 ALLOCATE (gbasis%am(n1, 0:n2 - 1))
818 IF (
ASSOCIATED(basis%cm))
THEN
819 n1 =
SIZE(basis%cm, 1)
820 n2 =
SIZE(basis%cm, 2)
821 n3 =
SIZE(basis%cm, 3)
822 ALLOCATE (gbasis%cm(n1, n2, 0:n3 - 1))
825 IF (
ASSOCIATED(basis%as))
THEN
826 n1 =
SIZE(basis%as, 1)
827 n2 =
SIZE(basis%as, 2)
828 ALLOCATE (gbasis%as(n1, 0:n2 - 1))
831 IF (
ASSOCIATED(basis%ns))
THEN
832 n1 =
SIZE(basis%ns, 1)
833 n2 =
SIZE(basis%ns, 2)
834 ALLOCATE (gbasis%ns(n1, 0:n2 - 1))
837 gbasis%eps_eig = basis%eps_eig
838 gbasis%geometrical = basis%geometrical
839 gbasis%aval = basis%aval
840 gbasis%cval = basis%cval
841 gbasis%start(0:
lmat) = basis%start(0:
lmat)
845 NULLIFY (gbasis%grid)
850 cpabort(
"The number of radial grid points must be greater than zero.")
854 gbasis%grid%rad(:) = r(:)
855 gbasis%grid%rad2(:) = r(:)*r(:)
856 gbasis%grid%wr(:) = rab(:)*gbasis%grid%rad2(:)
860 m = maxval(gbasis%nbas)
861 ALLOCATE (gbasis%bf(nr, m, 0:
lmat))
862 ALLOCATE (gbasis%dbf(nr, m, 0:
lmat))
863 ALLOCATE (gbasis%ddbf(nr, m, 0:
lmat))
868 SELECT CASE (gbasis%basis_type)
871 DO i = 1, gbasis%nbas(l)
874 rk = gbasis%grid%rad(k)
875 ear = exp(-al*gbasis%grid%rad(k)**2)
876 gbasis%bf(k, i, l) = rk**l*ear
877 gbasis%dbf(k, i, l) = (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear
878 gbasis%ddbf(k, i, l) = (real(l*(l - 1),
dp)*rk**(l - 2) - &
879 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))*ear
885 DO i = 1, gbasis%nprim(l)
888 rk = gbasis%grid%rad(k)
889 ear = exp(-al*gbasis%grid%rad(k)**2)
890 DO j = 1, gbasis%nbas(l)
891 gbasis%bf(k, j, l) = gbasis%bf(k, j, l) + rk**l*ear*gbasis%cm(i, j, l)
892 gbasis%dbf(k, j, l) = gbasis%dbf(k, j, l) &
893 + (real(l,
dp)*rk**(l - 1) - 2._dp*al*rk**(l + 1))*ear*gbasis%cm(i, j, l)
894 gbasis%ddbf(k, j, l) = gbasis%ddbf(k, j, l) + &
895 (real(l*(l - 1),
dp)*rk**(l - 2) - 2._dp*al*real(2*l + 1,
dp)*rk**(l) + 4._dp*al*rk**(l + 2))* &
896 ear*gbasis%cm(i, j, l)
903 DO i = 1, gbasis%nbas(l)
906 pf = (2._dp*al)**nl*sqrt(2._dp*al/
fac(2*nl))
908 rk = gbasis%grid%rad(k)
909 ear = rk**(nl - 1)*exp(-al*rk)
910 gbasis%bf(k, i, l) = pf*ear
911 gbasis%dbf(k, i, l) = pf*(real(nl - 1,
dp)/rk - al)*ear
912 gbasis%ddbf(k, i, l) = pf*(real((nl - 2)*(nl - 1),
dp)/rk/rk &
913 - al*real(2*(nl - 1),
dp)/rk + al*al)*ear
919 cpabort(
"Numerical basis set type not yet implemented.")
921 cpabort(
"Unknown basis set type specified. Check the code!")
933 IF (
ASSOCIATED(basis%am))
THEN
934 DEALLOCATE (basis%am)
936 IF (
ASSOCIATED(basis%cm))
THEN
937 DEALLOCATE (basis%cm)
939 IF (
ASSOCIATED(basis%as))
THEN
940 DEALLOCATE (basis%as)
942 IF (
ASSOCIATED(basis%ns))
THEN
943 DEALLOCATE (basis%ns)
945 IF (
ASSOCIATED(basis%bf))
THEN
946 DEALLOCATE (basis%bf)
948 IF (
ASSOCIATED(basis%dbf))
THEN
949 DEALLOCATE (basis%dbf)
951 IF (
ASSOCIATED(basis%ddbf))
THEN
952 DEALLOCATE (basis%ddbf)
967 cpassert(.NOT.
ASSOCIATED(
atom))
971 NULLIFY (
atom%zmp_section)
972 NULLIFY (
atom%xc_section)
974 atom%do_zmp = .false.
975 atom%doread = .false.
976 atom%read_vxc = .false.
978 atom%hfx_pot%scale_coulomb = 0.0_dp
979 atom%hfx_pot%scale_longrange = 0.0_dp
980 atom%hfx_pot%omega = 0.0_dp
991 cpassert(
ASSOCIATED(
atom))
994 NULLIFY (
atom%integrals)
995 IF (
ASSOCIATED(
atom%state))
THEN
996 DEALLOCATE (
atom%state)
998 IF (
ASSOCIATED(
atom%orbitals))
THEN
1028 SUBROUTINE set_atom(atom, basis, state, integrals, orbitals, potential, zcore, pp_calc, do_zmp, doread, &
1029 read_vxc, method_type, relativistic, coulomb_integral_type, exchange_integral_type, fmat)
1036 INTEGER,
INTENT(IN),
OPTIONAL :: zcore
1037 LOGICAL,
INTENT(IN),
OPTIONAL :: pp_calc, do_zmp, doread, read_vxc
1038 INTEGER,
INTENT(IN),
OPTIONAL :: method_type, relativistic, &
1039 coulomb_integral_type, &
1040 exchange_integral_type
1043 cpassert(
ASSOCIATED(
atom))
1045 IF (
PRESENT(basis))
atom%basis => basis
1046 IF (
PRESENT(state))
atom%state => state
1047 IF (
PRESENT(integrals))
atom%integrals => integrals
1048 IF (
PRESENT(orbitals))
atom%orbitals => orbitals
1049 IF (
PRESENT(potential))
atom%potential => potential
1050 IF (
PRESENT(zcore))
atom%zcore = zcore
1051 IF (
PRESENT(pp_calc))
atom%pp_calc = pp_calc
1053 IF (
PRESENT(do_zmp))
atom%do_zmp = do_zmp
1054 IF (
PRESENT(doread))
atom%doread = doread
1055 IF (
PRESENT(read_vxc))
atom%read_vxc = read_vxc
1057 IF (
PRESENT(method_type))
atom%method_type = method_type
1058 IF (
PRESENT(relativistic))
atom%relativistic = relativistic
1059 IF (
PRESENT(coulomb_integral_type))
atom%coulomb_integral_type = coulomb_integral_type
1060 IF (
PRESENT(exchange_integral_type))
atom%exchange_integral_type = exchange_integral_type
1062 IF (
PRESENT(fmat))
THEN
1077 INTEGER,
INTENT(IN) :: mbas, mo
1079 cpassert(.NOT.
ASSOCIATED(orbs))
1083 ALLOCATE (orbs%wfn(mbas, mo, 0:
lmat), orbs%wfna(mbas, mo, 0:
lmat), orbs%wfnb(mbas, mo, 0:
lmat))
1088 ALLOCATE (orbs%pmat(mbas, mbas, 0:
lmat), orbs%pmata(mbas, mbas, 0:
lmat), orbs%pmatb(mbas, mbas, 0:
lmat))
1093 ALLOCATE (orbs%ener(mo, 0:
lmat), orbs%enera(mo, 0:
lmat), orbs%enerb(mo, 0:
lmat))
1098 ALLOCATE (orbs%refene(mo, 0:
lmat, 2), orbs%refchg(mo, 0:
lmat, 2), orbs%refnod(mo, 0:
lmat, 2))
1102 ALLOCATE (orbs%wrefene(mo, 0:
lmat, 2), orbs%wrefchg(mo, 0:
lmat, 2), orbs%wrefnod(mo, 0:
lmat, 2))
1103 orbs%wrefene = 0._dp
1104 orbs%wrefchg = 0._dp
1105 orbs%wrefnod = 0._dp
1106 ALLOCATE (orbs%crefene(mo, 0:
lmat, 2), orbs%crefchg(mo, 0:
lmat, 2), orbs%crefnod(mo, 0:
lmat, 2))
1107 orbs%crefene = 0._dp
1108 orbs%crefchg = 0._dp
1109 orbs%crefnod = 0._dp
1110 ALLOCATE (orbs%rcmax(mo, 0:
lmat, 2))
1112 ALLOCATE (orbs%wpsir0(mo, 2), orbs%tpsir0(mo, 2))
1115 ALLOCATE (orbs%reftype(mo, 0:
lmat, 2))
1127 cpassert(
ASSOCIATED(orbs))
1129 IF (
ASSOCIATED(orbs%wfn))
THEN
1130 DEALLOCATE (orbs%wfn, orbs%wfna, orbs%wfnb)
1132 IF (
ASSOCIATED(orbs%pmat))
THEN
1133 DEALLOCATE (orbs%pmat, orbs%pmata, orbs%pmatb)
1135 IF (
ASSOCIATED(orbs%ener))
THEN
1136 DEALLOCATE (orbs%ener, orbs%enera, orbs%enerb)
1138 IF (
ASSOCIATED(orbs%refene))
THEN
1139 DEALLOCATE (orbs%refene)
1141 IF (
ASSOCIATED(orbs%refchg))
THEN
1142 DEALLOCATE (orbs%refchg)
1144 IF (
ASSOCIATED(orbs%refnod))
THEN
1145 DEALLOCATE (orbs%refnod)
1147 IF (
ASSOCIATED(orbs%wrefene))
THEN
1148 DEALLOCATE (orbs%wrefene)
1150 IF (
ASSOCIATED(orbs%wrefchg))
THEN
1151 DEALLOCATE (orbs%wrefchg)
1153 IF (
ASSOCIATED(orbs%wrefnod))
THEN
1154 DEALLOCATE (orbs%wrefnod)
1156 IF (
ASSOCIATED(orbs%crefene))
THEN
1157 DEALLOCATE (orbs%crefene)
1159 IF (
ASSOCIATED(orbs%crefchg))
THEN
1160 DEALLOCATE (orbs%crefchg)
1162 IF (
ASSOCIATED(orbs%crefnod))
THEN
1163 DEALLOCATE (orbs%crefnod)
1165 IF (
ASSOCIATED(orbs%rcmax))
THEN
1166 DEALLOCATE (orbs%rcmax)
1168 IF (
ASSOCIATED(orbs%wpsir0))
THEN
1169 DEALLOCATE (orbs%wpsir0)
1171 IF (
ASSOCIATED(orbs%tpsir0))
THEN
1172 DEALLOCATE (orbs%tpsir0)
1174 IF (
ASSOCIATED(orbs%reftype))
THEN
1175 DEALLOCATE (orbs%reftype)
1191 REAL(kind=
dp),
INTENT(OUT) :: hf_frac
1192 LOGICAL,
INTENT(OUT) :: do_hfx
1195 INTEGER,
INTENT(IN) :: extype
1197 INTEGER :: i, j, nr, nu, pot_type
1198 REAL(kind=
dp) :: scale_coulomb, scale_longrange
1199 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: abscissa, weights
1203 IF (
ASSOCIATED(
atom%xc_section))
THEN
1204 xc_section =>
atom%xc_section
1209 atom%hfx_pot%scale_longrange = 0.0_dp
1210 atom%hfx_pot%scale_coulomb = 1.0_dp
1223 SELECT CASE (pot_type)
1225 cpwarn(
"Potential not implemented, use Coulomb instead!")
1227 atom%hfx_pot%scale_longrange = 0.0_dp
1228 atom%hfx_pot%scale_coulomb = scale_coulomb
1230 atom%hfx_pot%scale_coulomb = 0.0_dp
1231 atom%hfx_pot%scale_longrange = scale_longrange
1233 atom%hfx_pot%scale_coulomb = 1.0_dp
1234 atom%hfx_pot%scale_longrange = -1.0_dp
1236 atom%hfx_pot%scale_coulomb = scale_coulomb
1237 atom%hfx_pot%scale_longrange = scale_longrange
1243 cpabort(
"Only numerical and semi-analytic lrHF exchange available!")
1246 IF (
atom%hfx_pot%scale_longrange /= 0.0_dp .AND. extype ==
do_numeric .AND. .NOT.
ALLOCATED(
atom%hfx_pot%kernel))
THEN
1249 IF (
atom%hfx_pot%do_gh)
THEN
1253 ALLOCATE (weights(
atom%hfx_pot%nr_gh), abscissa(
atom%hfx_pot%nr_gh))
1254 CALL get_gauss_hermite_weights(abscissa, weights,
atom%hfx_pot%nr_gh)
1256 nr =
atom%basis%grid%nr
1257 ALLOCATE (
atom%hfx_pot%kernel(nr,
atom%hfx_pot%nr_gh, 0:
atom%state%maxl_calc +
atom%state%maxl_occ))
1258 atom%hfx_pot%kernel = 0.0_dp
1259 DO nu = 0,
atom%state%maxl_calc +
atom%state%maxl_occ
1260 DO i = 1,
atom%hfx_pot%nr_gh
1263 *abscissa(i)*
atom%basis%grid%rad(j), nu)*sqrt(weights(i))
1271 nr =
atom%basis%grid%nr
1272 ALLOCATE (
atom%hfx_pot%kernel(nr, nr, 0:
atom%state%maxl_calc +
atom%state%maxl_occ))
1273 atom%hfx_pot%kernel = 0.0_dp
1274 DO nu = 0,
atom%state%maxl_calc +
atom%state%maxl_occ
1278 *
atom%basis%grid%rad(i)*
atom%basis%grid%rad(j), nu)
1285 NULLIFY (xc_section)
1297 SUBROUTINE get_gauss_hermite_weights(abscissa, weights, nn)
1298 REAL(kind=
dp),
DIMENSION(:),
INTENT(OUT) :: abscissa, weights
1299 INTEGER,
INTENT(IN) :: nn
1301 INTEGER :: counter, ii, info, liwork, lwork
1302 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: iwork
1303 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: diag, subdiag, work
1304 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: eigenvec
1308 ALLOCATE (work(1), iwork(1), diag(2*nn), subdiag(2*nn - 1), eigenvec(2*nn, 2*nn))
1313 subdiag(ii) = sqrt(real(ii, kind=
dp)/2.0_dp)
1317 CALL dstevd(
'V', 2*nn, diag, subdiag, eigenvec, 2*nn, work, lwork, iwork, liwork, info)
1320 cpabort(
'Finding size of working matrices failed!')
1324 lwork = int(work(1))
1326 DEALLOCATE (work, iwork)
1327 ALLOCATE (work(lwork), iwork(liwork))
1330 CALL dstevd(
'V', 2*nn, diag, subdiag, eigenvec, 2*nn, work, lwork, iwork, liwork, info)
1333 cpabort(
'Eigenvalue decomposition failed!')
1336 DEALLOCATE (work, iwork, subdiag)
1342 IF (diag(ii) > 0.0_dp)
THEN
1343 counter = counter + 1
1344 abscissa(counter) = diag(ii)
1345 weights(counter) =
rootpi*eigenvec(1, ii)**2
1348 IF (counter /= nn)
THEN
1349 cpabort(
'Have not found enough or too many zeros!')
1352 END SUBROUTINE get_gauss_hermite_weights
1362 INTEGER,
DIMENSION(0:lmat),
INTENT(IN) :: n
1363 INTEGER,
INTENT(IN),
OPTIONAL :: lmax
1368 IF (
PRESENT(lmax))
THEN
1374 cpassert(.NOT.
ASSOCIATED(opmat))
1379 ALLOCATE (opmat%op(m, m, 0:lm))
1391 cpassert(
ASSOCIATED(opmat))
1394 DEALLOCATE (opmat%op)
1411 cpassert(.NOT.
ASSOCIATED(opgrid))
1419 ALLOCATE (opgrid%op(nr))
1431 cpassert(
ASSOCIATED(opgrid))
1433 NULLIFY (opgrid%grid)
1434 DEALLOCATE (opgrid%op)
1450 INTEGER,
INTENT(IN) :: zval
1451 REAL(
dp),
INTENT(OUT) :: cval, aval
1452 INTEGER,
DIMENSION(0:lmat),
INTENT(OUT) :: ngto, ival
1465 cval = 2.14774520_dp
1466 aval = 0.04850670_dp
1469 cval = 2.08932430_dp
1470 aval = 0.02031060_dp
1473 cval = 2.09753060_dp
1474 aval = 0.03207070_dp
1477 cval = 2.10343410_dp
1478 aval = 0.03591970_dp
1482 cval = 2.10662820_dp
1483 aval = 0.05292410_dp
1487 cval = 2.13743840_dp
1488 aval = 0.06291970_dp
1492 cval = 2.08687310_dp
1493 aval = 0.08350860_dp
1497 cval = 2.12318180_dp
1498 aval = 0.09899170_dp
1502 cval = 2.13164810_dp
1503 aval = 0.11485350_dp
1507 cval = 2.11413310_dp
1508 aval = 0.00922630_dp
1513 cval = 2.12183620_dp
1514 aval = 0.01215850_dp
1519 cval = 2.06073230_dp
1520 aval = 0.01449350_dp
1525 cval = 2.08563660_dp
1526 aval = 0.01861460_dp
1531 cval = 2.04879270_dp
1532 aval = 0.02147790_dp
1537 cval = 2.06216660_dp
1538 aval = 0.01978920_dp
1543 cval = 2.04628670_dp
1544 aval = 0.02451470_dp
1549 cval = 2.08675200_dp
1550 aval = 0.02635040_dp
1555 cval = 2.02715220_dp
1556 aval = 0.01822040_dp
1561 cval = 2.01465650_dp
1562 aval = 0.01646570_dp
1567 cval = 2.01605240_dp
1568 aval = 0.01254190_dp
1574 cval = 2.01800000_dp
1575 aval = 0.01195490_dp
1582 cval = 1.98803560_dp
1583 aval = 0.02492140_dp
1590 cval = 1.98984000_dp
1591 aval = 0.02568400_dp
1598 cval = 2.01694380_dp
1599 aval = 0.02664480_dp
1606 cval = 2.01824090_dp
1607 aval = 0.01355000_dp
1614 cval = 1.98359400_dp
1615 aval = 0.01702210_dp
1622 cval = 1.96797340_dp
1623 aval = 0.02163180_dp
1630 cval = 1.98955180_dp
1631 aval = 0.02304480_dp
1638 cval = 1.98074320_dp
1639 aval = 0.02754320_dp
1646 cval = 2.00551070_dp
1647 aval = 0.02005530_dp
1654 cval = 2.00000030_dp
1655 aval = 0.02003000_dp
1662 cval = 2.00609100_dp
1663 aval = 0.02055620_dp
1670 cval = 2.00701000_dp
1671 aval = 0.02230400_dp
1678 cval = 2.01508710_dp
1679 aval = 0.02685790_dp
1686 cval = 2.01960430_dp
1687 aval = 0.02960430_dp
1694 cval = 2.00031000_dp
1695 aval = 0.00768400_dp
1703 cval = 1.99563960_dp
1704 aval = 0.01401940_dp
1711 cval = 1.98971210_dp
1712 aval = 0.01558470_dp
1718 cval = 1.97976190_dp
1719 aval = 0.01705520_dp
1725 cval = 1.97989290_dp
1726 aval = 0.01527040_dp
1732 cval = 1.97909240_dp
1733 aval = 0.01879720_dp
1739 cval = 1.98508430_dp
1740 aval = 0.01497550_dp
1747 cval = 1.98515010_dp
1748 aval = 0.01856670_dp
1755 cval = 1.98502970_dp
1756 aval = 0.01487000_dp
1763 cval = 1.97672850_dp
1764 aval = 0.01762500_dp
1772 cval = 1.97862730_dp
1773 aval = 0.01863310_dp
1780 cval = 1.97990020_dp
1781 aval = 0.01347150_dp
1788 cval = 1.97979410_dp
1789 aval = 0.00890265_dp
1796 cval = 1.98001000_dp
1797 aval = 0.00895215_dp
1804 cval = 1.97979980_dp
1805 aval = 0.01490290_dp
1812 cval = 1.98009310_dp
1813 aval = 0.01490390_dp
1820 cval = 1.97794750_dp
1821 aval = 0.01425880_dp
1829 cval = 1.97784450_dp
1830 aval = 0.01430130_dp
1838 cval = 1.97784450_dp
1839 aval = 0.00499318_dp
1847 cval = 1.97764820_dp
1848 aval = 0.00500392_dp
1856 cval = 1.97765150_dp
1857 aval = 0.00557083_dp
1865 cval = 1.97768750_dp
1866 aval = 0.00547531_dp
1876 cval = 1.96986600_dp
1877 aval = 0.00813143_dp
1887 cval = 1.97765720_dp
1888 aval = 0.00489201_dp
1898 cval = 1.97768120_dp
1899 aval = 0.00499000_dp
1909 cval = 1.97745700_dp
1910 aval = 0.00615587_dp
1920 cval = 1.97570240_dp
1921 aval = 0.00769959_dp
1931 cval = 1.97629350_dp
1932 aval = 0.00706610_dp
1942 cval = 1.96900000_dp
1943 aval = 0.01019150_dp
1953 cval = 1.97350000_dp
1954 aval = 0.01334320_dp
1964 cval = 1.97493000_dp
1965 aval = 0.01331360_dp
1974 cval = 1.97597670_dp
1975 aval = 0.01434040_dp
1985 cval = 1.97809240_dp
1986 aval = 0.01529430_dp
1996 cval = 1.97644360_dp
1997 aval = 0.01312770_dp
2007 cval = 1.96998000_dp
2008 aval = 0.01745150_dp
2018 cval = 1.97223830_dp
2019 aval = 0.01639750_dp
2029 cval = 1.97462110_dp
2030 aval = 0.01603680_dp
2040 cval = 1.97756000_dp
2041 aval = 0.02030570_dp
2051 cval = 1.97645760_dp
2052 aval = 0.02057180_dp
2062 cval = 1.97725820_dp
2063 aval = 0.02058210_dp
2073 cval = 1.97749380_dp
2074 aval = 0.02219380_dp
2084 cval = 1.97946280_dp
2085 aval = 0.02216280_dp
2095 cval = 1.97852130_dp
2096 aval = 0.02168500_dp
2106 cval = 1.98045190_dp
2107 aval = 0.02177860_dp
2117 cval = 1.97000000_dp
2118 aval = 0.02275000_dp
2128 cval = 1.97713580_dp
2129 aval = 0.02317030_dp
2139 cval = 1.97537880_dp
2140 aval = 0.02672860_dp
2150 cval = 1.97545360_dp
2151 aval = 0.02745360_dp
2161 cval = 1.97338370_dp
2162 aval = 0.02616310_dp
2172 cval = 1.97294240_dp
2173 aval = 0.02429220_dp
2183 cval = 1.98000000_dp
2184 aval = 0.01400000_dp
2194 cpabort(
"No geometrical basis set data are available for the selected atom number.")
2207 SUBROUTINE read_basis_set(element_symbol, basis, basis_set_name, basis_set_file, &
2210 CHARACTER(LEN=*),
INTENT(IN) :: element_symbol
2212 CHARACTER(LEN=*),
INTENT(IN) :: basis_set_name, basis_set_file
2215 INTEGER,
PARAMETER :: maxpri = 40, maxset = 20
2217 CHARACTER(len=20*default_string_length) :: line_att
2218 CHARACTER(LEN=240) :: line
2219 CHARACTER(LEN=242) :: line2
2220 CHARACTER(LEN=LEN(basis_set_name)) :: bsname
2221 CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
2222 CHARACTER(LEN=LEN(element_symbol)) :: symbol
2223 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
2224 INTEGER :: i, ii, ipgf, irep, iset, ishell, j, k, &
2225 lshell, nj, nmin, ns, nset, strlen1, &
2227 INTEGER,
DIMENSION(maxpri, maxset) :: l
2228 INTEGER,
DIMENSION(maxset) :: lmax, lmin, n, npgf, nshell
2229 LOGICAL :: found, is_ok, match, read_from_input
2230 REAL(
dp) :: expzet, gcca, prefac, zeta
2231 REAL(
dp),
DIMENSION(maxpri, maxpri, maxset) :: gcc
2232 REAL(
dp),
DIMENSION(maxpri, maxset) :: zet
2236 bsname = basis_set_name
2237 symbol = element_symbol
2249 read_from_input = .false.
2251 IF (read_from_input)
THEN
2258 CALL val_get(val, c_val=line_att)
2259 READ (line_att, *) nset
2260 cpassert(nset <= maxset)
2264 CALL val_get(val, c_val=line_att)
2265 READ (line_att, *) n(iset)
2267 READ (line_att, *) lmin(iset)
2269 READ (line_att, *) lmax(iset)
2271 READ (line_att, *) npgf(iset)
2273 cpassert(npgf(iset) <= maxpri)
2275 DO lshell = lmin(iset), lmax(iset)
2276 nmin = n(iset) + lshell - lmin(iset)
2277 READ (line_att, *) ishell
2279 nshell(iset) = nshell(iset) + ishell
2281 l(nshell(iset) - ishell + i, iset) = lshell
2284 cpassert(len_trim(line_att) == 0)
2285 DO ipgf = 1, npgf(iset)
2288 CALL val_get(val, c_val=line_att)
2289 READ (line_att, *) zet(ipgf, iset), (gcc(ipgf, ishell, iset), ishell=1, nshell(iset))
2306 line2 =
" "//line//
" "
2307 symbol2 =
" "//trim(symbol)//
" "
2308 bsname2 =
" "//trim(bsname)//
" "
2309 strlen1 = len_trim(symbol2) + 1
2310 strlen2 = len_trim(bsname2) + 1
2312 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
2313 (index(line2, bsname2(:strlen2)) > 0)) match = .true.
2318 cpassert(nset <= maxset)
2324 cpassert(npgf(iset) <= maxpri)
2326 DO lshell = lmin(iset), lmax(iset)
2327 nmin = n(iset) + lshell - lmin(iset)
2329 nshell(iset) = nshell(iset) + ishell
2331 l(nshell(iset) - ishell + i, iset) = lshell
2334 DO ipgf = 1, npgf(iset)
2336 DO ishell = 1, nshell(iset)
2347 cpabort(
"End of file reached and the requested basis set was not found.")
2360 DO j = lmin(i), min(lmax(i),
lmat)
2361 basis%nprim(j) = basis%nprim(j) + npgf(i)
2365 IF (k <=
lmat) basis%nbas(k) = basis%nbas(k) + 1
2369 nj = maxval(basis%nprim)
2370 ns = maxval(basis%nbas)
2371 ALLOCATE (basis%am(nj, 0:
lmat))
2373 ALLOCATE (basis%cm(nj, ns, 0:
lmat))
2380 IF (j >= lmin(i) .AND. j <= lmax(i))
THEN
2381 DO ipgf = 1, npgf(i)
2382 basis%am(nj + ipgf, j) = zet(ipgf, i)
2384 DO ii = 1, nshell(i)
2385 IF (l(ii, i) == j)
THEN
2387 DO ipgf = 1, npgf(i)
2388 basis%cm(nj + ipgf, ns, j) = gcc(ipgf, ii, i)
2399 expzet = 0.25_dp*real(2*j + 3,
dp)
2400 prefac = sqrt(
rootpi/2._dp**(j + 2)*
dfac(2*j + 1))
2401 DO ipgf = 1, basis%nprim(j)
2402 DO ii = 1, basis%nbas(j)
2403 gcca = basis%cm(ipgf, ii, j)
2404 zeta = 2._dp*basis%am(ipgf, j)
2405 basis%cm(ipgf, ii, j) = zeta**expzet*gcca/prefac
2410 END SUBROUTINE read_basis_set
2419 TYPE(section_vals_type),
POINTER :: opt_section
2421 INTEGER :: miter, ndiis
2422 REAL(kind=dp) :: damp, eps_diis, eps_scf
2424 CALL section_vals_val_get(opt_section,
"MAX_ITER", i_val=miter)
2425 CALL section_vals_val_get(opt_section,
"EPS_SCF", r_val=eps_scf)
2426 CALL section_vals_val_get(opt_section,
"N_DIIS", i_val=ndiis)
2427 CALL section_vals_val_get(opt_section,
"EPS_DIIS", r_val=eps_diis)
2428 CALL section_vals_val_get(opt_section,
"DAMPING", r_val=damp)
2430 optimization%max_iter = miter
2431 optimization%eps_scf = eps_scf
2432 optimization%n_diis = ndiis
2433 optimization%eps_diis = eps_diis
2434 optimization%damping = damp
2445 TYPE(section_vals_type),
POINTER :: potential_section
2446 INTEGER,
INTENT(IN) :: zval
2448 CHARACTER(LEN=default_string_length) :: pseudo_fn, pseudo_name
2450 REAL(dp),
DIMENSION(:),
POINTER :: convals
2451 TYPE(section_vals_type),
POINTER :: ecp_potential_section, &
2452 gth_potential_section
2455 CALL section_vals_val_get(potential_section,
"PSEUDO_TYPE", i_val=potential%ppot_type)
2457 SELECT CASE (potential%ppot_type)
2459 CALL section_vals_val_get(potential_section,
"POTENTIAL_FILE_NAME", c_val=pseudo_fn)
2460 CALL section_vals_val_get(potential_section,
"POTENTIAL_NAME", c_val=pseudo_name)
2461 gth_potential_section => section_vals_get_subs_vals(potential_section,
"GTH_POTENTIAL")
2462 CALL read_gth_potential(ptable(zval)%symbol, potential%gth_pot, &
2463 pseudo_name, pseudo_fn, gth_potential_section)
2465 CALL section_vals_val_get(potential_section,
"POTENTIAL_FILE_NAME", c_val=pseudo_fn)
2466 CALL section_vals_val_get(potential_section,
"POTENTIAL_NAME", c_val=pseudo_name)
2467 ecp_potential_section => section_vals_get_subs_vals(potential_section,
"ECP")
2469 pseudo_name, pseudo_fn, ecp_potential_section)
2471 CALL section_vals_val_get(potential_section,
"POTENTIAL_FILE_NAME", c_val=pseudo_fn)
2472 CALL section_vals_val_get(potential_section,
"POTENTIAL_NAME", c_val=pseudo_name)
2473 CALL atom_read_upf(potential%upf_pot, pseudo_fn)
2474 potential%upf_pot%pname = pseudo_name
2476 cpabort(
"Pseudopotential type SGP is not implemented.")
2480 cpabort(
"Invalid pseudopotential type selected. Check the code!")
2483 potential%ppot_type = no_pseudo
2488 CALL section_vals_val_get(potential_section,
"CONFINEMENT_TYPE", i_val=ic)
2489 potential%conf_type = ic
2490 IF (potential%conf_type == no_conf)
THEN
2491 potential%acon = 0.0_dp
2492 potential%rcon = 4.0_dp
2493 potential%scon = 2.0_dp
2494 potential%confinement = .false.
2495 ELSE IF (potential%conf_type == poly_conf)
THEN
2496 CALL section_vals_val_get(potential_section,
"CONFINEMENT", r_vals=convals)
2497 IF (
SIZE(convals) >= 1)
THEN
2498 IF (convals(1) > 0.0_dp)
THEN
2499 potential%confinement = .true.
2500 potential%acon = convals(1)
2501 IF (
SIZE(convals) >= 2)
THEN
2502 potential%rcon = convals(2)
2504 potential%rcon = 4.0_dp
2506 IF (
SIZE(convals) >= 3)
THEN
2507 potential%scon = convals(3)
2509 potential%scon = 2.0_dp
2512 potential%confinement = .false.
2515 potential%confinement = .false.
2517 ELSE IF (potential%conf_type == barrier_conf)
THEN
2518 potential%acon = 200.0_dp
2519 potential%rcon = 4.0_dp
2520 potential%scon = 12.0_dp
2521 potential%confinement = .true.
2522 CALL section_vals_val_get(potential_section,
"CONFINEMENT", r_vals=convals)
2523 IF (
SIZE(convals) >= 1)
THEN
2524 IF (convals(1) > 0.0_dp)
THEN
2525 potential%acon = convals(1)
2526 IF (
SIZE(convals) >= 2)
THEN
2527 potential%rcon = convals(2)
2529 IF (
SIZE(convals) >= 3)
THEN
2530 potential%scon = convals(3)
2533 potential%confinement = .false.
2546 potential%confinement = .false.
2548 CALL atom_release_upf(potential%upf_pot)
2559 SUBROUTINE read_gth_potential(element_symbol, potential, pseudo_name, pseudo_file, &
2562 CHARACTER(LEN=*),
INTENT(IN) :: element_symbol
2564 CHARACTER(LEN=*),
INTENT(IN) :: pseudo_name, pseudo_file
2565 TYPE(section_vals_type),
POINTER :: potential_section
2567 CHARACTER(LEN=240) :: line
2568 CHARACTER(LEN=242) :: line2
2569 CHARACTER(len=5*default_string_length) :: line_att
2570 CHARACTER(LEN=LEN(element_symbol)) :: symbol
2571 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
2572 CHARACTER(LEN=LEN(pseudo_name)) :: apname
2573 CHARACTER(LEN=LEN(pseudo_name)+2) :: apname2
2574 INTEGER :: i, ic, ipot, j, l, nlmax, strlen1, &
2576 INTEGER,
DIMENSION(0:lmat) :: elec_conf
2577 LOGICAL :: found, is_ok, match, read_from_input
2578 TYPE(cp_sll_val_type),
POINTER ::
list
2579 TYPE(val_type),
POINTER :: val
2583 apname = pseudo_name
2584 symbol = element_symbol
2586 potential%symbol = symbol
2587 potential%pname = apname
2589 potential%rc = 0._dp
2591 potential%cl = 0._dp
2593 potential%rcnl = 0._dp
2594 potential%hnl = 0._dp
2595 potential%soc = .false.
2596 potential%knl = 0._dp
2598 potential%lpotextended = .false.
2599 potential%lsdpot = .false.
2600 potential%nlcc = .false.
2601 potential%nexp_lpot = 0
2602 potential%nexp_lsd = 0
2603 potential%nexp_nlcc = 0
2605 read_from_input = .false.
2606 CALL section_vals_get(potential_section, explicit=read_from_input)
2607 IF (read_from_input)
THEN
2608 CALL section_vals_list_get(potential_section,
"_DEFAULT_KEYWORD_",
list=
list)
2609 CALL uppercase(symbol)
2610 CALL uppercase(apname)
2613 is_ok = cp_sll_val_next(
list, val)
2615 CALL val_get(val, c_val=line_att)
2616 READ (line_att, *) elec_conf(l)
2617 CALL remove_word(line_att)
2618 DO WHILE (len_trim(line_att) /= 0)
2620 READ (line_att, *) elec_conf(l)
2621 CALL remove_word(line_att)
2623 potential%econf(0:
lmat) = elec_conf(0:
lmat)
2624 potential%zion = real(sum(elec_conf), dp)
2626 is_ok = cp_sll_val_next(
list, val)
2628 CALL val_get(val, c_val=line_att)
2629 READ (line_att, *) potential%rc
2630 CALL remove_word(line_att)
2632 READ (line_att, *) potential%ncl
2633 CALL remove_word(line_att)
2634 DO i = 1, potential%ncl
2635 READ (line_att, *) potential%cl(i)
2636 CALL remove_word(line_att)
2640 is_ok = cp_sll_val_next(
list, val)
2642 CALL val_get(val, c_val=line_att)
2643 IF (index(line_att,
"LPOT") /= 0)
THEN
2644 potential%lpotextended = .true.
2645 CALL remove_word(line_att)
2646 READ (line_att, *) potential%nexp_lpot
2647 DO ipot = 1, potential%nexp_lpot
2648 is_ok = cp_sll_val_next(
list, val)
2650 CALL val_get(val, c_val=line_att)
2651 READ (line_att, *) potential%alpha_lpot(ipot)
2652 CALL remove_word(line_att)
2653 READ (line_att, *) potential%nct_lpot(ipot)
2654 CALL remove_word(line_att)
2655 DO ic = 1, potential%nct_lpot(ipot)
2656 READ (line_att, *) potential%cval_lpot(ic, ipot)
2657 CALL remove_word(line_att)
2660 ELSE IF (index(line_att,
"NLCC") /= 0)
THEN
2661 potential%nlcc = .true.
2662 CALL remove_word(line_att)
2663 READ (line_att, *) potential%nexp_nlcc
2664 DO ipot = 1, potential%nexp_nlcc
2665 is_ok = cp_sll_val_next(
list, val)
2667 CALL val_get(val, c_val=line_att)
2668 READ (line_att, *) potential%alpha_nlcc(ipot)
2669 CALL remove_word(line_att)
2670 READ (line_att, *) potential%nct_nlcc(ipot)
2671 CALL remove_word(line_att)
2672 DO ic = 1, potential%nct_nlcc(ipot)
2673 READ (line_att, *) potential%cval_nlcc(ic, ipot)
2675 potential%cval_nlcc(ic, ipot) = potential%cval_nlcc(ic, ipot)/(4.0_dp*pi)
2676 CALL remove_word(line_att)
2679 ELSE IF (index(line_att,
"LSD") /= 0)
THEN
2680 potential%lsdpot = .true.
2681 CALL remove_word(line_att)
2682 READ (line_att, *) potential%nexp_lsd
2683 DO ipot = 1, potential%nexp_lsd
2684 is_ok = cp_sll_val_next(
list, val)
2686 CALL val_get(val, c_val=line_att)
2687 READ (line_att, *) potential%alpha_lsd(ipot)
2688 CALL remove_word(line_att)
2689 READ (line_att, *) potential%nct_lsd(ipot)
2690 CALL remove_word(line_att)
2691 DO ic = 1, potential%nct_lsd(ipot)
2692 READ (line_att, *) potential%cval_lsd(ic, ipot)
2693 CALL remove_word(line_att)
2701 READ (line_att, *) nlmax
2702 CALL remove_word(line_att)
2703 IF (index(line_att,
"SOC") /= 0) potential%soc = .true.
2707 is_ok = cp_sll_val_next(
list, val)
2709 CALL val_get(val, c_val=line_att)
2710 READ (line_att, *) potential%rcnl(l)
2711 CALL remove_word(line_att)
2712 READ (line_att, *) potential%nl(l)
2713 CALL remove_word(line_att)
2714 DO i = 1, potential%nl(l)
2716 READ (line_att, *) potential%hnl(1, 1, l)
2717 CALL remove_word(line_att)
2719 cpassert(len_trim(line_att) == 0)
2720 is_ok = cp_sll_val_next(
list, val)
2722 CALL val_get(val, c_val=line_att)
2723 READ (line_att, *) potential%hnl(i, i, l)
2724 CALL remove_word(line_att)
2726 DO j = i + 1, potential%nl(l)
2727 READ (line_att, *) potential%hnl(i, j, l)
2728 potential%hnl(j, i, l) = potential%hnl(i, j, l)
2729 CALL remove_word(line_att)
2732 IF (potential%soc .AND. l /= 0)
THEN
2733 is_ok = cp_sll_val_next(
list, val)
2735 CALL val_get(val, c_val=line_att)
2736 DO i = 1, potential%nl(l)
2738 READ (line_att, *) potential%knl(1, 1, l)
2739 CALL remove_word(line_att)
2741 cpassert(len_trim(line_att) == 0)
2742 is_ok = cp_sll_val_next(
list, val)
2744 CALL val_get(val, c_val=line_att)
2745 READ (line_att, *) potential%knl(i, i, l)
2746 CALL remove_word(line_att)
2748 DO j = i + 1, potential%nl(l)
2749 READ (line_att, *) potential%knl(i, j, l)
2750 potential%knl(j, i, l) = potential%knl(i, j, l)
2751 CALL remove_word(line_att)
2755 cpassert(len_trim(line_att) == 0)
2760 TYPE(cp_parser_type) :: parser
2761 CALL parser_create(parser, pseudo_file)
2764 CALL parser_search_string(parser, trim(apname), .true., found, line)
2766 CALL uppercase(symbol)
2767 CALL uppercase(apname)
2770 CALL uppercase(line)
2771 line2 =
" "//line//
" "
2772 symbol2 =
" "//trim(symbol)//
" "
2773 apname2 =
" "//trim(apname)//
" "
2774 strlen1 = len_trim(symbol2) + 1
2775 strlen2 = len_trim(apname2) + 1
2777 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
2778 (index(line2, apname2(:strlen2)) > 0)) match = .true.
2783 CALL parser_get_object(parser, elec_conf(l), newline=.true.)
2784 DO WHILE (parser_test_next_token(parser) ==
"INT")
2786 CALL parser_get_object(parser, elec_conf(l))
2788 potential%econf(0:
lmat) = elec_conf(0:
lmat)
2789 potential%zion = real(sum(elec_conf), dp)
2791 CALL parser_get_object(parser, potential%rc, newline=.true.)
2793 CALL parser_get_object(parser, potential%ncl)
2794 DO i = 1, potential%ncl
2795 CALL parser_get_object(parser, potential%cl(i))
2799 CALL parser_get_next_line(parser, 1)
2800 IF (parser_test_next_token(parser) ==
"INT")
THEN
2802 ELSE IF (parser_test_next_token(parser) ==
"STR")
THEN
2803 CALL parser_get_object(parser, line)
2804 IF (index(line,
"LPOT") /= 0)
THEN
2806 potential%lpotextended = .true.
2807 CALL parser_get_object(parser, potential%nexp_lpot)
2808 DO ipot = 1, potential%nexp_lpot
2809 CALL parser_get_object(parser, potential%alpha_lpot(ipot), newline=.true.)
2810 CALL parser_get_object(parser, potential%nct_lpot(ipot))
2811 DO ic = 1, potential%nct_lpot(ipot)
2812 CALL parser_get_object(parser, potential%cval_lpot(ic, ipot))
2815 ELSE IF (index(line,
"NLCC") /= 0)
THEN
2817 potential%nlcc = .true.
2818 CALL parser_get_object(parser, potential%nexp_nlcc)
2819 DO ipot = 1, potential%nexp_nlcc
2820 CALL parser_get_object(parser, potential%alpha_nlcc(ipot), newline=.true.)
2821 CALL parser_get_object(parser, potential%nct_nlcc(ipot))
2822 DO ic = 1, potential%nct_nlcc(ipot)
2823 CALL parser_get_object(parser, potential%cval_nlcc(ic, ipot))
2825 potential%cval_nlcc(ic, ipot) = potential%cval_nlcc(ic, ipot)/(4.0_dp*pi)
2828 ELSE IF (index(line,
"LSD") /= 0)
THEN
2830 potential%lsdpot = .true.
2831 CALL parser_get_object(parser, potential%nexp_lsd)
2832 DO ipot = 1, potential%nexp_lsd
2833 CALL parser_get_object(parser, potential%alpha_lsd(ipot), newline=.true.)
2834 CALL parser_get_object(parser, potential%nct_lsd(ipot))
2835 DO ic = 1, potential%nct_lsd(ipot)
2836 CALL parser_get_object(parser, potential%cval_lsd(ic, ipot))
2840 cpabort(
"Parsing of extended potential type failed.")
2843 cpabort(
"Invalid input token found.")
2847 CALL parser_get_object(parser, nlmax)
2849 IF (parser_test_next_token(parser) ==
"STR")
THEN
2850 CALL parser_get_object(parser, line)
2851 IF (index(line,
"SOC") /= 0) potential%soc = .true.
2855 CALL parser_get_object(parser, potential%rcnl(l), newline=.true.)
2856 CALL parser_get_object(parser, potential%nl(l))
2857 DO i = 1, potential%nl(l)
2859 CALL parser_get_object(parser, potential%hnl(i, i, l))
2861 CALL parser_get_object(parser, potential%hnl(i, i, l), newline=.true.)
2863 DO j = i + 1, potential%nl(l)
2864 CALL parser_get_object(parser, potential%hnl(i, j, l))
2865 potential%hnl(j, i, l) = potential%hnl(i, j, l)
2868 IF (potential%soc .AND. l /= 0)
THEN
2869 DO i = 1, potential%nl(l)
2870 CALL parser_get_object(parser, potential%knl(i, i, l), newline=.true.)
2871 DO j = i + 1, potential%nl(l)
2872 CALL parser_get_object(parser, potential%knl(i, j, l))
2873 potential%knl(j, i, l) = potential%knl(i, j, l)
2883 cpabort(
"End of file reached unexpectedly")
2888 CALL parser_release(parser)
2892 END SUBROUTINE read_gth_potential
2902 SUBROUTINE read_ecp_potential_file(element_symbol, potential, pseudo_name, pseudo_file, &
2903 potential_section, potential_found)
2905 CHARACTER(LEN=*),
INTENT(IN) :: element_symbol
2907 CHARACTER(LEN=*),
INTENT(IN) :: pseudo_name, pseudo_file
2908 TYPE(section_vals_type),
POINTER :: potential_section
2909 LOGICAL,
INTENT(OUT),
OPTIONAL :: potential_found
2911 CHARACTER(LEN=240) :: line
2912 CHARACTER(len=5*default_string_length) :: line_att
2913 CHARACTER(LEN=LEN(element_symbol)+1) :: symbol
2914 CHARACTER(LEN=LEN(pseudo_name)) :: apname
2915 INTEGER :: i, ic, l, ncore, nel
2916 LOGICAL :: found, is_ok, read_from_input
2917 TYPE(cp_sll_val_type),
POINTER :: list
2918 TYPE(val_type),
POINTER :: val
2920 apname = pseudo_name
2921 symbol = element_symbol
2922 IF (
PRESENT(potential_found)) potential_found = .false.
2923 CALL get_ptable_info(symbol, number=ncore)
2925 potential%symbol = symbol
2926 potential%pname = apname
2932 potential%aloc = 0.0_dp
2933 potential%bloc = 0.0_dp
2936 potential%apot = 0.0_dp
2937 potential%bpot = 0.0_dp
2939 read_from_input = .false.
2940 CALL section_vals_get(potential_section, explicit=read_from_input)
2941 IF (read_from_input)
THEN
2942 CALL section_vals_list_get(potential_section,
"_DEFAULT_KEYWORD_",
list=
list)
2944 is_ok = cp_sll_val_next(
list, val)
2946 CALL val_get(val, c_val=line_att)
2947 CALL remove_word(line_att)
2948 CALL remove_word(line_att)
2950 READ (line_att, *) nel
2951 potential%zion = real(ncore - nel, kind=dp)
2953 is_ok = cp_sll_val_next(
list, val)
2955 CALL val_get(val, c_val=line_att)
2957 IF (.NOT. cp_sll_val_next(
list, val))
EXIT
2958 CALL val_get(val, c_val=line_att)
2959 IF (index(line_att, element_symbol) == 0)
THEN
2960 potential%nloc = potential%nloc + 1
2962 READ (line_att, *) potential%nrloc(ic), potential%bloc(ic), potential%aloc(ic)
2969 CALL val_get(val, c_val=line_att)
2970 IF (index(line_att, element_symbol) == 0)
THEN
2971 potential%npot(l) = potential%npot(l) + 1
2972 ic = potential%npot(l)
2973 READ (line_att, *) potential%nrpot(ic, l), potential%bpot(ic, l), potential%apot(ic, l)
2975 potential%lmax = potential%lmax + 1
2978 IF (.NOT. cp_sll_val_next(
list, val))
EXIT
2983 TYPE(cp_parser_type) :: parser
2984 CALL parser_create(parser, pseudo_file)
2987 CALL parser_search_string(parser, trim(apname), .true., found, line)
2990 CALL parser_get_object(parser, line, newline=.true.)
2991 IF (trim(line) == element_symbol)
THEN
2992 IF (
PRESENT(potential_found)) potential_found = .true.
2993 CALL parser_get_object(parser, line, lower_to_upper=.true.)
2994 cpassert(trim(line) ==
"NELEC")
2996 CALL parser_get_object(parser, nel)
2997 potential%zion = real(ncore - nel, kind=dp)
2999 CALL parser_get_object(parser, line, newline=.true.)
3002 CALL parser_read_line(parser, 1)
3003 IF (parser_test_next_token(parser) ==
"STR")
EXIT
3004 potential%nloc = potential%nloc + 1
3006 CALL parser_get_object(parser, potential%nrloc(ic))
3007 CALL parser_get_object(parser, potential%bloc(ic))
3008 CALL parser_get_object(parser, potential%aloc(ic))
3012 CALL parser_get_object(parser, symbol)
3013 IF (symbol == element_symbol)
THEN
3015 potential%lmax = potential%lmax + 1
3017 CALL parser_read_line(parser, 1)
3018 IF (parser_test_next_token(parser) ==
"STR")
EXIT
3019 potential%npot(l) = potential%npot(l) + 1
3020 ic = potential%npot(l)
3021 CALL parser_get_object(parser, potential%nrpot(ic, l))
3022 CALL parser_get_object(parser, potential%bpot(ic, l))
3023 CALL parser_get_object(parser, potential%apot(ic, l))
3030 ELSE IF (line ==
"END")
THEN
3031 IF (
PRESENT(potential_found))
THEN
3032 CALL parser_release(parser)
3035 cpabort(
"Element not found in ECP library")
3039 IF (
PRESENT(potential_found))
THEN
3040 CALL parser_release(parser)
3043 cpabort(
"ECP type not found in library")
3048 CALL parser_release(parser)
3052 IF (read_from_input .AND.
PRESENT(potential_found)) potential_found = .true.
3055 potential%econf(0:3) = ptable(ncore)%e_conv(0:3)
3058 cpabort(
"Unknown Core State")
3061 potential%econf(0:3) = potential%econf(0:3) - ptable(2)%e_conv(0:3)
3063 potential%econf(0:3) = potential%econf(0:3) - ptable(10)%e_conv(0:3)
3065 potential%econf(0:3) = potential%econf(0:3) - ptable(18)%e_conv(0:3)
3067 potential%econf(0:3) = potential%econf(0:3) - ptable(18)%e_conv(0:3)
3068 potential%econf(2) = potential%econf(2) - 10
3070 potential%econf(0:3) = potential%econf(0:3) - ptable(36)%e_conv(0:3)
3072 potential%econf(0:3) = potential%econf(0:3) - ptable(36)%e_conv(0:3)
3073 potential%econf(2) = potential%econf(2) - 10
3075 potential%econf(0:3) = potential%econf(0:3) - ptable(54)%e_conv(0:3)
3077 potential%econf(0:3) = potential%econf(0:3) - ptable(36)%e_conv(0:3)
3078 potential%econf(2) = potential%econf(2) - 10
3079 potential%econf(3) = potential%econf(3) - 14
3081 potential%econf(0:3) = potential%econf(0:3) - ptable(54)%e_conv(0:3)
3082 potential%econf(3) = potential%econf(3) - 14
3084 potential%econf(0:3) = potential%econf(0:3) - ptable(54)%e_conv(0:3)
3085 potential%econf(2) = potential%econf(2) - 10
3086 potential%econf(3) = potential%econf(3) - 14
3089 cpassert(all(potential%econf >= 0))
3091 END SUBROUTINE read_ecp_potential_file
3101 SUBROUTINE read_ecp_potential_files(element_symbol, potential, pseudo_name, pseudo_files, &
3104 CHARACTER(LEN=*),
INTENT(IN) :: element_symbol
3106 CHARACTER(LEN=*),
INTENT(IN) :: pseudo_name
3107 CHARACTER(LEN=*),
DIMENSION(:),
INTENT(IN) :: pseudo_files
3108 TYPE(section_vals_type),
POINTER :: potential_section
3110 CHARACTER(LEN=:),
ALLOCATABLE :: file_list
3112 LOGICAL :: potential_found
3114 DO i = 1,
SIZE(pseudo_files)
3115 CALL read_ecp_potential_file(element_symbol, potential, pseudo_name, pseudo_files(i), &
3116 potential_section, potential_found)
3117 IF (potential_found)
RETURN
3121 DO i = 1,
SIZE(pseudo_files)
3122 file_list = trim(file_list)//
"<"//trim(pseudo_files(i))//
"> "
3124 CALL cp_abort(__location__, &
3125 "The requested ECP potential <"//trim(pseudo_name)// &
3126 "> for element <"//trim(element_symbol)// &
3127 "> was not found in the potential files "//trim(file_list))
3129 END SUBROUTINE read_ecp_potential_files
3137 TYPE(grid_atom_type) :: grid1, grid2
3141 REAL(kind=dp) :: dr, dw
3144 IF (grid1%nr == grid2%nr)
THEN
3146 dr = abs(grid1%rad(i) - grid2%rad(i))
3147 dw = abs(grid1%wr(i) - grid2%wr(i))
3148 IF (dr + dw > 1.0e-12_dp)
THEN
Define the atom type and its sub types.
integer, parameter, public num_basis
subroutine, public read_atom_opt_section(optimization, opt_section)
...
subroutine, public create_atom_type(atom)
...
integer, parameter, public cgto_basis
integer, parameter, public gto_basis
integer, parameter, public sto_basis
subroutine, public release_atom_type(atom)
...
subroutine, public release_opmat(opmat)
...
subroutine, public release_atom_potential(potential)
...
subroutine, public init_atom_basis(basis, basis_section, zval, btyp)
Initialize the basis for the atomic code.
logical function, public atom_compare_grids(grid1, grid2)
...
subroutine, public release_opgrid(opgrid)
...
subroutine, public release_atom_orbs(orbs)
...
integer, parameter, public lmat
subroutine, public set_atom(atom, basis, state, integrals, orbitals, potential, zcore, pp_calc, do_zmp, doread, read_vxc, method_type, relativistic, coulomb_integral_type, exchange_integral_type, fmat)
...
subroutine, public init_atom_potential(potential, potential_section, zval)
...
subroutine, public release_atom_basis(basis)
...
subroutine, public create_atom_orbs(orbs, mbas, mo)
...
subroutine, public init_atom_basis_default_pp(basis)
...
subroutine, public create_opmat(opmat, n, lmax)
...
subroutine, public atom_basis_gridrep(basis, gbasis, r, rab)
...
subroutine, public setup_hf_section(hf_frac, do_hfx, atom, xc_section, extype)
...
subroutine, public clementi_geobas(zval, cval, aval, ngto, ival)
...
subroutine, public create_opgrid(opgrid, grid)
...
Routines that process Quantum Espresso UPF files.
subroutine, public atom_read_upf(pot, upf_filename, read_header)
...
pure subroutine, public atom_release_upf(upfpot)
...
Calculates Bessel functions.
elemental impure real(kind=dp) function, public bessel0(x, l)
...
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public limpanuparb2011
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_read_line(parser, nline, at_end)
Read the next line from a logical unit "unit" (I/O node only). Skip (nline-1) lines and skip also all...
subroutine, public parser_get_next_line(parser, nline, at_end)
Read the next input line and broadcast the input information. Skip (nline-1) lines and skip also all ...
character(len=3) function, public parser_test_next_token(parser, string_length)
Test next input object.
subroutine, public parser_search_string(parser, string, ignore_case, found, line, begin_line, search_from_begin_of_file)
Search a string pattern in a file defined by its logical unit number "unit". A case sensitive search ...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_release(parser)
releases the parser
subroutine, public parser_create(parser, file_name, unit_nr, para_env, end_section_label, separator_chars, comment_char, continuation_char, quote_char, section_char, parse_white_lines, initial_variables, apply_preprocessing)
Start a parser run. Initial variables allow to @SET stuff before opening the file.
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), dimension(-1:2 *maxfac+1), parameter, public dfac
real(kind=dp), parameter, public rootpi
real(kind=dp), dimension(0:maxfac), parameter, public fac
Periodic Table related data definitions.
type(atom), dimension(0:nelem), public ptable
subroutine, public get_ptable_info(symbol, number, amass, ielement, covalent_radius, metallic_radius, vdw_radius, found)
Pass information about the kind given the element symbol.
subroutine, public deallocate_grid_atom(grid_atom)
Deallocate a Gaussian-type orbital (GTO) basis set data set.
subroutine, public allocate_grid_atom(grid_atom)
Initialize components of the grid_atom_type structure.
subroutine, public create_grid_atom(grid_atom, nr, na, llmax, ll, quadrature)
...
Utilities for string manipulations.
character(len=1), parameter, public newline
subroutine, public remove_word(string)
remove a word from a string (words are separated by white spaces)
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.
Provides all information about a basis set.
Provides all information about a pseudopotential.
Provides info about hartree-fock exchange (For now, we only support potentials that can be represente...
Information on optimization procedure.
Holds atomic orbitals and energies.
Provides all information on states and occupation.
Provides all information about an atomic kind.