迷途远山林，迟暮堪自讼。----高启《姑苏杂咏·天平山》

---

前言

本文首先简单介绍一下微磁数据后处理软件 semargl 的功能界面，并结合 mumax3 官网上的后处理程序生成的数据，演示如何使用 semargl 获取频谱关系，特定频率模式的空间分布以及自旋波的色散关系。最后，本文还对前文笔记05中第三节的内容做了两个小小的补充：①优化多进程读取磁化数据文件的代码。②新增获取特定频率模式的空间分布的功能。

---

一、semargl 软件使用简介

1.semargl 软件简介

semargl 是 M. Dvornik 于2013年开发的一款专用于处理微磁软件输出数据的开源软件，主要功能就是对所选择的磁化分量进行快速傅里叶变换（FFT），可用于获取频谱关系，特定频率模式的空间分布，自旋波的色散关系等。虽然该软件的年龄有些老，而且有些数据保存的功能也有缺陷，但是它读取数据文件的速度较快，并且图像化操作也是比较方便的，所以可以作为一个对比分析工具。
semargl 的源码见GitHub 地址，编译过的 Windows 端程序见下载地址，其中 SSE 和 AVX 两个版本的程序经过测试并没有看出啥区别，任选一个使用就行了。注意，semargl 软件在每次计算前都需要重新启动以清除缓存数据，并且对输入的磁化文件有些要求：
①文件路径中不能有中文，并且磁化文件的后缀名必须是“.ovf”。（由于OOMMF输出的磁化文件的后缀名是.omf，所以需要修改一下）
②磁化文件中保存数据的方式为b4或者b8，不能是text文本格式。
③对于mumax3输出的磁化文件来说，在按时间顺序保存的最后一个磁化文件中需要用户手动修改一下：在表示模拟总时间的标签行 “# Desc: Total simulation time: **** s” 的前面补充一行表示每个磁化文件的保存时间间隔的标签行“# Desc: Stage simulation time: **** s”。至于OOMMF输出的磁化文件则是默认有这一标签行的，不用手动添加。

1.“File and info”标签页面
直接双击 .exe 文件打开 semargl 软件，默认进入“File and info”标签页面，如下所示：

该页面主要用于用户指定磁化文件的存放路径，待分析的磁化分量。在该页面中：
默认勾选框“Extract the path to the data from the static file”，接着直接点击“Static file”按钮选择基态磁化文件 $m(r,t_0)$ ，也就是按时间顺序保存的第一个磁化文件“xxx000.ovf”，选中之后，该页面根据磁化文件中的标签行相关信息来自动刷新，接着用户可以检查一下“Simulations details”信息框中的内容是否无误。

在选项框“Component of the magnetization to be analysed”中选择待分析的磁化分量Mx/y/z，至于后面两项Minp（计算平面分量）和div(M)/4Pi则不常用。

在选项框“Pre-processing method”中选择不减去（Bypass）基态磁化状态或者减去（Substrat static magnetization），即 $m(r,t)-m(r,t_0)$，通常都需要减去之后得到动态磁化状态用于后续FFT计算。

选项框“Normalize by magnetization of saturation”可以选择是否归一化待分析的磁化分量。

设置完以上内容之后，点击“GO”按钮开始读取磁化分量数据，等待页面左下角进度条停止即可。

2.“Settings”标签页面
读取完磁化文件之后，点击进入“Settings”标签页面，如下所示：

在该页面可以设置FFT所需的窗函数，至于剩下的选项框则用于设置输出文本数据和图像数据的格式，通常保持默认即可。

3.“Probe Configuration”标签页面
点击进入“Probe Configuration”标签页面，在该页面中可以设置空间采样范围。如下所示：

右侧的“Skin depth emulator”和“Optical spot emulator”选项框是用于模拟真实的实验条件，一般情况下保持默认即可。

