I am getting the following error:

IndexError: Values to set must be a 1D array

The lines of code causing the error are these ones:

self.Vect = VectorFunctionSpace(mesh, ‘P’, 1)
q_function = Function(self.Vect)
q_function.vector()[:] = q

Isn’t it possible to use .vector()[:] with 2D arrays?

please read Read before posting: How do I get my question answered? . The code you provide is not reproducible, and not correctly formatted.

In general, the answer to

Isn’t it possible to use .vector()[:] with 2D arrays?

is no. It would be good to clarify why you would want it to be a 2D array. Is it because you interpret the first column of the 2D array is the first component of `q`

, and the second column is the second component?

Yes, I want the first column of the 2D array to be the first component of `q`

, and the second column to be the second component. Here is a MWE:

```
import numpy as np
from fenics import *
width = 1000 # Width of the mesh (um)
height = 1000 # Height of the mesh (um)
c_0 = 1.4e-05
dt = 0.1
dt_environment = 1
t = 0.0
mesh = RectangleMesh(Point(0, 0), Point(width,height), 10, 10)
W_D1 = FunctionSpace(mesh, 'P', 1)
Vect = VectorFunctionSpace(mesh, 'P', 1)
v = TestFunction(W_D1)
u = Function(W_D1)
u_old = Function(W_D1)
u_initial= f"{c_0}"
u_aux = Expression(u_initial,degree=1, L = width)
u_0 = interpolate(u_aux, W_D1)
u.assign(u_0)
u_old.assign(u_0)
q = np.zeros([mesh.num_vertices(),2])
for vertex in vertices(mesh):
q[vertex_to_dof_map(W_D1)[vertex.index()],:] = np.random.uniform()
q_function = Function(Vect)
# q_function.vector()[:] = q
q_function = project(q,Vect)
F_D1 = (
inner(u - u_old, v) * dx
+ dt * q_function * grad(v) * dx
)
# Solver
for time in np.arange(0, dt_environment, dt):
set_log_level(30)
solver_parameters = {
"newton_solver": {
"linear_solver": "gmres",
"preconditioner": "ilu",
"maximum_iterations": 1000,
"relaxation_parameter": 0.7,
"relative_tolerance": 1e-7,
"absolute_tolerance": 1e-8
}
}
solve(F_D1 == 0, u, solver_parameters=solver_parameters)
u.vector()[:] = np.maximum(u.vector().get_local(), 0)
t += dt
u_old.assign(u)
```

dokken
August 25, 2024, 5:26pm
4
Consider the following approach:

```
from fenics import *
mesh = RectangleMesh(Point(0, 0), Point(width, height), 10, 10)
Vect = VectorFunctionSpace(mesh, "P", 1)
u = Function(Vect)
u1 = u.sub(1, deepcopy=True)
print(u1.vector().get_local())
# Modify u1
u1.vector()[:] = 2
# Create assigner from u1 to u
assigner = FunctionAssigner(Vect.sub(1), u1.function_space())
assigner.assign(u.sub(1), u1)
# Assign to u
print(u.vector().get_local())
File("u.pvd") << u
```