/
torque.go
191 lines (176 loc) · 6.94 KB
/
torque.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
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
package engine
import (
"reflect"
"github.com/mumax/3/cuda"
"github.com/mumax/3/data"
"github.com/mumax/3/util"
)
var (
Alpha = NewScalarParam("alpha", "", "Landau-Lifshitz damping constant")
Xi = NewScalarParam("xi", "", "Non-adiabaticity of spin-transfer-torque")
Pol = NewScalarParam("Pol", "", "Electrical current polarization")
Lambda = NewScalarParam("Lambda", "", "Slonczewski Λ parameter")
EpsilonPrime = NewScalarParam("EpsilonPrime", "", "Slonczewski secondairy STT term ε'")
FrozenSpins = NewScalarParam("frozenspins", "", "Defines spins that should be fixed") // 1 - frozen, 0 - free. TODO: check if it only contains 0/1 values
FreeLayerThickness = NewScalarParam("FreeLayerThickness", "m", "Slonczewski free layer thickness (if set to zero (default), then the thickness will be deduced from the mesh size)")
SOThickness = NewScalarParam("SOThickness", "m", "SOT freelayer thickness") // 自己加的
HFLoverHDL = NewScalarParam("HFLoverHDL", "", "SOT HFL/HDL") // 自己加的
SPINHALL = NewScalarParam("SPINHALL", "", "SOT spin hall angle") // 自己加的
FixedLayer = NewExcitation("FixedLayer", "", "Slonczewski fixed layer polarization")
SOTFixedLayer = NewExcitation("SOTFixedLayer", "", "SOT sigma") // 自己加的
Torque = NewVectorField("torque", "T", "Total torque/γ0", SetTorque)
LLTorque = NewVectorField("LLtorque", "T", "Landau-Lifshitz torque/γ0", SetLLTorque)
STTorque = NewVectorField("STTorque", "T", "Spin-transfer torque/γ0", AddSTTorque)
SOTorque = NewVectorField("SOTorque", "T", "Spin-orbit torque/γ0", AddSOTorque) // 自己加的
J = NewExcitation("J", "A/m2", "Electrical current density")
JSOT = NewExcitation("JSOT", "A/m2", "SOT Electrical current density") // 自己加的
MaxTorque = NewScalarValue("maxTorque", "T", "Maximum torque/γ0, over all cells", GetMaxTorque)
GammaLL float64 = 1.7595e11 // Gyromagnetic ratio of spins, in rad/Ts
Precess = true
DisableZhangLiTorque = false
DisableSlonczewskiTorque = false
DisableSOTorque = false // 自己加的
fixedLayerPosition = FIXEDLAYER_TOP // instructs mumax3 how free and fixed layers are stacked along +z direction
)
func init() {
Pol.setUniform([]float64{1}) // default spin polarization
Lambda.Set(1) // sensible default value (?).
DeclVar("GammaLL", &GammaLL, "Gyromagnetic ratio in rad/Ts")
DeclVar("DisableZhangLiTorque", &DisableZhangLiTorque, "Disables Zhang-Li torque (default=false)")
DeclVar("DisableSlonczewskiTorque", &DisableSlonczewskiTorque, "Disables Slonczewski torque (default=false)")
DeclVar("DoPrecess", &Precess, "Enables LL precession (default=true)")
DeclLValue("FixedLayerPosition", &flposition{}, "Position of the fixed layer: FIXEDLAYER_TOP, FIXEDLAYER_BOTTOM (default=FIXEDLAYER_TOP)")
DeclROnly("FIXEDLAYER_TOP", FIXEDLAYER_TOP, "FixedLayerPosition = FIXEDLAYER_TOP instructs mumax3 that fixed layer is on top of the free layer")
DeclROnly("FIXEDLAYER_BOTTOM", FIXEDLAYER_BOTTOM, "FixedLayerPosition = FIXEDLAYER_BOTTOM instructs mumax3 that fixed layer is underneath of the free layer")
}
// Sets dst to the current total torque
func SetTorque(dst *data.Slice) {
SetLLTorque(dst)
AddSTTorque(dst)
AddSOTorque(dst) // 自己加的
FreezeSpins(dst)
}
// Sets dst to the current Landau-Lifshitz torque
func SetLLTorque(dst *data.Slice) {
SetEffectiveField(dst) // calc and store B_eff
alpha := Alpha.MSlice()
defer alpha.Recycle()
if Precess {
cuda.LLTorque(dst, M.Buffer(), dst, alpha) // overwrite dst with torque
} else {
cuda.LLNoPrecess(dst, M.Buffer(), dst)
}
}
// Adds the current spin transfer torque to dst
func AddSTTorque(dst *data.Slice) {
if J.isZero() {
return
}
util.AssertMsg(!Pol.isZero(), "spin polarization should not be 0")
jspin, rec := J.Slice()
if rec {
defer cuda.Recycle(jspin)
}
fl, rec := FixedLayer.Slice()
if rec {
defer cuda.Recycle(fl)
}
if !DisableZhangLiTorque {
msat := Msat.MSlice()
defer msat.Recycle()
j := J.MSlice()
defer j.Recycle()
alpha := Alpha.MSlice()
defer alpha.Recycle()
xi := Xi.MSlice()
defer xi.Recycle()
pol := Pol.MSlice()
defer pol.Recycle()
cuda.AddZhangLiTorque(dst, M.Buffer(), msat, j, alpha, xi, pol, Mesh())
}
if !DisableSlonczewskiTorque && !FixedLayer.isZero() {
msat := Msat.MSlice()
defer msat.Recycle()
j := J.MSlice()
defer j.Recycle()
fixedP := FixedLayer.MSlice()
defer fixedP.Recycle()
alpha := Alpha.MSlice()
defer alpha.Recycle()
pol := Pol.MSlice()
defer pol.Recycle()
lambda := Lambda.MSlice()
defer lambda.Recycle()
epsPrime := EpsilonPrime.MSlice()
defer epsPrime.Recycle()
thickness := FreeLayerThickness.MSlice()
defer thickness.Recycle()
cuda.AddSlonczewskiTorque2(dst, M.Buffer(),
msat, j, fixedP, alpha, pol, lambda, epsPrime,
thickness,
CurrentSignFromFixedLayerPosition[fixedLayerPosition],
Mesh())
}
}
func FreezeSpins(dst *data.Slice) {
if !FrozenSpins.isZero() {
cuda.ZeroMask(dst, FrozenSpins.gpuLUT1(), regions.Gpu())
}
}
func GetMaxTorque() float64 {
torque := ValueOf(Torque)
defer cuda.Recycle(torque)
return cuda.MaxVecNorm(torque)
}
type FixedLayerPosition int
const (
FIXEDLAYER_TOP FixedLayerPosition = iota + 1
FIXEDLAYER_BOTTOM
)
var (
CurrentSignFromFixedLayerPosition = map[FixedLayerPosition]float64{
FIXEDLAYER_TOP: 1.0,
FIXEDLAYER_BOTTOM: -1.0,
}
)
type flposition struct{}
func (*flposition) Eval() interface{} { return fixedLayerPosition }
func (*flposition) SetValue(v interface{}) {
drainOutput()
fixedLayerPosition = v.(FixedLayerPosition)
}
func (*flposition) Type() reflect.Type { return reflect.TypeOf(FixedLayerPosition(FIXEDLAYER_TOP)) }
// 自己加的
func AddSOTorque(dst *data.Slice) {
if JSOT.isZero() {
return
}
util.AssertMsg(!SPINHALL.isZero(), "spin hall angle should not be 0")
jsotspin, sotrec := J.Slice()
if sotrec {
defer cuda.Recycle(jsotspin)
}
sotfl, sotrec := SOTFixedLayer.Slice()
if sotrec {
defer cuda.Recycle(sotfl)
}
if !DisableSOTorque {
msat := Msat.MSlice()
defer msat.Recycle()
jsot := JSOT.MSlice()
defer jsot.Recycle()
alpha := Alpha.MSlice()
defer alpha.Recycle()
sotfixedP := SOTFixedLayer.MSlice()
defer sotfixedP.Recycle()
spinhall := SPINHALL.MSlice()
defer spinhall.Recycle()
hfloverhdl := HFLoverHDL.MSlice()
defer hfloverhdl.Recycle()
sothickness := SOThickness.MSlice()
defer sothickness.Recycle()
cuda.AddSOTorque(dst, M.Buffer(),
msat, jsot, sotfixedP, alpha, spinhall, hfloverhdl,
sothickness, Mesh())
}
}