#!/usr/bin/env python
#coding=utf8
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import numpy as np
from math import *
alfa=21.0#Угол сужения
beta=11.5#Угол расширения
rkr=1.1#Радиус критического сечения
R0=2*rkr
r1=0.5*rkr# Радиус округления сужающейся части сопла
r2=0.8*rkr# Радиус округления расширяющейся части сопла
ye=rkr+r2
L=1.2*R0#Длина прямого участка сопла Лаваля
x0=0
y0=R0
xa=L
ya=y0
xc1=xa
yc1=ya-r1
xc=xa+r1*cos(radians(90-alfa))
yc=yc1+r1*sin(radians(90-alfa))
yd=ye-r2*sin(radians(90-alfa))
xd=xc+(yc-yd)/tan(radians(alfa))
xc2=xd+r2*sin(radians(alfa))
xe=xc2
xf=xe+r2*cos(radians(90-beta))
yf=ye-r2*sin(radians(90-beta))
def R(x):
if x0<=x<=xa:
return ya*1e-4
if xa<=x<=xc:
return (yc1+sqrt(r1**2-(x-xc1)**2))*1e-4
if xc<=x<=xd:
return (-tan(radians(alfa))*(x-xc)+yc)*1e-4
if xd<=x<=xf:
return (ye-sqrt(r2**2-(x-xe)**2))*1e-4
if xf<=x:
return (yf+tan(radians(beta))*(x-xf))*1e-4
y=[R(x+xe) for x in np.arange(-5,5,0.01)]
x=np.arange(-5,5,0.01)
plt.figure()
plt.title('Профиль сопла ')
plt.axis([-5.0, 5.0, 0.0, 3.0*10**-4])
plt.plot(x,y,'r')
plt.grid(True)
plt.figure()
plt.title('Изменение площади проходного \n сечения сопла вдоль его продольной оси ')
yy=[pi*R(x+xe)**2 for x in np.arange(-5,5,0.01)]
plt.plot(x,yy,'r')
plt.grid(True)
plt.show()
def lamda(z):
m=round(Q(z),2)
if z>= 0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,1.5)
return x[0]
if z<0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,0.5)
return x[0]
from sympy import *
import time
x = symbols('x',real=True)
start = time.time()
start = time.time()
d=solve( 1.5774409656148784068*x *(1.0-0.16666666666666666667*x**2)**2.5-0.25,x)
stop = time.time()
print ("Время работы решателя:",round(stop-start,3))
print(round(d[0],3))
print(round(d[1],3))
#!/usr/bin/env python
#coding=utf8
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import numpy as np
from math import *
from scipy.optimize import *
import time
start = time.time()
alfa=21.0#Угол сужения
beta=11.5#Угол расширения
rkr=1.1#Радиус критического сечения
R0=2*rkr
r1=0.5*rkr#Радиус округления сужающейся части сопла
r2=0.8*rkr#Радиус округления расширяющейся части сопла
ye=rkr+r2
L=1.2*R0#Длина прямого участка сопла Лаваля
x0=0
y0=R0
xa=L
ya=y0
xc1=xa
yc1=ya-r1
xc=xa+r1*cos(radians(90-alfa))
yc=yc1+r1*sin(radians(90-alfa))
yd=ye-r2*sin(radians(90-alfa))
xd=xc+(yc-yd)/tan(radians(alfa))
xc2=xd+r2*sin(radians(alfa))
xe=xc2
xf=xe+r2*cos(radians(90-beta))
yf=ye-r2*sin(radians(90-beta))
def R(x):
if x0<=x<=xa:
return ya*1e-4
if xa<=x<=xc:
return (yc1+sqrt(r1**2-(x-xc1)**2))*1e-4
if xc<=x<=xd:
return (-tan(radians(alfa))*(x-xc)+yc)*1e-4
if xd<=x<=xf:
return (ye-sqrt(r2**2-(x-xe)**2))*1e-4
if xf<=x:
return (yf+tan(radians(beta))*(x-xf))*1e-4
def S(x):
return pi*R(x+xe)**2
xg=2*xe
n=150
RG=287 #Газовая постоянная
Tt=611 #Температура торможения
k=1.4
def tau(x):
return 1-(1/6)*x**2
def eps(x):
return (1-(1/6)*x**2)**2.5
def pip(x):
return 1-(1/6)*x**2
def Q(z):
return S(0)/S(z)
x=[x0-xe+(i/n)*(xg-x0) for i in np.arange(0,n,1)]
def lamda(z):
m=round(Q(z),2)
if z>= 0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,1.5)
return x[0]
if z<0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,0.5)
return x[0]
plt.title('Газодинамическая функция приведенной скорости')
y=[lamda(z) for z in x]
stop = time.time()
print ("Время работы программы:",round(stop-start,3))
plt.ylabel('?(xi)-отношение скорости газа к скорости звука' )
plt.plot(0, 1, color='b', linestyle=' ', marker='o', label='Сужение сопла')
plt.xlabel('xi -координата вдоль оси сопла')
plt.plot(x,y,'r')
plt.legend(loc='best')
plt.grid(True)
plt.show()
#!/usr/bin/env python
#coding=utf8
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import numpy as np
from math import *
from scipy.optimize import *
import time
start = time.time()
alfa=21.0#Угол сужения
beta=11.5#Угол расширения
rkr=1.1#Радиус критического сечения
R0=2*rkr
r1=0.5*rkr#Радиус округления сужающейся части сопла
r2=0.8*rkr#Радиус округления расширяющейся части сопла
ye=rkr+r2
L=1.2*R0#Длина прямого участка сопла Лаваля
x0=0
y0=R0
xa=L
ya=y0
xc1=xa
yc1=ya-r1
xc=xa+r1*cos(radians(90-alfa))
yc=yc1+r1*sin(radians(90-alfa))
yd=ye-r2*sin(radians(90-alfa))
xd=xc+(yc-yd)/tan(radians(alfa))
xc2=xd+r2*sin(radians(alfa))
xe=xc2
xf=xe+r2*cos(radians(90-beta))
yf=ye-r2*sin(radians(90-beta))
def R(x):
if x0<=x<=xa:
return ya*1e-4
if xa<=x<=xc:
return (yc1+sqrt(r1**2-(x-xc1)**2))*1e-4
if xc<=x<=xd:
return (-tan(radians(alfa))*(x-xc)+yc)*1e-4
if xd<=x<=xf:
return (ye-sqrt(r2**2-(x-xe)**2))*1e-4
if xf<=x:
return (yf+tan(radians(beta))*(x-xf))*1e-4
def S(x):
return pi*R(x+xe)**2
xg=2*xe
n=150
RG=287 #Газовая постоянная
Tt=611 #Температура торможения
k=1.4
def tau(x):
return 1-(1/6)*x**2
def eps(x):
return (1-(1/6)*x**2)**2.5
def pip(x):
return 1-(1/6)*x**2
def Q(z):
return S(0)/S(z)
x=[x0-xe+(i/n)*(xg-x0) for i in np.arange(0,n,1)]
def lamda(z):
m=round(Q(z),2)
if z>= 0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,1.5)
return x[0]
if z<0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,0.5)
return x[0]
plt.title('Газодинамическая функция температуры')
t0=293
y=[ t0*tau(lamda(z)) for z in x]
stop = time.time()
print (" Время работы программы:",round(stop-start,3))
plt.ylabel('T(xi)-температура газового потока град. С' )
plt.xlabel('xi -координата вдоль оси сопла')
plt.plot(x,y,'r')
plt.grid(True)
plt.show()
#!/usr/bin/env python
#coding=utf8
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import numpy as np
from math import *
from scipy.optimize import *
import time
start = time.time()
alfa=21.0#Угол сужения
beta=11.5#Угол расширения
rkr=1.1#Радиус критического сечения
R0=2*rkr
r1=0.5*rkr#Радиус округления сужающейся части сопла
r2=0.8*rkr#Радиус округления расширяющейся части сопла
ye=rkr+r2
L=1.2*R0#Длина прямого участка сопла Лаваля
x0=0
y0=R0
xa=L
ya=y0
xc1=xa
yc1=ya-r1
xc=xa+r1*cos(radians(90-alfa))
yc=yc1+r1*sin(radians(90-alfa))
yd=ye-r2*sin(radians(90-alfa))
xd=xc+(yc-yd)/tan(radians(alfa))
xc2=xd+r2*sin(radians(alfa))
xe=xc2
xf=xe+r2*cos(radians(90-beta))
yf=ye-r2*sin(radians(90-beta))
def R(x):
if x0<=x<=xa:
return ya*1e-4
if xa<=x<=xc:
return (yc1+sqrt(r1**2-(x-xc1)**2))*1e-4
if xc<=x<=xd:
return (-tan(radians(alfa))*(x-xc)+yc)*1e-4
if xd<=x<=xf:
return (ye-sqrt(r2**2-(x-xe)**2))*1e-4
if xf<=x:
return (yf+tan(radians(beta))*(x-xf))*1e-4
def S(x):
return pi*R(x+xe)**2
xg=2*xe
n=150
RG=287 #Газовая постоянная
Tt=611 #Температура торможения
k=1.4
def tau(x):
return 1-(1/6)*x**2
def eps(x):
return (1-(1/6)*x**2)**2.5
def pip(x):
return 1-(1/6)*x**2
def Q(z):
return S(0)/S(z)
x=[x0-xe+(i/n)*(xg-x0) for i in np.arange(0,n,1)]
def lamda(z):
m=round(Q(z),2)
if z>= 0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,1.5)
return x[0]
if z<0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,0.5)
return x[0]
plt.title('Газодинамическая функция давления')
p0=3.6
y=[ 3.6*pip(lamda(z)) for z in x]
stop = time.time()
print ("Время работы программы:",round(stop-start,3))
plt.ylabel('P(xi)-давление газового потока в МПа' )
plt.xlabel('xi -координата вдоль оси сопла')
plt.plot(x,y,'r')
plt.grid(True)
plt.show()
#!/usr/bin/env python
#coding=utf8
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import numpy as np
from math import *
from scipy.optimize import *
import time
start = time.time()
alfa=21.0#Угол сужения
beta=11.5#Угол расширения
rkr=1.1#Радиус критического сечения
R0=2*rkr
r1=0.5*rkr#Радиус округления сужающейся части сопла
r2=0.8*rkr#Радиус округления расширяющейся части сопла
ye=rkr+r2
L=1.2*R0#Длина прямого участка сопла Лаваля
x0=0
y0=R0
xa=L
ya=y0
xc1=xa
yc1=ya-r1
xc=xa+r1*cos(radians(90-alfa))
yc=yc1+r1*sin(radians(90-alfa))
yd=ye-r2*sin(radians(90-alfa))
xd=xc+(yc-yd)/tan(radians(alfa))
xc2=xd+r2*sin(radians(alfa))
xe=xc2
xf=xe+r2*cos(radians(90-beta))
yf=ye-r2*sin(radians(90-beta))
def R(x):
if x0<=x<=xa:
return ya*1e-4
if xa<=x<=xc:
return (yc1+sqrt(r1**2-(x-xc1)**2))*1e-4
if xc<=x<=xd:
return (-tan(radians(alfa))*(x-xc)+yc)*1e-4
if xd<=x<=xf:
return (ye-sqrt(r2**2-(x-xe)**2))*1e-4
if xf<=x:
return (yf+tan(radians(beta))*(x-xf))*1e-4
def S(x):
return pi*R(x+xe)**2
xg=2*xe
n=150
RG=287 #Газовая постоянная
Tt=611 #Температура торможения
k=1.4
def tau(x):
return 1-(1/6)*x**2
def eps(x):
return (1-(1/6)*x**2)**2.5
def pip(x):
return 1-(1/6)*x**2
def Q(z):
return S(0)/S(z)
x=[x0-xe+(i/n)*(xg-x0) for i in np.arange(0,n,1)]
def lamda(z):
m=round(Q(z),2)
if z>= 0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,1.5)
return x[0]
if z<0:
x=fsolve(lambda x:x*( (1-(1/6)*x**2)**2.5)/((1-(1/6))**2.5)-m,0.5)
return x[0]
plt.title('Газодинамическая функция относительной плотности газа')
y=[ eps(lamda(z)) for z in x]
stop = time.time()
print ("Время работы программы:",round(stop-start,3))
plt.ylabel('Относительная плотность' )
plt.xlabel('xi -координата вдоль оси сопла')
plt.plot(x,y,'r')
plt.grid(True)
plt.show()
К сожалению, не доступен сервер mySQL