We present a new hybrid approach to computing penetration depth (PD) for general polygonal models. Our approach exploits both local and global approaches to PD computation and can compute error-bounded PD approximations for both deep and shallow penetrations. We use a two-step formulation: the first step corresponds to a global approximation approach that samples the configuration space with bounded error using support vector machines; the second step corresponds to a local optimization that performs a projection operation refining the penetration depth. We have implemented this hybrid algorithm on a standard PC platform and tested its performance with various benchmarks. The experimental results show that our algorithm offers significant benefits over previously developed local-only and global-only methods used to compute the PD.