2D example

The 2D code works similarly to the 1D version, with a few different function calls.

using FaADE

This example gives a 1D solution to a 2D problem.

We'll solve the 2D heat equation.

\[ \frac{\partial u}{\partial t} = K\frac\Delta u\]

with boundary conditions

\[ u(0,t) = 0, \qquad \partial_x u(1,t) = 0\]

and initial condition

\[ u(x,0) = \exp\left(\frac{-(x-0.5)^2}{0.02}\right)\]

We first need to create a domain to solve the PDE using Grid2D

𝒟x = [0.0,1.0]
𝒟y = [0.0,1.0]
nx = ny = 41
grid = Grid2D(𝒟x,𝒟y,nx,ny)

Grid2D{Float64, CartesianMetric, Matrix{Float64}}([0.0 0.0 … 0.0 0.0; 0.025 0.025 … 0.025 0.025; … ; 0.975 0.975 … 0.975 0.975; 1.0 1.0 … 1.0 1.0], [0.0 0.025 … 0.975 1.0; 0.0 0.025 … 0.975 1.0; … ; 0.0 0.025 … 0.975 1.0; 0.0 0.025 … 0.975 1.0], 0.025, 0.025, 41, 41, [1.0;;], [1.0;;], [1.0;;], [1.0;;], [1.0;;])

The initial condition

u₀(x,y) = exp(-((x-0.5)^2 + (y-0.5)^2) / 0.02)

u₀ (generic function with 1 method)

Set the FD order to 2,

order = 2

The boundary conditions are defined by creating SimultanousApproximationTerm objects, which will then be fed to the PDE structure

BoundaryLeft    = SAT_Dirichlet((X,t)->0.0,grid.Δx, Left, order)
BoundaryRight   = SAT_Neumann((X,t)->0.0, grid.Δx, Right, order)
BoundaryUp      = SAT_Periodic(grid.Δy, Up, order)
BoundaryDown    = SAT_Periodic(grid.Δy, Down, order)

BCs = (BoundaryLeft,BoundaryRight,BoundaryUp,BoundaryDown)

(SAT_Dirichlet{FaADE.Helpers.NodeType{:Left, 1}, :Cartesian, Float64, Vector{Float64}, Main.var"#1#2", typeof(eachcol)}(FaADE.Helpers.NodeType{:Left, 1}(), 1, 2, Main.var"#1#2"(), 6400.0, [3200.0, -1600.0], 0.025, 1.0, [-2.0], eachcol, 0.0, :Cartesian), SAT_Neumann{FaADE.Helpers.NodeType{:Right, 1}, :Cartesian, Float64, Vector{Float64}, Main.var"#3#4", typeof(eachcol)}(FaADE.Helpers.NodeType{:Right, 1}(), 1, 2, Main.var"#3#4"(), -80.0, [-40.0, 40.0], 0.025, 1.0, eachcol, 0.0, :Cartesian), SAT_Periodic{FaADE.Helpers.NodeType{:Right, 2}, :Cartesian, Float64, Vector{Float64}, FaADE.SATs.var"#τ₀#1"{Float64, Float64}, typeof(eachrow)}(FaADE.Helpers.BoundaryConditionType{:Periodic}(), FaADE.Helpers.NodeType{:Right, 2}(), 2, 2, [-80.0, 40.0], [-40.0, 80.0], [-40.0, 40.0], [-40.0, 40.0], 0.025, 40.0, -20.0, FaADE.SATs.var"#τ₀#1"{Float64, Float64}(0.025, 0.5), eachrow, 0.0, :Cartesian), SAT_Periodic{FaADE.Helpers.NodeType{:Left, 2}, :Cartesian, Float64, Vector{Float64}, FaADE.SATs.var"#τ₀#1"{Float64, Float64}, typeof(eachrow)}(FaADE.Helpers.BoundaryConditionType{:Periodic}(), FaADE.Helpers.NodeType{:Left, 2}(), 2, 2, [-80.0, 40.0], [-40.0, 80.0], [-40.0, 40.0], [-40.0, 40.0], 0.025, 40.0, -20.0, FaADE.SATs.var"#τ₀#1"{Float64, Float64}(0.025, 0.5), eachrow, 0.0, :Cartesian))

