code: def S(a,b):
      
      if b.contracted:
          return sum(cb*S(pb,a) for (cb,pb) in b)
      elif a.contracted:
          return sum(ca*S(b,pa) for (ca,pa) in a)
      return a.norm*b.norm*overlap(a.exponent,a.powers,
                                   a.origin,b.exponent,b.powers,b.origin)
docs: Simple interface to the overlap function.
      >>> from pyquante2 import pgbf,cgbf
      >>> s = pgbf(1)
      >>> isclose(S(s,s),1.0)
      True
      >>> sc = cgbf(exps=[1],coefs=[1])
      >>> isclose(S(sc,sc),1.0)
      True
  
      >>> sc = cgbf(exps=[1],coefs=[1])
      >>> isclose(S(sc,s),1.0)
      True
      >>> isclose(S(s,sc),1.0)
      True