def getNthPolyPrime( n, poly ):
    result = getNthPrime( n )

    for i in arange( 1, poly ):
        result = getNthPrime( result )

    return result

def OLDgetPartitionNumber( n ):
    if n < 0:
        return 0

    if n < 2:
        return 1

    result = mpmathify( 0 )

    for k in arange( 1, n + 1 ):
        #n1 = n - k * ( 3 * k - 1 ) / 2
        n1 = fsub( n, fdiv( fmul( k, fsub( fmul( 3, k ), 1 ) ), 2 ) )
        #n2 = n - k * ( 3 * k + 1 ) / 2
        n2 = fsub( n, fdiv( fmul( k, fadd( fmul( 3, k ), 1 ) ), 2 ) )

        result = fadd( result, fmul( power( -1, fadd( k, 1 ) ), fadd( getPartitionNumber( n1 ), getPartitionNumber( n2 ) ) ) )

        if n1 <= 0:
            break

    #old = NOT_QUITE_AS_OLDgetPartitionNumber( n )
    #
    #if ( old != result ):
    #    raise ValueError( "It's broke." )

    return result

def CalibrationHa(z,oldLum,it):
	newLum = [0. for x in range((it.shape[0]))]
	oldpara = np.array([[0. for i in range(3)] for i in range(it.shape[0])])
	oldLF = np.array([0. for i in range(it.shape[0])])
	interpolatingSize =100
	for i in range(oldpara.shape[0]):
	        oldpara[i][0] = 1.37e-3
	        oldpara[i][2] = -1.35
	        if z[it[i]] <1.3 :oldpara[i][1] = 5.1e41*np.power((1+z[it[i]]),3.1)
	        if z[it[i]] >=1.3:oldpara[i][1] =6.8e42
	        oldLF[i] = LuminosityFunction(oldpara[i],oldLum[i])
	rangea = oldLF.copy()
	rangea.sort()
	rangea = rangea[rangea>0]
	grid = np.array([[0. for i in range(interpolatingSize)] for i in range(2)])
	mpgrid = mpmath.arange(mpmath.log10(rangea.min()),mpmath.log10(rangea.max()),(np.log10(rangea.max())-np.log10(rangea.min()))/interpolatingSize)
	grid[0] = np.ogrid[np.log10(rangea.min()):np.log10(oldLF.max()):interpolatingSize*1j]
	grid[1] = np.ogrid[0:6:interpolatingSize*1j]
	newLF = np.array([[0.for i in range(interpolatingSize)] for i in range(interpolatingSize)])
	for i in range(interpolatingSize):
		for j in range(interpolatingSize):
			newLF[i][j]= optimize.brentq(rootfinding,1e30,1e50,args=(luminosityParameterHa(grid[1][i]),mpgrid[j]))
	newLum = interpolate.interpn([grid[0],grid[1]],newLF,np.array([np.log10(oldLF),z[it]]).T,bounds_error=False,fill_value=0.0)
	newLum = np.power(10,newLum)
	return newLum

def getNthDelannoyNumber( n ):
    result = 0

    for k in arange( 0, fadd( real( n ), 1 ) ):
        result = fadd( result, fmul( binomial( n, k ), binomial( fadd( n, k ), k ) ) )

    return result

def tetrate( i, j ):
    result = i

    for x in arange( 1, j ):
        result = power( result, i )

    return result

def tetrateRight( i, j ):
    result = i

    for x in arange( 1, j ):
        result = power( i, result )

    return result

def getNthFibonorial( n ):
    result = 1

    for i in arange( 2, real( n ) ):
        result = fmul( result, fib( i ) )

    return result

def getNthAperyNumber( n ):
    result = 0

    for k in arange( 0, real( n ) + 1 ):
        result = fadd( result, fmul( power( binomial( n, k ), 2 ),
                                     power( binomial( fadd( n, k ), k ), 2 ) ) )

    return result

def getNextPrimes( p, k, func=getNextPrimeCandidateForAny ):
    result = [ ]

    for i in arange( 0, k ):
        p = getNextPrime( p, func )
        result.append( p )

    return result

def getNthMotzkinNumber( n ):
    result = 0

    for j in arange( 0, floor( fdiv( real( n ), 3 ) ) + 1 ):
        result = fadd( result, fprod( [ power( -1, j ), binomial( fadd( n, 1 ), j ),
                                      binomial( fsub( fmul( 2, n ), fmul( 3, j ) ), n ) ] ) )

    return fdiv( result, fadd( n, 1 ) )

def getNthPolytopeNumber( n, d ):
    result = real_int( n )
    m = n + 1

    for i in arange( 1, d - 1 ):
        result = fmul( result, m )
        m += 1

    return fdiv( result, fac( d - 1 ) )

def rollMultipleDice( expression, times ):
    result = [ ]

    dice = parseDiceExpression( expression )

    for i in arange( 0, times ):
        result.append( evaluateDiceExpression( dice ) )

    return result

def getNthMerten( n ):
    if real( n ) == 1:
        return 1

    result = 0

    for i in arange( 1, n + 1 ):
        result = fadd( result, getMobius( i ) )

    return result

 def PlotPos(self, fignum=None, color='g', dt=0.01, lw=2):
     tVect = arange(0, self.duration, dt)
     if tVect[-1] < self.duration:
         tVect = np.append(tVect, self.duration)
         
     xVect = [self.EvalPos(t) for t in tVect]
     if fignum is not None:
         plt.figure(fignum)
     plt.plot(tVect, xVect, color=color, linewidth=lw)
     plt.show(False)

def getNthPrimorial( n ):
    if real_int( n ) == 0:
        return 1

    result = 2

    for i in arange( 1, n ):
        result = fmul( result, getNthPrime( i + 1 ) )

    return result

def getPreviousPrimes( p, k, func = getPreviousPrimeCandidateForAny ):
    if p < 1000000000 and getNthPrime( p ) < k:
        raiseError( 'There is no previous prime to 2.' )

    result = [ ]

    for i in arange( 0, k ):
        p = getPreviousPrime( p, func )
        result.append( p )

    return result

def getNthReversalAddition( n, k ):
    next = int( real_int( n ) )
    yield next

    previous = next

    for i in arange( k ):
        if isPalindrome( next ):
            break

        next = fadd( reverseDigits( next ), next )
        yield next

def getNthSylvester( n ):
    if real( n ) == 1:
        return 2
    elif n == 2:
        return 3
    else:
        list = [ 2, 3 ]

        for i in arange( 2, n ):
            list.append( fprod( list ) + 1 )

    return list[ -1 ]

def exponentiatePolynomial( n, k ):
    if not isinstance( n, list ):
        n = [ n ]

    if isinstance( k, list ):
        return [ multiplyPolynomials( n, i ) for i in k ]
    else:
        result = n

        for i in arange( 0, k - 1 ):
            result = multiplyPolynomials( result, n )

        return result

def getNthAlternatingFactorial( n ):
    result = 0

    negative = False

    for i in arange( real( n ), 0, -1 ):
        if negative:
            result = fadd( result, fneg( fac( i ) ) )
            negative = False
        else:
            result = fadd( result, fac( i ) )
            negative = True

    return result