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

Interface problem using discontinuous Galerkin Method

0 votes

I have been studying multi-phase materials with one diffusive system.

That mean my domain has two subdomains. The solution of the system is discontinuous.
I am tryng to use DG method. The equation of the system are discribe below:

$-\nabla \cdot \beta \nabla u =f \space in \space \space \Omega_1 \cup \Omega_2, $

$u=g \space on \space \partial \Omega $

$[u]= a \space on \space \Gamma_I$

$[\beta \partial_n u]= bn$

$ A(u,v)= \int_{\Omega }^{}{\beta \nabla u \cdot \nabla v \space dx} - \int_{\Gamma } \left( [u] \cdot \lbrace \beta \nabla v \rbrace + [v] \cdot \lbrace \beta \nabla u \rbrace \right)ds$
$+\sum\limits_{e\in \Gamma }{ \frac{\eta}{h_e} \int [u] \cdot [v] }dS $

$ L(v)= \int_{\Omega}{f\space v \space dx}+ \int_{\Gamma_I}{b \lbrace v \rbrace }ds -\int_{\Gamma_I} an \cdot \lbrace D \nabla u \rbrace ds $

$- \int_{\partial \Omega}{gn \cdot \lbrace \beta \nabla v \rbrace}ds +\sum\limits_{e\in \partial \Omega }{ \frac{\eta}{h_e} \int g \space [v] }dS $

$Domain = \Omega = \Omega_1 \cup \Omega_2$

$\Gamma_I = \text{Intern interface between }\space \Omega_1 and \space \Omega_2$

$\partial \Omega = \text{boundary for Dirichilet condition}$

$\Gamma = \text{all interfaces}$

The strong formulation and problem and figures was taken from this ddf in internet (Paper with similar problem

But i would like to use my formulation.

f = Expression('0.0', degree= 1)
u0 = Expression('100.0', degree= 1)
u02 = Expression('0.0', degree= 1)
ds = Measure('ds', domain=mesh, subdomain_data=boundaries)
dS = Measure('dS', domain=mesh, subdomain_data=boundaries)
dx = Measure('dx', domain=mesh, subdomain_data=subdomains)
#Define normal vector and mesh size
n = FacetNormal(mesh)
h = CellSize(mesh)
h_avg = (h('+') + h('-'))/2.0
u = TrialFunction(DG)
v = TestFunction(DG)
a = dot(D*grad(v), grad(u))*dx\
- dot(jump(u, n), avg(D*grad(v)))*dS(0)\
- dot(jump(v, n), avg(D*grad(u)))*dS(0)\
+ alpha/h_avg*dot(jump(v, n), jump(u, n))*dS(0)\
- dot(jump(u, n), avg(D*grad(v)))*dS(3)\
- dot(jump(v, n), avg(D*grad(u)))*dS(3)\
- lamb/h_avg*dot(jump(v, n), jump(u, n))*dS(3)\ 
a+=- dot(D*grad(v), u*n)*ds(1) - dot(v*n, D*grad(u))*ds(1)\
- dot(D*grad(v), u*n)*ds(2) - dot(v*n, D*grad(u))*ds(2)\
+ (gamma/h)*v*u*ds(1) + (gamma/h)*v*u*ds(2)
L = v*f*dx- u0*dot(D*grad(v), n)*ds(1) + (gamma/h)*u0*v*ds(1) \
-u02*dot(D*grad(v), n)*ds(2) + (gamma/h)*u02*v*ds(2)\
+u0*avg(v)*dS(3)\
+dot(u0*n,avg(D*grad(v)))*dS(3)

Can someone helpme with this particle this integral?
+dot(u0n,avg(Dgrad(v)))*dS(3)

ir doesnt work

error:

Calling FFC just-in-time (JIT) compiler, this may take some time.
Discontinuous type FacetNormal must be restricted.
Traceback (most recent call last):
File "/media/leocosta/C863-4EBF/fotos/fick/test_disc.py", line 149, in
solve(a==L , C)
File "/usr/lib/python2.7/dist-packages/dolfin/fem/solving.py", line 300, in solve
_solve_varproblem(*args, **kwargs)
File "/usr/lib/python2.7/dist-packages/dolfin/fem/solving.py", line 325, in _solve_varproblem
form_compiler_parameters=form_compiler_parameters)
File "/usr/lib/python2.7/dist-packages/dolfin/fem/solving.py", line 81, in __init__
L = Form(L, form_compiler_parameters=form_compiler_parameters)
File "/usr/lib/python2.7/dist-packages/dolfin/fem/form.py", line 89, in __init__
mpi_comm=mesh.mpi_comm())
File "/usr/lib/python2.7/dist-packages/dolfin/compilemodules/jit.py", line 68, in mpi_jit
return local_jit(*args, **kwargs)
File "/usr/lib/python2.7/dist-packages/dolfin/compilemodules/jit.py", line 133, in jit
"ffc.jit failed with message:\n%s" % (tb_text,))
File "/usr/lib/python2.7/dist-packages/dolfin/cpp/common.py", line 2973, in dolfin_error
return _common.dolfin_error(location, task, reason)
RuntimeError:

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at