在选项框“Probe type”中可以选择“NONE”方式表示计算的是磁化文件中默认的整个长方体形状的空间范围。“Rectangular ”方式则可以在下方的6个文本框中手动输入自定义的长方体空间范围，用户手动指定空间范围之后需要点击左边的“GENERATE”按钮来确定。在下方的视图中还会显示用户所指定的空间范围。

4.“Time domain and propagation analysies”标签页面
点击进入“Time domain and propagation analysies”标签页面，如下所示：

该页面可以显示和保存自旋波时域和频域的传播图像。

5.“mmDisplay”标签页面
点击进入“mmDisplay”标签页面，如下所示：

该页面主要用于显示和保存时间空间域$m(r,t)$和频率空间域$\tilde{m}(r,w)$图像。

6.“Frequnency domain analysies”标签页面
点击进入“Frequnency domain analysies”标签页面，如下所示：

在该页面中可设置空间区域范围内的单元格FFT的3种方法，具体参考文章《Micromagnetic Simulations in Magnonics》中8.3公式，或者前文笔记05中相关描述。在“Method”选项框中：

第一项“Temporal FFT -> Averaging (100% spatial BW) ”对应公式8.3(b)，即首先对所选空间范围中所有单元格进行一维FFT，得到每一个单元格的频谱图之后，再取平均值转化为一个频谱图。

第二项“Averaging -> Temporal FFT (excludes AS modes) ”对应公式8.3(a)，即首先把所选空间范围中所有单元格在时域上取平均，转化为一个单元格的数据，之后再对这一个单元格的数据进行一维FFT，得到一个频谱图。

第三项“Full DFT > Averaging (50% spatial BW)”则对应公式8.3©，是用于计算二维FFT的选项，该方法专门用于计算自旋波的色散关系。

至于第四项“Custom DFT -> Averaging”则不清楚有啥用。

选项框“Averaging Method”中的选项则会自动变化，保持默认即可。

选项框“Mode”中可以设置一维FFT结果数据的保存方式，一般都是“Amplitude and Phase”选项。

设置完以上内容之后，除了二维FFT即计算色散关系之外，点击“CALCULATE”按钮开始一维FFT计算，等待页面左下角进度条停止即可。当一维FFT计算完之后，点击右侧的“SAVE”按钮即可保存频谱关系（输出数据有三列：第一列为频率值，第二列为幅度，第三列为相位）！若一维FFT方法选中第一项“Temporal FFT -> Averaging (100% spatial BW) ”的话，计算完之后，可以点击进入“mmDisplay”页面查看特定频率模式的空间分布（即频率空间域图像）！

7.“Dispersion analysies”标签页面

点击进入“Dispersion analysies”标签页面，如下所示：

在该页面中可设置计算自旋波色散关系时相关选项，注意：计算色散关系之前，需先在“Frequnency domain analysies”页面中选中第三项“Full DFT > Averaging (50% spatial BW)”，但不点击“CALCULATE”按钮。
右侧的“Probe”选项框用于模拟真实实验条件，一般保持默认即可。

在“Averaging”和“Averaging method”选项框中设置波矢方向，即若波矢是X方向，则Y和Z方向的多余单元格的数据应该平均为一个单元格：“REAL SPACE”表示在FFT之前平均，“COMPLEX AMPLITUDE”则表示FFT之后平均。

设置完以上选项后，点击“GET SLICE”按钮开始计算二维FFT，等待页面左下角进度条停止，点击“VIEW”按钮即可刷新页面观察色散曲线。但“SAVE TXT”保存数据文件功能不可用，只能“SAVE IMAGE”保存图片文件。 图片中一个像素点对应二维FFT之后的矩阵中的一个数值，横轴范围为$[-\frac{1}{2dl},\frac{1}{2dl}]$，纵轴范围为$[-\frac{1}{2dt},\frac{1}{2dt}]$，其中dt和dl分别是采样时间和波矢方向的单元格尺寸。

2.准备演示软件操作所需的数据

