mesh-adaptation / movement

Mesh movement methods for finite element problems solved using Firedrake
MIT License
6 stars 4 forks source link

Monge-Ampere boundary movement for non-rectangular domains #13

Open jwallwork23 opened 2 years ago

jwallwork23 commented 2 years ago

Consider solving Poisson on a mesh of a non-rectangular, polygonal domain.

If all straight line boundary segments are uniquely tagged then the Monge-Ampere mover will fix the corners and only allow boundary movement tangential to each segment. However, this does not avoid the situation in the following plot, where the tangential movement passes a corner.

image

It isn't obvious how to fix this, but here are a few possibilities:

jwallwork23 commented 3 weeks ago

I think this should be fixed now? Presumably the image was generated with an old version of Movement.

jwallwork23 commented 2 weeks ago

@ddundo will try to replicate this issue.

jwallwork23 commented 2 weeks ago

We discussed in the meeting and determined that the DirichletBC for the EquationBC doesn't guarantee that nodes don't get moved beyond the corners (it doesn't know that the values have to do with mesh entities).

jwallwork23 commented 2 weeks ago

It's possible that there's no way to guarantee this effect doesn't manifest, i.e., problem is ill-posed.

ddundo commented 2 weeks ago

I tried to replicate this using this code: https://gist.github.com/ddundo/bef044545395ad94e0d28553a80d8578

It all looks good to me! I also tried other monitor functions that focus adaptation more towards the boundary. But please let me know if I missed/did something wrong :)

image

ddundo commented 2 weeks ago

Edit: sorry, boundary segments aren't marked in the gmsh code below.


Hmm just saw that there is some normal movement on this example:

image

Same code as in the gist above, but mesh created like this:

from math import cos, sin

import gmsh

gmsh.initialize()
gmsh.model.add("heptagon")

r = 1  # Radius of the circumscribed circle
n_sides = 7
lc = 0.1  # Characteristic length

points = []
for i in range(n_sides):
    angle = 2 * 3.14159 * i / n_sides
    x = r * cos(angle)
    y = r * sin(angle)
    points.append(gmsh.model.geo.addPoint(x, y, 0, lc))

lines = []
for i in range(n_sides):
    lines.append(gmsh.model.geo.addLine(points[i], points[(i + 1) % n_sides]))

loop = gmsh.model.geo.addCurveLoop(lines)
surface = gmsh.model.geo.addPlaneSurface([loop])

gmsh.model.geo.synchronize()
gmsh.model.mesh.generate(2)

gmsh.write("heptagon.msh")
gmsh.finalize()
jwallwork23 commented 2 weeks ago

Thanks for this @ddundo. Which branch are you using here?

ddundo commented 2 weeks ago

Thanks for this @ddundo. Which branch are you using here?

On 112_fix-specific-segment :)

jwallwork23 commented 2 weeks ago

Hm it looks like the mesh doesn't have any boundary tags. It needs to have Physical ID tags on boundary segments for that to work properly.

Added a warning on the branch cfdc231.

ddundo commented 2 weeks ago

Yeah I forgot about that in the second example, sorry. But the first one does have boundaries marked - and it works well I think