67 CHARACTER(len=*),
PARAMETER :: routinen =
'atom_pseudo_opt'
69 CHARACTER(LEN=2) :: elem
70 CHARACTER(LEN=default_string_length), &
71 DIMENSION(:),
POINTER :: tmpstringlist
72 INTEGER :: ads, do_eric, do_erie, handle, i, im, &
73 in, iw, k, l, maxl, mb, method, mo, &
74 n_meth, n_rep, nr_gh, reltyp, zcore, &
76 INTEGER,
DIMENSION(0:lmat) :: maxn
77 INTEGER,
DIMENSION(:),
POINTER :: cn
78 LOGICAL :: do_gh, eri_c, eri_e, pp_calc
79 REAL(kind=
dp) :: ne, nm
80 REAL(kind=
dp),
DIMENSION(0:lmat, 10) :: pocc
85 TYPE(
atom_p_type),
DIMENSION(:, :),
POINTER :: atom_info, atom_refs
90 opt_section, potential_section, &
91 powell_section, xc_section
93 CALL timeset(routinen, handle)
100 IF (
ptable(i)%symbol == elem)
THEN
105 IF (zz /= 1) zval = zz
108 ALLOCATE (ae_basis, pp_basis)
110 NULLIFY (ae_basis%grid)
112 NULLIFY (pp_basis%grid)
119 IF (iw > 0)
CALL atom_print_info(zval,
"Atomic Energy Calculation", iw)
129 NULLIFY (potential_section)
131 ALLOCATE (ae_pot, p_pot)
134 IF (.NOT. p_pot%confinement .AND. .NOT. ae_pot%confinement)
THEN
136 p_pot%confinement = .true.
141 p_pot%rcon = (2._dp*
ptable(zval)%covalent_radius*
bohr)**2
142 ae_pot%confinement = .true.
147 ae_pot%rcon = (2._dp*
ptable(zval)%covalent_radius*
bohr)**2
164 DO in = 1, min(
SIZE(cn), 4)
165 maxn(in - 1) = cn(in)
168 maxn(in) = min(maxn(in), ae_basis%nbas(in))
182 ALLOCATE (ae_int, pp_int)
184 ALLOCATE (atom_info(n_rep, n_meth), atom_refs(n_rep, n_meth))
188 WRITE (iw,
'(/," ",79("*"))')
189 WRITE (iw,
'(" ",26("*"),A,25("*"))')
" Calculate Reference States "
190 WRITE (iw,
'(" ",79("*"))')
197 NULLIFY (atom_info(in, im)%atom, atom_refs(in, im)%atom)
201 atom_info(in, im)%atom%optimization = optimization
202 atom_refs(in, im)%atom%optimization = optimization
204 atom_info(in, im)%atom%z = zval
205 atom_refs(in, im)%atom%z = zval
207 atom_info(in, im)%atom%xc_section => xc_section
208 atom_refs(in, im)%atom%xc_section => xc_section
210 ALLOCATE (state, statepp)
214 c_vals=tmpstringlist)
216 pp_calc = index(tmpstringlist(1),
"CORE") /= 0
217 cpassert(.NOT. pp_calc)
224 state%maxl_calc = max(maxl, state%maxl_occ)
225 state%maxl_calc = min(
lmat, state%maxl_calc)
227 DO k = 0, state%maxl_calc
229 IF (state%maxn_occ(k) == 0) ads = 1
230 state%maxn_calc(k) = max(maxn(k), state%maxn_occ(k) + ads)
231 state%maxn_calc(k) = min(state%maxn_calc(k), ae_basis%nbas(k))
234 CALL set_atom(atom_refs(in, im)%atom, zcore=zval, pp_calc=.false.)
236 IF (state%multiplicity /= -1)
THEN
241 nm = real((2*l + 1), kind=
dp)
243 ne = state%occupation(l, k)
244 IF (ne == 0._dp)
THEN
246 ELSEIF (ne == 2._dp*nm)
THEN
247 state%occa(l, k) = nm
248 state%occb(l, k) = nm
249 ELSEIF (state%multiplicity == -2)
THEN
250 state%occa(l, k) = min(ne, nm)
251 state%occb(l, k) = max(0._dp, ne - nm)
253 state%occa(l, k) = 0.5_dp*(ne + state%multiplicity - 1._dp)
254 state%occb(l, k) = ne - state%occa(l, k)
263 zcore = zval - nint(sum(statepp%core))
264 CALL set_atom(atom_info(in, im)%atom, zcore=zcore, pp_calc=.true.)
266 statepp%occ = state%occ - statepp%core
267 statepp%occupation = 0._dp
271 IF (statepp%occ(l, i) /= 0._dp)
THEN
273 statepp%occupation(l, k) = state%occ(l, i)
274 IF (state%multiplicity /= -1)
THEN
275 statepp%occa(l, k) = state%occa(l, i) - statepp%core(l, i)/2
276 statepp%occb(l, k) = state%occb(l, i) - statepp%core(l, i)/2
284 statepp%maxl_calc = state%maxl_calc
285 statepp%maxn_calc = 0
287 DO k = 0, statepp%maxl_calc
288 statepp%maxn_calc(k) = state%maxn_calc(k) - maxn(k)
289 statepp%maxn_calc(k) = min(statepp%maxn_calc(k), pp_basis%nbas(k))
291 statepp%multiplicity = state%multiplicity
295 CALL set_atom(atom_info(in, im)%atom, method_type=method)
296 CALL set_atom(atom_refs(in, im)%atom, method_type=method, relativistic=reltyp)
300 CALL atom_int_setup(pp_int, pp_basis, potential=p_pot, eri_coulomb=eri_c, eri_exchange=eri_e)
302 NULLIFY (pp_int%tzora, pp_int%hdkh)
306 CALL set_atom(atom_info(in, im)%atom, basis=pp_basis, integrals=pp_int, potential=p_pot)
307 statepp%maxn_calc(:) = min(statepp%maxn_calc(:), pp_basis%nbas(:))
308 cpassert(all(state%maxn_calc(:) >= state%maxn_occ))
313 eri_coulomb=eri_c, eri_exchange=eri_e)
319 CALL set_atom(atom_refs(in, im)%atom, basis=ae_basis, integrals=ae_int, potential=ae_pot)
320 state%maxn_calc(:) = min(state%maxn_calc(:), ae_basis%nbas(:))
321 cpassert(all(state%maxn_calc(:) >= state%maxn_occ))
323 CALL set_atom(atom_info(in, im)%atom, coulomb_integral_type=do_eric, &
324 exchange_integral_type=do_erie)
325 CALL set_atom(atom_refs(in, im)%atom, coulomb_integral_type=do_eric, &
326 exchange_integral_type=do_erie)
327 atom_info(in, im)%atom%hfx_pot%do_gh = do_gh
328 atom_info(in, im)%atom%hfx_pot%nr_gh = nr_gh
329 atom_refs(in, im)%atom%hfx_pot%do_gh = do_gh
330 atom_refs(in, im)%atom%hfx_pot%nr_gh = nr_gh
332 CALL set_atom(atom_info(in, im)%atom, state=statepp)
334 mo = maxval(statepp%maxn_calc)
335 mb = maxval(atom_info(in, im)%atom%basis%nbas)
337 CALL set_atom(atom_info(in, im)%atom, orbitals=orbitals)
339 CALL set_atom(atom_refs(in, im)%atom, state=state)
341 mo = maxval(state%maxn_calc)
342 mb = maxval(atom_refs(in, im)%atom%basis%nbas)
344 CALL set_atom(atom_refs(in, im)%atom, orbitals=orbitals)
346 IF (
atom_consistent_method(atom_refs(in, im)%atom%method_type, atom_refs(in, im)%atom%state%multiplicity))
THEN
361 WRITE (iw,
'(/," ",79("*"))')
362 WRITE (iw,
'(" ",21("*"),A,21("*"))')
" Optimize Pseudopotential Parameters "
363 WRITE (iw,
'(" ",79("*"))')
404 DEALLOCATE (atom_info, atom_refs)
406 DEALLOCATE (ae_pot, p_pot, ae_basis, pp_basis, ae_int, pp_int)
408 CALL timestop(handle)