为了方便演示 semargl 的使用方法，这里使用 mumax3 官网上的微磁数据后处理程序作为对比示例程序，只需关注其中第4和第5小节中的内容，即获取斯格明子的激发谱，两个频率模式的空间分布，和自旋波的色散关系。
在电脑上运行完这个程序之后，会看到生成了3个文件夹，其中的：skyrmion.out（共407个文件，其中磁化文件.ovf有201个，保存时间间隔为1e-11 s） 和 spinwaves.out （共805个文件，其中磁化文件.ovf有401个，保存时间间隔为2.5e-11 s）文件夹中保存了本文演示所需要的数据，如下所示：

对于mumax3输出的磁化文件来说，在按时间顺序保存的最后一个磁化文件中需要用户手动修改一下：在 skyrmion.out 文件夹中，对最后一个磁化文件即“m000200.ovf”，手动在语句“# Desc: Total simulation time: 2e-09 s”的前面添加语句“# Desc: Stage simulation time: 1e-11 s”。
在 spinwaves.out 文件夹中，对最后一个磁化文件即“m000400.ovf”，手动在语句“# Desc: Total simulation time: 1e-08 s”的前面添加语句“# Desc: Stage simulation time: 2.5e-11 s”。

3.使用 semargl 获取频谱关系

示例后处理程序片段如下：

# FAST FOURIER TRANSFORM
dm     = table["mz"] - table["mz"][0]   # average magnetization deviaton
spectr = np.abs(np.fft.fft(dm))         # the absolute value of the FFT of dm
freq   = np.linspace(0, 1/dt, len(dm))  # the frequencies for this FFT

# PLOT THE SPECTRUM
plt.plot(freq/1e9, spectr)
plt.xlim(0,fmax/1e9)
plt.ylabel("Spectrum (a.u.)")
plt.xlabel("Frequency (GHz)")
plt.show()

从而可以看出，斯格明子的激发谱是通过对整个磁体系的平均磁化z分量即$m_z$进行一维FFT获取的。那么使用 semargl 的操作步骤如下：①打开 semargl 软件，在“File and info”页面中点击“Static file”按钮，进入 skyrmion.out 文件夹，选中“m000000.ovf”，之后选择分析"Mz"，不归一化，最后点击”GO“按钮，等待右下角读条结束如下所示：

②接着点击进入“Probe Configuration”页面，这里只能选择采样整个空间范围的单元格来近似整个磁体系的平均磁化z分量即$m_z$，所以保持默认即可。

③接着点击进入“Frequnency domain analysies”页面，由于是先在时域平均再FFT，所以在“Method”中选择第二项“Averaging -> Temporal FFT (excludes AS modes) ”，最后点击“CALCULATE”按钮，等待右下角读条结束，在“Mode”中选择“Amplitude and Phase”，点击“SAVA”按钮保存结果文件.dat，使用 Origin 软件绘制频谱曲线，和示例程序对比如下：

虽然二者曲线形状类似，但是可看出 semargl 的幅度整体下降，推测可能原因是由于示例程序的磁体系是一个圆盘，而 semargl 在时域平均整个长方体空间范围的单元格的过程中，引入了圆盘形状之外的多余的无磁性单元格参与计算，导致结果偏低。

4.使用 semargl 获取特定频率模式的空间分布

示例后处理程序片段如下：

# Select the z component and the (only) layer z=0
mz = m[:,2,0,:,:]

# Apply the FFT for every cell
mz_fft = np.fft.fft(mz, axis=0)

# Select the the two modes
mode1 = mz_fft[mode1_idx]
mode2 = mz_fft[mode2_idx]

# Plot the result
plt.figure(figsize=(10,4))
plt.subplot(1,3,1)
plt.title("$m_z$")
plt.imshow(mz[0])
plt.subplot(1,3,2)
plt.title("Mode 1")
plt.imshow(np.abs(mode1)**2)
plt.subplot(1,3,3)
plt.title("Mode 2")
plt.imshow(np.abs(mode2)**2)
plt.show()

