11
$\begingroup$

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]

enter image description here

How to solve this equation over the Black Sea shape region?

 mask = GeoGraphics[{GeoStyling[RGBColor["Aqua"]], 
 Polygon[Entity["Ocean", "BlackSea"]]}, GeoBackground -> None]

enter image description here

Improve this question
$\endgroup$
6
  • 3
    $\begingroup$ Ω = DiscretizeRegion[ Entity["Ocean", "BlackSea"]["Polygon"] /. GeoPosition -> Identity] $\endgroup$
    – Edmund
    Nov 30, 2016 at 11:21
  • $\begingroup$ Somehow the answer seems trivially obvious. Could you say more about the problem you're facing? $\endgroup$
    – Michael E2
    Nov 30, 2016 at 11:23
  • 1
    $\begingroup$ @Edmund That, like my Ω = Cases[mask, _Polygon, Infinity][[1]], doesn't work. I mean, it produces a region, but does not produce a solution to the equations. $\endgroup$
    – corey979
    Nov 30, 2016 at 11:29
  • $\begingroup$ Might need {x, y} ∈ Ω in Plot3D instead of "ElementMesh". $\endgroup$
    – Edmund
    Nov 30, 2016 at 11:32
  • 1
    $\begingroup$ @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. $\endgroup$
    – Edmund
    Nov 30, 2016 at 11:38

2 Answers 2

Reset to default
8
$\begingroup$

If you use something like this:

Needs["NDSolve`FEM`"];
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]
Improve this answer
$\endgroup$
0
9
$\begingroup$

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}}]

enter image description here

Here is the solution for {t,0,4}:

enter image description here

EDIT: Another gif with 100 frames and PlotPoints -> 50

enter image description here

Improve this answer
$\endgroup$
3
  • 2
    $\begingroup$ Looks like the wavelength is too short for the granularity of the discretization ... $\endgroup$
    – Szabolcs
    Nov 30, 2016 at 20:20
  • $\begingroup$ @Szabolcs Yes, an accurate solution needs some tweaking :) $\endgroup$
    – Gyebro
    Nov 30, 2016 at 20:22
  • 6
    $\begingroup$ 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. $\endgroup$
    – Hugh
    Nov 30, 2016 at 21:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.