I recently asked a question about this topic, but too much was unclear to me at the time which yielded poor response. (You can find the previous question here). Anyway, I am now trying to solve the linear scalar advection equation analytically with solid wall boundary conditions (i.e. reflection at the boundaries). Something does not quite work in my code: I have a wave moving leftwards toward the left boundary, but suddenly the wave just diminshes by some factor. Then the wave proceeds to the boundary and is reflected as should be. See the results below.

enter image description here

enter image description here

Now I know what causes this problem actually, but I do not now how to resolve it. In my solution I let everything on the blue line travel to the left (propagate to the left) and then reflect on the left boundary. Thus there are waves of -0.5 amplitude going to the right "colliding" with waves going to the left. I solve this by just superpositioning these (two) waves and dividing by two. This works fine for when the colliding waves have the same amplitude, but when they are not of the same amplitude the total wave diminishes like in the picture. But how do I go about conserving amplitude or area with all these reflecting waves? What am I doing wrong?

You may test my code in MATLAB:

L = 10;
H = 1;
g = 9.61;
w = 0.4;
eps = 0.4;

lambda1 = -sqrt(g*h0(L,H,w,eps));
lambda2 = sqrt(g*h0(L,H,w,eps));

N = 100;
dx = L/N;
xgrid = 0:dx:L;

dt = 0.01;
time = 0:dt:2;
T = L/(-2*lambda1);

for t = time(1:100)
    sol = zeros(1,N+1);
    i = 0;

    if 0 <= t && t < T
        for x = xgrid
            i = i + 1;
            if -lambda1*t <= x && x <= L + lambda1*t
                sol(i) = z10(x-lambda1*t,L,H,w,eps);
            elseif 0 <= x && x < -lambda1*t
                xarg = x-lambda1*t;
                sol(i) = (1/2)*(zr(x,t,L,H,w,eps,lambda1) + z10(xarg,L,H,w,eps));
            elseif L+lambda1*t < x && x <= L
                sol(i) = z10(L,L,H,w,eps);
    elseif T <= t && t <= 2*T
        for x = xgrid
            i = i + 1;
            if L+lambda1*t <= x && x <= -lambda1*t
                sol(i) = zr(x,t,L,H,w,eps,lambda1);
            elseif 0 <= x && x < L+lambda1*t
                sol(i) = (1/2)*(zr(x,t,L,H,w,eps,lambda1) + z10(x-lambda1*t,L,H,w,eps));
            elseif -lambda1*t <= x && x <= L
                sol(i) = z10(L,L,H,w,eps);
    axis([0 L -0.7 -0.5])

axis([0 L -0.7 -0.5])
title(['Analytical solution for linearized system at t=',num2str(t)])

function val = h0(L,H,w,eps)
val = H + (eps*w*sqrt(pi))*erf(L/(2*w))/L;
function val = zr(x,t,L,H,w,eps,lambda1)
% Right boundary at x=0
tr = t + x/lambda1;
xarg = -lambda1*tr;
val = z10(xarg,L,H,w,eps);
function val = z10(x,L,H,w,eps)
% Initial data
val = -(H + eps*exp(-(x-L/2)^2/(w^2)))/2;

The code "works" for t < 3, but I have not implemented boundary conditions at the other end so it will not show correctly for times larger than that. I hope I make the question easy to answer!

Best regards/


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.