从而可以看出，需要获取频谱图中两个峰值频率处斯格明子模式的分布。在 semargl 中操作步骤如下：①打开 semargl 软件，在“File and info”页面中点击“Static file”按钮，进入 skyrmion.out 文件夹，选中“m000000.ovf”，之后选择分析"Mz"，不归一化，最后点击”GO“按钮，等待右下角读条结束。

②接着点击进入“Probe Configuration”页面，对于获取特定频率模式的空间分布的情况，必须对整个空间范围内的所有单元格进行一维FFT，所以本页面保持默认即可。

③接着点击进入“Frequnency domain analysies”页面，在“Method”中选择第一项“Temporal FFT -> Averaging “，最后点击“CALCULATE”按钮，等待右下角读条结束，点击进入“mmDisplay”页面查看特定频率模式的空间分布，和示例程序对比如下：

相比于示例程序中一次计算只能获取几个频率模式的分布，使用 semargl 可以一次性获取并查看频谱图中所有频率模式的分布，图像化操作也十分方便。

5.使用 semargl 获取自旋波的色散关系

示例后处理程序片段如下：

# Stack all snapshots of the magnetization on top of each other
m = np.stack([fields[key] for key in sorted(fields.keys())])

# Select the x component
mx = m[:,0,0,0,:]

# Apply the two dimensional FFT
mx_fft = np.fft.fft2(mx)
mx_fft = np.fft.fftshift(mx_fft)

plt.figure(figsize=(10,6))

plt.subplot(2,1,1)

# Show the intensity plot of the 2D FFT
extent = [ -(2*np.pi)/(2*dx), (2*np.pi)/(2*dx), -1/(2*dt), 1/(2*dt)] # extent of k values and frequencies
plt.imshow(np.abs(mx_fft)**2, extent=extent, aspect='auto', origin='lower', cmap="bwr")

# Plot the analytical derived dispersion relation 
k = np.linspace(-2e8,2e8,1000)
freq_theory = A*gamma*k**2 /(np.pi*Ms) + gamma*Bz /(2*np.pi)
plt.plot(k,freq_theory,'r--',lw=1)
plt.axhline(gamma*Bz/(2*np.pi),c='g',ls='--',lw=1)

plt.xlim([-2e8,2e8])
plt.ylim([0,fmax])
plt.ylabel("$f$ (Hz)")
plt.xlabel("$k$ (1/m)")

示例程序中的磁体系是一维的纳米线，波矢为X方向，并且选择分析的是x分量$m_x$进行二维FFT得到色散关系。在 semargl 中操作步骤如下：①打开 semargl 软件，在“File and info”页面中点击“Static file”按钮，进入 spinwaves.out 文件夹，选中“m000000.ovf”，之后选择分析"Mx"，不归一化，最后点击”GO“按钮，等待右下角读条结束。

②接着点击进入“Probe Configuration”页面，空间采样范围为整个磁体系，所以本页面保持默认即可。

③点击进入“Frequnency domain analysies”页面，在“Method”中选择第三项“Full DFT > Averaging (50% spatial BW)“，但不点击“CALCULATE”按钮。

④点击进入“Dispersion analysies”页面，在“Averaging”和“Averaging method”选项框中设置Y和Z方向的多余单元格平均为一维的方法，但该示例中Y和Z方向本来就是一维的，所以“Averaging method”保持默认即可，点击“GET SLICE”按钮开始计算二维FFT，等待页面左下角进度条停止，点击“VIEW”按钮即可刷新页面观察色散曲线，和示例程序对比如下：

由于“SAVE TXT”保存数据文件功能不可用，所以只能“SAVE IMAGE”保存图片文件，十分不方便用户画图，希望有大佬可以基于源码改进一下这个软件。

二、笔记05第三节内容的补充

1.优化多进程读取磁化数据文件的代码

