/** # Batch Reactors Collection of ideal reactor functions for the operator splitting solution of reactive systems ([Cuoci et al. 2013](#cuoci2013numerical)). The temperature and chemical species concentrations are considered perfectly mixed inside each cell. Neglecting the convective and diffusive transport in the conservation equations, the equations reduce to the ODE system of a batch reactor: $$ \begin{cases} \dfrac{d\omega_i}{dt} = \dfrac{1}{\rho} \sum \limits_{j=1}^{NR} \nu_{i,j} r_j MW_i \\ \dfrac{d T}{dt} = \dfrac{1}{\rho c_p} \sum \limits_{j=1}^{NR} r_j \Delta H_{R_j} \end{cases} $$ where $\nu$ is the matrix of the stoichiometric coefficients, while $r_j$ is the rate of the $j$-th reaction ([Fogler 2020](#fogler2020elements)). The term $\Delta H_R$ is the enthalpy change associated with the chemical reactions. The reaction rate $r_j$ is usually computed as: $$ r_j = k_j \prod \limits_{i=1}^{NS} c_i^{\nu_{i,j}} $$ and the kinetic constant $k_j$ is computed from a modified Arrhenius equation: $$ k_j = AT^n \exp^{-E_a/(RT)} $$ This procedure is managed by the OpenSMOKE++ library ([opensmoke.h](opensmoke.h)). */ #include "radiation.h" /** ## User Data Structure that gathers arguments that can be used inside the batch reactor system functions. */ typedef struct { double rho; double cp; double P; double T; } UserDataODE; /** ## *batch_isothermal_constantpressure()*: solve chemical species reactions * *y*: vector (length = NS) containing the initial values of the system * *dt*: simulation time step * *dy*: right-hand-side of the batch system of equations * *args*: structure with additional arguments to be passed to the system */ void batch_isothermal_constantpressure (const double * y, const double dt, double * dy, void * args) { /** Unpack data for the ODE system. */ UserDataODE data = *(UserDataODE *)args; double rho = data.rho; double cp = data.cp; double Pressure = data.P; double Temperature = data.T; /** Set temperature and pressure of the system, for the calculation of the kinetic constants. */ OpenSMOKE_GasProp_SetTemperature (Temperature); OpenSMOKE_GasProp_SetPressure (Pressure); /** We create a vector with the mass fractions and we convert it to mole fractions. */ double massfracs[NGS], molefracs[NGS]; for (int jj=0; jj