A numerical method is presented to study cavity and bubble dynamics. The liquid phase is assumed to be inviscid and incom-pressible and separated from the gas or vacuum phase by a free surface. On the free surface the stress tensor reduces to a spatially constant pressure. The flow in the bulk of the liquid is computed using a second-order-in-time projection method. The interface is advected and reconstructed using a Volume-of-Fluid (VOF) method. Setting the pressure on the free surface to the prescribed value involves a modified stencil on nodes close to the interface. This modifed stencil contains interpolated pressures on branches that are cut by the interface. Capillary effects are taken into account by adding the Laplace law pressure increment to these prescribed pressures. The curvature that appears in the Laplace law is computed using the height-function method. The VOF advection and momentum advection schemes both require an extension of the velocity in a two-layer wide ghost cell region on the grid across the free surface. This ghost layer is computed in two stages. In the preliminary stage a first-order velocity extrapolation of the liquid velocity field to the ghost layers is performed. In the second stage the ghost layer velocities are projected on the space of divergence free velocities using an auxiliary projection step. The whole procedure is implemented in a free code developed with the help of Gretar Tryggvason and Yue (Stanley) Ling and is available at http://parissimulator.sf.net. Tests are perfomed on radial flows with spherical symmetry except for boundary conditions far from the bubble in a cubic box. In such a geometry the flow is predicted by solutions of the Rayleigh-Plesset equation. Good comparison to the Rayleigh-Plesset solution for a single bubble with low and moderate amplitude oscillations is shown. Perspectives for parallel simulations involving very large numbers of bubbles are given.