# Solving advection equation - periodic conditions - using roll python function [closed]

The original post was on stackoverflow : I transfert it here.

I have to solve numerically the advection equation with periodic boundaries conditions : u(t,0) = u(t,L) with L the length of system to solve.

I start also with u(0,x) = uexacte(0,x) = sin(2*pi*x/L)

Here the main part of the code with loop time (we use here FTCS scheme) :

V=1
L=1
# analytical solution --------------------------
def uexacte(t,x):
return sin(2*pi*(x-V*t)/L)

# 1. Centre FTCS (Forward Time Centered Space)
cfl = 0.25
nx = 10
tend = 1

#
dx = L/(nx-1.)
dt = cfl*dx/V
nt = int(tend/dt)+1
print "CFL=%5.2f tend=%4.1f --> %i iterations en temps"%(cfl,tend,nt)

# Arrays
x = linspace(0,L,nx)
u0 = uexacte(0,x)

# Starting solution
t=0.0 ; u=copy(u0)

# Time loop
for i in range(1,nt):
# FTCS
#u[1:nx-1] = u[1:nx-1] - cfl/2*(u[2:nx] - u[0:nx-2])
# Using roll
u = u + - cfl/2*(roll(u,-1)- roll(u,1))
# Update time
t = t+dt


I don't understand the solution given by teacher who uses the python function roll in this way :

# Using roll
u = u - cfl/2*(roll(u,-1)- roll(u,1))


One says that with the using of roll, we are sure to respect the periodic boundary conditions but I don't understand why ?

Indeed, my first approach was to do :

u[0] = u[nx-1]
u[1:nx-1] = u[1:nx-1] - cfl/2*(u[2:nx] - u[0:nx-2])


but this doesn't work and I don't know how to implement this periodic conditions in this way (without using roll function).

If someone could explain this matter and the trick with roll function, this would be nice to tell it.

UPDATE 1 :

I tried with classical approach (simple recurrence formula ) like this :

# Time loop
for i in range(1,nt):
# FTCS
u[1:nx-1] = u[1:nx-1] - cfl/2*(u[2:nx] - u[0:nx-2])

# Try to impose periodic boundary conditions but without success
u[0] = u[0] - cfl/2*(u[0] - u[nx-1])

# Update time
t = t+dt


Indeed, the result is bad (values for each side are not the same). I could impose the theorical values at each step but in practise, we don't know always the analytical solution.

What is the trick to impose this periodic boundary condition on the numerical solution at each step of time ?

The problem with periodic boundary conditions is formulated as :

UPDATE 2 :

The solution given by LonelyProf below works fine. Thanks

• I have a general query, unrelated to my answer (which addressed only the original question, concerning periodic boundaries). You seem to be using the FTCS scheme for the advection equation. Isn't this unstable? See en.m.wikipedia.org/wiki/FTCS_scheme – user28077 Dec 3 '18 at 21:38
• yes, you're right but I think if I choose a little CFL, I have, at least at the beginning, the numerical solution closed to analytical solution – youpilat13 Dec 3 '18 at 21:39
• With the updates, this question has definitely moved into the realm of "debug my code" and doesn't seem to be applicable to others any longer. – David Ketcheson Dec 11 '18 at 11:43

Similar issues arise in other languages, so this is not Python-specific.

The NumPy function roll applies a cyclic shift to the array given as its first argument.

• roll(u,-1) returns an array with elements u[1],u[2],u[3]...,u[nx-1],u[0]
• roll(u,1) returns an array with elements u[nx-1],u[0],u[1]...,u[nx-3],u[nx-2].

This is exactly what you want on the right hand side of your assignment, accessing the elements on either side of every element of u, using periodic boundaries.

You could do this by array slices, but then you would need to cope with the end effects by hand. So you can write a single statement to update all of u[1:nx-2] in terms of u[0:nx-3] and u[2:nx-1], and write separate statements to update u[0] in terms of u[nx-1] and u[1], and to update u[nx-1] in terms of u[nx-2] and u[0].

Currently, you are making mistakes in the way you try to do this by hand (using an illegal index nx, for instance, or simply writing the wrong indices on the right).

Another issue is that using a function such as roll enables you to do the update effectively simultaneously for all the elements of u. In some applications, perhaps yours, it may be an error (and is certainly different) to update any of the elements, and then use the updated value on the right-hand side to update other elements. You can avoid this by defining temporary variables to hold the updated values, if you are sticking with the "by hand" approach. Or define an array du to hold the change in u, calculated from the current values, and actually perform the update u=u+du in a separate statement, all at once. Using the roll function (or similar) allows the code to be more compact and simpler.

utemp1=u[0] - cfl/2*(u[nx-1] - u[1]) utemp2=u[nx-1] - cfl/2*(u[nx-2] - u[0]) u[1:nx-2] = u[1:nx-2] - cfl/2*(u[0:nx-3] - u[2:nx-1]) u[0] = utemp1 u[nx-1] = utemp2 
• -@LonelyProf I am near to the solution but what do you mean by So you can write a single statement to update all of u[1:nx-2] in terms of u[0:nx-3] and u[2:nx-1], and write separate statements to update u[0] in terms of u[nx-1] and u[1], and to update u[nx-1] in terms of u[nx-2] and u[0]. ?, I mean what's the traduction in terms of python code ? – youpilat13 Dec 3 '18 at 21:35