lululxvi / deepxde

A library for scientific machine learning and physics-informed learning
https://deepxde.readthedocs.io
GNU Lesser General Public License v2.1
2.64k stars 741 forks source link

Solution of the Helmholtz equation for a plane wave #649

Open alabaykazakh opened 2 years ago

alabaykazakh commented 2 years ago

Dear Lu Lu,

Thank you for such a great package. I really consider deepXDE as a future working tool for my research interests. Here is my specific problem: image

  1. First I tried to train the network for a plane wave propagating in z-direction through the single domain (uniform eps1 distribution). The results are quite satisfactory: image Then I tried to get more periods by enlarging the domain (increasing the z range) or using larger eps. Unfortunately, I always obtain the non-converging result. How could I solve this problem? Should I modify the number of training points?

  2. Then I tried to solve PDE for the domain consisting of two materials (see pic). image I have read the similar issues and implemented tf.where for eps(x,z) to split the domain into the right and left parts, corresponding to different materials. However, the network could not solve this problem. Could you please check my implementation? Here is the part of my code where I tried to split the domain:

    
    def solution(x):
    kz = np.zeros_like(x[:,1:2])
    kz = np.where(x[:,1:2] > 2.0, k0 * np.sqrt(eps2), k0 * np.sqrt(eps1))
    return E0*np.cos(kz*x[:,1:2]) 

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)

mask = tf.where(x[:, 1:2] > 2, x[:, 1:2] *0 + eps2 * k0 ** 2, x[:, 1:2] * 0 + eps1 * k0 ** 2) 
return -dy_xx - dy_yy  - mask * y    

k0 = np.pi # wavenumber in vacuum (spatial frequency) E0 = 1 # amplitude of wave kx = 0 # spatial frequency in x-direction w = 0 # frequency t = 0 # time eps1 = 1 # permittivity of vacuum eps2 = 4 # permittivity of dielectric

Here is my code where I defined boundary conditions:

x = [-1,1] , z = [0,3]

geom = dde.geometry.Rectangle(xmin=[-1, 0], xmax=[1, 3])

Starting point of wave propagation (z=0)

bc_l = dde.icbc.DirichletBC( geom, lambda x: np.ones((len(x[:, 0:1]),1)) * E0, lambda x, on_boundary: on_boundary and np.isclose(x[1], 0), )

Periodic boundary condition (down)

bc_d = dde.icbc.DirichletBC( geom, lambda x: solution(x), lambda x, on_boundary: on_boundary and np.isclose(x[0], -1), )

Periodic boundary condition (up)

bc_u = dde.icbc.DirichletBC( geom, lambda x: solution(x), lambda x, on_boundary: on_boundary and np.isclose(x[0], 1), )


Other possible solutions from my side to add _eps(x,z)_ to my problem: 
a) use _eps(x,z)_ as third input. However, the problem is that I need a third dimension in this case, as I know.
b) implement subdomain inside the current domain and define extra boundary conditions for it. However, I plan to work with multilayer structures in the future, so the implementation of such many boundary conditions would be cumbersome.
I would be glad to get your opinion and feedback on my problem. Thank you in advance.
lululxvi commented 2 years ago
  1. Yes, try more training points, a larger network, and more iterations, etc.
  2. Not sure what is the problem. But for PBC, you may need to also add the first-order derivative.
tuzu123 commented 2 years ago

Hi @alabaykazakh I also tried to implement the multilayer problem for the wave equation in time domain, but the result did not converge. Maybe your solution b) will work (though the BC is cumbersome), but I don't know how to implement the subdomain now. Do you implement the subdomain solution?

Hoping for your reply and wish to keep in touch with you. Thanks in advance.

alabaykazakh commented 2 years ago

Hi @alabaykazakh I also tried to implement the multilayer problem for the wave equation in time domain, but the result did not converge. Maybe your solution b) will work (though the BC is cumbersome), but I don't know how to implement the subdomain now. Do you implement the subdomain solution?

Hoping for your reply and wish to keep in touch with you. Thanks in advance.

I rewrote my code and completely got rid of the time domain. The domain is now represented by two rectangles united using | operation. It seems to work for me.

