Instead of treating the polytrope as a free boundary problem I have defined an interpolating polynominal, accurate to within 0.01 percent over n=[0,5) which is used to set the size of the domain for a given n
the polytrope module will be used as an initial guess to the solver. A skeleton of this has been imported from https://github.com/tboudreaux/FEMPolytrope
This module will need major updates still to handle 3D, proper boundary conditions, and to incorporate it with the rest of our meshing scheme
