Create a network geometry

[DeepaMahm]

Discussed in https://github.com/lululxvi/deepxde/discussions/1500

^{Originally posted by **DeepaMahm** September 27, 2023} Dear All, I would like to create a network geometry (3D) in Deepxde to set up a diffusion problem. I looked for the geometry creation options given in the documentation and I think I may need to use `PointCloud`. I've used `vedo` in the past to create pointcloud ``` import networkx as nx from vedo import * G = nx.gnm_random_graph(n=30, m=55, seed=1) nxpos = nx.spring_layout(G, dim=3) nxpts = [nxpos[pt] for pt in sorted(nxpos)] # node values vals = list(range(10, 10+len(nxpts))) nx_lines= [] for i, j in G.edges(): nx_lines.append([nxpos[i], nxpos[j]]) nx_pts = Points(nxpts, r=12) nx_edg = Lines(nx_lines).lw(2) ``` Basically, I have the coordinates of the nodes in the geometry and `Lines` in `vedo` is used to define the edge connectivity i.e. the edge between two nodes. `Points` in `vedo` would result in a pointcloud which I can pass to the class `Geometry` in `Deepxde`. I am a bit unsure if this is the right way tot proceed. And it is not clear to me how the lines have to be defined in `Deepxde`. Suggestions on how to proceed would be of great help. Thanks so much for the amazing library.

[DeepaMahm]

Hi Prof. @lululxvi ,

If I interpolate and obtain the coordinates of points lying on the edge formed between 2 nodes, can these corrdinates be specified in the pointcloud geometry class to get the final network geometry?

Could someone please give me suggestions? I am not sure how to define a line object.

[praksharma]

Hey, PINNs are meshless. This means you don't need lines, surfaces , planes etc. Just pass your the nodal coordinates as anchors. Here is an example,

data = dde.data.PDE(
    geom,
    pde,
    [observed],    # BCs
    num_domain = 0,
    num_boundary = 0,
    anchors=nodes  # collocation points
)

You just plug your collocation points i.e. the nodal coordinates as nodes and BC coordinates and its value as [observed]. You can define observed as follows:

observed = dde.icbc.PointSetBC(boundary_nodal_coordinates, boundary_solution)

[DeepaMahm]

Hi @praksharma ,

Thanks a lot. If I may ask another question,

To solve 1D diffusion problem, I tried to define Dirichlet boundary (value of 5) on the left and Neumann boundary on the right (zero gradient)

import deepxde as dde
import numpy as np

def pde(x, y):
    """
    :param x: 2-dimensional vector; x[:, 0:1] = x-coordinate; x[:, 1:] = t-coordinate
    :param y: network output; i.e., solution y(x,t)
    :return:
    """
    dc_t = dde.grad.jacobian(y, x, i=0, j=1)
    dc_xx = dde.grad.hessian(y, x, i=0, j=0)
    d = 1
    return (
        dc_t
        - d * dc_xx
    )

def boundary_l(x, on_boundary):
    """
    Because of rounding-off errors, dde.utils.isclose(x[0], 0)
    :param x:
    :param on_boundary: all boundary points
    :return:
    """
    return on_boundary and dde.utils.isclose(x[0], 0)

def boundary_r(x, on_boundary):
    """
    Because of rounding-off errors, dde.utils.isclose(x[0], 0)
    :param x:
    :param on_boundary: all boundary points
    :return:
    """
    return on_boundary and dde.utils.isclose(x[0], 1)

