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