How to use the splipy.Volume function in Splipy

"""
    if len(surfaces) == 1: # probably gives input as a list-like single variable
        surfaces = surfaces[0]
    if len(surfaces) == 2:
        surf1 = surfaces[0].clone()
        surf2 = surfaces[1].clone()
        Surface.make_splines_identical(surf1, surf2)
        (n1, n2, d) = surf1.controlpoints.shape  # d = dimension + rational

        controlpoints = np.zeros((n1, n2, 2, d))
        controlpoints[:, :, 0, :] = surf1.controlpoints
        controlpoints[:, :, 1, :] = surf2.controlpoints

        # Volume constructor orders control points in a different way, so we
        # create it from scratch here
        result = Volume(surf1.bases[0], surf1.bases[1], BSplineBasis(2), controlpoints,
                         rational=surf1.rational, raw=True)

        return result
    elif len(surfaces) == 6:
        if any([surf.rational for surf in surfaces]):
            raise RuntimeError('edge_surfaces not supported for rational splines')

        # coons patch (https://en.wikipedia.org/wiki/Coons_patch)
        umin = surfaces[0]
        umax = surfaces[1]
        vmin = surfaces[2]
        vmax = surfaces[3]
        wmin = surfaces[4]
        wmax = surfaces[5]
        vol1 = edge_surfaces(umin,umax)
        vol2 = edge_surfaces(vmin,vmax)

def cube(size=1, lower_left=(0,0,0)):
    """  Create a cube with parmetric origin at *(0,0,0)*.

    :param float size: Size(s), either a single scalar or a tuple of scalars per axis
    :param array-like lower_left: local origin, the lower bottom left corner of the cube
    :return: A linear parametrized box
    :rtype: Volume
    """
    result = Volume()
    result.scale(size)
    result += lower_left
    return result

return result
    elif len(surfaces) == 6:
        if any([surf.rational for surf in surfaces]):
            raise RuntimeError('edge_surfaces not supported for rational splines')

        # coons patch (https://en.wikipedia.org/wiki/Coons_patch)
        umin = surfaces[0]
        umax = surfaces[1]
        vmin = surfaces[2]
        vmax = surfaces[3]
        wmin = surfaces[4]
        wmax = surfaces[5]
        vol1 = edge_surfaces(umin,umax)
        vol2 = edge_surfaces(vmin,vmax)
        vol3 = edge_surfaces(wmin,wmax)
        vol4 = Volume(controlpoints=vol1.corners(order='F'), rational=vol1.rational)
        vol1.swap(0, 2)
        vol1.swap(1, 2)
        vol2.swap(1, 2)
        vol4.swap(1, 2)
        Volume.make_splines_identical(vol1, vol2)
        Volume.make_splines_identical(vol1, vol3)
        Volume.make_splines_identical(vol1, vol4)
        Volume.make_splines_identical(vol2, vol3)
        Volume.make_splines_identical(vol2, vol4)
        Volume.make_splines_identical(vol3, vol4)
        result                =   vol1.clone()
        result.controlpoints +=   vol2.controlpoints
        result.controlpoints +=   vol3.controlpoints
        result.controlpoints -= 2*vol4.controlpoints
        return result
    else:

:return: Interpolated volume
    :rtype: Volume
    """
    vol_shape = [b.num_functions() for b in bases]
    dim = x.shape[-1]
    if len(x.shape) == 2:
        x = x.reshape(vol_shape + [dim])
    if u is None:
        u = [b.greville() for b in bases]
    N_all = [b(t) for b,t in zip(bases, u)]
    N_all.reverse()
    cp = x
    for N in N_all:
        cp = np.tensordot(np.linalg.inv(N), cp, axes=(1,2))

    return Volume(bases[0], bases[1], bases[2], cp.transpose(2,1,0,3).reshape((np.prod(vol_shape),dim)))

def get_cm1_mesh(self):
        # Create the C^{-1} mesh
        nx, ny, nz = self.n
        Xm1 = self.raw.get_discontinuous_all()
        b1 = BSplineBasis(2, sorted(list(range(self.n[0]+1))*2))
        b2 = BSplineBasis(2, sorted(list(range(self.n[1]+1))*2))
        b3 = BSplineBasis(2, sorted(list(range(self.n[2]+1))*2))
        discont_vol = Volume(b1, b2, b3, Xm1)
        return discont_vol

# does not work with rational surfaces, so we'll just manually try
        # and add some inner controlpoints
        cp   = np.zeros((5,5,5,4))
        cp[ :, :, 0,:] = wmin[:,:,:]
        cp[ :, :,-1,:] = wmax[:,:,:]
        cp[ :, 0, :,:] = vmin[:,:,:]
        cp[ :,-1, :,:] = vmax[:,:,:]
        cp[ 0, :, :,:] = umin[:,:,:]
        cp[-1, :, :,:] = umax[:,:,:]
        inner = np.linspace(-.5,.5, 3)
        Y, X, Z = np.meshgrid(inner,inner,inner)
        cp[1:4,1:4,1:4,0] = X
        cp[1:4,1:4,1:4,1] = Y
        cp[1:4,1:4,1:4,2] = Z
        cp[1:4,1:4,1:4,3] = 1
        ball = Volume(b,b,b,cp,rational=True, raw=True)
        return r*ball + center
    else:
        raise ValueError('invalid type argument')

def write(self, obj, n=None):
        if isinstance(obj, SplineModel):
            if obj.pardim == 3: # volume model
                for surface in obj.boundary():
                    self.write_surface(surface.obj,n)
            elif obj.pardim == 2: # surface model
                for surface in obj:
                    self.write_surface(surface, n)

        elif isinstance(obj, Volume):
            for surface in obj.faces():
                self.write_surface(surface, n)

        elif isinstance(obj, Surface):
            self.write_surface(obj, n)

        else:
            raise ValueError('Unsopported object for STL format')

def controlpoints(spline):
    """ Return controlpoints according to nutils ordering """
    n = len(spline)
    dim = spline.dimension
    if isinstance(spline, Curve):
        return np.reshape(spline[:,:]                  , (n, dim), order='F')
    elif isinstance(spline, Surface):
        return np.reshape(spline[:,:,:].swapaxes(0,1)  , (n, dim), order='F')
    elif isinstance(spline, Volume):
        return np.reshape(spline[:,:,:,:].swapaxes(0,2), (n, dim), order='F')
    raise RuntimeError('Non-spline argument detected')

def get_mixed_cont_mesh(self):
        # Create mixed discontinuity mesh: C^0, C^0, C^{-1}
        nx, ny, nz = self.n
        Xz = self.raw.get_discontinuous_z()
        b1 = BSplineBasis(2, sorted(list(range(self.n[0]+1))+[0,self.n[0]]))
        b2 = BSplineBasis(2, sorted(list(range(self.n[1]+1))+[0,self.n[1]]))
        b3 = BSplineBasis(2, sorted(list(range(self.n[2]+1))*2))
        true_vol = Volume(b1, b2, b3, Xz, raw=True)
        return true_vol