Generate a small random disturbance on a flat surface for NDSlove

Question

I want to generate a small random perturbation on a flat surface defined on a square domain, e.g. 1 + 0.05 rand(x,y), which can be used as initial condition with periodic boundary conditions in NDSolve. In other word, rand(x,y) is a special pseudo-random function distributed in the interval (-1,1) which has periodicity between the sides of the square-domain [0,L]*[0,L], that is, rand(0,y)=rand(L,y) and rand(x,0)=rand(x,L). In fact, this question is very similar to this question.

However, after I went through this post in details, all comment, voted answers and even "accepted" answer as well as the answer given by the author, I found this question is still unsolved.

Based on acl's comment and Vitaliy's answer, I have the following working example with two different initial conditions.

Off[NDSolve::mxsst];
L = 4*π;
c = -(1/20);
tmax = 1000;
ini1 = Interpolation@Flatten[Table[{{x, y}, 1 + c*RandomReal[{-1, 1}]}, {x, 0, L + 1}, {y, 0, L + 1}], 1]; (*initial condition given by acl*)
cond3D = First[h /. NDSolve[{D[h[x, y, t], t] + 
Div[h[x, y, t]^3*Grad[Laplacian[h[x, y, t], {x, y}], {x, y}], {x, y}] + 
Div[h[x, y, t]^3*Grad[h[x, y, t], {x, y}], {x, y}] == 0,
h[0, y, t] == h[L, y, t],
h[x, 0, t] == h[x, L, t],
h[x, y, 0] == ini1[x, y]
 },
h, {x, 0, L}, {y, 0, L}, {t, 0, tmax},
Method -> {"BDF", "MaxDifferenceOrder" -> 1}, 
MaxStepFraction -> 1/50]]

I have the following error:

NDSolve::ibcinc: Warning: boundary and initial conditions are inconsistent.

If I use the other initial condition given by @Vitaliy just with a slight modification:

ini2 = 1 + c*BSplineFunction[RandomReal[1, {30, 30, 1}], SplineClosed -> True][##] &;

where, because I want a small perturbation so a small parameter c. Just replace the initial condition as

h[x, y, 0] == ini2[x, y]

The NDSlove then eats up the memory and spit a lot a error related to data.

I realized that Vitaliy's answer is wonderful in SplineClosed -> True that makes I can use it with periodic boundary conditions in NDSolve. But another confused issue is when I try, ini2[9, 0] === ini2[9, L] // FullSimplify using my ini2, MMA gives False.I know @Vitaliy produced a disturbed surface in the domain[0, 1]* [0, 1], but how to get such a surface in a square domain, [0, L]*[0, L].

How can I generate a random disturbance consisting with periodic boundary condition meanwhile. If anyone can give such an initial condition that has something like the SplineClosed->True feature that can be fed in NDSolve, that would be truly awesome! I admit I am being very silly when I try to use RandomReal[] in this way.

Answer 1

Every function of the following form satisfies your boundary conditions: $C_{i,j}=\left(a_{i,j}\cos(2iπ x)+b_{i,j}\sin(2iπ x\right)\left(c_{i,j}\cos(2jπ y)+d_{i,j}\sin(2jπ y\right)$

Therefore $f=∑_{i,j} C_{i,j}$ also satisfies, for any random coefficients $a$, $b$, $c$, $d$. Here's a possible implementation:

component[i_, j_] := (a[i, j] Cos[2 i π x / lL] + 
                      b[i, j] Sin[2 i π x / lL]) (c[i, j] Cos[2 j π y / lL] + 
                      d[i, j] Sin[2 j π y / lL])
fn[imax_, jmax_] := Sum[component[i, j], {i, 0, imax}, {j, 0, jmax}]
a[i_, j_] := RandomReal[{-1, 1}]
b[i_, j_] := RandomReal[{-1, 1}]
c[i_, j_] := RandomReal[{-1, 1}]
d[i_, j_] := RandomReal[{-1, 1}]
lL=4 π;
f = fn[5, 5]/5;
Plot3D[Evaluate@f, {x, 0, lL}, {y, 0, lL}]

It works with your code, as long as you set the initial condition to be h[x,y,0]==f.