Averaging quantities over the domain efficiently

Question

I have a code that solves for u and p. Suppose I want to average u*p over all elements, what I do currently is

import numpy as np
...
u, p = u_.split(deepcopy=True)[:]  # get u and p from a mixed function space
u_array = u.vector().array()
p_array = p.vector().array()
avgUP = np.average(u_array*p_array)

I wish to run this code in parallel.
Is there some way to compute avgUP with fenics (without invoking numpy)?

Answer 1

Try:

m = u.vector()*v.vector()
avgUP = m.sum()/m.size()

Comments

Thank you! This is such a simple and elegant solution. wonder why I haven't thought about this before :D

---

suppose I have to average curl(u)*c or a component of u times c,

cu = curl(u).vector()*c.vector()   # doesn't work
I = Identity(3)
uz =  dot(u, I[:, 2]) # is a dot object and has no attribute 'vector' 
Ew = uz.vector()*c.vector()  # fails
Ew_avg = Ew.sum()/Ew.size()

what should I do?

---

uz is a UFL expression and does not have a vector attribute. You option for now is to first project it to a Function.

Ew = project(uz*cu)

It would be nice to be able to interpolate it to a Function as project is less acuaret and involves solving a linear system.

Dependent on what you want you might also assemble an average value:

vol = assemble(Constant(1)*dx, mesh=mesh)
Ew = assemble(uz*dx)/vol