I am attempting to solve for waves on a water surface starting with a two dimensional solution. The equations are that the water must satisfy Laplace's equation everywhere with a time dependent boundary condition on the top surface. The equations are
$$\begin{align*} &\nabla^2\phi=0\\ &\left(\frac{\partial\phi}{\partial y}+\frac1{g}\frac{\partial^2\phi}{\partial t^2}\right)_{y=y_2}=0 \end{align*}$$
where $x$ is horizontal and $y$ vertical; also, $g$ is the acceleration due to gravity. The second equation applies at the surface. The first term in this equation is the vertical gradient so I think I can use a Neumann condition to enter this equation.
So I start by making a grid.
Needs["NDSolve`FEM`"];
x2 = 4; y2 = 1;
reg = ImplicitRegion[0 <= x <= x2 && 0 <= y <= y2, {x, y}];
mesh = ToElementMesh[reg, "BoundaryMeshGenerator" -> {"Continuation"},
MaxCellMeasure -> .002, "MaxBoundaryCellMeasure" -> 0.01];
Show[mesh["Wireframe"], Frame -> True, PlotRange -> All]
I now make some initial conditions. I have worked out an initial velocity for the top surface and solve to give the initial value of the potential function.
a = 0.40825787026798765;
b = 0.1689925306573793;
g = 9.81;
solIC = NDSolveValue[{
Laplacian[ϕ0[x, y], {x, y}] ==
NeumannValue[-b (E^-(x - 1)^2 - a), 0 <= x <= x2 && y == y2],
DirichletCondition[ϕ0[x, y] == 0, x == 0 && y == 0]
}, ϕ0, {x, y} ∈ mesh];
The velocity on the surface and throughout the water looks reasonable.
Plot[Evaluate[(D[solIC[x, y], y]) /. y -> y2], {x, 0, x2},
AspectRatio -> 1/4]
ClearAll[f];
f[x_, y_] := Evaluate[Grad[solIC[x, y], {x, y}, "Cartesian"]];
StreamPlot[f[x, y], {x, 0, x2}, {y, 0, y2}, AspectRatio -> Automatic]
Now I attempt to set up the time dependent problem.
sol = NDSolveValue[{
Laplacian[ϕ[x, y, t], {x, y}] ==
NeumannValue[1/g Derivative[0, 0, 2][ϕ][x, y, t],
0 <= x <= x2 && y == y2],
DirichletCondition[ϕ[x, y, t] == 0, x == 0 && y == 0],
ϕ[x, y, 0] == solIC[x, y]
}, ϕ,
{x, y} ∈ mesh, {t, 0, 0.1}
]
This does not work and gives me the error
CoefficientArrays::poly: ϕ$1664+ϕ$1665-NeumannValue[0.101937 ϕ$1666,0<=x<=4&&y==1] is not a polynomial. >>
Now I am stuck. This does not help. Any suggestions?
Suggestions from comments
If I put in a time derivative in the equation this gives the same error.
If I take the time derivative out of the
NeumannValue
and put a time derivative into the equation then it solves but is meaningless. Does this suggest we can't have time derivatives in the boundary conditions?
Can you suggest how we might set up an alternative equation who's solution might go into the boundary condition?
sol = NDSolveValue[{
Laplacian[ϕ[x, y, t], {x, y}] == 0.1 D[ϕ[x, y, t], t] +
NeumannValue[ϕ[x, y, t], 0 <= x <= x2 && y == y2],
DirichletCondition[ϕ[x, y, t] == 0, x == 0 && y == 0],
ϕ[x, y, 0] == solIC[x, y]
}, ϕ,
{x, y} ∈ mesh, {t, 0, 0.1}
];
Plot3D[sol[x, y2, t], {x, 0, x2}, {t, 0, 0.1}]
So the problem could be the derivative in the NeumannValue
.
Derivative[0, 0, 1][\[Phi]][x, y, t]
orD[\[Phi][x, y, t], {t, 1}]
in the time dependent equation? Also is this a wave equation? (second order time deriv? - but maybe the title is miss leading me here.) $\endgroup$