I’m trying to use a basic advection/diffusion equation with MethodOfLines.jl. I have an initial condition which is a Gaussian shaped density, and I’m struggling to get periodic boundary conditions to work.
My initial state looks like this, so perhaps the left boundary is trying to take the value of the right? I have tried swapping the order in bcs but no luck.
Any help appreciated!
I get
Warning: The system contains interface boundaries, which are not compatible with system transformation. The system will not be transformed. Please post an issue if you need this feature.
using ModelingToolkit
using DifferentialEquations
using MethodOfLines
using Revise
using DomainSets
@variables x t
@variables c(..)
∂t = Differential(t)
∂x = Differential(x)
# advection-diffusion constants
R = 1.0
D = 0.1
v = 1.0
cₒ = 1.0
x_min = 0.0
x_max = 1.0
t_min = 0.0
t_max = 1.0
c_normal(x) = exp(-(10*x)^2)
# advection-diffusion equation
eqn = R * ∂t(c(x, t)) ~ D * ∂x(c(x, t))^2 - v * ∂x(c(x, t))
# initial and boundary conditions
bcs = [c(x,t_min) ~ c_normal(x),
c(x_min,t) ~ c(x_max,t)]
# space and time domains
dom = [x ∈ Interval(x_min, x_max), t ∈ Interval(t_min, t_max)]
# define PDE system
@named sys = PDESystem(eqn, bcs, dom, [x, t], [c(x, t)])
# convert the PDE into an ODE problem
prob = discretize(sys, MOLFiniteDifference([x => 100], t))#, advection_scheme=WENOScheme()))
# solve the problem
sol = solve(prob, Rosenbrock23(), saveat=0.2)
c = sol[c(x,t)]
using Plots
plot(xlim=(0,100),ylim=(0,10))
n = size(c, 2)
anim = @animate for t ∈ 1:n
plot!(c[:,t])
end
gif(anim, "diff.gif", fps=1)
(edit: I have also tried adding the BC on dc/dx (although this is not specified in the Periodic BCs section of MethodOfLines) and the result is not correct still)