### creating a Fenics Function array

694
views
0
10 months ago by
I want to create a Fenics Function array Phi_s:
#define mesh
self.mesh = BoxMesh(Point(0, 0, 0), Point(self.L, self.W, self.W), 10, 3, 3)
#define function space
self.V = VectorFunctionSpace(self.mesh, 'P', 1)

self.temp = Function(self.V)

for i in range(0,self.MAX_ITER+1):
for j in range(0,self.MAX_ITER+1):
self.Phi_s = self.temp
self.Phi_s.vector()[:] = self.phi_s[i,j,:] # correct syntax but not desired (phi_s is a compatible np array)
self.Phi_s.vector()[i,j,:] = self.phi_s[i,j,:] # desired but incorrect syntax

This gives me the following error:
self.Phi_s.vector()[i,j,:] = self.phi_s[i,j,:]
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1284, in __setitem__
indices = self._check_indices(indices)
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1209, in _check_indices
raise TypeError("expected an int or a list or numpy array of "\
TypeError: expected an int or a list or numpy array of integers or a boolean numpy array as indices.
Aborted

I can see the fatality of this approach because a Fenics VectorFunction can only accept vector coefficients (I think).

To resolve this issue, I define a TensorFunctionSpace instead:

#define mesh
self.mesh = BoxMesh(Point(0, 0, 0), Point(self.L, self.W, self.W), 10, 3, 3)
#define function space
self.V = TensorFunctionSpace(self.mesh, 'P', 1)

self.temp = Function(self.V)

for i in range(0,self.MAX_ITER+1):
for j in range(0,self.MAX_ITER+1):
self.Phi_s = self.temp
self.Phi_s.vector()[:] = self.phi_s[i,j,:] # incorrect syntax and not desired (phi_s is a compatible np array)
self.Phi_s.vector()[i,j,:] = self.phi_s[i,j,:] # desired but incorrect syntax
self.Phi_s.tensor()[i,j,:] = self.phi_s[i,j,:] # desired but incorrect syntax
These result in the following errors:
self.Phi_s.vector()[:] = self.phi_s[i,j,:]
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1194, in __setslice__
self[indices] = values
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1303, in __setitem__
raise IndexError("expected same size of indices and values")
IndexError: expected same size of indices and values
Aborted

self.Phi_s.vector()[i,j,:] = self.phi_s[i,j,:]
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1284, in __setitem__
indices = self._check_indices(indices)
File "/usr/local/home/bhattah/temp/pwd/lib/python2.7/site-packages/dolfin/cpp/la.py", line 1209, in _check_indices
raise TypeError("expected an int or a list or numpy array of "\
TypeError: expected an int or a list or numpy array of integers or a boolean numpy array as indices.
Aborted

self.Phi_s.tensor()[i,j,:] = self.phi_s[i,j,:]
AttributeError: 'Function' object has no attribute 'tensor'
Aborted

How do I create a Fenics Function array? What are the right keywords for a google search here? Thanks!
Community: FEniCS Project

0
10 months ago by
A minimal working example would be helpful.

If you want to set elements of a vector --- as you appear to desire to accomplish by reading your first code snippet --- you can do something like this:

from fenics import *

mesh = UnitSquareMesh(2,2)
V    = VectorFunctionSpace(mesh, 'P', 1)
Q    = FunctionSpace(mesh, 'P', 1)

u    = Function(V)
v    = Function(V)
f    = Function(Q)

f.vector()[:] = 42.0

u_x, u_y  = u.split()

assu_x    =  FunctionAssigner(u_x.function_space(), Q)
assu_y    =  FunctionAssigner(u_y.function_space(), Q)

assu_x.assign(u_x, f)
assu_y.assign(u_y, f)

u_x, u_y  = u.split(True)

from numpy import unique

print unique(u_x.vector().array())
print unique(u_y.vector().array())
​

This will set some elements, and set them efficiently; the "FunctionAssigner" object has been optimized in this regard.

Note also how you can copy and paste my code into an interpreter and reproduce results.  This is what is meant by "minimal working example".

Thanks for the advice but it doesn't quite tackle the problem although my subscript notation might have been misleading.

I am looking for a way to create an array (or list) of fenics Functions. It is mentioned in
https://fenicsproject.org/qa/8182/indexing-part-of-functionspace
on how to accomplish this. So here we go:
from fenics import *
from numpy import unique

m = 4 # declare number of elements in the Function list
mesh = UnitSquareMesh(2,2)
Vm    = VectorFunctionSpace(mesh, 'P', 1, m) # define list of VectorFunctions

u    = Function(Vm) # define list of Functions u

u[0].vector()[:] = 42.0 # assign coefficients for the first Function in the list u[0:m]

print(unique(u.vector().array())) # print the coefficient vector of Function u​

This script throws the following error:

Traceback (most recent call last):
File "test1.py", line 13, in <module>
u[0].vector()[:] = 42.0
AttributeError: 'Indexed' object has no attribute 'vector'
Aborted

Which I can get rid off by modifying the erroneous line to be

u.vector()[:] = 42.0

instead. But this doesn't give me access to the coefficient vector to each Function of the list. I want to be able to assign a value to every individual Function vector separately. I hope this makes my question clear.

fenics 2017.1.0 py27_4 conda-forge  4.3.22
written 10 months ago by Ashish Bhatt
I see, you want to avoid deep copy and set the elements of a mixed function directly.  See here:

written 10 months ago by pf4d
Note that you can set the elements of a Function to values of an array, and then use the FunctionAssigner as I've indicated above to set your individual components as well.  That is, "assu_x" will assign to only the first, "assu_y" will only assign to the second, etc if more were defined.
written 10 months ago by pf4d
I know how to assign values to a Function's array and its components' arrays, I am looking to create an array (or list) of Functions instead. Like so:
from fenics import *

m = 4 # declare number of elements in the Function list
mesh = UnitSquareMesh(2,2)

for i in xrange(m):
V[i]    = VectorFunctionSpace(mesh, 'P', 1, i) # define a VectorFunction

u[i]    = Function(V[i]) # define a Function​

But this, of course, is not the right syntax.

written 10 months ago by Ashish Bhatt
OK, this is an entirely different question...  again.  So I'll answer it:

from fenics import *

mesh = UnitIntervalMesh(3)
P1   = FiniteElement("CG", interval, 1)
M4   = MixedElement([P1]*4)
Q    = FunctionSpace(mesh, M4)
f    = Function(Q)

# your "list of Functions" :
m1,m2,m3,m4 = f
​