Migration of magma within an active volcano produces a deformation signature at the Earth's surface. The internal structure of a volcano and specific movements of the magma control the actual deformation that is observed. Relatively simple models that simulate magma injection as a pressurized body embedded in a homogeneous elastic half-space (e.g., Mogi) can predict the characteristic radially-symmetric deformation patterns that are commonly observed for episodes of volcano inflation or deflation. Inverse methods, based on half-space models, can precisely and efficiently estimate the non-linear parameters that describe the geometry (position and shape) of the deformation source, as well as the linear parameter that describes the strength (pressure) of the deformation source. However, although such models can accurately predict the observed deformation, actual volcanoes have internal structures that are not compatible with the elastic half-space assumptions inherent to Mogi-type models. This incompatibility translates to errors in source parameter estimations. Alternatively, Finite Element Models (FEMs) can simulate a pressurized body embedded in a problem domain having an arbitrary distribution of material properties that better corresponds to the internal structure of an active volcano. FEMs can be used in inverse methods for estimating linear deformation source parameters, such as the source pressure. However, perturbations of the non-linear parameters that describe the geometry of the source require automated re-meshing of the problem domain - a significant obstacle to implementing FEM-based nonlinear inverse methods in volcano deformation studies. I present a parametric executable (C++ source code), which automatically generates FEMs that simulate a pressurized ellipsoid embedded in an axisymmetric problem domain, having an a priori distribution of material properties. I demonstrate this executable by analyzing Interferometric Synthetic Aperture Radar (InSAR) deformation data of the 1997 eruption of Okmok Volcano, Alaska as an example. This executable facilitates an inverse analysis that estimates the non-linear parameters that describe the depth and radius of the spherical source, as well as the linear strength parameter that best accounts for the InSAR data. The strong radial symmetry and high signal-to-noise ratio of the InSAR data, along with known seismic tomography data, provide robust constraints for estimated parameters and sensitivity analyses.