Using ForwardDiff to compute the jacobian and gradient of a function that includes another function

Specific Domains Numerics
1

hi guys I have a situation like following:

function update_step(kernelPCA::KernelPCA,
                    withRespectTo::Type{WithRespectToData},
                    abstractApplication::Type{Denoising},
                    X₀::AbstractVecOrMat{T},
                    tₘₐₓ::Int,
                    cost_function,
                    m::Int,
                    eigenvecsVec::AbstractVecOrMat{T};
                    ηᵢₙᵢₜ = 0.2,
                    ρ=0.5) where {T<:Real}
    Xstep = Array{Float64,2}(undef, 2, tₘₐₓ)
    Xstep[:, 1] = X₀ 
    η = ηᵢₙᵢₜ
    for t in 2:tₘₐₓ
        xTemp = Xstep[:,t-1]
        jac_dev_temp = jacobian_derivative(withRespectTo, kernelPCA.ker, kernelPCA, abstractApplication, xTemp, tₘₐₓ, m, eigenvecsVec)

        jac_dev_temp_FD = ForwardDiff.jacobian(cost_function, xTemp)
        grad_dev_temp_FD = ForwardDiff.gradient(cost_function, xTemp)
	.
	.
	.
        
        Xstep[:, t] = Xstep[:, t-1] - η .* Diagonal(Xstep[:, t-1]) * jac_dev_temp
    end
    return Xstep
end





function cost_function(kernelPCA::KernelPCA,
                        abstractApplication::Type{Denoising},
                        X₀::AbstractVecOrMat{T},
                        m::Int,
                        eigenvecsVec::AbstractVecOrMat{T}) where T <: Real
    dims = size(kernelPCA.X)
    cost_fun = 0.0
    for j in 1:dims[2]
        gammaJ = construct_gamma(Denoising, kernelPCA, kernelPCA.X[:,j], j, m, eigenvecsVec)
        cost_fun += gammaJ * kernelPCA.ker(kernelPCA.X[:,j], X₀)
    end
    cost_fun = -cost_fun
    cost_fun += (1/2) * kernelPCA.ker(X₀, X₀)
end




function construct_gamma(::Type{Denoising},
                        kernelPCA::KernelPCA,
                        kernel::SquaredExponentialKernel,
                        X₀::AbstractVecOrMat{T₃},
                        j::Int,
                        m::Int,
                        eigenvecsVec::AbstractVecOrMat{T₃}) where T₃ <: Real
    n = size(kernelPCA.X, 2) 
    tempSum = 0.0
    for k in 1:m
        for i in 1:n
            tempSum += eigenvecsVec[k, i] * eigenvecsVec[k, j] * kernel(X₀, kernelPCA.X[:, i])
        end
    end
    return tempSum
end

and as you can see within function update_step I calculate the Jacobian of cost_function by using the function jacobian_derivative which contains the already calculated formula of the Jacobian for that function. I was wandering if it is possible to use ForwardDiff to do this calculation automatically. Please bare in mind that cost_function has within its body the call for another function namely construct_gamma.

Cheers.

Ergnoor