# Assign() cannot be recognized in dolfin_adjoint

I want to use Picard iteration to solve the stationary optimal control NS equations.
However I enconter a problem:once I import dolfin_adjoint repository,the function assign() cannot be recognized.

# Here is my code:

``````import matplotlib.pyplot as plt
from dolfin import *
import numpy as np
from dolfin_adjoint import *
#Load mesh and subdomains
mesh = Mesh("data/dolfin_fine.xml")
sub_domains = MeshFunction("size_t", mesh, "data/dolfin_fine_subdomains.xml")

#coefficient
Re=0.01
max_iteration=100

#Define function spaces
V_h = VectorElement("Lagrange", mesh.ufl_cell(), 2)
Q_h = FiniteElement("Lagrange", mesh.ufl_cell(), 1)
TH = V_h * Q_h
W = FunctionSpace(mesh, TH)
g=Function(W.sub(0).collapse(),name="control")
P=VectorFunctionSpace(mesh,"Lagrange",2,2)

#No-slip boundary condition for velocity
#x1 = 0, x1 = 1 and around the dolphin
noslip = Constant((0, 0))
bc0 = DirichletBC(W.sub(0), g, sub_domains, 0)

#Inflow boundary condition for velocity
#x0 = 1
inflow = Expression(("-sin(x[1]*pi)", "0.0"), degree=2)
bc1 = DirichletBC(W.sub(0), inflow, sub_domains, 1)

#Collect boundary conditions
bcs = [bc0, bc1]

#Define variational problem
(u, p) = TrialFunctions(W)
(v, q) = TestFunctions(W)
f = Constant((0, 0))
u_n_minus_1=interpolate(f,P)
L = inner(f, v)*dx
w = Function(W,name="State")
ds=ds(subdomain_data=sub_domains)

n=0
for n in range(max_iteration):
#Update current iteration
print("#################",n+1)
# Compute solution
solve(a == L, w, bcs)

#Split the mixed solution using deepcopy
#(needed for further computation on coefficient vector)
(u, p) = w.split(True)

#Compute error at vertices
vertex_values_u_n_minus_1=u_n_minus_1.compute_vertex_values(mesh)
vertex_values_u=u.compute_vertex_values(mesh)
error_max = np.max(np.abs(vertex_values_u - vertex_values_u_n_minus_1))

print(error_max)
if error_max<1e-10:
break
#Update previous solution
u_n_minus_1.assign(u)
n+=1

#Split the mixed solution using a shallow copy
(u, p) = w.split()
alpha=Constant(10)

J = assemble(1./2*inner(grad(u), grad(u))*dx + alpha/2*inner(g, g)*ds(0))
m = Control(g)
Jhat = ReducedFunctional(J, m)

g_opt = minimize(Jhat, options={"disp": True})
plot(g_opt, title="Optimised boundary")

g.assign(g_opt)
A,b=assemble_system(a,L,bcs)
solve(A,w.vector(),b)

plot(w.sub(0),title="Velocity")
plot(w.sub(1),title="Pressure")

#Display plots
plt.show()
``````

# Here is my error:

``````
Traceback (most recent call last):
File "/home/rsrg/PycharmProjects/pde_project/ns_iterative_bug_algorithm.py", line 80, in <module>
u_n_minus_1.assign(u)
File "/home/rsrg/.local/lib/python3.6/site-packages/fenics_adjoint/types/function.py", line 70, in assign
ret = super(Function, self).assign(other, *args, **kwargs)
File "/usr/lib/python3/dist-packages/dolfin/function/function.py", line 408, in assign
self._cpp_object._assign(rhs._cpp_object)

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
***
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to initialize vector of degrees of freedom for function.
*** Reason:  Cannot re-initialize a non-empty vector. Consider creating a new function.
*** Where:   This error was encountered inside Function.cpp.
*** Process: 0
***
*** DOLFIN version: 2019.1.0
*** Git changeset:  unknown
*** -------------------------------------------------------------------------
``````

