Averaged volume of a scalar variable and its gradient

General
1

Hi Everyone,
Hope you are doing great. I am looking for obtaining the averaged volume of a scalar parameter as well as the same for the gradient of that scalar parameter, at each gauss quadrature point of a 2D domain. What I am trying to do is given below:

from dolfin import *
mesh = RectangleMesh(Point(0., 0.), Point(L, H), nFEx, nFEy)



lmbda = Constant(24.0)             # GPa = nN/nm^2
mu = Constant(19.4)           	   # GPa = nN/nm^2
kappa = Constant(0.0878)           # nJ/m
kappaT = as_tensor([[2*kappa, 0], [0, kappa/2.]])

ss = as_vector(( cos(thta), sin(thta)))
mm = as_vector((-sin(thta), cos(thta)))


BoundaryU = Expression(('gam*x[1] + 1e-3*t'), gam = 0.07, t = 0, degree = 1)

# Function Spaces
V_u = VectorElement('CG', mesh.ufl_cell(), 1) # Vector finite element
V_eta = FiniteElement('CG', mesh.ufl_cell(), 1) # Scalar finite element
V = FunctionSpace(mesh, MixedElement([V_u, V_eta]))

### Define the space for the gradient of the scalar parameter
V_geta = VectorFunctionSpace(mesh, 'DG', 0)
##########################################################
U_ = TestFunction(V)
(u_, eta_) = split(U_)

U_tr = TrialFunction(V)

dU = Function(V)
(du, deta) = split(dU)

Uold = Function(V)
(uold, eta_old) = split(Uold)

delta = Identity(2)
i, j, k, l, p, q, r, s = indices(8)

# Interpolation Function
phi_eta = deta**2*(3-2*deta)

## Gradient of eta
grad_eta = as_vector(deta.dx(i), i)

### Obtaining the volume average of the variables
volume_element = Constant(1.0)
volume = assemble(volume_element * dx(mesh))
dx = Measure("dx", domain=mesh)

etaAVG = np.zeros((num_steps+1, 1))
gradetaAVG = np.zeros((num_steps+1, 2))

for n in range(num_steps):
    BoundaryU.t = t
    solve(res == 0, dU, bc, J=Gain, solver_parameters={"newton_solver": {"linear_solver": "cg", "preconditioner": "ilu", "relative_tolerance": 1e-6, "absolute_tolerance": 1e-7}},
          form_compiler_parameters={"cpp_optimize": True, "quadrature_degree": 2})
    
    ### project the order parameter into its function space
    etaScalar = project(deta, FunctionSpace(mesh, "CG", 1))
    ### Compute the volume averaged of the order parameter, etaAVG
    etaAVG[n+1,0] = assemble(etaScalar*dx) / volume
    print(assemble(etaScalar*dx) / volume)

    grad_etaScalar = project(grad(dU.split()[1]), V_geta)
    for i in range(2):
        grad_component_i = grad_etaScalar[i]
        print(f"Component {i} of gradient of etaScalar:")
        print(assemble(grad_component_i * dx)/volume)

    Uold.assign(dU)
    t += Dt

However, the values for the components of the gradient of the order parameters are zero.
Your help regarding this would be highly appreciated.