Search code examples
juliapde

No Method Matching Neural PDE error in julia


I am trying to run the code that solves physics informed partial differential equation particularly poison equation on JULIA. The code shows this, 'no method matching' error whenever, I try to run it. Please help me resolve it.

using NeuralPDE, Flux, ModelingToolkit, GalacticOptim, Optim, DiffEqFlux

@parameters x y
@variables u(..)
Dxx = Differential(x)^2
Dyy = Differential(y)^2

# 2D PDE
eq  = Dxx(u(x,y)) + Dyy(u(x,y)) ~ -sin(pi*x)*sin(pi*y)

# Boundary conditions
bcs = [u(0,y) ~ 0.f0, u(1,y) ~ -sin(pi*1)*sin(pi*y),
       u(x,0) ~ 0.f0, u(x,1) ~ -sin(pi*x)*sin(pi*1)]
# Space and time domains
domains = [x ∈ IntervalDomain(0.0,1.0),
           y ∈ IntervalDomain(0.0,1.0)]

# Neural network
dim = 2 # number of dimensions
chain = FastChain(FastDense(dim,16,Flux.σ),FastDense(16,16,Flux.σ),FastDense(16,1))

# Discretization
dx = 0.05
discretization = PhysicsInformedNN(chain,GridTraining(dx))

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

#Optimizer
opt = Optim.BFGS()

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

res = GalacticOptim.solve(prob, opt, cb = cb, maxiters=1000)
phi = discretization.phi

using Plots

xs,ys = [domain.domain.lower:dx/10:domain.domain.upper for domain in domains]
analytic_sol_func(x,y) = (sin(pi*x)*sin(pi*y))/(2pi^2)

u_predict = reshape([first(phi([x,y],res.minimizer)) for x in xs for y in ys],(length(xs),length(ys)))
u_real = reshape([analytic_sol_func(x,y) for x in xs for y in ys], (length(xs),length(ys)))
diff_u = abs.(u_predict .- u_real)

p1 = plot(xs, ys, u_real, linetype=:contourf,title = "analytic");
p2 = plot(xs, ys, u_predict, linetype=:contourf,title = "predict");
p3 = plot(xs, ys, diff_u,linetype=:contourf,title = "error");
plot(p1,p2,p3)

The code uses uses Neural networks to solve the partial differential equation. This is the image of the error I am countering: This is the image of the error I am countering.

I am trying this for the first time, So i really don't understand what type of output would should be generated, But by basic understanding there must be different types of prediction and error plots. i have attached the images of the plots listed on the julia site: i have attached the images of the plots listed on the julia site.


Solution

  • You are not using a current documentation. Out of curiosity, could you share the page that you happened to find?

    If you follow the latest documentation, the documentation for v5.3 (which is surely the version that you are on), then it should work. See: https://docs.sciml.ai/NeuralPDE/v5.3/.

    https://docs.sciml.ai/NeuralPDE/v5.3/tutorials/pdesystem/

    using NeuralPDE, Lux, Optimization, OptimizationOptimJL
    import ModelingToolkit: Interval
    
    @parameters x y
    @variables u(..)
    Dxx = Differential(x)^2
    Dyy = Differential(y)^2
    
    # 2D PDE
    eq  = Dxx(u(x,y)) + Dyy(u(x,y)) ~ -sin(pi*x)*sin(pi*y)
    
    # Boundary conditions
    bcs = [u(0,y) ~ 0.0, u(1,y) ~ 0.0,
           u(x,0) ~ 0.0, u(x,1) ~ 0.0]
    # Space and time domains
    domains = [x ∈ Interval(0.0,1.0),
               y ∈ Interval(0.0,1.0)]
    
    # Neural network
    dim = 2 # number of dimensions
    chain = Lux.Chain(Dense(dim,16,Lux.σ),Dense(16,16,Lux.σ),Dense(16,1))
    
    # Discretization
    dx = 0.05
    discretization = PhysicsInformedNN(chain,GridTraining(dx))
    
    @named pde_system = PDESystem(eq,bcs,domains,[x,y],[u(x, y)])
    prob = discretize(pde_system,discretization)
    
    #Optimizer
    opt = OptimizationOptimJL.BFGS()
    
    #Callback function
    callback = function (p,l)
        println("Current loss is: $l")
        return false
    end
    
    res = Optimization.solve(prob, opt, callback = callback, maxiters=1000)
    phi = discretization.phi
    
    using Plots
    
    xs,ys = [infimum(d.domain):dx/10:supremum(d.domain) for d in domains]
    analytic_sol_func(x,y) = (sin(pi*x)*sin(pi*y))/(2pi^2)
    
    u_predict = reshape([first(phi([x,y],res.u)) for x in xs for y in ys],(length(xs),length(ys)))
    u_real = reshape([analytic_sol_func(x,y) for x in xs for y in ys], (length(xs),length(ys)))
    diff_u = abs.(u_predict .- u_real)
    
    p1 = plot(xs, ys, u_real, linetype=:contourf,title = "analytic");
    p2 = plot(xs, ys, u_predict, linetype=:contourf,title = "predict");
    p3 = plot(xs, ys, diff_u,linetype=:contourf,title = "error");
    plot(p1,p2,p3)
    

    Note: you should fully remove GalacticOptim from your system (]rm GalacticOptim) since that's a deprecated package.