Hard magnet tutorial#

Introduction#

  • In this notebook we explore hard magnet properties such as Hc as function of temperature for Fe16N2.

  • We query databases to get temperature-dependent inputs for micromagnetic simulations from DFT and spin dynamics simulations

  • We run hysteresis simulations and compute derived quantities.

Requirements:

The MODA diagram is provided at the bottom of the notebook.

[1]:
%config InlineBackend.figure_format = "retina"

import math

import mammos_analysis
import mammos_dft
import mammos_entity as me
import mammos_mumag
import mammos_spindynamics
import mammos_units as u
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from matplotlib import colormaps
[2]:
# Allow convenient conversions between A/m and T
u.set_enabled_equivalencies(u.magnetic_flux_field());

DFT data: magnetization and anisotropy at zero Kelvin#

The first step loads spontaneous magnetization Ms_0 and the uniaxial anisotropy constant K1_0 from a database of DFT calculations (at T=0K).

We can use the print_info flag to trigger printing of crystallographic information.

[3]:
results_dft = mammos_dft.db.get_micromagnetic_properties("Fe16N2", print_info=True)
Found material in database.
Chemical Formula: Fe16N2
Space group name: P42/mnm
Space group number: 136
Cell length a: 8.78 Angstrom
Cell length b: 8.78 Angstrom
Cell length c: 12.12 Angstrom
Cell angle alpha: 90.0 deg
Cell angle beta: 90.0 deg
Cell angle gamma: 90.0 deg
Cell volume: 933.42 Angstrom3
ICSD_label:
OQMD_label:

[4]:
results_dft
[4]:
MicromagneticProperties(Ms_0=SpontaneousMagnetization(value=1671126.902464901, unit=A / m), K1_0=UniaxialAnisotropyConstant(value=1100000.0, unit=J / m3))
[5]:
results_dft.Ms_0
[5]:
SpontaneousMagnetization(value=1671126.902464901, unit=A / m)
[6]:
results_dft.K1_0
[6]:
UniaxialAnisotropyConstant(value=1100000.0, unit=J / m3)

Temperature-dependent magnetization data from spindynamics database lookup#

In the second step we use a spin dynamics calculation database to load data for the temperature-dependent magnetization.

[7]:
results_spindynamics = mammos_spindynamics.db.get_spontaneous_magnetization("Fe16N2")

We can visualize the pre-computed data using .plot.

[8]:
results_spindynamics.plot();
../../_images/examples_workflows_hard-magnet-tutorial_11_0.png

We can access T and Ms and get mammos_entity.Entity objects:

[9]:
results_spindynamics.T
[9]:
ThermodynamicTemperature(value=
[  10.   20.   40.   60.   80.  100.  120.  140.  160.  180.  200.  220.
  240.  250.  260.  270.  290.  310.  330.  350.  370.  390.  410.  430.
  450.  470.  490.  510.  530.  550.  570.  590.  610.  630.  650.  670.
  690.  710.  730.  750.  770.  790.  810.  830.  850.  870.  890.  910.
  930.  950.  970.  990. 1010. 1030. 1050. 1070. 1090. 1110. 1130. 1150.
 1170. 1190. 1210. 1230. 1250. 1270. 1290. 1310. 1330. 1350. 1370. 1390.
 1410. 1430. 1450. 1470. 1490. 1510. 1530. 1550. 1570. 1600.],
 unit=K)
