add_ecpder Subroutine

public subroutine add_ecpder(basis, coord, denab, de)

Arguments

Type IntentOptional Attributes Name
type(basis_set), intent(in) :: basis
real(kind=real64), intent(in), contiguous :: coord(:,:)
real(kind=dp), intent(inout) :: denab(:)
real(kind=dp), intent(inout) :: de(:,:)

Calls

proc~~add_ecpder~~CallsGraph proc~add_ecpder add_ecpder interface~compute_first_derivs compute_first_derivs proc~add_ecpder->interface~compute_first_derivs interface~free_integrator free_integrator proc~add_ecpder->interface~free_integrator interface~free_result free_result proc~add_ecpder->interface~free_result interface~init_integrator init_integrator proc~add_ecpder->interface~init_integrator interface~init_integrator_instance init_integrator_instance proc~add_ecpder->interface~init_integrator_instance interface~set_ecp_basis set_ecp_basis proc~add_ecpder->interface~set_ecp_basis none~bf_to_shell basis_set%bf_to_shell proc~add_ecpder->none~bf_to_shell

Called by

proc~~add_ecpder~~CalledByGraph proc~add_ecpder add_ecpder proc~grad_1e_ecp grad_1e_ecp proc~grad_1e_ecp->proc~add_ecpder proc~hf_gradient hf_gradient proc~hf_gradient->proc~grad_1e_ecp proc~tdhf_1e_grad tdhf_1e_grad proc~tdhf_1e_grad->proc~grad_1e_ecp proc~tdhf_gradient tdhf_gradient proc~tdhf_gradient->proc~tdhf_1e_grad proc~tdhf_gradient_c tdhf_gradient_C proc~tdhf_gradient_c->proc~tdhf_gradient

Source Code

    subroutine add_ecpder(basis, coord, denab, de)

        real(real64), contiguous, intent(in) :: coord(:,:)
        type(basis_set), intent(in) :: basis
        REAL(kind=dp), INTENT(INOUT) :: denab(:)
        REAL(kind=dp), intent(INOUT) :: de(:,:)

        type(ecp_result) :: result_ptr
        type(c_ptr) :: integrator
        real(c_double), pointer :: libecp_res(:)
        real(real64), allocatable :: res_x(:), res_y(:), res_z(:)
        real(real64), allocatable :: deloc(:,:)
        integer :: i, j, c, n, natm, prim
        integer :: tri_size, full_size
        integer(c_int) :: driv_order

        if (.not.(basis%ecp_params%is_ecp)) then
            return
        end if

        driv_order = 1

        tri_size = basis%nbf * (basis%nbf + 1) / 2
        full_size = basis%nbf * basis%nbf

        allocate(res_x(full_size))
        allocate(res_y(full_size))
        allocate(res_z(full_size))

        natm = size(coord, dim=2)

        allocate(deloc(3, natm))
        deloc = 0

        call set_integrator(integrator, basis, coord, driv_order)

        result_ptr = compute_first_derivs(integrator)

        call c_f_pointer(result_ptr%data, libecp_res, [result_ptr%size])


        do n = 1, natm
            res_x = libecp_res(full_size * (3 * n - 3)+1:full_size * (3 * n - 2))
            res_y = libecp_res(full_size * (3 * n - 2)+1:full_size * (3 * n - 1))
            res_z = libecp_res(full_size * (3 * n - 1)+1:full_size * (3 * n))
            call transform_matrix(basis, res_x)
            call transform_matrix(basis, res_y)
            call transform_matrix(basis, res_z)

            do j = 1, basis%nbf
                do i = 1, j
                    c = j * (j - 1) / 2 + i

                    if (i == j) then
                        prim = 1
                    else
                        prim = 2
                    end if

                    deloc(1, n) = deloc(1, n) + prim * res_x((i - 1) * basis%nbf + j) * denab(c)
                    deloc(2, n) = deloc(2, n) + prim * res_y((i - 1) * basis%nbf + j) * denab(c)
                    deloc(3, n) = deloc(3, n) + prim * res_z((i - 1) * basis%nbf + j) * denab(c)
                end do
            end do

        end do

        de(:, 1:natm) = de(:, 1:natm) + deloc(:, 1:natm)

        result_ptr%data = c_null_ptr
        result_ptr%size = 0
        nullify(libecp_res)

        call free_integrator(integrator)
        call free_result(result_ptr)

    end subroutine add_ecpder