# Vacuum domain:
x_min1  =  -1.0
x_max1 =   1.0
z_min1  =   0.0
z_max1  =   2.0

# Dielectric domain:
x_min2 = -1.0
x_max2 =  1.0
z_min2 =   2.0
z_max2 =   3.0

# Vacuum geometry:
geom1 = dde.geometry.Rectangle(xmin=[x_min1, z_min1], xmax=[x_max1, z_max1])
# Dielectric geometry:
geom2 = dde.geometry.Rectangle(xmin=[x_min2, z_min2], xmax=[x_max2, z_max2])
# Final domain:
geom_final = geom1 | geom2

The eps properties you can define inside def pde(x, y). Feel free to ask more questions. I am open for any discussion.

tuzu123 commented 2 years ago

Thanks for your quick reply!

  1. Is the eps properties still define inside def pde(x,y) using tf.where() or any other forms?
  2. How to define the connection condition in the subdomain interface z = 2? My wave equation is for the elastic plane wave and I don't know how to convert to the frequency domain at present because of my poor physical knowledge.
alabaykazakh commented 2 years ago
  1. Yes, tf.where() works well to define required properties in different parts of the domain.
  2. For such a simple case, I didn't define a special condition at z = 2 and could get the solution below. image

At the moment I am stuck trying to solve the Helmholz equation for a more complicated case "dielectric hanging in a vacuum".

I think the dielectric domain requires extra boundary conditions, which, in my opinion, should be:

Did someone already manage to correctly define boundary conditions for the inner domain or solve a similar problem? I would be glad to get some pieces of advice, @lululxvi, please.

lululxvi commented 2 years ago

Yes, it is doable. Check the paper https://doi.org/10.1364/OE.384875

alabaykazakh commented 2 years ago

Dear @lululxvi, thank you again for deepXDE and the paper.

As I saw from your paper on metamaterials design, you also implemented a plane wave to illuminate an array. Could you please give me a hint on how I could implement oblique plane wave illumination (for example 45°)? I think I should somehow give the network information about propagation direction in the derivative part of PeriodicBC. However, since I tried first the classic PeriodicBC (see pic), where no changes on the border,

pbc1 = dde.icbc.PeriodicBC(geom_final, 0, per_boundary, derivative_order=0, component=0)
pbc2 = dde.icbc.PeriodicBC(geom_final, 0, per_boundary, derivative_order=1, component=0)

the network returns the symmetric output.

lululxvi commented 2 years ago

Not sure what you are doing. PeriodicBC enforces the Periodic BC, not directly related to the plane wave.

lshil00 commented 2 years ago
  1. Yes, tf.where() works well to define required properties in different parts of the domain.
  2. For such a simple case, I didn't define a special condition at z = 2 and could get the solution below. image

At the moment I am stuck trying to solve the Helmholz equation for a more complicated case "dielectric hanging in a vacuum".

I think the dielectric domain requires extra boundary conditions, which, in my opinion, should be:

Did someone already manage to correctly define boundary conditions for the inner domain or solve a similar problem? I would be glad to get some pieces of advice, @lululxvi, please.

sorry to bother you, but every time when I use tf.where() in my pde function, L-BFGS will stop much earlier even though it did not converge. Have you met this problem before? Thanks for your help!

alabaykazakh commented 2 years ago
  1. Yes, tf.where() works well to define required properties in different parts of the domain.
  2. For such a simple case, I didn't define a special condition at z = 2 and could get the solution below. image

At the moment I am stuck trying to solve the Helmholz equation for a more complicated case "dielectric hanging in a vacuum". I think the dielectric domain requires extra boundary conditions, which, in my opinion, should be:

Did someone already manage to correctly define boundary conditions for the inner domain or solve a similar problem? I would be glad to get some pieces of advice, @lululxvi, please.

sorry to bother you, but every time when I use tf.where() in my pde function, L-BFGS will stop much earlier even though it did not converge. Have you met this problem before? Thanks for your help!

tf.where() just defines the distribution of material properties over the domain, so it probably should not affect an early stop of the optimizer. Can you show the corresponding piece of your code, namely how the L-BFGS and tf.where() parts are implemented?

