Get a transform matrix map design variable to field

Dear experienced fenics users,

I am trying to get a transformation matrix A which computes f = A*d, which f is a scaler field and d is a vector of design variables.

from fenics import *
import numpy as np

mesh = UnitSquareMesh(5, 5)
V = VectorElement("CG", mesh.ufl_cell(), 2)
Vh = FunctionSpace(mesh, V)

#actual expression is f = Expression(("d0*x[0]*x[1] + d1*(x[0]+x[1])"), degree=2)
#trying to get matrix A from differentiated u_D, such that f = A*d 
u_D = Expression(("x[0]*x[1]", "x[0]+x[1]"), degree=2)
u = Function(Vh)

u.interpolate(u_D)
print(np.shape(u.vector().get_local())) #However I get a single column vector instead of two column matrix

I have some questions regarding how fenics discrete vector or mixed function space. If I have a vector field expression in 2D, I will expect the discrete vector is of the shape of (dof_dimension * vector_dimension). However, as the following example shows, I’m getting a single column vector instead of a matrix.

The output is (242,), I think it is of dimension (2*dof, ). So does it make sense if I just make it shape of (dof, 2) by cut it in the middle.

Hi, some additional thoughts,

so if I have a variational problem formulated in a vector function space and I assembled the system getting a linear system A*x = b.

Did I get A of shape (dof, dof) and b of shape (dof, vector_dimension)
or I get A of shape(dof, dofvector_space) and b of shape (dofvector_space, )

I found a workaround,

A1 = u.sub(0).vector()
A2 = u.sub(1).vector()

Concatenate these two vectors gave me the transformation matrix A.

But it will be nice if someone could tell me proper ways to get such a matrix, it will be troublesome when design parameters gets more.