Boundary Conditions for geometries defined by discrete points

[lm0bo]

Hello all,

I am working with a geometry where the walls are defined by a set of discrete points rather than a continuous surface representation. I need to correctly apply Dirichlet boundary conditions using the coordinates of these wall points, but I am encountering some challenges in doing so.

Thank you in advance!

[kyouma]

Hello. I would try PointSetBC instead of DirichletBC, or DirichletBC with interpolation inside func and on_boundary parameters.

[lm0bo]

Thank you for your response! I have the coordinates of the points on the wall of a complex 3D geometry, and I want to impose a Dirichlet boundary condition with a velocity of 0 at those points.

I have defined the wall as follows:
geom_wall = dde.geometry.PointCloud(wall_coords, boundary_points=wall_coords)

And the boundary condition:
bc_u_wall = dde.icbc.PointSetBC(wall_coords, np.zeros((wall_coords.shape[0], 1)), component=0)

However, the solution obtained from the network using PointSetBC for the boundary conditions is not accurate.

[kyouma]

Must your equation be solved inside the 3D geometry (volume) or on the wall (surface)? Probably you should use dde.geometry.PointCloud(INNER_VOLUME_POINTS, boundary_points=WALL_POINTS) and dde.icbc.PointSetBC(WALL_POINTS,...)?

Also, soft boundary conditions (conditions set by loss functions) cannot provide accurate solution, there will always be some residue. You can control its size using loss weights in model.compile(). For choosing the weights see the "An Expert's Guide to Training Physics-informed Neural Networks" paper. You could also try to find complex transforms to obtain hard constrains on the boundary.

[LordPontiag]

Hello all,

I am working with a geometry where the walls are defined by a set of discrete points rather than a continuous surface representation. I need to correctly apply Dirichlet boundary conditions using the coordinates of these wall points, but I am encountering some challenges in doing so.

Thank you in advance!

Hello dear @lm0bo , as @kyouma mentioned using PointSetBC is a perfect alternative for DirichletBC. How ever you can use following code that I wrote and used and it works correctly.

PINN_walls =[]
# Define boundary conditions
def boundary_wall(X, on_boundary):
    # Convert X to a NumPy array for consistency
    X = np.array(X)
    tolerance = 1e-4  # Define a small tolerance for floating-point comparison
    wall_cond = any(np.all(np.isclose(wall_coords_2D, X, atol=tolerance), axis=1))
    # Combine with the on_boundary flag
    on_wall = wall_cond & on_boundary
    if on_wall:
        PINN_walls.append(X)
        print(f"Point:",X,"is recognized on_Wall")
    return on_wall

and then call this function:
// here I called it for u and v velocity components, you can extend it to your case

bc_wall_u = dde.DirichletBC(geom, lambda X: 0., boundary_wall, component=0)  # u = 0 on wall
bc_wall_v = dde.DirichletBC(geom, lambda X: 0., boundary_wall, component=1)  # v = 0 on wall

[lm0bo]

Thank you all so much for your advice @LordPontiag @kyouma!
Applying hard constraints on the boundary has significantly improved the results. However, I’m still unsure how to handle the case where the wall is defined by a set of points. Any suggestions on the best way to implement this?

[LordPontiag]

Thank you all so much for your advice @LordPontiag @kyouma! Applying hard constraints on the boundary has significantly improved the results. However, I’m still unsure how to handle the case where the wall is defined by a set of points. Any suggestions on the best way to implement this?

You're most welcome!
I cannot understand your question. both ways work for implementing wall.

[lm0bo]

I apologize if my question wasn't clear. I’m not sure how to apply hard constraints when the wall is defined by a set of discrete points.

Thank you in advance!

[LordPontiag]

I apologize if my question wasn't clear. I’m not sure how to apply hard constraints when the wall is defined by a set of discrete points.

Thank you in advance!

