def fp(L):
x,A, L, w, t, k =symbols('x A L w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate(dx2,(x,0,L)))/(A**2*w**2*cos(t*w)**2)
def fh(h):
x,A, h, w, t, k =symbols(' x A h w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate (dx2,(x,0,h)))/(A**2*w**2*cos(t*w)**2)
def pl(L):
x,A, J, w, t, k,L =symbols('x A J w t k L ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(x))**2
return factor(integrate (dx2,(x,0,L)))/(A**2*sin(t*w)**2)
from sympy import *
from numpy import arange,pi
import matplotlib.pyplot as plt
import time
start = time.time()
def fp(L):
x,A, L, w, t, k =symbols('x A L w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate(dx2,(x,0,L)))/(A**2*w**2*cos(t*w)**2)
def fh(h):
x,A, h, w, t, k =symbols(' x A h w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate (dx2,(x,0,h)))/(A**2*w**2*cos(t*w)**2)
def pl(L):
x,A, J, w, t, k,L =symbols('x A J w t k L ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(x))**2
return factor(integrate (dx2,(x,0,L)))/(A**2*sin(t*w)**2)
""" Установка численных значений символьных переменных """
L,h,m,m0,E,J,q =symbols(' L,h m0 m E J q ')
k1=q.subs({q:0.16})
k2=pl(L).subs({L:1})
k3=E.subs({E:196e9})*J.subs({J:2.3e-08})
k4=m0.subs({m0:1.005})
k5=m.subs({m:0.98})
k6=fp(L).subs({L:1})
"""Решения для частот чувствительностей к уровню и траекторий конца ттрубки """
x=arange(0.0,1.0,0.01)
f=k1*((k2+k3)/(k4*k6+k5*fh(h)))**0.5
y=[k1*((k2+k3)/(k4*k6+k5*fh(h).subs({h:w})))**0.5 for w in arange(0.0,1.0,0.01)]
s=f.diff(h)
y1=[s.subs({h:w}) for w in arange(0.0,1.0,0.01)]
k3=E.subs({E:196e9})*J.subs({J:1.7e-08})
f1=k1*((k2+k3)/(k4*k6+k5*fh(h)))**0.5
y2=[k1*((k2+k3)/(k4*k6+k5*fh(h).subs({h:w})))**0.5 for w in arange(0.0,1.0,0.01)]
s1=f1.diff(h)
y3=[s1.subs({h:w}) for w in arange(0.0,1.0,0.01)]
k7=fh(h).subs({h:0.7})
k3=E.subs({E:196e9})*J.subs({J:2.3e-08})
f1=k1*((k2+k3)/(k4*k6+k5*k7))**0.5
k3=E.subs({E:196e9})*J.subs({J:1.7e-08})
f2=k1*((k2+k3)/(k4*k6+k5*k7))**0.5
x1=[sin(f1*t) for t in arange(0.,2*pi,0.01)]
y4=[sin(f2*t) for t in arange(0.,2*pi,0.01)]
""" Построение графиков"""
plt.subplot(221)
plt.plot(x, y, label='f1-zox')
plt.plot(x, y2,label='f2- zoy')
plt.xlabel(' h')
plt.ylabel(' f 1,f2')
plt.legend(loc='best')
plt.grid(True)
plt.subplot(222)
plt.plot(x, y1,label='s=df/dh zox')
plt.plot(x, y3,label='s1=df1/dh - zoy')
plt.ylabel('s')
plt.xlabel(' h ')
plt.legend(loc='best')
plt.grid(True)
plt.subplot(223)
plt.plot(x1, y4,label='X,Y для h=0.7')
plt.ylabel(' Y')
plt.xlabel('X')
plt.grid(True)
plt.legend(loc='best')
stop = time.time()
print ("Время :",round(stop-start,3))
plt.show()
from sympy import *
from numpy import arange,pi
import matplotlib.pyplot as plt
import time
start = time.time()
def fp(L):
x,A, L, w, t, k =symbols('x A L w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate(dx2,(x,0,L)))/(A**2*w**2*cos(t*w)**2)
def fh(h):
x,A, h, w, t, k =symbols(' x A h w t k ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(t))**2
return factor(integrate (dx2,(x,0,h)))/(A**2*w**2*cos(t*w)**2)
def pl(L):
x,A, J, w, t, k,L =symbols('x A J w t k L ')
z=(sin(k*x)-sinh(k*x)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*x)-cosh(k *x))).subs({k:1.875})
dx=z*A*sin(w*t)
dx2=(dx.diff(x))**2
return factor(integrate (dx2,(x,0,L)))/(A**2*sin(t*w)**2)
L,h,m,m0,E,J,q =symbols(' L,h m0 m E J q ')
f=q*((pl(L)+E*J)/(m0*fp(L)+m*fh(h)))**0.5
s=f.diff(h)
print(fp(L))
print(fh(h))
print(pl(L))
print(s)
stop = time.time()
print ("Время :",round(stop-start,3))
from numpy import (pi,cos,cosh,sin,sinh,arange)
import matplotlib.pyplot as plt
import matplotlib as mpl
mpl.rcParams['font.family'] = 'fantasy'
mpl.rcParams['font.fantasy'] = 'Comic Sans MS, Arial'
import time
start = time.time()
L=1# длинна консольной трубки в м.
d=1e-3# толщина стенки в м.
e=0.8# эксцентриситет сечения б.р.
a0=20e-3# большая внутренняя полуось в м.
a=a0+d# большая наружная полуось в м.
b0=a0*e# малая внутренняя полуось в м.
b=b0+d#малая наружная полуось в м.
E=196e9# модуль Юнга материала в н/м2.
rm=7.9e3# массовая плотность материала трубки в кг/м3.
rg=1e3#массовая плотность жидкости в кг/м3.
J=(pi/4)*(a**3*b-a0**3*b0)# статический момент инерции по малой оси в м4.
J1=(pi/4)*(a*b**3-a0*b0**3)#статический момент инерции по большй оси в м4.
m0=pi*(a*b-a0*b0)*rm# масса единицы длинны трубки в кг/м3
m=pi*(a0*b0)*rg#масса единицы длинны жидкости в кг/м3
q=1/(2*pi)# коэффициент для перевода круговой частоты в циклическую
def fp(L):# функция кинитической энергии пустой консольной трубки L
return (1.83012983499975*L*sin(1.875*L)**2 + 1.83012983499975*L*cos(1.875*L)**2 - 0.830129834999752*L*sinh(1.875*L)**2 + 0.830129834999752*L*cosh(1.875*L)**2 - 0.869882798951083*sin(1.875*L)**2 + 0.442735911999868*sin(1.875*L)*cos(1.875*L) + 1.73976559790217*sin(1.875*L)*sinh(1.875*L) - 1.9521384906664*sin(1.875*L)*cosh(1.875*L) - 0.885471823999735*cos(1.875*L)*sinh(1.875*L) + 0.976069245333201*sinh(1.875*L)*cosh(1.875*L) - 0.869882798951083*cosh(1.875*L)**2 + 0.869882798951083)
def fh(h):# функция кинитической энергии жидкости уровня h
return (1.83012983499975*h*sin(1.875*h)**2 + 1.14978095631541e-15*h*sin(1.875*h)*sinh(1.875*h) + 1.83012983499975*h*cos(1.875*h)**2 - 1.07791964654569e-15*h*cos(1.875*h)*cosh(1.875*h) - 0.830129834999752*h*sinh(1.875*h)**2 + 0.830129834999752*h*cosh(1.875*h)**2 - 0.869882798951083*sin(1.875*h)**2 + 0.442735911999868*sin(1.875*h)*cos(1.875*h) + 1.73976559790217*sin(1.875*h)*sinh(1.875*h) - 1.9521384906664*sin(1.875*h)*cosh(1.875*h) - 0.885471823999735*cos(1.875*h)*sinh(1.875*h) + 0.976069245333201*sinh(1.875*h)*cosh(1.875*h) - 0.869882798951083*cosh(1.875*h)**2 + 0.869882798951083)
def pl(L):# функция потенциальной энергии консольной трубки жёсткостью E*J/L**3
return (6.434050201171*L*sin(1.875*L)**2 + 1.48843545191282e-15*L*sin(1.875*L)*sinh(1.875*L) + 6.434050201171*L*cos(1.875*L)**2 - 2.79081647233653e-15*L*cos(1.875*L)*cosh(1.875*L) + 2.918425201171*L*sinh(1.875*L)**2 - 2.918425201171*L*cosh(1.875*L)**2 + 3.0581817150624*sin(1.875*L)**2 - 1.55649344062453*sin(1.875*L)*cos(1.875*L) + 3.11298688124907*sin(1.875*L)*cosh(1.875*L) - 6.86298688124907*cos(1.875*L)*sinh(1.875*L) + 6.1163634301248*cos(1.875*L)*cosh(1.875*L) - 3.0581817150624*sinh(1.875*L)**2 + 3.43149344062453*sinh(1.875*L)*cosh(1.875*L) - 6.1163634301248)
def f(m0,q,E,J,h,m):# чувствительность уровнимера s=df/dh
return -0.5*m0*q*((E*J + 6.434050201171*L*sin(1.875*L)**2 + 6.434050201171*L*cos(1.875*L)**2 + 2.918425201171*L*sinh(1.875*L)**2 - 2.918425201171*L*cosh(1.875*L)**2 + 3.0581817150624*sin(1.875*L)**2 - 1.55649344062453*sin(1.875*L)*cos(1.875*L) + 3.11298688124907*sin(1.875*L)*cosh(1.875*L) - 6.86298688124907*cos(1.875*L)*sinh(1.875*L) + 6.1163634301248*cos(1.875*L)*cosh(1.875*L) - 3.0581817150624*sinh(1.875*L)**2 + 3.43149344062453*sinh(1.875*L)*cosh(1.875*L) - 6.1163634301248)/(m*(1.83012983499975*L*sin(1.875*L)**2 + 1.83012983499975*L*cos(1.875*L)**2 - 0.830129834999752*L*sinh(1.875*L)**2 + 0.830129834999752*L*cosh(1.875*L)**2 - 0.869882798951083*sin(1.875*L)**2 + 0.442735911999868*sin(1.875*L)*cos(1.875*L) + 1.73976559790217*sin(1.875*L)*sinh(1.875*L) - 1.9521384906664*sin(1.875*L)*cosh(1.875*L) - 0.885471823999735*cos(1.875*L)*sinh(1.875*L) - 0.869882798951083*sinh(1.875*L)**2 + 0.976069245333201*sinh(1.875*L)*cosh(1.875*L)) + m0*(1.83012983499975*h*sin(1.875*h)**2 + 1.83012983499975*h*cos(1.875*h)**2 - 0.830129834999752*h*sinh(1.875*h)**2 + 0.830129834999752*h*cosh(1.875*h)**2 - 0.869882798951083*sin(1.875*h)**2 + 0.442735911999868*sin(1.875*h)*cos(1.875*h) + 1.73976559790217*sin(1.875*h)*sinh(1.875*h) - 1.9521384906664*sin(1.875*h)*cosh(1.875*h) - 0.885471823999735*cos(1.875*h)*sinh(1.875*h) - 0.869882798951083*sinh(1.875*h)**2 + 0.976069245333201*sinh(1.875*h)*cosh(1.875*h))))**0.5*(1.0*sin(1.875*h)**2 - 3.26206049606656*sin(1.875*h)*cos(1.875*h) - 2.0*sin(1.875*h)*sinh(1.875*h) + 3.26206049606656*sin(1.875*h)*cosh(1.875*h) + 2.6602596699995*cos(1.875*h)**2 + 3.26206049606656*cos(1.875*h)*sinh(1.875*h) - 5.32051933999901*cos(1.875*h)*cosh(1.875*h) + 1.0*sinh(1.875*h)**2 - 3.26206049606656*sinh(1.875*h)*cosh(1.875*h) + 2.6602596699995*cosh(1.875*h)**2)/(m*(1.83012983499975*L*sin(1.875*L)**2 + 1.83012983499975*L*cos(1.875*L)**2 - 0.830129834999752*L*sinh(1.875*L)**2 + 0.830129834999752*L*cosh(1.875*L)**2 - 0.869882798951083*sin(1.875*L)**2 + 0.442735911999868*sin(1.875*L)*cos(1.875*L) + 1.73976559790217*sin(1.875*L)*sinh(1.875*L) - 1.9521384906664*sin(1.875*L)*cosh(1.875*L) - 0.885471823999735*cos(1.875*L)*sinh(1.875*L) - 0.869882798951083*sinh(1.875*L)**2 + 0.976069245333201*sinh(1.875*L)*cosh(1.875*L)) + m0*(1.83012983499975*h*sin(1.875*h)**2 + 1.83012983499975*h*cos(1.875*h)**2 - 0.830129834999752*h*sinh(1.875*h)**2 + 0.830129834999752*h*cosh(1.875*h)**2 - 0.869882798951083*sin(1.875*h)**2 + 0.442735911999868*sin(1.875*h)*cos(1.875*h) + 1.73976559790217*sin(1.875*h)*sinh(1.875*h) - 1.9521384906664*sin(1.875*h)*cosh(1.875*h) - 0.885471823999735*cos(1.875*h)*sinh(1.875*h) - 0.869882798951083*sinh(1.875*h)**2 + 0.976069245333201*sinh(1.875*h)*cosh(1.875*h)))
"""Графический анализ результатов"""
x=arange(0.0,1.0,0.01)
y1=[q*((pl(L)+E*J)/(m0*fp(L)+m*fh(h)))**0.5 for h in arange(0.0,1.0,0.01)] #числовая функция частоты трубки в сечении ZOX
J=(pi/4)*(a*b**3-a0*b0**3)#момент инерции в плоскости ZOY
y2=[q*((pl(L)+E*J)/(m0*fp(L)+m*fh(h)))**0.5 for h in arange(0.0,1.0,0.01)]#числовая функция частоты трубки в сечении ZOX
J=(pi/4)*(a**3*b-a0**3*b0)#момент инерции в плоскости ZOX
y3=[f(m0,q,E,J,h,m) for h in arange(0.0,1.0,0.01)]#числовая функция чувствительности трубки в сечении ZOX
J=(pi/4)*(a*b**3-a0*b0**3)#момент инерции в плоскости ZOY
y4=[f(m0,q,E,J,h,m) for h in arange(0.0,1.0,0.01)]#числовая функция чувствительности трубки в сечении ZOX
J=(pi/4)*(a**3*b-a0**3*b0)#момент инерции в плоскости ZOX
h=0.7
f1=q*((pl(L)+E*J)/(m0*fp(L)+m*fh(h)))**0.5 # частота колебаний для h=0.7 в плоскости ZOX
J=(pi/4)*(a*b**3-a0*b0**3)#момент инерции в плоскости ZOY
f2=q*((pl(L)+E*J)/(m0*fp(L)+m*fh(h)))**0.5 # частота колебаний для h=0.7 в плоскости ZOY
x1=[sin(f1*t) for t in arange(0.,2*pi,0.01)]#траектория конца трубки по Х
y5=[sin(f2*t) for t in arange(0.,2*pi,0.01)]#траектория конца трубки по Y
""" Построение графиков"""
plt.subplot(221)
plt.plot(x, y1, label='f-ZOX')
plt.plot(x, y2,label='f1-ZOY')
plt.xlabel(' h')
plt.ylabel(' f 1,f2')
plt.legend(loc='best')
plt.grid(True)
plt.subplot(222)
plt.plot(x, y3,label='s1=df1/dh - ZOX')
plt.plot(x, y4,label='s2=df2/dh - ZOY')
plt.ylabel('s')
plt.xlabel(' h ')
plt.legend(loc='best')
plt.grid(True)
plt.subplot(223)
plt.plot(x1, y5,label='X,Y для h=0.7')
plt.ylabel(' Y')
plt.xlabel('X')
plt.grid(True)
plt.legend(loc='best')
stop = time.time()
print ("Время :",round(stop-start,3))
plt.show()
#!/usr/bin/env python
#coding=utf8
import pylab
import numpy
from numpy import (zeros,arange,ones,pi,sin,cos,sinh,cosh)
from matplotlib.colors import LinearSegmentedColormap
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D
def trubka():
a=10;b=10;k=q[0]
v = arange(0, 2.05*pi, 0.05*pi)
u= zeros([len(v),1])
for i in arange(0,len(v)):
u[i,0]=[sin(k*w)-sinh(k*w)+((cosh(k) -cos(k))/(sin(k)-sinh(k)))*(cos(k*w)-cosh(k *w) )for w in arange(0, 2.05*pi, 0.05*pi)][i]
x=a*u*cos(v)
y=b*u*sin(v)
z=u*ones(len(v))
return x,y,z
x,y,z=trubka()
fig = pylab.figure()
axes = Axes3D(fig)
axes.plot_surface(x, y, z, rstride=4, cstride=4, cmap = cm.jet)
pylab.show()
К сожалению, не доступен сервер mySQL