(git:b17b328)
Loading...
Searching...
No Matches
qs_environment.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \par History
10!> - Merged with the Quickstep MODULE method_specification (17.01.2002,MK)
11!> - USE statements cleaned, added
12!> (25.09.2002,MK)
13!> - Added more LSD structure (01.2003,Joost VandeVondele)
14!> - New molecule data types introduced (Sep. 2003,MK)
15!> - Cleaning; getting rid of pnode (02.10.2003,MK)
16!> - Sub-system setup added (08.10.2003,MK)
17!> \author MK (18.05.2000)
18! **************************************************************************************************
30 USE bibliography, ONLY: iannuzzi2006,&
32 cite_reference,&
34 USE cell_types, ONLY: cell_type
44 USE cp_control_utils, ONLY: &
52 USE cp_output_handling, ONLY: cp_p_file,&
62 USE ec_environment, ONLY: ec_env_create,&
81 USE gamma, ONLY: init_md_ftable
84 USE header, ONLY: dftb_header,&
85 qs_header,&
86 se_header,&
91 USE input_constants, ONLY: &
106 USE kinds, ONLY: default_string_length,&
107 dp
111 USE kpoint_types, ONLY: get_kpoint_info,&
119 USE machine, ONLY: m_flush
120 USE mathconstants, ONLY: pi
125 USE mp2_setup, ONLY: read_mp2_section
126 USE mp2_types, ONLY: mp2_env_create,&
135 USE pw_env_types, ONLY: pw_env_type
152 USE qs_gcp_types, ONLY: qs_gcp_type
153 USE qs_gcp_utils, ONLY: qs_gcp_env_set,&
166 USE qs_kind_types, ONLY: &
170 USE qs_ks_types, ONLY: qs_ks_env_create,&
174 USE qs_mo_types, ONLY: allocate_mo_set,&
177 USE qs_rho0_methods, ONLY: init_rho0
182 USE qs_subsys_types, ONLY: qs_subsys_get,&
205 USE tblite_interface, ONLY: tb_get_basis,&
207 tb_init_wf,&
210 USE xtb_parameters, ONLY: init_xtb_basis,&
216#include "./base/base_uses.f90"
217
218 IMPLICIT NONE
219
220 PRIVATE
221
222 ! *** Global parameters ***
223 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_environment'
224
225 ! *** Public subroutines ***
226 PUBLIC :: qs_init
227
228CONTAINS
229
230! **************************************************************************************************
231!> \brief Read the input and the database files for the setup of the
232!> QUICKSTEP environment.
233!> \param qs_env ...
234!> \param para_env ...
235!> \param root_section ...
236!> \param globenv ...
237!> \param cp_subsys ...
238!> \param kpoint_env ...
239!> \param cell ...
240!> \param cell_ref ...
241!> \param qmmm ...
242!> \param qmmm_env_qm ...
243!> \param force_env_section ...
244!> \param subsys_section ...
245!> \param use_motion_section ...
246!> \param silent ...
247!> \author Creation (22.05.2000,MK)
248! **************************************************************************************************
249 SUBROUTINE qs_init(qs_env, para_env, root_section, globenv, cp_subsys, kpoint_env, cell, cell_ref, &
250 qmmm, qmmm_env_qm, force_env_section, subsys_section, &
251 use_motion_section, silent)
252
253 TYPE(qs_environment_type), POINTER :: qs_env
254 TYPE(mp_para_env_type), POINTER :: para_env
255 TYPE(section_vals_type), OPTIONAL, POINTER :: root_section
256 TYPE(global_environment_type), OPTIONAL, POINTER :: globenv
257 TYPE(cp_subsys_type), OPTIONAL, POINTER :: cp_subsys
258 TYPE(kpoint_type), OPTIONAL, POINTER :: kpoint_env
259 TYPE(cell_type), OPTIONAL, POINTER :: cell, cell_ref
260 LOGICAL, INTENT(IN), OPTIONAL :: qmmm
261 TYPE(qmmm_env_qm_type), OPTIONAL, POINTER :: qmmm_env_qm
262 TYPE(section_vals_type), POINTER :: force_env_section, subsys_section
263 LOGICAL, INTENT(IN) :: use_motion_section
264 LOGICAL, INTENT(IN), OPTIONAL :: silent
265
266 CHARACTER(LEN=default_string_length) :: basis_type
267 INTEGER :: ikind, method_id, nelectron_total, &
268 nkind, nkp_grid(3)
269 LOGICAL :: do_admm_rpa, do_ec_hfx, do_et, do_exx, do_hfx, do_kpoints, is_identical, is_semi, &
270 mp2_present, my_qmmm, qmmm_decoupl, same_except_frac, use_ref_cell
271 REAL(kind=dp), DIMENSION(:, :), POINTER :: rtmat
272 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
273 TYPE(cell_type), POINTER :: my_cell, my_cell_ref
274 TYPE(cp_blacs_env_type), POINTER :: blacs_env
275 TYPE(dft_control_type), POINTER :: dft_control
276 TYPE(distribution_1d_type), POINTER :: local_particles
277 TYPE(energy_correction_type), POINTER :: ec_env
278 TYPE(excited_energy_type), POINTER :: exstate_env
279 TYPE(harris_type), POINTER :: harris_env
280 TYPE(kpoint_type), POINTER :: kpoints
281 TYPE(lri_environment_type), POINTER :: lri_env
282 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
283 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
284 TYPE(qs_ks_env_type), POINTER :: ks_env
285 TYPE(qs_subsys_type), POINTER :: subsys
286 TYPE(qs_wf_history_type), POINTER :: wf_history
287 TYPE(rel_control_type), POINTER :: rel_control
288 TYPE(scf_control_type), POINTER :: scf_control
289 TYPE(section_vals_type), POINTER :: dft_section, ec_hfx_section, ec_section, &
290 et_coupling_section, hfx_section, kpoint_section, mp2_section, rpa_hfx_section, &
291 transport_section
292
293 NULLIFY (my_cell, my_cell_ref, atomic_kind_set, particle_set, &
294 qs_kind_set, kpoint_section, dft_section, ec_section, &
295 subsys, ks_env, dft_control, blacs_env)
296
297 CALL set_qs_env(qs_env, input=force_env_section)
298 IF (.NOT. ASSOCIATED(subsys_section)) THEN
299 subsys_section => section_vals_get_subs_vals(force_env_section, "SUBSYS")
300 END IF
301
302 ! QMMM
303 my_qmmm = .false.
304 IF (PRESENT(qmmm)) my_qmmm = qmmm
305 qmmm_decoupl = .false.
306 IF (PRESENT(qmmm_env_qm)) THEN
307 IF (qmmm_env_qm%qmmm_coupl_type == do_qmmm_gauss .OR. &
308 qmmm_env_qm%qmmm_coupl_type == do_qmmm_swave) THEN
309 ! For GAUSS/SWAVE methods there could be a DDAPC decoupling requested
310 qmmm_decoupl = my_qmmm .AND. qmmm_env_qm%periodic .AND. qmmm_env_qm%multipole
311 END IF
312 qs_env%qmmm_env_qm => qmmm_env_qm
313 END IF
314 CALL set_qs_env(qs_env=qs_env, qmmm=my_qmmm)
315
316 ! Possibly initialize arrays for SE
317 CALL section_vals_val_get(force_env_section, "DFT%QS%METHOD", i_val=method_id)
318 SELECT CASE (method_id)
321 CALL init_se_intd_array()
322 is_semi = .true.
324 is_semi = .true.
325 CASE DEFAULT
326 is_semi = .false.
327 END SELECT
328
329 ALLOCATE (subsys)
330 CALL qs_subsys_create(subsys, para_env, &
331 force_env_section=force_env_section, &
332 subsys_section=subsys_section, &
333 use_motion_section=use_motion_section, &
334 root_section=root_section, &
335 cp_subsys=cp_subsys, cell=cell, cell_ref=cell_ref, &
336 elkind=is_semi, silent=silent)
337
338 ALLOCATE (ks_env)
339 CALL qs_ks_env_create(ks_env)
340 CALL set_ks_env(ks_env, subsys=subsys)
341 CALL set_qs_env(qs_env, ks_env=ks_env)
342
343 CALL qs_subsys_get(subsys, &
344 cell=my_cell, &
345 cell_ref=my_cell_ref, &
346 use_ref_cell=use_ref_cell, &
347 atomic_kind_set=atomic_kind_set, &
348 qs_kind_set=qs_kind_set, &
349 particle_set=particle_set)
350
351 CALL set_ks_env(ks_env, para_env=para_env)
352 IF (PRESENT(globenv)) THEN
353 CALL cp_blacs_env_create(blacs_env, para_env, globenv%blacs_grid_layout, &
354 globenv%blacs_repeatable)
355 ELSE
356 CALL cp_blacs_env_create(blacs_env, para_env)
357 END IF
358 CALL set_ks_env(ks_env, blacs_env=blacs_env)
359 CALL cp_blacs_env_release(blacs_env)
360
361 ! *** Setup the grids for the G-space Interpolation if any
362 CALL cp_ddapc_ewald_create(qs_env%cp_ddapc_ewald, qmmm_decoupl, my_cell, &
363 force_env_section, subsys_section, para_env)
364
365 ! kpoints
366 IF (PRESENT(kpoint_env)) THEN
367 kpoints => kpoint_env
368 CALL set_qs_env(qs_env=qs_env, kpoints=kpoints)
369 CALL kpoint_initialize(kpoints, particle_set, my_cell)
370 ELSE
371 NULLIFY (kpoints)
372 CALL kpoint_create(kpoints)
373 CALL set_qs_env(qs_env=qs_env, kpoints=kpoints)
374 kpoint_section => section_vals_get_subs_vals(qs_env%input, "DFT%KPOINTS")
375 CALL read_kpoint_section(kpoints, kpoint_section, my_cell%hmat)
376 CALL kpoint_initialize(kpoints, particle_set, my_cell)
377 dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
378 CALL write_kpoint_info(kpoints, dft_section)
379 END IF
380
381 CALL qs_init_subsys(qs_env, para_env, subsys, my_cell, my_cell_ref, use_ref_cell, &
382 subsys_section, silent=silent)
383
384 CALL get_qs_env(qs_env, dft_control=dft_control)
385 IF (method_id == do_method_lrigpw .OR. dft_control%qs_control%lri_optbas) THEN
386 CALL get_qs_env(qs_env=qs_env, lri_env=lri_env)
387 CALL lri_env_basis("LRI", qs_env, lri_env, qs_kind_set)
388 ELSE IF (method_id == do_method_rigpw) THEN
389 CALL cp_warn(__location__, "Experimental code: "// &
390 "RIGPW should only be used for testing.")
391 CALL get_qs_env(qs_env=qs_env, lri_env=lri_env)
392 CALL lri_env_basis("RI", qs_env, lri_env, qs_kind_set)
393 END IF
394
395 IF (my_qmmm .AND. PRESENT(qmmm_env_qm) .AND. .NOT. dft_control%qs_control%commensurate_mgrids) THEN
396 IF (qmmm_env_qm%qmmm_coupl_type == do_qmmm_gauss .OR. qmmm_env_qm%qmmm_coupl_type == do_qmmm_swave) THEN
397 CALL cp_abort(__location__, "QM/MM with coupling GAUSS or S-WAVE requires "// &
398 "keyword FORCE_EVAL/DFT/MGRID/COMMENSURATE to be enabled.")
399 END IF
400 END IF
401
402 ! more kpoint stuff
403 CALL get_qs_env(qs_env=qs_env, do_kpoints=do_kpoints, blacs_env=blacs_env)
404 IF (do_kpoints) THEN
405 CALL kpoint_env_initialize(kpoints, para_env, blacs_env, with_aux_fit=dft_control%do_admm)
406 CALL kpoint_initialize_mos(kpoints, qs_env%mos)
407 CALL get_qs_env(qs_env=qs_env, wf_history=wf_history)
408 CALL wfi_create_for_kp(wf_history)
409 END IF
410 ! basis set symmetry rotations
411 IF (do_kpoints) THEN
412 CALL qs_basis_rotation(qs_env, kpoints)
413 END IF
414
415 do_hfx = .false.
416 hfx_section => section_vals_get_subs_vals(qs_env%input, "DFT%XC%HF")
417 CALL section_vals_get(hfx_section, explicit=do_hfx)
418 CALL get_qs_env(qs_env, dft_control=dft_control, scf_control=scf_control, nelectron_total=nelectron_total)
419 IF (do_hfx) THEN
420 ! Retrieve particle_set and atomic_kind_set (needed for both kinds of initialization)
421 nkp_grid = 1
422 IF (do_kpoints) CALL get_kpoint_info(kpoints, nkp_grid=nkp_grid)
423 IF (dft_control%do_admm) THEN
424 basis_type = 'AUX_FIT'
425 ELSE
426 basis_type = 'ORB'
427 END IF
428 CALL hfx_create(qs_env%x_data, para_env, hfx_section, atomic_kind_set, &
429 qs_kind_set, particle_set, dft_control, my_cell, orb_basis=basis_type, &
430 nelectron_total=nelectron_total, nkp_grid=nkp_grid)
431 END IF
432
433 mp2_section => section_vals_get_subs_vals(qs_env%input, "DFT%XC%WF_CORRELATION")
434 CALL section_vals_get(mp2_section, explicit=mp2_present)
435 IF (mp2_present) THEN
436 cpassert(ASSOCIATED(qs_env%mp2_env))
437 CALL read_mp2_section(qs_env%input, qs_env%mp2_env)
438 ! create the EXX section if necessary
439 do_exx = .false.
440 rpa_hfx_section => section_vals_get_subs_vals(qs_env%input, "DFT%XC%WF_CORRELATION%RI_RPA%HF")
441 CALL section_vals_get(rpa_hfx_section, explicit=do_exx)
442 IF (do_exx) THEN
443
444 ! do_exx in call of hfx_create decides whether to go without ADMM (do_exx=.TRUE.) or with
445 ! ADMM (do_exx=.FALSE.)
446 CALL section_vals_val_get(mp2_section, "RI_RPA%ADMM", l_val=do_admm_rpa)
447
448 ! Reuse the HFX integrals from the qs_env if applicable
449 qs_env%mp2_env%ri_rpa%reuse_hfx = .true.
450 IF (.NOT. do_hfx) qs_env%mp2_env%ri_rpa%reuse_hfx = .false.
451 CALL compare_hfx_sections(hfx_section, rpa_hfx_section, is_identical, same_except_frac)
452 IF (.NOT. (is_identical .OR. same_except_frac)) qs_env%mp2_env%ri_rpa%reuse_hfx = .false.
453 IF (dft_control%do_admm .AND. .NOT. do_admm_rpa) qs_env%mp2_env%ri_rpa%reuse_hfx = .false.
454
455 IF (.NOT. qs_env%mp2_env%ri_rpa%reuse_hfx) THEN
456 IF (do_admm_rpa) THEN
457 basis_type = 'AUX_FIT'
458 ELSE
459 basis_type = 'ORB'
460 END IF
461 CALL hfx_create(qs_env%mp2_env%ri_rpa%x_data, para_env, rpa_hfx_section, atomic_kind_set, &
462 qs_kind_set, particle_set, dft_control, my_cell, orb_basis=basis_type, &
463 nelectron_total=nelectron_total)
464 ELSE
465 qs_env%mp2_env%ri_rpa%x_data => qs_env%x_data
466 END IF
467 END IF
468 END IF
469
470 IF (dft_control%qs_control%do_kg) THEN
471 CALL cite_reference(iannuzzi2006)
472 CALL kg_env_create(qs_env, qs_env%kg_env, qs_kind_set, qs_env%input)
473 END IF
474
475 dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
476 CALL section_vals_val_get(dft_section, "EXCITED_STATES%_SECTION_PARAMETERS_", &
477 l_val=qs_env%excited_state)
478 NULLIFY (exstate_env)
479 CALL exstate_create(exstate_env, qs_env%excited_state, dft_section)
480 CALL set_qs_env(qs_env, exstate_env=exstate_env)
481
482 et_coupling_section => section_vals_get_subs_vals(qs_env%input, &
483 "PROPERTIES%ET_COUPLING")
484 CALL section_vals_get(et_coupling_section, explicit=do_et)
485 IF (do_et) CALL et_coupling_create(qs_env%et_coupling)
486
487 transport_section => section_vals_get_subs_vals(qs_env%input, "DFT%TRANSPORT")
488 CALL section_vals_get(transport_section, explicit=qs_env%do_transport)
489 IF (qs_env%do_transport) THEN
490 CALL transport_env_create(qs_env)
491 END IF
492
493 CALL get_qs_env(qs_env, harris_env=harris_env)
494 IF (qs_env%harris_method) THEN
495 ! initialize the Harris input density and potential integrals
496 CALL get_qs_env(qs_env, local_particles=local_particles)
497 CALL harris_rhoin_init(harris_env%rhoin, "RHOIN", qs_kind_set, atomic_kind_set, &
498 local_particles, dft_control%nspins)
499 ! Print information of the HARRIS section
500 CALL harris_write_input(harris_env)
501 END IF
502
503 NULLIFY (ec_env)
504 dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
505 CALL section_vals_val_get(dft_section, "ENERGY_CORRECTION%_SECTION_PARAMETERS_", &
506 l_val=qs_env%energy_correction)
507 ec_section => section_vals_get_subs_vals(qs_env%input, "DFT%ENERGY_CORRECTION")
508 CALL ec_env_create(qs_env, ec_env, dft_section, ec_section)
509 CALL set_qs_env(qs_env, ec_env=ec_env)
510
511 IF (qs_env%energy_correction) THEN
512 ! Energy correction with Hartree-Fock exchange
513 ec_hfx_section => section_vals_get_subs_vals(ec_section, "XC%HF")
514 CALL section_vals_get(ec_hfx_section, explicit=do_ec_hfx)
515
516 IF (ec_env%do_ec_hfx) THEN
517
518 ! Hybrid functionals require same basis
519 IF (ec_env%basis_inconsistent) THEN
520 CALL cp_abort(__location__, &
521 "Energy correction methods with hybrid functionals: "// &
522 "correction and ground state need to use the same basis. "// &
523 "Checked by comparing basis set names only.")
524 END IF
525
526 ! Similar to RPA_HFX we can check if HFX integrals from the qs_env can be reused
527 IF (ec_env%do_ec_admm .AND. .NOT. dft_control%do_admm) THEN
528 CALL cp_abort(__location__, "Need an ADMM input section for ADMM EC to work")
529 END IF
530
531 ec_env%reuse_hfx = .true.
532 IF (.NOT. do_hfx) ec_env%reuse_hfx = .false.
533 CALL compare_hfx_sections(hfx_section, ec_hfx_section, is_identical, same_except_frac)
534 IF (.NOT. (is_identical .OR. same_except_frac)) ec_env%reuse_hfx = .false.
535 IF (dft_control%do_admm .AND. .NOT. ec_env%do_ec_admm) ec_env%reuse_hfx = .false.
536
537 IF (.NOT. ec_env%reuse_hfx) THEN
538 IF (ec_env%do_ec_admm) THEN
539 basis_type = 'AUX_FIT'
540 ELSE
541 basis_type = 'ORB'
542 END IF
543 CALL hfx_create(ec_env%x_data, para_env, ec_hfx_section, atomic_kind_set, &
544 qs_kind_set, particle_set, dft_control, my_cell, orb_basis=basis_type, &
545 nelectron_total=nelectron_total)
546 ELSE
547 ec_env%x_data => qs_env%x_data
548 END IF
549 END IF
550
551 ! Print information of the EC section
552 CALL ec_write_input(ec_env)
553
554 END IF
555
556 IF (dft_control%qs_control%do_almo_scf) THEN
557 CALL almo_scf_env_create(qs_env)
558 END IF
559
560 ! see if we have atomic relativistic corrections
561 CALL get_qs_env(qs_env, rel_control=rel_control)
562 IF (rel_control%rel_method /= rel_none) THEN
563 IF (rel_control%rel_transformation == rel_trans_atom) THEN
564 nkind = SIZE(atomic_kind_set)
565 DO ikind = 1, nkind
566 NULLIFY (rtmat)
567 CALL calculate_atomic_relkin(atomic_kind_set(ikind), qs_kind_set(ikind), rel_control, rtmat)
568 IF (ASSOCIATED(rtmat)) CALL set_qs_kind(qs_kind_set(ikind), reltmat=rtmat)
569 END DO
570 END IF
571 END IF
572
573 END SUBROUTINE qs_init
574
575! **************************************************************************************************
576!> \brief Initialize the qs environment (subsys)
577!> \param qs_env ...
578!> \param para_env ...
579!> \param subsys ...
580!> \param cell ...
581!> \param cell_ref ...
582!> \param use_ref_cell ...
583!> \param subsys_section ...
584!> \param silent ...
585!> \author Creation (22.05.2000,MK)
586! **************************************************************************************************
587 SUBROUTINE qs_init_subsys(qs_env, para_env, subsys, cell, cell_ref, use_ref_cell, subsys_section, &
588 silent)
589
590 TYPE(qs_environment_type), POINTER :: qs_env
591 TYPE(mp_para_env_type), POINTER :: para_env
592 TYPE(qs_subsys_type), POINTER :: subsys
593 TYPE(cell_type), POINTER :: cell, cell_ref
594 LOGICAL, INTENT(in) :: use_ref_cell
595 TYPE(section_vals_type), POINTER :: subsys_section
596 LOGICAL, INTENT(in), OPTIONAL :: silent
597
598 CHARACTER(len=*), PARAMETER :: routinen = 'qs_init_subsys'
599
600 CHARACTER(len=2) :: element_symbol
601 INTEGER :: gfn_type, handle, ikind, ispin, iw, lmax_sphere, maxl, maxlgto, maxlgto_lri, &
602 maxlgto_nuc, maxlppl, maxlppnl, method_id, multiplicity, my_ival, n_ao, n_mo_add, natom, &
603 nelectron, ngauss, nkind, output_unit, sort_basis, tnadd_method
604 INTEGER, DIMENSION(2) :: n_mo, nelectron_spin
605 INTEGER, DIMENSION(5) :: occ
606 LOGICAL :: all_potential_present, be_silent, cneo_potential_present, do_kpoints, do_ri_hfx, &
607 do_ri_mp2, do_ri_rpa, do_ri_sos_mp2, do_rpa_ri_exx, do_wfc_im_time, e1terms, &
608 has_unit_metric, lribas, mp2_present, orb_gradient
609 REAL(kind=dp) :: alpha, ccore, ewald_rcut, fxx, maxocc, &
610 rcut, total_zeff_corr, verlet_skin, &
611 zeff_correction
612 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
613 TYPE(cp_logger_type), POINTER :: logger
614 TYPE(dft_control_type), POINTER :: dft_control
615 TYPE(dftb_control_type), POINTER :: dftb_control
616 TYPE(distribution_1d_type), POINTER :: local_molecules, local_particles
617 TYPE(ewald_environment_type), POINTER :: ewald_env
618 TYPE(ewald_pw_type), POINTER :: ewald_pw
619 TYPE(fist_nonbond_env_type), POINTER :: se_nonbond_env
620 TYPE(gapw_control_type), POINTER :: gapw_control
621 TYPE(gto_basis_set_type), POINTER :: aux_fit_basis, lri_aux_basis, &
622 rhoin_basis, ri_aux_basis_set, &
623 ri_hfx_basis, ri_xas_basis, &
624 tmp_basis_set
625 TYPE(harris_type), POINTER :: harris_env
626 TYPE(local_rho_type), POINTER :: local_rho_set
627 TYPE(lri_environment_type), POINTER :: lri_env
628 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos, mos_last_converged
629 TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
630 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
631 TYPE(mp2_type), POINTER :: mp2_env
632 TYPE(nddo_mpole_type), POINTER :: se_nddo_mpole
633 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
634 TYPE(pw_env_type), POINTER :: pw_env
635 TYPE(qs_control_type), POINTER :: qs_control
636 TYPE(qs_dftb_pairpot_type), DIMENSION(:, :), &
637 POINTER :: dftb_potential
638 TYPE(qs_dispersion_type), POINTER :: dispersion_env
639 TYPE(qs_energy_type), POINTER :: energy
640 TYPE(qs_force_type), DIMENSION(:), POINTER :: force
641 TYPE(qs_gcp_type), POINTER :: gcp_env
642 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
643 TYPE(qs_kind_type), POINTER :: qs_kind
644 TYPE(qs_ks_env_type), POINTER :: ks_env
645 TYPE(qs_wf_history_type), POINTER :: wf_history
646 TYPE(rho0_mpole_type), POINTER :: rho0_mpole
647 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set
648 TYPE(scf_control_type), POINTER :: scf_control
649 TYPE(se_taper_type), POINTER :: se_taper
650 TYPE(section_vals_type), POINTER :: dft_section, et_coupling_section, et_ddapc_section, &
651 ewald_section, harris_section, lri_section, mp2_section, nl_section, poisson_section, &
652 pp_section, print_section, qs_section, rixs_section, se_section, tddfpt_section, &
653 xc_section
654 TYPE(semi_empirical_control_type), POINTER :: se_control
655 TYPE(semi_empirical_si_type), POINTER :: se_store_int_env
656 TYPE(xtb_control_type), POINTER :: xtb_control
657
658 CALL timeset(routinen, handle)
659 NULLIFY (logger)
660 logger => cp_get_default_logger()
661 output_unit = cp_logger_get_default_io_unit(logger)
662
663 be_silent = .false.
664 IF (PRESENT(silent)) be_silent = silent
665
666 CALL cite_reference(cp2kqs2020)
667
668 ! Initialise the Quickstep environment
669 NULLIFY (mos, se_taper)
670 NULLIFY (dft_control)
671 NULLIFY (energy)
672 NULLIFY (force)
673 NULLIFY (local_molecules)
674 NULLIFY (local_particles)
675 NULLIFY (scf_control)
676 NULLIFY (dft_section)
677 NULLIFY (et_coupling_section)
678 NULLIFY (ks_env)
679 NULLIFY (mos_last_converged)
680 dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
681 qs_section => section_vals_get_subs_vals(dft_section, "QS")
682 et_coupling_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%ET_COUPLING")
683 ! reimplemented TDDFPT
684 tddfpt_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%TDDFPT")
685 rixs_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%RIXS")
686
687 CALL qs_subsys_get(subsys, particle_set=particle_set, &
688 qs_kind_set=qs_kind_set, &
689 atomic_kind_set=atomic_kind_set, &
690 molecule_set=molecule_set, &
691 molecule_kind_set=molecule_kind_set)
692
693 ! *** Read the input section with the DFT control parameters ***
694 CALL read_dft_control(dft_control, dft_section)
695
696 ! set periodicity flag
697 dft_control%qs_control%periodicity = sum(cell%perd)
698
699 ! Read the input section with the Quickstep control parameters
700 CALL read_qs_section(dft_control%qs_control, qs_section)
701
702 ! *** Print the Quickstep program banner (copyright and version number) ***
703 IF (.NOT. be_silent) THEN
704 iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%PROGRAM_BANNER", extension=".Log")
705 CALL section_vals_val_get(qs_section, "METHOD", i_val=method_id)
706 SELECT CASE (method_id)
707 CASE DEFAULT
708 CALL qs_header(iw)
711 CALL se_header(iw)
712 CASE (do_method_dftb)
713 CALL dftb_header(iw)
714 CASE (do_method_xtb)
715 IF (dft_control%qs_control%xtb_control%do_tblite) THEN
716 CALL tblite_header(iw, dft_control%qs_control%xtb_control%tblite_method)
717 ELSE
718 gfn_type = dft_control%qs_control%xtb_control%gfn_type
719 CALL xtb_header(iw, gfn_type)
720 END IF
721 END SELECT
722 CALL cp_print_key_finished_output(iw, logger, dft_section, &
723 "PRINT%PROGRAM_BANNER")
724 END IF
725
726 IF (dft_control%do_sccs .AND. dft_control%qs_control%gapw) THEN
727 cpabort("SCCS is not yet implemented with GAPW")
728 END IF
729 CALL get_qs_env(qs_env=qs_env, do_kpoints=do_kpoints)
730 IF (do_kpoints) THEN
731 ! reset some of the settings for wfn extrapolation for kpoints
732 SELECT CASE (dft_control%qs_control%wf_interpolation_method_nr)
734 dft_control%qs_control%wf_interpolation_method_nr = wfi_use_guess_method_nr
735 END SELECT
736 END IF
737
738 ! ******* check if any kind of electron transfer calculation has to be performed
739 CALL section_vals_val_get(et_coupling_section, "TYPE_OF_CONSTRAINT", i_val=my_ival)
740 dft_control%qs_control%et_coupling_calc = .false.
741 IF (my_ival == do_et_ddapc) THEN
742 et_ddapc_section => section_vals_get_subs_vals(et_coupling_section, "DDAPC_RESTRAINT_A")
743 dft_control%qs_control%et_coupling_calc = .true.
744 dft_control%qs_control%ddapc_restraint = .true.
745 CALL read_ddapc_section(dft_control%qs_control, ddapc_restraint_section=et_ddapc_section)
746 END IF
747
748 CALL read_mgrid_section(dft_control%qs_control, dft_section)
749
750 ! reimplemented TDDFPT
751 CALL read_tddfpt2_control(dft_control%tddfpt2_control, tddfpt_section, dft_control%qs_control)
752
753 ! rixs
754 CALL section_vals_get(rixs_section, explicit=qs_env%do_rixs)
755 IF (qs_env%do_rixs) THEN
756 CALL read_rixs_control(dft_control%rixs_control, rixs_section, dft_control%qs_control)
757 END IF
758
759 ! Create relativistic control section
760 block
761 TYPE(rel_control_type), POINTER :: rel_control
762 ALLOCATE (rel_control)
763 CALL rel_c_create(rel_control)
764 CALL rel_c_read_parameters(rel_control, dft_section)
765 CALL set_qs_env(qs_env, rel_control=rel_control)
766 END block
767
768 ! *** Read DFTB parameter files ***
769 IF (dft_control%qs_control%method_id == do_method_dftb) THEN
770 NULLIFY (ewald_env, ewald_pw, dftb_potential)
771 dftb_control => dft_control%qs_control%dftb_control
772 CALL qs_dftb_param_init(atomic_kind_set, qs_kind_set, dftb_control, dftb_potential, &
773 subsys_section=subsys_section, para_env=para_env)
774 CALL set_qs_env(qs_env, dftb_potential=dftb_potential)
775 ! check for Ewald
776 IF (dftb_control%do_ewald) THEN
777 ALLOCATE (ewald_env)
778 CALL ewald_env_create(ewald_env, para_env)
779 poisson_section => section_vals_get_subs_vals(dft_section, "POISSON")
780 CALL ewald_env_set(ewald_env, poisson_section=poisson_section)
781 ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD")
782 print_section => section_vals_get_subs_vals(qs_env%input, "PRINT%GRID_INFORMATION")
783 CALL get_qs_kind_set(qs_kind_set, basis_rcut=ewald_rcut)
784 CALL read_ewald_section_tb(ewald_env, ewald_section, cell_ref%hmat)
785 ALLOCATE (ewald_pw)
786 CALL ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section=print_section)
787 CALL set_qs_env(qs_env, ewald_env=ewald_env, ewald_pw=ewald_pw)
788 END IF
789 ELSEIF (dft_control%qs_control%method_id == do_method_xtb) THEN
790 ! *** Read xTB parameter file ***
791 xtb_control => dft_control%qs_control%xtb_control
792 CALL get_qs_env(qs_env, nkind=nkind)
793 IF (xtb_control%do_tblite) THEN
794 ! put geometry to tblite
795 CALL tb_init_geometry(qs_env, qs_env%tb_tblite)
796 ! select tblite method
797 CALL tb_set_calculator(qs_env%tb_tblite, xtb_control%tblite_method)
798 !set up wave function
799 CALL tb_init_wf(qs_env%tb_tblite)
800 !get basis set
801 DO ikind = 1, nkind
802 qs_kind => qs_kind_set(ikind)
803 ! Setup proper xTB parameters
804 cpassert(.NOT. ASSOCIATED(qs_kind%xtb_parameter))
805 CALL allocate_xtb_atom_param(qs_kind%xtb_parameter)
806 ! Set default parameters
807 CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
808
809 NULLIFY (tmp_basis_set)
810 CALL tb_get_basis(qs_env%tb_tblite, tmp_basis_set, element_symbol, qs_kind%xtb_parameter, occ)
811 CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, "ORB")
812 CALL set_xtb_atom_param(qs_kind%xtb_parameter, occupation=occ)
813
814 !setting the potential for the computation
815 zeff_correction = 0.0_dp
816 CALL init_potential(qs_kind%all_potential, itype="BARE", &
817 zeff=real(sum(occ), dp), zeff_correction=zeff_correction)
818 END DO
819 ELSE
820 NULLIFY (ewald_env, ewald_pw)
821 DO ikind = 1, nkind
822 qs_kind => qs_kind_set(ikind)
823 ! Setup proper xTB parameters
824 cpassert(.NOT. ASSOCIATED(qs_kind%xtb_parameter))
825 CALL allocate_xtb_atom_param(qs_kind%xtb_parameter)
826 ! Set default parameters
827 gfn_type = dft_control%qs_control%xtb_control%gfn_type
828 CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
829 CALL xtb_parameters_init(qs_kind%xtb_parameter, gfn_type, element_symbol, &
830 xtb_control%parameter_file_path, xtb_control%parameter_file_name, &
831 para_env)
832 ! set dependent parameters
833 CALL xtb_parameters_set(qs_kind%xtb_parameter)
834 ! Generate basis set
835 NULLIFY (tmp_basis_set)
836 IF (qs_kind%xtb_parameter%z == 1) THEN
837 ! special case hydrogen
838 ngauss = xtb_control%h_sto_ng
839 ELSE
840 ngauss = xtb_control%sto_ng
841 END IF
842 IF (qs_kind%xtb_parameter%defined) THEN
843 CALL init_xtb_basis(qs_kind%xtb_parameter, tmp_basis_set, ngauss)
844 CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, "ORB")
845 ELSE
846 CALL set_qs_kind(qs_kind, ghost=.true.)
847 IF (ASSOCIATED(qs_kind%all_potential)) THEN
848 DEALLOCATE (qs_kind%all_potential%elec_conf)
849 DEALLOCATE (qs_kind%all_potential)
850 END IF
851 END IF
852 ! potential
853 IF (qs_kind%xtb_parameter%defined) THEN
854 zeff_correction = 0.0_dp
855 CALL init_potential(qs_kind%all_potential, itype="BARE", &
856 zeff=qs_kind%xtb_parameter%zeff, zeff_correction=zeff_correction)
857 CALL get_potential(qs_kind%all_potential, alpha_core_charge=alpha)
858 ccore = qs_kind%xtb_parameter%zeff*sqrt((alpha/pi)**3)
859 CALL set_potential(qs_kind%all_potential, ccore_charge=ccore)
860 qs_kind%xtb_parameter%zeff = qs_kind%xtb_parameter%zeff - zeff_correction
861 END IF
862 END DO
863 !
864 ! set repulsive potential range
865 !
866 ALLOCATE (xtb_control%rcpair(nkind, nkind))
867 CALL xtb_pp_radius(qs_kind_set, xtb_control%rcpair, xtb_control%eps_pair, xtb_control%kf)
868 ! check for Ewald
869 IF (xtb_control%do_ewald) THEN
870 ALLOCATE (ewald_env)
871 CALL ewald_env_create(ewald_env, para_env)
872 poisson_section => section_vals_get_subs_vals(dft_section, "POISSON")
873 CALL ewald_env_set(ewald_env, poisson_section=poisson_section)
874 ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD")
875 print_section => section_vals_get_subs_vals(qs_env%input, "PRINT%GRID_INFORMATION")
876 IF (gfn_type == 0) THEN
877 CALL read_ewald_section_tb(ewald_env, ewald_section, cell_ref%hmat, &
878 silent=silent, pset="EEQ")
879 ELSE
880 CALL read_ewald_section_tb(ewald_env, ewald_section, cell_ref%hmat, &
881 silent=silent)
882 END IF
883 ALLOCATE (ewald_pw)
884 CALL ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section=print_section)
885 CALL set_qs_env(qs_env, ewald_env=ewald_env, ewald_pw=ewald_pw)
886 END IF
887 END IF
888 END IF
889 ! lri or ri env initialization
890 lri_section => section_vals_get_subs_vals(qs_section, "LRIGPW")
891 IF (dft_control%qs_control%method_id == do_method_lrigpw .OR. &
892 dft_control%qs_control%lri_optbas .OR. &
893 dft_control%qs_control%method_id == do_method_rigpw) THEN
894 CALL lri_env_init(lri_env, lri_section)
895 CALL set_qs_env(qs_env, lri_env=lri_env)
896 END IF
897
898 ! *** Check basis and fill in missing parts ***
899 CALL check_qs_kind_set(qs_kind_set, dft_control, subsys_section=subsys_section)
900
901 ! *** Check that no all-electron potential is present if GPW or GAPW_XC
902 CALL get_qs_kind_set(qs_kind_set, all_potential_present=all_potential_present)
903 IF ((dft_control%qs_control%method_id == do_method_gpw) .OR. &
904 (dft_control%qs_control%method_id == do_method_gapw_xc) .OR. &
905 (dft_control%qs_control%method_id == do_method_ofgpw)) THEN
906 IF (all_potential_present) THEN
907 cpabort("All-electron calculations with GPW, GAPW_XC, and OFGPW are not implemented")
908 END IF
909 END IF
910
911 ! *** Check that no cneo potential is present if not GAPW
912 CALL get_qs_kind_set(qs_kind_set, cneo_potential_present=cneo_potential_present)
913 IF (cneo_potential_present .AND. &
914 dft_control%qs_control%method_id /= do_method_gapw) THEN
915 cpabort("CNEO calculations require GAPW method")
916 END IF
917
918 ! DFT+U
919 CALL get_qs_kind_set(qs_kind_set, dft_plus_u_atom_present=dft_control%dft_plus_u)
920
921 IF (dft_control%do_admm) THEN
922 ! Check if ADMM basis is available
923 CALL get_qs_env(qs_env, nkind=nkind)
924 DO ikind = 1, nkind
925 NULLIFY (aux_fit_basis)
926 qs_kind => qs_kind_set(ikind)
927 CALL get_qs_kind(qs_kind, basis_set=aux_fit_basis, basis_type="AUX_FIT")
928 IF (.NOT. (ASSOCIATED(aux_fit_basis))) THEN
929 ! AUX_FIT basis set is not available
930 cpabort("AUX_FIT basis set is not defined. ")
931 END IF
932 END DO
933 END IF
934
935 lribas = .false.
936 e1terms = .false.
937 IF (dft_control%qs_control%method_id == do_method_lrigpw) THEN
938 lribas = .true.
939 CALL get_qs_env(qs_env, lri_env=lri_env)
940 e1terms = lri_env%exact_1c_terms
941 END IF
942 IF (dft_control%qs_control%do_kg) THEN
943 CALL section_vals_val_get(dft_section, "KG_METHOD%TNADD_METHOD", i_val=tnadd_method)
944 IF (tnadd_method == kg_tnadd_embed_ri) lribas = .true.
945 END IF
946 IF (lribas) THEN
947 ! Check if LRI_AUX basis is available, auto-generate if needed
948 CALL get_qs_env(qs_env, nkind=nkind)
949 DO ikind = 1, nkind
950 NULLIFY (lri_aux_basis)
951 qs_kind => qs_kind_set(ikind)
952 CALL get_qs_kind(qs_kind, basis_set=lri_aux_basis, basis_type="LRI_AUX")
953 IF (.NOT. (ASSOCIATED(lri_aux_basis))) THEN
954 ! LRI_AUX basis set is not yet loaded
955 CALL cp_warn(__location__, "Automatic Generation of LRI_AUX basis. "// &
956 "This is experimental code.")
957 ! Generate a default basis
958 CALL create_lri_aux_basis_set(lri_aux_basis, qs_kind, dft_control%auto_basis_lri_aux, e1terms)
959 CALL add_basis_set_to_container(qs_kind%basis_sets, lri_aux_basis, "LRI_AUX")
960 END IF
961 END DO
962 END IF
963
964 CALL section_vals_val_get(qs_env%input, "DFT%XC%HF%RI%_SECTION_PARAMETERS_", l_val=do_ri_hfx)
965 CALL section_vals_val_get(qs_env%input, "DFT%XC%WF_CORRELATION%RI_RPA%HF%RI%_SECTION_PARAMETERS_", &
966 l_val=do_rpa_ri_exx)
967 IF (do_ri_hfx .OR. do_rpa_ri_exx) THEN
968 CALL get_qs_env(qs_env, nkind=nkind)
969 CALL section_vals_val_get(qs_env%input, "DFT%SORT_BASIS", i_val=sort_basis)
970 DO ikind = 1, nkind
971 NULLIFY (ri_hfx_basis)
972 qs_kind => qs_kind_set(ikind)
973 CALL get_qs_kind(qs_kind=qs_kind, basis_set=ri_hfx_basis, &
974 basis_type="RI_HFX")
975 IF (.NOT. (ASSOCIATED(ri_hfx_basis))) THEN
976 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto)
977 IF (dft_control%do_admm) THEN
978 CALL create_ri_aux_basis_set(ri_hfx_basis, qs_kind, dft_control%auto_basis_ri_hfx, &
979 basis_type="AUX_FIT", basis_sort=sort_basis)
980 ELSE
981 CALL create_ri_aux_basis_set(ri_hfx_basis, qs_kind, dft_control%auto_basis_ri_hfx, &
982 basis_sort=sort_basis)
983 END IF
984 CALL add_basis_set_to_container(qs_kind%basis_sets, ri_hfx_basis, "RI_HFX")
985 END IF
986 END DO
987 END IF
988
989 IF (dft_control%qs_control%method_id == do_method_rigpw) THEN
990 ! Check if RI_HXC basis is available, auto-generate if needed
991 CALL get_qs_env(qs_env, nkind=nkind)
992 DO ikind = 1, nkind
993 NULLIFY (ri_hfx_basis)
994 qs_kind => qs_kind_set(ikind)
995 CALL get_qs_kind(qs_kind, basis_set=ri_hfx_basis, basis_type="RI_HXC")
996 IF (.NOT. (ASSOCIATED(ri_hfx_basis))) THEN
997 ! Generate a default basis
998 CALL create_ri_aux_basis_set(ri_hfx_basis, qs_kind, dft_control%auto_basis_ri_hxc)
999 CALL add_basis_set_to_container(qs_kind%basis_sets, ri_hfx_basis, "RI_HXC")
1000 END IF
1001 END DO
1002 END IF
1003
1004 ! Harris method
1005 NULLIFY (harris_env)
1006 CALL section_vals_val_get(dft_section, "HARRIS_METHOD%_SECTION_PARAMETERS_", &
1007 l_val=qs_env%harris_method)
1008 harris_section => section_vals_get_subs_vals(dft_section, "HARRIS_METHOD")
1009 CALL harris_env_create(qs_env, harris_env, harris_section)
1010 CALL set_qs_env(qs_env, harris_env=harris_env)
1011 !
1012 IF (qs_env%harris_method) THEN
1013 CALL get_qs_env(qs_env, nkind=nkind)
1014 ! Check if RI_HXC basis is available, auto-generate if needed
1015 DO ikind = 1, nkind
1016 NULLIFY (tmp_basis_set)
1017 qs_kind => qs_kind_set(ikind)
1018 CALL get_qs_kind(qs_kind, basis_set=rhoin_basis, basis_type="RHOIN")
1019 IF (.NOT. (ASSOCIATED(rhoin_basis))) THEN
1020 ! Generate a default basis
1021 CALL create_ri_aux_basis_set(tmp_basis_set, qs_kind, dft_control%auto_basis_ri_hxc)
1022 IF (qs_env%harris_env%density_source == hden_atomic) THEN
1023 CALL create_primitive_basis_set(tmp_basis_set, rhoin_basis, lmax=0)
1024 CALL deallocate_gto_basis_set(tmp_basis_set)
1025 ELSE
1026 rhoin_basis => tmp_basis_set
1027 END IF
1028 CALL add_basis_set_to_container(qs_kind%basis_sets, rhoin_basis, "RHOIN")
1029 END IF
1030 END DO
1031 END IF
1032
1033 mp2_section => section_vals_get_subs_vals(qs_env%input, "DFT%XC%WF_CORRELATION")
1034 CALL section_vals_get(mp2_section, explicit=mp2_present)
1035 IF (mp2_present) THEN
1036
1037 ! basis should be sorted for imaginary time RPA/GW
1038 CALL section_vals_val_get(qs_env%input, "DFT%SORT_BASIS", i_val=sort_basis)
1039 CALL section_vals_val_get(qs_env%input, "DFT%XC%WF_CORRELATION%LOW_SCALING%_SECTION_PARAMETERS_", &
1040 l_val=do_wfc_im_time)
1041
1042 IF (do_wfc_im_time .AND. sort_basis /= basis_sort_zet) THEN
1043 CALL cp_warn(__location__, &
1044 "Low-scaling RPA requires SORT_BASIS EXP keyword (in DFT input section) for good performance")
1045 END IF
1046
1047 ! Check if RI_AUX basis (for MP2/RPA) is given, auto-generate if not
1048 CALL mp2_env_create(qs_env%mp2_env)
1049 CALL get_qs_env(qs_env, mp2_env=mp2_env, nkind=nkind)
1050 CALL section_vals_val_get(qs_env%input, "DFT%XC%WF_CORRELATION%RI_MP2%_SECTION_PARAMETERS_", l_val=do_ri_mp2)
1051 CALL section_vals_val_get(qs_env%input, "DFT%XC%WF_CORRELATION%RI_SOS_MP2%_SECTION_PARAMETERS_", l_val=do_ri_sos_mp2)
1052 CALL section_vals_val_get(qs_env%input, "DFT%XC%WF_CORRELATION%RI_RPA%_SECTION_PARAMETERS_", l_val=do_ri_rpa)
1053 IF (do_ri_mp2 .OR. do_ri_sos_mp2 .OR. do_ri_rpa) THEN
1054 DO ikind = 1, nkind
1055 NULLIFY (ri_aux_basis_set)
1056 qs_kind => qs_kind_set(ikind)
1057 CALL get_qs_kind(qs_kind=qs_kind, basis_set=ri_aux_basis_set, &
1058 basis_type="RI_AUX")
1059 IF (.NOT. (ASSOCIATED(ri_aux_basis_set))) THEN
1060 ! RI_AUX basis set is not yet loaded
1061 ! Generate a default basis
1062 CALL create_ri_aux_basis_set(ri_aux_basis_set, qs_kind, dft_control%auto_basis_ri_aux, basis_sort=sort_basis)
1063 CALL add_basis_set_to_container(qs_kind%basis_sets, ri_aux_basis_set, "RI_AUX")
1064 ! Add a flag, which allows to check if the basis was generated
1065 ! when applying ERI_METHOD OS to mp2, ri-rpa, gw etc
1066 qs_env%mp2_env%ri_aux_auto_generated = .true.
1067 END IF
1068 END DO
1069 END IF
1070
1071 END IF
1072
1073 IF (dft_control%do_xas_tdp_calculation .OR. qs_env%do_rixs) THEN
1074 ! Check if RI_XAS basis is given, auto-generate if not
1075 CALL get_qs_env(qs_env, nkind=nkind)
1076 DO ikind = 1, nkind
1077 NULLIFY (ri_xas_basis)
1078 qs_kind => qs_kind_set(ikind)
1079 CALL get_qs_kind(qs_kind, basis_set=ri_xas_basis, basis_type="RI_XAS")
1080 IF (.NOT. ASSOCIATED(ri_xas_basis)) THEN
1081 ! Generate a default basis
1082 CALL create_ri_aux_basis_set(ri_xas_basis, qs_kind, dft_control%auto_basis_ri_xas)
1083 CALL add_basis_set_to_container(qs_kind%basis_sets, ri_xas_basis, "RI_XAS")
1084 END IF
1085 END DO
1086 END IF
1087
1088 ! *** Initialize the spherical harmonics and ***
1089 ! *** the orbital transformation matrices ***
1090 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto, maxlppl=maxlppl, maxlppnl=maxlppnl)
1091
1092 ! CNEO nuclear basis contributes to GAPW rho0
1093 IF (cneo_potential_present) THEN
1094 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto_nuc, basis_type="NUC")
1095 maxlgto = max(maxlgto, maxlgto_nuc)
1096 END IF
1097 lmax_sphere = dft_control%qs_control%gapw_control%lmax_sphere
1098 IF (lmax_sphere .LT. 0) THEN
1099 lmax_sphere = 2*maxlgto
1100 dft_control%qs_control%gapw_control%lmax_sphere = lmax_sphere
1101 END IF
1102 IF (dft_control%qs_control%method_id == do_method_lrigpw .OR. dft_control%qs_control%lri_optbas) THEN
1103 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto_lri, basis_type="LRI_AUX")
1104 !take maxlgto from lri basis if larger (usually)
1105 maxlgto = max(maxlgto, maxlgto_lri)
1106 ELSE IF (dft_control%qs_control%method_id == do_method_rigpw) THEN
1107 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto_lri, basis_type="RI_HXC")
1108 maxlgto = max(maxlgto, maxlgto_lri)
1109 END IF
1110 IF (dft_control%do_xas_tdp_calculation .OR. qs_env%do_rixs) THEN
1111 !done as a precaution
1112 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto_lri, basis_type="RI_XAS")
1113 maxlgto = max(maxlgto, maxlgto_lri)
1114 END IF
1115 maxl = max(2*maxlgto, maxlppl, maxlppnl, lmax_sphere) + 1
1116
1117 CALL init_orbital_pointers(maxl)
1118
1119 CALL init_spherical_harmonics(maxl, 0)
1120
1121 ! *** Initialise the qs_kind_set ***
1122 CALL init_qs_kind_set(qs_kind_set)
1123
1124 ! *** Initialise GAPW soft basis and projectors
1125 IF (dft_control%qs_control%method_id == do_method_gapw .OR. &
1126 dft_control%qs_control%method_id == do_method_gapw_xc) THEN
1127 qs_control => dft_control%qs_control
1128 CALL init_gapw_basis_set(qs_kind_set, qs_control, qs_env%input)
1129 END IF
1130
1131 ! *** Initialise CNEO nuclear soft basis
1132 IF (cneo_potential_present) THEN
1133 CALL init_cneo_basis_set(qs_kind_set, qs_control)
1134 END IF
1135
1136 ! *** Initialize the pretabulation for the calculation of the ***
1137 ! *** incomplete Gamma function F_n(t) after McMurchie-Davidson ***
1138 CALL get_qs_kind_set(qs_kind_set, maxlgto=maxlgto)
1139 maxl = max(3*maxlgto + 1, 0)
1140 CALL init_md_ftable(maxl)
1141
1142 ! *** Initialize the atomic interaction radii ***
1143 CALL init_interaction_radii(dft_control%qs_control, qs_kind_set)
1144 !
1145 IF (dft_control%qs_control%method_id == do_method_xtb) THEN
1146 IF (.NOT. dft_control%qs_control%xtb_control%do_tblite) THEN
1147 ! cutoff radius
1148 DO ikind = 1, nkind
1149 qs_kind => qs_kind_set(ikind)
1150 IF (qs_kind%xtb_parameter%defined) THEN
1151 CALL get_qs_kind(qs_kind, basis_set=tmp_basis_set)
1152 rcut = xtb_control%coulomb_sr_cut
1153 fxx = 2.0_dp*xtb_control%coulomb_sr_eps*qs_kind%xtb_parameter%eta**2
1154 fxx = 0.80_dp*(1.0_dp/fxx)**0.3333_dp
1155 rcut = min(rcut, xtb_control%coulomb_sr_cut)
1156 qs_kind%xtb_parameter%rcut = min(rcut, fxx)
1157 ELSE
1158 qs_kind%xtb_parameter%rcut = 0.0_dp
1159 END IF
1160 END DO
1161 END IF
1162 END IF
1163
1164 IF (.NOT. be_silent) THEN
1165 CALL write_pgf_orb_radii("orb", atomic_kind_set, qs_kind_set, subsys_section)
1166 CALL write_pgf_orb_radii("aux", atomic_kind_set, qs_kind_set, subsys_section)
1167 CALL write_pgf_orb_radii("lri", atomic_kind_set, qs_kind_set, subsys_section)
1168 CALL write_pgf_orb_radii("nuc", atomic_kind_set, qs_kind_set, subsys_section)
1169 CALL write_core_charge_radii(atomic_kind_set, qs_kind_set, subsys_section)
1170 CALL write_ppl_radii(atomic_kind_set, qs_kind_set, subsys_section)
1171 CALL write_ppnl_radii(atomic_kind_set, qs_kind_set, subsys_section)
1172 CALL write_paw_radii(atomic_kind_set, qs_kind_set, subsys_section)
1173 END IF
1174
1175 ! *** Distribute molecules and atoms using the new data structures ***
1176 CALL distribute_molecules_1d(atomic_kind_set=atomic_kind_set, &
1177 particle_set=particle_set, &
1178 local_particles=local_particles, &
1179 molecule_kind_set=molecule_kind_set, &
1180 molecule_set=molecule_set, &
1181 local_molecules=local_molecules, &
1182 force_env_section=qs_env%input)
1183
1184 ! *** SCF parameters ***
1185 ALLOCATE (scf_control)
1186 ! set (non)-self consistency
1187 IF (dft_control%qs_control%dftb) THEN
1188 scf_control%non_selfconsistent = .NOT. dft_control%qs_control%dftb_control%self_consistent
1189 END IF
1190 IF (dft_control%qs_control%xtb) THEN
1191 IF (dft_control%qs_control%xtb_control%do_tblite) THEN
1192 scf_control%non_selfconsistent = .false.
1193 ELSE
1194 scf_control%non_selfconsistent = (dft_control%qs_control%xtb_control%gfn_type == 0)
1195 END IF
1196 END IF
1197 IF (qs_env%harris_method) THEN
1198 scf_control%non_selfconsistent = .true.
1199 END IF
1200 CALL scf_c_create(scf_control)
1201 CALL scf_c_read_parameters(scf_control, dft_section)
1202 ! *** Allocate the data structure for Quickstep energies ***
1203 CALL allocate_qs_energy(energy)
1204
1205 ! check for orthogonal basis
1206 has_unit_metric = .false.
1207 IF (dft_control%qs_control%semi_empirical) THEN
1208 IF (dft_control%qs_control%se_control%orthogonal_basis) has_unit_metric = .true.
1209 END IF
1210 IF (dft_control%qs_control%dftb) THEN
1211 IF (dft_control%qs_control%dftb_control%orthogonal_basis) has_unit_metric = .true.
1212 END IF
1213 CALL set_qs_env(qs_env, has_unit_metric=has_unit_metric)
1214
1215 ! *** Activate the interpolation ***
1216 CALL wfi_create(wf_history, &
1217 interpolation_method_nr= &
1218 dft_control%qs_control%wf_interpolation_method_nr, &
1219 extrapolation_order=dft_control%qs_control%wf_extrapolation_order, &
1220 has_unit_metric=has_unit_metric)
1221
1222 ! *** Set the current Quickstep environment ***
1223 CALL set_qs_env(qs_env=qs_env, &
1224 scf_control=scf_control, &
1225 wf_history=wf_history)
1226
1227 CALL qs_subsys_set(subsys, &
1228 cell_ref=cell_ref, &
1229 use_ref_cell=use_ref_cell, &
1230 energy=energy, &
1231 force=force)
1232
1233 CALL get_qs_env(qs_env, ks_env=ks_env)
1234 CALL set_ks_env(ks_env, dft_control=dft_control)
1235
1236 CALL qs_subsys_set(subsys, local_molecules=local_molecules, &
1237 local_particles=local_particles, cell=cell)
1238
1239 CALL distribution_1d_release(local_particles)
1240 CALL distribution_1d_release(local_molecules)
1241 CALL wfi_release(wf_history)
1242
1243 CALL get_qs_env(qs_env=qs_env, &
1244 atomic_kind_set=atomic_kind_set, &
1245 dft_control=dft_control, &
1246 scf_control=scf_control)
1247
1248 ! decide what conditions need mo_derivs
1249 ! right now, this only appears to be OT
1250 IF (dft_control%qs_control%do_ls_scf .OR. &
1251 dft_control%qs_control%do_almo_scf) THEN
1252 CALL set_qs_env(qs_env=qs_env, requires_mo_derivs=.false.)
1253 ELSE
1254 IF (scf_control%use_ot) THEN
1255 CALL set_qs_env(qs_env=qs_env, requires_mo_derivs=.true.)
1256 ELSE
1257 CALL set_qs_env(qs_env=qs_env, requires_mo_derivs=.false.)
1258 END IF
1259 END IF
1260
1261 ! XXXXXXX this is backwards XXXXXXXX
1262 dft_control%smear = scf_control%smear%do_smear
1263
1264 ! Periodic efield needs equal occupation and orbital gradients
1265 IF (.NOT. (dft_control%qs_control%dftb .OR. dft_control%qs_control%xtb)) THEN
1266 IF (dft_control%apply_period_efield) THEN
1267 CALL get_qs_env(qs_env=qs_env, requires_mo_derivs=orb_gradient)
1268 IF (.NOT. orb_gradient) THEN
1269 CALL cp_abort(__location__, "Periodic Efield needs orbital gradient and direct optimization."// &
1270 " Use the OT optimization method.")
1271 END IF
1272 IF (dft_control%smear) THEN
1273 CALL cp_abort(__location__, "Periodic Efield needs equal occupation numbers."// &
1274 " Smearing option is not possible.")
1275 END IF
1276 END IF
1277 END IF
1278
1279 ! Initialize the GAPW local densities and potentials
1280 IF (dft_control%qs_control%method_id == do_method_gapw .OR. &
1281 dft_control%qs_control%method_id == do_method_gapw_xc) THEN
1282 ! *** Allocate and initialize the set of atomic densities ***
1283 NULLIFY (rho_atom_set)
1284 gapw_control => dft_control%qs_control%gapw_control
1285 CALL init_rho_atom(rho_atom_set, atomic_kind_set, qs_kind_set, dft_control, para_env)
1286 CALL set_qs_env(qs_env=qs_env, rho_atom_set=rho_atom_set)
1287 IF (dft_control%qs_control%method_id /= do_method_gapw_xc) THEN
1288 CALL get_qs_env(qs_env=qs_env, local_rho_set=local_rho_set, natom=natom)
1289 ! *** Allocate and initialize the compensation density rho0 ***
1290 CALL init_rho0(local_rho_set, qs_env, gapw_control)
1291 ! *** Allocate and Initialize the local coulomb term ***
1292 CALL init_coulomb_local(qs_env%hartree_local, natom)
1293 END IF
1294 ! NLCC
1295 CALL init_gapw_nlcc(qs_kind_set)
1296 ELSE IF (dft_control%qs_control%method_id == do_method_lrigpw) THEN
1297 ! allocate local ri environment
1298 ! nothing to do here?
1299 ELSE IF (dft_control%qs_control%method_id == do_method_rigpw) THEN
1300 ! allocate ri environment
1301 ! nothing to do here?
1302 ELSE IF (dft_control%qs_control%semi_empirical) THEN
1303 NULLIFY (se_store_int_env, se_nddo_mpole, se_nonbond_env)
1304 natom = SIZE(particle_set)
1305 se_section => section_vals_get_subs_vals(qs_section, "SE")
1306 se_control => dft_control%qs_control%se_control
1307
1308 ! Make the cutoff radii choice a bit smarter
1309 CALL se_cutoff_compatible(se_control, se_section, cell, output_unit)
1310
1311 SELECT CASE (dft_control%qs_control%method_id)
1312 CASE DEFAULT
1315 ! Neighbor lists have to be MAX(interaction range, orbital range)
1316 ! set new kind radius
1317 CALL init_se_nlradius(se_control, atomic_kind_set, qs_kind_set, subsys_section)
1318 END SELECT
1319 ! Initialize to zero the max multipole to treat in the EWALD scheme..
1320 se_control%max_multipole = do_multipole_none
1321 ! check for Ewald
1322 IF (se_control%do_ewald .OR. se_control%do_ewald_gks) THEN
1323 ALLOCATE (ewald_env)
1324 CALL ewald_env_create(ewald_env, para_env)
1325 poisson_section => section_vals_get_subs_vals(dft_section, "POISSON")
1326 CALL ewald_env_set(ewald_env, poisson_section=poisson_section)
1327 ewald_section => section_vals_get_subs_vals(poisson_section, "EWALD")
1328 print_section => section_vals_get_subs_vals(qs_env%input, &
1329 "PRINT%GRID_INFORMATION")
1330 CALL read_ewald_section(ewald_env, ewald_section)
1331 ! Create ewald grids
1332 ALLOCATE (ewald_pw)
1333 CALL ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, &
1334 print_section=print_section)
1335 ! Initialize ewald grids
1336 CALL ewald_pw_grid_update(ewald_pw, ewald_env, cell%hmat)
1337 ! Setup the nonbond environment (real space part of Ewald)
1338 CALL ewald_env_get(ewald_env, rcut=ewald_rcut)
1339 ! Setup the maximum level of multipoles to be treated in the periodic SE scheme
1340 IF (se_control%do_ewald) THEN
1341 CALL ewald_env_get(ewald_env, max_multipole=se_control%max_multipole)
1342 END IF
1343 CALL section_vals_val_get(se_section, "NEIGHBOR_LISTS%VERLET_SKIN", &
1344 r_val=verlet_skin)
1345 ALLOCATE (se_nonbond_env)
1346 CALL fist_nonbond_env_create(se_nonbond_env, atomic_kind_set, do_nonbonded=.true., &
1347 do_electrostatics=.true., verlet_skin=verlet_skin, ewald_rcut=ewald_rcut, &
1348 ei_scale14=0.0_dp, vdw_scale14=0.0_dp, shift_cutoff=.false.)
1349 ! Create and Setup NDDO multipole environment
1350 CALL nddo_mpole_setup(se_nddo_mpole, natom)
1351 CALL set_qs_env(qs_env, ewald_env=ewald_env, ewald_pw=ewald_pw, &
1352 se_nonbond_env=se_nonbond_env, se_nddo_mpole=se_nddo_mpole)
1353 ! Handle the residual integral part 1/R^3
1354 CALL semi_empirical_expns3_setup(qs_kind_set, se_control, &
1355 dft_control%qs_control%method_id)
1356 END IF
1357 ! Taper function
1358 CALL se_taper_create(se_taper, se_control%integral_screening, se_control%do_ewald, &
1359 se_control%taper_cou, se_control%range_cou, &
1360 se_control%taper_exc, se_control%range_exc, &
1361 se_control%taper_scr, se_control%range_scr, &
1362 se_control%taper_lrc, se_control%range_lrc)
1363 CALL set_qs_env(qs_env, se_taper=se_taper)
1364 ! Store integral environment
1365 CALL semi_empirical_si_create(se_store_int_env, se_section)
1366 CALL set_qs_env(qs_env, se_store_int_env=se_store_int_env)
1367 END IF
1368
1369 ! Initialize possible dispersion parameters
1370 IF (dft_control%qs_control%method_id == do_method_gpw .OR. &
1371 dft_control%qs_control%method_id == do_method_gapw .OR. &
1372 dft_control%qs_control%method_id == do_method_gapw_xc .OR. &
1373 dft_control%qs_control%method_id == do_method_lrigpw .OR. &
1374 dft_control%qs_control%method_id == do_method_rigpw .OR. &
1375 dft_control%qs_control%method_id == do_method_ofgpw) THEN
1376 ALLOCATE (dispersion_env)
1377 NULLIFY (xc_section)
1378 xc_section => section_vals_get_subs_vals(dft_section, "XC")
1379 CALL qs_dispersion_env_set(dispersion_env, xc_section)
1380 IF (dispersion_env%type == xc_vdw_fun_pairpot) THEN
1381 NULLIFY (pp_section)
1382 pp_section => section_vals_get_subs_vals(xc_section, "VDW_POTENTIAL%PAIR_POTENTIAL")
1383 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, pp_section, para_env)
1384 ELSE IF (dispersion_env%type == xc_vdw_fun_nonloc) THEN
1385 NULLIFY (nl_section)
1386 nl_section => section_vals_get_subs_vals(xc_section, "VDW_POTENTIAL%NON_LOCAL")
1387 CALL qs_dispersion_nonloc_init(dispersion_env, para_env)
1388 END IF
1389 CALL set_qs_env(qs_env, dispersion_env=dispersion_env)
1390 ELSE IF (dft_control%qs_control%method_id == do_method_dftb) THEN
1391 ALLOCATE (dispersion_env)
1392 ! set general defaults
1393 dispersion_env%doabc = .false.
1394 dispersion_env%c9cnst = .false.
1395 dispersion_env%lrc = .false.
1396 dispersion_env%srb = .false.
1397 dispersion_env%verbose = .false.
1398 NULLIFY (dispersion_env%c6ab, dispersion_env%maxci, dispersion_env%r0ab, dispersion_env%rcov, &
1399 dispersion_env%r2r4, dispersion_env%cn, dispersion_env%cnkind, dispersion_env%cnlist, &
1400 dispersion_env%d3_exclude_pair)
1401 NULLIFY (dispersion_env%q_mesh, dispersion_env%kernel, dispersion_env%d2phi_dk2, &
1402 dispersion_env%d2y_dx2, dispersion_env%dftd_section)
1403 NULLIFY (dispersion_env%sab_vdw, dispersion_env%sab_cn)
1404 IF (dftb_control%dispersion .AND. dftb_control%dispersion_type == dispersion_d3) THEN
1405 dispersion_env%type = xc_vdw_fun_pairpot
1406 dispersion_env%pp_type = vdw_pairpot_dftd3
1407 dispersion_env%eps_cn = dftb_control%epscn
1408 dispersion_env%s6 = dftb_control%sd3(1)
1409 dispersion_env%sr6 = dftb_control%sd3(2)
1410 dispersion_env%s8 = dftb_control%sd3(3)
1411 dispersion_env%domol = .false.
1412 dispersion_env%kgc8 = 0._dp
1413 dispersion_env%rc_disp = dftb_control%rcdisp
1414 dispersion_env%exp_pre = 0._dp
1415 dispersion_env%scaling = 0._dp
1416 dispersion_env%nd3_exclude_pair = 0
1417 dispersion_env%parameter_file_name = dftb_control%dispersion_parameter_file
1418 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, para_env=para_env)
1419 ELSEIF (dftb_control%dispersion .AND. dftb_control%dispersion_type == dispersion_d3bj) THEN
1420 dispersion_env%type = xc_vdw_fun_pairpot
1421 dispersion_env%pp_type = vdw_pairpot_dftd3bj
1422 dispersion_env%eps_cn = dftb_control%epscn
1423 dispersion_env%s6 = dftb_control%sd3bj(1)
1424 dispersion_env%a1 = dftb_control%sd3bj(2)
1425 dispersion_env%s8 = dftb_control%sd3bj(3)
1426 dispersion_env%a2 = dftb_control%sd3bj(4)
1427 dispersion_env%domol = .false.
1428 dispersion_env%kgc8 = 0._dp
1429 dispersion_env%rc_disp = dftb_control%rcdisp
1430 dispersion_env%exp_pre = 0._dp
1431 dispersion_env%scaling = 0._dp
1432 dispersion_env%nd3_exclude_pair = 0
1433 dispersion_env%parameter_file_name = dftb_control%dispersion_parameter_file
1434 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, para_env=para_env)
1435 ELSEIF (dftb_control%dispersion .AND. dftb_control%dispersion_type == dispersion_d2) THEN
1436 dispersion_env%type = xc_vdw_fun_pairpot
1437 dispersion_env%pp_type = vdw_pairpot_dftd2
1438 dispersion_env%exp_pre = dftb_control%exp_pre
1439 dispersion_env%scaling = dftb_control%scaling
1440 dispersion_env%parameter_file_name = dftb_control%dispersion_parameter_file
1441 dispersion_env%rc_disp = dftb_control%rcdisp
1442 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, para_env=para_env)
1443 ELSE
1444 dispersion_env%type = xc_vdw_fun_none
1445 END IF
1446 CALL set_qs_env(qs_env, dispersion_env=dispersion_env)
1447 ELSE IF (dft_control%qs_control%method_id == do_method_xtb) THEN
1448 IF (.NOT. (dft_control%qs_control%xtb_control%do_tblite)) THEN
1449 ALLOCATE (dispersion_env)
1450 ! set general defaults
1451 dispersion_env%doabc = .false.
1452 dispersion_env%c9cnst = .false.
1453 dispersion_env%lrc = .false.
1454 dispersion_env%srb = .false.
1455 dispersion_env%verbose = .false.
1456 NULLIFY (dispersion_env%c6ab, dispersion_env%maxci, &
1457 dispersion_env%r0ab, dispersion_env%rcov, &
1458 dispersion_env%r2r4, dispersion_env%cn, &
1459 dispersion_env%cnkind, dispersion_env%cnlist, &
1460 dispersion_env%d3_exclude_pair)
1461 NULLIFY (dispersion_env%q_mesh, dispersion_env%kernel, dispersion_env%d2phi_dk2, &
1462 dispersion_env%d2y_dx2, dispersion_env%dftd_section)
1463 NULLIFY (dispersion_env%sab_vdw, dispersion_env%sab_cn)
1464 dispersion_env%type = xc_vdw_fun_pairpot
1465 dispersion_env%eps_cn = xtb_control%epscn
1466 dispersion_env%s6 = xtb_control%s6
1467 dispersion_env%s8 = xtb_control%s8
1468 dispersion_env%a1 = xtb_control%a1
1469 dispersion_env%a2 = xtb_control%a2
1470 dispersion_env%domol = .false.
1471 dispersion_env%kgc8 = 0._dp
1472 dispersion_env%rc_disp = xtb_control%rcdisp
1473 dispersion_env%rc_d4 = xtb_control%rcdisp
1474 dispersion_env%exp_pre = 0._dp
1475 dispersion_env%scaling = 0._dp
1476 dispersion_env%nd3_exclude_pair = 0
1477 dispersion_env%parameter_file_name = xtb_control%dispersion_parameter_file
1478 !
1479 SELECT CASE (xtb_control%vdw_type)
1481 dispersion_env%pp_type = vdw_pairpot_dftd3bj
1482 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, para_env=para_env)
1483 IF (xtb_control%vdw_type == xtb_vdw_type_none) dispersion_env%type = xc_vdw_fun_none
1484 CASE (xtb_vdw_type_d4)
1485 dispersion_env%pp_type = vdw_pairpot_dftd4
1486 dispersion_env%ref_functional = "none"
1487 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, &
1488 dispersion_env, para_env=para_env)
1489 dispersion_env%cnfun = 2
1490 CASE DEFAULT
1491 cpabort("vdw type")
1492 END SELECT
1493 CALL set_qs_env(qs_env, dispersion_env=dispersion_env)
1494 END IF
1495 ELSE IF (dft_control%qs_control%semi_empirical) THEN
1496 ALLOCATE (dispersion_env)
1497 ! set general defaults
1498 dispersion_env%doabc = .false.
1499 dispersion_env%c9cnst = .false.
1500 dispersion_env%lrc = .false.
1501 dispersion_env%srb = .false.
1502 dispersion_env%verbose = .false.
1503 NULLIFY (dispersion_env%c6ab, dispersion_env%maxci, dispersion_env%r0ab, dispersion_env%rcov, &
1504 dispersion_env%r2r4, dispersion_env%cn, dispersion_env%cnkind, dispersion_env%cnlist, &
1505 dispersion_env%d3_exclude_pair)
1506 NULLIFY (dispersion_env%q_mesh, dispersion_env%kernel, dispersion_env%d2phi_dk2, &
1507 dispersion_env%d2y_dx2, dispersion_env%dftd_section)
1508 NULLIFY (dispersion_env%sab_vdw, dispersion_env%sab_cn)
1509 IF (se_control%dispersion) THEN
1510 dispersion_env%type = xc_vdw_fun_pairpot
1511 dispersion_env%pp_type = vdw_pairpot_dftd3
1512 dispersion_env%eps_cn = se_control%epscn
1513 dispersion_env%s6 = se_control%sd3(1)
1514 dispersion_env%sr6 = se_control%sd3(2)
1515 dispersion_env%s8 = se_control%sd3(3)
1516 dispersion_env%domol = .false.
1517 dispersion_env%kgc8 = 0._dp
1518 dispersion_env%rc_disp = se_control%rcdisp
1519 dispersion_env%exp_pre = 0._dp
1520 dispersion_env%scaling = 0._dp
1521 dispersion_env%nd3_exclude_pair = 0
1522 dispersion_env%parameter_file_name = se_control%dispersion_parameter_file
1523 CALL qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, para_env=para_env)
1524 ELSE
1525 dispersion_env%type = xc_vdw_fun_none
1526 END IF
1527 CALL set_qs_env(qs_env, dispersion_env=dispersion_env)
1528 END IF
1529
1530 ! Initialize possible geomertical counterpoise correction potential
1531 IF (dft_control%qs_control%method_id == do_method_gpw .OR. &
1532 dft_control%qs_control%method_id == do_method_gapw .OR. &
1533 dft_control%qs_control%method_id == do_method_gapw_xc .OR. &
1534 dft_control%qs_control%method_id == do_method_lrigpw .OR. &
1535 dft_control%qs_control%method_id == do_method_rigpw .OR. &
1536 dft_control%qs_control%method_id == do_method_ofgpw) THEN
1537 ALLOCATE (gcp_env)
1538 NULLIFY (xc_section)
1539 xc_section => section_vals_get_subs_vals(dft_section, "XC")
1540 CALL qs_gcp_env_set(gcp_env, xc_section)
1541 CALL qs_gcp_init(qs_env, gcp_env)
1542 CALL set_qs_env(qs_env, gcp_env=gcp_env)
1543 END IF
1544
1545 ! *** Allocate the MO data types ***
1546 CALL get_qs_kind_set(qs_kind_set, nsgf=n_ao, nelectron=nelectron)
1547
1548 ! the total number of electrons
1549 nelectron = nelectron - dft_control%charge
1550
1551 IF (dft_control%multiplicity == 0) THEN
1552 IF (modulo(nelectron, 2) == 0) THEN
1553 dft_control%multiplicity = 1
1554 ELSE
1555 dft_control%multiplicity = 2
1556 END IF
1557 END IF
1558
1559 multiplicity = dft_control%multiplicity
1560
1561 IF ((dft_control%nspins < 1) .OR. (dft_control%nspins > 2)) THEN
1562 cpabort("nspins should be 1 or 2 for the time being ...")
1563 END IF
1564
1565 IF ((modulo(nelectron, 2) /= 0) .AND. (dft_control%nspins == 1)) THEN
1566 IF (.NOT. dft_control%qs_control%ofgpw .AND. .NOT. dft_control%smear) THEN
1567 cpabort("Use the LSD option for an odd number of electrons")
1568 END IF
1569 END IF
1570
1571 ! The transition potential method to calculate XAS needs LSD
1572 IF (dft_control%do_xas_calculation) THEN
1573 IF (dft_control%nspins == 1) THEN
1574 cpabort("Use the LSD option for XAS with transition potential")
1575 END IF
1576 END IF
1577
1578 ! assigning the number of states per spin initial version, not yet very
1579 ! general. Should work for an even number of electrons and a single
1580 ! additional electron this set of options that requires full matrices,
1581 ! however, makes things a bit ugly right now.... we try to make a
1582 ! distinction between the number of electrons per spin and the number of
1583 ! MOs per spin this should allow the use of fractional occupations later
1584 ! on
1585 IF (dft_control%qs_control%ofgpw) THEN
1586
1587 IF (dft_control%nspins == 1) THEN
1588 maxocc = nelectron
1589 nelectron_spin(1) = nelectron
1590 nelectron_spin(2) = 0
1591 n_mo(1) = 1
1592 n_mo(2) = 0
1593 ELSE
1594 IF (modulo(nelectron + multiplicity - 1, 2) /= 0) THEN
1595 cpabort("LSD: try to use a different multiplicity")
1596 END IF
1597 nelectron_spin(1) = (nelectron + multiplicity - 1)/2
1598 nelectron_spin(2) = (nelectron - multiplicity + 1)/2
1599 IF (nelectron_spin(1) < 0) THEN
1600 cpabort("LSD: too few electrons for this multiplicity")
1601 END IF
1602 maxocc = maxval(nelectron_spin)
1603 n_mo(1) = min(nelectron_spin(1), 1)
1604 n_mo(2) = min(nelectron_spin(2), 1)
1605 END IF
1606
1607 ELSE
1608
1609 IF (dft_control%nspins == 1) THEN
1610 maxocc = 2.0_dp
1611 nelectron_spin(1) = nelectron
1612 nelectron_spin(2) = 0
1613 IF (modulo(nelectron, 2) == 0) THEN
1614 n_mo(1) = nelectron/2
1615 ELSE
1616 n_mo(1) = int(nelectron/2._dp) + 1
1617 END IF
1618 n_mo(2) = 0
1619 ELSE
1620 maxocc = 1.0_dp
1621
1622 ! The simplist spin distribution is written here. Special cases will
1623 ! need additional user input
1624 IF (modulo(nelectron + multiplicity - 1, 2) /= 0) THEN
1625 cpabort("LSD: try to use a different multiplicity")
1626 END IF
1627
1628 nelectron_spin(1) = (nelectron + multiplicity - 1)/2
1629 nelectron_spin(2) = (nelectron - multiplicity + 1)/2
1630
1631 IF (nelectron_spin(2) < 0) THEN
1632 cpabort("LSD: too few electrons for this multiplicity")
1633 END IF
1634
1635 n_mo(1) = nelectron_spin(1)
1636 n_mo(2) = nelectron_spin(2)
1637
1638 END IF
1639
1640 END IF
1641
1642 ! Read the total_zeff_corr here [SGh]
1643 CALL get_qs_kind_set(qs_kind_set, total_zeff_corr=total_zeff_corr)
1644 ! store it in qs_env
1645 qs_env%total_zeff_corr = total_zeff_corr
1646
1647 ! store the number of electrons once an for all
1648 CALL qs_subsys_set(subsys, &
1649 nelectron_total=nelectron, &
1650 nelectron_spin=nelectron_spin)
1651
1652 ! Check and set number of added (unoccupied) MOs
1653 IF (dft_control%nspins == 2) THEN
1654 IF (scf_control%added_mos(2) < 0) THEN
1655 n_mo_add = n_ao - n_mo(2) ! use all available MOs
1656 ELSEIF (scf_control%added_mos(2) > 0) THEN
1657 n_mo_add = scf_control%added_mos(2)
1658 ELSE
1659 n_mo_add = scf_control%added_mos(1)
1660 END IF
1661 IF (n_mo_add > n_ao - n_mo(2)) THEN
1662 cpwarn("More ADDED_MOs requested for beta spin than available.")
1663 END IF
1664 scf_control%added_mos(2) = min(n_mo_add, n_ao - n_mo(2))
1665 n_mo(2) = n_mo(2) + scf_control%added_mos(2)
1666 END IF
1667
1668 ! proceed alpha orbitals after the beta orbitals; this is essential to avoid
1669 ! reduction in the number of available unoccupied molecular orbitals.
1670 ! E.g. n_ao = 10, nelectrons = 10, multiplicity = 3 implies n_mo(1) = 6, n_mo(2) = 4;
1671 ! added_mos(1:2) = (6,undef) should increase the number of molecular orbitals as
1672 ! n_mo(1) = min(n_ao, n_mo(1) + added_mos(1)) = 10, n_mo(2) = 10.
1673 ! However, if we try to proceed alpha orbitals first, this leads us n_mo(1:2) = (10,8)
1674 ! due to the following assignment instruction above:
1675 ! IF (scf_control%added_mos(2) > 0) THEN ... ELSE; n_mo_add = scf_control%added_mos(1); END IF
1676 IF (scf_control%added_mos(1) < 0) THEN
1677 scf_control%added_mos(1) = n_ao - n_mo(1) ! use all available MOs
1678 ELSEIF (scf_control%added_mos(1) > n_ao - n_mo(1)) THEN
1679 CALL cp_warn(__location__, &
1680 "More added MOs requested than available. "// &
1681 "The full set of unoccupied MOs will be used. "// &
1682 "Use 'ADDED_MOS -1' to always use all available MOs "// &
1683 "and to get rid of this warning.")
1684 END IF
1685 scf_control%added_mos(1) = min(scf_control%added_mos(1), n_ao - n_mo(1))
1686 n_mo(1) = n_mo(1) + scf_control%added_mos(1)
1687
1688 IF (dft_control%nspins == 2) THEN
1689 IF (n_mo(2) > n_mo(1)) &
1690 CALL cp_warn(__location__, &
1691 "More beta than alpha MOs requested. "// &
1692 "The number of beta MOs will be reduced to the number alpha MOs.")
1693 n_mo(2) = min(n_mo(1), n_mo(2))
1694 cpassert(n_mo(1) >= nelectron_spin(1))
1695 cpassert(n_mo(2) >= nelectron_spin(2))
1696 END IF
1697
1698 ! kpoints
1699 CALL get_qs_env(qs_env=qs_env, do_kpoints=do_kpoints)
1700 IF (do_kpoints .AND. dft_control%nspins == 2) THEN
1701 ! we need equal number of calculated states
1702 IF (n_mo(2) /= n_mo(1)) &
1703 CALL cp_warn(__location__, &
1704 "Kpoints: Different number of MOs requested. "// &
1705 "The number of beta MOs will be set to the number alpha MOs.")
1706 n_mo(2) = n_mo(1)
1707 cpassert(n_mo(1) >= nelectron_spin(1))
1708 cpassert(n_mo(2) >= nelectron_spin(2))
1709 END IF
1710
1711 ! Compatibility checks for smearing
1712 IF (scf_control%smear%do_smear) THEN
1713 IF (scf_control%added_mos(1) == 0) THEN
1714 cpabort("Extra MOs (ADDED_MOS) are required for smearing")
1715 END IF
1716 END IF
1717
1718 ! *** Some options require that all MOs are computed ... ***
1719 IF (btest(cp_print_key_should_output(logger%iter_info, dft_section, &
1720 "PRINT%MO/CARTESIAN"), &
1721 cp_p_file) .OR. &
1722 (scf_control%level_shift /= 0.0_dp) .OR. &
1723 (scf_control%diagonalization%eps_jacobi /= 0.0_dp) .OR. &
1724 (dft_control%roks .AND. (.NOT. scf_control%use_ot))) THEN
1725 n_mo(:) = n_ao
1726 END IF
1727
1728 ! Compatibility checks for ROKS
1729 IF (dft_control%roks .AND. (.NOT. scf_control%use_ot)) THEN
1730 IF (scf_control%roks_scheme == general_roks) THEN
1731 cpwarn("General ROKS scheme is not yet tested!")
1732 END IF
1733 IF (scf_control%smear%do_smear) THEN
1734 CALL cp_abort(__location__, &
1735 "The options ROKS and SMEAR are not compatible. "// &
1736 "Try UKS instead of ROKS")
1737 END IF
1738 END IF
1739 IF (dft_control%low_spin_roks) THEN
1740 SELECT CASE (dft_control%qs_control%method_id)
1741 CASE DEFAULT
1743 CALL cp_abort(__location__, &
1744 "xTB/DFTB methods are not compatible with low spin ROKS.")
1747 CALL cp_abort(__location__, &
1748 "SE methods are not compatible with low spin ROKS.")
1749 END SELECT
1750 END IF
1751
1752 ! in principle the restricted calculation could be performed
1753 ! using just one set of MOs and special casing most of the code
1754 ! right now we'll just take care of what is effectively an additional constraint
1755 ! at as few places as possible, just duplicating the beta orbitals
1756 IF (dft_control%restricted .AND. (output_unit > 0)) THEN
1757 ! it is really not yet tested till the end ! Joost
1758 WRITE (output_unit, *) ""
1759 WRITE (output_unit, *) " **************************************"
1760 WRITE (output_unit, *) " restricted calculation cutting corners"
1761 WRITE (output_unit, *) " experimental feature, check code "
1762 WRITE (output_unit, *) " **************************************"
1763 END IF
1764
1765 ! no point in allocating these things here ?
1766 IF (dft_control%qs_control%do_ls_scf) THEN
1767 NULLIFY (mos)
1768 ELSE
1769 ALLOCATE (mos(dft_control%nspins))
1770 DO ispin = 1, dft_control%nspins
1771 CALL allocate_mo_set(mo_set=mos(ispin), &
1772 nao=n_ao, &
1773 nmo=n_mo(ispin), &
1774 nelectron=nelectron_spin(ispin), &
1775 n_el_f=real(nelectron_spin(ispin), dp), &
1776 maxocc=maxocc, &
1777 flexible_electron_count=dft_control%relax_multiplicity)
1778 END DO
1779 END IF
1780
1781 CALL set_qs_env(qs_env, mos=mos)
1782
1783 ! allocate mos when switch_surf_dip is triggered [SGh]
1784 IF (dft_control%switch_surf_dip) THEN
1785 ALLOCATE (mos_last_converged(dft_control%nspins))
1786 DO ispin = 1, dft_control%nspins
1787 CALL allocate_mo_set(mo_set=mos_last_converged(ispin), &
1788 nao=n_ao, &
1789 nmo=n_mo(ispin), &
1790 nelectron=nelectron_spin(ispin), &
1791 n_el_f=real(nelectron_spin(ispin), dp), &
1792 maxocc=maxocc, &
1793 flexible_electron_count=dft_control%relax_multiplicity)
1794 END DO
1795 CALL set_qs_env(qs_env, mos_last_converged=mos_last_converged)
1796 END IF
1797
1798 IF (.NOT. be_silent) THEN
1799 ! Print the DFT control parameters
1800 CALL write_dft_control(dft_control, dft_section)
1801
1802 ! Print the vdW control parameters
1803 IF (dft_control%qs_control%method_id == do_method_gpw .OR. &
1804 dft_control%qs_control%method_id == do_method_gapw .OR. &
1805 dft_control%qs_control%method_id == do_method_gapw_xc .OR. &
1806 dft_control%qs_control%method_id == do_method_lrigpw .OR. &
1807 dft_control%qs_control%method_id == do_method_rigpw .OR. &
1808 dft_control%qs_control%method_id == do_method_dftb .OR. &
1809 (dft_control%qs_control%method_id == do_method_xtb .AND. &
1810 (.NOT. dft_control%qs_control%xtb_control%do_tblite)) .OR. &
1811 dft_control%qs_control%method_id == do_method_ofgpw) THEN
1812 CALL get_qs_env(qs_env, dispersion_env=dispersion_env)
1813 CALL qs_write_dispersion(qs_env, dispersion_env)
1814 END IF
1815
1816 ! Print the Quickstep control parameters
1817 CALL write_qs_control(dft_control%qs_control, dft_section)
1818
1819 ! Print the ADMM control parameters
1820 IF (dft_control%do_admm) THEN
1821 CALL write_admm_control(dft_control%admm_control, dft_section)
1822 END IF
1823
1824 ! Print XES/XAS control parameters
1825 IF (dft_control%do_xas_calculation) THEN
1826 CALL cite_reference(iannuzzi2007)
1827 !CALL write_xas_control(dft_control%xas_control,dft_section)
1828 END IF
1829
1830 ! Print the unnormalized basis set information (input data)
1831 CALL write_gto_basis_sets(qs_kind_set, subsys_section)
1832
1833 ! Print the atomic kind set
1834 CALL write_qs_kind_set(qs_kind_set, subsys_section)
1835
1836 ! Print the molecule kind set
1837 CALL write_molecule_kind_set(molecule_kind_set, subsys_section)
1838
1839 ! Print the total number of kinds, atoms, basis functions etc.
1840 CALL write_total_numbers(qs_kind_set, particle_set, qs_env%input)
1841
1842 ! Print the atomic coordinates
1843 CALL write_qs_particle_coordinates(particle_set, qs_kind_set, subsys_section, label="QUICKSTEP")
1844
1845 ! Print the interatomic distances
1846 CALL write_particle_distances(particle_set, cell, subsys_section)
1847
1848 ! Print the requested structure data
1849 CALL write_structure_data(particle_set, cell, subsys_section)
1850
1851 ! Print symmetry information
1852 CALL write_symmetry(particle_set, cell, subsys_section)
1853
1854 ! Print the SCF parameters
1855 IF ((.NOT. dft_control%qs_control%do_ls_scf) .AND. &
1856 (.NOT. dft_control%qs_control%do_almo_scf)) THEN
1857 CALL scf_c_write_parameters(scf_control, dft_section)
1858 END IF
1859 END IF
1860
1861 ! Sets up pw_env, qs_charges, mpools ...
1862 CALL qs_env_setup(qs_env)
1863
1864 ! Allocate and initialise rho0 soft on the global grid
1865 IF (dft_control%qs_control%method_id == do_method_gapw) THEN
1866 CALL get_qs_env(qs_env=qs_env, pw_env=pw_env, rho0_mpole=rho0_mpole)
1867 CALL rho0_s_grid_create(pw_env, rho0_mpole)
1868 END IF
1869
1870 IF (output_unit > 0) CALL m_flush(output_unit)
1871 CALL timestop(handle)
1872
1873 END SUBROUTINE qs_init_subsys
1874
1875! **************************************************************************************************
1876!> \brief Write the total number of kinds, atoms, etc. to the logical unit
1877!> number lunit.
1878!> \param qs_kind_set ...
1879!> \param particle_set ...
1880!> \param force_env_section ...
1881!> \author Creation (06.10.2000)
1882! **************************************************************************************************
1883 SUBROUTINE write_total_numbers(qs_kind_set, particle_set, force_env_section)
1884
1885 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1886 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1887 TYPE(section_vals_type), POINTER :: force_env_section
1888
1889 INTEGER :: maxlgto, maxlppl, maxlppnl, natom, &
1890 natom_q, ncgf, nkind, nkind_q, npgf, &
1891 nset, nsgf, nshell, output_unit
1892 TYPE(cp_logger_type), POINTER :: logger
1893
1894 NULLIFY (logger)
1895 logger => cp_get_default_logger()
1896 output_unit = cp_print_key_unit_nr(logger, force_env_section, "PRINT%TOTAL_NUMBERS", &
1897 extension=".Log")
1898
1899 IF (output_unit > 0) THEN
1900 natom = SIZE(particle_set)
1901 nkind = SIZE(qs_kind_set)
1902
1903 CALL get_qs_kind_set(qs_kind_set, &
1904 maxlgto=maxlgto, &
1905 ncgf=ncgf, &
1906 npgf=npgf, &
1907 nset=nset, &
1908 nsgf=nsgf, &
1909 nshell=nshell, &
1910 maxlppl=maxlppl, &
1911 maxlppnl=maxlppnl)
1912
1913 WRITE (unit=output_unit, fmt="(/,/,T2,A)") &
1914 "TOTAL NUMBERS AND MAXIMUM NUMBERS"
1915
1916 IF (nset + npgf + ncgf > 0) THEN
1917 WRITE (unit=output_unit, fmt="(/,T3,A,(T30,A,T71,I10))") &
1918 "Total number of", &
1919 "- Atomic kinds: ", nkind, &
1920 "- Atoms: ", natom, &
1921 "- Shell sets: ", nset, &
1922 "- Shells: ", nshell, &
1923 "- Primitive Cartesian functions: ", npgf, &
1924 "- Cartesian basis functions: ", ncgf, &
1925 "- Spherical basis functions: ", nsgf
1926 ELSE IF (nshell + nsgf > 0) THEN
1927 WRITE (unit=output_unit, fmt="(/,T3,A,(T30,A,T71,I10))") &
1928 "Total number of", &
1929 "- Atomic kinds: ", nkind, &
1930 "- Atoms: ", natom, &
1931 "- Shells: ", nshell, &
1932 "- Spherical basis functions: ", nsgf
1933 ELSE
1934 WRITE (unit=output_unit, fmt="(/,T3,A,(T30,A,T71,I10))") &
1935 "Total number of", &
1936 "- Atomic kinds: ", nkind, &
1937 "- Atoms: ", natom
1938 END IF
1939
1940 IF ((maxlppl > -1) .AND. (maxlppnl > -1)) THEN
1941 WRITE (unit=output_unit, fmt="(/,T3,A,(T30,A,T75,I6))") &
1942 "Maximum angular momentum of the", &
1943 "- Orbital basis functions: ", maxlgto, &
1944 "- Local part of the GTH pseudopotential: ", maxlppl, &
1945 "- Non-local part of the GTH pseudopotential: ", maxlppnl
1946 ELSEIF (maxlppl > -1) THEN
1947 WRITE (unit=output_unit, fmt="(/,T3,A,(T30,A,T75,I6))") &
1948 "Maximum angular momentum of the", &
1949 "- Orbital basis functions: ", maxlgto, &
1950 "- Local part of the GTH pseudopotential: ", maxlppl
1951 ELSE
1952 WRITE (unit=output_unit, fmt="(/,T3,A,T75,I6)") &
1953 "Maximum angular momentum of the orbital basis functions: ", maxlgto
1954 END IF
1955
1956 ! LRI_AUX BASIS
1957 CALL get_qs_kind_set(qs_kind_set, &
1958 maxlgto=maxlgto, &
1959 ncgf=ncgf, &
1960 npgf=npgf, &
1961 nset=nset, &
1962 nsgf=nsgf, &
1963 nshell=nshell, &
1964 basis_type="LRI_AUX")
1965 IF (nset + npgf + ncgf > 0) THEN
1966 WRITE (unit=output_unit, fmt="(/,T3,A,/,T3,A,(T30,A,T71,I10))") &
1967 "LRI_AUX Basis: ", &
1968 "Total number of", &
1969 "- Shell sets: ", nset, &
1970 "- Shells: ", nshell, &
1971 "- Primitive Cartesian functions: ", npgf, &
1972 "- Cartesian basis functions: ", ncgf, &
1973 "- Spherical basis functions: ", nsgf
1974 WRITE (unit=output_unit, fmt="(T30,A,T75,I6)") &
1975 " Maximum angular momentum ", maxlgto
1976 END IF
1977
1978 ! RI_HXC BASIS
1979 CALL get_qs_kind_set(qs_kind_set, &
1980 maxlgto=maxlgto, &
1981 ncgf=ncgf, &
1982 npgf=npgf, &
1983 nset=nset, &
1984 nsgf=nsgf, &
1985 nshell=nshell, &
1986 basis_type="RI_HXC")
1987 IF (nset + npgf + ncgf > 0) THEN
1988 WRITE (unit=output_unit, fmt="(/,T3,A,/,T3,A,(T30,A,T71,I10))") &
1989 "RI_HXC Basis: ", &
1990 "Total number of", &
1991 "- Shell sets: ", nset, &
1992 "- Shells: ", nshell, &
1993 "- Primitive Cartesian functions: ", npgf, &
1994 "- Cartesian basis functions: ", ncgf, &
1995 "- Spherical basis functions: ", nsgf
1996 WRITE (unit=output_unit, fmt="(T30,A,T75,I6)") &
1997 " Maximum angular momentum ", maxlgto
1998 END IF
1999
2000 ! AUX_FIT BASIS
2001 CALL get_qs_kind_set(qs_kind_set, &
2002 maxlgto=maxlgto, &
2003 ncgf=ncgf, &
2004 npgf=npgf, &
2005 nset=nset, &
2006 nsgf=nsgf, &
2007 nshell=nshell, &
2008 basis_type="AUX_FIT")
2009 IF (nset + npgf + ncgf > 0) THEN
2010 WRITE (unit=output_unit, fmt="(/,T3,A,/,T3,A,(T30,A,T71,I10))") &
2011 "AUX_FIT ADMM-Basis: ", &
2012 "Total number of", &
2013 "- Shell sets: ", nset, &
2014 "- Shells: ", nshell, &
2015 "- Primitive Cartesian functions: ", npgf, &
2016 "- Cartesian basis functions: ", ncgf, &
2017 "- Spherical basis functions: ", nsgf
2018 WRITE (unit=output_unit, fmt="(T30,A,T75,I6)") &
2019 " Maximum angular momentum ", maxlgto
2020 END IF
2021
2022 ! NUCLEAR BASIS
2023 CALL get_qs_kind_set(qs_kind_set, &
2024 nkind_q=nkind_q, &
2025 natom_q=natom_q, &
2026 maxlgto=maxlgto, &
2027 ncgf=ncgf, &
2028 npgf=npgf, &
2029 nset=nset, &
2030 nsgf=nsgf, &
2031 nshell=nshell, &
2032 basis_type="NUC")
2033 IF (nset + npgf + ncgf > 0) THEN
2034 WRITE (unit=output_unit, fmt="(/,T3,A,/,T3,A,(T30,A,T71,I10))") &
2035 "Nuclear Basis: ", &
2036 "Total number of", &
2037 "- Quantum atomic kinds: ", nkind_q, &
2038 "- Quantum atoms: ", natom_q, &
2039 "- Shell sets: ", nset, &
2040 "- Shells: ", nshell, &
2041 "- Primitive Cartesian functions: ", npgf, &
2042 "- Cartesian basis functions: ", ncgf, &
2043 "- Spherical basis functions: ", nsgf
2044 WRITE (unit=output_unit, fmt="(T30,A,T75,I6)") &
2045 " Maximum angular momentum ", maxlgto
2046 END IF
2047
2048 END IF
2049 CALL cp_print_key_finished_output(output_unit, logger, force_env_section, &
2050 "PRINT%TOTAL_NUMBERS")
2051
2052 END SUBROUTINE write_total_numbers
2053
2054END MODULE qs_environment
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
almo_scf_env methods
subroutine, public almo_scf_env_create(qs_env)
Creation and basic initialization of the almo environment.
calculate the orbitals for a given atomic kind type
subroutine, public calculate_atomic_relkin(atomic_kind, qs_kind, rel_control, rtmat)
...
Define the atomic kind types and their sub types.
Automatic generation of auxiliary basis sets of different kind.
Definition auto_basis.F:15
subroutine, public create_lri_aux_basis_set(lri_aux_basis_set, qs_kind, basis_cntrl, exact_1c_terms, tda_kernel)
Create a LRI_AUX basis set using some heuristics.
Definition auto_basis.F:233
subroutine, public create_ri_aux_basis_set(ri_aux_basis_set, qs_kind, basis_cntrl, basis_type, basis_sort)
Create a RI_AUX basis set using some heuristics.
Definition auto_basis.F:57
subroutine, public add_basis_set_to_container(container, basis_set, basis_set_type)
...
integer, parameter, public basis_sort_zet
subroutine, public deallocate_gto_basis_set(gto_basis_set)
...
subroutine, public create_primitive_basis_set(basis_set, pbasis, lmax)
...
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public cp2kqs2020
integer, save, public iannuzzi2007
integer, save, public iannuzzi2006
Handles all functions related to the CELL.
Definition cell_types.F:15
methods related to the blacs parallel environment
subroutine, public cp_blacs_env_release(blacs_env)
releases the given blacs_env
subroutine, public cp_blacs_env_create(blacs_env, para_env, blacs_grid_layout, blacs_repeatable, row_major, grid_2d)
allocates and initializes a type that represent a blacs context
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
Utilities to set up the control types.
subroutine, public write_qs_control(qs_control, dft_section)
Purpose: Write the QS control parameters to the output unit.
subroutine, public read_rixs_control(rixs_control, rixs_section, qs_control)
Reads the input and stores in the rixs_control_type.
subroutine, public read_tddfpt2_control(t_control, t_section, qs_control)
Read TDDFPT-related input parameters.
subroutine, public read_dft_control(dft_control, dft_section)
...
subroutine, public read_qs_section(qs_control, qs_section)
...
subroutine, public write_admm_control(admm_control, dft_section)
Write the ADMM control parameters to the output unit.
subroutine, public write_dft_control(dft_control, dft_section)
Write the DFT control parameters to the output unit.
subroutine, public read_ddapc_section(qs_control, qs_section, ddapc_restraint_section)
reads the input parameters needed for ddapc.
subroutine, public read_mgrid_section(qs_control, dft_section)
...
contains information regarding the decoupling/recoupling method of Bloechl
subroutine, public cp_ddapc_ewald_create(cp_ddapc_ewald, qmmm_decoupl, qm_cell, force_env_section, subsys_section, para_env)
...
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
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...
types that represent a subsys, i.e. a part of the system
Work with symmetry.
Definition cp_symmetry.F:13
subroutine, public write_symmetry(particle_set, cell, input_section)
Write symmetry information to output.
Definition cp_symmetry.F:61
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
subroutine, public distribution_1d_release(distribution_1d)
releases the given distribution_1d
Distribution methods for atoms, particles, or molecules.
subroutine, public distribute_molecules_1d(atomic_kind_set, particle_set, local_particles, molecule_kind_set, molecule_set, local_molecules, force_env_section, prev_molecule_kind_set, prev_local_molecules)
Distribute molecules and particles.
Types needed for a for a Energy Correction.
Energy correction environment setup and handling.
subroutine, public ec_write_input(ec_env)
Print out the energy correction input section.
subroutine, public ec_env_create(qs_env, ec_env, dft_section, ec_section)
Allocates and intitializes ec_env.
Definition and initialisation of the et_coupling data type.
subroutine, public et_coupling_create(et_coupling)
...
subroutine, public ewald_env_set(ewald_env, ewald_type, alpha, epsilon, eps_pol, gmax, ns_max, precs, o_spline, para_env, poisson_section, interaction_cutoffs, cell_hmat)
Purpose: Set the EWALD environment.
subroutine, public ewald_env_create(ewald_env, para_env)
allocates and intitializes a ewald_env
subroutine, public read_ewald_section(ewald_env, ewald_section)
Purpose: read the EWALD section.
subroutine, public read_ewald_section_tb(ewald_env, ewald_section, hmat, silent, pset)
Purpose: read the EWALD section for TB methods.
subroutine, public ewald_env_get(ewald_env, ewald_type, alpha, eps_pol, epsilon, gmax, ns_max, o_spline, group, para_env, poisson_section, precs, rcut, do_multipoles, max_multipole, do_ipol, max_ipol_iter, interaction_cutoffs, cell_hmat)
Purpose: Get the EWALD environment.
subroutine, public ewald_pw_grid_update(ewald_pw, ewald_env, cell_hmat)
Rescales pw_grids for given box, if necessary.
subroutine, public ewald_pw_create(ewald_pw, ewald_env, cell, cell_ref, print_section)
creates the structure ewald_pw_type
Types for excited states potential energies.
subroutine, public exstate_create(ex_env, excited_state, dft_section)
Allocates and intitializes exstate_env.
Definition of the atomic potential types.
subroutine, public fist_nonbond_env_create(fist_nonbond_env, atomic_kind_set, potparm14, potparm, do_nonbonded, do_electrostatics, verlet_skin, ewald_rcut, ei_scale14, vdw_scale14, shift_cutoff)
allocates and intitializes a fist_nonbond_env
Calculation of the incomplete Gamma function F_n(t) for multi-center integrals over Cartesian Gaussia...
Definition gamma.F:15
subroutine, public init_md_ftable(nmax)
Initialize a table of F_n(t) values in the range 0 <= t <= 12 with a stepsize of 0....
Definition gamma.F:540
Define type storing the global information of a run. Keep the amount of stored data small....
subroutine, public init_coulomb_local(hartree_local, natom)
...
subroutine, public dftb_header(iw)
...
Definition header.F:193
subroutine, public qs_header(iw)
...
Definition header.F:310
subroutine, public se_header(iw)
...
Definition header.F:286
subroutine, public xtb_header(iw, gfn_type)
...
Definition header.F:217
subroutine, public tblite_header(iw, tb_type)
...
Definition header.F:252
Types and set/get functions for HFX.
Definition hfx_types.F:15
subroutine, public hfx_create(x_data, para_env, hfx_section, atomic_kind_set, qs_kind_set, particle_set, dft_control, cell, orb_basis, ri_basis, nelectron_total, nkp_grid)
This routine allocates and initializes all types in hfx_data
Definition hfx_types.F:595
subroutine, public compare_hfx_sections(hfx_section1, hfx_section2, is_identical, same_except_frac)
Compares the non-technical parts of two HFX input section and check whether they are the same Ignore ...
Definition hfx_types.F:2957
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_method_ofgpw
integer, parameter, public xc_vdw_fun_none
integer, parameter, public dispersion_d3
integer, parameter, public do_method_rigpw
integer, parameter, public do_method_gpw
integer, parameter, public vdw_pairpot_dftd3
integer, parameter, public do_method_pdg
integer, parameter, public wfi_linear_wf_method_nr
integer, parameter, public wfi_linear_ps_method_nr
integer, parameter, public do_method_pnnl
integer, parameter, public xc_vdw_fun_nonloc
integer, parameter, public vdw_pairpot_dftd4
integer, parameter, public xtb_vdw_type_d3
integer, parameter, public wfi_use_guess_method_nr
integer, parameter, public kg_tnadd_embed_ri
integer, parameter, public dispersion_d3bj
integer, parameter, public do_method_rm1
integer, parameter, public do_qmmm_swave
integer, parameter, public do_method_pm3
integer, parameter, public do_method_mndo
integer, parameter, public xtb_vdw_type_d4
integer, parameter, public do_method_gapw
integer, parameter, public do_method_mndod
integer, parameter, public rel_trans_atom
integer, parameter, public xtb_vdw_type_none
integer, parameter, public do_method_am1
integer, parameter, public do_method_dftb
integer, parameter, public do_qmmm_gauss
integer, parameter, public wfi_ps_method_nr
integer, parameter, public rel_none
integer, parameter, public general_roks
integer, parameter, public vdw_pairpot_dftd2
integer, parameter, public do_method_lrigpw
integer, parameter, public dispersion_d2
integer, parameter, public wfi_aspc_nr
integer, parameter, public do_method_xtb
integer, parameter, public do_method_pm6fm
integer, parameter, public xc_vdw_fun_pairpot
integer, parameter, public do_et_ddapc
integer, parameter, public vdw_pairpot_dftd3bj
integer, parameter, public do_method_gapw_xc
integer, parameter, public do_method_pm6
integer, parameter, public hden_atomic
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Routines for a Kim-Gordon-like partitioning into molecular subunits.
subroutine, public kg_env_create(qs_env, kg_env, qs_kind_set, input)
Allocates and intitializes kg_env.
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Routines needed for kpoint calculation.
subroutine, public kpoint_initialize_mos(kpoint, mos, added_mos, for_aux_fit)
Initialize a set of MOs and density matrix for each kpoint (kpoint group)
subroutine, public kpoint_initialize(kpoint, particle_set, cell)
Generate the kpoints and initialize the kpoint environment.
subroutine, public kpoint_env_initialize(kpoint, para_env, blacs_env, with_aux_fit)
Initialize the kpoint environment.
Types and basic routines needed for a kpoint calculation.
subroutine, public read_kpoint_section(kpoint, kpoint_section, a_vec)
Read the kpoint input section.
subroutine, public write_kpoint_info(kpoint, dft_section)
Write information on the kpoints to output.
subroutine, public kpoint_create(kpoint)
Create a kpoint environment.
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym)
Retrieve information from a kpoint environment.
initializes the environment for lri lri : local resolution of the identity
subroutine, public lri_env_init(lri_env, lri_section)
initializes the lri env
subroutine, public lri_env_basis(ri_type, qs_env, lri_env, qs_kind_set)
initializes the lri env
contains the types and subroutines for dealing with the lri_env lri : local resolution of the identit...
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:131
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
Interface to the message passing library MPI.
Define the molecule kind structure types and the corresponding functionality.
subroutine, public write_molecule_kind_set(molecule_kind_set, subsys_section)
Write a moleculeatomic kind set data set to the output unit.
Define the data structure for the molecule information.
Types needed for MP2 calculations.
Definition mp2_setup.F:14
subroutine, public read_mp2_section(input, mp2_env)
...
Definition mp2_setup.F:62
Types needed for MP2 calculations.
Definition mp2_types.F:14
subroutine, public mp2_env_create(mp2_env)
...
Definition mp2_types.F:471
Multipole structure: for multipole (fixed and induced) in FF based MD.
integer, parameter, public do_multipole_none
Provides Cartesian and spherical orbital pointers and indices.
subroutine, public init_orbital_pointers(maxl)
Initialize or update the orbital pointers.
Calculation of the spherical harmonics and the corresponding orbital transformation matrices.
subroutine, public init_spherical_harmonics(maxl, output_unit)
Initialize or update the orbital transformation matrices.
Define methods related to particle_type.
subroutine, public write_qs_particle_coordinates(particle_set, qs_kind_set, subsys_section, label)
Write the atomic coordinates to the output unit.
subroutine, public write_structure_data(particle_set, cell, input_section)
Write structure data requested by a separate structure data input section to the output unit....
subroutine, public write_particle_distances(particle_set, cell, subsys_section)
Write the matrix of the particle distances to the output unit.
Define the data structure for the particle information.
container for various plainwaves related things
subroutine, public qs_basis_rotation(qs_env, kpoints)
Construct basis set rotation matrices.
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.
Calculation of non local dispersion functionals Some routines adapted from: Copyright (C) 2001-2009 Q...
subroutine, public qs_dispersion_nonloc_init(dispersion_env, para_env)
...
Calculation of dispersion using pair potentials.
subroutine, public qs_dispersion_pairpot_init(atomic_kind_set, qs_kind_set, dispersion_env, pp_section, para_env)
...
Definition of disperson types for DFT calculations.
Set disperson types for DFT calculations.
subroutine, public qs_dispersion_env_set(dispersion_env, xc_section)
...
subroutine, public qs_write_dispersion(qs_env, dispersion_env, ounit)
...
subroutine, public allocate_qs_energy(qs_energy)
Allocate and/or initialise a Quickstep energy data structure.
qs_environment methods that use many other modules
subroutine, public qs_env_setup(qs_env)
initializes various components of the qs_env, that need only atomic_kind_set, cell,...
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, rhoz_cneo_set, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, harris_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, wanniercentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Set the QUICKSTEP environment.
subroutine, public qs_init(qs_env, para_env, root_section, globenv, cp_subsys, kpoint_env, cell, cell_ref, qmmm, qmmm_env_qm, force_env_section, subsys_section, use_motion_section, silent)
Read the input and the database files for the setup of the QUICKSTEP environment.
Definition of gCP types for DFT calculations.
Set disperson types for DFT calculations.
subroutine, public qs_gcp_env_set(gcp_env, xc_section)
...
subroutine, public qs_gcp_init(qs_env, gcp_env)
...
Types needed for a for a Harris model calculation.
subroutine, public harris_rhoin_init(rhoin, basis_type, qs_kind_set, atomic_kind_set, local_particles, nspin)
...
Harris method environment setup and handling.
subroutine, public harris_write_input(harris_env)
Print out the Harris method input section.
subroutine, public harris_env_create(qs_env, harris_env, harris_section)
Allocates and intitializes harris_env.
Calculate the interaction radii for the operator matrix calculation.
subroutine, public write_pgf_orb_radii(basis, atomic_kind_set, qs_kind_set, subsys_section)
Write the orbital basis function radii to the output unit.
subroutine, public write_paw_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the one center projector.
subroutine, public write_ppnl_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the projector functions of the Goedecker pseudopotential (GTH, non-local part) to ...
subroutine, public init_se_nlradius(se_control, atomic_kind_set, qs_kind_set, subsys_section)
...
subroutine, public write_ppl_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the exponential functions of the Goedecker pseudopotential (GTH,...
subroutine, public init_interaction_radii(qs_control, qs_kind_set)
Initialize all the atomic kind radii for a given threshold value.
subroutine, public write_core_charge_radii(atomic_kind_set, qs_kind_set, subsys_section)
Write the radii of the core charge distributions to the output unit.
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, cneo_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, 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_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_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
subroutine, public init_gapw_nlcc(qs_kind_set)
...
subroutine, public init_cneo_basis_set(qs_kind_set, qs_control)
...
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.
subroutine, public write_qs_kind_set(qs_kind_set, subsys_section)
Write an atomic kind set data set to the output unit.
subroutine, public init_gapw_basis_set(qs_kind_set, qs_control, force_env_section, modify_qs_control)
...
subroutine, public init_qs_kind_set(qs_kind_set)
Initialise an atomic kind set data set.
subroutine, public check_qs_kind_set(qs_kind_set, dft_control, subsys_section)
...
subroutine, public write_gto_basis_sets(qs_kind_set, subsys_section)
Write all the GTO basis sets of an atomic kind set to the output unit (for the printing of the unnorm...
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
subroutine, public qs_ks_env_create(ks_env)
Allocates a new instance of ks_env.
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
subroutine, public rho0_s_grid_create(pw_env, rho0_mpole)
...
subroutine, public init_rho0(local_rho_set, qs_env, gapw_control, zcore)
...
subroutine, public init_rho_atom(rho_atom_set, atomic_kind_set, qs_kind_set, dft_control, para_env)
...
Routines that work on qs_subsys_type.
subroutine, public qs_subsys_create(subsys, para_env, root_section, force_env_section, subsys_section, use_motion_section, cp_subsys, cell, cell_ref, elkind, silent)
Creates a qs_subsys. Optionally an existsing cp_subsys is used.
types that represent a quickstep subsys
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
subroutine, public qs_subsys_set(subsys, cp_subsys, local_particles, local_molecules, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, nelectron_total, nelectron_spin)
...
Storage of past states of the qs_env. Methods to interpolate (or actually normally extrapolate) the n...
subroutine, public wfi_create_for_kp(wf_history)
...
subroutine, public wfi_create(wf_history, interpolation_method_nr, extrapolation_order, has_unit_metric)
...
interpolate the wavefunctions to speed up the convergence when doing MD
subroutine, public wfi_release(wf_history)
releases a wf_history of a wavefunction (see doc/ReferenceCounting.html)
parameters that control a relativistic calculation
subroutine, public rel_c_create(rel_control)
allocates and initializes an rel control object with the default values
subroutine, public rel_c_read_parameters(rel_control, dft_section)
reads the parameters of the relativistic section into the given rel_control
parameters that control an scf iteration
subroutine, public scf_c_read_parameters(scf_control, inp_section)
reads the parameters of the scf section into the given scf_control
subroutine, public scf_c_write_parameters(scf_control, dft_section)
writes out the scf parameters
subroutine, public scf_c_create(scf_control)
allocates and initializes an scf control object with the default values
Methods for handling the 1/R^3 residual integral part.
subroutine, public semi_empirical_expns3_setup(qs_kind_set, se_control, method_id)
Setup the quantity necessary to handle the slowly convergent residual integral term 1/R^3.
Arrays of parameters used in the semi-empirical calculations \References Everywhere in this module TC...
subroutine, public init_se_intd_array()
Initialize all arrays used for the evaluation of the integrals.
Setup and Methods for semi-empirical multipole types.
subroutine, public nddo_mpole_setup(nddo_mpole, natom)
Setup NDDO multipole type.
Definition of the semi empirical multipole integral expansions types.
Type to store integrals for semi-empirical calculations.
subroutine, public semi_empirical_si_create(store_int_env, se_section, compression)
Allocate semi-empirical store integrals type.
Definition of the semi empirical parameter types.
subroutine, public se_taper_create(se_taper, integral_screening, do_ewald, taper_cou, range_cou, taper_exc, range_exc, taper_scr, range_scr, taper_lrc, range_lrc)
Creates the taper type used in SE calculations.
Working with the semi empirical parameter types.
subroutine, public se_cutoff_compatible(se_control, se_section, cell, output_unit)
Reset cutoffs trying to be somehow a bit smarter.
interface to tblite
subroutine, public tb_set_calculator(tb, typ)
...
subroutine, public tb_init_wf(tb)
initialize wavefunction ...
subroutine, public tb_init_geometry(qs_env, tb)
intialize geometry objects ...
subroutine, public tb_get_basis(tb, gto_basis_set, element_symbol, param, occ)
...
routines for DFT+NEGF calculations (coupling with the quantum transport code OMEN)
Definition transport.F:19
subroutine, public transport_env_create(qs_env)
creates the transport environment
Definition transport.F:124
Read xTB parameters.
subroutine, public xtb_parameters_set(param)
Read atom parameters for xTB Hamiltonian from input file.
subroutine, public xtb_parameters_init(param, gfn_type, element_symbol, parameter_file_path, parameter_file_name, para_env)
...
subroutine, public init_xtb_basis(param, gto_basis_set, ngauss)
...
xTB (repulsive) pair potentials Reference: Stefan Grimme, Christoph Bannwarth, Philip Shushkov JCTC 1...
subroutine, public xtb_pp_radius(qs_kind_set, ppradius, eps_pair, kfparam)
...
Definition of the xTB parameter types.
Definition xtb_types.F:20
subroutine, public allocate_xtb_atom_param(xtb_parameter)
...
Definition xtb_types.F:99
subroutine, public set_xtb_atom_param(xtb_parameter, aname, typ, defined, z, zeff, natorb, lmax, nao, lao, rcut, rcov, kx, eta, xgamma, alpha, zneff, nshell, nval, lval, kpoly, kappa, hen, zeta, xi, kappa0, alpg, electronegativity, occupation, chmax, en, kqat2, kcn, kq)
...
Definition xtb_types.F:299
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
type of a logger, at the moment it contains just a print level starting at which level it should be l...
represents a system: atoms, molecules, their pos,vel,...
structure to store local (to a processor) ordered lists of integers.
Contains information on the energy correction functional for KG.
Contains information on the excited states energy.
contains the initially parsed file and the initial parallel environment
Contains information about kpoints.
stores all the informations relevant to an mpi environment
contained for different pw related things
Contains information on the Harris method.
Provides all information about a quickstep kind.
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
keeps track of the previous wavefunctions and can extrapolate them for the next step of md
contains the parameters needed by a relativistic calculation
Global Multipolar NDDO information type.
Taper type use in semi-empirical calculations.