def testSerializationLinearBoudaryBoundingBox(self):
        """Uses Linear grid for tests"""
        from pysgpp import Grid
        
        factory = Grid.createLinearBoundaryGrid(2)
        self.failIfEqual(factory, None)

        gen = factory.createGridGenerator()
        gen.regular(3)

        boundingBox = factory.getBoundingBox()
        tempBound = boundingBox.getBoundary(0)
        tempBound.leftBoundary = 0.0
        tempBound.rightBoundary = 100.0
        tempBound.bDirichletLeft = False;
        tempBound.bDirichletRight = False;
        boundingBox.setBoundary(0, tempBound)

        str = factory.serialize()
        self.assert_(len(str) > 0)
        
        newfac = Grid.unserialize(str)
        self.failIfEqual(newfac, None)
        
        self.assertEqual(factory.getStorage().size(), newfac.getStorage().size())
            
        boundingBox = newfac.getBoundingBox()
        tempBound = boundingBox.getBoundary(0)
        self.assertEqual(0.0, tempBound.leftBoundary)
        self.assertEqual(100.0, tempBound.rightBoundary)
        self.assertEqual(False, tempBound.bDirichletLeft)
        self.assertEqual(False, tempBound.bDirichletRight)

def createGrid(grid, dim, deg=1, addTruncatedBorder=False):
    # create new grid
    gridType = grid.getType()

    if gridType in [Poly, PolyBoundary]:
        deg = max(deg, grid.getDegree())

    # print gridType, deg
    if deg > 1:
        if gridType in [LinearBoundary, PolyBoundary]:
            return Grid.createPolyBoundaryGrid(dim, deg)
        elif gridType == LinearL0Boundary:
            raise NotImplementedError("there is no full boundary polynomial grid")
        elif gridType in [Linear, Poly]:
            return Grid.createPolyGrid(dim, deg)
        else:
            raise Exception('unknown grid type %s' % gridType)
    else:
        if gridType == Linear:
            return Grid.createLinearGrid(dim)
        elif gridType == LinearBoundary:
            return Grid.createLinearBoundaryGrid(dim)
        elif gridType == LinearL0Boundary:
            return Grid.createLinearBoundaryGrid(dim, 0)
        else:
            raise Exception('unknown grid type %s' % gridType)

    def testRefinement2d_two(self):
        from pysgpp import Grid, DataVector, SurplusRefinementFunctor
        factory = Grid.createLinearBoundaryGrid(2)
        storage = factory.getStorage()
        
        gen = factory.createGridGenerator()
        gen.regular(0)
        
        alpha = DataVector(4)
    
        for i in xrange(len(alpha)):
            alpha[i] = 0.0

        alpha[0] = 1.0
        func = SurplusRefinementFunctor(alpha)
            
        gen.refine(func)
        
        alpha2 = DataVector(8)
    
        for i in xrange(len(alpha2)):
            alpha2[i] = 0.0

        alpha2[4] = 1.0
        func = SurplusRefinementFunctor(alpha2)
            
        gen.refine(func)    
        self.failUnlessEqual(storage.size(), 13)        

 def testOperationB(self):
     from pysgpp import Grid, DataVector, DataMatrix
     factory = Grid.createLinearBoundaryGrid(1)
     gen = factory.createGridGenerator()
     gen.regular(2)
     
     alpha = DataVector(factory.getStorage().size())
     p = DataMatrix(1,1)
     beta = DataVector(1)
     
     
     alpha.setAll(0.0)
     p.set(0,0,0.25)
     beta[0] = 1.0
     
     opb = factory.createOperationB()
     opb.mult(beta, p, alpha)
     
     self.failUnlessAlmostEqual(alpha[0], 0.75)
     self.failUnlessAlmostEqual(alpha[1], 0.25)
     self.failUnlessAlmostEqual(alpha[2], 0.5)
     self.failUnlessAlmostEqual(alpha[3], 1.0)
     self.failUnlessAlmostEqual(alpha[4], 0.0)
     
     alpha.setAll(0.0)
     alpha[2] = 1.0
     
     p.set(0,0, 0.25)
     
     beta[0] = 0.0
     
     opb.multTranspose(alpha, p, beta)
     self.failUnlessAlmostEqual(beta[0], 0.5)

