CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   Main CFD Forum (https://www.cfd-online.com/Forums/main/)
-   -   The Blasius equation in Python (https://www.cfd-online.com/Forums/main/87947-blasius-equation-python.html)

javiergra24 May 3, 2011 14:49

The Blasius equation in Python
 
Hi everybody

I'm writing a script in Python to solve the Blasius equation but it does not work, numerical results does not match with data I've seen in fluid mechanics books.

Please anyone could help me? Thank you very much

Script is

import sys, pylab, numpy
from pylab import *
from numpy import *
import matplotlib.pyplot as plt


deta=0.0001
total=10


e = []
f = []
g = []
h = []


fvec=0; gvec=0; hvec=0.3219

etavec=0.0

while etavec <= 10:

fvec = fvec + gvec * deta
gvec = gvec + hvec * deta
hvec = hvec -(1/2) * fvec * hvec * deta
etavec = etavec + deta
e.append(etavec)
f.append(fvec)
g.append(gvec)
h.append(hvec)
print etavec,fvec,gvec

numpy.savetxt("blasius.dat", transpose((e,f,g,h)),delimiter=' ')

pylab.figure()
pylab.plot(e,g)
pylab.plot(e,h)
xticklines = getp(gca(), 'xticklines')
yticklines = getp(gca(), 'yticklines')
xgridlines = getp(gca(), 'xgridlines')
ygridlines = getp(gca(), 'ygridlines')
xticklabels = getp(gca(), 'xticklabels')
yticklabels = getp(gca(), 'yticklabels')
#lines1 = pylab.plot(e,f)
lines2 = pylab.plot(e,g)
lines3 = pylab.plot(e,h)
#pylab.setp(lines1, color='b', linewidth=3.0)
pylab.setp(lines2, color='r', linewidth=3.0)
pylab.setp(lines3, color='k', linewidth=3.0)
setp(xticklines, 'linewidth', 3)
setp(yticklines, 'linewidth', 3)
setp(xgridlines, 'linestyle', '--')
setp(ygridlines, 'linestyle', '--')
setp(yticklabels, fontsize='xx-large')
setp(xticklabels, fontsize='xx-large')
xlabel('$\eta$',fontsize=30)
ylabel('$f(\eta) \quad \partial_{\eta} f(\eta) \quad \partial_{\eta \eta}f(\eta)$',fontsize=30)
axis([0.0,10,-0.9,20.9])
#legend( (lines1, lines2, lines3), ('$\f(\eta)$', '$\partial_{\eta} f(\eta)$', '$\partial_{\eta \eta}f(\eta)$'))


grid(True)
pylab.show()


All times are GMT -4. The time now is 00:04.