Line data    Source code

       1              : !--------------------------------------------------------------------------------------------------!
       2              : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3              : !   Copyright 2000-2026 CP2K developers group <https://cp2k.org>                                   !
       4              : !                                                                                                  !
       5              : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6              : !--------------------------------------------------------------------------------------------------!
       7              : 
       8            2 : PROGRAM qs_fod_unittest
       9            2 :    USE kinds,                           ONLY: dp
      10              :    USE qs_fod,                          ONLY: fod_weights
      11              : 
      12              :    IMPLICIT NONE
      13              :    REAL(dp), DIMENSION(4) :: occ, eig, weight, expected
      14              : 
      15            2 :    eig = [-1.0_dp, -0.1_dp, 0.0_dp, 0.2_dp]
      16            2 :    occ = [2.0_dp, 1.6_dp, 1.0_dp, 0.2_dp]
      17            2 :    expected = [0.0_dp, 0.4_dp, 1.0_dp, 0.2_dp]
      18            2 :    weight = fod_weights(occ, eig, 0.0_dp, 2.0_dp)
      19           10 :    IF (MAXVAL(ABS(weight - expected)) > 1.0E-14_dp) ERROR STOP "Restricted FOD weights"
      20           10 :    weight = fod_weights(0.5_dp*occ, eig, 0.0_dp, 1.0_dp)
      21           10 :    IF (MAXVAL(ABS(2.0_dp*weight - expected)) > 1.0E-14_dp) ERROR STOP "Unrestricted FOD weights"
      22              : 
      23              :    ! Fractional spin populations: the electron-count/index split is not valid.
      24            2 :    eig = [0.1_dp, 0.2_dp, 0.3_dp, 0.4_dp]
      25            2 :    occ = [0.4_dp, 0.3_dp, 0.2_dp, 0.1_dp]
      26            2 :    weight = fod_weights(occ, eig, 0.0_dp, 1.0_dp)
      27           10 :    IF (MAXVAL(ABS(weight - occ)) > 1.0E-14_dp) ERROR STOP "FOD must use the chemical potential"
      28            2 :    weight = fod_weights(occ(4:1:-1), eig(4:1:-1), 0.0_dp, 1.0_dp)
      29           10 :    IF (MAXVAL(ABS(weight - occ(4:1:-1))) > 1.0E-14_dp) ERROR STOP "FOD orbital permutation"
      30              : 
      31              :    ! Integer occupations vanish; a shift of the energy zero changes nothing.
      32            2 :    eig = [-2.0_dp, -1.0_dp, 1.0_dp, 2.0_dp]
      33            2 :    occ = [2.0_dp, 2.0_dp, 0.0_dp, 0.0_dp]
      34           10 :    weight = fod_weights(occ, eig + 3.0_dp, 3.0_dp, 2.0_dp)
      35           10 :    IF (ANY(weight /= 0.0_dp)) ERROR STOP "Integer occupations should have no FOD"
      36            2 : END PROGRAM qs_fod_unittest