I have a set of two simultaneous ODEs of a heat exchanger system with the boundary conditions x(0)=753K and y(4)=573K. Now I need to incorporate a reaction term into one of the equations whose rate changes with the temp and the extent. How to solve?

I cannot use ODE45 since its not an initial value problem and discretisation gives a profile does not comply to the expected result since an overshoot is observed. Is there any other way to solve this? The reaction term is in the power law model.