def checkPositivity(grid, alpha):
    # define a full grid of maxlevel of the grid
    gs = grid.getStorage()
    fullGrid = Grid.createLinearGrid(gs.dim())
    fullGrid.createGridGenerator().full(gs.getMaxLevel())
    fullHashGridStorage = fullGrid.getStorage()
    A = DataMatrix(fullHashGridStorage.size(), fullHashGridStorage.dim())
    p = DataVector(gs.dim())
    for i in xrange(fullHashGridStorage.size()):
        fullHashGridStorage.get(i).getCoords(p)
        A.setRow(i, p)

    res = evalSGFunctionMulti(grid, alpha, A)
    ymin, ymax, cnt = 0, -1e10, 0
    for i, yi in enumerate(res.array()):
        if yi < 0. and abs(yi) > 1e-13:
            cnt += 1
            ymin = min(ymin, yi)
            ymax = max(ymax, yi)
            A.getRow(i, p)
            print "  %s = %g" % (p, yi)
    if cnt > 0:
        print "warning: function is not positive"
        print "%i/%i: [%g, %g]" % (cnt, fullHashGridStorage.size(), ymin, ymax)
    return cnt == 0

    def testOperationTest_test(self):
        from pysgpp import Grid, DataVector, DataMatrix

        factory = Grid.createLinearBoundaryGrid(1)
        gen = factory.createGridGenerator()
        gen.regular(1)
        
        alpha = DataVector(factory.getStorage().size())        
        
        data = DataMatrix(1,1)
        data.setAll(0.25)
        classes = DataVector(1)
        classes.setAll(1.0)

        testOP = factory.createOperationTest()

        alpha[0] = 0.0
        alpha[1] = 0.0
        alpha[2] = 1.0
        
        c = testOP.test(alpha, data, classes)
        self.failUnless(c > 0.0)
        
        alpha[0] = 0.0
        alpha[1] = 0.0
        alpha[2] = -1.0
        c = testOP.test(alpha, data, classes)
        self.failUnless(c == 0.0)

def eval_fullGrid(level, dim, border=True):
    if border:
        grid = Grid.createLinearBoundaryGrid(dim)
    else:
        grid = Grid.createLinearGrid(dim)

    grid.createGridGenerator().full(level)
    gs = grid.getStorage()
    ans = DataMatrix(gs.size(), dim)
    p = DataVector(dim)

    for i in xrange(gs.size()):
        gs.get(i).getCoords(p)
        ans.setRow(i, p)

    return ans

 def test34(self):
     from pysgpp import Grid, DataVector, FullGrid, FullGridSet
     
     dim = 2
     level = 9
     function = buildParableBoundary(dim)
     grid = Grid.createLinearBoundaryGrid(dim)
     testFG(self, grid, level, function)

 def testHierarchisationDModLinear(self):
     from pysgpp import Grid
     
     dim = 3
     level = 5
     function = buildParable(dim)
     grid = Grid.createModLinearGrid(dim)
     testHierarchisationDehierarchisation(self, grid, level, function)

    def testSerializationLinearBoudary(self):
        """Uses Linear grid for tests"""
        from pysgpp import Grid
        
        factory = Grid.createLinearBoundaryGrid(2)
        self.failIfEqual(factory, None)

        gen = factory.createGridGenerator()
        gen.regular(3)

        str = factory.serialize()
        self.assert_(len(str) > 0)
        
        newfac = Grid.unserialize(str)
        self.failIfEqual(newfac, None)
        
        self.assertEqual(factory.getStorage().size(), newfac.getStorage().size())        

 def testHierarchisationDBoundary(self):
     from pysgpp import Grid
     
     dim = 3
     level = 5
     function = buildParableBoundary(dim)
     grid = Grid.createLinearBoundaryGrid(dim)
     testHierarchisationDehierarchisation(self, grid, level, function)

 def test44(self):
     from pysgpp import Grid, DataVector, FullGrid, FullGridSet
     
     dim = 2
     level = 8
     function = buildParableBoundary(dim)
     grid = Grid.createSquareRootGrid(dim)
     testFG(self, grid, level, function)  

 def testCreation(self):
     """Uses Linear grid for tests"""
     from pysgpp import Grid
     
     factory = Grid.createLinearGrid(2)
     self.failIfEqual(factory, None)
     
     storage = factory.getStorage()
     self.failIfEqual(storage, None)

 def setUp(self):
     self.grid = Grid.createLinearGrid(2)  # a simple 2D grid
     self.grid.createGridGenerator().regular(3)  # max level 3 => 17 points
     self.HashGridStorage = self.grid.getStorage()
     alpha = DataVector(self.grid.getSize())
     alpha.setAll(1.0)
     for i in [9, 10, 11, 12]:
         alpha[i] = 0.0
     coarseningFunctor = SurplusCoarseningFunctor(alpha, 4, 0.5)
     self.grid.createGridGenerator().coarsen(coarseningFunctor, alpha)

    def testHatRegulardD_two(self):
        from pysgpp import Grid
        
        factory = Grid.createLinearBoundaryGrid(3)

        m = generateLaplaceMatrix(factory, 4)
        m_ref = readReferenceMatrix(self, factory.getStorage(), 'data/C_laplace_phi_li_hut_l0_rand_dim_3_nopsgrid_297_float.dat.gz')

        # compare
        compareStiffnessMatrices(self, m, m_ref)    

    def testHatRegular1D_two(self):
        from pysgpp import Grid
        
        factory = Grid.createLinearTrapezoidBoundaryGrid(1)

        m = generateLaplaceMatrix(factory, 5)
        m_ref = readReferenceMatrix(self, factory.getStorage(), 'data/C_laplace_phi_li_hut_trapezrand_dim_1_nopsgrid_33_float.dat.gz')

        # compare
        compareStiffnessMatrices(self, m, m_ref) 

    def testHatRegular1D(self):
        from pysgpp import Grid
        
        factory = Grid.createModLinearGrid(1)

        m = generateLaplaceMatrix(factory, 5)
        m_ref = readReferenceMatrix(self, factory.getStorage(), 'data/C_laplace_phi_li_ausgeklappt_dim_1_nopsgrid_31_float.dat.gz')

        # compare
        compareStiffnessMatrices(self, m, m_ref)

