docs: def sigmoid(f=Ellipsis, mu=0, sigma=1, scale=1, invert=False, normalize=False):
      '''
      sigmoid() yields a potential function that is equivalent to the integral of gaussian(), i.e.,
        the error function, but scaled to match gaussian().
      sigmoid(f) is equivalent to compose(sigmoid(), f).
  
      All options that are accepted by the gaussian() function are accepted by sigmoid() with the same
      default values and are handled in an equivalent manner with the exception of the invert option;
      when a sigmoid is inverted, the function approaches its maximum value at -inf and approaches 0
      at inf.
  
      Note that because sigmoid() explicitly matches gaussian(), the base formula used is as follows:
        f(x) = scale * sigma * sqrt(pi/2) * erf((x - mu) / (sqrt(2) * sigma))
        k*sig*Sqrt[Pi/2] Erf[(x - mu)/sig/Sqrt[2]]
      '''
      f = to_potential(f)
      F = erf((f - mu) / (sigma * np.sqrt(2.0)))
      if invert: F = 1 - F
      F = np.sqrt(np.pi / 2) * scale * F
      if normalize: F = F / (np.sqrt(2.0*np.pi) * sigma)
      return F