28 atom_basis_type, atom_gthpot_type, atom_integrals, atom_optimization_type, atom_orbitals, &
56 #include "./base/base_uses.f90"
62 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'atom_energy'
80 TYPE(section_vals_type),
POINTER :: atom_section
82 CHARACTER(len=*),
PARAMETER :: routinen =
'atom_energy_opt'
84 CHARACTER(LEN=2) :: elem
85 CHARACTER(LEN=default_string_length) :: filename
86 CHARACTER(LEN=default_string_length), &
87 DIMENSION(:),
POINTER :: tmpstringlist
88 INTEGER :: do_eric, do_erie, handle, i, im, in, iw, k, maxl, mb, method, mo, n_meth, n_rep, &
89 nder, nr_gh, num_gau, num_gto, num_pol, reltyp, zcore, zval, zz
90 INTEGER,
DIMENSION(0:lmat) :: maxn
91 INTEGER,
DIMENSION(:),
POINTER :: cn
92 LOGICAL :: dm, do_gh, do_zmp, doread, eri_c, eri_e, &
93 had_ae, had_pp, lcomp, lcond, pp_calc, &
95 REAL(kind=
dp) :: crad, delta, lambda
96 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: ext_density, ext_vxc
97 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: pocc
98 TYPE(atom_basis_type),
POINTER :: ae_basis, pp_basis
99 TYPE(atom_integrals),
POINTER :: ae_int, pp_int
100 TYPE(atom_optimization_type) :: optimization
101 TYPE(atom_orbitals),
POINTER :: orbitals
102 TYPE(atom_p_type),
DIMENSION(:, :),
POINTER :: atom_info
103 TYPE(atom_potential_type),
POINTER :: ae_pot, p_pot
104 TYPE(atom_state),
POINTER :: state
105 TYPE(cp_logger_type),
POINTER :: logger
106 TYPE(section_vals_type),
POINTER :: admm_section, basis_section, external_vxc_section, &
107 method_section, opt_section, potential_section, powell_section, sgp_section, xc_section, &
108 zmp_restart_section, zmp_section
110 CALL timeset(routinen, handle)
117 IF (
ptable(i)%symbol == elem)
THEN
122 IF (zz /= 1) zval = zz
125 ALLOCATE (ae_basis, pp_basis)
127 NULLIFY (ae_basis%grid)
129 NULLIFY (pp_basis%grid)
136 IF (iw > 0)
CALL atom_print_info(zval,
"Atomic Energy Calculation", iw)
146 NULLIFY (potential_section)
148 ALLOCATE (ae_pot, p_pot)
166 DO in = 1, min(
SIZE(cn), 4)
167 maxn(in - 1) = cn(in)
170 maxn(in) = min(maxn(in), ae_basis%nbas(in))
171 maxn(in) = min(maxn(in), pp_basis%nbas(in))
188 ALLOCATE (ae_int, pp_int)
190 ALLOCATE (atom_info(n_rep, n_meth))
195 NULLIFY (atom_info(in, im)%atom)
198 atom_info(in, im)%atom%optimization = optimization
200 atom_info(in, im)%atom%z = zval
202 atom_info(in, im)%atom%xc_section => xc_section
211 atom_info(in, im)%atom%do_zmp = do_zmp
213 atom_info(in, im)%atom%ext_file = filename
215 r_val=atom_info(in, im)%atom%zmpgrid_tol)
217 atom_info(in, im)%atom%lambda = lambda
219 atom_info(in, im)%atom%dm = dm
223 atom_info(in, im)%atom%doread = doread
225 atom_info(in, im)%atom%zmp_restart_file = filename
230 atom_info(in, im)%atom%read_vxc = read_vxc
232 atom_info(in, im)%atom%ext_vxc_file = filename
234 r_val=atom_info(in, im)%atom%zmpvxcgrid_tol)
240 c_vals=tmpstringlist)
248 state%maxl_calc = max(maxl, state%maxl_occ)
249 state%maxl_calc = min(
lmat, state%maxl_calc)
251 DO k = 0, state%maxl_calc
252 state%maxn_calc(k) = max(maxn(k), state%maxn_occ(k))
256 pp_calc = any(index(tmpstringlist(1:),
"CORE") /= 0)
261 zcore = zval - nint(sum(state%core))
262 CALL set_atom(atom_info(in, im)%atom, zcore=zcore, pp_calc=.true.)
265 CALL set_atom(atom_info(in, im)%atom, zcore=zval, pp_calc=.false.)
270 CALL set_atom(atom_info(in, im)%atom, method_type=method, relativistic=reltyp)
290 potential=p_pot, eri_coulomb=eri_c, eri_exchange=eri_e)
294 NULLIFY (pp_int%tzora, pp_int%hdkh)
296 CALL set_atom(atom_info(in, im)%atom, basis=pp_basis, integrals=pp_int, potential=p_pot)
297 state%maxn_calc(:) = min(state%maxn_calc(:), pp_basis%nbas(:))
298 cpassert(all(state%maxn_calc(:) >= state%maxn_occ))
303 eri_coulomb=eri_c, eri_exchange=eri_e)
309 CALL set_atom(atom_info(in, im)%atom, basis=ae_basis, integrals=ae_int, potential=ae_pot)
310 state%maxn_calc(:) = min(state%maxn_calc(:), ae_basis%nbas(:))
311 cpassert(all(state%maxn_calc(:) >= state%maxn_occ))
315 CALL set_atom(atom_info(in, im)%atom, state=state)
317 CALL set_atom(atom_info(in, im)%atom, coulomb_integral_type=do_eric, &
318 exchange_integral_type=do_erie)
319 atom_info(in, im)%atom%hfx_pot%do_gh = do_gh
320 atom_info(in, im)%atom%hfx_pot%nr_gh = nr_gh
323 mo = maxval(state%maxn_calc)
324 mb = maxval(atom_info(in, im)%atom%basis%nbas)
326 CALL set_atom(atom_info(in, im)%atom, orbitals=orbitals)
335 IF (atom_info(in, im)%atom%do_zmp)
THEN
338 ALLOCATE (ext_density(atom_info(in, im)%atom%basis%grid%nr))
342 atom=atom_info(in, im)%atom, &
344 DEALLOCATE (ext_density)
347 IF (atom_info(in, im)%atom%read_vxc)
THEN
348 ALLOCATE (ext_vxc(atom_info(in, im)%atom%basis%grid%nr))
352 atom=atom_info(in, im)%atom, &
369 CALL atom_fit_density(atom_info(in, im)%atom, num_gto, 0, iw, powell_section=powell_section)
379 CALL atom_fit_kgpot(atom_info(in, im)%atom, num_gau, num_pol, iw, powell_section=powell_section)
388 CALL atom_response_basis(atom_info(in, im)%atom, delta, nder, iw)
394 file_position=
"REWIND")
396 CALL atom_write_upf(atom_info(in, im)%atom, iw)
405 iw =
cp_print_key_unit_nr(logger, atom_section,
"PRINT%GEOMETRICAL_RESPONSE_BASIS", extension=
".log")
414 CALL section_vals_val_get(atom_section,
"PRINT%ANALYZE_BASIS%OVERLAP_CONDITION_NUMBER", l_val=lcond)
432 CALL atom_admm(atom_info, admm_section, iw)
437 iw =
cp_print_key_unit_nr(logger, atom_section,
"PRINT%SEPARABLE_GAUSSIAN_PSEUDO", extension=
".log")
466 DEALLOCATE (atom_info)
468 DEALLOCATE (ae_pot, p_pot, ae_basis, pp_basis, ae_int, pp_int)
470 CALL timestop(handle)
484 SUBROUTINE atom_response_basis(atom, delta, nder, iw)
486 TYPE(atom_type),
POINTER ::
atom
487 REAL(kind=
dp),
INTENT(IN) :: delta
488 INTEGER,
INTENT(IN) :: nder, iw
490 INTEGER :: i, ider, j, k, l, lhomo, m, n, nhomo, &
492 REAL(kind=
dp) :: dene, emax, expzet, fhomo, o, prefac, &
494 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: amat
495 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: rbasis
496 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :, :) :: wfn
497 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: ovlp
498 TYPE(atom_state),
POINTER :: state
500 WRITE (iw,
'(/," ",79("*"),/,T30,A,A/," ",79("*"))')
"RESPONSE BASIS for ",
ptable(
atom%z)%symbol
503 ovlp =>
atom%integrals%ovlp
509 DO l = 0, state%maxl_occ
510 DO i = 1, state%maxn_occ(l)
511 IF (
atom%orbitals%ener(i, l) > emax)
THEN
514 emax =
atom%orbitals%ener(i, l)
515 fhomo = state%occupation(l, i)
520 s1 =
SIZE(
atom%orbitals%wfn, 1)
521 s2 =
SIZE(
atom%orbitals%wfn, 2)
522 ALLOCATE (wfn(s1, s2, 0:
lmat, -nder:nder))
523 s2 = maxval(state%maxn_occ) + nder
524 ALLOCATE (rbasis(s1, s2, 0:
lmat))
527 DO ider = -nder, nder
528 dene = real(ider, kind=
dp)*delta
529 cpassert(fhomo > abs(dene))
530 state%occupation(lhomo, nhomo) = fhomo + dene
532 wfn(:, :, :, ider) =
atom%orbitals%wfn
533 state%occupation(lhomo, nhomo) = fhomo
536 DO l = 0, state%maxl_occ
538 DO i = 1, max(state%maxn_occ(l), 1)
539 rbasis(:, i, l) = wfn(:, i, l, 0)
543 i = max(state%maxn_occ(l), 1)
546 rbasis(:, i + 1, l) = 0.5_dp*(wfn(:, i, l, 1) - wfn(:, i, l, -1))/delta
548 rbasis(:, i + 2, l) = 0.25_dp*(wfn(:, i, l, 2) - 2._dp*wfn(:, i, l, 0) + wfn(:, i, l, -2))/delta**2
550 rbasis(:, i + 3, l) = 0.125_dp*(wfn(:, i, l, 3) - 3._dp*wfn(:, i, l, 1) &
551 + 3._dp*wfn(:, i, l, -1) - wfn(:, i, l, -3))/delta**3
558 n = state%maxn_occ(l) + nder
559 m =
atom%basis%nbas(l)
562 o = dot_product(rbasis(1:m, j, l), reshape(matmul(ovlp(1:m, 1:m, l), rbasis(1:m, i:i, l)), (/m/)))
563 rbasis(1:m, i, l) = rbasis(1:m, i, l) - o*rbasis(1:m, j, l)
565 o = dot_product(rbasis(1:m, i, l), reshape(matmul(ovlp(1:m, 1:m, l), rbasis(1:m, i:i, l)), (/m/)))
566 rbasis(1:m, i, l) = rbasis(1:m, i, l)/sqrt(o)
570 ALLOCATE (amat(n, n))
571 amat(1:n, 1:n) = matmul(transpose(rbasis(1:m, 1:n, l)), matmul(ovlp(1:m, 1:m, l), rbasis(1:m, 1:n, l)))
573 amat(i, i) = amat(i, i) - 1._dp
575 IF (maxval(abs(amat)) > 1.e-12)
THEN
576 WRITE (iw,
'(/,A,G20.10)')
" Orthogonality error ", maxval(abs(amat))
581 WRITE (iw,
'(/,A,T30,I3)')
" Angular momentum :", l
583 WRITE (iw,
'(/,A,I0,A,I0,A)')
" Exponent Coef.(Quickstep Normalization), first ", &
584 n - nder,
" valence ", nder,
" response"
585 expzet = 0.25_dp*real(2*l + 3,
dp)
586 prefac = sqrt(sqrt(
pi)/2._dp**(l + 2)*
dfac(2*l + 1))
588 zeta = (2._dp*
atom%basis%am(i, l))**expzet
589 WRITE (iw,
'(4X,F20.10,4X,15ES20.6)')
atom%basis%am(i, l), ((prefac*rbasis(i, k, l)/zeta), k=1, n)
594 DEALLOCATE (wfn, rbasis)
596 WRITE (iw,
'(" ",79("*"))')
598 END SUBROUTINE atom_response_basis
607 SUBROUTINE atom_write_upf(atom, iw)
609 TYPE(atom_type),
POINTER ::
atom
610 INTEGER,
INTENT(IN) :: iw
612 CHARACTER(LEN=default_string_length) :: string
613 INTEGER :: i, ibeta, j, k, l, lmax, nbeta, nr, &
616 REAL(kind=
dp) :: pf, rl, rmax
617 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: beta, corden, dens, ef, locpot, rp
618 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: dij
619 TYPE(atom_gthpot_type),
POINTER :: pot
621 IF (.NOT.
atom%pp_calc)
RETURN
622 IF (
atom%potential%ppot_type /= gth_pseudo)
RETURN
623 pot =>
atom%potential%gth_pot
624 cpassert(.NOT. pot%lsdpot)
626 WRITE (iw,
'(A)')
'<UPF version="2.0.1">'
628 WRITE (iw,
'(T4,A)')
'<PP_INFO>'
629 WRITE (iw,
'(T8,A)')
'Converted from CP2K GTH format'
630 WRITE (iw,
'(T8,A)')
'<PP_INPUTFILE>'
632 WRITE (iw,
'(T8,A)')
'</PP_INPUTFILE>'
633 WRITE (iw,
'(T4,A)')
'</PP_INFO>'
634 WRITE (iw,
'(T4,A)')
'<PP_HEADER'
635 WRITE (iw,
'(T8,A)')
'generated="Generated in analytical, separable form"'
636 WRITE (iw,
'(T8,A)')
'author="Goedecker/Hartwigsen/Hutter/Teter"'
637 WRITE (iw,
'(T8,A)')
'date="Phys.Rev.B58, 3641 (1998); B54, 1703 (1996)"'
638 WRITE (iw,
'(T8,A)')
'comment="This file generated by CP2K ATOM code"'
639 CALL compose(string,
"element", cval=pot%symbol)
640 WRITE (iw,
'(T8,A)') trim(string)
641 WRITE (iw,
'(T8,A)')
'pseudo_type="NC"'
642 WRITE (iw,
'(T8,A)')
'relativistic="no"'
643 WRITE (iw,
'(T8,A)')
'is_ultrasoft="F"'
644 WRITE (iw,
'(T8,A)')
'is_paw="F"'
645 WRITE (iw,
'(T8,A)')
'is_coulomb="F"'
646 WRITE (iw,
'(T8,A)')
'has_so="F"'
647 WRITE (iw,
'(T8,A)')
'has_wfc="F"'
648 WRITE (iw,
'(T8,A)')
'has_gipaw="F"'
649 WRITE (iw,
'(T8,A)')
'paw_as_gipaw="F"'
651 WRITE (iw,
'(T8,A)')
'core_correction="T"'
653 WRITE (iw,
'(T8,A)')
'core_correction="F"'
655 WRITE (iw,
'(T8,A)')
'functional="DFT"'
656 CALL compose(string,
"z_valence", rval=pot%zion)
657 WRITE (iw,
'(T8,A)') trim(string)
658 CALL compose(string,
"total_psenergy", rval=2._dp*
atom%energy%etot)
659 WRITE (iw,
'(T8,A)') trim(string)
660 WRITE (iw,
'(T8,A)')
'wfc_cutoff="0.0E+00"'
661 WRITE (iw,
'(T8,A)')
'rho_cutoff="0.0E+00"'
664 IF (pot%nl(l) > 0) lmax = l
666 CALL compose(string,
"l_max", ival=lmax)
667 WRITE (iw,
'(T8,A)') trim(string)
668 WRITE (iw,
'(T8,A)')
'l_max_rho="0"'
669 WRITE (iw,
'(T8,A)')
'l_local="-3"'
670 nr =
atom%basis%grid%nr
671 CALL compose(string,
"mesh_size", ival=nr)
672 WRITE (iw,
'(T8,A)') trim(string)
673 nwfc = sum(
atom%state%maxn_occ)
674 CALL compose(string,
"number_of_wfc", ival=nwfc)
675 WRITE (iw,
'(T8,A)') trim(string)
677 CALL compose(string,
"number_of_proj", ival=nbeta)
678 WRITE (iw,
'(T8,A)') trim(string)//
'/>'
681 up =
atom%basis%grid%rad(1) <
atom%basis%grid%rad(nr)
682 WRITE (iw,
'(T4,A)')
'<PP_MESH'
683 WRITE (iw,
'(T8,A)')
'dx="1.E+00"'
684 CALL compose(string,
"mesh", ival=nr)
685 WRITE (iw,
'(T8,A)') trim(string)
686 WRITE (iw,
'(T8,A)')
'xmin="1.E+00"'
687 rmax = maxval(
atom%basis%grid%rad)
688 CALL compose(string,
"rmax", rval=rmax)
689 WRITE (iw,
'(T8,A)') trim(string)
690 WRITE (iw,
'(T8,A)')
'zmesh="1.E+00">'
691 WRITE (iw,
'(T8,A)')
'<PP_R type="real"'
692 CALL compose(string,
"size", ival=nr)
693 WRITE (iw,
'(T12,A)') trim(string)
694 WRITE (iw,
'(T12,A)')
'columns="4">'
696 WRITE (iw,
'(T12,4ES25.12E3)') (
atom%basis%grid%rad(i), i=1, nr)
698 WRITE (iw,
'(T12,4ES25.12E3)') (
atom%basis%grid%rad(i), i=nr, 1, -1)
700 WRITE (iw,
'(T8,A)')
'</PP_R>'
701 WRITE (iw,
'(T8,A)')
'<PP_RAB type="real"'
702 CALL compose(string,
"size", ival=nr)
703 WRITE (iw,
'(T12,A)') trim(string)
704 WRITE (iw,
'(T12,A)')
'columns="4">'
706 WRITE (iw,
'(T12,4ES25.12E3)') (
atom%basis%grid%wr(i)/
atom%basis%grid%rad2(i), i=1, nr)
708 WRITE (iw,
'(T12,4ES25.12E3)') (
atom%basis%grid%wr(i)/
atom%basis%grid%rad2(i), i=nr, 1, -1)
710 WRITE (iw,
'(T8,A)')
'</PP_RAB>'
711 WRITE (iw,
'(T4,A)')
'</PP_MESH>'
715 WRITE (iw,
'(T4,A)')
'<PP_NLCC type="real"'
716 CALL compose(string,
"size", ival=nr)
717 WRITE (iw,
'(T8,A)') trim(string)
718 WRITE (iw,
'(T8,A)')
'columns="4">'
719 ALLOCATE (corden(nr))
723 WRITE (iw,
'(T8,4ES25.12E3)') (corden(i), i=1, nr)
725 WRITE (iw,
'(T8,4ES25.12E3)') (corden(i), i=nr, 1, -1)
728 WRITE (iw,
'(T4,A)')
'</PP_NLCC>'
732 WRITE (iw,
'(T4,A)')
'<PP_LOCAL type="real"'
733 CALL compose(string,
"size", ival=nr)
734 WRITE (iw,
'(T8,A)') trim(string)
735 WRITE (iw,
'(T8,A)')
'columns="4">'
736 ALLOCATE (locpot(nr))
740 WRITE (iw,
'(T8,4ES25.12E3)') (2.0_dp*locpot(i), i=1, nr)
742 WRITE (iw,
'(T8,4ES25.12E3)') (2.0_dp*locpot(i), i=nr, 1, -1)
745 WRITE (iw,
'(T4,A)')
'</PP_LOCAL>'
748 WRITE (iw,
'(T4,A)')
'<PP_NONLOCAL>'
749 ALLOCATE (rp(nr), ef(nr), beta(nr))
752 IF (pot%nl(l) == 0) cycle
754 rp(:) =
atom%basis%grid%rad
755 ef(:) = exp(-0.5_dp*rp*rp/(rl*rl))
757 pf = rl**(l + 0.5_dp*(4._dp*i - 1._dp))
759 pf = sqrt(2._dp)/(pf*sqrt(
gamma1(j)))
760 beta(:) = pf*rp**(l + 2*i - 2)*ef
762 CALL compose(string,
"<PP_BETA", counter=ibeta)
763 WRITE (iw,
'(T8,A)') trim(string)
764 CALL compose(string,
"angular_momentum", ival=l)
765 WRITE (iw,
'(T12,A)') trim(string)
766 WRITE (iw,
'(T12,A)')
'type="real"'
767 CALL compose(string,
"size", ival=nr)
768 WRITE (iw,
'(T12,A)') trim(string)
769 WRITE (iw,
'(T12,A)')
'columns="4">'
772 WRITE (iw,
'(T12,4ES25.12E3)') (2._dp*beta(j), j=1, nr)
774 WRITE (iw,
'(T12,4ES25.12E3)') (2._dp*beta(j), j=nr, 1, -1)
776 CALL compose(string,
"</PP_BETA", counter=ibeta, isfinal=.true.)
777 WRITE (iw,
'(T8,A)') trim(string)
780 DEALLOCATE (rp, ef, beta)
782 ALLOCATE (dij(nbeta, nbeta))
785 IF (pot%nl(l) == 0) cycle
786 ibeta = sum(pot%nl(0:l - 1)) + 1
787 i = ibeta + pot%nl(l) - 1
788 dij(ibeta:i, ibeta:i) = pot%hnl(1:pot%nl(l), 1:pot%nl(l), l)
790 WRITE (iw,
'(T8,A)')
'<PP_DIJ type="real"'
791 WRITE (iw,
'(T12,A)')
'columns="4">'
792 WRITE (iw,
'(T12,4ES25.12E3)') ((0.5_dp*dij(i, j), j=1, nbeta), i=1, nbeta)
793 WRITE (iw,
'(T8,A)')
'</PP_DIJ>'
795 WRITE (iw,
'(T4,A)')
'</PP_NONLOCAL>'
799 WRITE (iw,
'(T4,A)')
'<PP_PSWFC>'
803 IF (abs(
atom%state%occupation(l, i)) == 0._dp) cycle
805 CALL compose(string,
"<PP_CHI", counter=nwfn)
806 WRITE (iw,
'(T8,A)') trim(string)
807 CALL compose(string,
"l", ival=l)
808 WRITE (iw,
'(T12,A)') trim(string)
809 CALL compose(string,
"occupation", rval=
atom%state%occupation(l, i))
810 WRITE (iw,
'(T12,A)') trim(string)
811 CALL compose(string,
"pseudo_energy", rval=2._dp*
atom%orbitals%ener(i, l))
812 WRITE (iw,
'(T12,A)') trim(string)
813 WRITE (iw,
'(T12,A)')
'columns="4">'
815 DO k = 1,
atom%basis%nbas(l)
816 beta(:) = beta(:) +
atom%orbitals%wfn(k, i, l)*
atom%basis%bf(:, k, l)
818 beta(:) = beta*
atom%basis%grid%rad
820 WRITE (iw,
'(T12,4ES25.12E3)') (beta(j)*
atom%basis%grid%rad(j), j=1, nr)
822 WRITE (iw,
'(T12,4ES25.12E3)') (beta(j)*
atom%basis%grid%rad(j), j=nr, 1, -1)
824 CALL compose(string,
"</PP_CHI", counter=nwfn, isfinal=.true.)
825 WRITE (iw,
'(T8,A)') trim(string)
828 WRITE (iw,
'(T4,A)')
'</PP_PSWFC>'
833 WRITE (iw,
'(T4,A)')
'<PP_RHOATOM type="real"'
834 CALL compose(string,
"size", ival=nr)
835 WRITE (iw,
'(T8,A)') trim(string)
836 WRITE (iw,
'(T8,A)')
'columns="4">'
838 "RHO",
atom%basis%grid%rad)
840 WRITE (iw,
'(T8,4ES25.12E3)') (4._dp*
pi*dens(j)*
atom%basis%grid%rad2(j), j=1, nr)
842 WRITE (iw,
'(T8,4ES25.12E3)') (4._dp*
pi*dens(j)*
atom%basis%grid%rad2(j), j=nr, 1, -1)
844 WRITE (iw,
'(T4,A)')
'</PP_RHOATOM>'
847 WRITE (iw,
'(A)')
'</UPF>'
849 END SUBROUTINE atom_write_upf
863 SUBROUTINE compose(string, tag, counter, rval, ival, cval, isfinal)
864 CHARACTER(len=*),
INTENT(OUT) :: string
865 CHARACTER(len=*),
INTENT(IN) :: tag
866 INTEGER,
INTENT(IN),
OPTIONAL :: counter
867 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: rval
868 INTEGER,
INTENT(IN),
OPTIONAL :: ival
869 CHARACTER(len=*),
INTENT(IN),
OPTIONAL :: cval
870 LOGICAL,
INTENT(IN),
OPTIONAL :: isfinal
872 CHARACTER(LEN=default_string_length) :: str
875 IF (
PRESENT(counter))
THEN
876 WRITE (str,
"(I12)") counter
877 ELSEIF (
PRESENT(rval))
THEN
878 WRITE (str,
"(G18.8)") rval
879 ELSEIF (
PRESENT(ival))
THEN
880 WRITE (str,
"(I12)") ival
881 ELSEIF (
PRESENT(cval))
THEN
882 WRITE (str,
"(A)") trim(adjustl(cval))
884 WRITE (str,
"(A)")
""
887 IF (
PRESENT(isfinal)) fin = isfinal
888 IF (
PRESENT(counter))
THEN
890 WRITE (string,
"(A,A1,A,A1)") trim(adjustl(tag)),
'.', trim(adjustl(str)),
'>'
892 WRITE (string,
"(A,A1,A)") trim(adjustl(tag)),
'.', trim(adjustl(str))
896 WRITE (string,
"(A,A2,A,A2)") trim(adjustl(tag)),
'="', trim(adjustl(str)),
'>"'
898 WRITE (string,
"(A,A2,A,A1)") trim(adjustl(tag)),
'="', trim(adjustl(str)),
'"'
901 END SUBROUTINE compose
static GRID_HOST_DEVICE orbital up(const int i, const orbital a)
Increase i'th component of given orbital angular momentum.
subroutine, public atom_admm(atom_info, admm_section, iw)
Analysis of ADMM approximation to exact exchange.
subroutine, public calculate_atom(atom, iw, noguess, converged)
General routine to perform electronic structure atomic calculations.
subroutine, public atom_energy_opt(atom_section)
Compute the atomic energy.
routines that fit parameters for /from atomic calculations
subroutine, public atom_fit_density(atom, num_gto, norder, iunit, powell_section, results)
Fit the atomic electron density using a geometrical Gaussian basis set.
subroutine, public atom_fit_kgpot(atom, num_gau, num_pol, iunit, powell_section, results)
...
subroutine, public atom_grb_construction(atom_info, atom_section, iw)
Construct geometrical response basis set.
Calculate the atomic operator matrices.
subroutine, public atom_ppint_release(integrals)
Release memory allocated for atomic integrals (core electrons).
subroutine, public atom_int_setup(integrals, basis, potential, eri_coulomb, eri_exchange, all_nu)
Set up atomic integrals.
subroutine, public atom_relint_setup(integrals, basis, reltyp, zcore, alpha)
...
subroutine, public atom_relint_release(integrals)
Release memory allocated for atomic integrals (relativistic effects).
subroutine, public atom_ppint_setup(integrals, basis, potential)
...
subroutine, public atom_int_release(integrals)
Release memory allocated for atomic integrals (valence electrons).
Routines that print various information about an atomic kind.
subroutine, public atom_print_basis(atom_basis, iw, title)
Print atomic basis set.
subroutine, public atom_write_pseudo_param(gthpot, iunit)
Print GTH pseudo-potential parameters.
subroutine, public atom_print_method(atom, iw)
Print information about the electronic structure method in use.
subroutine, public atom_print_orbitals(atom, iw)
Print atomic orbitals.
subroutine, public atom_print_potential(potential, iw)
Print information about the pseudo-potential.
subroutine, public atom_print_info(zval, info, iw)
Print an information string related to the atomic kind.
subroutine, public atom_sgp_construction(atom_info, input_section, iw)
...
Define the atom type and its sub types.
subroutine, public read_atom_opt_section(optimization, opt_section)
...
subroutine, public create_atom_type(atom)
...
subroutine, public release_atom_type(atom)
...
subroutine, public release_atom_potential(potential)
...
subroutine, public init_atom_basis(basis, basis_section, zval, btyp)
Initialize the basis for the atomic code.
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)
...
Some basic routines for atomic calculations.
subroutine, public atom_condnumber(basis, crad, iw)
Print condition numbers of the given atomic basis set.
pure subroutine, public atom_local_potential(locpot, gthpot, rr)
...
subroutine, public atom_read_external_vxc(vxc, atom, iw)
ZMP subroutine to read external v_xc in the atomic code.
pure logical function, public atom_consistent_method(method, multiplicity)
Check that the atomic multiplicity is consistent with the electronic structure method.
subroutine, public atom_core_density(corden, potential, typ, rr)
...
pure integer function, dimension(0:lmat), public get_maxn_occ(occupation)
Return the maximum principal quantum number of occupied orbitals.
subroutine, public atom_density(density, pmat, basis, maxl, typ, rr)
Map the electron density on an atomic radial grid.
subroutine, public atom_read_external_density(density, atom, iw)
ZMP subroutine to read external density from linear grid of density matrix.
subroutine, public atom_completeness(basis, zv, iw)
Estimate completeness of the given atomic basis set.
subroutine, public atom_write_zmp_restart(atom)
ZMP subroutine to write external restart file.
subroutine, public atom_set_occupation(ostring, occupation, wfnocc, multiplicity)
Set occupation of atomic orbitals.
pure integer function, public get_maxl_occ(occupation)
Return the maximum orbital quantum number of occupied orbitals.
routines that build the integrals of the Vxc potential calculated for the atomic code
subroutine, public calculate_atom_ext_vxc(vxc, atom, lprint, xcmat)
ZMP subroutine for the scf external density from the external v_xc.
subroutine, public calculate_atom_zmp(ext_density, atom, lprint, xcmat)
ZMP subroutine for the calculation of the effective constraint potential in the atomic code.
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
Definition of mathematical constants and functions.
real(kind=dp), dimension(0:maxfac), parameter, public gamma1
real(kind=dp), parameter, public pi
real(kind=dp), dimension(-1:2 *maxfac+1), parameter, public dfac
Periodic Table related data definitions.
type(atom), dimension(0:nelem), public ptable
integer, parameter, public nelem
Definition of physical constants:
real(kind=dp), parameter, public bohr