g = 9.81
@parameters x z t
@variables φ(..) φ̃(..) η(..)
Dt = Differential(t)
Dx = Differential(x)
Dz = Differential(z)
Dxx = Differential(x)^2
Dzz = Differential(z)^2
eqs = [Dxx(φ(t, x, z)) + Dzz(φ(t, x, z)) ~ 0,
Dt(φ̃(t,x)) ~ -g*η(t,x),
Dt(η(t,x)) ~ Dz(φ(t, x, 1.))
]
bcs = [
φ(0, x, z) ~ 0,
φ̃(0., x) ~ 0.,
η(0., x) ~ cos(2*π*x),
φ(t, x, 1.) ~ φ̃(t, x),
Dx(φ(t, 0., z)) ~ 0.0,
Dx(φ(t, 1., z)) ~ 0.0,
Dz(φ(t, x, 0.)) ~ 0.0,
Dx(φ̃(t, 0.)) ~ 0.0,
Dx(φ̃(t, 1.)) ~ 0.0,
Dx(η(t, 0.)) ~ 0.0,
Dx(η(t, 1.)) ~ 0.0,
]
domains = [x ∈ Interval(0.0, 1.0),
z ∈ Interval(0.0, 1.0),
t ∈ Interval(0.0, 10.0)]
@named pdesys = PDESystem(eqs, bcs, domains, [t, x, z],
[φ(t, x, z), φ̃(t, x), η(t, x)])
dx = 0.1
dz = 0.1
order = 2
discretization = MOLFiniteDifference([x=>dx, z=>dz], t,
approx_order=order,
grid_align=center_align)
println("Discretization:")
prob = discretize(pdesys, discretization)