[10]:
results_spindynamics.Ms
[10]:
SpontaneousMagnetization(value=
[1639272.04060506 1639375.49131807 1629109.99748864 1618804.71492344
 1608308.44642653 1597692.8117223  1587021.47278799 1576206.89440489
 1565368.44278034 1554346.96297122 1543229.99019625 1531906.11599527
 1520542.4530585  1514956.11455598 1509250.40984613 1503385.5501932
 1491950.26753205 1479981.81581153 1467981.53310241 1456148.36308352
 1444132.16488009 1431606.67085875 1418874.2754114  1405958.85177949
 1393528.85072402 1380183.70874576 1366997.7217106  1353159.19940876
 1339583.28276302 1325235.46464328 1311126.37893819 1296293.13824202
 1281237.08062553 1265600.10746675 1250249.61320554 1234525.10482806
 1217774.04706764 1200091.93289013 1183277.21315247 1165443.90177902
 1146377.13959661 1128058.40564674 1108681.2913253  1086869.10637455
 1065781.07641488 1044048.46893568 1020732.26977272  995935.929639
  969261.5611768   942356.41804711  911663.38727184  878965.00421361
  844937.67738056  808204.71651495  766450.41719479  720597.87809002
  664272.9437298   601223.71302394  523094.55146013  407754.96420144
  171823.67656201   49003.80697799   30979.50967284   24024.43865972
   20817.46655642   17053.4521523    15780.21260756   14578.59278722
   13249.6490124    12533.45176849   11634.22634002   11244.29672944
   10520.14173837   10026.76141479    9803.94449446    9469.71911397
    9358.3106538     8809.22610014    8554.57819119    8284.01478793
    8164.64858061    7933.87391313],
 unit=A / m)

We can get also the data in the form of a pandas.DataFrame, which only contains the values in SI units:

[11]:
results_spindynamics.dataframe.head()
[11]:
T Ms
0 10.0 1.639272e+06
1 20.0 1.639375e+06
2 40.0 1.629110e+06
3 60.0 1.618805e+06
4 80.0 1.608308e+06

Calculate micromagnetic intrinsic properties using Kuz’min formula#

  • We use Kuz’min equations to compute Ms(T), A(T), K1(T)

  • Kuz’min, M.D., Skokov, K.P., Diop, L.B. et al. Exchange stiffness of ferromagnets. Eur. Phys. J. Plus 135, 301 (2020). https://doi.org/10.1140/epjp/s13360-020-00294-y

  • Additional details about inputs and outputs are available in the API reference

[12]:
results_kuzmin = mammos_analysis.kuzmin_properties(
    T=results_spindynamics.T,
    Ms=results_spindynamics.Ms,
    K1_0=results_dft.K1_0,
)

The plot method of the returned object can be used to visualize temperature-dependence of all three quantities. The temperature range matches that of the fit data:

[13]:
results_kuzmin.plot();
../../_images/examples_workflows_hard-magnet-tutorial_20_0.png
[14]:
results_kuzmin
[14]:
KuzminResult(Ms=Ms(T), A=A(T), Tc=CurieTemperature(value=1172.0457204795637, unit=K), s=<Quantity 1.68813857>, K1=K1(T))
  • The attributes Ms, A and K1 provide fit results as function of temperature. They each have a plot method.

  • Tc is the fitted Curie temperature.

  • s is a fit parameter in the Kuzmin equation.

To visually assess the accuracy of the fit, we can combine the plot methods of results_kuzmin.Ms and results_spindynamics:

[15]:
ax = results_kuzmin.Ms.plot(label="Kuzmin fit")
results_spindynamics.plot(ax=ax, label="Spin dynamics");
../../_images/examples_workflows_hard-magnet-tutorial_23_0.png

To get inputs for the micromagnetic simulation at a specific temperature we call the three attributes Ms, A and K1. We can pass a mammos_entity.Entity, an astropy.units.Quantity or a number:

[16]:
temperature = me.T(300)
temperature
[16]:
ThermodynamicTemperature(value=300.0, unit=K)
[17]:
results_kuzmin.Ms(temperature)  # Evaluation with Entity
[17]:
SpontaneousMagnetization(value=1524422.5569099565, unit=A / m)
[18]:
results_kuzmin.A(300 * u.K)  # Evaluation with Quantity
[18]:
ExchangeStiffnessConstant(value=6.90010473546342e-12, unit=J / m)
[19]:
results_kuzmin.K1(300)  # Evaluation with number
[19]:
UniaxialAnisotropyConstant(value=884617.8319031609, unit=J / m3)

