Fe3O4 Nanosphere LLG Simulation Accuracy Help

[MichaelSherburne]

Hello,

I am fairly close in being a bit more certain that I can get Ubermag to replicate FMR values in this journal article regarding the LLG FMR of a Fe3O4 nanosphere: I Nature Report. I am starting this as a thread to better document changes. I followed their supplementary article section S4 to the tee for the simulation inputs.

I was able to get it to be fairly consistent with the paper's results and to first principles if I multiply gamma0 by 2pi. I believe it may be due to the conversion from Hz/Oe -> m/As. For a while, I was doing 2.8MHz/Oe / 79.577 to convert over to m/As, but it appears from trial and error that I need to multiply this by 2pi to work consistently across various Hdc settings. I tried changing other values but none of them would relatively accurately match for other Hdc values after calibrating the 1000Oe Hdc case. Would multiplying gamma0 by 2pi in my case be physically valid or is it a weird coincidence that turns out to be a correction factor?

I also tried multiplying the frequency of the x-axis of the FMR output plot by 2pi instead of the 2pi correction for gamma0, this also appears to work, however looking in Ubermag's output, it appears the FFT function is in Hz, so multiplying this by 2pi does not make sense to get a frequency output or maybe it is okay?

At Hdc = 1000Oe, one should get ~3GHz FMR. At Hdc = 2000Oe, one should get ~6GHz. Going towards 500Oe, the error between first principles and Ubermag becomes higher, however there is a note in the journal article that mentioned this: "For Hdc > 500Oe, the Zeeman energy sufficiently overcomes the thermal energy of spin relaxation in the particles at room temperature. Therefore, our numerical simulation carried out at T = 0K can represent the general feature of the precession in a single spherical Fe3O4 nanoparticle..." Ubermag is doing its simulations at 0K, therefore I am thinking higher error here compared to a simplified first principles model could be due to effects that occur at <500Oe Hdc. Therefore, I could trust Ubermag's output?

There is going to be some error in general since this gets fairly complicated to make fully accurate due to orbital momentums of Fe3O4, whether you have more or less Fe2+ and Fe3+ ions. This can change the Lande factor and thus the gamma0 can go from 2.8MHz/Oe to another value. Getting in the approximate ballpark is what counts here.

I have a feeling a correction factor of 2pi is needed to gamma0s given in literature as MHz/Oe along with a division of 79.577 or it is indeed multiplying the x-axis by 2pi to get the proper frequency in Hz in the FMR plot? Please let me know which one makes sense or if something else is going on? There is certainly a 2pi correction that needs to occur somewhere and I want to make sure it is in the correct spot.

The code I wrote is below:

#Import Settings
import discretisedfield as df
import micromagneticmodel as mm
import oommfc as mc

#For Post-Process Plotting
import matplotlib.pyplot as plt
import numpy as np

#Set Simulation Parameters
UnitCell = 3.75e-9 #2.375e-9 # Unit Cell Size (m) 

#Set Material Parameters for Fe3O4 Nanoparticles
# Sphere 7.5 nm radius
r = 7.5e-9
MsUser = 102.24 #emu/g lowering this will make FMR greatly higher, but barely change with Hdc
K1 = 1.36e4/2 #(1.36e4/1) #J/m^3 (Lee 2021 paper has it divided by 2 for an average)
Aex = 13.2e-12 #J/m
Hint = 0 #Internal Magnetic Field of Fe3O4 Nanoparticle (Oe)
#Lee's paper used num modeling to get Hint, appears they did not use it in input
Hdc = 1000 #Applied DC Field (Oe)
alphaUser = 0.246
gammaUser = (2.8E6/79.577471546)*(2*np.pi) #5 #*2*np.pi #Hz/Oe -> m/As
Hac = (5*79.577471546) #HfieldAC <-This does not affect freq peak loc unless 0
gammaUser #check to make sure the number makes sense (2.21e5 with these settings)

#Process Variables
Ms = MsUser*5.18*1000 #emu/g * g/cm^3 -> 1 emu/cm^3 -> 1000 A/m <- major effect on peak
Hk = (2*K1)/Ms
Hstatic = (Hint+Hdc+Hk)*79.5774715459 #A/m
print(Ms)
print(Hk)

UserCell =  (UnitCell,UnitCell,UnitCell)  # mesh discretisation (m)

