Hi Qi,
a product of linear functionals as bilinearform cannot be achieved directly this way. One reason is that the resulting matrix would be completely dense.
However, you can use a mixed formulation by introducing a new unknown representing one integral
[tex]
\mu = \int_{\Omega}p dx,
[/tex]
where the additional unknown is a single scalar.
This can be implemented in NGSolve via
Code:
Q = H1(mesh, order=1)
N = NumberSpace(mesh)
fes = Q*N
(p,mu), (q,lam) = fes.TnT()
a = BilinearForm(fes, symmetric=True)
a += (p - mu)*lam*dx
a += mu*q*dx
Then you obtain a sparse matrix, where only one row and column is completely filled.
Best,
Michael