Integrate passive scalars into the EOS

[msbc]

I want to extend the General EOS capabilities to include passive scalars (see #359 ). For example, one could have a composition dependant EOS and use the passive scalars to track composition. I'd like to gather input from others to help me implement this feature myself.

The primary issue here is that the EOS is called within the Riemann Solver. Therefore the face reconstructed passive scalars need to be passed to the Riemann Solver, and then passed to the EOS calls within. Currently the reconstruction and the flux calculations for the hydro and scalar variables are completely separated with different function calls and different tasks. Ideally, these would be unified in some way.

The two main ways I see of doing this is either:

1. Extending the hydro athena array to include the scalars, so that the zeroth scalar index is NHYDRO (when MHD is disabled).
2. Keep the hydro and scalar arrays separated but combine the two separate flux/reconstruction functions/tasks.

A possible alternative would be to add scalar reconstruction to the hydro steps but to only use this computation for the Riemann problem and to disregard it and recompute it in the scalar task. This adds unnecessary computation but is probably least likely to break something.