Run micromagnetic simulation to compute hysteresis loop#

  • We now compute a hysteresis loop (using a finite-element micromagnetic simulation) with the material parameters we have obtained.

  • We simulate a 20x20x20 nm cube for which a pre-defined mesh is available.

  • Additional documentation of this step is available this notebook.

[20]:
results_hysteresis = mammos_mumag.hysteresis.run(
    mesh_filepath=mammos_mumag.mesh.CUBE_20_nm,
    Ms=results_kuzmin.Ms(temperature),
    A=results_kuzmin.A(temperature),
    K1=results_kuzmin.K1(temperature),
    hstart=(1.5 * u.T).to(u.A / u.m),
    hfinal=(-1.5 * u.T).to(u.A / u.m),
    hnsteps=30,
)

The returned results_hysteresis object provides a plot method to visualize the computed data. mammos_mumag.hysteresis only computes half a hysteresis loop, going from hstart to hfinal. To show a full loop this function mirrors the computed data and plots it twice:

[21]:
results_hysteresis.plot(marker=".");  # blue: simulation output, orange: mirrored data
../../_images/examples_workflows_hard-magnet-tutorial_32_0.png

The result object provides access to H and M:

[22]:
results_hysteresis.H
[22]:
ExternalMagneticField(value=
[ 1.19366207e+06  1.11408460e+06  1.03450713e+06  9.54929658e+05
  8.75352187e+05  7.95774715e+05  7.16197244e+05  6.36619772e+05
  5.57042301e+05  4.77464829e+05  3.97887358e+05  3.18309886e+05
  2.38732415e+05  1.59154943e+05  7.95774715e+04 -1.54610297e-10
 -7.95774715e+04 -1.59154943e+05 -2.38732415e+05 -3.18309886e+05
 -3.97887358e+05 -4.77464829e+05 -5.57042301e+05 -6.36619772e+05
 -7.16197244e+05],
 unit=A / m)
[23]:
results_hysteresis.M
[23]:
SpontaneousMagnetization(value=
[ 1520754.32089275  1520522.44169762  1520267.63344469  1519986.74508936
  1519676.05406227  1519331.15952891  1518946.78905411  1518516.61539684
  1518032.97097948  1517486.49476462  1516865.63479532  1516156.06490823
  1515339.72008483  1514393.6087704   1513288.01242433  1511984.00604483
  1510429.5214326   1508553.485197    1506255.68170452  1503391.93986854
  1499742.04803147  1494950.56315574  1488379.26535328  1478620.2951192
 -1518946.80143828],
 unit=A / m)

The dataframe property generates a dataframe in the SI units.

[24]:
results_hysteresis.dataframe.head()
[24]:
configuration_type H M energy_density
0 1 1.193662e+06 1.520754e+06 -2.725976e+06
1 1 1.114085e+06 1.520522e+06 -2.573912e+06
2 1 1.034507e+06 1.520268e+06 -2.421872e+06
3 1 9.549297e+05 1.519987e+06 -2.269859e+06
4 1 8.753522e+05 1.519676e+06 -2.117876e+06

We can generate a table in alternate units:

[25]:
df = pd.DataFrame(
    {
        "mu0_H": results_hysteresis.H.q.to(u.T),
        "J": results_hysteresis.M.q.to(u.T),
    },
)
df.head()
[25]:
mu0_H J
0 1.5 1.911036
1 1.4 1.910745
2 1.3 1.910425
3 1.2 1.910072
4 1.1 1.909681

Plotting of magnetization configurations#

Simulation stores specific magnetization field configurations:

