Chapter 6 - Convection-Diffusion Problems

Exa 6.1 Page 106

from __future__ import division
from numpy import zeros,exp
#given convective diffusive eqn--> -u*(dc/dx)+D*(d2C/dx2)=0
C_ini=0             #at x=0
C_end=1             #at x=1


#using central difference method for both diffusion and convective term 
#-u*(C(i+1)-C(i-1))/(2*delta_x) + D*(C(i+1)+C(i-1)-2*C(i))/delta_x**2 = 0
delta_x=1/50
#on solving the given eqns and by using the given boundary eqns we have
Pe=50                              #given 
Pe_local=50*delta_x                #u/D=50 as l=1
alpha=Pe_local-2                   #co-eff of C(i+1)
Beta=Pe_local+2                    #co-eff of C(i-1)
#multipling with -2*delta_x**2/D we get
#-(Pe_local+2)*C(i-1) + 4*C(i) + (Pe_local-2)*C(i+1)=0
#solving eqns using TDMA method
a=[0]
for i in range(2,50):
    a.append(-Beta)            #sub diagonal assignment

b=[]
for j in range(1,50):
    b.append(4)           #main diagonal assignment

c=[]
for k in range(1,49):
    c.append(alpha)            #super diagonal assignment

d=zeros(49)    
d[0]=Beta*C_ini
d[-1]=-alpha*C_end
for l in range(2,49):
    d[l-1]=0
              #given values assignment
i=1#
n=49#
beta1=[b[i-1]]                #initial b is equal to beta since a1=0
gamma1=[d[i-1]/beta1[i-1]]      #since c7=0
m=i+1#
for j in range(m,n+1):
    beta1.append(b[j-1]-a[j-1]*c[j-2]/beta1[j-2])
    gamma1.append((d[j-1]-a[j-1]*gamma1[j-2])/beta1[j-1])

x=zeros(n)
x[-1]=gamma1[n-1]               #since c7=0
n1=n-i#
for k in range(1,n1+1):
    j=n-k#
    x[j-1]=gamma1[j-1]-c[j-1]*x[j]/beta1[j-1]

print "the values of conc. using CDS method for x=.84 to 1 are"
for i in range(42,50):
    print x[i-1]
    
#part (ii) using CDS and UDS method
#multipling with -delta_x**2/D we get
#-(Pe_local+1)*C(i-1) + (Pe_local+2)*C(i)-C(i+1)=0
BEta=Pe_local+2
Gamma=Pe_local+1
a=[0]
for i in range(2,50):
    a.append(-Gamma)            #sub diagonal assignment
b=[]
for j in range(1,50):
    b.append(BEta)           #main diagonal assignment
c=[]
for k in range(1,49):
    c.append(-1)            #super diagonal assignment

d=zeros(49)    
d[0]=Gamma*C_ini
d[-1]=C_end
for l in range(2,49):
    d[l-1]=0            #given values assignment
i=1#
n=49#
beta1=[b[i-1]]                #initial b is equal to beta since a1=0
gamma1=[d[i-1]/beta1[i-1]]      #since c7=0
m=i+1#
for j in range(m,n+1):
    beta1.append(b[j-1]-a[j-1]*c[j-2]/beta1[j-2])
    gamma1.append((d[j-1]-a[j-1]*gamma1[j-2])/beta1[j-1])
x=zeros(n)
x[-1]=gamma1[n-1]               #since c7=0
n1=n-i#
for k in range(1,n1+1):
    j=n-k# 
    x[j-1]=gamma1[j-1]-c[j-1]*x[j]/beta1[j-1]

print "the values of conc. using CDS/UDS method for x=.84 to 1 are"
for i in range(42,50):
    print x[i-1]

print "the analytical soln is"
C_an=zeros(50)
for i in range(42,50):
     C_an[i-1]=C_ini+((exp(Pe*.02*i)-1)*(C_end-C_ini)/(exp(Pe)-1))
     print C_an[i-1]    

the values of conc. using CDS method for x=.84 to 1 are
0.000152415790276
0.000457247370828
0.00137174211248
0.00411522633745
0.0123456790123
0.037037037037
0.111111111111
0.333333333333
the values of conc. using CDS/UDS method for x=.84 to 1 are
0.00390625
0.0078125
0.015625
0.03125
0.0625
0.125
0.25
0.5
the analytical soln is
0.000335462627903
0.000911881965555
0.00247875217667
0.00673794699909
0.0183156388887
0.0497870683679
0.135335283237
0.367879441171

Exa 6.2 Page 112

