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
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?
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
@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?
