title: solving system of coupled odes with odeint
body: <p>I'm using a system of ode's to model coffee bean roasting for a class assignment. The equations are below.</p>
  
  <p><img src="https://i.stack.imgur.com/TDGOi.png" alt="Systems of Equations"></p>
  
  <p>The parameters (other than X_b and T_b) are all constants. </p>
  
  <p>When I try to use odeint to solve this system, it gives a constant T_b and X_b profile (which conceptually doesn't make sense). </p>
  
  <p>Below is the code I'm using</p>
  
  <pre><code>from scipy.integrate import odeint
  import numpy as np
  import matplotlib.pyplot as plt
  
  # Write function for bean temperature T_b differential equation
  def deriv(X,t):
      T_b, X_b = X 
      dX_b = (-4.32*10**9*X_b**2)/(l_b**2)*np.exp(-9889/T_b)
      dT_b = ((h_gb*A_gb*(T_gi - T_b))+(m_b*A_arh*np.exp(-H_a/R_g/T_b))+ 
      (m_b*lam*dX_b))/(m_b*(1.099+0.0070*T_b+5*X_b)*1000)
      return [dT_b, dX_b]
  
  # Establish initial conditions
  t = 0 #seconds
  T_b = 298 # degrees K
  X_b = 0.1 # mass fraction of moisture
  
  # Set time step
  dt = 1 # second
  
  # Establish location to store data
  history = [[t,T_b, X_b]]
  
  # Use odeint to solve DE
  while t &lt; 600:
    T_b, X_b = odeint(deriv, [T_b, X_b], [t+dt])[-1]
    t += dt
    history.append([t,T_b, X_b])
  
  # Plot Results
  def plot_history(history, labels):
    """Plots a simulation history."""
    history = np.array(history)
    t = history[:,0]
    n = len(labels) - 1
    plt.figure(figsize=(8,1.95*n))
    for k in range(0,n):
      plt.subplot(n, 1, k+1)
      plt.plot(t, history[:,k+1])
      plt.title(labels[k+1])
      plt.xlabel(labels[0])
      plt.grid()
    plt.tight_layout()
  
  
  plot_history(history, ['t (s)','Bean Temperature $T_b$ (K)', 'Bean Moisture Content $X_b$'])
  plt.show()
  </code></pre>
  
  <p>Do you have any ideas why the integration step isn't working?</p>
  
  <p>Thank You!!</p>