forked from mumax/3
/
minimizer.go
162 lines (132 loc) · 3.28 KB
/
minimizer.go
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
package engine
// Minimize follows the steepest descent method as per Exl et al., JAP 115, 17D118 (2014).
import (
"github.com/mumax/3/cuda"
"github.com/mumax/3/data"
)
var (
DmSamples int = 10 // number of dm to keep for convergence check
StopMaxDm float64 = 1e-6 // stop minimizer if sampled dm is smaller than this
)
func init() {
DeclFunc("Minimize", Minimize, "Use steepest conjugate gradient method to minimize the total energy")
DeclVar("MinimizerStop", &StopMaxDm, "Stopping max dM for Minimize")
DeclVar("MinimizerSamples", &DmSamples, "Number of max dM to collect for Minimize convergence check.")
}
// fixed length FIFO. Items can be added but not removed
type fifoRing struct {
count int
tail int // index to put next item. Will loop to 0 after exceeding length
data []float64
}
func FifoRing(length int) fifoRing {
return fifoRing{data: make([]float64, length)}
}
func (r *fifoRing) Add(item float64) {
r.data[r.tail] = item
r.count++
r.tail = (r.tail + 1) % len(r.data)
if r.count > len(r.data) {
r.count = len(r.data)
}
}
func (r *fifoRing) Max() float64 {
max := r.data[0]
for i := 1; i < r.count; i++ {
if r.data[i] > max {
max = r.data[i]
}
}
return max
}
type Minimizer struct {
k *data.Slice // torque saved to calculate time step
lastDm fifoRing
h float32
}
func (mini *Minimizer) Step() {
m := M.Buffer()
size := m.Size()
if mini.k == nil {
mini.k = cuda.Buffer(3, size)
torqueFn(mini.k)
}
k := mini.k
h := mini.h
// save original magnetization
m0 := cuda.Buffer(3, size)
defer cuda.Recycle(m0)
data.Copy(m0, m)
// make descent
cuda.Minimize(m, m0, k, h)
// calculate new torque for next step
k0 := cuda.Buffer(3, size)
defer cuda.Recycle(k0)
data.Copy(k0, k)
torqueFn(k)
setMaxTorque(k) // report to user
// just to make the following readable
dm := m0
dk := k0
// calculate step difference of m and k
cuda.Madd2(dm, m, m0, 1., -1.)
cuda.Madd2(dk, k, k0, -1., 1.) // reversed due to LLNoPrecess sign
// get maxdiff and add to list
max_dm := cuda.MaxVecNorm(dm)
mini.lastDm.Add(max_dm)
setLastErr(mini.lastDm.Max()) // report maxDm to user as LastErr
// adjust next time step
var nom, div float32
if NSteps%2 == 0 {
nom = cuda.Dot(dm, dm)
div = cuda.Dot(dm, dk)
} else {
nom = cuda.Dot(dm, dk)
div = cuda.Dot(dk, dk)
}
if div != 0. {
mini.h = nom / div
} else { // in case of division by zero
mini.h = 1e-4
}
M.normalize()
// as a convention, time does not advance during relax
NSteps++
}
func (mini *Minimizer) Free() {
mini.k.Free()
}
func Minimize() {
Refer("exl2014")
SanityCheck()
// Save the settings we are changing...
prevType := solvertype
prevFixDt := FixDt
prevPrecess := Precess
t0 := Time
relaxing = true // disable temperature noise
// ...to restore them later
defer func() {
SetSolver(prevType)
FixDt = prevFixDt
Precess = prevPrecess
Time = t0
relaxing = false
}()
Precess = false // disable precession for torque calculation
// remove previous stepper
if stepper != nil {
stepper.Free()
}
// set stepper to the minimizer
mini := Minimizer{
h: 1e-4,
k: nil,
lastDm: FifoRing(DmSamples)}
stepper = &mini
cond := func() bool {
return (mini.lastDm.count < DmSamples || mini.lastDm.Max() > StopMaxDm)
}
RunWhile(cond)
pause = true
}