def createGrid(dim, level, borderType, isFull=False):
    from pysgpp.extensions.datadriven.learner.Types import BorderTypes
    if borderType == BorderTypes.NONE:
        grid = Grid.createLinearGrid(dim)
    elif borderType == BorderTypes.TRAPEZOIDBOUNDARY:
        grid = Grid.createLinearTrapezoidBoundaryGrid(dim)
    elif borderType == BorderTypes.COMPLETEBOUNDARY:
        grid = Grid.createLinearBoundaryGrid(dim, 0)
    else:
        raise Exception('Unknown border type')

    # create regular grid of level accLevel
    gridGen = grid.createGridGenerator()
    if isFull:
        gridGen.full(level)
    else:
        gridGen.regular(level)

    return grid

    def general_test(self, d, l, bb, xs):

        test_desc = "dim=%d, level=%d, len(x)=%s" % (d, l, len(xs))

        print test_desc

        self.grid = Grid.createLinearGrid(d)
        self.grid_gen = self.grid.createGridGenerator()
        self.grid_gen.regular(l)

        alpha = DataVector([self.get_random_alpha() for i in xrange(self.grid.getSize())])

        bb_ = BoundingBox(d)

        for d_k in xrange(d):
            dimbb = DimensionBoundary()
            dimbb.leftBoundary = bb[d_k][0]
            dimbb.rightBoundary = bb[d_k][1]
            bb_.setBoundary(d_k, dimbb)

        # Calculate the expected value without the bounding box

        expected_normal = [self.calc_exp_value_normal(x, d, bb, alpha) for x in xs]
        #expected_transposed = [self.calc_exp_value_transposed(x, d, bb, alpha) for x in xs]

        # Now set the bounding box

        self.grid.getStorage().setBoundingBox(bb_)

        dm = DataMatrix(len(xs), d)
        for k, x in enumerate(xs):
            dv = DataVector(x)
            dm.setRow(k, dv)

        multEval = createOperationMultipleEval(self.grid, dm)

        actual_normal = DataVector(len(xs))
        #actual_transposed = DataVector(len(xs))

        multEval.mult(alpha, actual_normal)
        #multEval.mult(alpha, actual_transposed)

        actual_normal_list = []
        for k in xrange(len(xs)):
            actual_normal_list.append(actual_normal.__getitem__(k))

        #actual_transposed_list = []
        #for k in xrange(len(xs)):
        #    actual_transposed_list.append(actual_transposed.__getitem__(k))

        self.assertAlmostEqual(actual_normal_list, expected_normal)
        #self.assertAlmostEqual(actual_tranposed_list, expected_tranposed)

        del self.grid

def makePositive(grid, alpha):
    """
    insert full grid points if they are negative and the father
    node is part of the sparse grid
    @param grid:
    @param alpha:
    """
    # copy old sg function
    jgrid = copyGrid(grid)

    # evaluate the sparse grid function at all full grid points
    level = grid.getStorage().getMaxLevel()
    fg = Grid.createLinearGrid(grid.getStorage().dim())
    fg.createGridGenerator().full(level)

    # copy the old grid and use it as reference
    jgs = jgrid.getStorage()
    fgs = fg.getStorage()

    # run over all results and check where the function value
    # is lower than zero
    cnt = 1
    while True:
        print "run %i: full grid size = %i" % (cnt, fgs.size())
        gps = []

        # insert those fg points, which are not yet positive
        values = computeNodalValues(fg, grid, alpha)
        for i in xrange(len(values)):
            gp = fgs.get(i)
            if values[i] < 0 and not jgs.has_key(gp):
                gps += insertPoint(jgrid, gp)
                gps += insertHierarchicalAncestors(jgrid, gp)

        jgrid.getStorage().recalcLeafProperty()

        # 1. compute nodal values for new grid points
        jnodalValues = computeNodalValues(jgrid, grid, alpha)
        # 2. set the new ones to zero
        jgs = jgrid.getStorage()
        for gp in gps:
            jnodalValues[jgs.seq(gp)] = 0.
        # 3. hierarchize
        jalpha = hierarchize(jgrid, jnodalValues)
        # stop loop if no points have been added
        if len(gps) == 0:
            break
        # 4. reset values for next loop
        grid = copyGrid(jgrid)
        alpha = DataVector(jalpha)
        cnt += 1

    return jgrid, jalpha