Create a 1D mesh in a 2D space

Hi all,

I have a function defined on a 2D mesh (r,z). Now I would like to create a 1D mesh (ultimately multiple) in the r direction and use this 2D-mesh function in th 1D problem.

I think I should try to create a 1D mesh with multiple coordinates but I haven’t yet found examples to do so in dolfinx.

Any pointers? Thanks in advance!

Rem

Actually I found I could modify the geometry of the mesh with:

mesh_1d = dolfinx.mesh.create_interval(MPI.COMM_WORLD, 10, [0, 1])


mesh_1d.geometry.x[:, 1] = mesh_1d.geometry.x[:, 0]
mesh_1d.geometry.x[:, 0] = 2

mesh_1d.geometry.x
array([[2. , 0. , 0. ],
       [2. , 0.1, 0. ],
       [2. , 0.2, 0. ],
       [2. , 0.3, 0. ],
       [2. , 0.4, 0. ],
       [2. , 0.5, 0. ],
       [2. , 0.6, 0. ],
       [2. , 0.7, 0. ],
       [2. , 0.8, 0. ],
       [2. , 0.9, 0. ],
       [2. , 1. , 0. ]])

But then when naively using the Function from the 2D mesh in the 1D problem, I get:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[32], line 27
     22 p0 = 133.3  # Pa
     23 bcs_pressure = [
     24     dolfinx.fem.dirichletbc(dolfinx.fem.Constant(mesh_1d, p0), dofs_bot, V_p),
     25 ]
---> 27 problem = NonlinearProblem(F_pressure, p, bcs=bcs_pressure)
     28 solver = NewtonSolver(MPI.COMM_WORLD, problem)
     30 # solver.rtol = 1e-7

File ~/anaconda3/envs/plasma-centrifuge-env/lib/python3.13/site-packages/dolfinx/fem/petsc.py:922, in NonlinearProblem.__init__(self, F, u, bcs, J, form_compiler_options, jit_options)
    894 def __init__(
    895     self,
    896     F: ufl.form.Form,
   (...)    901     jit_options: typing.Optional[dict] = None,
    902 ):
    903     """Initialize solver for solving a non-linear problem using Newton's method`.
    904 
    905     Args:
   (...)    920         problem = LinearProblem(F, u, [bc0, bc1])
    921     """
--> 922     self._L = _create_form(
    923         F, form_compiler_options=form_compiler_options, jit_options=jit_options
    924     )
...
    220 if len(gdims) != 1:
--> 221     raise ValueError("Found domains with different geometric dimensions.")
    223 return domains

ValueError: Found domains with different geometric dimensions.

So I assume I need to do some interpolation or something.

Do you want these meshes to align with the facet lines of the 2D mesh? If so, I would use create submesh over those interval(s).

This doesn’t change mesh.geometry.dim, which means that these coordinates will not be moved further into assembly. Dolfinx stores 3 coordinates for any mesh, to make it easier/faster for generated code and the core in general.

This one is a bit harder to grasp, as I don’t have the code to reproduce this bit;)

I can restrict it to facets of the 2D mesh yes. Would it raise an error if I try to make a 1D submesh if it doesn’t match?

Ok so this not a good route then?

Working on a MWE

No, there would just be no cells in the submesh, ie the entity list you send into create_submesh would be empty.

I wouldn’t recommend it. For intervals, it is quite easy to make custom meshes, see for instance

1 Like

Ok so made this little example:

import dolfinx
from mpi4py import MPI

import numpy as np


mesh = dolfinx.mesh.create_unit_square(
    MPI.COMM_WORLD, 10, 10, dolfinx.cpp.mesh.CellType.triangle
)


entities = dolfinx.mesh.locate_entities(
    mesh, mesh.topology.dim - 1, lambda x: np.isclose(x[0], 0.5)
)

mesh_1d, entity_map, vmap, _ = dolfinx.mesh.create_submesh(
    mesh, dim=1, entities=entities
)

V_2d = dolfinx.fem.functionspace(mesh, ("CG", 1))
V_1d = dolfinx.fem.functionspace(mesh_1d, ("CG", 1))


u_2d = dolfinx.fem.Function(V_2d)
u_2d.interpolate(lambda x: x[0] + x[1])

u_1d = dolfinx.fem.Function(V_1d)
u_1d.interpolate(u_2d)  # can't interpolate

It seems I can’t interpolate from triangles to intervals

Traceback (most recent call last):
  File "/home/remidm/plasma-centrifuge/mwe.py", line 28, in <module>
    u_1d.interpolate(u_2d)
    ~~~~~~~~~~~~~~~~^^^^^^
  File "/home/remidm/anaconda3/envs/plasma-centrifuge-env/lib/python3.13/site-packages/dolfinx/fem/function.py", line 459, in interpolate
    _interpolate(u0)
    ~~~~~~~~~~~~^^^^
  File "/home/remidm/anaconda3/envs/plasma-centrifuge-env/lib/python3.13/functools.py", line 934, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
           ~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^
  File "/home/remidm/anaconda3/envs/plasma-centrifuge-env/lib/python3.13/site-packages/dolfinx/fem/function.py", line 445, in _
    self._cpp_object.interpolate(u0._cpp_object, cells0, cells1)  # type: ignore
    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
RuntimeError: Cannot interpolate between elements defined on different cell types.

Using this non-matching meshes interpolation function seems to work.
Taken from GJK error in interpolation between non matching second ordered 3D meshes - #6 by chenyongxin

def nmm_interpolate(f_out: dolfinx.fem.Function, f_in: dolfinx.fem.Function):
    """Non Matching Mesh Interpolate: interpolate one function (f_in) from one mesh into
    another function (f_out) with a mismatching mesh

    args:
        f_out: function to interpolate into
        f_in: function to interpolate from

    notes:
    https://fenicsproject.discourse.group/t/gjk-error-in-interpolation-between-non-matching-second-ordered-3d-meshes/16086/6
    """

    dim = f_out.function_space.mesh.topology.dim
    index_map = f_out.function_space.mesh.topology.index_map(dim)
    ncells = index_map.size_local + index_map.num_ghosts
    cells = np.arange(ncells, dtype=np.int32)
    interpolation_data = dolfinx.fem.create_interpolation_data(
        f_out.function_space, f_in.function_space, cells, padding=1e-11
    )
    f_out.interpolate_nonmatching(f_in, cells, interpolation_data=interpolation_data)

Interpolating from a 1D submesh into 2D is open for quite a lot of interpretation, as what should happen wherever there isn’t a 1D mesh matching the 2D mesh?

As you state in your second post, non-matching interpolation will work, and this is because all cells that do not have a cell with a facet «close enough» to the 1D mesh, is assign 0 to all dofs.

What is the particular use-case for moving the 1D to 2D?:slight_smile:

My use case is:

I have this function v defined on a 2D space (z is the horizontal direction, r is the vertical direction):

In this image, I would like to have several vertical lines on which I solve

\frac{dp}{dr} = p \frac{v^2}{r} \\ p = p_0 , \ \text{on} \ r=R_0

@dokken I see there may be some confusion: I want to interpolate from 2D to 1D not from 1D to 2D, as per my previous reply

Just curious, but if you define several physical lines using gmsh, wouldn’t it work?