dx1/dt = -b2*x1 + b2*x3; dx2/dt = a2*x4 - a2*x2; dx3/dt = (ka*kc*b1*x1 - ka*kc*b1*x3 + ka*b1*x4*x1 - ka*b1*x4*x3 - n*kc*x3 + kc*b1*x1*x3 - kc*b1*x3^2) / (ka*kc + ka*x4 + kc*x3); dx4/dt = (-ka*n*x4 - ka*kc*a1*x4 + ka*kc*a1*x2 - ka*a1*x4^2 + ka*a1*x4*x2 - kc*a1*x4*x3 + kc*a1*x2*x3) / (ka*kc + ka*x4 + kc*x3); y1 = x4