A closer look at
shape_divergence reveals, that it just calls
cv.dNdx[basefunc, q_point] and
sum(cv.dNdx[basefunc, q_point]). These values are a combination of
cv.dNdξ and the inverse of the jacobian. The first is a precomputed value (only depends on the gauss point location) and the latter is computed in
reinit! for each cell.
So ideally, you’d like to subtype
CellValues and extend it with
dNdξ² which saves the second order derivative information at a given gauss point.
you’d need to specify a constructor, something like: JuAFEM.jl/cell_values.jl at master · KristofferC/JuAFEM.jl · GitHub
JuAFEM.jl/cell_values.jl at master · KristofferC/JuAFEM.jl · GitHub this line would need to be changed to something like
dNdξ²_temp, dNdξ_temp, N_temp = hessian(ξ -> value(func_interpol, basefunc, ξ), ξ, :all)
and further implement a
reinit! that updates the
cv.dNdx² in the same fashion
cv.dNdx is updated here: JuAFEM.jl/cell_values.jl at master · KristofferC/JuAFEM.jl · GitHub