diff --git a/ackwork/ackwork_test.go b/ackwork/ackwork_test.go
index 05885dae23dff7a4d161d90d4015ee495e5d9ca8..30e1c58f35e22c819763d6eddd9e107dda1894d8 100644
--- a/ackwork/ackwork_test.go
+++ b/ackwork/ackwork_test.go
@@ -14,6 +14,20 @@ func TestErr(t *testing.T) {
}
}
+var sTypo = `
+ini_temp = 1
+final_temp 0
+x_ini 0.8,0.8,0.8
+n_step = 100000
+n_output = 5000
+f_outName testReal`
+
+func TestTypo(t *testing.T) { // Make sure a typo in input does provoke an error
+ if err := realmain(strings.NewReader(sTypo)); err == nil {
+ t.Fatal("line with typo should provoke error\n", err, "\nInput", sTypo)
+ }
+}
+
var s1d = `
ini_temp 1
final_temp 0.05
@@ -30,6 +44,23 @@ n_step 10000
n_output 500
f_outName testanneal3d`
+// Five dimensions works easily. 10 is too hard.
+var s5d = `
+ini_temp 0.15
+final_temp 0.02
+x_ini 7,7,7,7,7
+n_step 100000
+n_output 10000
+f_outName testanneal5d`
+
+var s10d = `
+ini_temp 0.1
+final_temp 0.02
+x_ini 7,7,7,7,7,7,7,7,7,7
+n_step 10000000
+n_output 10000
+f_outName testanneal10d`
+
var sCold = `
ini_temp 1e-10
final_temp 1e-12
@@ -46,12 +77,41 @@ n_step 10000
n_output 500
f_outName testhot`
-var plottable_test = []string {s1d, s3d, sHot, sCold}
-func TestPlottable (t *testing.T) {
+var sRealAnneal = `
+ini_temp 0.5
+final_temp 0.05
+x_ini 7,1,7
+n_step 100000
+n_output 5000
+f_outName testreal`
+
+// These values seem to work well. There is something of a phase
+// transition near temperature 0.09. Good values for an optimisation
+// seem to be 0.17 and 0.01 (initial and final temperatures), but this
+// needs a few steps (500000).
+
+// constant temp..
+// 0.75 jump around and nicely visits all locations
+// 0.1 jumps around, acceptance is too high and really uncontrolled.
+// 0.09375 visit 0, 1, 2, 3, 4, 5, 6, .. Can fully explore local minima and sometimes get to next
+// 0.0875 visits 6.7, 6.8, .. 7.2
+// 0.075 only see points at 6.8, 6.9, 7.0, 7.1
+// 0.05 only see points at 6.8, 6.9, 7.0, 7.1
+var sPhaseTrans = `
+ini_temp 0.18
+final_temp 0.01
+x_ini 7,1,7
+n_step 100000
+n_output 50000
+f_outName testphase`
+
+var plottable_test = []string{
+ s1d, s3d, sHot, sCold, sRealAnneal, sPhaseTrans, s5d, s10d}
+
+func TestPlottable(t *testing.T) {
for _, s := range plottable_test {
- if err := realmain(strings.NewReader(sHot)); err != nil {
+ if err := realmain(strings.NewReader(s)); err != nil {
t.Fatal("plottable failed with\n", err, "\nInput", s)
}
}
}
-
diff --git a/ackwork/dorun.go b/ackwork/dorun.go
index 7463c440d3298cbd529c226f5a3db0c2013358c6..4f0879debd7931097485fee1adfca6801f7dbd66 100644
--- a/ackwork/dorun.go
+++ b/ackwork/dorun.go
@@ -3,9 +3,9 @@ package ackwork
import (
"fmt"
+ "io"
"math"
"math/rand"
- "io"
"os"
"sync"
@@ -40,7 +40,6 @@ func setupRun(mcPrm *mcPrm, cprm *cprm) error {
}
cprm.rand = rand.New(rand.NewSource(getSeed()))
-
if mcPrm.iniTmp != mcPrm.fnlTmp {
var coolrate float64
cprm.coolme = true
@@ -49,6 +48,10 @@ func setupRun(mcPrm *mcPrm, cprm *cprm) error {
}
nStep := float64(mcPrm.nStep)
fnlTmp, iniTmp := float64(mcPrm.fnlTmp), float64(mcPrm.iniTmp)
+ if fnlTmp < 0.01 {
+ breaker()
+ fnlTmp = math.SmallestNonzeroFloat32
+ }
coolrate = -1 / nStep
coolrate *= math.Log(fnlTmp / iniTmp)
cprm.coolMult = math.Exp(-coolrate)
@@ -57,22 +60,22 @@ func setupRun(mcPrm *mcPrm, cprm *cprm) error {
cprm.nEvery = uint32(math.Ceil(float64(mcPrm.nStep) / float64(mcPrm.nOutput)))
// change the next line to handle the case of multiple runs
ntmp := mcPrm.fOutName + ".csv"
- if cprm.fOut, err = os.Create (ntmp); err != nil {
+ if cprm.fOut, err = os.Create(ntmp); err != nil {
return err
}
-
+
return nil
}
// newx gets a candidate x slice
func newx(xold []float32, xT []float32, rand *rand.Rand, xDlta float32) {
for i, x := range xold {
- t := 2. * rand.Float32() - 1
+ t := 2.*rand.Float32() - 1
t *= xDlta
xT[i] = x + t
}
}
-
+func breaker() {}
func doRun(mcPrm *mcPrm) error {
var cprm cprm
if err := setupRun(mcPrm, &cprm); err != nil {
@@ -81,27 +84,48 @@ func doRun(mcPrm *mcPrm) error {
if cprm.fOut != nil {
defer cprm.fOut.Close()
}
- var nAcc int // number of accepted moves
- x := make([]float32, len(mcPrm.xIni))
- xT := make([]float32, len(mcPrm.xIni))
- nout := uint32(1)
-
+ const accRateIdeal = 0.05 // Arbitrary target of 5% acceptance
+ const nRunAccAdj = 200 // every n steps, check acceptance rate
+ const maxxDlta = 2. // max step size
+ var nAcc int // Counter of number of accepted moves
+ x := make([]float32, len(mcPrm.xIni)) // current position
+ xT := make([]float32, len(mcPrm.xIni)) // trial position
+ nout := uint32(1) // counter for printing output
+ runAcc := accRateIdeal // Running acceptance rate, exponentially weighted
copy(x, mcPrm.xIni)
- fOld := ackley.Ackley(x)
+ fOld := ackley.Ackley(x) // Initial function value
tmprtr := float64(mcPrm.iniTmp)
- var n uint32
- for ; n < mcPrm.nStep; n++ {
+ nRunAcc := nRunAccAdj // Every nRunAccAdj, try adjusting the step size.
+ const runMult = 0.99
+ xDlta := mcPrm.xDlta // Adaptable step size
+ for n := uint32(0); n < mcPrm.nStep; n++ {
var acc bool
nout--
- if nout == 0 {
- fmt.Fprintf (cprm.fOut, "%d,%.4g,%.5g", n, tmprtr,fOld)
- for _,xx := range x {
- fmt.Fprintf(cprm.fOut, ",%.2g", xx)}
+ if nout == 0 { // Do we want to print out results on this step ?
+ fmt.Fprintf(cprm.fOut, "%d,%.4g,%.5g", n, tmprtr, fOld)
+ for _, xx := range x {
+ fmt.Fprintf(cprm.fOut, ",%.2g", xx)
+ }
fmt.Fprintln(cprm.fOut)
nout = cprm.nEvery
}
-
- newx(x, xT, cprm.rand, mcPrm.xDlta)
+ nRunAcc--
+ if nRunAcc == 0 { // Do we want to try adjusting the step size ?
+ const step = "step"
+ if runAcc < accRateIdeal-0.01 { // acceptance rate too low
+ xDlta *= 0.9
+ if mcPrm.verbose {
+ fmt.Println(step, xDlta, "decrease")
+ }
+ } else if (xDlta < maxxDlta) && (runAcc > accRateIdeal+0.01) {
+ xDlta *= 1.1
+ if mcPrm.verbose {
+ fmt.Println(step, xDlta, "increase")
+ }
+ }
+ nRunAcc = nRunAccAdj
+ }
+ newx(x, xT, cprm.rand, xDlta)
fTrial := ackley.Ackley(xT)
if fOld < fTrial {
@@ -116,12 +140,16 @@ func doRun(mcPrm *mcPrm) error {
}
if acc == true {
nAcc++
+ runAcc = runMult*runAcc + (1.0 - runMult)
copy(x, xT)
fOld = fTrial
+ } else {
+ runAcc = runMult * runAcc
}
if cprm.coolme {
tmprtr *= cprm.coolMult
}
}
+ fmt.Println("n accepted:", nAcc, "of", mcPrm.nStep+1)
return nil
}
diff --git a/ackwork/plotres.r b/ackwork/plotres.r
new file mode 100644
index 0000000000000000000000000000000000000000..bada2a0c0c1472000c4e8d17e695bca2e79b160b
--- /dev/null
+++ b/ackwork/plotres.r
@@ -0,0 +1,10 @@
+# Obviously this only makes sense for local testing of some of the results
+# We can plot out the column f_val as a function of temperature, x-coordinates, ...
+
+names <- c("step", "temperature", "f_val", "x1", "x2", "x3", "x4", "x5", "x6", "x7", "x8", "x9", "x10")
+df <- read.csv("testanneal5d.csv", col.names = names)
+
+
+plot (xlab = "step", ylab = "f val", x = df$step, y = df$f_val,
+ type = "p")
+
diff --git a/ackwork/prm_file.go b/ackwork/rdprm.go
similarity index 100%
rename from ackwork/prm_file.go
rename to ackwork/rdprm.go