Issue with symbolic solver - disparity between methods/results

I have a system of six ODEs that I wish to solve symbolically and then find the numeric solutions over a range of operating conditions (two parameters, a, b).

The slow method:

- Iteratively (two for loops for range(a) and range(b)) - Prior to iteration, create matlabFunction for each equation (1-6) - Iteratively replace a and b into the matlab functions to give eq1* - eq6* (with only the 6 variables to be solved as unknowns) - Use assume(variable>=0) so that solver knows all variable solutions are non-negative. - Iteratively solve the set of equations using [solutions]=solve(eqs*,variables,'Real', true);

-- This method gives the correct solutions for each pair of (a,b).

The fast method:

- Solve the system of equations without replacing (a,b) - [sols]=solve(eqs,variables) gives symbolic solutions with (a,b) in each solution - Create matlabFunctions for each solution and iteratively replace (a,b) to get the final numerical solutions

-- This method gives incorrect results. In some cases when the solutions should be zero, there are negative results, for example.

I have similar systems of three, four and five ODEs, and the Fast Method works well in all cases. I know things become challenging when the number of equations > 4 or 5.

Is this merely a case of Matlab struggling to get precise solutions for such a complex set of equations when there are two unknowns in the solution and then subsequent loss of precision when those values are entered into the solutions?

It takes < 2 minutes to get the solutions using the Fast Method for a relatively large parameter space (50 x 50), whereas > 2 hours using the slow method.

I am looking to anyone who could confirm why this inaccuracy may occur and also any suggestions for speeding up this process given that I cannot see any other way other than iteratively.

Regards Matt