I found three problems when I used deepxde to reproduce other people's papers.

[wongsihan]

Dear Dr. lu.
Hello, thank you very much for the deepxde framework, it works very well. But recently I had some problems with it when I was reproducing someone else's paper (Physics Informed Neural Networks for Electromagnetic Analysis - IEEE Transactions on Magnetics). There are three examples in this paper, and I encountered the problem when reproducing the third example: Multi-domain Magnetostatic. Its purpose is to give the initial current and find out its magnetic vector potential distribution. I ran into three problems.

1. When I set the vacuum permeability u0 to a very small value (0.000001256), the training loss often returns as nan, resulting in an incorrect output magnetic vector potential A. The code will be given below. However, when I give a larger vacuum permeability for u0, it returns the correct output magnetic vector potential, specifically, the vacuum permeability returns the correct distribution from 300 to 0.00003, but less than this order of magnitude will return nan during the first or second stage of the optimizer iteration, as well as the wrong A distribution. The code is as follows:

"""Backend: pytorch"""

import deepxde as dde
import matplotlib.pyplot as plt
import numpy as np
import torch,os,datetime,sys

def pde(x, y):
    dy_xx = dde.grad.hessian(y, x, i=0, j=0)
    dy_yy = dde.grad.hessian(y, x, i=1, j=1)
    dy1_x = dde.grad.jacobian(y, x, i=0, j=0)
    dy2_y = dde.grad.jacobian(y, x, i=0, j=1)

    dy1_x_dum = dy1_x.cpu().detach().numpy()
    dy2_y_dum = dy2_y.cpu().detach().numpy()
    device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")
    values = torch.Tensor(len(y), 1, device=device)

    u0 = 0.000001256
    for ind, ele in enumerate(x):
        ele_dum=ele.cpu().detach().numpy()
        if np.linalg.norm(ele_dum) <= 0.1:
            values[ind] = -dy_xx[ind] - dy_yy[ind] - 1
        else:
            u = u0 * ((5000 / (1 + 0.05 * (np.square(dy1_x_dum[ind]) + np.square(dy2_y_dum[ind])))) + 200)
            values[ind] = torch.from_numpy(-(1/u) * (dy1_x_dum[ind] + dy2_y_dum[ind]))
    return values

def boundary(x, on_boundary):
    return on_boundary and np.isclose(np.linalg.norm(x),1)

def func1(x):
    values = np.zeros((len(x),1))
    return values


geom1 = dde.geometry.geometry_2d.Disk((0,0), 0.1)
geom2 = dde.geometry.geometry_2d.Disk((0,0), 1)
geom = geom1 | geom2
bc = dde.icbc.DirichletBC(geom, func1, boundary)

data = dde.data.PDE(geom, pde, bc, num_domain=2800, num_boundary=500, num_test=3000)
net = dde.nn.FNN([2] + [50] * 4 + [1], "tanh", "Glorot uniform")
model = dde.Model(data, net)

model.compile("adam", lr=0.001)
model.train(epochs=15000)
model.compile("L-BFGS")
time = datetime.datetime.now()
time1_str = datetime.datetime.strftime(time,'%Y-%m-%d-%H-%M-%S')
filename = sys.argv[0][sys.argv[0].rfind(os.sep) + 1:]
losshistory, train_state = model.train(model_save_path="/root/data/wsh/deepxde/examples/pinn_forward/model/model-{time}-{filename}.ckpt".format(time=time1_str,filename=filename))

dde.saveplot(losshistory, train_state, issave=True, isplot=True,  loss_fname="{filename}-loss.dat".format(filename=filename), train_fname="{filename}-train.dat".format(filename=filename),test_fname="{filename}-test.dat".format(filename=filename))

When the value of loss is not nan, the output is close to that in the paper, as shown in Figure 1, but it still has problems with the order of magnitude, which is 100 times different from that in the original paper：

Fig. 1
When the vacuum permeability is too small, resulting in a return loss of nan, its output is erroneous, as shown in Figure 2, and the order of magnitude and distribution do not match the original.