def norm_fun(point):
    x, y, z = point

    if x**2 + y**2 + z**2 <= r**2:
        return Ms
    else:
        return 0

#Create Mesh
mesh = df.Mesh(p1=(-r, -r, -r), p2=(r, r, r), cell=UserCell)
m_init = df.Field(mesh=mesh, dim=3, value=(0, 0, 1), norm=norm_fun)

m_init.plane('x').mpl()
m_init.norm.k3d.nonzero()
system = mm.System(name='nanosphere',T=0) #Specify finite temperature here
#Default is 0K

#J/m for Aex  #mm.Exchange(A=Aex)
#Ground State Energy Setup
system.energy = (mm.Exchange(A=Aex) +
                 mm.Demag() +
                 mm.CubicAnisotropy(K=K1, u1=(1, 0, 0), u2=(0, 1, 0)) +
                 mm.Zeeman(H=(0,0,Hstatic)))  # put DC field here
system.energy

system.m = m_init
system.dynamics = mm.Precession(gamma0=gammaUser) + mm.Damping(alpha=alphaUser)
system.dynamics
md = mc.MinDriver()
md.drive(system)
system.m.plane('y').mpl()
system.energy += mm.Zeeman(H=(Hac,0,0), wave='sinc', f=10e9, t0= 10e-9, name='ac')
td = mc.TimeDriver()
td.drive(system, t=30e-9, n=500,n_threads = 20)
system.table.data.plot('t', 'mx')
fft = system.table.rfft()
susceptibility =  (fft.data['ft_mx']*Ms)/((fft.data['ft_Bx_ac']/1e3)/(mm.consts.mu0))  
from scipy.signal import find_peaks
lst = -1*np.imag(susceptibility)
peaks, _ = find_peaks(lst, height=5) #change height to get the exact peak
print(peaks)      
freqpeak = np.array(fft.data['f']/1e9)
print(freqpeak[peaks]*1e9) #Frequency Peak of X''
#Create Susceptibility Plot (Frequency is 100% in Hz from FFT Function

fig, ax1 = plt.subplots() 

ax1.set_xlabel('Frequency (GHz)') 
ax1.set_ylabel('Real (X)', color = 'red') 
ax1.plot((fft.data['f']/1e9),np.real(susceptibility), color = 'red') 
#x1.plot(((fft.data['f']/1e9)*2*np.pi),np.real(susceptibility), color = 'red') 
ax1.tick_params(axis ='y', labelcolor = 'red') 
ax2 = ax1.twinx() 

ax2.set_ylabel('Imaginary (X)', color = 'blue') 
ax2.plot((fft.data['f']/1e9),-1*np.imag(susceptibility), color = 'blue') 
#x2.plot(((fft.data['f']/1e9)*2*np.pi),-1*np.imag(susceptibility), color = 'blue') 
ax2.tick_params(axis ='y', labelcolor = 'blue') 
plt.xlim((0,15))
plt.show()

#FFT Power
#plt.plot((fft.data['f']/1e9),np.abs(susceptibility)**2)
plt.plot(((fft.data['f']/1e9)),np.abs(susceptibility)**2)
plt.title('FFT Power')
plt.xlabel('Frequency (GHz)')
plt.ylabel('FFT Power (a.u.)')
plt.xlim((0,10))
lst = np.abs(susceptibility)**2
peaks, _ = find_peaks(lst, height=40)
freqpeak = np.array((fft.data['f']/1e9))
print(freqpeak[peaks]*1e9) #Frequency Peak of X''

#Calculate the SAR
SAR = ((mm.consts.mu0*2*np.pi*fft.data['f'])/2)*(Hac)*-1*np.imag(susceptibility)

#plt.plot((fft.data['f']/1e9), SAR)
plt.plot(((fft.data['f']/1e9)), SAR)
plt.title('SAR of Fe3O4 Nanoparticle')
plt.xlabel('Frequency (GHz)')
plt.ylabel('SAR (W/g)')
plt.xlim((0,10))

[MichaelSherburne]

@marijanbeg @lang-m after a few discussions, is the fftr function outputting everything in Hz or potentially needs a 2pi correction factor to show in Hz and what definition of Fourier is it using? We are thinking this is the most likely place for the error to occur.

[marijanbeg]

Hi @MichaelSherburne, is this issue resolved?

[MichaelSherburne]

It is!