#python #scipy #equation
#питон #сципи #уравнение
Вопрос:
Мне нужно использовать метод Эйлера в следующих уравнениях и получить реальное значение температуры, которое в этом упражнении 146.076550012783
.
В моей попытке, хотя я делаю что-то не так, и ответ неправильный,
import numpy as np import matplotlib.pyplot as plt from scipy.integrate import solve_ivp from scipy.integrate import odeint # PARÁMETROS DEL SISTEMA # Reactor FA0 = 80 # lb mol A/h (Flujo molar de la especie A) FB0 = 1000 # lb mol B/h (Flujo molar de la especie B) FM0 = 100 # lb mol M/h (Flujo molar de la especie M) V = 500*(0.133681) # ft3 (Volumen del reactor) dHrx = -36000 # BTU/lb mol A (Entalpía de reacción) CPA = 35 # BTU/lb mol ºF (Capacidad calorífica de la especie A) CPB = 18 # BTU/lb mol ºF (Capacidad calorífica de la especie B) CPC = 46 # BTU/lb mol ºF (Capacidad calorífica de la especie C) CPM = 19.5 # BTU/lb mol ºF (Capacidad calorífica de la especie M) T0 = 75 459.67 # F (Temperatura de entrada) C𝐴0 = 0.932 # lb mol/ft3 (Concentración molar de entrada de la especie A) CB0 = 3.45 # lb mol/ft3 (Concentración molar de entrada de la especie B) CM0 = 1.54 # lb mol/ft3 (Concentración molar de entrada de la especie M) R = 1.986 # BTU/lb mol R (Constante de los gases ideales) CAi = 0.14 # lb mol/ft3 (Concentración molar de la especie A en el reactor a tiempo cero) CBi = 2.93 # lb mol/ft3 (Concentración molar de la especie B en el reactor a tiempo cero) CCi = 0 # lb mol/ft3 (Concentración molar de la especie C en el reactor a tiempo cero) CMi = 0 # lb mol/ft3 (Concentración molar de la especie M en el reactor a tiempo cero) Ti = 160 459.67 # °F (Temperatura inicial del reactor) F = (FA0/CA0) (FB0/CB0) (FM0/CM0) # ft3/h (Flujo molar total) # Intercambiador de calor FR = 1000 # lb mol R/h (Flujo molar de agua) UA = 16000 # BTU/h °F (Coeficiente de transferencia de calor) CPR = 18 # BTU/lb mol °F (Capacidad calorífica del agua) TR = 60 459.67 # °F (Temperatura de entrada del agua) # Constante de velocidad de reacción def k(T): return 16.96*(10**12)*np.exp(-3200/(R*T)) # Velocidad de reacción de primer orden con respecto a la especie A def RA(T,CA): return -k(T)*CA # Velocidad de enfriamiento del reactor por el intercambiador de calor def Q(T): return FR*CPR*(TR-T)*(1 - np.exp((-UA)/(FR*CPR))) # ECUACIONES DIFERENCIALES def balances(x,t): CA,CB,CC,CM,T = x dCA = FA0/V - (CA*F)/V RA(T,CA) dCB = FB0/V - (CB*F)/V RA(T,CA) dCC = (-CC*F)/V - RA(T,CA) dCM = FM0 - (CM*F)/V dT = ((Q(T)/V) (((T0 - T) * ((FA0*CPA) (FB0*CPB) (FM0*CPM)))/V) (dHrx*RA(T,CA)))/((CA*CPA) (CB*CPB) (CC*CPC) (CM*CPM)) return [dCA, dCB, dCC, dCM, dT] #Condiciones iniciales CI = [CAi, CBi, CCi, CMi, Ti] #Tiempo en horas t_inicial = 0.000 t_final = 4.000 t = np.linspace(t_inicial, t_final, 4000) sol = odeint(balances,CI,t) CA,CB,CC,CM,T = sol.transpose() plt.figure(1) plt.plot(t,CA) plt.show plt.figure(2) plt.plot(t,CB) plt.figure(3) plt.plot(t,CC) plt.figure(4) plt.plot(t,CM) plt.figure(5) plt.plot(t,T - 459.67) TF = T - 459.67