kinnala / scikit-fem

Simple finite element assemblers
https://scikit-fem.readthedocs.io
BSD 3-Clause "New" or "Revised" License
470 stars 76 forks source link

1D mesh in 2D, 3D space #1074

Closed dad616610 closed 7 months ago

dad616610 commented 7 months ago

Is it possible to generate a mesh of 1D elements in a higher space?

I'm trying to create a mesh, similar to this one, using 1D elements. image

Can I do this with build-in mesh generation? Can I run a simulation if I import such a mesh using meshio?

kinnala commented 7 months ago

There isn't a "proper" support for stuff like that but I have successfully made such simulations work by creating triangles whose edges correspond to the "1D elements". Then I simply apply FacetBasis to integrate over a subset of the edges and finally remove those third extra nodes in each triangle from the linear system by slicing.

So it depends how much effort you are willing to put into it. There are probably better libraries for this purpose though.

kinnala commented 7 months ago

In case you want to try, I made an example where I project $$f(x, y)=\sin(x + y)$$ onto the finite element mesh on a skeleton mesh which is the dark blue part in this picture: Figure_1 I end up with: Figure_2 Code is here:


import numpy as np
from skfem import *

m1 = MeshTri().refined(2)
m2 = MeshTri().refined(2).translated((-1., -1.))

m = (m1 + m2).with_boundaries({
    'skeleton': lambda x: (x[0] == 0) + (x[1] == 0)
})

m.draw(boundaries=True).show()

fbasis = FacetBasis(m, ElementTriP1(), facets='skeleton')

@BilinearForm
def mass(u, v, _):
    return u * v

I = fbasis.get_dofs(facets='skeleton').all()
MII = mass.assemble(fbasis)[I][:, I]

@LinearForm
def proj(v, w):
    return np.sin(w.x[0] + w.x[1])

y = solve(MII, proj.assemble(fbasis)[I])

Y = fbasis.zeros()
Y[I] = y

m.plot(Y, shading='gouraud', colorbar=True).show()

``
kinnala commented 7 months ago

If you want to use external mesh for this purpose, then I suggest you load the mesh with meshio and then try to turn this 1D mesh into 2D mesh by adding one layer of nodes, while somehow preserving the tagging so that you are later able to find the correct edges to integrate over.

dad616610 commented 7 months ago

Thank you for such a detailed answer! I'll play with the code around to see if it fits my needs