We study the seismic inverse problem associated with the time-harmonic wave equations for the reconstruction of subsurface media. The reconstruction is conducted using the full waveform inversion (FWI) method and relies on iterative minimization algorithm, which we adapt for large scale situation. In particular the inverse problem shows a conditional Lipschitz stability when assuming piecewise constant representation of the parameters. We give the analytical lower and upper bounds for the stability constant and provide quantitative numerical estimates to demonstrate the sharpness of these bounds in the geophysical context. We further study the convergence of the minimization problem and are able to numerically estimate the size of the basin of attraction, depending on the frequency. From these stability and convergence results, we design a multi-level algorithm with simultaneous progression in frequency and scale. Eventually we carry out numerical experiments for acoustic and elastic parameters reconstruction assuming no-prior information in the initial models, in two and three dimensions.