Using a scalar function as component of a vector function

dolfinx
1

Hello fenics community,

I am currently porting my code from fenics2019 to fenicsx and I’m having some issues. I can probably explain myself better if I explain what my fenics2019 code did: I solve a scalar equation and then use that solution as initial condition for one of the components of a vector function (using a mixed element space)

import dolfin as fem

# solve the scalar equation (I understand how to do this in fenicsx)
Vs = fem.FunctionSpace(mesh, 'CG', 1)

u, v = fem.TrialFunction(Vs), fem.TestFunction(Vs)

a = (fem.dot(fem.grad(u), fem.grad(v)) + u*v)*fem.dx + u*v*fem.ds

eq_mag = fem.Function(Vs)

fem.solve(a == v*fem.dx, eq_mag)

# the idea is to solve a vector equation using the obtained scalar
# solution as initial condition for the first component of the
# vector solution and zero for the second component
P = fem.FiniteElement('P', mesh.ufl_cell(), 1)

V = fem.FunctionSpace(mesh, fem.MixedElement((P, P)))

u, v = fem.TestFunctions(V)

m = fem.Function(V)

# from here on I don't know how to translate this to fenicsx
vertex2dof = fem.vertex_to_dof_map(V)

m_local   = m.vector().get_local()
meq_local = eq_mag.vector().get_local()

for i in range(mesh.coordinates().shape[0]):
    m_local[vertex2dof[2*i+0]] = meq_local[vertex2dof[2*i]//2]
    m_local[vertex2dof[2*i+1]] = 0.

m.vector().set_local(m_local)
m.vector().apply('insert')

# at this point I can split m into its components and solve my equations
mx, my = fem.split(m)

So my question is how to do this. Also, I can’t find the equivalent of vertex_to_dof_map which I would like to know even if it isn’t needed for this particular case.

Thank you,
Miguel

2

Perhaps not as efficient as directly copying DoFs to the subvector, but consider interpolation:

import dolfinx
import ufl
from mpi4py import MPI

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 32, 32)
Vs = dolfinx.fem.FunctionSpace(mesh, ('CG', 1))

P = ufl.FiniteElement('P', mesh.ufl_cell(), 1)
V = dolfinx.fem.FunctionSpace(mesh, ufl.MixedElement((P, P)))

eq_mag = dolfinx.fem.Function(Vs)
eq_mag.interpolate(lambda x: x[0] * x[1])

m = dolfinx.fem.Function(V)
m.sub(0).interpolate(eq_mag)


import pyvista, febug
plotter = pyvista.Plotter(shape=(1, 3))
plotter.subplot(0, 0)
febug.plot_function(eq_mag, plotter=plotter)
plotter.add_title("eq_mag")

plotter.subplot(0, 1)
febug.plot_function(m.sub(0).collapse(), plotter=plotter)
plotter.add_title("m.sub(0)")

plotter.subplot(0, 2)
febug.plot_function(m.sub(1).collapse(), plotter=plotter)
plotter.add_title("m.sub(1)")

plotter.show()

image
image877×432 21.9 KB

1 Like