NeuralPDE.jl icon indicating copy to clipboard operation
NeuralPDE.jl copied to clipboard

[WIP] Generate set of points for each variable

Open KirillZubov opened this issue 2 years ago • 2 comments

https://github.com/SciML/NeuralPDE.jl/issues/469

KirillZubov avatar Dec 29 '21 15:12 KirillZubov

@parameters t, x
@variables u(..)
∂t, ∂x = Differential(t), Differential(x)
eq = ∂t(u(t, x)) + ∂x(u(t, x)) ~ 0
domains = [t ∈ Interval(0.0, 1.0), x ∈ Interval(0.0, 1.0)]
bcs = [u(t, 0) ~ u(t, 1), u(0, x) ~ sin(2π * x)]
dx = 0.05
chain = FastChain(FastDense(2, 16, Flux.tanh), FastDense(16, 16, Flux.tanh), FastDense(16, 1))
initθ = Float64.(DiffEqFlux.initial_params(chain))
strategy = NeuralPDE.GridTraining(dx)
# strategy = NeuralPDE.QuasiRandomTraining(100)
# quadrature_strategy = NeuralPDE.QuadratureTraining(quadrature_alg=CubatureJLh(),
#                                                     reltol=1e-3,abstol=1e-3,
#                                                     maxiters =50, batch=100)
discretization = NeuralPDE.PhysicsInformedNN(chain,strategy ; init_params = initθ)
@named pde_system = PDESystem(eq, bcs, domains, [t, x], [u(t, x)])
prob = NeuralPDE.discretize(pde_system, discretization)
sym_prob = NeuralPDE.symbolic_discretize(pde_system, discretization)

prob.f(initθ, nothing)

pde_inner_loss_functions = prob.f.f.loss_function.pde_loss_function.pde_loss_functions.contents
bcs_inner_loss_functions = prob.f.f.loss_function.bcs_loss_function.bc_loss_functions.contents

s1 = prob.f.f.loss_function.bcs_loss_function.bc_loss_functions.contents[1].train_set
s2 = prob.f.f.loss_function.bcs_loss_function.bc_loss_functions.contents[2].train_set

cb_ = function (p,l)
    println("loss: ", l )
    println("pde_losses: ", map(l_ -> l_(p), pde_inner_loss_functions))
    println("bcs_losses: ", map(l_ -> l_(p), bcs_inner_loss_functions))
    return false
end
opt = Optim.BFGS()
res = GalacticOptim.solve(prob, opt; cb = cb_, maxiters = 1200)
phi = discretization.phi
xs = 0.0:0.01:1.0
sol0 = [phi([0, x], res.u)[1] for x in xs]
sol1 = [phi([0.5, x], res.u)[1] for x in xs]
sol2 = [phi([0.75, x], res.u)[1] for x in xs]
plot(xs, sol0, label = "t=0")
plot!(xs, sol1,label = "t=0.5")

periodic

KirillZubov avatar Dec 29 '21 15:12 KirillZubov

@parameters t, x
@variables u(..)
∂t, ∂x = Differential(t), Differential(x)
eq = ∂t(u(t, x)) + ∂x(u(t, x)) ~ 0
domains = [t ∈ Interval(0.0, 1.0), x ∈ Interval(0.0, 1.0)]
bcs = [u(t, 0) ~ sin(-2π * t),
       u(t, 1) ~ u(t, 0),
       u(1, 1) * u(t, 1) ~ + u(t, 0) - u(t, x)]
dx = 0.05
chain = FastChain(FastDense(2, 16, Flux.tanh), FastDense(16, 16, Flux.tanh), FastDense(16, 1))
initθ = Float64.(DiffEqFlux.initial_params(chain))
# strategy = NeuralPDE.GridTraining(dx)
strategy = NeuralPDE.QuasiRandomTraining(100)
discretization = NeuralPDE.PhysicsInformedNN(chain,strategy ; init_params = initθ)
@named pde_system = PDESystem(eq, bcs, domains, [t, x], [u(t, x)])
# prob = NeuralPDE.discretize(pde_system, discretization)
sym_prob = NeuralPDE.symbolic_discretize(pde_system, discretization)
(Expr[:((cord, var"##θ#476", phi, derivative, integral, u, p)->begin
          begin
              let (t, x) = (cord[[1], :], cord[[2], :])
                  begin
                      cord1 = vcat(t, x)
                  end
                  (+).(derivative(phi, u, cord1, [[0.0, 6.0554544523933395e-6]], 1, var"##θ#476"), derivative(phi, u, cord1, [[6.0554544523933395e-6, 0.0]], 1, var"##θ#476")) .- 0
              end
          end
      end)], Expr[:((cord, var"##θ#476", phi, derivative, integral, u, p)->begin
          begin
              let (t, x) = (cord[[1], :], cord[[2], :])
                  begin
                      cord1 = vcat(t, fill(0, size(cord[[1], :])))
                  end
                  u(cord1, var"##θ#476", phi) .- (sin).((*).(-6.283185307179586, t))
              end
          end
      end), :((cord, var"##θ#476", phi, derivative, integral, u, p)->begin
          begin
              let (t, x) = (cord[[1], :], cord[[2], :])
                  begin
                      cord1 = vcat(t, fill(1, size(cord[[1], :])))
                      cord2 = vcat(t, fill(0, size(cord[[1], :])))
                  end
                  u(cord1, var"##θ#476", phi) .- u(cord2, var"##θ#476", phi)
              end
          end
      end), :((cord, var"##θ#476", phi, derivative, integral, u, p)->begin
          begin
              let (t, x) = (cord[[1], :], cord[[2], :])
                  begin
                      cord1 = vcat(t, fill(1, size(cord[[1], :])))
                      cord2 = vcat(fill(1, size(cord[[1], :])), fill(1, size(cord[[1], :])))
                      cord3 = vcat(t, fill(0, size(cord[[1], :])))
                      cord4 = vcat(t, x)
                  end
                  (*).(u(cord2, var"##θ#476", phi), u(cord1, var"##θ#476", phi)) .- (+).((*).(-1, u(cord4, var"##θ#476", phi)), u(cord3, var"##θ#476", phi))
              end
          end
      end)])

KirillZubov avatar Dec 29 '21 15:12 KirillZubov