alabaykazakh commented 1 year ago

Dear @lululxvi , Talking about the same problem described above, I want to make the geometry of the scatterer more complex. Instead of a simple rectangle, I would like to describe it as a polygon with four vertices. This way, eps(x,z) becomes a function of 9 input parameters [eps2, x1, z1, x2, z2, x3, z3, x4, z4] defining the material and geometry distribution for the scatterer. Since I have to solve a parametric PDE now, I started with an investigation of possible approaches. I have read the issues to related problems opened by other users. To my best knowledge, I have two options:

I. Increase the dimensionality of the input dataset from 2D to 11D (using ND Hypercube). This is the actual approach I use at the moment and it looks working. However, the higher the dimensionality of the input, the smaller the concentration of points in every single layer (x,z). Increasing the total number of collocation points helped to some extent till I reached computational limits (~80000 points). Q1: Is an increase in the total number of collocation points according to dimension the best way to improve accuracy prediction, in general? Should I also do additional steps (e.g., modify a network size)? Q2: In the case huge amount of training points, I would like to split them into batches and feed them to FNN portion-wise. This way I can avoid exceeding the limits of GPU memory. However, I did not find a straightforward application of batches in DeepXDE. I tried to define batch_size in deepxde.model.train(), but the inputs tensor size remains the same during training. Some users also advised using PDEPointResampler(), but I think it just shuffles the points inside the same single huge batch during training steps. Could you please suggest how I reduce the batch size of my network?

II. DeepONet/ Physics-informed DeepONet. In such a case, we stay in the 2D domain and combine DeepONet with PINNs to train on parameter variations. However, as far as I know, DeepONet is not unsupervised anymore and requires input-output pairs to get trained. Meanwhile, physics-informed DeepONet is less supervised but still is data-driven. My intention is to find an architecture that could learn operators for the given problem of light scattering, i.e. mapping of eps(x,z) to E(x,z). Q3: Is there any possibility to use DeepONet without ground truth data, so to keep it purely physics-informed? Otherwise, PINN loses its advantage in unsupervised training.

Thank you!

lululxvi commented 1 year ago

Q1: Is an increase in the total number of collocation points according to dimension the best way to improve accuracy prediction, in general? Should I also do additional steps (e.g., modify a network size)?

Yes, you may also need a larger network.

Q2: Could you please suggest how I reduce the batch size of my network?

Use PDEPointResampler

II. Q3: Is there any possibility to use DeepONet without ground truth data, so to keep it purely physics-informed? Otherwise, PINN loses its advantage in unsupervised training.

You can use physics-informed DeepONet. DeepXDE supports it, but we don't have public examples yet.

Wang946992 commented 1 year ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

alabaykazakh commented 1 year ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

Hi,

Please check out our recent publication https://ieeexplore.ieee.org/abstract/document/10221390. We have solved the scattered light computation problem for a hanging absorber of different materials using the DeepXDE framework.

lululxvi commented 1 year ago

@alabaykazakh Feel free to add your paper to https://deepxde.readthedocs.io/en/latest/user/research.html

akokoro commented 1 week ago

Yes, it is doable. Check the paper https://doi.org/10.1364/OE.384875

Dear @lululxvi I have read the paper you recommended. It describes the reverse design of metamaterials. There are a lot of boundary condition Settings that I'm interested in. Where can I see the code for this paper? I want to learn more about the conditions for this electromagnetic case.

akokoro commented 1 week ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

Hi,

Please check out our recent publication https://ieeexplore.ieee.org/abstract/document/10221390. We have solved the scattered light computation problem for a hanging absorber of different materials using the DeepXDE framework.

@alabaykazakh,Hi I have to say I've learned a lot from your work. I have finished studying your paper. But there's a problem. In the above two comments you said that the suspended target boundary should also have a definition of the boundary condition, but I saw in the paper that your setup is just a distribution of two different dielectric constants, and there is no definition of the boundary condition.

alabaykazakh commented 1 week ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

Hi, Please check out our recent publication https://ieeexplore.ieee.org/abstract/document/10221390. We have solved the scattered light computation problem for a hanging absorber of different materials using the DeepXDE framework.

