Quoting Bartosz Sawicki <sawickib@xxxxxxxxxxxxx>:
On 30/01/09 01:04 AM, Kristian Oelgaard wrote:
Quoting Bartosz Sawicki <sawickib@xxxxxxxxxxxxx>:
I have a problem how to define form for integration over internal facets
inside VectorElement("Discontinuous Lagrange", "tetrahedron", 0).
Having vector function "j" which is constant over element, I want to
integrate over all internal facets term dot(j,n), where n is facet
normal vector.
Consider form:
=====
element = FiniteElement("Discontinuous Lagrange", "tetrahedron", 0)
vectorelement = VectorElement("Discontinuous Lagrange", "tetrahedron", 0)
v = TestFunction(element)
e = TrialFunction(element)
j = Function(vectorelement)
n = FacetNormal("tetrahedron")
a = v*e*dx
L = v*dot(j, n)*dS
=====
Unfortunately it doesn't work. FFC communicates that "Integrand must be
restricted ('+') or ('-') in an interior facet integral". That is
confusing for me.
Functions on interior facets can potentially be double valued (in your case
they
really are). So you need to tell FFC which value you want.
Maybe you want something like:
L = v*jump(j, n)*dS = v*(dot(j('+'), n('+')) + dot(j('-'), n('-')))*dS
The discontinuous operators are explained in more detail in the FFC
manual.
Thanks for your answer.
I've read FCC manual carefully, but there is "..., one _may_ restrict
expressions to the positive or negative side of the facet". What sugest
that I don't have to use restrictions.
It seems that, when using dS operator we _have to_ restrict every
function. Even your example doesn't compile, because v in not restricted.
Embarrassing... but yes, you have to restrict all functions when computing
interior facet integrals.
Something what compiles correctly is:
a = v*e*dx
L = v('+')*dot(j('+'),n('+'))*dS + v('-')*dot(j('-'),n('-'))*dS
but I'm not sure if it is what I need. Moreover when I've tried to a
solve the problem:
I don't know what you want to compute, so I can't help you here.
Function e;
VariationalProblem problem(a,L);
problem.solve(e);
I get runtime error "*** Error: Coefficient 1 has not been defined."
So something is definitely wrong with this form. :(
Did you initialise the normal vector? and the functin j in you linear form?
Have a look at the DG demos demo/pde/dg/
Have you got any other idea, how to describe separate integration for
every cell?
The linear form that you have compiling should do just that.
You can shorten it a bit by doing:
L = mult(v, dot(j, n))('+')*dS + mult(v, dot(j, n))('-')*dS
which is equivalent, I tested it :)
Kristian
regrs.
BArtek
Kristian
Maybe interior integration operator works differently?
So, would it be possible to use form language to describe such simple
algorithm?
iterate for every cell c:
e_i = 0
integrate for every facet in cell:
e_i += dot(j,n)
e(c) = e_i
regrds.
BArtek
_______________________________________________
FFC-dev mailing list
FFC-dev@xxxxxxxxxx
http://www.fenics.org/mailman/listinfo/ffc-dev