from __future__ import division
from numpy import zeros,exp
#given convective diffusive eqn--> -u*(dc/dx)+D*(d2C/dx2)=0
C_ini=0             #at x=0
C_end=1             #at x=1

#using central difference method for both diffusion and convective term 
#-u*(C(i+1)-C(i-1))/(2*delta_x) + D*(C(i+1)+C(i-1)-2*C(i))/delta_x**2 = 0
delta_x=1/50
#on solving the given eqns and by using the given boundary eqns we have
Pe=500                              #given 
Pe_local=500*delta_x                #u/D=50 as l=1
alpha=Pe_local-2                   #co-eff of C(i+1)
Beta=Pe_local+2                    #co-eff of C(i-1)
#multipling with -2*delta_x**2/D we get
#-(Pe_local+2)*C(i-1) + 4*C(i) + (Pe_local-2)*C(i+1)=0
#solving eqns using TDMA method
a=[0]
for i in range(2,50):
    a.append(-Beta)            #sub diagonal assignment
b=[]
for j in range(1,50):
    b.append(4)           #main diagonal assignment

c=[]
for k in range(1,49):
    c.append(alpha)#            #super diagonal assignment

d=zeros(49)
d[0]=Beta*C_ini
d[-1]=-alpha*C_end
for l in range(2,49):
    d[l-1]=0            #given values assignment
i=1#
n=49#
beta1=[b[i-1]]                #initial b is equal to beta since a1=0
gamma1 = [d[i-1]/beta1[i-1]]      #since c7=0
m=i+1#
for j in range(m,n+1):
    beta1.append(b[j-1]-a[j-1]*c[j-2]/beta1[j-2])
    gamma1.append((d[j-1]-a[j-1]*gamma1[j-2])/beta1[j-1])

x=zeros(n)
x[-1]=gamma1[n-1]               #since c7=0
n1=n-i#
for k in range(1,n1+1):
    j=n-k# 
    x[j-1]=gamma1[j-1]-c[j-1]*x[j]/beta1[j-1]

print "the values of conc. using CDS method for x=.84 to 1 are"
for i in range(42,50):
    print x[i-1]
    
#part (ii) using CDS and UDS method
#multipling with -delta_x**2/D we get
#-(Pe_local+1)*C(i-1) + (Pe_local+2)*C(i)-C(i+1)=0
BEta=Pe_local+2
Gamma=Pe_local+1
a=[0]
for i in range(2,50):
    a.append(-Gamma)            #sub diagonal assignment
b=[]
for j in range(1,50):
    b.append(BEta)           #main diagonal assignment
c=[]
for k in range(1,49):
    c.append(-1)            #super diagonal assignment

d=zeros(49)
d[0]=Gamma*C_ini
d[-1]=C_end
for l in range(2,49):
    d[l-1]=0            #given values assignment
i=1#
n=49#
beta1 = [b[i-1]]                #initial b is equal to beta since a1=0
gamma1 = [d[i-1]/beta1[i-1]]      #since c7=0
m=i+1#
for j in range(m,n+1):
    beta1.append(b[j-1]-a[j-1]*c[j-2]/beta1[j-2])
    gamma1.append((d[j-1]-a[j-1]*gamma1[j-2])/beta1[j-1])
x=zeros(n)
x[-1]=gamma1[n-1]               #since c7=0
n1=n-i#
for k in range(1,n1+1):
    j=n-k# 
    x[j-1]=gamma1[j-1]-c[j-1]*x[j]/beta1[j-1]
print "the values of conc. using CDS/UDS method for x=.84 to 1 are"
for i in range(42,50):
    print x[i-1]

print "the analytical soln is"
C_an=zeros(50)
for i in range(42,50):
    C_an[i-1]=C_ini+((exp(Pe*.02*i)-1)*(C_end-C_ini)/(exp(Pe)-1))
    print C_an[i-1]
        

the values of conc. using CDS method for x=.84 to 1 are
0.0390184408035
-0.0585276651261
0.0877914937683
-0.131687244573
0.197530862939
-0.296296298329
0.444444443573
-0.666666669281
the values of conc. using CDS/UDS method for x=.84 to 1 are
4.6650738021e-09
5.13158118231e-08
5.64473930054e-07
6.20921323059e-06
6.83013455365e-05
0.000751314800902
0.00826446280992
0.0909090909091
the analytical soln is
1.80485138785e-35
3.97544973591e-31
8.7565107627e-27
1.92874984796e-22
4.24835425529e-18
9.35762296884e-14
2.06115362244e-09
4.53999297625e-05