Ask Your Question

Revision history [back]

Black box numerical optimization of a convex function

I would like to numerically optimize a convex function without using any of its internals, that is, a python function. Is there an implementation of this in sage? I tried sage.numerical.optimize.find_local_minimum but it fails internally with a

TypeError: unable to find a common ring for all elements

My function to optimize is in fact the second largest eigenvalue of a matrix that is constructed from one parameter.

Black box numerical optimization of a convex function

I would like to numerically optimize a convex function without using any of its internals, that is, a python function. Is there an implementation of this in sage? I tried sage.numerical.optimize.find_local_minimum but it fails internally with a

TypeError: unable to find a common ring for all elements

My function to optimize is in fact the second largest eigenvalue of a matrix that is constructed from one parameter.parameter. Here it is in all gory detail:

def binaryStrings(k):
    for i in range(2**k):
        s = str(bin(i))[2:]
        yield "0"^(k-len(s)) + s;

def hammingDistanceInt (a,b):
    return bin(a^^b).count("1")

def cubeAdjacency (k):
    # Note: We would like to use QQ instead of CDF, but this will
    # result in inexact computations failing due to rounding later
    # (Yes it does!)
    # https://groups.google.com/forum/#!topic/sage-support/3IIbu0OBJX4
    return matrix(CDF, 2^k, 2^k, lambda i, j: 1 if hammingDistanceInt(i,j)==1 else 0);

def hemmeckeFiberAdjacency (k, p1):
    p2 = (1 - 4*k*p1)/2;
    M = p1 * block_diagonal_matrix (cubeAdjacency(k),cubeAdjacency(k));
    M[0,2^k] = p2
    M[2^k,0] = p2
    rowsums = M * vector ([1]*2^(k+1))
    M += diagonal_matrix(vector([1]*2^(k+1))-rowsums);
    return M

def slem(k,p):
    e = map(abs, hemmeckeFiberAdjacency(k,p).eigenvalues())
    e.sort();
    return e[-2];

def slem4(p): return slem(4,p)
plot(slem4, [0,1/16]) # Works fine
sage.numerical.optimize.find_local_minimum(slem4, [0,1/16]) # Type error