# Cavity Flow CFD Boundary conditions and strange waves

So I have a PDE that I use to describe how material flows through a volume(2D or 3D). $$\frac{\partial C}{\partial t} + \vec{u} \cdot \nabla C = (D' + D )\nabla^2C$$ Now using finite differences I get to the following: $$C^{i+1}_{j,k} =C^i_{j,k}+\Delta t \biggr[ (D+D')\biggr[\frac{C^i_{j+1,k}-2C^i_{j,k}+C^i_{j-1,k}}{(\Delta x)^2}+\frac{C^i_{j,k+1}-2C^i_{j,k}+C^i_{j,k-1}}{(\Delta y)^2}\biggr]-u_x \frac{C^i_{j,k} -C^i_{j-1,k}}{\Delta x} - u_y \frac{C^i_{j,k} -C^i_{j,k-1}}{\Delta y}\biggr]$$ where D is diffusion coefficient and superscript is temporal, subscript is spatial. (Just imagine a grid). We are given $u_x$ and $u_y$ through some already validated code. Now, we use a boundary condition where $C=0$ at the edges. I know this is a silly boundary condition but it's irrelevant because the problem I'm about to show you persists even if we use small times.(where the concentration doesn't reach the edge)

[Look in the comments...It's telling me I can't post 2 links... Where do I have the second?]

here is the code. It's relatively simple.

import numpy as np
import matplotlib.pyplot as plt
import time
%matplotlib inline
from scipy import ndimage
import scipy.io as sio
b = a['VelField']
N =50
M =300
ux = b[0][0][0]
uy = b[0][0][1]
#ux = np.zeros((50,50))
#uy = np.ones((50,50))
dx = 2/N
dy =2/N
D= .001
Dprime =0
dt = .0001
C = np.zeros((N,N))
#C[np.floor((.25+1)/dx):np.floor((.75+1)/dx),np.floor((.25+1)/dy)\
#  :np.floor((.75+1)/dy)]=1
C[np.floor((.25+1)/dy)\
:np.floor((.75+1)/dy),.1/dx:np.floor((.25)/dx)]=1
f, axarr = plt.subplots(2, sharex=True)
axarr[0].contourf(C)
Chist= C
Ctemp = C
for T in range(M):
for j in range(N):
for k in range(N):
if k==49  or j ==49 or k ==0 or j ==0:
Ctemp[j,k]=0
else:

Ctemp[j,k]= C[j,k]+dt*((D+Dprime)*(((C[j+1,k]-2*C[j,k]+\
C[j-1,k])/(dx**2))+((C[j,k+1]-2*C[j,k]+C[j,k-1])/(dy**2)))-\
ux[j,k]*(C[j,k]-C[j-1,k])/(dx)-uy[j,k]*(C[j,k]-C[j,k-1])/dy)
#Chist.append(list(C))
C = Ctemp
axarr[0].set_title('Original') axarr[2].contourf(C)
axarr[2].set_title('{} seconds later'.format(dt*M))


So the question is why do I have these waves and if anyone knows effective boundary conditions I would be grateful. Thanks in advance.

• I can't add the picture until I'm 10 reputation. If you up-vote twice. I will put the image. – Aviv Moshe Albeg Aug 25 '16 at 19:47
• With a picture of the waves, it's possible to give a physical explanation to it -- but odds are good it's a numerical problem. We will have to wait and see what the problem actually is. – tpg2114 Aug 25 '16 at 20:41
• I have 2 possible explanations. (a) You're writing dy=2/N, but bot 2 and N are integers, so dy might be 0 in that case. Try something like dy=2.0/N. The same goes for dx. (b) The problem might also be with your input velocity field ux & uy. Try (a), if that doesn't work, take a look at the velocity field. – Hayk Hakobyan Aug 25 '16 at 23:09
• I'm voting to close this question as off-topic because it is about debugging specific code. – ACuriousMind Aug 26 '16 at 0:02
• Are you using a value of $\frac{u\Delta t}{\Delta x}>1$? – Chester Miller Aug 26 '16 at 3:20