|
#!/usr/bin/env python2 |
|
""" bump2 optimization problem |
|
from Keane, https://www.southampton.ac.uk/~ajk/bump2.html 1994 |
|
MAXimize |
|
{ abs ( sum i=1 to n { cos^4 ( x_i ) } |
|
- 2 prod i=1 to n { cos^2 ( x_i ) } |
|
) } |
|
/ { sqrt { sum i=1 to n { i x_i^2 } } } |
|
|
|
constraints: |
|
0 < x_i < 10 |
|
av x_i < 7.5 |
|
prod x_i > .75 |
|
|
|
n=20: 0.76 ~ 1000 n fevals (derivs ?) |
|
|
|
What's here: |
|
bump2(x) = - Keane^2 -- optimizers don't like sharp corners, sqrt |
|
random_bumpx( nrand, dim ) -- random x with product .75 etc. |
|
|
|
run this.py dim=4 nrand=1e4: .622 .621 ... |
|
run this.py dim=20 nrand=1e6: .725 .725 ... |
|
nrand=1e7: .745 .740 ... 7 minutes |
|
|
|
Why ? Get an overview of rough terrain, *before* zooming in |
|
and fun with constraints: |
|
notice that if xi < xj, swapping them decreases the denominator |
|
while leaving the numerator unchanged; so, constrain x0 >= x1 >= x2 ... |
|
|
|
Vanaret et al. "Certified Global Minima for a Benchmark of Difficult Optimization Problems" 2014 |
|
code not on the web ? |
|
x = np.r_[ 3.065318, 1.531047, 0.405617, 0.393987 ] |
|
f2 = abs( bump2( x )) |
|
print "bump2( Vanaret %s ): %.3g = sqrt %.3g" % ( # .622 = sqrt .387 |
|
x, np.sqrt(f2), f2 ) |
|
""" |
|
|
|
|
|
from __future__ import division |
|
import numpy as np |
|
from numpy import cos, pi, prod, sin, sqrt, sum |
|
|
|
#............................................................................... |
|
def bump2( x ): |
|
""" Keane's bump function squared, - to minimize |
|
https://www.southampton.ac.uk/~ajk/bump2.html 1994 |
|
""" |
|
x = np.asarray(x) |
|
n = len(x) |
|
cos2 = cos(x) ** 2 |
|
num2 = (sum( cos2**2 ) - 2 * prod( cos2 )) ** 2 |
|
denom2 = np.arange( 1., n+1 ) .dot( x**2 ) |
|
return - num2 / denom2 # - to minimize |
|
|
|
def random_bumpx( nrand, n ): |
|
""" -> nrand x n random for bump2: prod .75, trim, sort down """ |
|
X = np.random.uniform( high=2, size=[int(nrand), n] ).astype( np.float32 ) |
|
prods = X.prod( axis=1 ) |
|
scale = (.75 / prods) ** (1/n) |
|
X *= scale[:,np.newaxis] # scale rows |
|
rowavs = X.mean( axis=1 ) |
|
rowmax = X.max( axis=1 ) |
|
X = X[ (rowavs < 7.5) & (rowmax < 10) ] # trim < 1 % |
|
X = np.sort( X, axis=1 )[:,::-1] # sort down |
|
return X |
|
|
|
__version__ = "2018-08-24 Aug denis-bz-py t-online de" |
|
|
|
#............................................................................... |
|
if __name__ == "__main__": # standalone test |
|
import sys |
|
|
|
np.set_printoptions( threshold=20, edgeitems=10, linewidth=140, |
|
formatter = dict( float = lambda x: "%.3g" % x )) # float arrays %.3g |
|
print "\n", 80 * "=" |
|
|
|
dim = 4 |
|
nrand = 1e4 |
|
seed = 3 |
|
tag = "" # top few > tag.nptxt |
|
|
|
# to change these params, run this.py a=1 b=None c='"str"' ... in sh or ipython |
|
for arg in sys.argv[1:]: |
|
exec( arg ) # name not defined error: forgot '' ? |
|
|
|
np.random.seed( seed ) |
|
if not tag: |
|
tag = "%dd-%g-bump" % (dim, nrand) |
|
|
|
#........................................................................... |
|
Xrand = random_bumpx( nrand, dim ) |
|
info = "bump2: dim %d nrand %g - %d %.2g bytes " % ( |
|
dim, nrand, nrand - len(Xrand), Xrand.nbytes ) |
|
print info |
|
print "average Xrand:", Xrand.mean( axis=0 ) |
|
|
|
#........................................................................... |
|
f2 = np.fabs( map( bump2, Xrand )) # [bump2(row0), bump2(row1) ...] |
|
f = np.sqrt( f2 ) |
|
jsort = f.argsort()[::-1] [:10] |
|
|
|
print "# j f x " |
|
for j in jsort: |
|
x = Xrand[j] |
|
print "%8d %-8.3g %s " % ( |
|
j, f[j], x ) |
|
|
|
# save the top few, plot-bump.py -- |
|
fx = np.c_[ f[jsort], Xrand[jsort] ] |
|
out = tag + ".nptxt" |
|
print "\nsaving", out |
|
np.savetxt( out, fx, fmt="%-4.3g", delimiter="\t", |
|
header = "f x " + info ) |
|
|
