How to construct a matrix of gradients depending on nodal locations in gridap

Hi,
I want to calculate the gradients of residuals(displacement) with respect to another vector. based on mathematical calculations, as residuals are a vector this gradient should be a matrix. This matrix should has value when we are calculating the gradient of residuals wrt variable at the same node and it should be zero in other nodes. I couldnt find a way to directly get the node number and compare with other nodes to make others zero. so i wrote a function to logically do this as below

function vector_to_matrixut(vecu)
    mat = zeros(sr, ur)
    for i in 1:sr
        mat[i, (i-1)*2+1] = vecu[(i-1)*2+1]
        mat[i, 2*i] = vecu[2*i]
    end
    return mat
end

the input parameter is the global vector of the gradients. as below:

function dσfun(ε_in,s_in)
    σ_elas = C_mat⊙ε_in
    if tr(ε_in) >= 0
        dσ = (2 * s_in)*σ_elas
    elseif  tr(ε_in) < 0
        dσ = (2 * s_in)*P_dev ⊙ σ_elas
    end
    return  dσ
end
function  drus(fem_params,pth, uh_in,sh_in, vApp)
    uApp1(x) = VectorValue(0.0,0.0)
    uApp2(x) = VectorValue(0.0,0.0)
    uApp3(x) = VectorValue(0.0,-vApp)
    U_Disp = TrialFESpace(V0_Disp ,[uApp1 ,uApp2 ,uApp3])
    dA_Disp(v,pth,uh_in,sh_in) =  ((p->Em(p))∘pth) * (ε(v) ⊙ (dσfun∘(ε(uh_in), sh_in)))
    return  assemble_vector(v -> ∫(dA_Disp(v,pth,uh_in,sh_in))fem_params.dΩ ,U_Disp)
end

Is there any way that gridap do this calculation and considering nodal number without an explicit function?
The issue that i am facing here is that due to boundary conditions the number of that global vector is not matching with the number of nodes so some values of that vector can not be used.
let me know if you have any questions to clarify
thank you