*** fenics-support@googlegroups.com


*** Remember to include the error message listed below and, if possible,
*** include a minimal running example to reproduce the error.


*** -------------------------------------------------------------------------
*** Error: Unable to perform just-in-time compilation of form.
*** Reason: ffc.jit failed with message:
Traceback (most recent call last):
File "/usr/lib/python2.7/dist-packages/dolfin/compilemodules/jit.py", line 128, in jit
result = ffc.jit(ufl_object, parameters=p)
File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 198, in jit
module = jit_build(ufl_object, module_name, parameters)
File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 120, in jit_build
generate=jit_generate)
File "/usr/lib/python2.7/dist-packages/dijitso/jit.py", line 160, in jit
header, source, dependencies = generate(jitable, name, signature, params["generator"])
File "/usr/lib/python2.7/dist-packages/ffc/jitcompiler.py", line 66, in jit_generate
prefix=module_name, parameters=parameters, jit=True)
File "/usr/lib/python2.7/dist-packages/ffc/compiler.py", line 141, in compile_form
prefix, parameters, jit)
File "/usr/lib/python2.7/dist-packages/ffc/compiler.py", line 183, in compile_ufl_objects
analysis = analyze_ufl_objects(ufl_objects, kind, parameters)
File "/usr/lib/python2.7/dist-packages/ffc/analysis.py", line 89, in analyze_ufl_objects
for form in forms)
File "/usr/lib/python2.7/dist-packages/ffc/analysis.py", line 89, in
for form in forms)
File "/usr/lib/python2.7/dist-packages/ffc/analysis.py", line 166, in _analyze_form
form_data = compute_form_data(form)
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/compute_form_data.py", line 309, in compute_form_data
form = apply_restrictions(form)
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/apply_restrictions.py", line 179, in apply_restrictions
only_integral_type=integral_types)
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/map_integrands.py", line 58, in map_integrand_dags
form, only_integral_type)
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/map_integrands.py", line 39, in map_integrands
for itg in form.integrals()]
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/map_integrands.py", line 46, in map_integrands
return itg.reconstruct(function(itg.integrand()))
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/map_integrands.py", line 57, in
return map_integrands(lambda expr: map_expr_dag(function, expr, compress),
File "/usr/lib/python2.7/dist-packages/ufl/corealg/map_dag.py", line 37, in map_expr_dag
result, = map_expr_dags(function, [expression], compress=compress)
File "/usr/lib/python2.7/dist-packages/ufl/corealg/map_dag.py", line 84, in map_expr_dags
r = handlersv._ufl_typecode_
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/apply_restrictions.py", line 166, in facet_normal
return self._opposite(o)
File "/usr/lib/python2.7/dist-packages/ufl/algorithms/apply_restrictions.py", line 71, in _opposite
error("Discontinuous type (null) must be restricted." 300._ufl_class_.name)
File "/usr/lib/python2.7/dist-packages/ufl/log.py", line 171, in error
raise self._exception_type(self._format_raw(*message))
UFLException: Discontinuous type FacetNormal must be restricted.
.
*** Where: This error was encountered inside jit.py.
*** Process: 0


*** DOLFIN version: 2016.2.0
*** Git changeset: unknown
*** -------------------------------------------------------------------------

asked May 3, 2017 by LeoCosta FEniCS User (1,190 points)
edited May 5, 2017 by LeoCosta

1 Answer

0 votes

When integrating over facets, you have to specify on which side of the facet you want to integrate. See this question with the minimal working example provided in the answer.

answered May 3, 2017 by jmmal FEniCS User (5,890 points)

I substitute for "+dot(u0n,avg(Dgrad(v('-'))))*ds(3) " but it doest work too.

It does the another error:

raise self._exception_type(self._format_raw(*message))
Exception: Expression is restricted twice: (Grad(Argument(FunctionSpace(Mesh(VectorElement(FiniteElement('Lagrange', triangle, 1), dim=2), 0), FiniteElement('Discontinuous Lagrange', triangle, 1)), 0, None)),)
.
*** Where: This error was encountered inside jit.py.
*** Process: 0


*** DOLFIN version: 2016.2.0
*** Git changeset: unknown

Look at the definition of avg: it already includes restrictions and so does v('-').

...