# Solving a wave equation over the Black Sea shaped region

There is a code for solving a wave equation over an arbitrarily shaped region.

Ω =
RegionDifference[RegionUnion[Disk[], Rectangle[{0, 0}, {2, 2}]],
Disk[{1/4, 1/4}, 1/5]];

uifWave = NDSolveValue[{D[u[t, x, y], t, t] - Inactive[Laplacian][u[t, x, y], {x, y}] == 0,
u[0, x, y] == E^(-5*((x - 3/2)^2 + (y - 3/2)^2)),
Derivative[1, 0, 0][u][0, x, y] == 0,
DirichletCondition[u[t, x, y] == 0,
True]},
u, {t, 0, 2 π}, {x, y} ∈ Ω] //
Quiet;

framesWEQ =
Table[Plot3D[
uifWave[t, x, y], {x, y} ∈ uifWave["ElementMesh"],
PlotRange -> {-1, 1}, Boxed -> False, Axes -> False,
Mesh -> None], {t, 0, 2 π, 2 π/50}];

Manipulate[framesWEQ[[i]], {{i, 16, "time"}, 1, Length[framesWEQ], 1},
SaveDefinitions -> True] How to solve this equation over the Black Sea shape region?

 mask = GeoGraphics[{GeoStyling[RGBColor["Aqua"]],
Polygon[Entity["Ocean", "BlackSea"]]}, GeoBackground -> None] • Ω = DiscretizeRegion[ Entity["Ocean", "BlackSea"]["Polygon"] /. GeoPosition -> Identity] Nov 30, 2016 at 11:21
• Somehow the answer seems trivially obvious. Could you say more about the problem you're facing? Nov 30, 2016 at 11:23
• @Edmund That, like my Ω = Cases[mask, _Polygon, Infinity][], doesn't work. I mean, it produces a region, but does not produce a solution to the equations. Nov 30, 2016 at 11:29
• Might need {x, y} ∈ Ω in Plot3D instead of "ElementMesh". Nov 30, 2016 at 11:32
• @corey979 I get a solution for uifWave that plots but I think the wave is very small over the area. May need to increase the amplitude in u. Nov 30, 2016 at 11:38

If you use something like this:

Needs["NDSolveFEM"];
dg = DiscretizeGraphics[
Entity["Ocean", "BlackSea"]["Polygon"] /.
GeoPosition -> Identity];
mesh = ToElementMesh[dg, MaxCellMeasure -> 0.005 (*, "MeshOrder"->1 *)]


You will be able to refine the mesh to your needs. This will then do the time integration and monitor it's progress:

Monitor[uifWave =
NDSolveValue[{D[u[t, x, y], t, t] -
Inactive[Laplacian][u[t, x, y], {x, y}] == 0,
u[0, x, y] == Exp[-10 ((x - 44)^2 + (y - 31)^2)],
Derivative[1, 0, 0][u][0, x, y] == 0,
DirichletCondition[u[t, x, y] == 0, True]},
u, {t, 0, 2}, {x, y} \[Element] mesh,
EvaluationMonitor :> (monitor = Row[{"t = ", CForm[t]}])], monitor]


For completeness:

framesWEQ =
Table[Plot3D[uifWave[t, x, y], {x, y} \[Element] dg,
PlotRange -> {-1, 1}, Boxed -> False, Axes -> False,
Mesh -> None], {t, 0, 2, 2/50}];
Manipulate[framesWEQ[[i]], {{i, 16, "time"}, 1, Length[framesWEQ], 1},
SaveDefinitions -> True]


It works with Edmund's suggestion, only the initial state is changed to

u[0, x, y] == Exp[-10 ((x - 44)^2 + (y - 31)^2)].

Coordinates {44, 31} are inside the Black Sea.

For the sake of simplicity I've restricted the region to the western third of the sea.

Ω =  DiscretizeRegion[Entity["Ocean", "BlackSea"]["Polygon"] /. GeoPosition -> Identity,
{{40, 50}, {27, 33}}] Here is the solution for {t,0,4}: EDIT: Another gif with 100 frames and PlotPoints -> 50 • Looks like the wavelength is too short for the granularity of the discretization ... Nov 30, 2016 at 20:20
• @Szabolcs Yes, an accurate solution needs some tweaking :) Nov 30, 2016 at 20:22
• Just being picky; these are membrane waves. We need to solve for water surface waves which we don't know how to do yet see here.
– Hugh
Nov 30, 2016 at 21:29