  Studies of neutron stars are extremely timely given the recent detection of
gravitational waves from a binary neutron star merger GW170817, and an
International Space Station payload NICER currently in operation that aims to
determine radii of neutron stars to a precision better than 5%. In many cases,
neutron star solutions are constructed numerically due to the complexity of the
field equations with realistic equations of state. However, in order to relate
observables like the neutron star mass and radius to interior quantities like
central density and pressure, it would be useful to provide an accurate,
analytic modeling of a neutron star interior. One such solution for static and
isolated neutron stars is the Tolman VII solution characterized only by two
parameters (e.g. mass and radius), though its agreement with numerical
solutions is not perfect. We here introduce an improved analytic model based on
the Tolman VII solution by introducing an additional parameter to make the
analytic density profile agree better with the numerically obtained one. This
additional parameter can be fitted in terms of the stellar mass, radius and
central density in an equation-of-state-insensitive way. In most cases, we find
that the new model more accurately describes realistic profiles than the
original Tolman VII solution by a factor of 2-5. Our results are first-step
calculations towards constructing analytic interior solutions for more
realistic neutron stars under rotation or tidal deformation.