@alabaykazakh,Hi I have to say I've learned a lot from your work. I have finished studying your paper. But there's a problem. In the above two comments you said that the suspended target boundary should also have a definition of the boundary condition, but I saw in the paper that your setup is just a distribution of two different dielectric constants, and there is no definition of the boundary condition.

Dear @akokoro, glad to hear about your interest in our paper. I guess you mean the internal boundary conditions/interface boundary conditions. You're right, indeed we did not implement them in our paper since the material contrast between the air and materials in the EUV spectral range is quite small. In reality, sharp material transition or abrupt changes in material properties can lead to bad regularity and numerical inaccuracies, even for numerical solvers. To handle the discontinuity on the interface and improve the numerical convergence, different techniques are applied such as smoothing schemes or the aforementioned interface boundary conditions.

akokoro commented 1 week ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

Hi, Please check out our recent publication https://ieeexplore.ieee.org/abstract/document/10221390. We have solved the scattered light computation problem for a hanging absorber of different materials using the DeepXDE framework.

@alabaykazakh,Hi I have to say I've learned a lot from your work. I have finished studying your paper. But there's a problem. In the above two comments you said that the suspended target boundary should also have a definition of the boundary condition, but I saw in the paper that your setup is just a distribution of two different dielectric constants, and there is no definition of the boundary condition.

Dear @akokoro, glad to hear about your interest in our paper. I guess you mean the internal boundary conditions/interface boundary conditions. You're right, indeed we did not implement them in our paper since the material contrast between the air and materials in the EUV spectral range is quite small. In reality, sharp material transition or abrupt changes in material properties can lead to bad regularity and numerical inaccuracies, even for numerical solvers. To handle the discontinuity on the interface and improve the numerical convergence, different techniques are applied such as smoothing schemes or the aforementioned interface boundary conditions.

@alabaykazakh ,Hi Thank you for your clarification. I understand that in your paper you show a coordinate point of 91800. You should not be using the deepxde framework, you should be writing your own PINN network. Is it convenient to learn your code? Especially PML and periodic boundary conditions.

alabaykazakh commented 1 week ago

@alabaykazakh Hi, In"dielectric hanging in a vacuum promblem".Now, have you solved this problem about how define partial differential equations? Thanks

Hi, Please check out our recent publication https://ieeexplore.ieee.org/abstract/document/10221390. We have solved the scattered light computation problem for a hanging absorber of different materials using the DeepXDE framework.

@alabaykazakh,Hi I have to say I've learned a lot from your work. I have finished studying your paper. But there's a problem. In the above two comments you said that the suspended target boundary should also have a definition of the boundary condition, but I saw in the paper that your setup is just a distribution of two different dielectric constants, and there is no definition of the boundary condition.

Dear @akokoro, glad to hear about your interest in our paper. I guess you mean the internal boundary conditions/interface boundary conditions. You're right, indeed we did not implement them in our paper since the material contrast between the air and materials in the EUV spectral range is quite small. In reality, sharp material transition or abrupt changes in material properties can lead to bad regularity and numerical inaccuracies, even for numerical solvers. To handle the discontinuity on the interface and improve the numerical convergence, different techniques are applied such as smoothing schemes or the aforementioned interface boundary conditions.

@alabaykazakh ,Hi Thank you for your clarification. I understand that in your paper you show a coordinate point of 91800. You should not be using the deepxde framework, you should be writing your own PINN network. Is it convenient to learn your code? Especially PML and periodic boundary conditions.

91800 is a number of collocation points we sampled across the computational domain using the Sobol sampling sequence. As mentioned in my comments above, we used the DeepXDE framework to solve a given forward 2D scattering problem. We utilized PML and periodic boundary conditions from Physics-informed neural networks with hard constraints for inverse design. Please feel free to check out other recent publications on using convolution-based PINNs for 3D modeling of light diffraction from more complex nano-optical devices.: 3D mask simulation and lithographic imaging using physics-informed neural networks 3D EUV mask simulator based on physics-informed neural networks: effects of polarization and illumination