When attempting to call as_vector() with too few arguments (in this case, 1 instead of 2), the program runs if we only solving for a single field, but breaks when trying to solve for multiple functions simultaneously. Code below.
=====================================================
from dolfin import *
mesh = UnitSquareMesh(8,8)
S = FunctionSpace(mesh, "BDM", 1)
V = FunctionSpace(mesh, "DG", 0)
W = S*V
u0 = interpolate(Constant((-1.0,1.0)),S)
h0 = interpolate(Constant("1.0"),V)
u = TrialFunction(S)
v = TestFunction(S)
a = dot(as_vector([-u[1]]),v)*dx # works
#a = dot(as_vector([-u[1],0.0]),v)*dx # also works
L = dot(u0,v)*dx
u = Function(S)
solve(a==L,u)
print "First part compiles"
plot(u)
interactive()
u,h = TrialFunctions(W)
v,phi = TestFunctions(W)
a = dot(as_vector([-u[1]]),v)*dx + h*phi*dx # fails
#a = dot(as_vector([-u[1],0.0]),v)*dx + h*phi*dx # works
L = dot(u0,v)*dx + h0*phi*dx
w = Function(W)
solve(a==L,w)
u,h = w.split(deepcopy=True)
print "Second part compiles"
plot(u)
interactive()
=====================================================
atm112@ae-amcrae:~/Dropbox/FEniCS/scratch$ python asvector2.py
Solving linear variational problem.
First part compiles
Calling FFC just-in-time (JIT) compiler, this may take some time.
Traceback (most recent call last):
File "asvector2.py", line 27, in <module>
solve(a==L,w)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/fem/solving.py", line 268, in solve
_solve_varproblem(*args, **kwargs)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/fem/solving.py", line 292, in _solve_varproblem
form_compiler_parameters=form_compiler_parameters)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/fem/solving.py", line 80, in __init__
a = Form(a, form_compiler_parameters=form_compiler_parameters)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/fem/form.py", line 56, in __init__
common_cell)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/compilemodules/jit.py", line 66, in mpi_jit
return local_jit(*args, **kwargs)
File "/home/atm112/local/fenics/src/dolfin/current/local/lib/python2.7/site-packages/dolfin/compilemodules/jit.py", line 154, in jit
return jit_compile(form, parameters=p, common_cell=common_cell)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/jitcompiler.py", line 77, in jit
return jit_form(ufl_object, parameters, common_cell)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/jitcompiler.py", line 197, in jit_form
parameters=parameters)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/compiler.py", line 156, in compile_form
ir = compute_ir(analysis, parameters)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/representation.py", line 96, in compute_ir
for (i, fd) in enumerate(form_datas)]
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/representation.py", line 216, in _compute_integral_ir
parameters)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturerepresentation.py", line 87, in compute_integral_ir
itg_data.domain_type, form_data.cell)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturerepresentation.py", line 145, in _transform_integrals_by_type
terms = _transform_integrals(transformer, integrals_dict, domain_type)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturerepresentation.py", line 363, in _transform_integrals
terms = transformer.generate_terms(integral.integrand(), domain_type)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturetransformerbase.py", line 739, in generate_terms
terms = self.visit(integrand)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ufl/algorithms/transformer.py", line 101, in visit
r = h(o, *map(self.visit, o.operands()))
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ufl/algorithms/transformer.py", line 105, in visit
r = h(o)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturetransformerbase.py", line 549, in index_sum
ops.append(self.visit(summand))
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ufl/algorithms/transformer.py", line 101, in visit
r = h(o, *map(self.visit, o.operands()))
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ufl/algorithms/transformer.py", line 105, in visit
r = h(o)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturetransformerbase.py", line 512, in indexed
code = self.visit(indexed_expr)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ufl/algorithms/transformer.py", line 105, in visit
r = h(o)
File "/home/atm112/local/fenics/lib/python2.7/site-packages/ffc/quadrature/quadraturetransformerbase.py", line 714, in list_tensor
op = o.operands()[c0]
IndexError: tuple index out of range
Aside: the output seems to be crazy, at least compared to the related program where we move the rotation onto the RHS. Not sure if bug or genuine...
======= ======= ======= ======= =
from dolfin import *
mesh = UnitSquareMesh(8,8)
S = FunctionSpace(mesh, "BDM", 1)
V = FunctionSpace(mesh, "DG", 0)
W = S*V
u0 = interpolate( Constant( (-1.0,1. 0)),S) Constant( "1.0"), V)
h0 = interpolate(
u = TrialFunction(S) vector( [-u0[1] ]),v)*dx # works vector( [-u0[1] ,0.0]), v)*dx # also works
v = TestFunction(S)
a = dot(u,v)*dx
#L = dot(as_
L = dot(as_
u = Function(S)
solve(a==L,u)
print "First part compiles"
plot(u)
interactive()
u,h = TrialFunctions(W) vector( [-u0[1] ]),v)*dx + h0*phi*dx # fails vector( [-u0[1] ,0.0]), v)*dx + h0*phi*dx # works deepcopy= True)
v,phi = TestFunctions(W)
a = dot(u,v)*dx + h*phi*dx
#L = dot(as_
L = dot(as_
w = Function(W)
solve(a==L,w)
u,h = w.split(
print "Second part compiles"
plot(u)
interactive()