code: def plot(self, xmin=-1, xmax=1, center=0,
               resolution_outside=20, resolution_inside=200):
          
          if self.eps == 0:
              return [xmin, center, xmax], [0, 0, xmax]
          else:
              n = float(resolution_inside)/self.eps
              x = np.concatenate((
                  np.linspace(xmin, center-self.eps, resolution_outside+1),
                  np.linspace(center-self.eps, center+self.eps, n+1),
                  np.linspace(center+self.eps, xmax, resolution_outside+1)))
              y = self(x)
              return x, y
docs: Return arrays x, y for plotting the Heaviside function
          H(x-`center`) on [`xmin`, `xmax`]. For the exact
          Heaviside function,
          ``x = [xmin, center, xmax]; y = [0, 0, 1]``,
          while for the smoothed version, the ``x`` array
          is computed on basis of the `eps` parameter, with
          `resolution_outside` intervals on each side of the smoothed
          region and `resolution_inside` intervals in the smoothed region.