LCOV - code coverage report
Current view: top level - src - qs_band_structure.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:42dac4a) Lines: 92.4 % 184 170
Test Date: 2025-07-25 12:55:17 Functions: 100.0 % 4 4

            Line data    Source code
       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              : !> \brief Calculation of band structures
      10              : !> \par History
      11              : !>       2015.06 created [JGH]
      12              : !> \author JGH
      13              : ! **************************************************************************************************
      14              : MODULE qs_band_structure
      15              :    USE cell_types,                      ONLY: cell_type
      16              :    USE cp_blacs_env,                    ONLY: cp_blacs_env_type
      17              :    USE cp_control_types,                ONLY: dft_control_type
      18              :    USE cp_dbcsr_api,                    ONLY: dbcsr_p_type
      19              :    USE cp_files,                        ONLY: close_file,&
      20              :                                               open_file
      21              :    USE cp_log_handling,                 ONLY: cp_logger_get_default_io_unit
      22              :    USE cp_parser_methods,               ONLY: read_float_object
      23              :    USE input_section_types,             ONLY: section_vals_get,&
      24              :                                               section_vals_get_subs_vals,&
      25              :                                               section_vals_type,&
      26              :                                               section_vals_val_get
      27              :    USE kinds,                           ONLY: default_string_length,&
      28              :                                               dp,&
      29              :                                               max_line_length
      30              :    USE kpoint_methods,                  ONLY: kpoint_env_initialize,&
      31              :                                               kpoint_init_cell_index,&
      32              :                                               kpoint_initialize_mo_set,&
      33              :                                               kpoint_initialize_mos
      34              :    USE kpoint_types,                    ONLY: get_kpoint_info,&
      35              :                                               kpoint_create,&
      36              :                                               kpoint_env_type,&
      37              :                                               kpoint_release,&
      38              :                                               kpoint_sym_create,&
      39              :                                               kpoint_type
      40              :    USE machine,                         ONLY: m_walltime
      41              :    USE mathconstants,                   ONLY: twopi
      42              :    USE message_passing,                 ONLY: mp_para_env_type
      43              :    USE physcon,                         ONLY: angstrom,&
      44              :                                               evolt
      45              :    USE qs_environment_types,            ONLY: get_qs_env,&
      46              :                                               qs_env_release,&
      47              :                                               qs_environment_type
      48              :    USE qs_gamma2kp,                     ONLY: create_kp_from_gamma
      49              :    USE qs_mo_types,                     ONLY: get_mo_set
      50              :    USE qs_neighbor_list_types,          ONLY: neighbor_list_set_p_type
      51              :    USE qs_scf_diagonalization,          ONLY: do_general_diag_kp
      52              :    USE qs_scf_types,                    ONLY: qs_scf_env_type
      53              :    USE scf_control_types,               ONLY: scf_control_type
      54              :    USE string_utilities,                ONLY: uppercase
      55              : #include "./base/base_uses.f90"
      56              : 
      57              :    IMPLICIT NONE
      58              : 
      59              :    PRIVATE
      60              : 
      61              :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_band_structure'
      62              : 
      63              :    PUBLIC :: calculate_band_structure, calculate_kp_orbitals, calculate_kpoints_for_bs
      64              : 
      65              : ! **************************************************************************************************
      66              : 
      67              : CONTAINS
      68              : 
      69              : ! **************************************************************************************************
      70              : !> \brief Main routine for band structure calculation
      71              : !> \param qs_env ...
      72              : !> \author JGH
      73              : ! **************************************************************************************************
      74        41442 :    SUBROUTINE calculate_band_structure(qs_env)
      75              :       TYPE(qs_environment_type), POINTER                 :: qs_env
      76              : 
      77              :       LOGICAL                                            :: do_kpoints, explicit
      78              :       TYPE(section_vals_type), POINTER                   :: bs_input
      79              : 
      80        20721 :       bs_input => section_vals_get_subs_vals(qs_env%input, "DFT%PRINT%BAND_STRUCTURE")
      81        20721 :       CALL section_vals_get(bs_input, explicit=explicit)
      82        20721 :       IF (explicit) THEN
      83            8 :          CALL get_qs_env(qs_env, do_kpoints=do_kpoints)
      84            8 :          IF (do_kpoints) THEN
      85            6 :             CALL do_calculate_band_structure(qs_env)
      86              :          ELSE
      87              :             BLOCK
      88              :                TYPE(qs_environment_type), POINTER :: qs_env_kp
      89            2 :                CALL create_kp_from_gamma(qs_env, qs_env_kp)
      90            2 :                CALL do_calculate_band_structure(qs_env_kp)
      91            2 :                CALL qs_env_release(qs_env_kp)
      92            2 :                DEALLOCATE (qs_env_kp)
      93              :             END BLOCK
      94              :          END IF
      95              :       END IF
      96              : 
      97        20721 :    END SUBROUTINE calculate_band_structure
      98              : 
      99              : ! **************************************************************************************************
     100              : !> \brief band structure calculation
     101              : !> \param qs_env ...
     102              : !> \author JGH
     103              : ! **************************************************************************************************
     104           16 :    SUBROUTINE do_calculate_band_structure(qs_env)
     105              :       TYPE(qs_environment_type), POINTER                 :: qs_env
     106              : 
     107              :       CHARACTER(LEN=default_string_length)               :: filename, ustr
     108              :       CHARACTER(LEN=default_string_length), &
     109            8 :          DIMENSION(:), POINTER                           :: spname, strptr
     110              :       CHARACTER(LEN=max_line_length)                     :: error_message
     111              :       INTEGER                                            :: bs_data_unit, i, i_rep, ik, ikk, ikpgr, &
     112              :                                                             imo, ip, ispin, n_ptr, n_rep, nadd, &
     113              :                                                             nkp, nmo, npline, npoints, nspins, &
     114              :                                                             unit_nr
     115              :       INTEGER, DIMENSION(2)                              :: kp_range
     116              :       LOGICAL                                            :: explicit, io_default, my_kpgrp
     117              :       REAL(KIND=dp)                                      :: t1, t2
     118              :       REAL(KIND=dp), DIMENSION(3)                        :: kpptr
     119            8 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: eigenvalues, eigval, occnum, &
     120            8 :                                                             occupation_numbers, wkp
     121            8 :       REAL(kind=dp), DIMENSION(:, :), POINTER            :: kpgeneral, kspecial, xkp
     122              :       TYPE(cell_type), POINTER                           :: cell
     123              :       TYPE(dft_control_type), POINTER                    :: dft_control
     124              :       TYPE(kpoint_env_type), POINTER                     :: kp
     125              :       TYPE(kpoint_type), POINTER                         :: kpoint
     126              :       TYPE(mp_para_env_type), POINTER                    :: para_env
     127              :       TYPE(section_vals_type), POINTER                   :: bs_input, kpset
     128              : 
     129           16 :       bs_input => section_vals_get_subs_vals(qs_env%input, "DFT%PRINT%BAND_STRUCTURE")
     130            8 :       CALL section_vals_get(bs_input, explicit=explicit)
     131            8 :       IF (explicit) THEN
     132            8 :          CALL section_vals_val_get(bs_input, "FILE_NAME", c_val=filename)
     133            8 :          CALL section_vals_val_get(bs_input, "ADDED_MOS", i_val=nadd)
     134            8 :          unit_nr = cp_logger_get_default_io_unit()
     135            8 :          CALL get_qs_env(qs_env=qs_env, para_env=para_env)
     136            8 :          CALL get_qs_env(qs_env, cell=cell)
     137            8 :          kpset => section_vals_get_subs_vals(bs_input, "KPOINT_SET")
     138            8 :          CALL section_vals_get(kpset, n_repetition=n_rep)
     139            8 :          IF (unit_nr > 0) THEN
     140            4 :             WRITE (unit_nr, FMT="(/,T2,A)") "KPOINTS| Band Structure Calculation"
     141            4 :             WRITE (unit_nr, FMT="(T2,A,T71,I10)") "KPOINTS| Number of k-point sets", n_rep
     142            4 :             IF (nadd /= 0) THEN
     143            1 :                WRITE (unit_nr, FMT="(T2,A,T71,I10)") "KPOINTS| Number of added MOs/bands", nadd
     144              :             END IF
     145              :          END IF
     146            8 :          IF (filename == "") THEN
     147              :             ! use standard output file
     148            2 :             bs_data_unit = unit_nr
     149            2 :             io_default = .TRUE.
     150              :          ELSE
     151            6 :             io_default = .FALSE.
     152            6 :             IF (para_env%is_source()) THEN
     153              :                CALL open_file(filename, unit_number=bs_data_unit, file_status="UNKNOWN", file_action="WRITE", &
     154            3 :                               file_position="APPEND")
     155              :             ELSE
     156            3 :                bs_data_unit = -1
     157              :             END IF
     158              :          END IF
     159           18 :          DO i_rep = 1, n_rep
     160           10 :             t1 = m_walltime()
     161           10 :             CALL section_vals_val_get(kpset, "NPOINTS", i_rep_section=i_rep, i_val=npline)
     162           10 :             CALL section_vals_val_get(kpset, "UNITS", i_rep_section=i_rep, c_val=ustr)
     163           10 :             CALL uppercase(ustr)
     164           10 :             CALL section_vals_val_get(kpset, "SPECIAL_POINT", i_rep_section=i_rep, n_rep_val=n_ptr)
     165           30 :             ALLOCATE (kspecial(3, n_ptr))
     166           30 :             ALLOCATE (spname(n_ptr))
     167           30 :             DO ip = 1, n_ptr
     168           20 :                CALL section_vals_val_get(kpset, "SPECIAL_POINT", i_rep_section=i_rep, i_rep_val=ip, c_vals=strptr)
     169           20 :                IF (SIZE(strptr(:), 1) == 4) THEN
     170            6 :                   spname(ip) = strptr(1)
     171           24 :                   DO i = 1, 3
     172           18 :                      CALL read_float_object(strptr(i + 1), kpptr(i), error_message)
     173           24 :                      IF (LEN_TRIM(error_message) > 0) CPABORT(TRIM(error_message))
     174              :                   END DO
     175           14 :                ELSE IF (SIZE(strptr(:), 1) == 3) THEN
     176           14 :                   spname(ip) = "not specified"
     177           56 :                   DO i = 1, 3
     178           42 :                      CALL read_float_object(strptr(i), kpptr(i), error_message)
     179           56 :                      IF (LEN_TRIM(error_message) > 0) CPABORT(TRIM(error_message))
     180              :                   END DO
     181              :                ELSE
     182            0 :                   CPABORT("Input SPECIAL_POINT invalid")
     183              :                END IF
     184           10 :                SELECT CASE (ustr)
     185              :                CASE ("B_VECTOR")
     186           48 :                   kspecial(1:3, ip) = kpptr(1:3)
     187              :                CASE ("CART_ANGSTROM")
     188              :                   kspecial(1:3, ip) = (kpptr(1)*cell%hmat(1, 1:3) + &
     189              :                                        kpptr(2)*cell%hmat(2, 1:3) + &
     190            0 :                                        kpptr(3)*cell%hmat(3, 1:3))/twopi*angstrom
     191              :                CASE ("CART_BOHR")
     192              :                   kspecial(1:3, ip) = (kpptr(1)*cell%hmat(1, 1:3) + &
     193              :                                        kpptr(2)*cell%hmat(2, 1:3) + &
     194           56 :                                        kpptr(3)*cell%hmat(3, 1:3))/twopi
     195              :                CASE DEFAULT
     196           20 :                   CPABORT("Unknown unit <"//TRIM(ustr)//"> specified for k-point definition")
     197              :                END SELECT
     198              :             END DO
     199           10 :             npoints = (n_ptr - 1)*npline + 1
     200           10 :             CPASSERT(npoints >= 1)
     201              : 
     202              :             ! Initialize environment and calculate MOs
     203           30 :             ALLOCATE (kpgeneral(3, npoints))
     204           70 :             kpgeneral(1:3, 1) = kspecial(1:3, 1)
     205           10 :             ikk = 1
     206           20 :             DO ik = 2, n_ptr
     207          100 :                DO ip = 1, npline
     208           80 :                   ikk = ikk + 1
     209              :                   kpgeneral(1:3, ikk) = kspecial(1:3, ik - 1) + &
     210              :                                         REAL(ip, KIND=dp)/REAL(npline, KIND=dp)* &
     211          570 :                                         (kspecial(1:3, ik) - kspecial(1:3, ik - 1))
     212              :                END DO
     213              :             END DO
     214           10 :             NULLIFY (kpoint)
     215           10 :             CALL calculate_kp_orbitals(qs_env, kpoint, "GENERAL", nadd, kpgeneral=kpgeneral)
     216           10 :             DEALLOCATE (kpgeneral)
     217              : 
     218           10 :             CALL get_qs_env(qs_env, dft_control=dft_control)
     219           10 :             nspins = dft_control%nspins
     220           10 :             kp => kpoint%kp_env(1)%kpoint_env
     221           10 :             CALL get_mo_set(kp%mos(1, 1), nmo=nmo)
     222           40 :             ALLOCATE (eigval(nmo), occnum(nmo))
     223           10 :             CALL get_kpoint_info(kpoint, nkp=nkp, kp_range=kp_range, xkp=xkp, wkp=wkp)
     224              : 
     225           10 :             IF (unit_nr > 0) THEN
     226              :                WRITE (UNIT=unit_nr, FMT="(T2,A,I4,T71,I10)") &
     227            5 :                   "KPOINTS| Number of k-points in set ", i_rep, npoints
     228              :                WRITE (UNIT=unit_nr, FMT="(T2,A)") &
     229            5 :                   "KPOINTS| In units of b-vector [2pi/Bohr]"
     230           15 :                DO ip = 1, n_ptr
     231              :                   WRITE (UNIT=unit_nr, FMT="(T2,A,I5,1X,A11,3(1X,F12.6))") &
     232           15 :                      "KPOINTS| Special point ", ip, ADJUSTL(TRIM(spname(ip))), kspecial(1:3, ip)
     233              :                END DO
     234              :             END IF
     235           10 :             IF (bs_data_unit > 0 .AND. (bs_data_unit /= unit_nr)) THEN
     236              :                WRITE (UNIT=bs_data_unit, FMT="(4(A,I0),A)") &
     237            4 :                   "# Set ", i_rep, ": ", n_ptr, " special points, ", npoints, " k-points, ", nmo, " bands"
     238           13 :                DO ip = 1, n_ptr
     239              :                   WRITE (UNIT=bs_data_unit, FMT="(A,I0,T20,T24,3(1X,F14.8),2X,A)") &
     240           13 :                      "#  Special point ", ip, kspecial(1:3, ip), ADJUSTL(TRIM(spname(ip)))
     241              :                END DO
     242              :             END IF
     243              : 
     244          100 :             DO ik = 1, nkp
     245           90 :                my_kpgrp = (ik >= kp_range(1) .AND. ik <= kp_range(2))
     246          190 :                DO ispin = 1, nspins
     247           90 :                   IF (my_kpgrp) THEN
     248           90 :                      ikpgr = ik - kp_range(1) + 1
     249           90 :                      kp => kpoint%kp_env(ikpgr)%kpoint_env
     250           90 :                      CALL get_mo_set(kp%mos(1, ispin), eigenvalues=eigenvalues, occupation_numbers=occupation_numbers)
     251         3370 :                      eigval(1:nmo) = eigenvalues(1:nmo)
     252         3370 :                      occnum(1:nmo) = occupation_numbers(1:nmo)
     253              :                   ELSE
     254            0 :                      eigval(1:nmo) = 0.0_dp
     255            0 :                      occnum(1:nmo) = 0.0_dp
     256              :                   END IF
     257         3370 :                   CALL kpoint%para_env_inter_kp%sum(eigval)
     258         3370 :                   CALL kpoint%para_env_inter_kp%sum(occnum)
     259          180 :                   IF (bs_data_unit > 0) THEN
     260              :                      WRITE (UNIT=bs_data_unit, FMT="(A,I0,T15,A,I0,A,T24,3(1X,F14.8),3X,F14.8)") &
     261           45 :                         "#  Point ", ik, "  Spin ", ispin, ":", xkp(1:3, ik), wkp(ik)
     262              :                      WRITE (UNIT=bs_data_unit, FMT="(A)") &
     263           45 :                         "#   Band    Energy [eV]     Occupation"
     264          865 :                      DO imo = 1, nmo
     265              :                         WRITE (UNIT=bs_data_unit, FMT="(T2,I7,2(1X,F14.8))") &
     266          865 :                            imo, eigval(imo)*evolt, occnum(imo)
     267              :                      END DO
     268              :                   END IF
     269              :                END DO
     270              :             END DO
     271              : 
     272           10 :             DEALLOCATE (kspecial, spname)
     273           10 :             DEALLOCATE (eigval, occnum)
     274           10 :             CALL kpoint_release(kpoint)
     275           10 :             t2 = m_walltime()
     276           48 :             IF (unit_nr > 0) THEN
     277            5 :                WRITE (UNIT=unit_nr, FMT="(T2,A,T67,F14.3)") "KPOINTS| Time for k-point line ", t2 - t1
     278              :             END IF
     279              : 
     280              :          END DO
     281              : 
     282              :          ! Close output files
     283            8 :          IF (.NOT. io_default) THEN
     284            6 :             IF (para_env%is_source()) CALL close_file(bs_data_unit)
     285              :          END IF
     286              : 
     287              :       END IF
     288              : 
     289            8 :    END SUBROUTINE do_calculate_band_structure
     290              : 
     291              : ! **************************************************************************************************
     292              : !> \brief diagonalize KS matrices at a set of kpoints
     293              : !> \param qs_env ...
     294              : !> \param kpoint ...
     295              : !> \param scheme ...
     296              : !> \param nadd ...
     297              : !> \param mp_grid ...
     298              : !> \param kpgeneral ...
     299              : !> \param group_size_ext ...
     300              : !> \author JGH
     301              : ! **************************************************************************************************
     302           10 :    SUBROUTINE calculate_kp_orbitals(qs_env, kpoint, scheme, nadd, mp_grid, kpgeneral, group_size_ext)
     303              :       TYPE(qs_environment_type), POINTER                 :: qs_env
     304              :       TYPE(kpoint_type), POINTER                         :: kpoint
     305              :       CHARACTER(LEN=*), INTENT(IN)                       :: scheme
     306              :       INTEGER, INTENT(IN)                                :: nadd
     307              :       INTEGER, DIMENSION(3), INTENT(IN), OPTIONAL        :: mp_grid
     308              :       REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), &
     309              :          OPTIONAL                                        :: kpgeneral
     310              :       INTEGER, INTENT(IN), OPTIONAL                      :: group_size_ext
     311              : 
     312              :       LOGICAL                                            :: diis_step
     313              :       TYPE(cp_blacs_env_type), POINTER                   :: blacs_env
     314           10 :       TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER       :: matrix_ks, matrix_s
     315              :       TYPE(dft_control_type), POINTER                    :: dft_control
     316              :       TYPE(mp_para_env_type), POINTER                    :: para_env
     317              :       TYPE(neighbor_list_set_p_type), DIMENSION(:), &
     318           10 :          POINTER                                         :: sab_nl
     319              :       TYPE(qs_scf_env_type), POINTER                     :: scf_env
     320              :       TYPE(scf_control_type), POINTER                    :: scf_control
     321              : 
     322           10 :       CALL calculate_kpoints_for_bs(kpoint, scheme, group_size_ext, mp_grid, kpgeneral)
     323              : 
     324           10 :       CALL get_qs_env(qs_env=qs_env, para_env=para_env, blacs_env=blacs_env)
     325           10 :       CALL kpoint_env_initialize(kpoint, para_env, blacs_env)
     326              : 
     327           10 :       CALL kpoint_initialize_mos(kpoint, qs_env%mos, nadd)
     328           10 :       CALL kpoint_initialize_mo_set(kpoint)
     329              : 
     330           10 :       CALL get_qs_env(qs_env, sab_kp=sab_nl, dft_control=dft_control)
     331           10 :       CALL kpoint_init_cell_index(kpoint, sab_nl, para_env, dft_control)
     332              : 
     333              :       CALL get_qs_env(qs_env, matrix_ks_kp=matrix_ks, matrix_s_kp=matrix_s, &
     334           10 :                       scf_env=scf_env, scf_control=scf_control)
     335           10 :       CALL do_general_diag_kp(matrix_ks, matrix_s, kpoint, scf_env, scf_control, .FALSE., diis_step)
     336              : 
     337           10 :    END SUBROUTINE calculate_kp_orbitals
     338              : 
     339              : ! **************************************************************************************************
     340              : !> \brief ...
     341              : !> \param kpoint ...
     342              : !> \param scheme ...
     343              : !> \param group_size_ext ...
     344              : !> \param mp_grid ...
     345              : !> \param kpgeneral ...
     346              : ! **************************************************************************************************
     347           42 :    SUBROUTINE calculate_kpoints_for_bs(kpoint, scheme, group_size_ext, mp_grid, kpgeneral)
     348              : 
     349              :       TYPE(kpoint_type), POINTER                         :: kpoint
     350              :       CHARACTER(LEN=*), INTENT(IN)                       :: scheme
     351              :       INTEGER, INTENT(IN), OPTIONAL                      :: group_size_ext
     352              :       INTEGER, DIMENSION(3), INTENT(IN), OPTIONAL        :: mp_grid
     353              :       REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), &
     354              :          OPTIONAL                                        :: kpgeneral
     355              : 
     356              :       INTEGER                                            :: i, ix, iy, iz, npoints
     357              : 
     358           58 :       CPASSERT(.NOT. ASSOCIATED(kpoint))
     359              : 
     360           58 :       CALL kpoint_create(kpoint)
     361              : 
     362           58 :       kpoint%kp_scheme = scheme
     363           58 :       kpoint%symmetry = .FALSE.
     364           58 :       kpoint%verbose = .FALSE.
     365           58 :       kpoint%full_grid = .FALSE.
     366           58 :       kpoint%use_real_wfn = .FALSE.
     367           58 :       kpoint%eps_geo = 1.e-6_dp
     368           58 :       IF (PRESENT(group_size_ext)) THEN
     369           48 :          kpoint%parallel_group_size = group_size_ext
     370              :       ELSE
     371           10 :          kpoint%parallel_group_size = -1
     372              :       END IF
     373            0 :       SELECT CASE (scheme)
     374              :       CASE ("GAMMA")
     375            0 :          kpoint%nkp = 1
     376            0 :          ALLOCATE (kpoint%xkp(3, 1), kpoint%wkp(1))
     377            0 :          kpoint%xkp(1:3, 1) = 0.0_dp
     378            0 :          kpoint%wkp(1) = 1.0_dp
     379            0 :          kpoint%symmetry = .TRUE.
     380            0 :          ALLOCATE (kpoint%kp_sym(1))
     381            0 :          NULLIFY (kpoint%kp_sym(1)%kpoint_sym)
     382            0 :          CALL kpoint_sym_create(kpoint%kp_sym(1)%kpoint_sym)
     383              :       CASE ("MONKHORST-PACK")
     384           16 :          CPASSERT(PRESENT(mp_grid))
     385           16 :          npoints = mp_grid(1)*mp_grid(2)*mp_grid(3)
     386           64 :          kpoint%nkp_grid(1:3) = mp_grid(1:3)
     387           16 :          kpoint%full_grid = .TRUE.
     388           16 :          kpoint%nkp = npoints
     389           80 :          ALLOCATE (kpoint%xkp(3, npoints), kpoint%wkp(npoints))
     390           32 :          kpoint%wkp(:) = 1._dp/REAL(npoints, KIND=dp)
     391              :          i = 0
     392           32 :          DO ix = 1, mp_grid(1)
     393           48 :             DO iy = 1, mp_grid(2)
     394           48 :                DO iz = 1, mp_grid(3)
     395           16 :                   i = i + 1
     396           16 :                   kpoint%xkp(1, i) = REAL(2*ix - mp_grid(1) - 1, KIND=dp)/(2._dp*REAL(mp_grid(1), KIND=dp))
     397           16 :                   kpoint%xkp(2, i) = REAL(2*iy - mp_grid(2) - 1, KIND=dp)/(2._dp*REAL(mp_grid(2), KIND=dp))
     398           32 :                   kpoint%xkp(3, i) = REAL(2*iz - mp_grid(3) - 1, KIND=dp)/(2._dp*REAL(mp_grid(3), KIND=dp))
     399              :                END DO
     400              :             END DO
     401              :          END DO
     402              :          ! default: no symmetry settings
     403           64 :          ALLOCATE (kpoint%kp_sym(kpoint%nkp))
     404           32 :          DO i = 1, kpoint%nkp
     405           16 :             NULLIFY (kpoint%kp_sym(i)%kpoint_sym)
     406           32 :             CALL kpoint_sym_create(kpoint%kp_sym(i)%kpoint_sym)
     407              :          END DO
     408              :       CASE ("MACDONALD")
     409            0 :          CPABORT("MACDONALD not implemented")
     410              :       CASE ("GENERAL")
     411           42 :          CPASSERT(PRESENT(kpgeneral))
     412           42 :          npoints = SIZE(kpgeneral, 2)
     413           42 :          kpoint%nkp = npoints
     414          210 :          ALLOCATE (kpoint%xkp(3, npoints), kpoint%wkp(npoints))
     415          372 :          kpoint%wkp(:) = 1._dp/REAL(npoints, KIND=dp)
     416         1362 :          kpoint%xkp(1:3, 1:npoints) = kpgeneral(1:3, 1:npoints)
     417              :          ! default: no symmetry settings
     418          456 :          ALLOCATE (kpoint%kp_sym(kpoint%nkp))
     419          372 :          DO i = 1, kpoint%nkp
     420          330 :             NULLIFY (kpoint%kp_sym(i)%kpoint_sym)
     421          372 :             CALL kpoint_sym_create(kpoint%kp_sym(i)%kpoint_sym)
     422              :          END DO
     423              :       CASE DEFAULT
     424           58 :          CPABORT("Unknown kpoint scheme requested")
     425              :       END SELECT
     426              : 
     427           58 :    END SUBROUTINE calculate_kpoints_for_bs
     428              : 
     429              : END MODULE qs_band_structure
        

Generated by: LCOV version 2.0-1