dolfin team mailing list archive

[Melanie Jahny <question155482@xxxxxxxxxxxxxxxxxxxxx>]

Question #155482 on DOLFIN changed:
https://answers.launchpad.net/dolfin/+question/155482

    Status: Answered => Open

Melanie Jahny is still having a problem:
Thanks B. Emek Abali for your advices.
The part with the cell_domains works I think.
But I've got two problems left:

1. If I use dS for the interior facet domains I get the following error,
which I cannot interpret:

Code:

bound_mem = (inner(conc*vel,normal)*eta)*dS(0)
vector_bound = assemble(bound_mem, interior_facet_domains = boundary_parts)

Error:
Calling FFC just-in-time (JIT) compiler, this may take some time.
Form argument must be restricted.
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/usr/local/lib/python2.6/site-packages/dolfin/fem/assemble.py", line 100, in assemble
    common_cell=common_cell)
  File "/usr/local/lib/python2.6/site-packages/dolfin/fem/form.py", line 34, in __init__
    (self._compiled_form, module, self.form_data) = jit(form, form_compiler_parameters, common_cell)
  File "/usr/local/lib/python2.6/site-packages/dolfin/compilemodules/jit.py", line 44, in mpi_jit
    return local_jit(*args, **kwargs)
  File "/usr/local/lib/python2.6/site-packages/dolfin/compilemodules/jit.py", line 103, in jit
    return jit_compile(form, parameters=p, common_cell=common_cell)
  File "/usr/local/lib/python2.6/site-packages/ffc/jitcompiler.py", line 64, in jit
    return jit_form(object, parameters, common_cell)
  File "/usr/local/lib/python2.6/site-packages/ffc/jitcompiler.py", line 122, in jit_form
    compile_form(preprocessed_form, prefix=jit_object.signature(), parameters=parameters)
  File "/usr/local/lib/python2.6/site-packages/ffc/compiler.py", line 140, in compile_form
    ir = compute_ir(analysis, parameters)
  File "/usr/local/lib/python2.6/site-packages/ffc/representation.py", line 66, in compute_ir
    irs = [_compute_integral_ir(f, i, parameters) for (i, f) in enumerate(forms)]
  File "/usr/local/lib/python2.6/site-packages/ffc/representation.py", line 183, in _compute_integral_ir
    parameters)
  File "/usr/local/lib/python2.6/site-packages/ffc/quadrature/quadraturerepresentation.py", line 119, in compute_integral_ir
    terms[i][j] = _transform_integrals(transformer, integrals_dict, domain_type)
  File "/usr/local/lib/python2.6/site-packages/ffc/quadrature/quadraturerepresentation.py", line 278, in _transform_integrals
    integrand = propagate_restrictions(integrand)
  File "/usr/local/lib/python2.6/site-packages/ufl/algorithms/propagate_restrictions.py", line 85, in propagate_restrictions
    return RestrictionPropagator().visit(expression)
  File "/usr/local/lib/python2.6/site-packages/ufl/algorithms/transformations.py", line 129, in visit
    r = h(o, *map(self.visit, o.operands()))
  File "/usr/local/lib/python2.6/site-packages/ufl/algorithms/transformations.py", line 133, in visit
    r = h(o)
  File "/usr/local/lib/python2.6/site-packages/ufl/algorithms/propagate_restrictions.py", line 43, in form_argument
    ufl_assert(self.current_restriction is not None, "Form argument must be restricted.")
  File "/usr/local/lib/python2.6/site-packages/ufl/assertions.py", line 20, in ufl_assert
    if not condition: error(*message)
  File "/usr/local/lib/python2.6/site-packages/ufl/log.py", line 124, in error
    raise UFLException(self._format_raw(*message))
ufl.log.UFLException: Form argument must be restricted.


2. If I try to apply Dirichlet Boundary Conditions to the linear problem, I also get an error:

Applying boundary conditions to linear system.
[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Matrix is missing diagonal entry in row 7594!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.1.0, Patch 3, Fri Jun  4 15:34:52 CDT 2010
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Unknown Name on a linux-gnu named nb1-58 by melanie Tue May  3 11:31:05 2011
[0]PETSC ERROR: Libraries linked from /opt/petsc-3.1-p3/lib
[0]PETSC ERROR: Configure run at Wed Jul 14 19:53:49 2010
[0]PETSC ERROR: Configure options --with-shared --with-dynamic --with-c++-support --with-boost=1 --with-boost-dir=/usr/include/boost --with-blas-lapack-lib="[-llapack,-lcblas,-lf77blas,-latlas]" --with-umfpack=1 --with-umfpack-dir=/opt/SuiteSparse --with-parmetis --with-parmetis-dir=/usr/local --with-netcdf -with-netcdf-dir=/usr --with-debugging=no --prefix=/opt/petsc-3.1-p3
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: MatZeroRows_SeqAIJ() line 1407 in src/mat/impls/aij/seq/aij.c
[0]PETSC ERROR: MatZeroRows() line 5006 in src/mat/interface/matrix.c
[0]PETSC ERROR: MatZeroRowsIS() line 5063 in src/mat/interface/matrix.c


Here is the program-code:

mesh = Rectangle(0,0,150,78,150, 78)
mesh.order()

V = FunctionSpace(mesh, "Lagrange", 1)
conc = TrialFunction(V) # solution
conc_prev = Function(V)
eta = TestFunction(V)

#Defining the subdomain

class Channel(SubDomain): 
   def inside(self,x,on_boundary):
        return True if x[1] >= 50.0 else False

subdomains = MeshFunction('uint', mesh, 2)
sub_channel = Channel()
sub_channel.mark(subdomains,0)

velocity = Constant((0.0, -0.2))
# Defining the equation
F_time = inner((conc-conc_prev),eta)*dx(0) 
F_Diff = (inner(grad(conc), grad(eta))*dx(0))  
F_Adv = inner(inner(velocity, grad(conc)),eta)*dx(0)      

# Defining the boundary conditions
boundary_parts = MeshFunction("uint", mesh, 1)

class RobinBoundary_Mem(SubDomain):
    def inside(self, x, on_boundary):
        return x[1] == 50.0 and x[0] > 20.0 and x[0] < 70.0  


Gamma_Mem = RobinBoundary_Mem()
Gamma_Mem.mark(boundary_parts,0)

normal = FacetNormal(mesh)
vel = Constant((0.0, -0.02))

bound_mem = (inner(conc*vel,normal)*eta)*dS(0)
vector_bound = assemble(bound_mem, interior_facet_domains = boundary_parts) # ERROR !!!

#######################################

F = 0.1*(F_Adv + F_Diff + bound_mem) + F_time
a = lhs(F) 
L = rhs(F)

#######################################
# applying dirichlet boundary condition

c_const = Constant(0.0)
def boundary_up(x,  on_boundary):
    return x[1] == 78.0 and on_boundary

bc_up = DirichletBC(V, c_const, boundary_up)


A = assemble(a, cell_domains=subdomains, interior_facet_domains = boundary_parts)
b = assemble(L, cell_domains=subdomains, interior_facet_domains = boundary_parts)

bc_up.apply(A,b)  #  BOUNDARY CONDITION ERROR


 Could you help me to interpet the error messages? Do you know how to solve
the problems?
Thanks a lot in advance. Melanie

-- 
You received this question notification because you are a member of
DOLFIN Team, which is an answer contact for DOLFIN.