import numpy as np
import scipy as sp
import matplotlib.pyplot as plt
from scipy.optimize import fsolve
# Загрузка данных
increase = np.loadtxt("increase.txt", delimiter=" ", dtype=np.float)
decrease = np.loadtxt("decrease.txt", delimiter=" ", dtype=np.float)
dataIncrease = [[], []]
dataDecrease = [[], []]
for element in increase:
dataIncrease[0].append(element[0])
dataIncrease[1].append(element[1] * 133.3)
for element in decrease:
dataDecrease[0].append(element[0])
dataDecrease[1].append(element[1] * 133.3)
dataIncrease = np.array(dataIncrease)
dataDecrease = np.array(dataDecrease)
delta_P = 6.66 # Погрешность измерения давления (в Паскалях)
delta_T = 0.1 # Погрешность измерения температуры
# Аппроксимация
approxDegree = 2 # Степень многочлена, которым аппроксимируем. Два замечательно подходит
fp, residuals, rank, sv, rcond = sp.polyfit(dataIncrease[0], dataIncrease[1], approxDegree, full=True)
f = sp.poly1d(fp)
print(f) # Выведет многочлен со всеми коэффициентами
# Вывод графика
fig, ax = plt.subplots()
ax.plot(dataIncrease[0], dataIncrease[1], marker='.', markersize=5, markerfacecolor='#2980b9', color='#2980b9', ls='None')
ax.plot(dataDecrease[0], dataDecrease[1], marker='.', markersize=5, markerfacecolor='#e74c3c', color='#e74c3c', ls='None')
ax.errorbar(dataIncrease[0], dataIncrease[1], xerr=delta_T, yerr=delta_P, linestyle='None', color='#2980b9') # Кресты
ax.errorbar(dataDecrease[0], dataDecrease[1], xerr=delta_T, yerr=delta_P, linestyle='None', color='#e74c3c') # погрешностей
ax.plot(dataIncrease[0], f(dataIncrease[0]), color='#000000', linewidth=2) # Рисуем аппроксимационную кривую
ax.set_xlabel('Температура, K')
ax.set_ylabel('Разность давлений, Па')
ax.text(300.5, 7000, 'Повышение', color='#2980b9')
ax.text(298.5, 8000, 'Понижение', color='#e74c3c')
ax.grid()
plt.show()
# Вычисление итогового значения
def d(x):
""" Производная многочлена, который вернул нам SciPy """
return 13.712 * x - 3876
# Значение при 25 градусах Цельсия
print("T = " + str(dataIncrease[0][5]))
print("L = " + str(((8.31 * (dataIncrease[0][5] ** 2)) / dataIncrease[1][5]) * d(dataIncrease[0][5])))
