NeuralPDE to simulate Burguer 1D - not a satisfactorial solution

Specific Domains Machine Learning
1

I’m writting a research in which I heavily used NeuralPDE.jl, and I would like to showcase solution vs analytical, or classical simulations vs provided by the package.

In one of these instances, I’m getting an answer I don’t understand. Because, I believe it’s wrong. The solver is not obeying my initial condition.

eq: ∂u/∂t + u(x,t)(∂u/∂x) = 0,
bc: u(0,t) = (1/(2√π))exp((-1/2)(x-3)^2)

f(x) = (1/(2*√π))*exp((-1/2)*(x-3)^2)

  using NeuralPDE, Flux, ModelingToolkit, GalacticOptim, Optim, DiffEqFlux
  import ModelingToolkit: Interval, infimum, supremum

  @parameters t, x
  @variables u(..)
  Dt = Differential(t)
  Dx = Differential(x)

  #2D PDE
  eq  = Dt(u(t,x)) + u(t,x)*Dx(u(t,x)) ~ 0

  # Initial and boundary conditions
  bcs = [u(0,x) ~ f(x),
         u(0,x) - f(x) ~ 0]

  # Space and time domains
  domains = [t ∈ Interval(0.0,20.0),
             x ∈ Interval(0.0,10.0)]
  # Discretization
  nx=100;
  dx = 10/(nx)
  endTime = 20
  nt = 100          
  dt = endTime/nt

The neural network:

  # Neural network
  dim = 2 # number of dimensions
  chain = FastChain(FastDense(dim,25,Flux.σ),FastDense(25,25,Flux.σ),FastDense(25,1))
  # Initial parameters of Neural network
  initθ = Float64.(DiffEqFlux.initial_params(chain))

  # Discretization
  discretization = PhysicsInformedNN(chain,GridTraining(dx),init_params=initθ)

  @named pde_system = PDESystem(eq,bcs,domains,[x,t],[u(x, t)])
  prob = discretize(pde_system,discretization)

  cb = function (p,l)
      println("Current loss is: $l")
      return false
  end

  opt = Optim.BFGS()
  res = GalacticOptim.solve(prob,opt; cb = cb, maxiters=5000)
  phi = discretization.phi

Finally, the plot,

  using Plots; pyplot();
 
  xs,ts = [infimum(d.domain):dx:supremum(d.domain) for (d,dx) in zip(domains,[dx,dt])]
  u_predict = reshape([Array(phi([x, t], res.minimizer))[1] for x in xs for t in ts], length(xs), length(ts))
  
  plot(xs,ts,u_predict',st=:surface, title="Burguer equation, PINN", xlabel="X", ylabel="Time", zlabel="U")

which gives me:

Figure_31
Figure_311854×986 403 KB

And, the initial condition is pretty different from what I specified:

p1=plot(xs, u_predict[1,:], title="Predicted")
p2=plot(xs, map(x -> f(x), xs), title="Analytical")
plot(p1,p2)

By the Euler Explicit method,

nxx=100;
delta_xx = 10/(nxx - 1)
xx = range(0, stop=delta_xx*(nxx-1), length=nxx)

endTime = 20   
nt = 1000          
delta_t = endTime/nt  
t = range(0, stop=endTime, length=nt) 

f(x) = (1/(2*√π))*exp((-1/2)*(x-3)^2)
v_zero = f.(xx)  
 
v=zeros((nxx,nt+1))
v[:,1]=copy(v_zero)

for n in 1:nt       
    v[:,n+1] = copy(v[:,n]) 
    for i in 2:nxx 
        v[i,n+1] = v[i,n] - v[i,n] * delta_t/delta_xx * (v[i,n] - v[i-1,n])
    end
end

using Plots; pyplot()

xxs = collect(xx)
ts = collect(t)

plot(collect(xx),collect(t),v'[1:1000,1:100],st=:surface, title="Burguer equation", xlabel="X", ylabel="Time", zlabel="V")
2 Likes
2

Initial condition comparison:

Init_1
Init_11854×986 55.2 KB

3

Euler Explicit solution plot

Burguer_Euler_1
Burguer_Euler_11854×986 358 KB

4

In case anyone needs to verify the Euler Equation,
pic-selected-220126-1016-23

5

Don’t know anything about PDEs, but some nice surface plots you got there :wink:

6

Thank you @nilshg, our suites provided by Julia are truly amazing :slight_smile:.

EDIT:
I’m editing this comment, because I’m unable to respond due to the maximum replies achieved by a newcomer in a day.

But, the results are still spurious using ADAM beforehand and, also, when utilizing only ADAM.

Result ADAM → BFGS:

res = GalacticOptim.solve(prob,ADAM(10^-3); cb = cb, maxiters=300)
prob = remake(prob,u0=res.minimizer)
opt = Optim.BFGS()
res = GalacticOptim.solve(prob,opt; cb = cb, maxiters=500)
phi = discretization.phi

Figure_21
Figure_211854×986 298 KB

7

It’s not recommended to start with BFGS. That will easily hit local minima. The NeuralPDE paper has a whole section dedicated to that:

Normally 300 iterations of ADAM helps achieve a better minima.

2 Likes