Create a network geometry

[DeepaMahm]

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.

[lululxvi]

Your approach seems good. You don't need to define lines.

[DeepaMahm]

Hi Prof. @lululxvi ,

Thanks so much for the clarification and sorry for the late reply.
Kindly excuse me for another question.

I tried to set up the problem for 1D diffusion with Dirichlet and Neumann boundary conditions.

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 don't know how to fix this error.

Could you please have a look at 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)

Alternatively, I also tried

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

But I still have some errors and am not sure how to fix them.

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)

When I want to try for 2D geometry like the below (geometry generated by specifying point cloud)

I am unsure how to specify the points for the two ends for the geometry like above. On one end , I would like to specify the Dirichlet and Neumann on the other end.

For 1D, I have the function:

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)

[lululxvi]

Delete metrics=["l2 relative error"]

[DeepaMahm]

Thanks a lot Prof.,

I tried with the change you had suggested.

Please see the code below:

"""
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)
    losshistory, train_state = model.train(iterations=10)

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

Unfortunately, I still couldn't run it successfully

Traceback (most recent call last):
File "I:\xx\xxx\simdeep_diff1d.py", line 82, in
dde.saveplot(losshistory, train_state, issave=True, isplot=True)
File "I:\xxx\venv\Lib\site-packages\deepxde\utils\external.py", line 175, in saveplot
save_best_state(train_state, train_fname, test_fname)
File "I:\xxx\venv\Lib\site-packages\deepxde\utils\external.py", line 344, in save_best_state
test = np.hstack((train_state.X_test, best_y))
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "<array_function internals>", line 200, in hstack
File "I:\xxxx\venv\Lib\site-packages\numpy\core\shape_base.py", line 370, in hstack
return _nx.concatenate(arrs, 1, dtype=dtype, casting=casting)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "<array_function internals>", line 200, in concatenate
ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 2 dimension(s) and the array at index 1 has 1 dimension(s)

Could you please have a look?

layer_size = [2] + [32] * 3 + [1]
I have [2] (for time and space) and [1] for concentration.

Could you please also look at my initial and boundary conditions?
I am not sure if any of the syntax is wrong.

[lululxvi]

Strange. Not sure why. You can just delete the last line.

[DeepaMahm]

Prof.,
Thank you so much for the advise.

I excluded the last line dde.saveplot(losshistory, train_state, issave=True, isplot=True).

Unfortunately, we are still facing issues. The training stops after the first iteration and the loss is inf.

I:\xxxx\venv\Lib\site-packages\keras\src\initializers\initializers.py:120: UserWarning: The initializer GlorotUniform is unseeded and being called multiple times, which will return identical values each time (even if the initializer is unseeded). Please update your code to provide a seed to the initializer, or avoid using the same initializer instance more than once.
  warnings.warn(
WARNING:tensorflow:AutoGraph could not transform <function <lambda> at 0x000002A75DBAE520> and will run it as-is.
Cause: could not parse the source code of <function <lambda> at 0x000002A75DBAE520>: 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 <lambda> at 0x000002A75DBAE7A0> and will run it as-is.
Cause: could not parse the source code of <function <lambda> at 0x000002A75DBAE7A0>: 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
Step      Train loss                                  Test loss                                   Test metric
0         [2.71e-02, 2.46e+01, nan, 6.92e-01]         [2.69e-02, 2.46e+01, nan, 6.92e-01]         []  

Best model at step 0:
  train loss: inf
  test loss: inf
  test metric: 

'train' took 7.742125 s

I notice again the problem is with the lambda function defined in the boundary/initial condition.

I would like to share the python file (1D diffusion) for your kind reference. Could you please have a look it when it is possible for you?

[lululxvi]

There is nan in the loss. You need to fix this.