# NeuralPDE to simulate Burguer 1D - not a satisfactorial solution

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:

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

Initial condition comparison:

Euler Explicit solution plot

In case anyone needs to verify the Euler Equation,

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

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

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.

``````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
``````

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