LCOV - code coverage report
Current view: top level - src/xc - xc_fxc_kernel.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:42dac4a) Lines: 79.8 % 119 95
Test Date: 2025-07-25 12:55:17 Functions: 100.0 % 2 2

            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 Exchange and Correlation kernel functionals
      10              : !> \author JGH
      11              : ! **************************************************************************************************
      12              : MODULE xc_fxc_kernel
      13              :    USE input_section_types,             ONLY: section_vals_type,&
      14              :                                               section_vals_val_get
      15              :    USE kinds,                           ONLY: dp
      16              :    USE pw_methods,                      ONLY: pw_axpy,&
      17              :                                               pw_copy,&
      18              :                                               pw_scale,&
      19              :                                               pw_zero
      20              :    USE pw_pool_types,                   ONLY: pw_pool_type
      21              :    USE pw_types,                        ONLY: pw_c1d_gs_type,&
      22              :                                               pw_r3d_rs_type
      23              :    USE xc_b97_fxc,                      ONLY: b97_fcc_eval,&
      24              :                                               b97_fxc_eval
      25              :    USE xc_input_constants,              ONLY: xc_deriv_pw
      26              :    USE xc_pade,                         ONLY: pade_fxc_eval,&
      27              :                                               pade_init
      28              :    USE xc_perdew_wang,                  ONLY: perdew_wang_fxc_calc
      29              :    USE xc_rho_cflags_types,             ONLY: xc_rho_cflags_setall,&
      30              :                                               xc_rho_cflags_type
      31              :    USE xc_util,                         ONLY: xc_pw_gradient
      32              :    USE xc_xalpha,                       ONLY: xalpha_fxc_eval
      33              : #include "../base/base_uses.f90"
      34              : 
      35              :    IMPLICIT NONE
      36              :    PRIVATE
      37              :    PUBLIC :: calc_fxc_kernel
      38              : 
      39              :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'xc_fxc_kernel'
      40              : 
      41              : CONTAINS
      42              : 
      43              : ! **************************************************************************************************
      44              : !> \brief Exchange and Correlation kernel functional calculations
      45              : !> \param fxc_rspace ...
      46              : !> \param rho_r the value of the density in the real space
      47              : !> \param rho_g value of the density in the g space (needs to be associated
      48              : !>        only for gradient corrections)
      49              : !> \param tau_r value of the kinetic density tau on the grid (can be null,
      50              : !>        used only with meta functionals)
      51              : !> \param xc_kernel which functional to calculate, and how to do it
      52              : !> \param triplet ...
      53              : !> \param pw_pool the pool for the grids
      54              : !> \author JGH
      55              : ! **************************************************************************************************
      56           12 :    SUBROUTINE calc_fxc_kernel(fxc_rspace, rho_r, rho_g, tau_r, xc_kernel, triplet, pw_pool)
      57              :       TYPE(pw_r3d_rs_type), DIMENSION(:)                 :: fxc_rspace, rho_r
      58              :       TYPE(pw_c1d_gs_type), DIMENSION(:)                 :: rho_g
      59              :       TYPE(pw_r3d_rs_type), DIMENSION(:)                 :: tau_r
      60              :       TYPE(section_vals_type), POINTER                   :: xc_kernel
      61              :       LOGICAL, INTENT(IN)                                :: triplet
      62              :       TYPE(pw_pool_type), POINTER                        :: pw_pool
      63              : 
      64              :       CHARACTER(len=*), PARAMETER                        :: routineN = 'calc_fxc_kernel'
      65              :       REAL(KIND=dp), PARAMETER                           :: eps_rho = 1.E-10_dp
      66              : 
      67              :       CHARACTER(len=20)                                  :: fxc_name
      68              :       INTEGER                                            :: handle, i, idir, j, k, nspins
      69              :       INTEGER, DIMENSION(2, 3)                           :: bo
      70              :       LOGICAL                                            :: lsd
      71              :       REAL(KIND=dp)                                      :: scalec, scalex
      72              :       REAL(KIND=dp), DIMENSION(3)                        :: ccaa, ccab, cxaa, g_ab
      73           12 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: rvals
      74              :       TYPE(pw_c1d_gs_type)                               :: rhog, tmpg
      75              :       TYPE(pw_r3d_rs_type)                               :: fxa, fxb, norm_drhoa, norm_drhob, rhoa, &
      76              :                                                             rhob
      77           48 :       TYPE(pw_r3d_rs_type), DIMENSION(3)                 :: drhoa
      78              :       TYPE(xc_rho_cflags_type)                           :: needs
      79              : 
      80            0 :       CPASSERT(ASSOCIATED(xc_kernel))
      81           12 :       CPASSERT(ASSOCIATED(pw_pool))
      82              : 
      83           12 :       CALL timeset(routineN, handle)
      84              : 
      85           12 :       nspins = SIZE(rho_r)
      86           12 :       lsd = (nspins == 2)
      87           12 :       IF (triplet) THEN
      88            0 :          CPASSERT(nspins == 1)
      89              :       END IF
      90              : 
      91           12 :       CALL section_vals_val_get(xc_kernel, "_SECTION_PARAMETERS_", c_val=fxc_name)
      92           12 :       CALL section_vals_val_get(xc_kernel, "SCALE_X", r_val=scalex)
      93           12 :       CALL section_vals_val_get(xc_kernel, "SCALE_C", r_val=scalec)
      94              : 
      95           12 :       CALL xc_rho_cflags_setall(needs, .FALSE.)
      96           12 :       CALL fxc_kernel_info(fxc_name, needs, lsd)
      97              : 
      98           12 :       CALL pw_pool%create_pw(rhoa)
      99           12 :       CALL pw_pool%create_pw(rhob)
     100           12 :       IF (lsd) THEN
     101            0 :          CALL pw_copy(rho_r(1), rhoa)
     102            0 :          CALL pw_copy(rho_r(2), rhob)
     103           12 :       ELSE IF (triplet) THEN
     104            0 :          CALL pw_copy(rho_r(1), rhoa)
     105            0 :          CALL pw_copy(rho_r(1), rhob)
     106              :       ELSE
     107           12 :          CALL pw_copy(rho_r(1), rhoa)
     108           12 :          CALL pw_copy(rho_r(1), rhob)
     109           12 :          CALL pw_scale(rhoa, 0.5_dp)
     110           12 :          CALL pw_scale(rhob, 0.5_dp)
     111              :       END IF
     112           12 :       IF (needs%norm_drho) THEN
     113              :          ! deriv rho
     114           16 :          DO idir = 1, 3
     115           16 :             CALL pw_pool%create_pw(drhoa(idir))
     116              :          END DO
     117            4 :          CALL pw_pool%create_pw(norm_drhoa)
     118            4 :          CALL pw_pool%create_pw(norm_drhob)
     119            4 :          CALL pw_pool%create_pw(rhog)
     120            4 :          CALL pw_pool%create_pw(tmpg)
     121            4 :          IF (lsd) THEN
     122            0 :             CALL pw_copy(rho_g(1), rhog)
     123            4 :          ELSE IF (triplet) THEN
     124            0 :             CALL pw_copy(rho_g(1), rhog)
     125              :          ELSE
     126            4 :             CALL pw_copy(rho_g(1), rhog)
     127            4 :             CALL pw_scale(rhog, 0.5_dp)
     128              :          END IF
     129            4 :          CALL xc_pw_gradient(rhoa, rhog, tmpg, drhoa(:), xc_deriv_pw)
     130           40 :          bo(1:2, 1:3) = rhoa%pw_grid%bounds_local(1:2, 1:3)
     131            4 : !$OMP    PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(bo,norm_drhoa,drhoa)
     132              :          DO k = bo(1, 3), bo(2, 3)
     133              :             DO j = bo(1, 2), bo(2, 2)
     134              :                DO i = bo(1, 1), bo(2, 1)
     135              :                   norm_drhoa%array(i, j, k) = SQRT(drhoa(1)%array(i, j, k)**2 + &
     136              :                                                    drhoa(2)%array(i, j, k)**2 + &
     137              :                                                    drhoa(3)%array(i, j, k)**2)
     138              :                END DO
     139              :             END DO
     140              :          END DO
     141            4 :          IF (lsd) THEN
     142            0 :             CALL pw_copy(rho_g(2), rhog)
     143            0 :             CALL xc_pw_gradient(rhob, rhog, tmpg, drhoa(:), xc_deriv_pw)
     144            0 :             bo(1:2, 1:3) = rhob%pw_grid%bounds_local(1:2, 1:3)
     145            0 : !$OMP       PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(bo,norm_drhob,drhoa)
     146              :             DO k = bo(1, 3), bo(2, 3)
     147              :                DO j = bo(1, 2), bo(2, 2)
     148              :                   DO i = bo(1, 1), bo(2, 1)
     149              :                      norm_drhob%array(i, j, k) = SQRT(drhoa(1)%array(i, j, k)**2 + &
     150              :                                                       drhoa(2)%array(i, j, k)**2 + &
     151              :                                                       drhoa(3)%array(i, j, k)**2)
     152              :                   END DO
     153              :                END DO
     154              :             END DO
     155              :          ELSE
     156       365733 :             norm_drhob%array(:, :, :) = norm_drhoa%array(:, :, :)
     157              :          END IF
     158            4 :          CALL pw_pool%give_back_pw(rhog)
     159            4 :          CALL pw_pool%give_back_pw(tmpg)
     160              :       END IF
     161           12 :       IF (needs%tau) THEN
     162              :          MARK_USED(tau_r)
     163            0 :          CPABORT("Meta functionals not available.")
     164              :       END IF
     165              : 
     166           14 :       SELECT CASE (TRIM(fxc_name))
     167              :       CASE ("PADEFXC")
     168            2 :          IF (scalec == scalex) THEN
     169            2 :             CALL pade_init(eps_rho)
     170            2 :             CALL pade_fxc_eval(rhoa, rhob, fxc_rspace(1), fxc_rspace(2), fxc_rspace(3))
     171            2 :             IF (scalex /= 1.0_dp) THEN
     172            0 :                CALL pw_scale(fxc_rspace(1), scalex)
     173            0 :                CALL pw_scale(fxc_rspace(2), scalex)
     174            0 :                CALL pw_scale(fxc_rspace(3), scalex)
     175              :             END IF
     176              :          ELSE
     177            0 :             CPABORT("PADE Fxc Kernel functional needs SCALE_X==SCALE_C")
     178              :          END IF
     179              :       CASE ("LDAFXC")
     180            6 :          CALL pw_zero(fxc_rspace(1))
     181            6 :          CALL pw_zero(fxc_rspace(2))
     182            6 :          CALL pw_zero(fxc_rspace(3))
     183            6 :          CALL xalpha_fxc_eval(rhoa, rhob, fxc_rspace(1), fxc_rspace(3), scalex, eps_rho)
     184              :          CALL perdew_wang_fxc_calc(rhoa, rhob, fxc_rspace(1), fxc_rspace(2), fxc_rspace(3), &
     185            6 :                                    scalec, eps_rho)
     186              :       CASE ("GGAFXC")
     187              :          ! get parameter
     188            4 :          CALL section_vals_val_get(xc_kernel, "GAMMA", r_vals=rvals)
     189           16 :          g_ab(1:3) = rvals(1:3)
     190            4 :          CALL section_vals_val_get(xc_kernel, "C_XAA", r_vals=rvals)
     191           16 :          cxaa(1:3) = rvals(1:3)
     192            4 :          CALL section_vals_val_get(xc_kernel, "C_CAA", r_vals=rvals)
     193           16 :          ccaa(1:3) = rvals(1:3)
     194            4 :          CALL section_vals_val_get(xc_kernel, "C_CAB", r_vals=rvals)
     195           16 :          ccab(1:3) = rvals(1:3)
     196              :          ! correlation
     197            4 :          CALL pw_zero(fxc_rspace(1))
     198            4 :          CALL pw_zero(fxc_rspace(2))
     199            4 :          CALL pw_zero(fxc_rspace(3))
     200              :          CALL perdew_wang_fxc_calc(rhoa, rhob, fxc_rspace(1), fxc_rspace(2), fxc_rspace(3), &
     201            4 :                                    scalec, eps_rho)
     202            4 :          CALL b97_fxc_eval(rhoa, norm_drhoa, fxc_rspace(1), g_ab(1), ccaa, eps_rho)
     203            4 :          CALL b97_fxc_eval(rhob, norm_drhob, fxc_rspace(3), g_ab(3), ccaa, eps_rho)
     204            4 :          CALL b97_fcc_eval(rhoa, rhob, norm_drhoa, norm_drhob, fxc_rspace(2), g_ab(2), ccab, eps_rho)
     205              :          ! exchange
     206            4 :          CALL pw_pool%create_pw(fxa)
     207            4 :          CALL pw_pool%create_pw(fxb)
     208            4 :          CALL pw_zero(fxa)
     209            4 :          CALL pw_zero(fxb)
     210            4 :          CALL xalpha_fxc_eval(rhoa, rhob, fxa, fxb, scalex, eps_rho)
     211            4 :          CALL b97_fxc_eval(rhoa, norm_drhoa, fxa, g_ab(1), cxaa, eps_rho)
     212            4 :          CALL b97_fxc_eval(rhob, norm_drhob, fxb, g_ab(1), cxaa, eps_rho)
     213            4 :          CALL pw_axpy(fxa, fxc_rspace(1))
     214            4 :          CALL pw_axpy(fxb, fxc_rspace(3))
     215            4 :          CALL pw_pool%give_back_pw(fxa)
     216            4 :          CALL pw_pool%give_back_pw(fxb)
     217              :       CASE ("NONE")
     218            0 :          CALL pw_zero(fxc_rspace(1))
     219            0 :          CALL pw_zero(fxc_rspace(2))
     220            0 :          CALL pw_zero(fxc_rspace(3))
     221              :       CASE default
     222           12 :          CPABORT("Fxc Kernel functional is defined incorrectly")
     223              :       END SELECT
     224              : 
     225           12 :       CALL pw_pool%give_back_pw(rhoa)
     226           12 :       CALL pw_pool%give_back_pw(rhob)
     227           12 :       IF (needs%norm_drho) THEN
     228            4 :          CALL pw_pool%give_back_pw(norm_drhoa)
     229            4 :          CALL pw_pool%give_back_pw(norm_drhob)
     230           16 :          DO idir = 1, 3
     231           16 :             CALL pw_pool%give_back_pw(drhoa(idir))
     232              :          END DO
     233              :       END IF
     234              : 
     235           12 :       CALL timestop(handle)
     236              : 
     237           12 :    END SUBROUTINE calc_fxc_kernel
     238              : 
     239              : ! **************************************************************************************************
     240              : !> \brief ...
     241              : !> \param fxc_name ...
     242              : !> \param needs ...
     243              : !> \param lsd ...
     244              : ! **************************************************************************************************
     245           12 :    SUBROUTINE fxc_kernel_info(fxc_name, needs, lsd)
     246              :       CHARACTER(len=20), INTENT(IN)                      :: fxc_name
     247              :       TYPE(xc_rho_cflags_type), INTENT(INOUT)            :: needs
     248              :       LOGICAL, INTENT(IN)                                :: lsd
     249              : 
     250           20 :       SELECT CASE (TRIM(fxc_name))
     251              :       CASE ("PADEFXC", "LDAFXC")
     252            8 :          IF (lsd) THEN
     253            0 :             needs%rho_spin = .TRUE.
     254              :          ELSE
     255            8 :             needs%rho = .TRUE.
     256              :          END IF
     257              :       CASE ("GGAFXC")
     258            4 :          IF (lsd) THEN
     259            0 :             needs%rho_spin = .TRUE.
     260            0 :             needs%norm_drho_spin = .TRUE.
     261            0 :             needs%norm_drho = .TRUE.
     262              :          ELSE
     263            4 :             needs%rho = .TRUE.
     264            4 :             needs%norm_drho = .TRUE.
     265              :          END IF
     266              :       CASE ("NONE")
     267              :       CASE default
     268           12 :          CPABORT("Fxc Kernel functional is defined incorrectly")
     269              :       END SELECT
     270              : 
     271           12 :    END SUBROUTINE fxc_kernel_info
     272              : 
     273              : END MODULE xc_fxc_kernel
     274              : 
        

Generated by: LCOV version 2.0-1