There is no problem with being discrete for the points, make sure to have enough rational points at your boundary, then append their properties to them. The Neural Network will recognize the near region trend and will behave with that region like a solid boundary. for instance if your point sets are wall, it will recognize all wall correctly. You should just provide it with a number of points which add a weight to that area. For example for a 1m length wall, at least having 20 points with increments of 20cm might be a reasonable choice. I'm not sure about this estimation, but if it make sense, it might work!.

After that in your data add those properties of the points. like this:

data = dde.data.PDE(geom,
                   pde,
                   [Wall_U, Wall_V],
                   anchors = wall_coords_2D[::20],
                   num_test = 1000
                   ) 

In this example I directly set the wall coordinates with increments of 20 points as the training points, and their properties as boundary conditions.

[kyouma]

how to apply hard constraints when the wall is defined by a set of discrete points

Perhaps (or more possibly: for sure) my suggestion is not the best option, but you could try using a transformation $u(x) \cdot f(x) + g(x)$, where $u(x)$ is the network output; a (perhaps "continuous" is the best option) $f(x)=0$ for $x \in points_{BC}$; and (perhaps not necessary "continuous") $g(x)=u_{BC}(x)$ for $x \in points_{BC}$. The $f(x)$ may look like $\prod_{x_i \in points_{BC}} (x - x_i)$.

I have not tested this code, but maybe it may help you. I have tried to write a backend-agnostic approach.

• BC_COORDS is a NxD Numpy array of N boundary points, each one is a vector of D input coordinates $(x_1, ..., x_N)$;
• BC_VALUES is a NxC Numpy array of N values of $u_{BC}(x)$ for each boundary point; C is the output size (e.g. for a vector-function).

See the examples:

• Poisson equation in 1D with hard boundary conditions
• Diffusion equation with hard initial and boundary conditions

and papers:

• Exact imposition of boundary conditions with distance functions in physics-informed deep neural networks

and for periodic conditions:

• Physics-informed neural networks with hard constraints for inverse design
• A method for representing periodic functions and enforcing exactly periodic boundary conditions with deep neural networks
• the time-periodic case can be found here: An expert's guide to training physics-informed neural networks

import scipy.interpolate
INTERPOLATOR = scipy.interpolate.NearestNDInterpolator(BC_COORDS, BC_VALUES)  # Or LinearNDInterpolator with fill_value!=NaN


def output_transform(coords, values):  # coords ~ MxD backend tensor; values ~ MxC backend tensor
    bc_coords = dde.backend.as_tensor(BC_COORDS, dtype=coords.dtype)  # You can also store this outside of this function for speed, but check for the dtype and device placement
    differences = dde.backend.expand_dims(coords, 1) - dde.backend.expand_dims(bc_coords, 0)  # Calculate pairwise differences of coordinates (be careful: very memory-intensive; more efficient ways must exist; you can try to find one for your backend)
    distances = dde.backend.sum(dde.backend.abs(differences), dim=-1, keepdims=False)  # You get MxN pairwise distances between input points and BC points (any type: must be 0 at x in BC_COORDS and >0 otherwise)
    
    normalized_distances = dde.backend.tanh(distances)  # Normalize to [0; 1] to avoid overflow and rounding errors during multiplication later
    # normalized_distances = dde.backend.sigmoid(distances) * 2 - 1  # Another way (I don't know which one is faster)
    remoteness_from_bc_points = dde.backend.prod(normalized_distances, dim=-1, keepdims=True)  # Mx1 array of elements: 0 when at any BC point, !=0 otherwise
    
    # OR
    
    remoteness_from_bc_points = dde.backend.min(distances, dim=-1, keepdims=True) # An alternative to "normalized_distances = ..." and "remoteness_from_bc_points = ..." (not sure which one is better)
    
    zeroed_network_output = remoteness_from_bc_points * values  # Zero the network outputs at BC points
    condition_values = INTERPOLATOR(dde.backend.to_numpy(coords))  # Any function that produces exactly BC_VALUES at corresponding BC_POINTS
    condition_values = dde.backend.as_tensor(condition_values, dtype=zeroed_network_output.dtype)
    
    return zeroed_network_output + condition_values


net.apply_output_transform(output_transform)

# Create model and train it