if __name__ == '__main__':
    geom = dde.geometry.Interval(0, 1)
    timedomain = dde.geometry.TimeDomain(0, 1)
    geomtime = dde.geometry.GeometryXTime(geom, timedomain)

    bc_l = dde.icbc.DirichletBC(geom, lambda x: 5, lambda _, on_boundary: boundary_l)
    bc_r = dde.icbc.NeumannBC(geom, lambda x: 0, lambda _, on_boundary: boundary_r)

    ic = dde.icbc.IC(geomtime, lambda x: 1, lambda _, on_initial: on_initial)

    data = dde.data.TimePDE(
        geomtime,
        pde,
        [bc_l, bc_r, ic],
        num_domain=40,  # no. of training residual points sampled inside the domain
        num_boundary=20,  # no. of training points sampled on the boundary
        num_initial=10,  # initial residual points for initial conditions
        # solution=func, # reference solution to compute error of the dde solution, can be ignored if no exact soln
        num_test=10000,  # points for testing the residual
    )

    layer_size = [2] + [32] * 3 + [1]  # width 32, depth 3
    activation = "tanh"
    initializer = "Glorot uniform"
    net = dde.nn.FNN(layer_size, activation, initializer)

    model = dde.Model(data, net)
    model.compile("adam", lr=0.001, metrics=["l2 relative error"])
    losshistory, train_state = model.train(iterations=10)

    dde.saveplot(losshistory, train_state, issave=True, isplot=True)

I get the following error

    return X[self.on_boundary(X, self.geom.on_boundary(X))]
           ~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
IndexError: arrays used as indices must be of integer (or boolean) type

I'm not sure how to fix this error.

Could you please have a look the boundary conditions? I'm not sure what's wrong.

    bc_l = dde.icbc.DirichletBC(geom, lambda x: 5, lambda _, on_boundary: boundary_l)
    bc_r = dde.icbc.NeumannBC(geom, lambda x: 0, lambda _, on_boundary: boundary_r)

[praksharma]

The error is pretty obvious. You are using function within function. Try this.

bc_l = dde.icbc.DirichletBC(geom, lambda x: 5, boundary_l)
bc_r = dde.icbc.NeumannBC(geom, lambda x: 0, boundary_r)

[DeepaMahm]

Hi @praksharma , Thank you for the suggestion.

I tried and I get the following error. Could you please have a look?

lambda x, on: np.array([on_boundary(x[i], on[i]) for i in range(len(x))])
To silence this warning, decorate the function with @tf.autograph.experimental.do_not_convert
Traceback (most recent call last):
  File "I:\xxx\xxx\xxx_diff1d.py", line 79, in <module>
    losshistory, train_state = model.train(iterations=10)

I am referring to the example provided here https://github.com/lululxvi/deepxde/blob/b3921720188f0bf9183c3997dc60e5b26c17ae48/docs/demos/pinn_forward/eulerbeam.rst

[DeepaMahm]

Hi @praksharma ,

When you have a chance to look into this, could you please let me know?

[praksharma]

Sorry what was the error? Can you please paste the full error trace? It looks like you are using tensorflow.

[DeepaMahm]

Hello @praksharma ,

Thank you for your reply. Yes, I am using tensorflow.

Below is the complete code:

"""
Solving 1D concentration diffusion problem

References:
   Boundary condition:  https://github.com/lululxvi/deepxde/blob/b3921720188f0bf9183c3997dc60e5b26c17ae48/docs/demos/pinn_forward/eulerbeam.rst
   Initial condition: https://github.com/lululxvi/deepxde/blob/b3921720188f0bf9183c3997dc60e5b26c17ae48/docs/demos/pinn_forward/ode.2nd.rst#L11
"""
import deepxde as dde
import numpy as np

def pde(x, y):
    """
    :param x: 2-dimensional vector; x[:, 0:1] = x-coordinate; x[:, 1:] = t-coordinate
    :param y: network output; i.e., solution y(x,t)
    :return:
    """
    dc_t = dde.grad.jacobian(y, x, i=0, j=1)
    dc_xx = dde.grad.hessian(y, x, i=0, j=0)
    d = 1
    return (
        dc_t
        - d * dc_xx
    )

def boundary_l(x, on_boundary):
    """
    Because of rounding-off errors, dde.utils.isclose(x[0], 0)
    :param x:
    :param on_boundary: all boundary points
    :return:
    """
    return on_boundary and dde.utils.isclose(x[0], 0)  # définit le boundary x=0, left

