13 atom_basis_type, atom_gthpot_type, atom_potential_type, atom_state, atom_type,
cgto_basis, &
40 #include "./base/base_uses.f90"
46 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'atom_output'
64 INTEGER,
INTENT(IN) :: zval
65 CHARACTER(len=*),
INTENT(IN) :: info
66 INTEGER,
INTENT(IN) :: iw
68 WRITE (iw,
'(/," ",A,T40,A," [",A,"]",T62,"Atomic number:",T78,I3,/)') &
69 adjustl(trim(info)), trim(
ptable(zval)%name), trim(
ptable(zval)%symbol), zval
83 TYPE(atom_state) :: state
84 INTEGER,
INTENT(IN) :: iw
86 CHARACTER(LEN=1),
DIMENSION(0:7),
PARAMETER :: &
87 label = (/
"S",
"P",
"D",
"F",
"G",
"H",
"I",
"K"/)
89 INTEGER :: j, l, mc, mlc, mlo, mm(0:
lmat), mo
92 WRITE (iw,
'(/,T2,A)')
"Electronic structure"
93 WRITE (iw,
'(T5,A,T71,F10.2)')
"Total number of core electrons", sum(state%core)
94 WRITE (iw,
'(T5,A,T71,F10.2)')
"Total number of valence electrons", sum(state%occ)
95 WRITE (iw,
'(T5,A,T71,F10.2)')
"Total number of electrons", sum(state%occ + state%core)
96 SELECT CASE (state%multiplicity)
98 WRITE (iw,
'(T5,A,T68,A)')
"Multiplicity",
"not specified"
100 WRITE (iw,
'(T5,A,T72,A)')
"Multiplicity",
"high spin"
102 WRITE (iw,
'(T5,A,T73,A)')
"Multiplicity",
"low spin"
104 WRITE (iw,
'(T5,A,T74,A)')
"Multiplicity",
"singlet"
106 WRITE (iw,
'(T5,A,T74,A)')
"Multiplicity",
"doublet"
108 WRITE (iw,
'(T5,A,T74,A)')
"Multiplicity",
"triplet"
110 WRITE (iw,
'(T5,A,T74,A)')
"Multiplicity",
"quartet"
112 WRITE (iw,
'(T5,A,T74,A)')
"Multiplicity",
"quintet"
114 WRITE (iw,
'(T5,A,T75,A)')
"Multiplicity",
"sextet"
116 WRITE (iw,
'(T5,A,T75,A)')
"Multiplicity",
"septet"
124 IF (state%multiplicity == -1)
THEN
125 DO l = 0, max(mlo, mlc)
126 mo = state%maxn_occ(l)
127 IF (sum(state%core(l, :)) == 0)
THEN
128 WRITE (iw,
'(A5,T10,10F6.2)') label(l), (state%occ(l, j), j=1, mo)
131 cpassert(sum(state%occ(l, 1:mc)) == 0)
132 WRITE (iw, advance=
"no", fmt=
'(A5,T9,A1,10F6.2)') label(l),
"[", (state%core(l, j), j=1, mc)
133 WRITE (iw, fmt=
'(A1,F5.2,10F6.2)')
"]", (state%occ(l, j), j=mc + 1, mc + mo)
137 WRITE (iw,
'(T5,A)')
"Alpha Electrons"
138 DO l = 0, max(mlo, mlc)
139 mo = state%maxn_occ(l)
140 IF (sum(state%core(l, :)) == 0)
THEN
141 WRITE (iw,
'(A5,T10,10F6.2)') label(l), (state%occa(l, j), j=1, mo)
144 WRITE (iw, advance=
"no", fmt=
'(A5,T9,A1,10F6.2)') label(l),
"[", (0.5_dp*state%core(l, j), j=1, mc)
145 WRITE (iw, fmt=
'(A1,F5.2,10F6.2)')
"]", (state%occa(l, j), j=1, mo)
148 WRITE (iw,
'(T5,A)')
"Beta Electrons"
149 DO l = 0, max(mlo, mlc)
150 mo = state%maxn_occ(l)
151 IF (sum(state%core(l, :)) == 0)
THEN
152 WRITE (iw,
'(A5,T10,10F6.2)') label(l), (state%occb(l, j), j=1, mo)
155 WRITE (iw, advance=
"no", fmt=
'(A5,T9,A1,10F6.2)') label(l),
"[", (0.5_dp*state%core(l, j), j=1, mc)
156 WRITE (iw, fmt=
'(A1,F5.2,10F6.2)')
"]", (state%occb(l, j), j=1, mo)
175 TYPE(atom_type) ::
atom
176 INTEGER,
INTENT(IN) :: iw
179 REAL(kind=
dp) :: drho
181 WRITE (iw,
'(/,A,T36,A,T61,F20.12)')
" Energy components [Hartree]", &
182 " Total Energy ::",
atom%energy%etot
183 WRITE (iw,
'(T36,A,T61,F20.12)')
" Band Energy ::",
atom%energy%eband
184 WRITE (iw,
'(T36,A,T61,F20.12)')
" Kinetic Energy ::",
atom%energy%ekin
185 WRITE (iw,
'(T36,A,T61,F20.12)')
"Potential Energy ::",
atom%energy%epot
186 IF (
atom%energy%ekin /= 0.0_dp)
THEN
187 WRITE (iw,
'(T36,A,T61,F20.12)')
" Virial (-V/T) ::", -
atom%energy%epot/
atom%energy%ekin
189 WRITE (iw,
'(T36,A,T61,F20.12)')
" Core Energy ::",
atom%energy%ecore
190 IF (
atom%energy%exc /= 0._dp) &
191 WRITE (iw,
'(T36,A,T61,F20.12)')
" XC Energy ::",
atom%energy%exc
192 WRITE (iw,
'(T36,A,T61,F20.12)')
" Coulomb Energy ::",
atom%energy%ecoulomb
193 IF (
atom%energy%eexchange /= 0._dp) &
194 WRITE (iw,
'(T34,A,T61,F20.12)')
"HF Exchange Energy ::",
atom%energy%eexchange
195 IF (
atom%potential%ppot_type /= no_pseudo)
THEN
196 WRITE (iw,
'(T20,A,T61,F20.12)')
" Total Pseudopotential Energy ::",
atom%energy%epseudo
197 WRITE (iw,
'(T20,A,T61,F20.12)')
" Local Pseudopotential Energy ::",
atom%energy%eploc
198 IF (
atom%energy%elsd /= 0._dp) &
199 WRITE (iw,
'(T20,A,T61,F20.12)')
" Local Spin-potential Energy ::",
atom%energy%elsd
200 WRITE (iw,
'(T20,A,T61,F20.12)')
" Nonlocal Pseudopotential Energy ::",
atom%energy%epnl
202 IF (
atom%potential%confinement)
THEN
203 WRITE (iw,
'(T36,A,T61,F20.12)')
" Confinement ::",
atom%energy%econfinement
206 IF (
atom%state%multiplicity == -1)
THEN
207 WRITE (iw,
'(/,A,T20,A,T30,A,T36,A,T49,A,T71,A,/)')
" Orbital energies", &
208 "State",
"L",
"Occupation",
"Energy[a.u.]",
"Energy[eV]"
209 DO l = 0,
atom%state%maxl_calc
210 n =
atom%state%maxn_calc(l)
212 WRITE (iw,
'(T23,I2,T30,I1,T36,F10.3,T46,F15.6,T66,F15.6)') &
213 i, l,
atom%state%occupation(l, i),
atom%orbitals%ener(i, l),
atom%orbitals%ener(i, l)*
evolt
215 IF (n > 0)
WRITE (iw, *)
218 WRITE (iw,
'(/,A,T20,A,T30,A,T36,A,T42,A,T55,A,T71,A,/)')
" Orbital energies", &
219 "State",
"Spin",
"L",
"Occupation",
"Energy[a.u.]",
"Energy[eV]"
220 DO l = 0,
atom%state%maxl_calc
221 n =
atom%state%maxn_calc(l)
223 WRITE (iw,
'(T23,I2,T29,A,T36,I1,T42,F10.3,T52,F15.6,T68,F13.6)') &
224 i,
"alpha", l,
atom%state%occa(l, i),
atom%orbitals%enera(i, l),
atom%orbitals%enera(i, l)*
evolt
227 WRITE (iw,
'(T23,I2,T29,A,T36,I1,T42,F10.3,T52,F15.6,T68,F13.6)') &
228 i,
" beta", l,
atom%state%occb(l, i),
atom%orbitals%enerb(i, l),
atom%orbitals%enerb(i, l)*
evolt
230 IF (n > 0)
WRITE (iw, *)
235 WRITE (iw,
'(/,A,T66,F15.6)')
" Total Electron Density at R=0: ", drho
248 INTEGER,
INTENT(IN) :: iter
249 REAL(
dp),
INTENT(IN) :: deps
250 TYPE(atom_type),
INTENT(IN) ::
atom
251 INTEGER,
INTENT(IN) :: iw
254 WRITE (iw,
'(/," ",79("*"),/,T33,"Integral",T48,"Integral",/,T3,A,T16,A,T33,A,T46,A,T69,A/," ",79("*"))') &
255 "Iteration",
"Convergence",
"rho diff.",
"rho*v_xc[au]",
"Energy[au]"
257 WRITE (iw,
'(T3,I9,T15,G13.6,T30,G13.6,T46,G13.6,T61,F20.12)') iter, deps,
atom%rho_diff_integral, &
272 INTEGER,
INTENT(IN) :: iter
273 REAL(
dp),
INTENT(IN) :: deps, etot
274 INTEGER,
INTENT(IN) :: iw
277 WRITE (iw,
'(/," ",79("*"),/,T19,A,T38,A,T70,A,/," ",79("*"))') &
278 "Iteration",
"Convergence",
"Energy [au]"
280 WRITE (iw,
'(T20,i8,T34,G14.6,T61,F20.12)') iter, deps, etot
294 INTEGER,
INTENT(IN) :: iw
295 CHARACTER(len=*) :: title
299 WRITE (iw,
'(/,A)') trim(title)
303 WRITE (iw,
'(/," ",21("*"),A,22("*"))')
" Geometrical Gaussian Type Orbitals "
304 WRITE (iw,
'(A,F15.8,T41,A,F15.8)')
" Initial exponent: ",
atom_basis%aval, &
307 WRITE (iw,
'(/," ",21("*"),A,21("*"))')
" Uncontracted Gaussian Type Orbitals "
313 WRITE (iw,
'(/,T2,A,(T30,I5,T51,F30.8))') &
316 WRITE (iw,
'(/,T2,A,(T30,I5,T51,F30.8))') &
319 WRITE (iw,
'(/,T2,A,(T30,I5,T51,F30.8))') &
322 WRITE (iw,
'(/,T2,A,(T30,I5,T51,F30.8))') &
325 WRITE (iw,
'(/,T2,A,(T30,I5,T51,F30.8))') &
330 WRITE (iw,
'(" ",79("*"))')
332 WRITE (iw,
'(/," ",22("*"),A,22("*"))')
" Contracted Gaussian Type Orbitals "
335 IF (l == 0)
WRITE (iw,
'(A)')
" s Functions"
336 IF (l == 1)
WRITE (iw,
'(A)')
" p Functions"
337 IF (l == 2)
WRITE (iw,
'(A)')
" d Functions"
338 IF (l == 3)
WRITE (iw,
'(A)')
" f Functions"
339 IF (l >= 3)
WRITE (iw,
'(A)')
" x Functions"
341 WRITE (iw,
'(F15.6,5(T21,6F10.6,/))') &
346 WRITE (iw,
'(" ",79("*"))')
348 WRITE (iw,
'(/," ",28("*"),A,29("*"))')
" Slater Type Orbitals "
365 WRITE (iw,
'(" ",79("*"))')
387 REAL(kind=
dp),
DIMENSION(:, :, 0:),
OPTIONAL :: wfn
389 INTEGER :: i, im, iw, l
390 REAL(kind=
dp) :: expzet, prefac, zeta
392 CALL open_file(file_name=
"OPT_BASIS", file_status=
"UNKNOWN", file_action=
"WRITE", unit_number=iw)
396 WRITE (iw,
'(/," ",21("*"),A,22("*"))')
" Geometrical Gaussian Type Orbitals "
397 WRITE (iw,
'(A,F15.8,T41,A,F15.8)')
" Initial exponent: ",
atom_basis%aval, &
400 WRITE (iw,
'(T3,A)')
"BASIS_TYPE GAUSSIAN"
406 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
409 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
412 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
415 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
418 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
426 WRITE (iw,
'(T3,A)')
"BASIS_TYPE SLATER"
431 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
433 WRITE (iw,
'(T3,A,60I3)') &
436 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
438 WRITE (iw,
'(T3,A,60I3)') &
441 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
443 WRITE (iw,
'(T3,A,60I3)') &
446 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
448 WRITE (iw,
'(T3,A,60I3)') &
451 WRITE (iw,
'(T3,A,(T15,F20.8,:," \"))') &
453 WRITE (iw,
'(T3,A,60I3)') &
464 IF (
PRESENT(wfn))
THEN
469 WRITE (iw,
'(/,T3,A)')
"ORBITAL COEFFICENTS (Quickstep normalization)"
470 im = min(6,
SIZE(wfn, 2))
473 WRITE (iw,
'(T3,A,I3)')
"L Quantum Number:", l
475 expzet = 0.25_dp*real(2*l + 3,
dp)
476 prefac = sqrt(
rootpi/2._dp**(l + 2)*
dfac(2*l + 1))
479 WRITE (iw,
'(T5,F14.8,2x,6F12.8)')
atom_basis%am(i, l), wfn(i, 1:im, l)*prefac/zeta
504 TYPE(atom_type) ::
atom
505 INTEGER,
INTENT(IN) :: iw
507 CHARACTER(len=160) :: shortform
508 CHARACTER(LEN=20) :: tmpstr
509 CHARACTER(len=:),
ALLOCATABLE :: reference
510 INTEGER :: ifun, il, meth, myfun, reltyp
512 TYPE(section_vals_type),
POINTER :: xc_fun, xc_fun_section, xc_section
514 NULLIFY (xc_fun, xc_fun_section, xc_section)
516 meth =
atom%method_type
518 xc_section =>
atom%xc_section
539 IF (iw > 0)
WRITE (iw, fmt=
"(/,' METHOD | Restricted Kohn-Sham Calculation')")
541 IF (iw > 0)
WRITE (iw, fmt=
"(/,' METHOD | Unrestricted Kohn-Sham Calculation')")
543 IF (iw > 0)
WRITE (iw, fmt=
"(/,' METHOD | Restricted Hartree-Fock Calculation')")
545 IF (iw > 0)
WRITE (iw, fmt=
"(/,' METHOD | Unrestricted Hartree-Fock Calculation')")
547 IF (iw > 0)
WRITE (iw, fmt=
"(/,' METHOD | Restricted Open-Shell Kohn-Sham Calculation')")
551 IF (
atom%do_zmp)
THEN
552 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Method on atomic radial density')")
553 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Lambda : ',F5.1)")
atom%lambda
554 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Reading external density : ',A20)")
atom%ext_file
556 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | The file is in the form of a density matrix')")
558 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | The file is in the form of a linear density')")
560 IF (
atom%doread)
THEN
561 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Restarting calculation from ',A20,' file if present')")
atom%zmp_restart_file
563 ELSE IF (
atom%read_vxc)
THEN
564 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Calculating density from external V_xc')")
565 IF (iw > 0)
WRITE (iw, fmt=
"(' ZMP | Reading external v_xc file : ',A20)")
atom%ext_vxc_file
568 IF (
atom%pp_calc)
THEN
569 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Nonrelativistic Calculation')")
571 reltyp =
atom%relativistic
577 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Nonrelativistic Calculation')")
579 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using ZORA(MP)')")
581 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using scaled ZORA(MP)')")
583 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using Douglas-Kroll 0th order')")
584 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using kietic energy scaling')")
586 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using Douglas-Kroll 1st order')")
587 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using Foldy-Wouthuysen transformation')")
589 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using Douglas-Kroll 2nd order')")
591 IF (iw > 0)
WRITE (iw, fmt=
"(' METHOD | Relativistic Calculation using Douglas-Kroll 3rd order')")
599 IF (iw > 0)
WRITE (iw, fmt=
"(' FUNCTIONAL| ROUTINE=',a)") trim(tmpstr)
606 IF (.NOT.
ASSOCIATED(xc_fun))
EXIT
608 ALLOCATE (
CHARACTER(LEN=libxc_get_reference_length(xc_fun, lsd)) :: reference)
610 ALLOCATE (
CHARACTER(LEN=20*default_string_length) :: reference)
613 WRITE (iw, fmt=
"(' FUNCTIONAL| ',a,':')") &
614 trim(xc_fun%section%name)
615 DO il = 1, len_trim(reference), 67
616 WRITE (iw, fmt=
"(' FUNCTIONAL| ',a67)") reference(il:)
618 DEALLOCATE (reference)
622 IF (iw > 0)
WRITE (iw, fmt=
"(' FUNCTIONAL| NO EXCHANGE-CORRELATION FUNCTIONAL USED.')")
640 TYPE(atom_potential_type) :: potential
641 INTEGER,
INTENT(IN) :: iw
643 CHARACTER(len=60) :: pline
644 INTEGER :: i, j, k, l
646 SELECT CASE (potential%ppot_type)
648 WRITE (iw,
'(/," ",28("*"),A,27("*"))')
" All Electron Potential "
650 WRITE (iw,
'(/," ",29("*"),A,29("*"))')
" GTH Pseudopotential "
651 WRITE (iw,
'(T10,A,T76,F5.1)')
" Core Charge ", potential%gth_pot%zion
652 WRITE (iw,
'(T10,A,T66,F15.6)')
" Rc ", potential%gth_pot%rc
653 WRITE (pline,
'(5F12.6)') (potential%gth_pot%cl(i), i=1, potential%gth_pot%ncl)
654 WRITE (iw,
'(T10,A,T21,A60)')
" C1 C2 ... ", adjustr(pline)
655 IF (potential%gth_pot%lpotextended)
THEN
656 DO k = 1, potential%gth_pot%nexp_lpot
657 WRITE (iw,
'(T10,A,F10.6,T38,A,4F10.6)')
" LPot: rc=", potential%gth_pot%alpha_lpot(k), &
658 "CX=", (potential%gth_pot%cval_lpot(i, k), i=1, potential%gth_pot%nct_lpot(k))
661 IF (potential%gth_pot%nlcc)
THEN
662 DO k = 1, potential%gth_pot%nexp_nlcc
663 WRITE (iw,
'(T10,A,F10.6,T38,A,4F10.6)')
" LSDPot: rc=", potential%gth_pot%alpha_nlcc(k), &
664 "CX=", (potential%gth_pot%cval_nlcc(i, k)*4.0_dp*
pi, i=1, potential%gth_pot%nct_nlcc(k))
667 IF (potential%gth_pot%lsdpot)
THEN
668 DO k = 1, potential%gth_pot%nexp_lsd
669 WRITE (iw,
'(T10,A,F10.6,T38,A,4F10.6)')
" LSDPot: rc=", potential%gth_pot%alpha_lsd(k), &
670 "CX=", (potential%gth_pot%cval_lsd(i, k), i=1, potential%gth_pot%nct_lsd(k))
674 IF (potential%gth_pot%nl(l) > 0)
THEN
675 WRITE (iw,
'(T10,A,T76,I5)')
" Angular momentum ", l
676 WRITE (iw,
'(T10,A,T66,F15.6)')
" Rcnl ", potential%gth_pot%rcnl(l)
677 WRITE (iw,
'(T10,A,T76,I5)')
" Nl ", potential%gth_pot%nl(l)
678 WRITE (pline,
'(5F12.6)') (potential%gth_pot%hnl(1, j, l), j=1, potential%gth_pot%nl(l))
679 WRITE (iw,
'(T10,A,T21,A60)')
" Hnl ", adjustr(pline)
680 DO i = 2, potential%gth_pot%nl(l)
681 WRITE (pline,
'(T21,5F12.6)') (potential%gth_pot%hnl(i, j, l), j=i, potential%gth_pot%nl(l))
682 WRITE (iw,
'(T21,A60)') adjustr(pline)
687 WRITE (iw,
'(/," ",29("*"),A,29("*"))')
" UPF Pseudopotential "
688 DO k = 1, potential%upf_pot%maxinfo
689 WRITE (iw,
'(A80)') potential%upf_pot%info(k)
692 WRITE (iw,
'(/," ",29("*"),A,29("*"))')
" SGP Pseudopotential "
693 WRITE (iw,
'(T10,A,T76,F5.1)')
" Core Charge ", potential%sgp_pot%zion
695 WRITE (iw,
'(/," ",26("*"),A,27("*"))')
" Effective Core Potential "
696 WRITE (iw,
'(T10,A,T76,F5.1)')
" Core Charge ", potential%ecp_pot%zion
697 DO k = 1, potential%ecp_pot%nloc
699 WRITE (iw,
'(T10,A,T40,I3,T49,2F16.8)')
" Local Potential ", potential%ecp_pot%nrloc(k), &
700 potential%ecp_pot%bloc(k), potential%ecp_pot%aloc(k)
702 WRITE (iw,
'(T40,I3,T49,2F16.8)') potential%ecp_pot%nrloc(k), &
703 potential%ecp_pot%bloc(k), potential%ecp_pot%aloc(k)
706 DO l = 0, potential%ecp_pot%lmax
707 WRITE (iw,
'(T10,A,I3)')
" ECP l-value ", l
708 DO k = 1, potential%ecp_pot%npot(l)
709 WRITE (iw,
'(T40,I3,T49,2F16.8)') potential%ecp_pot%nrpot(k, l), &
710 potential%ecp_pot%bpot(k, l), potential%ecp_pot%apot(k, l)
716 IF (potential%confinement)
THEN
717 IF (potential%conf_type ==
poly_conf)
THEN
718 WRITE (iw,
'(/,T10,A,T51,F12.6," * (R /",F6.2,")**",F6.2)') &
719 " Confinement Potential ", potential%acon, potential%rcon, potential%scon
721 WRITE (iw,
'(/,T10,A)')
" Confinement Potential s*F[(r-ron)/w] "
722 WRITE (iw,
'(T57,A,F12.6,A)')
"s =", potential%acon,
" Ha"
723 WRITE (iw,
'(T57,A,F12.6,A)')
"w =", potential%rcon,
" Bohr"
724 WRITE (iw,
'(T57,A,F12.6,A)')
"ron =", potential%scon,
" Bohr"
729 WRITE (iw,
'(/,T10,A)')
" No Confinement Potential is applied "
731 WRITE (iw,
'(" ",79("*"))')
745 TYPE(atom_gthpot_type),
INTENT(INOUT) :: gthpot
746 INTEGER,
INTENT(IN),
OPTIONAL :: iunit
748 INTEGER :: i, iw, j, k, n
750 IF (
PRESENT(iunit))
THEN
753 CALL open_file(file_name=
"GTH-PARAMETER", file_status=
"UNKNOWN", file_action=
"WRITE", unit_number=iw)
755 WRITE (iw,
'(A2,A)') gthpot%symbol, adjustl(trim(gthpot%pname))
756 WRITE (iw,
'(4I5)') gthpot%econf(0:3)
757 WRITE (iw,
'(F20.14,I8,5F20.14)') gthpot%rc, gthpot%ncl, (gthpot%cl(i), i=1, gthpot%ncl)
758 IF (gthpot%lpotextended)
THEN
759 WRITE (iw,
'(A,I5)')
" LPOT", gthpot%nexp_lpot
760 DO i = 1, gthpot%nexp_lpot
761 WRITE (iw,
'(F20.14,I8,5F20.14)') gthpot%alpha_lpot(i), gthpot%nct_lpot(i), &
762 (gthpot%cval_lpot(j, i), j=1, gthpot%nct_lpot(i))
765 IF (gthpot%lsdpot)
THEN
766 WRITE (iw,
'(A,I5)')
" LSD ", gthpot%nexp_lsd
767 DO i = 1, gthpot%nexp_lsd
768 WRITE (iw,
'(F20.14,I8,5F20.14)') gthpot%alpha_lsd(i), gthpot%nct_lsd(i), &
769 (gthpot%cval_lsd(j, i), j=1, gthpot%nct_lsd(i))
772 IF (gthpot%nlcc)
THEN
773 WRITE (iw,
'(A,I5)')
" NLCC ", gthpot%nexp_nlcc
774 DO i = 1, gthpot%nexp_nlcc
775 WRITE (iw,
'(F20.14,I8,5F20.14)') gthpot%alpha_nlcc(i), gthpot%nct_nlcc(i), &
776 (gthpot%cval_nlcc(j, i)*4.0_dp*
pi, j=1, gthpot%nct_nlcc(i))
781 IF (gthpot%nl(i) > 0)
THEN
788 WRITE (iw,
'(F20.14,I8,5F20.14)') gthpot%rcnl(i), gthpot%nl(i), (gthpot%hnl(1, k, i), k=1, gthpot%nl(i))
789 SELECT CASE (gthpot%nl(i))
791 WRITE (iw,
'(T49,F20.14)') gthpot%hnl(2, 2, i)
793 WRITE (iw,
'(T49,2F20.14)') gthpot%hnl(2, 2, i), gthpot%hnl(2, 3, i)
794 WRITE (iw,
'(T69,F20.14)') gthpot%hnl(3, 3, i)
796 DO j = 2, gthpot%nl(i)
797 WRITE (iw,
'(T29,5F20.14)') (gthpot%hnl(j, k, i), k=j, gthpot%nl(i))
801 IF (.NOT.
PRESENT(iunit))
CALL close_file(unit_number=iw)
813 TYPE(atom_type),
POINTER ::
atom
814 INTEGER,
INTENT(IN) :: iw
816 SELECT CASE (
atom%method_type)
820 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfn,
"", iw)
822 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfna,
"Alpha", iw)
823 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfnb,
"Beta", iw)
825 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfn,
"", iw)
827 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfna,
"Alpha", iw)
828 CALL atom_print_orbitals_helper(
atom,
atom%orbitals%wfnb,
"Beta", iw)
844 SUBROUTINE atom_print_orbitals_helper(atom, wfn, description, iw)
845 TYPE(atom_type),
POINTER ::
atom
846 REAL(kind=
dp),
DIMENSION(:, :, 0:),
INTENT(INOUT) :: wfn
847 CHARACTER(len=*),
INTENT(IN) :: description
848 INTEGER,
INTENT(IN) :: iw
850 INTEGER :: b, l, maxl, nb, nv, v
852 WRITE (iw,
'(/,A,A,A)')
" Atomic orbital expansion coefficients [", description,
"]"
854 maxl =
atom%state%maxl_calc
857 nb =
atom%basis%nbas(l)
858 nv =
atom%state%maxn_calc(l)
859 IF (nb > 0 .AND. nv > 0)
THEN
860 nv = min(nv,
SIZE(wfn, 2))
862 WRITE (iw,
'(/," ORBITAL L = ",I1," State = ",I3)') l, v
864 WRITE (iw,
'(" ",ES23.15)') wfn(b, v, l)
869 END SUBROUTINE atom_print_orbitals_helper
Routines that print various information about an atomic kind.
subroutine, public atom_print_energies(atom, iw)
Print energy components.
subroutine, public atom_print_basis(atom_basis, iw, title)
Print atomic basis set.
subroutine, public atom_print_iteration(iter, deps, etot, iw)
Print convergence information.
subroutine, public atom_print_basis_file(atom_basis, wfn)
Print the optimized atomic basis set into a file.
subroutine, public atom_write_pseudo_param(gthpot, iunit)
Print GTH pseudo-potential parameters.
subroutine, public atom_print_state(state, iw)
Print information about electronic state.
subroutine, public atom_print_zmp_iteration(iter, deps, atom, iw)
Printing of the atomic iterations when ZMP is active.
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.
Define the atom type and its sub types.
integer, parameter, public num_basis
integer, parameter, public cgto_basis
integer, parameter, public gto_basis
integer, parameter, public sto_basis
integer, parameter, public lmat
Some basic routines for atomic calculations.
pure integer function, dimension(0:lmat), public get_maxn_occ(occupation)
Return the maximum principal quantum number of occupied orbitals.
subroutine, public get_rho0(atom, rho0)
Calculate the total electron density at R=0.
pure integer function, public get_maxl_occ(occupation)
Return the maximum orbital quantum number of occupied orbitals.
Utility routines to open and close files. Tracking of preconnections.
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
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
Periodic Table related data definitions.
type(atom), dimension(0:nelem), public ptable
Definition of physical constants:
real(kind=dp), parameter, public evolt
subroutine, public xc_functional_get_info(functional, lsd, reference, shortform, needs, max_deriv, print_warn)
get the information about the given functional
calculates a functional from libxc and its derivatives
logical function, public libxc_check_existence_in_libxc(libxc_params)
This function checks whether a functional name belongs to LibXC.
integer function, public libxc_get_reference_length(libxc_params, lsd)
This function returns the maximum length of the reference string for a given LibXC functional.