I wrote a fully working code to numerically solve a non-linear variation of the classic 1+2D waves equation. Here's a stripped down version of it (classic waves equation), with static initial conditions (initial velocities are all 0):


size = 30;
tmax = 20;

equDiff = D[phi[t,x,y],t,t] - D[phi[t,x,y],x,x] - D[phi[t,x,y],y,y]) == 0;

bones = Table[{x,y,RandomReal[{-1,1}]}, {x,-size,size,2}, {y,-size,size,2}];
flesh = Interpolation[Flatten[bones, 1], Method -> "Spline"];

Phi0[t_, x_, y_] := flesh[x,y] (* Random function of X and Y only. *)

boundariesCond = {
   phi[t,-size,y] == Phi0[t,-size,y],
   phi[t,size,y] == Phi0[t,size,y],
   phi[t,x,-size] == Phi0[t,x,-size],
   phi[t,x,size] == Phi0[t,x,size]

initialCond = {phi[0,x,y] == Phi0[0,x,y], (D[phi[t,x,y],t]/.t -> 0) == D[Phi0[t,x,y],t]/.t -> 0};

fieldSolution = NDSolve[
   Flatten@{equDiff, initialCond, boundariesCond}, phi,
   {t,0,tmax}, {x,-size,size}, {y,-size,size},
   Method -> {"MethodOfLines","SpatialDiscretization" -> {"TensorProductGrid","MaxPoints" -> 50,"MinPoints" -> 50,"DifferenceOrder" -> 2}}

  Evaluate[phi[t,x,y]/.fieldSolution/.t -> time],
  {x,-10,10}, {y,-10,10},
  PlotPoints -> ControlActive[20, 60],
  PlotRange -> {{-10,10}, {-10,10}, {-3,3}},
  Axes -> True,
  Mesh -> True,
  MeshFunctions -> {(#3&)},
  MeshStyle -> GrayLevel[0.25],
  ColorFunction -> "Rainbow",
  ImageSize -> 500,
  Method -> {"RotationControl" -> "Globe"},
  SphericalRegion -> True
 {time, 0, tmax, 0.1,ImageSize -> Large,Appearance -> {"Labeled", "Closed"},AppearanceElements -> {"InputField", "Slider"}},
 ControlPlacement -> Bottom

This code is working pretty well on my system.

My problem is to define smooth (and derivable) random initial conditions for both the initial field and its initial time derivative (or "velocities"), that would be consistent with the boundary conditions. Currently, the code above is using only static initial conditions, so all initial velocities are 0. This isn't very "natural" to me. So how would you define smooth and derivable random initial conditions for this system, that would be compatible with this kind of boundary conditions?

All the code should remain the same (unless you notice some mistakes or problems), except the bones and flesh functions that would need to change, so Phi0[t_, x_, y_] := flesh[x,y] would be replaced by Phi0[t_, x_, y_] := someNewRandomFunction[t,x,y].



Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.