Fig. 2

2. This example in the paper is two circles, the outer large circle is the solution domain and the inner small circle is the current source. When the large circle remains unchanged and the small circle is changed to a square, the output A distribution is wrong, but the returned loss is normal. Its code is almost the same as the previous problem, only the geometric structure and boundary conditions are modified. Its code is as follows:

"""Backend: pytorch"""
import deepxde as dde
import matplotlib.pyplot as plt
import numpy as np
import torch,os,datetime,sys
from decimal import Decimal

def pde(x, y):
    dy_xx = dde.grad.hessian(y, x, i=0, j=0)
    dy_yy = dde.grad.hessian(y, x, i=1, j=1)
    dy1_x = dde.grad.jacobian(y, x, i=0, j=0)
    dy2_x = dde.grad.jacobian(y, x, i=0, j=1)

    dy1_x_dum = dy1_x.cpu().detach().numpy()
    dy2_x_dum = dy2_x.cpu().detach().numpy()
    device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")
    values = torch.Tensor(len(y), 1, device=device)
    u0 = 0.000001256

    for ind, ele in enumerate(x):
        ele_dum=ele.cpu().detach().numpy()
        if ele_dum[0] >= -3 and ele_dum[0] <= 3 and ele_dum[1] >= -3 and ele_dum[1] <= -3:
            values[ind] = -dy_xx[ind] - dy_yy[ind] - 1
        else:
            u = u0 * ((5000 / (1 + 0.05 * (np.square(dy1_x_dum[ind]) + np.square(dy2_x_dum[ind])))) + 200)
            values[ind] = torch.from_numpy(-(1/u) * (dy1_x_dum[ind] + dy2_x_dum[ind]))

    return values


def boundary(x, on_boundary):
    return on_boundary and np.isclose(np.linalg.norm(x),10)

def func1(x):
    values = np.zeros((len(x),1))
    return values


geom1 = dde.geometry.geometry_2d.Rectangle((-3,-3), (3,3))
geom2 = dde.geometry.geometry_2d.Disk((0,0), 10)
geom = geom1 | geom2
bc = dde.icbc.DirichletBC(geom, func1, boundary)

data = dde.data.PDE(geom, pde, bc, num_domain=2000, num_boundary=500, num_test=3000)
net = dde.nn.FNN([2] + [50] * 4 + [1], "tanh", "Glorot uniform")
model = dde.Model(data, net)

model.compile("adam", lr=0.001)
model.train(epochs=15000)
model.compile("L-BFGS")
time = datetime.datetime.now()
time1_str = datetime.datetime.strftime(time,'%Y-%m-%d-%H-%M-%S')
filename = sys.argv[0][sys.argv[0].rfind(os.sep) + 1:]
losshistory, train_state = model.train(model_save_path="/root/data/wsh/deepxde/examples/pinn_forward/model/model-{time}-{filename}.ckpt".format(time=time1_str,filename=filename))

Its output is shown in Figure 3, which I think is wrong because just changing a circle to a square should not have this change, I think the correct output should be close to when it is a circle and should not have a negative number.

Fig.3

3. The above examples are two-dimensional planes, I want to predict and analyze the distribution in three-dimensional space, so I want to draw three-dimensional cylinders and I want to insert some one-dimensional currents, so please ask if deepxde has built-in functions to build three-dimensional cylinders, I just found cubes and spheres. Also, is it possible to insert one-dimensional currents into a three-dimensional cylinder for solving? When I use the 1D geometry function, I don't find the corresponding settings for the 3D coordinates.

Thank you again very much for the deepxde framework! And I wish you all the best!

[lululxvi]

One possible issue of your code: The following code will make some trouble for the back-propgragation.

    for ind, ele in enumerate(x):
        ele_dum=ele.cpu().detach().numpy()

[mengxia123]

Hello brother, I am currently reproducing the case you mentioned in the paper. Unfortunately, in the second case, I will not handle different media interfaces. Do you know if you have any good suggestions or could you take a look at your code. It's really presumptuous. Good luck