[26]:
results_hysteresis.plot(configuration_marks=True);
../../_images/examples_workflows_hard-magnet-tutorial_41_0.png
[27]:
results_hysteresis.configurations
[27]:
{1: PosixPath('/home/petrocch/repo/mammos/mammos/examples/hystloop/hystloop_0001.vtu'),
 2: PosixPath('/home/petrocch/repo/mammos/mammos/examples/hystloop/hystloop_0002.vtu')}
[28]:
results_hysteresis.plot_configuration(1);

Analyze hysteresis loop#

We can extract extrinsic properties with the extrinsic_properties function from the mammos_analysis package:

[29]:
extrinsic_properties = mammos_analysis.hysteresis.extrinsic_properties(
    results_hysteresis.H,
    results_hysteresis.M,
    demagnetization_coefficient=1 / 3,
)
[30]:
extrinsic_properties.Hc
[30]:
CoercivityHcExternal(value=675873.2267754405, unit=A / m)
[31]:
extrinsic_properties.Mr
[31]:
Remanence(value=1511984.0060448316, unit=A / m)
[32]:
extrinsic_properties.BHmax
[32]:
MaximumEnergyProduct(value=311501.3278528385, unit=J / m3)

We can combine the results_hysteresis.plot method with some custom code to show Hc and Mr in the hysteresis plot:

[33]:
ax = results_hysteresis.plot()
ax.scatter(0, extrinsic_properties.Mr.value, c="r", label="Mr")
ax.scatter(-extrinsic_properties.Hc.value, 0, c="g", label="Hc")
ax.axhline(0, c="k")  # Horizontal line at M=0
ax.axvline(0, c="k")  # Vertical line at H=0
ax.legend();
../../_images/examples_workflows_hard-magnet-tutorial_50_0.png

Compute Hc(T)#

We can leverage mammos to calculate Hc(T) for multiple values of T.

First, we run hysteresis simulations at 7 different temperatures and collect all simulation results:

[34]:
T = np.linspace(0, 1.1 * results_kuzmin.Tc.q, 7)

simulations = []
for temperature in T:
    print(f"Running simulation for T={temperature:.0f}")
    results_hysteresis = mammos_mumag.hysteresis.run(
        mesh_filepath=mammos_mumag.mesh.CUBE_20_nm,
        Ms=results_kuzmin.Ms(temperature),
        A=results_kuzmin.A(temperature),
        K1=results_kuzmin.K1(temperature),
        hstart=(1.5 * u.T).to(u.A / u.m),
        hfinal=(-1.5 * u.T).to(u.A / u.m),
        hnsteps=30,
    )
    simulations.append(results_hysteresis)
Running simulation for T=0 K
Running simulation for T=215 K
Running simulation for T=430 K
Running simulation for T=645 K
Running simulation for T=860 K
Running simulation for T=1074 K
Running simulation for T=1289 K

We can now use mammos_analysis.hysteresis as shown before to extract Hc for all simulations and visualize Hc(T):

[35]:
Hcs = []
for res in simulations:
    cf = mammos_analysis.hysteresis.extract_coercive_field(H=res.H, M=res.M).value
    if np.isnan(cf):  # Above Tc
        cf = 0
    Hcs.append(cf)
[36]:
plt.plot(T, Hcs, linestyle="-", marker="o")
plt.xlabel(me.T().axis_label)
plt.ylabel(me.Hc().axis_label);
../../_images/examples_workflows_hard-magnet-tutorial_55_0.png

We can also show the hysteresis loops of all simulations:

[37]:
colors = colormaps["plasma"].colors[:: math.ceil(256 / len(T))]

fix, ax = plt.subplots()
for temperature, sim, color in zip(T, simulations, colors, strict=False):
    if np.isnan(sim.M.q).all():  # no Ms above Tc
        continue
    sim.plot(ax=ax, label=f"{temperature:.0f}", color=color, duplicate_change_color=False)
ax.legend(loc="lower right");
../../_images/examples_workflows_hard-magnet-tutorial_57_0.png

MODA for the workflow#

MODA