The 2 input to the periodic boundary ensures it is along the y-axis.

Set the diffusion in $x$ and $y$ directions to 1

Kx = Ky = 1.0

1.0

Now we can create a Problem2D object to pass to the solver,

P = Problem2D(order,u₀,Kx,Ky,grid,BCs)

2 dimensional PDE Problem

Lastly before solving we define our time step and simulation time,

Δt = 0.01grid.Δx;
t_f = 200Δt;

Finally we call the solver (currently not working with Documenter.jl)

soln = solve(P,grid,Δt,t_f)

FaADE.solvers.solution{Float64, Matrix{Float64}, Grid2D{Float64, CartesianMetric, Matrix{Float64}}, Problem2D{Float64, 2, Float64, FaADE.Helpers.SourceTerm{Nothing}, Tuple{SAT_Dirichlet{FaADE.Helpers.NodeType{:Left, 1}, :Cartesian, Float64, Vector{Float64}, Main.var"#1#2", typeof(eachcol)}, SAT_Neumann{FaADE.Helpers.NodeType{:Right, 1}, :Cartesian, Float64, Vector{Float64}, Main.var"#3#4", typeof(eachcol)}, SAT_Periodic{FaADE.Helpers.NodeType{:Right, 2}, :Cartesian, Float64, Vector{Float64}, FaADE.SATs.var"#τ₀#1"{Float64, Float64}, typeof(eachrow)}, SAT_Periodic{FaADE.Helpers.NodeType{:Left, 2}, :Cartesian, Float64, Vector{Float64}, FaADE.SATs.var"#τ₀#1"{Float64, Float64}, typeof(eachrow)}}, Nothing}}([[1.3887943864964021e-11 4.698230849877748e-11 … 4.698230849877748e-11 1.3887943864964021e-11; 4.698230849877748e-11 1.5893910094516368e-10 … 1.5893910094516368e-10 4.698230849877748e-11; … ; 4.698230849877748e-11 1.5893910094516368e-10 … 1.5893910094516368e-10 4.698230849877748e-11; 1.3887943864964021e-11 4.698230849877748e-11 … 4.698230849877748e-11 1.3887943864964021e-11], [-2.3973166445963694e-9 1.1624443972996291e-9 … 1.162444446060902e-9 -2.3973166445963694e-9; 0.004249272771844638 0.004264705095316439 … 0.004264705093887542 0.004249272771844638; … ; 0.03751274636963689 0.03764898860648547 … 0.03764898860021413 0.03751274636963689; 0.037377449267445684 0.03751319980431178 … 0.03751319981099281 0.037377449267445684]], Grid2D{Float64, CartesianMetric, Matrix{Float64}}([0.0 0.0 … 0.0 0.0; 0.025 0.025 … 0.025 0.025; … ; 0.975 0.975 … 0.975 0.975; 1.0 1.0 … 1.0 1.0], [0.0 0.025 … 0.975 1.0; 0.0 0.025 … 0.975 1.0; … ; 0.0 0.025 … 0.975 1.0; 0.0 0.025 … 0.975 1.0], 0.025, 0.025, 41, 41, [1.0;;], [1.0;;], [1.0;;], [1.0;;], [1.0;;]), [0.00025, 0.00025], [0.0, 0.05], 2 dimensional PDE Problem, 0.0026295125322778664, Float64[])

The solver ourputs a solution data structure, with everything packaged in that we would need to reconstruct the problem from the final state if we wanted to restart.

No visualisation routines are written at the moment but we imported the Plots.jl package earlier so we'll use that

using Plots
surface(grid.gridx, grid.gridy, soln.u[2])

This page was generated using Literate.jl.