def boundary_r(x, on_boundary):
    """
    Because of rounding-off errors, dde.utils.isclose(x[0], 0)
    :param x:
    :param on_boundary: all boundary points
    :return:
    """
    return on_boundary and dde.utils.isclose(x[0], 1)  # définit le boundary x=L, right

if __name__ == '__main__':
    geom = dde.geometry.Interval(0, 1)
    timedomain = dde.geometry.TimeDomain(0, 1)
    geomtime = dde.geometry.GeometryXTime(geom, timedomain)

    # bc_l = dde.icbc.DirichletBC(geom, lambda x: 5, lambda _, on_boundary: boundary_l)
    # bc_r = dde.icbc.NeumannBC(geom, lambda x: 0, lambda _, on_boundary: boundary_r)

    bc_l = dde.icbc.DirichletBC(geomtime, lambda x: 5, boundary_l)
    bc_r = dde.icbc.NeumannBC(geom, lambda x: 0, boundary_r)

    ic = dde.icbc.IC(geomtime, lambda x: 1, lambda _, on_initial: on_initial)

    data = dde.data.TimePDE(
        geomtime,
        pde,
        [bc_l, bc_r, ic],
        num_domain=40,  # no. of training residual points sampled inside the domain
        num_boundary=20,  # no. of training points sampled on the boundary
        num_initial=10,  # initial residual points for initial conditions
        # solution=func, # reference solution to compute error of the dde solution, can be ignored if no exact soln
        num_test=10000,  # points for testing the residual
    )

    layer_size = [2] + [32] * 3 + [1]  # width 32, depth 3
    activation = "tanh"
    initializer = "Glorot uniform"
    net = dde.nn.FNN(layer_size, activation, initializer)

    model = dde.Model(data, net)
    model.compile("adam", lr=0.001, metrics=["l2 relative error"])
    losshistory, train_state = model.train(iterations=10)

    dde.saveplot(losshistory, train_state, issave=True, isplot=True)

Could you please try running?

Error:

WARNING:tensorflow:AutoGraph could not transform <function at 0x00000255A562E520> and will run it as-is. Cause: could not parse the source code of <function at 0x00000255A562E520>: no matching AST found among candidates:

coding=utf-8

lambda x, on: np.array([on_boundary(x[i], on[i]) for i in range(len(x))]) To silence this warning, decorate the function with @tf.autograph.experimental.do_not_convert WARNING:tensorflow:AutoGraph could not transform <function at 0x00000255A562E7A0> and will run it as-is. Cause: could not parse the source code of <function at 0x00000255A562E7A0>: no matching AST found among candidates:

coding=utf-8

lambda x, on: np.array([on_boundary(x[i], on[i]) for i in range(len(x))]) To silence this warning, decorate the function with @tf.autograph.experimental.do_not_convert Traceback (most recent call last): File "I:\c2s\simgraph\simdeep_diff1d.py", line 80, in losshistory, train_state = model.train(iterations=10) ^^^^^^^^^^^^^^^^^^^^^^^^^^ File "I:\IISC-Project\venv\Lib\site-packages\deepxde\utils\internal.py", line 22, in wrapper result = f(*args, **kwargs) ^^^^^^^^^^^^^^^^^^ File "I:\xxx\venv\Lib\site-packages\deepxde\model.py", line 628, in train self._test() File "I:\xxx\venv\Lib\site-packages\deepxde\model.py", line 837, in _test self.train_state.metrics_test = [ ^ File "I:\xxx\venv\Lib\site-packages\deepxde\model.py", line 838, in m(self.train_state.y_test, self.train_state.y_pred_test) File "I:\xxx\venv\Lib\site-packages\deepxde\metrics.py", line 12, in l2_relative_error return np.linalg.norm(y_true - y_pred) / np.linalg.norm(y_true)

TypeError: unsupported operand type(s) for -: 'NoneType' and 'float'