def get_harmonic_bonds(system):
    """Get a list of all harmonic bonds in the system"""
    bonds = []
    for f in system.getForces():
        if "HarmonicBond" in str(f.__class__):
            for i in range(f.getNumBonds()):
                bonds.append(f.getBondParameters(i))
    return bonds