Golf a Brain-Flak Integer Integers are tedious to represent in Brain-Flak. There are 8 operators:

()      Evaluates to 1, but does not push anything on any stack
[]      Evaluates to an indeterminate value for the purposes of this question
{}      Removes the top of the stack and evaluates to it
&lt;&gt;      Switches to or back from the alternate stack and evaluates to zero
(foo)   Pushes the value of the expression foo to the stack and evaluates to it
[foo]   Evaluates to the negation of foo
{foo}   Evaluates the expression foo until the top of the stack is zero
&lt;foo&gt;   Evaluates to zero but executes foo anyway


foo may consist of multiple operators, in which case they are evaluated and summed. For example (()()) pushes 2 to the stack (and evaluates to 2 too).

Obviously, the (()...()) mechanism is not useful in Code Golf as large numbers would take n*2+2 bytes to represent. Your challenge is therefore to write a program or function that will output in as few bytes as possible a Brain-Flak program that will push a given positive integer n to the active stack. This program must not make any assumptions about the existing contents of the stacks, so it must not leave the stacks exchanged or add or remove extra values from the stacks.

Although your program or function must be capable of returning a working Brain-Flak program for all inputs from 1 to 1,000,000 the winner will be the program or function that generates the smallest set of appropriate Brain-Flak programs for all 1061 prime numbers between 1,000 and 10,000. You should note the total size of your outputs for those 1061 inputs as part of your submission. Your program or function may accept the integer and return the (string) Brain-Flak program in any of the usual acceptable I/O formats. Ties will be broken using the size of your program or function.
Python 2, 59394 59244 58534 58416 58394 58250

Ok here is my solution.
import re
import math

cache = {0:&quot;&lt;()&gt;&quot;}

def find(x,i,j):
    return i*((x**2+x)/2)+(j+1)*((x**2-x)/2)

def solve(x, i, j):
    a = (i + j + 1)/2.
    b = (i - j - 1)/2.
    c = -x
    return (-b + math.sqrt(b**2 - 4*a*c))/(2*a)

def size(i,j=0):
    return 4*(i+j)+14

def polynomials(n):
    upperBound = int(4*math.log(n,2))
    i = 0
    answers = []
    while size(i) &lt; upperBound:
        for j in range(i):
            sol = int(solve(n, i-j, j)+.5)
            if find(sol, i-j, j) == n:
                answers.append((sol, i-j, j))
        i += 1
    return answers

def complement(character):
        dict = {&quot;(&quot;:&quot;)&quot;,&quot;{&quot;:&quot;}&quot;,&quot;&lt;&quot;:&quot;&gt;&quot;,&quot;[&quot;:&quot;]&quot;,&quot;)&quot;:&quot;(&quot;,&quot;}&quot;:&quot;{&quot;,&quot;&gt;&quot;:&quot;&lt;&quot;,&quot;]&quot;:&quot;[&quot;}
        return dict[character]

