The ISPH package solves large-scale multi-physics mesoscale flow problems using implicit smoothed particle hydrodynamics methods. It relies on an incremental pressure correction scheme (effective for flows with moderately small Reynolds numbers) and differential operator renormalizations to achieve second-order accuracy in both time and space. ISPH, built on LAMMPS and Trilinos, provides a massively parallel scalable simulation capability for incompressible, multiphase flow coupled with electrostatics. In particular, the code has been used to study pore-scale flows in a complex geometry discretized with 165M particles, electroosmotic flow and diffusive mixing of solute in micro channels. Results from these researches are published in 3 journal papers thus far. The code is available at https://github.com/kyungjoo-kim/implicit-sph.