经过测试，发现使用之前的代码在第一次多进程读取磁化数据文件时，磁盘占用100%导致读取速度反而下降，所以需要区分是否第一次读入，从而使用多进程还是单进程读取文件：

    if isFirstRead == 1:
        m_allFiles = readUseOneCore(Dir,StartStage,EndStage,export_x,export_y,export_z)
        print('已经使用单进程读完了所有文件！',flush=True)
    else:
        m_allFiles = readUseMultiCore(Dir,StartStage,EndStage,export_x,export_y,export_z,CoreNum)
        print('已经使用多进程读完了所有文件！',flush=True)

2.新增获取特定频率模式的空间分布的代码

新建文件“getModeProfile.py”，内容如下：

# -*- coding: utf-8 -*-
"""
Created on 22/11/8
获取指定频率点自旋波的空间模式分布
@author: YQYUN
"""

import numpy as np

#对矩阵（行代表时间，列代表空间单元格）进行一维傅里叶变换
#Mx,My,Mz：三者之一保存着待分析的矩阵
#export_x,export_y,export_z：待分析的磁化分量
#T：每个磁化文件的间隔时间（采样时间）
#frequencyIndexList：需要获取的频率点索引列表
def getSWModeProfile(Mx,My,Mz,export_x,export_y,export_z,T,frequencyIndexList):
    if export_z == 1:
        M_matrix = Mz
    elif export_y == 1:
        M_matrix = My
    else:
        M_matrix = Mx
    
    #得到M_matrix的行(M)和列(N)数
    #行是采样时间间隔:t0,t1,t2...tM-1，列是指定区域单元格（空间排列顺序为x->y->z）数量cell0,cell1...cellN-1
    M,N = M_matrix.shape
    #采样频率
    Fs = 1 / T
    #N个实数点进行FFT之后也有N个复数点，这里需要对每个单元格进行FFT，所以结果点数为N*M
    fftMatrix = np.zeros((N,M))
    #在时间（M行）-单元格（N列）矩阵中，按照列的方向来求每个单元格的FFT
    for i in range(0,N):
        #将FFT得到的M*N个数据点按照N行（表示单元格），M列（表示时间）来存放在fftMatrix中
        fftMatrix[i,:] = np.abs(np.fft.fft(M_matrix[:,i]))
        #对FFT后的复数结果的模值进行处理，得到正确的幅值:
        #除了第一个点（0Hz）外，对模值要先除以采样点数量，再乘以2
        fftMatrix[i,:] = fftMatrix[i,:] / M * 2
        #对于第一个点（0Hz），只需对模值除以采样点数量即可
        fftMatrix[i,0] = fftMatrix[i,0] / 2
    
    #只取前一半部分正的频率点，舍弃后一半部分负的频率点
    fftMatrix = fftMatrix[0:N,range(int(M / 2))]
    
    #归一化幅度
    #fftMatrix = fftMatrix / np.max(fftMatrix)
    
    #将结果保存到文件
    modeProfile = np.zeros((N,1))
    for frequencyIndex in frequencyIndexList:
        modeProfile[0:N,0] = fftMatrix[0:N,frequencyIndex]
        #根据每个频率点获取对应的频率值（GHz），频率分辨率为采样频率/采样点数
        frequency = Fs / M * (frequencyIndex-1) / 1e9
        #列之间用空格分隔
        fileName = f"f= {frequency} GHz,modeProfile.txt"
        np.savetxt(fileName,modeProfile,delimiter=' ')
        print(f"第 {frequencyIndex} 个频率点的模式分布已经写入 {fileName}")
    
    return;

在主程序中新加入调用上述函数的代码：

    if ModeProfile_Switch == 1:
        from getModeProfile import getSWModeProfile
        getSWModeProfile(Mx,My,Mz,export_x,export_y,export_z,Ts,frequencyIndexList)
    else:
        ModeProfile_Switch = ModeProfile_Switch 

改进后的整个程序见下载地址。

---

总结

关于微磁数据后处理的程序的介绍差不多已经结束了，觉得对您有帮助的可以帮忙点一个赞。