

In [6]:

    
%pylab inline
import pylab as pl
import math
import numpy as np
import matplotlib.pyplot as plt
import csv
# number of the grid points 

ymax = 31;

# material properties 

rho = 1000.0;
mu = 0.00102;
nu = mu / rho ; 

#dimensions of the computational domain
l = 0.5 ;
h = 1 ;
w = 1.0 ;

# define Reynolds number

Re = 2000 ;

# compute the average velocity

ua = (Re*nu)/(2.0*h);

# compute the pressure drop

dp = (12.0*mu*l*ua)/(h*h);

# compute the volume flow rate 

Q = (w*dp*h*h*h)/(12.0*mu*l);

# create the grid 

dy = h/(ymax-1.0)

y = np.zeros(ymax)

y[0] = 0.0

for i in range(1,ymax):
    y[i] = y[i-1]+dy
      
# compute the analytical solution
u =np.zeros(ymax)

#define the boundary conditions...define them inside the loop
u[0] = 0.0
u[ymax-1] = 0.0
const = dp/(2.0*mu*l )

for j in range(1,ymax-1):
    u[j] = const * y[j]*(h-y[j])
   

print ua

plt.plot(u,y,'x')
plt.ylabel('Position [m]')
plt.xlabel('Velocity [m/s]')
plt.title('1D laminar velocity profile in a channel ')
plt.show()









    



Populating the interactive namespace from numpy and matplotlib
0.00102