How to get the Function value on some points or the mesh boundary?

dolfin-adjoint
1

Hi everyone,
As a very beginner of fenics, i need someone to help me for this problerm :
While i solve an equation in this mesh using Neumann boundary condition. Now i want to get the Function value in some points or the mesh boundary
This is a simple example:

from fenics import *
import numpy as np
#create mesh
mesh = UnitSquareMesh(2, 2)
x = SpatialCoordinate(mesh)
#bulid FunctionSoace
V = FunctionSpace(mesh, “P”, 1)
#interpolation a function
u = project(sin(x[0]) + 1 , V)
print(u.compute_vertex_values(mesh))

Now, I want to get the u value on some points (not the mesh vertex) or the mesh boundary.

2

See: How to get nemrical solution value for a given point? - #2 by dokken

3

Thank you so much for giving me the example.

May I ask you one more favor, if don’t mind?
I want to know how to get the boundary points from a mesh.

4

What do you imply by boundary points? Do you mean those vertices that are part of a facet that is on the boundary?

5

Thanks for your kind help .
I want to rewrite my problem as flows:

from fenics import *
from mshr import *

mesh1 = UnitSquareMesh(2, 2)
domain =  Rectangle(Point(0, 0), Point(1, 1))
mesh2 = generate_mesh(domain, 30)

We can get mesh1 and mesh2 by running this code .How to get the (left) boundary points coordinate ?
This may be a simple problem in this example, but if the mesh boundary is moved using ALE.move(), what can i do to get the new boundary points coordinate ?
the code is :

from dolfin import *
from mshr import *
import numpy as np
import matplotlib.pyplot as plt


W = 2.0
H = 1.0
domain = Rectangle(Point(0, 0), Point(W, H))
mesh = generate_mesh(domain, 20)


class Right(SubDomain):
    def __init__(self):
        SubDomain.__init__(self)
        
    def inside(self, x, on_boundary):
        return on_boundary and near(x[0], W)
    
right = Right()
#bmesh = BoundaryMesh(mesh, "exterior", True)
bmf = MeshFunction("size_t", mesh, mesh.topology().dim() - 1, 0)
bmf.set_all(0)
right.mark(bmf, 1)

V= VectorFunctionSpace(mesh, "CG", 1)
x = SpatialCoordinate(mesh)
u = Function(V)
#boundary displacement
bcs = project(as_vector((0.1*x[1]*(1 - x[1]), 0)), V)
#bcs = Expression(("0.05*sin(4*(pi/L)*x[1])", "0.0"), L=1, degree=1)
bc = DirichletBC(V, bcs, bmf, 1)
bc.apply(u.vector())
plot(mesh)
plt.show()
ALE.move(mesh, u)
plot(mesh)
plt.show()

How to get the new mesh (right) boundary points coordinate ?

6

You can call the following:


# Find the facets that are marked with 1
facet_indices = np.flatnonzero(bmf.array() == 1)
# Create facet to vertex connectivity
mesh.init(mesh.topology().dim()-1, 0)  
# Get mesh nodes
coords = mesh.coordinates()
# Get facet to vertex connectivity
f_to_c = mesh.topology()(mesh.topology().dim()-1, 0) 
for facet in facet_indices:
    vertices = f_to_c(facet)
    print(coords[vertices])