This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

Computing integral over boundary

0 votes

Dear All,

assuming, I've already solved a pde with solution X. Now, I would like to compute the following integral:

enter image description here

The Problem is, that I can't enter grad( X (s+1 ) ) directly in Fenics (version 1.4.x). Is there an easy way to do this?

I've already tried to create an Expression and than to integrate, but the assemble-command doesn't work here (Error: "Expecting a completed form with domains at this point."):

class Source(Expression):  
    def eval(self, values, x):
        values[0] = inner(grad(X(x[0])),n)*inner(grad(X(x[0]+1)),n)  
f = Source()
int=f*ds(1)
assemble_int=assemble(int)

Here is a minimal code example:

# -*- coding: utf-8 -*-
from dolfin import *

omega=UnitSquareMesh(4,4)
V=FunctionSpace(omega,"CG",1)

# Marking the top of the boundary 
class DirichletBoundaryTop(SubDomain):
    def inside(self, x, on_boundary):
        return bool((x[1] == 1) and on_boundary)

dbt = DirichletBoundaryTop()
mf = MeshFunction("size_t", omega, omega.topology().dim() - 1)
mf.set_all(0)
dbt.mark(mf, 1)

# solving PDE for X
X_Top = Constant(1.0)
DirXtop = DirichletBC(V, X_Top, dbt)
bc_X=[DirXtop]
X = TrialFunction(V)
t = TestFunction(V)
a = inner(nabla_grad(X), nabla_grad(t))*dx
L = Constant(0.0)*t*dx
X = Function(V)
solve(a==L, X, bc_X)

# Compute integral
n=FacetNormal(omega)
ds = Measure('ds')[mf]
integral=inner(grad(X),n)*ds(1)      
integral_value=assemble(integral)
print integral_value

Thank you for your support in advance.

asked Mar 3, 2015 by puffin FEniCS Novice (440 points)
...