Please format your code by using ``` to encapsulate it. Additionally, Please include the full error message.

To solve your issue, you can use a function assigner (see for example: https://raw.githubusercontent.com/dolfin-adjoint/pyadjoint/master/tests/fenics_adjoint/test_functionassigner.py)
The problem with the code above is that you are creating multiple duplicate function-spaces.

I have followed your suggestion.However,assigner cannot be recognized in

``````g_opt = minimize(Jhat, options={"disp": True})
``````

# Here is the full code.

``````import matplotlib.pyplot as plt
from fenics import *
from dolfin_adjoint import *
import numpy as np

# Load mesh and subdomains
mesh = Mesh("data/dolfin_fine.xml")
sub_domains = MeshFunction("size_t", mesh, "data/dolfin_fine_subdomains.xml")

#coefficient
Re=0.01
max_iteration=100

# Define function spaces
V = VectorElement("Lagrange", mesh.ufl_cell(), 2)
Q = FiniteElement("Lagrange", mesh.ufl_cell(), 1)
TH = V * Q
W = FunctionSpace(mesh, TH)
P=VectorFunctionSpace(mesh,"Lagrange",2)

u_n_minus_1=Function(P)
u_n_minus_1=interpolate(Constant((0, 0)),P)

# No-slip boundary condition for velocity
# x1 = 0, x1 = 1 and around the dolphin
noslip = Constant((0, 0))
bc0 = DirichletBC(W.sub(0), noslip, sub_domains, 0)

# Inflow boundary condition for velocity
# x0 = 1
inflow = Expression(("-sin(x[1]*pi)", "0.0"), degree=2)
bc1 = DirichletBC(W.sub(0), inflow, sub_domains, 1)

# Collect boundary conditions
bcs = [bc0, bc1]

# Define variational problem
(u, p) = TrialFunctions(W)
(v, q) = TestFunctions(W)
f=Constant((0,0))

L = inner(f, v)*dx
w = Function(W)
g = interpolate(Constant((0,0)), W.sub(0).collapse())
ds= Measure("ds",domain=mesh,subdomain_data=sub_domains)

#Time-stepping
n=0
for n in range(max_iteration):
#Update current iteration
print("#################",n+1)
# Compute solution
solve(a == L, w, bcs)

# Split the mixed solution using deepcopy
# (needed for further computation on coefficient vector)
(u, p) = w.split(True)
# print("Norm of velocity coefficient vector: %.15g" % u.vector().norm("l2"))
# print("Norm of pressure coefficient vector: %.15g" % p.vector().norm("l2"))

#Compute error at vertices
vertex_values_u_n_minus_1=u_n_minus_1.compute_vertex_values(mesh)
vertex_values_u=u.compute_vertex_values(mesh)
error_max = np.max(np.abs(vertex_values_u - vertex_values_u_n_minus_1))

print(error_max)

if error_max<1e-10:
break
#Update previous solution
fa= FunctionAssigner(P,W.sub(0))
fa.assign(u_n_minus_1,u)
n+=1

# Split the mixed solution using a shallow copy
(u, p) = w.split()

# optimization
alpha=Constant(10)

J = assemble(1./2*inner(grad(u), grad(u))*dx + alpha/2*inner(g, g)*ds(0))
Jhat = ReducedFunctional(J, Control(g))

g_opt = minimize(Jhat, options={"disp": True})
plot(g_opt, title="Optimised boundary")

g.assign(g_opt)
A,b=assemble_system(a,L,bcs)
solve(A,w.vector(),b)

plot(w.sub(0),title="Velocity")
plot(w.sub(1),title="Pressure")
``````

# error

``````Traceback (most recent call last):
File "/home/rsrg/PycharmProjects/pde_project/ns_iterative_bug_algorithm.py", line 91, in <module>
g_opt = minimize(Jhat, options={"disp": True})
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/tape.py", line 46, in wrapper
return function(*args, **kwargs)
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/optimization/optimization.py", line 254, in minimize
opt = algorithm(rf_np, **kwargs)
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/optimization/optimization.py", line 135, in minimize_scipy_generic
res = scipy_minimize(J, m_global, method=method, **kwargs)
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/_minimize.py", line 618, in minimize
callback=callback, **options)
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/lbfgsb.py", line 308, in _minimize_lbfgsb
finite_diff_rel_step=finite_diff_rel_step)
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/optimize.py", line 262, in _prepare_scalar_function
finite_diff_rel_step, bounds, epsilon=epsilon)
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/_differentiable_functions.py", line 76, in __init__
self._update_fun()
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/_differentiable_functions.py", line 166, in _update_fun
self._update_fun_impl()
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/_differentiable_functions.py", line 73, in update_fun
self.f = fun_wrapped(self.x)
File "/home/rsrg/.local/lib/python3.6/site-packages/scipy/optimize/_differentiable_functions.py", line 70, in fun_wrapped
return fun(x, *args)
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/reduced_functional_numpy.py", line 36, in __call__
return self.rf.__call__(self.set_local(m_copies, m_array))
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/tape.py", line 46, in wrapper
return function(*args, **kwargs)
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/reduced_functional.py", line 135, in __call__
blocks[i].recompute()
File "/home/rsrg/.local/lib/python3.6/site-packages/pyadjoint/block.py", line 345, in recompute
prepared = self.prepare_recompute_component(inputs, relevant_outputs)
File "/home/rsrg/.local/lib/python3.6/site-packages/fenics_adjoint/types/function_assigner.py", line 91, in prepare_recompute_component
self.assigner.input_spaces.delist(inputs))
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
***
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to assign functions.
*** Reason:  The assigning Function is not in the assigning FunctionSpaces.
*** Where:   This error was encountered inside FunctionAssigner.cpp.
*** Process: 0
***
*** DOLFIN version: 2019.1.0
*** Git changeset:  unknown
*** -------------------------------------------------------------------------``````

It will be easier to work out what’s up if you can create a smaller example that gives this error (see Read before posting: How do I get my question answered?).

In particular, it would be helpful if you could create an example on a smaller (preferably built-in) mesh and remove as much code as possible to reduce the number of forms and solvers involved.

Next time, please follow @mscroggs and use the following code as a guideline for a minimal example:

``````from fenics import *
from dolfin_adjoint import *

mesh = UnitSquareMesh(2,2)

V = VectorElement("Lagrange", mesh.ufl_cell(), 2)
Q = FiniteElement("Lagrange", mesh.ufl_cell(), 1)
TH = V * Q
W = FunctionSpace(mesh, TH)
P = W.sub(0).collapse()
Z = W.sub(1).collapse()

u_n_minus_1=interpolate(Constant((0, 0)),P)
w = Function(W)
n=0
while n< 5:
w.vector()[:] = n
u,p = w.split(True)
print("-"*10, "n = {0:d}".format(n),"-"*10)
print("Pre assign", u_n_minus_1.vector().get_local())

fa= FunctionAssigner(P,W.sub(0))
fa.assign(u_n_minus_1,u)
print("Post assign", u_n_minus_1.vector().get_local())
n+=1

Jhat = ReducedFunctional(assemble(inner(u,u)*dx), Control(w))
Jhat(w)
``````

This tiny code, which is runnable for anyone, reproduces your problem.

This can be fixed by doing the following:

``````from fenics import *
from dolfin_adjoint import *

mesh = UnitSquareMesh(2,2)

V = VectorElement("Lagrange", mesh.ufl_cell(), 2)
Q = FiniteElement("Lagrange", mesh.ufl_cell(), 1)
TH = V * Q
W = FunctionSpace(mesh, TH)
P = W.sub(0).collapse()
Z = W.sub(1).collapse()

u_n_minus_1=interpolate(Constant((0, 0)),P)
w = Function(W)
n=0
fa= FunctionAssigner([P,Z], W)
while n< 5:
w.vector()[:] = n
u,p = w.split(True)
print("-"*10, "n = {0:d}".format(n),"-"*10)
print("Pre assign", u_n_minus_1.vector().get_local())
fa.assign([u_n_minus_1, Function(Z)],w)
print("Post assign", u_n_minus_1.vector().get_local())
n+=1

Jhat = ReducedFunctional(assemble(inner(u,u)*dx), Control(w))
Jhat(w)
``````

Note the following changes:

• The function assigner is created outside the for-loop and is reused.
• The function assigner operates from the full space to both the sub spaces (W->(P,Z). This matches the usage in the link from my previous reply.