52 #include "./base/base_uses.f90"
60 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_dftb_parameters'
62 REAL(dp),
PARAMETER :: slako_d0 = 1._dp
80 subsys_section, para_env)
81 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
82 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
83 TYPE(dftb_control_type),
INTENT(inout) :: dftb_control
84 TYPE(qs_dftb_pairpot_type),
DIMENSION(:, :), &
85 POINTER :: dftb_potential
86 TYPE(section_vals_type),
POINTER :: subsys_section
87 TYPE(mp_para_env_type),
POINTER :: para_env
89 CHARACTER(LEN=2) :: iel, jel
90 CHARACTER(LEN=6) :: cspline
91 CHARACTER(LEN=default_path_length) :: file_name
92 CHARACTER(LEN=default_path_length),
ALLOCATABLE, &
93 DIMENSION(:, :) :: sk_files
94 CHARACTER(LEN=default_string_length) :: iname, jname, name_a, name_b, skfn
95 INTEGER :: ikind, isp, jkind, k, l, l1, l2, llm, &
96 lmax, lmax_a, lmax_b, lp, m, n_urpoly, &
97 ngrd, nkind, output_unit, runit, &
99 LOGICAL :: at_end, found, ldum, search, sklist
100 REAL(dp) :: da, db, dgrd, dij, energy, eps_disp, ra, &
101 radmax, rb, rcdisp, rmax6, s_cut, xij, &
103 REAL(dp),
ALLOCATABLE,
DIMENSION(:, :) :: fmat, scoeff, smat, spxr
104 REAL(dp),
DIMENSION(0:3) :: eta, occupation, skself
105 REAL(dp),
DIMENSION(10) :: fwork, swork, uwork
106 REAL(dp),
DIMENSION(1:2) :: surr
107 REAL(dp),
DIMENSION(1:3) :: srep
108 TYPE(cp_logger_type),
POINTER :: logger
109 TYPE(qs_dftb_atom_type),
POINTER :: dftb_atom_a, dftb_atom_b
115 "PRINT%KINDS/BASIS_SET"),
cp_p_file))
THEN
117 "PRINT%KINDS", extension=
".Log")
118 IF (output_unit > 0)
THEN
119 WRITE (output_unit,
"(/,A)")
" DFTB| A set of relativistic DFTB "// &
120 "parameters for material sciences."
121 WRITE (output_unit,
"(A)")
" DFTB| J. Frenzel, N. Jardillier, A.F. Oliveira,"// &
122 " T. Heine, G. Seifert"
123 WRITE (output_unit,
"(A)")
" DFTB| TU Dresden, 2002-2007"
124 WRITE (output_unit,
"(/,A)")
" DFTB| Non-SCC parameters "
125 WRITE (output_unit,
"(A,T25,A)")
" DFTB| C,H :", &
126 " D. Porezag et al, PRB 51 12947 (1995)"
127 WRITE (output_unit,
"(A,T25,A)")
" DFTB| B,N :", &
128 " J. Widany et al, PRB 53 4443 (1996)"
129 WRITE (output_unit,
"(A,T25,A)")
" DFTB| Li,Na,K,Cl :", &
130 " S. Hazebroucq et al, JCP 123 134510 (2005)"
131 WRITE (output_unit,
"(A,T25,A)")
" DFTB| F :", &
132 " T. Heine et al, JCSoc-Perkins Trans 2 707 (1999)"
133 WRITE (output_unit,
"(A,T25,A)")
" DFTB| Mo,S :", &
134 " G. Seifert et al, PRL 85 146 (2000)"
135 WRITE (output_unit,
"(A,T25,A)")
" DFTB| P :", &
136 " G. Seifert et al, EPS 16 341 (2001)"
137 WRITE (output_unit,
"(A,T25,A)")
" DFTB| Sc,N,C :", &
138 " M. Krause et al, JCP 115 6596 (2001)"
144 sklist = (dftb_control%sk_file_list /=
"")
146 nkind =
SIZE(atomic_kind_set)
147 ALLOCATE (sk_files(nkind, nkind))
149 ALLOCATE (dftb_potential(nkind, nkind))
153 CALL get_atomic_kind(atomic_kind_set(ikind), name=iname, element_symbol=iel)
158 CALL get_atomic_kind(atomic_kind_set(jkind), name=jname, element_symbol=jel)
163 DO k = 1,
SIZE(dftb_control%sk_pair_list, 2)
164 name_a = trim(dftb_control%sk_pair_list(1, k))
165 name_b = trim(dftb_control%sk_pair_list(2, k))
168 IF ((iname == name_a .AND. jname == name_b))
THEN
169 sk_files(ikind, jkind) = trim(dftb_control%sk_file_path)//
"/"// &
170 trim(dftb_control%sk_pair_list(3, k))
175 IF (.NOT. found .AND. sklist)
THEN
176 file_name = trim(dftb_control%sk_file_path)//
"/"// &
177 trim(dftb_control%sk_file_list)
179 TYPE(cp_parser_type) :: parser
185 CALL parser_get_object(parser, name_a, lower_to_upper=.true.)
186 CALL parser_get_object(parser, name_b, lower_to_upper=.true.)
188 IF ((iname == name_a .AND. jname == name_b))
THEN
189 CALL parser_get_object(parser, skfn, string_length=8)
190 sk_files(ikind, jkind) = trim(dftb_control%sk_file_path)//
"/"// &
196 IF ((iel == name_a .AND. jel == name_b))
THEN
197 CALL parser_get_object(parser, skfn, string_length=8)
198 sk_files(ikind, jkind) = trim(dftb_control%sk_file_path)//
"/"// &
208 CALL cp_abort(__location__, &
209 "Failure in assigning KINDS <"//trim(iname)//
"> and <"//trim(jname)// &
210 "> to a DFTB interaction pair!")
218 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
219 IF (.NOT.
ASSOCIATED(dftb_atom_a))
THEN
221 CALL set_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
228 CALL get_qs_kind(qs_kind_set(jkind), dftb_parameter=dftb_atom_b)
230 IF (output_unit > 0)
THEN
231 WRITE (output_unit,
"(A,T30,A50)")
" DFTB| Reading parameter file ", &
232 adjustr(trim(sk_files(jkind, ikind)))
237 IF (para_env%is_source())
THEN
239 CALL open_file(file_name=sk_files(jkind, ikind), unit_number=runit)
241 READ (runit, fmt=*,
END=1, err=1) dgrd, ngrd
249 READ (runit, fmt=*,
END=1, err=1) skself(2:0:-1), energy, &
250 eta(2:0:-1), occupation(2:0:-1)
252 READ (runit, fmt=*,
END=1, err=1) uwork(1:10)
254 IF (dot_product(uwork(2:10), uwork(2:10)) >= 1.e-12_dp)
THEN
257 IF (abs(uwork(k)) >= 1.e-12_dp) n_urpoly = k
262 IF (n_urpoly < 2) n_urpoly = 0
265 CALL para_env%bcast(n_urpoly)
266 CALL para_env%bcast(uwork)
267 CALL para_env%bcast(ngrd)
268 CALL para_env%bcast(dgrd)
270 CALL para_env%bcast(skself)
271 CALL para_env%bcast(energy)
272 CALL para_env%bcast(eta)
273 CALL para_env%bcast(occupation)
275 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, &
276 z=z, zeff=sum(occupation), defined=.true., &
277 skself=skself, energy=energy, eta=eta, occupation=occupation)
280 ALLOCATE (fmat(ngrd, 10))
281 ALLOCATE (smat(ngrd, 10))
282 IF (para_env%is_source())
THEN
284 READ (runit, fmt=*,
END=1, err=1) fwork(1:10), swork(1:10)
285 fmat(k, 1:10) = fwork(1:10)
286 smat(k, 1:10) = swork(1:10)
289 CALL para_env%bcast(fmat)
290 CALL para_env%bcast(smat)
303 CALL get_qs_kind(qs_kind_set(ikind), lmax_dftb=lmax)
321 IF ((abs(skself(l)) > 0._dp) .OR. &
322 (sum(abs(fmat(ngrd/10:ngrd, lp))) > 0._dp)) lmax = l
328 CALL cp_abort(__location__,
"Maximum L allowed is d. "// &
329 "Use KIND/LMAX_DFTB to set smaller values if needed.")
332 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, &
333 lmax=lmax, natorb=(lmax + 1)**2)
336 IF (n_urpoly == 0)
THEN
337 IF (para_env%is_source())
THEN
341 READ (runit, fmt=
'(A6)',
END=1, err=1) cspline
342 IF (cspline ==
'Spline')
THEN
345 READ (runit, fmt=*,
END=1, err=1) spdim, s_cut
346 ALLOCATE (spxr(spdim, 2))
347 ALLOCATE (scoeff(spdim, 4))
349 READ (runit, fmt=*,
END=1, err=1) srep(1:3)
350 DO isp = 1, spdim - 1
352 READ (runit, fmt=*,
END=1, err=1) spxr(isp, 1:2), scoeff(isp, 1:4)
355 READ (runit, fmt=*,
END=1, err=1) spxr(spdim, 1:2), scoeff(spdim, 1:4), surr(1:2)
361 IF (para_env%is_source())
THEN
362 CALL close_file(unit_number=runit)
365 CALL para_env%bcast(spdim)
366 IF (spdim > 0 .AND. (.NOT. para_env%is_source()))
THEN
367 ALLOCATE (spxr(spdim, 2))
368 ALLOCATE (scoeff(spdim, 4))
371 CALL para_env%bcast(spxr)
372 CALL para_env%bcast(scoeff)
373 CALL para_env%bcast(surr)
374 CALL para_env%bcast(srep)
375 CALL para_env%bcast(s_cut)
380 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_a, lmax=lmax_a)
381 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_b, lmax=lmax_b)
383 DO l1 = 0, max(lmax_a, lmax_b)
384 DO l2 = 0, min(l1, lmax_a, lmax_b)
390 CALL qs_dftb_pairpot_create(dftb_potential(ikind, jkind), &
394 dftb_potential(ikind, jkind)%n_urpoly = n_urpoly
395 dftb_potential(ikind, jkind)%urep_cut = uwork(10)
396 dftb_potential(ikind, jkind)%urep(:) = 0._dp
397 dftb_potential(ikind, jkind)%urep(1) = uwork(10)
398 dftb_potential(ikind, jkind)%urep(2:n_urpoly) = uwork(2:n_urpoly)
401 dftb_potential(ikind, jkind)%dgrd = dgrd
402 CALL skreorder(fmat, lmax_a, lmax_b)
403 dftb_potential(ikind, jkind)%fmat(:, 1:llm) = fmat(:, 1:llm)
404 CALL skreorder(smat, lmax_a, lmax_b)
405 dftb_potential(ikind, jkind)%smat(:, 1:llm) = smat(:, 1:llm)
406 dftb_potential(ikind, jkind)%ngrdcut = ngrd + int(slako_d0/dgrd)
409 dftb_potential(ikind, jkind)%s_cut = s_cut
410 dftb_potential(ikind, jkind)%srep = srep
411 dftb_potential(ikind, jkind)%spxr = spxr
412 dftb_potential(ikind, jkind)%scoeff = scoeff
413 dftb_potential(ikind, jkind)%surr = surr
427 CALL get_atomic_kind(atomic_kind_set(ikind), z=z, name=iname)
428 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
430 IF (.NOT.
ASSOCIATED(dftb_atom_a))
THEN
431 CALL allocate_dftb_atom_param(dftb_atom_a)
432 CALL set_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
437 IF (ikind == jkind) cycle
438 CALL get_atomic_kind(atomic_kind_set(jkind), name=jname)
439 CALL get_qs_kind(qs_kind_set(jkind), dftb_parameter=dftb_atom_b)
441 IF (output_unit > 0)
THEN
442 WRITE (output_unit,
"(A,T30,A50)")
" DFTB| Reading parameter file ", &
443 adjustr(trim(sk_files(ikind, jkind)))
448 IF (para_env%is_source())
THEN
449 runit = get_unit_number()
450 CALL open_file(file_name=sk_files(ikind, jkind), unit_number=runit)
452 READ (runit, fmt=*,
END=1, err=1) dgrd, ngrd
459 IF (ikind == jkind)
THEN
461 READ (runit, fmt=*,
END=1, err=1) skself(2:0:-1), energy, &
462 eta(2:0:-1), occupation(2:0:-1)
465 READ (runit, fmt=*,
END=1, err=1) uwork(1:10)
467 IF (dot_product(uwork(2:10), uwork(2:10)) >= 1.e-12_dp)
THEN
470 IF (abs(uwork(k)) >= 1.e-12_dp) n_urpoly = k
475 IF (n_urpoly < 2) n_urpoly = 0
478 CALL para_env%bcast(n_urpoly)
479 CALL para_env%bcast(uwork)
480 CALL para_env%bcast(ngrd)
481 CALL para_env%bcast(dgrd)
484 ALLOCATE (fmat(ngrd, 10))
485 ALLOCATE (smat(ngrd, 10))
486 IF (para_env%is_source())
THEN
488 READ (runit, fmt=*,
END=1, err=1) fwork(1:10), swork(1:10)
489 fmat(k, 1:10) = fwork(1:10)
490 smat(k, 1:10) = swork(1:10)
493 CALL para_env%bcast(fmat)
494 CALL para_env%bcast(smat)
497 IF (n_urpoly == 0)
THEN
498 IF (para_env%is_source())
THEN
502 READ (runit, fmt=
'(A6)',
END=1, err=1) cspline
503 IF (cspline ==
'Spline')
THEN
506 READ (runit, fmt=*,
END=1, err=1) spdim, s_cut
507 ALLOCATE (spxr(spdim, 2))
508 ALLOCATE (scoeff(spdim, 4))
510 READ (runit, fmt=*,
END=1, err=1) srep(1:3)
511 DO isp = 1, spdim - 1
513 READ (runit, fmt=*,
END=1, err=1) spxr(isp, 1:2), scoeff(isp, 1:4)
516 READ (runit, fmt=*,
END=1, err=1) spxr(spdim, 1:2), scoeff(spdim, 1:4), surr(1:2)
522 IF (para_env%is_source())
THEN
523 CALL close_file(unit_number=runit)
526 CALL para_env%bcast(spdim)
527 IF (spdim > 0 .AND. (.NOT. para_env%is_source()))
THEN
528 ALLOCATE (spxr(spdim, 2))
529 ALLOCATE (scoeff(spdim, 4))
532 CALL para_env%bcast(spxr)
533 CALL para_env%bcast(scoeff)
534 CALL para_env%bcast(surr)
535 CALL para_env%bcast(srep)
536 CALL para_env%bcast(s_cut)
541 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_a, lmax=lmax_a)
542 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_b, lmax=lmax_b)
544 DO l1 = 0, max(lmax_a, lmax_b)
545 DO l2 = 0, min(l1, lmax_a, lmax_b)
551 CALL qs_dftb_pairpot_create(dftb_potential(ikind, jkind), &
555 dftb_potential(ikind, jkind)%n_urpoly = n_urpoly
556 dftb_potential(ikind, jkind)%urep_cut = uwork(10)
557 dftb_potential(ikind, jkind)%urep(:) = 0._dp
558 dftb_potential(ikind, jkind)%urep(1) = uwork(10)
559 dftb_potential(ikind, jkind)%urep(2:n_urpoly) = uwork(2:n_urpoly)
562 dftb_potential(ikind, jkind)%dgrd = dgrd
563 CALL skreorder(fmat, lmax_a, lmax_b)
564 dftb_potential(ikind, jkind)%fmat(:, 1:llm) = fmat(:, 1:llm)
565 CALL skreorder(smat, lmax_a, lmax_b)
566 dftb_potential(ikind, jkind)%smat(:, 1:llm) = smat(:, 1:llm)
567 dftb_potential(ikind, jkind)%ngrdcut = ngrd + int(slako_d0/dgrd)
570 dftb_potential(ikind, jkind)%s_cut = s_cut
571 dftb_potential(ikind, jkind)%srep = srep
572 dftb_potential(ikind, jkind)%spxr = spxr
573 dftb_potential(ikind, jkind)%scoeff = scoeff
574 dftb_potential(ikind, jkind)%surr = surr
587 DEALLOCATE (sk_files)
590 IF (dftb_control%dispersion)
THEN
592 IF (dftb_control%dispersion_type == dispersion_uff)
THEN
593 file_name = trim(dftb_control%sk_file_path)//
"/"// &
594 trim(dftb_control%uff_force_field)
596 TYPE(cp_parser_type) :: parser
598 CALL get_atomic_kind(atomic_kind_set(ikind), name=iname)
599 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
602 CALL parser_create(parser, file_name, para_env=para_env)
606 CALL parser_get_next_line(parser, 1, at_end)
608 CALL parser_get_object(parser, name_a)
610 IF (name_a(1:1) ==
'"') name_a(1:m) = name_a(2:m + 1)
611 IF (name_a(1:m) == trim(iname))
THEN
612 CALL parser_get_object(parser, rb)
613 CALL parser_get_object(parser, rb)
614 CALL parser_get_object(parser, ra)
615 CALL parser_get_object(parser, da)
619 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, name=iname, xi=ra, di=da)
623 CALL parser_release(parser)
632 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
633 radmax = (dftb_potential(ikind, ikind)%ngrdcut + 1)* &
634 dftb_potential(ikind, ikind)%dgrd*0.5_dp
635 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, cutoff=radmax)
638 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
639 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_a, cutoff=ra)
641 CALL get_qs_kind(qs_kind_set(jkind), dftb_parameter=dftb_atom_b)
642 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_b, cutoff=rb)
643 radmax = (dftb_potential(ikind, jkind)%ngrdcut + 1)* &
644 dftb_potential(ikind, jkind)%dgrd
645 IF (ra + rb < radmax)
THEN
646 ra = ra + (radmax - ra - rb)*0.5_dp
647 rb = rb + (radmax - ra - rb)*0.5_dp
648 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, cutoff=ra)
649 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_b, cutoff=rb)
656 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
657 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_a, zeff=zeff)
658 CALL set_potential(potential=qs_kind_set(ikind)%all_potential, &
659 zeff=zeff, zeff_correction=0.0_dp)
663 IF (dftb_control%dftb3_diagonal)
THEN
665 CALL get_qs_kind(qs_kind_set(ikind), dftb3_param=db)
666 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
667 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, dudq=db)
672 IF (dftb_control%dispersion)
THEN
673 IF (dftb_control%dispersion_type == dispersion_uff)
THEN
674 eps_disp = dftb_control%eps_disp
676 CALL get_qs_kind(qs_kind_set(ikind), dftb_parameter=dftb_atom_a)
677 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_a, xi=ra, di=da)
680 CALL get_qs_kind(qs_kind_set(jkind), dftb_parameter=dftb_atom_b)
681 CALL get_dftb_atom_param(dftb_parameter=dftb_atom_b, xi=rb, di=db)
684 dftb_potential(ikind, jkind)%xij = xij
685 dftb_potential(ikind, jkind)%dij = dij
686 dftb_potential(ikind, jkind)%x0ij = xij*(0.5_dp**(1.0_dp/6.0_dp))
687 dftb_potential(ikind, jkind)%a = dij*396.0_dp/25.0_dp
688 dftb_potential(ikind, jkind)%b = &
689 dij/(xij**5)*672.0_dp*2.0_dp**(5.0_dp/6.0_dp)/25.0_dp
690 dftb_potential(ikind, jkind)%c = &
691 -dij/(xij**10)*2.0_dp**(2.0_dp/3.0_dp)*552.0_dp/25.0_dp
692 rmax6 = ((8._dp*pi*dij/eps_disp)*xij**6)**0.25_dp
693 rcdisp = max(rcdisp, rmax6*0.5_dp)
695 CALL set_dftb_atom_param(dftb_parameter=dftb_atom_a, rcdisp=rcdisp)
730 SUBROUTINE skreorder(xmat, la, lb)
731 REAL(dp),
DIMENSION(:, :),
INTENT(INOUT) :: xmat
732 INTEGER,
INTENT(IN) :: la, lb
734 INTEGER :: i, l1, l2, llm, m
735 REAL(dp) :: skllm(0:3, 0:3, 0:3)
737 DO i = 1,
SIZE(xmat, 1)
739 skllm(0, 0, 0) = xmat(i, 10)
740 skllm(1, 0, 0) = xmat(i, 9)
741 skllm(2, 0, 0) = xmat(i, 8)
742 skllm(1, 1, 1) = xmat(i, 7)
743 skllm(1, 1, 0) = xmat(i, 6)
744 skllm(2, 1, 1) = xmat(i, 5)
745 skllm(2, 1, 0) = xmat(i, 4)
746 skllm(2, 2, 2) = xmat(i, 3)
747 skllm(2, 2, 1) = xmat(i, 2)
748 skllm(2, 2, 0) = xmat(i, 1)
750 DO l1 = 0, max(la, lb)
751 DO l2 = 0, min(l1, la, lb)
754 xmat(i, llm) = skllm(l1, l2, m)
760 END SUBROUTINE skreorder
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
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.
integer function, public get_unit_number(file_name)
Returns the first logical unit that is not preconnected.
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.
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,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
Utility routines to read data from files. Kept as close as possible to the old parser because.
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 ...
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.
Definition of the atomic potential types.
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
integer, parameter, public default_path_length
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
Interface to the message passing library MPI.
Definition of physical constants:
real(kind=dp), parameter, public kcalmol
real(kind=dp), parameter, public angstrom
logical function, public qmmm_ff_precond_only_qm(id1, id2, id3, id4, is_link)
This function handles the atom names and modifies the "_QM_" prefix, in order to find the parameters ...
subroutine, public qs_dftb_param_init(atomic_kind_set, qs_kind_set, dftb_control, dftb_potential, subsys_section, para_env)
...
Definition of the DFTB parameter types.
subroutine, public qs_dftb_pairpot_init(pairpot)
...
subroutine, public qs_dftb_pairpot_create(pairpot, ngrd, llm, spdim)
...
Working with the DFTB parameter types.
subroutine, public set_dftb_atom_param(dftb_parameter, name, typ, defined, z, zeff, natorb, lmax, skself, occupation, eta, energy, cutoff, xi, di, rcdisp, dudq)
...
subroutine, public get_dftb_atom_param(dftb_parameter, name, typ, defined, z, zeff, natorb, lmax, skself, occupation, eta, energy, cutoff, xi, di, rcdisp, dudq)
...
subroutine, public allocate_dftb_atom_param(dftb_parameter)
...
Define the quickstep kind type and their sub types.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_r3d_rs_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, U_of_dft_plus_u, J_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, J0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public set_qs_kind(qs_kind, paw_atom, ghost, floating, hard_radius, hard0_radius, covalent_radius, vdw_radius, lmax_rho0, zeff, no_optimize, dispersion, u_minus_j, reltmat, dftb_parameter, xtb_parameter, elec_conf, pao_basis_size)
Set the components of an atomic kind data set.
Utilities for string manipulations.
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.