forked from tsdev/spinw
-
Notifications
You must be signed in to change notification settings - Fork 15
/
genmagstr.m
517 lines (471 loc) · 18.5 KB
/
genmagstr.m
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
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
function genmagstr(obj, varargin)
% generates magnetic structure
%
% ### Syntax
%
% `genmagstr(obj,Name,Value)`
%
% ### Description
%
% `genmagstr(obj,Name,Value)` is a Swiss knife for generating magnetic
% structures. It allows the definition of magnetic structures using several
% different ways, depending on the `mode` parameter. The generated magnetic
% structure is stored in the [obj.mag_str] field. The magetic structure is
% stored as Fourier components with arbitrary number of wave vectors in the
% [spinw] object. However spin waves can be only calculated if the magnetic
% structure has a single propagation vector (plus a k=0 ferromagnetic
% component perpendicular to the incommensurate magnetization), we simply
% call it single-k magnetic structure. Thus `genmagstr` enables to input
% magnetic structures that comply with this restriction by defining a
% magnetic structure by the moment directions (`S`) in the crystallographic
% cell, a propagation vector (`km`) and a vector that defines the normal of
% the rotation of the spin spiral (`n`). The function converts these values
% into Fourier components to store. To solve spin wave dispersion of
% magnetic structures with multiple propagation vectors, a magnetic
% supercell has to be defined where the propagation vector can be
% approximated to zero.
%
% ### Examples
%
% The example creates the multi-k magnetic structure of USb given by the
% `FQ` Fourier components and plots the magnetic structure:
%
% ```
% >>USb = spinw
% >>USb.genlattice('lat_const',[6.203 6.203 6.203],'spgr','F m -3 m')
% >>USb.addatom('r',[0 0 0],'S',1)
% >>FQ = cat(3,[0;0;1+1i],[0;1+1i;0],[1+1i;0;0])>>
% >>k = [0 0 1;0 1 0;1 0 0];
% >>USb.genmagstr('mode','fourier','S',FQ,'nExt',[1 1 1],'k',k)
% >>plot(USb,'range',[1 1 1])
% >>snapnow
% ```
%
% ### Input Arguments
%
% `obj`
% : [spinw] object.
%
% ### Name-Value Pair Arguments
%
% `'mode'`
% : Mode that determines how the magnetic structure is generated:
% * `'random'` (optionally reads `k`, `n`, `nExt`)
% generates a random structure in the structural cell if no other
% arguments are specified here or previously in this spinw
% object. If `nExt` is specified all spins in the supercell are
% randomised. If `k` is specified a random helical structure with
% moments perpendicular to `n` (default value: `[0 0 1]`) with
% the specified `k` propagation vector is generated. (`n` is not
% otherwise used).
% * `'direct'` (reads `S`, optionally reads `k`, `nExt`)
% direct input of the magnetic structure using the
% parameters of the single-k magnetic structure.
% * `'tile'` (reads `S`, optionally reads `nExt`)
% Simply extends the existing or input structure
% (`S`) into a magnetic supercell by replicating it.
% If no structure is stored in the [spinw] object a random
% structure is generated automatically. If defined,
% `S` is used as starting structure for extension
% overwriting the stored structure. If the original
% structure is already extended with other size, only the
% moments in the crystallographic cell wil be replicated.
% Magnetic ordering wavevector `k` will be set to zero. To
% generate structure with non-zero k, use the `'helical'` or
% `'direct'` option.
% * `'helical'` (reads `S`, optionally reads `n`, `k`, `r0`, `nExt`, `epsilon`)
% generates helical structure in a single cell or in a
% supercell. In contrary to the `'extend'` option the
% magnetic structure is not generated by replication but
% by rotation of the moments using the following formula:
%
% $\mathbf{S}^{gen}_i(\mathbf{r}) = R(2 \pi \mathbf{k_m} \cdot \mathbf{r})\cdot \mathbf{S}_i$
%
% where $S_i$ has either a single moment or as many moments
% as the number of magnetic atoms in the crystallographic
% cell. In the first case $r$ denotes the atomic
% positions, while for the second case $r$ denotes the
% position of the origin of the cell.
% * `'rotate'` (optionally reads `S`, `phi`, `phid`, `n`, `nExt`)
% uniform rotation of all magnetic moments with a
% `phi` angle around the given `n` vector. If
% `phi=0`, all moments are rotated so, that the first
% moment is parallel to `n` vector in case of
% collinear structure or in case of planar structure
% `n` defines the normal of the plane of the magnetic
% moments.
% * `'func'` (reads `func`, `x0`)
% function that defines the parameters of the single-k
% magnetic structure: moment vectors, propagation vector
% and normal vector from arbitrary parameters in the
% following form:
% ```
% [S, k, n] = @(x)func(S0, x)
% ```
% where `S` is matrix with dimensions of $[3\times n_{magExt}]$. `k` is
% the propagation vector in a 3-element row vector. `n` is the
% normal vector of the spin rotation plane also 3-element row
% vector. The default value for `func` is `@gm_spherical3d`. For planar
% magnetic structure use `@gm_planar`. Only `func` and `x`
% have to be defined for this mode.
% * `'fourier'` (reads `S`, optionally reads `k`, `r0`, `nExt`, `epsilon`)
% same as `'helical'`, except the `S` option is taken as the
% Fourier components, thus if it contains real numbers, it will
% generate a sinusoidally modulated structure instead of
% a spiral.
%
% `'phi'`
% : Angle of rotation of the magnetic moments in radian. Default
% value is 0.
%
% `'phid'`
% : Angle of rotation of the magnetic moments in \\deg. Default
% value is 0.
%
% `'nExt'`
% : Size of the magnetic supercell in multiples of the
% crystallographic cell, dimensions are $[1\times 3]$. Default value is
% stored in `obj`. If `nExt` is a single number, then the size of the
% extended unit cell is automatically determined from the first
% magnetic ordering wavevector. E.g. if `nExt = 0.01`, then the number
% of unit cells is determined so, that in the extended unit cell,
% the magnetic ordering wave vector is `[0 0 0]`, within the given
% 0.01 rlu tolerance.
%
% `'k'`
% : Magnetic ordering wavevector in rlu, dimensions are $[n_K\times 3]$.
% Default value is defined in `obj`.
%
% `'n'`
% : Normal vector to the spin rotation plane for single-k magnetic
% structures, stored in a 3-element row vector. Default value `[0 0 1]`. The
% coordinate system of the vector is determined by `unit`.
%
% `'S'`
% : Vector values of the spins (expectation value), dimensions are $[3\times n_{spin} n_K]$.
% Every column defines the three $(S_x, S_y, S_z)$ components of
% the spin (magnetic moment) in the $xyz$ Descartes coodinate system or
% in lu. Default value is stored in `obj`.
%
% `'unit'`
% : Determines the coordinate system for `S` and `n` vectors using a
% string:
% * `'xyz'` Cartesian coordinate system, fixed to the lattice.
% Default value.
% * `'lu'` Lattice coordinate system (not necessarily
% Cartesian. The three coordinate vectors are
% parallel to the lattice vectors but normalized to
% unity.
%
% `'epsilon'`
% : The smalles value of incommensurability that is
% tolerated without warning in lattice units. Default is $10^{-5}$.
%
% `'func'`
% : Function handle that produces the magnetic moments, ordering wave
% vector and normal vector from the `x` parameters in the
% following form:
% ```
% [M, k, n] = @(x)func(M0,x)
% ```
% where `M` is a matrix with dimensions of $[3\times n_{magExt}]$, `k` is
% the propagation vector, `n` is the normal vector of the spin rotation
% plane. The default function is [gm_spherical3d]. For planar magnetic
% structure use [gm_planar].
%
% `'x0'`
% : Input parameters for `func` function, row vector with $n_X$ number of
% elements.
%
% `'norm'`
% : Set the length of the generated magnetic moments to be equal to
% the spin of the magnetic atoms. Default is `true`.
%
% `'r0'`
% : If `true` and only a single spin direction is given, the spin
% phases are determined by atomic position times k-vector, while
% if it is `false`, the first spin will have 0 phase. Default is
% `true`.
%
% ### Output Arguments
%
% The [obj.mag_str] field will contain the new magnetic structure.
%
% ### See Also
%
% [spinw] \| [spinw.anneal] \| [spinw.optmagstr] \| [gm_spherical3d] \| [gm_planar]
%
% *[rlu]: Reciprocal Lattice Units
% *[lu]: Lattice Units
%
if isempty(obj.matom.r)
error('spinw:genmagstr:NoMagAtom','There are no magnetic atoms (S>0) in the unit cell!')
end
if isempty(obj.mag_str.k)
k0 = [0 0 0];
else
k0 = obj.mag_str.k';
end
inpForm.fname = {'mode' 'nExt' 'k' 'n' };
inpForm.defval = {'tile' obj.mag_str.nExt [] nan(1,3)};
inpForm.size = {[1 -1] [1 -4] [-6 3] [-6 3] };
inpForm.soft = {false false true false };
inpForm.fname = [inpForm.fname {'func' 'x0' 'norm' 'r0' }];
inpForm.defval = [inpForm.defval {@gm_spherical3d [] true true }];
inpForm.size = [inpForm.size {[1 1] [1 -3] [1 1] [1 1]}];
inpForm.soft = [inpForm.soft {false true false false}];
inpForm.fname = [inpForm.fname {'S' 'phi' 'phid' 'epsilon' 'unit'}];
inpForm.defval = [inpForm.defval {[] 0 0 1e-5 'xyz'}];
inpForm.size = [inpForm.size {[3 -7 -6] [1 1] [1 1] [1 1] [1 -5]}];
inpForm.soft = [inpForm.soft {true true false false false }];
param = sw_readparam(inpForm, varargin{:});
% Error if S or k is provided but is empty. This means it has failed the
% input validation, but hasn't caused an error because soft=True
err_str = [];
if any(strcmp(varargin, 'S')) && isempty(param.S)
err_str = ["S"];
end
if any(strcmp(varargin, 'k')) && isempty(param.k)
err_str = [err_str "k"];
end
if length(err_str) > 0
error('spinw:genmagstr:WrongInput', 'Incorrect input size for ' + ...
join(err_str, ', '));
end
if strcmpi(param.mode, 'rotate') && isempty(obj.mag_str.F)
error('spinw:genmagstr:WrongInput', ['rotate mode requires a magnetic ' ...
'structure to be defined with another mode first'])
end
if isempty(param.k)
noK = true;
param.k = k0;
else
noK = false;
end
if prod(double(param.nExt)) == 0
error('spinw:genmagstr:WrongInput','''nExt'' has to be larger than 0!');
end
if strcmp(param.mode,'extend')
warning('spinw:genmagstr:DeprecationWarning',...
'extend mode is deprecated, please use tile mode instead');
param.mode = 'tile';
end
% input type for S, check whether it is complex type
cmplxS = ~isreal(param.S);
if strcmpi(param.mode, 'helical') && cmplxS
error('spinw:genmagstr:WrongInput', ...
['S must be real for helical mode. To specify complex basis ' ...
'vectors directly use fourier mode.'])
end
switch lower(param.unit)
case 'lu'
% convert the moments from lattice units to xyz
BV = obj.basisvector(true);
%param.S = BV*param.S;
if ~isempty(param.S)
param.S = mmat(BV,param.S);
end
param.n = (BV*param.n')';
case 'xyz'
% do nothing
otherwise
error('spinw:genmagstr:WrongInput','Option ''unit'' has to be either ''xyz'' or ''lu''!');
end
if isempty(param.S)
% use the complex Fourier components from the stored magnetic structure
param.S = obj.mag_str.F;
cmplxS = true;
end
% Magnetic ordering wavevector(s)
k = param.k;
% number of k-vectors
nK = size(k,1);
nExt = double(param.nExt);
% automatic determination of the size of the extended unit cell based on
% the given k-vectors if nExt is a single number
if numel(nExt) == 1
[~, nExtT] = rat(param.k(1,:),nExt);
if nK>1
for ii = 2:nK
[~, nExtT2] = rat(param.k(ii,:),nExt);
nExtT = lcm(nExtT,nExtT2);
end
end
nExt = nExtT;
end
mAtom = obj.matom;
nMagAtom = size(mAtom.r,2);
% number of magnetic atoms in the supercell
nMagExt = nMagAtom*prod(nExt);
% Create mAtom.Sext matrix.
mAtom = sw_extendlattice(nExt, mAtom);
% normalized axis of rotation, size (nK,3)
if isnan(param.n(1))
% default value
param.n = repmat([0 0 1],[nK 1]);
end
n = bsxfunsym(@rdivide,param.n,sqrt(sum(param.n.^2,2)));
if size(param.n,1) ~= nK
error('spinw:genmagstr:WrongInput',['The number of normal vectors has'...
' to be equal to the number of k-vectors!'])
end
% convert input into symbolic variables
if obj.symb
k = sym(k);
param.S = sym(param.S);
n = sym(n);
end
if ~cmplxS && ~strcmpi(param.mode,'fourier') && ~strcmpi(param.mode,'direct') && any(k(:))
param.S = param.S + 1i*cross(repmat(permute(n,[2 3 1]),[1 size(param.S,2) 1]),param.S);
end
switch param.mode
case 'tile'
% effectively tiles the magnetic supercell with the given magnetic
% moments if:
% -the new number of extended cells does not equal to the number of
% cells defined in obj
% -the number of spins stored in obj is not equal to the number
% of spins in the final structure
if nMagAtom ~= size(param.S,2) && nMagExt ~= size(param.S,2)
error('spinw:genmagstr:WrongInput', ['Incorrect input size for S, ' ...
'S must be provided for each magnetic atom']);
elseif any(double(obj.mag_str.nExt) - double(param.nExt)) || (size(param.S,2) ~= nMagExt)
S = repmat(param.S,[1 prod(nExt) 1]);
else
S = param.S;
end
% sum up all kvectors and keep the real part only
S = real(sum(S,3));
k = [0 0 0];
case 'random'
% Create random spin directions and use a single k-vector
S = randn(nMagExt,3);
S = bsxfun(@rdivide,S,sqrt(sum(S.^2,2)));
S = bsxfunsym(@times,S,mAtom.Sext')';
if noK
k = [0 0 0];
end
% keep the normal vector
S = S + 1i*cross(repmat(permute(n,[2 3 1]),[1 size(S,2) 1]),S);
case {'helical' 'fourier'}
S0 = param.S;
% Magnetic ordering wavevector in the extended unit cell.
kExt = k.*nExt;
% Warns about the non sufficient extension of the unit cell.
% we substitute random values for symbolic km
skExt = sw_sub1(kExt,'rand');
if any(abs(skExt(:)-round(skExt(:)))>param.epsilon) && prod(nExt) > 1
warning('spinw:genmagstr:UCExtNonSuff','In the extended unit cell k is still larger than epsilon!');
end
% number of spins in the input
nSpin = size(param.S,2);
if (nSpin~= nMagAtom) && (nSpin==1)
% there is only a single given spin, use the fractional atomic position
if param.r0
r = mAtom.RRext;
else
r = bsxfun(@minus,mAtom.RRext,mAtom.RRext(:,1));
end
% add phase
addPhase = true;
elseif nSpin == nMagAtom
% moments in the crystallographic unit cell are defined, use
% only unit cell position.
r = bsxfun(@rdivide,floor(bsxfun(@times,mAtom.RRext,nExt')),nExt');
% add phase
addPhase = true;
elseif nSpin == nMagExt
% no phase, magnetic structure is already in superstructure
addPhase = false;
else
error('spinw:genmagstr:WrongNumberSpin','Wrong number of input spins!');
end
if addPhase
% additional phase for each spin in the magnetic supercell
%phi = sum(bsxfun(@times,2*pi*kExt',r),1);
phi = sum(bsxfunsym(@times,2*pi*permute(kExt,[2 3 1]),r),1);
% add the extra phase for each spin in the unit cell
% TODO check
S = bsxfunsym(@times,S0(:,mod(0:(nMagExt-1),nSpin)+1,:),exp(-1i*phi));
else
S = S0;
end
case 'direct'
% direct input of real magnetic moments
S = param.S;
if size(S,2) == nMagAtom
% single unit cell
nExt = [1 1 1];
end
if size(S,2) ~= nMagExt
error('spinw:genmagstr:WrongSpinSize','Wrong dimensions of param.S matrix!');
end
case 'rotate'
if param.phi == 0
% use degrees for phi if given
param.phi = param.phid*pi/180;
end
S = param.S;
if ~isreal(param.phi)
error('spinw:genmagstr:ComplexPhi', ...
['Phi should be real! If you really intended phi to have ' ...
'a complex component, this undocumented feature has been ' ...
'removed. Please contact the spinw developers.'])
end
if param.phi == 0
% The starting vector, size (1,3):
incomm = mod(bsxfun(@times,k,nExt),1);
incomm = any(incomm(:));
if incomm
S1 = sw_nvect(S);
else
S1 = sw_nvect(real(S));
end
% Axis of rotation defined by the spin direction
nRot = cross(n,S1);
if all(nRot == 0)
error('spinw:genmagstr:InvalidRotation', ...
['Cannot automatically determine rotation as the ' ...
'first spin is parallel to n']);
end
% Angle of rotation.
phi = -atan2(norm(cross(S1,n)),dot(S1,n));
else
nRot = n;
% Angle of rotation.
phi = param.phi(1);
end
% Rotate the spins.
S = sw_rot(nRot,phi,S);
k = obj.mag_str.k';
case 'func'
S = mAtom.S;
S = repmat(S,[prod(nExt) 1]);
if obj.symbolic
[S, k, n] = param.func(sym(S), sym(param.x0));
else
[S, k, n] = param.func(S,param.x0);
end
if any(k)
S = S + 1i*cross(repmat(permute(n,[2 3 1]),[1 size(S,2) 1]),S);
end
otherwise
error('spinw:genmagstr:WrongMode','Wrong param.mode value!');
end
% normalize the magnetic moments
if param.norm
normS = sqrt(sum(real(S).^2,1))./repmat(mAtom.S,[1 prod(nExt)]);
normS(normS==0) = 1;
S = bsxfunsym(@rdivide,S,normS);
end
% simplify expressions
if obj.symbolic
S = simplify(sym(S));
k = simplify(sym(k));
end
mag_str.nExt = int32(nExt(:))';
mag_str.k = k';
mag_str.F = S;
obj.mag_str = mag_str;
spinw.validate(obj);
end