docs: def _rk4_integrate(self, x, y, uv0fun, uv1fun, dt):
          """Integrates positions x, y using velocity functions
             uv0fun, uv1fun. Returns dx and dy, the displacements."""
          u0, v0 = uv0fun(x, y)
          k1u = dt*u0
          k1v = dt*v0
          x11 = self._wrap_x(x + 0.5*k1u)
          y11 = self._wrap_y(y + 0.5*k1v)
          u11, v11 = uv1fun(x11, y11)
          k2u = dt*u11
          k2v = dt*v11
          x12 = self._wrap_x(x + 0.5*k2u)
          y12 = self._wrap_y(y + 0.5*k2v)
          u12, v12 = uv1fun(x12, y12)
          k3u = dt*u12
          k3v = dt*v12
          x13 = self._wrap_x(x + k3u)
          y13 = self._wrap_y(y + k3v)
          u13, v13 = uv1fun(x13, y13)
          k4u = dt*u13
          k4v = dt*v13
          
          # update
          dx = 6**-1*(k1u + 2*k2u + 2*k3u + k4u)
          dy = 6**-1*(k1v + 2*k2v + 2*k3v + k4v)
          return dx, dy