def findMatch(snippet, index):
        increment = 1 if snippet[index] in &quot;({&lt;[&quot; else -1
        stack = []
        if snippet[index] in &quot;(){}&lt;&gt;[]&quot;:
                stack.append(snippet[index])
        while len(stack) &gt; 0 and index + increment &lt; len(snippet):
                index += increment
                if snippet[index] in &quot;(){}&lt;&gt;[]&quot;:
                        if complement(snippet[index]) == stack[-1]:
                                stack = stack[:-1]
                        else:
                                stack.append(snippet[index])
        return index

def isPrime(n):
    return not [0 for x in range(2,int(n**.5)+1) if n%x==0] and n&gt;1

def getPrimeFactors(n):
    return [x for x in range(2,n/2) if n%x==0 and isPrime(x)]

def divHardcode(n,m):
    assert n%m == 0
    assert m != 1
    assert n != 1
    binary = bin(m)[3:]
    return (binary.count(&quot;1&quot;)+len(binary))*&quot;(&quot;+getBF(n/m)+&quot;)&quot;*binary.count(&quot;1&quot;)+binary.replace(&quot;1&quot;,&quot;){}{}&quot;).replace(&quot;0&quot;,&quot;){}&quot;)

def isTriangular(n):
    #Triangles must be between sqrt(2n) and cbrt(2n)
    if n &lt; 0: return isTriangular(-n)
    for x in range(int((2*n)**(1/3.)),int((2*n)**.5)+1):
        if (x**2+x) == 2*n:
            return True
    return False

def getTriangle(n):
    if n &lt; 0: return -getTriangle(-n)
    #Triangles must be between sqrt(2n) and cbrt(2n)
    for x in range(int((2*n)**(1/3.)),int((2*n)**.5)+1):
        if (x**2+x) == 2*n:
            return x
    #If we don't find one we made a mistake
    assert False

def getSimpleBF(n):
    if n in cache:return cache[n]
    if n &lt; 0:
        # There is room for better solutions here
        return &quot;[&quot;+getSimpleBF(-n)+&quot;]&quot;
    elif n == 0:
        return &quot;&quot;
    elif n &lt; 6:
        return &quot;()&quot;*n
    #Non-edge cases
    solutions = []
    factors = getPrimeFactors(n)
    if n &gt;= 78 and isTriangular(n):
        solutions.append(
           min([push(getTriangle(n))+&quot;{({}[()])}{}&quot;,&quot;&lt;&quot;+push(getTriangle(n)+1)+&quot;&gt;{({}[()])}{}&quot;],key=len)
        )
    polynomialSolutions = polynomials(n)
    for polynomial in polynomialSolutions:
        solutions.append(&quot;&lt;%s&gt;{%s({}[()])%s}{}&quot;%(push(polynomial[0]),&quot;({})&quot;*polynomial[1],&quot;({})&quot;*polynomial[2]))
        #Mod 3 tricks
    if n % 3 == 2:
       solutions.append((&quot;((%s)()){}{}&quot;)%getBF(n/3))
    elif n % 3 == 1:
       solutions.append((&quot;((%s)()()){}{}&quot;)%getBF(n/3-1))
    #Basic solutions
    if isPrime(n):
        solutions.append(getSimpleBF(n-1) + &quot;()&quot;)
    else:
        #TODO multithread
        solutions += map(lambda m:divHardcode(n,m),factors)
    return min(solutions,key=lambda x:len(unpack(x)))

def getBF(n):
    if n in cache: return cache[n]
    result = getSimpleBF(n)
    index = n - 1
    while index &gt; n-(len(result)/2):
        score = getSimpleBF(index)+getSimpleBF(n-index)
        if len(score) &lt; len(result):result = score
        index -= 1
    index = n + 1
    while index &lt; n+(len(result)/2):
        score = getSimpleBF(index)+getSimpleBF(n-index)
        if len(score) &lt; len(result):result = score
        index += 1
    cache[n] = result
    return result

def unpack(string):
    reMatch = re.match(&quot;\(*&lt;&quot;,string)
    if reMatch:
        location =reMatch.span()
        return string[location[1]:findMatch(string,location[1]-1)] +string[:location[1]-1] + string[findMatch(string,location[1]-1)+1:]
    return string

def push(n):
    return unpack(&quot;(&quot;+getBF(n)+&quot;)&quot;)

def kolmo(string):
    code = push(ord(string[-1]))
    stringVector = map(ord,string)
    for x,y in zip(stringVector[-1:0:-1],stringVector[-2::-1]):
        code = &quot;(&quot;+code+getBF(y-x)+&quot;)&quot;
    code = code.replace(&quot;&lt;()&gt;)&quot;,&quot;)&quot;)
    return code

def kolmo(stringVector):
    code = push(stringVector[-1])
    for x,y in zip(stringVector[-1:0:-1],stringVector[-2::-1]):
        code = &quot;(&quot;+code+getBF(y-x)+&quot;)&quot;
    code = code.replace(&quot;&lt;()&gt;)&quot;,&quot;)&quot;)
    return code


if __name__ == &quot;__main__&quot;:
    import primes
    sum = 0
    for prime in primes.nums:
        print push(prime)
        sum += len(push(prime))
    print sum

The relevant function is push(n). To call it simply call push on the integer you wish to represent.
Explanation
The main optimization done by the program is multiplication hard-coding.  The idea of multiplication hardcoding is pretty simple.  You push the a number and then pop and push it to create a new value.  For example to multiply by two you can use the following code ((n){}) where n code producing a specific number.  This works because both (n) and {} have a value of n.
This simple idea can be made more complex for larger numbers.  Take for instance 5 it was discovered a while ago that the best way to multiply by five was (((n)){}){}{}.  This code makes two copies of the n multiplies one by 4 and adds the two.  Using the same strategy I make every multiplication based on a number's binary representation.  I won't get into the details of how this works right now but I do this by chopping off the first one of the binary representation and replacing 0 with ){} and 1 with ){}{}.  It then makes sure that n is pushed the appropriate number of times and balances all the parentheses.  (If you want to know how this is done you can look at my code).  If you want to know why this works just ask me in a comment.  I don't think anyone actually reads all the updates to my post so I left the explanation out.
When the algorithm attempts to find a multiplication hardcode it attempts all of a numbers prime factors.  It ignores composite factors because at one point composite factors could always be expressed more concisely as its own prime factors it is not known if this is still true.
The other byte saving mechanism is a polynomial solution finder.  There are certain forms of polynomials that are easy to represent with decrementing loops.  These polynomials include, but are not limited to, polygonal numbers. This optimization finds polynomials that fit the form and creates the code that makes them.
Output
paste-bin