Tutorial

This tutorial walks through a complete mag4 study, from a CIF file to Heisenberg couplings and a critical-temperature estimate — every command is real and the outputs shown are what mag4 actually prints. The running example is NiO (rocksalt, S = 1), with La₂CuO₄ (S = ½, ring exchange) for the advanced sections.

The whole loop is three commands:

mag4-magnetic --cif NiO.cif --element Ni --spin 1      # 1 · GENERATE
# … run your DFT code in every NiO/configN/ …          # 2 · RUN
cd NiO && mag4-extract --code vasp --spin 1            # 3 · EXTRACT

Everything below is refinement of that loop.

Step 0 — Install

git clone https://gitlab.com/xrocquef/MAG4.git
cd MAG4
pip install -e .

Five commands are installed: mag4-magnetic, mag4-extract, mag4-magnon, mag4-angles, mag4-moments (see the CLI reference).

Step 1 — Inspect the structure and choose the J shells

Always start with a dry run: it prints the metal–metal distance shells, the couplings mag4 would extract at the chosen --cutoff, and the supercell it would build — without writing a single file.

mag4-magnetic --cif NiO.cif --element Ni --cutoff 4.3 --spin 1 --dry-run

  Coupling    d_mean (Å)  ...  Pairs
  J1              2.9543  ...  (nearest neighbour, 90° Ni–O–Ni)
  J2              4.1780  ...  (next-nearest, 180° Ni–O–Ni superexchange)
  => 2 inequivalent coupling(s): J1, J2

  Optimal supercell : 2×2×2

Two things happen here:

  • Shell selection. Only shells with d ≤ --cutoff are extracted. Raising the cutoff adds shells and may grow the supercell. Tune the cutoff until the J set matches your physics.

  • Supercell search. mag4 finds the smallest cell in which every coupling can be isolated by the four-state formulas — including the no-self-folding criterion (no coupling may connect an atom to its own periodic image, which would contaminate another shell). For NiO the compact cells fail this test and mag4 selects 2×2×2.

If shells at the same distance hide different superexchange paths, add --ligand O to split them by their M–L–M bridging-angle fingerprint, and use mag4-angles for a Goodenough–Kanamori overview of the angles.

Step 2 — Generate the DFT inputs

mag4-magnetic --cif NiO.cif --element Ni --cutoff 4.3 --spin 1 \
    --code vasp --functional PBE+U --ldaul 2 --ldauu 6 --ldauj 0

This writes one directory per symmetry-inequivalent spin configuration (crystal symmetry + time reversal deduplicate the raw uu/dd/ud/du states):

NiO/
├── four_state_report.dat     # formulas, energy equations, MAGMOMs — the contract
├── config1/  POSCAR  INCAR  KPOINTS
├── config2/  …
└── config4/  …

The report is the heart of the method — it states, before any DFT is run, exactly how J will be computed:

  EXTRACTION FORMULAS
  J1 (d=2.9543 Å) = ( E[#1] + E[#2] - E[#3] - E[#3] ) / 4
  J2 (d=4.1780 Å) = ( E[#1] + E[#4] - E[#3] - E[#3] ) / 8

  MAGNETIC-ORDER ENERGY EQUATIONS  (per config, to check the formulas)
  E[#1] = E0  +192·(J1·S²)  +96·(J2·S²)   [|Sz|=32]
  …

Useful options at this stage:

  • --kspacing 0.3 (default) sets the k-mesh à la VASP KSPACING; metals need ~0.15–0.2.

  • --ref-couplings "J1:AFM" uses an antiferromagnetic reference bath — strongly recommended for strong antiferromagnets (the SCF converges to the intended state far more readily). --ref-magmom lets you write the bath directly as the supercell MAGMOM string.

  • --code wien2k (or both) writes configN.struct + configN.inst

    • job.init/job.run instead of VASP inputs. job.init follows P. Blaha’s workflow — setrmt, then init_lapw -sp -prec 2n -noreduc -nohdlo -nokshift with a -numk mesh matching --kspacing, and an interactive sgroup origin-shift check. -nokshift is always on, so the k-mesh is unshifted: a shifted mesh would sample different k-points from one config to the next, and J is an energy difference. For PBE+U there is also an active init_orb block built from --ldaul/--ldauu/--ldauj — U and J are converted from eV to Ry (÷ 13.6057), which is the unit init_orb expects.

Step 3 — Run the DFT calculations

Run your DFT code in every configN/. mag4 is not involved — use your cluster, your queue, your settings. The only requirements: identical numerical settings for every config (k-mesh, smearing, U, precision — J is an energy difference), and converged SCF runs.

Step 4 — Extract the couplings

cd NiO
mag4-extract --code vasp --spin 1

mag4 walks the configN/ tree, verifies convergence, reads the total energies and applies the report’s formulas:

  Configurations read: 4  (vasp outputs)
  All 4 configurations converged. ✓

  Band gaps (per config)          # metallic configs are flagged loudly
  config1    gap =  3.021 eV
  …

  name  d_mean(A)    JS2(eV)    JS2(meV)     J(meV)     J(K)   status
  J1     2.9543    -0.001000    -1.0000     -1.000    -11.60   OK
  J2     4.1780     0.015200    15.2000     15.200    176.39   OK

  J values written to: j_values.dat
  Full output saved to: mag4_extract.out

Sign convention: J > 0 is antiferromagnetic (H = +Σ J S·S). Three layers of sanity checks come free:

  • convergence per config (unconverged runs are excluded, with reasons);

  • band gap / metallicity per config (a metallic config in an insulating series signals a wrong SCF solution);

  • magnetic state: --check-mag compares each config’s converged total moment against its intended Sz and excludes inconsistent ones; mag4-moments gives the per-site view when you need to know which atom flipped.

The complete screen output is always persisted to mag4_extract.out.

Step 5 — Magnons and the critical temperature

Feed the extraction straight into the spin-wave module — no hand-copying:

mag4-magnon --cif NiO.cif --element Ni --from-extract NiO/

INFO: Couplings read from mag4-extract output: J1=-1.0000 meV.S2, J2=15.2000 meV.S2
INFO: Spin S = 1 taken from the j_values.dat header.
…
  T_c^LT-MF  =    705.6 K   (mean field — always an overestimate)
  T_c^LT-RPA =    522.0 K   (Tyablikov RPA)

For NiO the RPA value lands within a few Kelvin of the experimental T_N = 525 K. Two guards keep the numbers honest:

  • Mermin–Wagner: if the significant couplings span < 3 dimensions, mag4 states plainly that the temperature is an energy scale of short-range correlations, not a Néel/Curie temperature. A coupling counts as zero only when the fit cannot resolve it from zero (|JS²| below its own error bar or the fit RMS) — never for physical smallness alone.

  • Flat-axis q₀: components of the ordering vector along directions with no resolvable dispersion are reported as undetermined, not as spurious incommensurate modulations.

Advanced 1 — Energy mapping (all orders, least squares)

The four-state method computes a few targeted configs per bond. The energy-mapping method instead enumerates every inequivalent collinear order of the supercell and fits them all at once — with error bars, and with the option of a four-spin ring exchange column:

mag4-magnetic --cif La2CuO4.cif --element Cu --method energy-mapping \
    --supercell 2 2 1 --cutoff 7 --spin 0.5 \
    --code vasp --functional PBE+U --ldaul 2 --ldauu 8 --ldauj 0 --jring
# … run DFT …
mag4-extract --code vasp --spin 0.5

Fitted parameters  (value ± 2σ)
  J1    =       130.5059 ± 0.0028 meV   (raw fit J1*S^2 = 32.6265 meV)
  J2    =         6.4116 ± 0.0026 meV
  J3    =        -0.0008 ± 0.0017 meV   (statistically zero — interlayer)
  Jring =        32.6926 ± 0.0170 meV   (raw fit Jring*S^4 = 2.0433 meV)
  RMS error =  0.0086 meV

The huge cuprate ring exchange (J_ring ≈ 0.25 J₁) is fitted alongside the pairwise couplings. See Theoretical background for the model.

Advanced 2 — The 16-state ring method

To measure only J_ring from a handful of targeted runs, the 16-state method flips one square plaquette over a fixed bath and projects with the four-spin parity χ = s₁s₂s₃s₄ (which cancels every ≤ 3-spin term exactly):

mag4-magnetic --cif La2CuO4.cif --element Cu --method 16-state \
    --ref-couplings "J1:AFM" --spin 0.5
# … run the ~6–9 configs, then:
mag4-extract --code vasp --spin 0.5     # → Jring·S⁴ (and Jring)

The report labels every configuration with its space group and its magnetic (Shubnikov) group — BNS number, type, and magnetic point group (grey G1' groups included). --check-degeneracy keeps the time-reversal-redundant configurations as extra runs and mag4-extract verifies their degeneracy to µeV — a direct test of the antiunitary symmetry and of the spin model’s completeness.

Advanced 3 — T_c including the ring exchange

A bilinear spin-wave model cannot contain J_ring directly, but at the LSWT level about the Néel state the ring operator reduces exactly to renormalised bilinear couplings (Coldea et al., PRL 86, 5377 (2001)). One flag applies the mapping, and the report documents it:

mag4-magnon --cif La2CuO4.cif --element Cu --from-extract . --jring

  Exchange parameters  (meV.S2)   [ring-renormalised EFFECTIVE couplings]
  Shell     d_mean (A)     JS2 (meV)
  J1            3.7817       28.5399     ← 32.6265 − 2·JringS⁴
  J2            5.3481       -0.4404     ←  1.6029 − 1·JringS⁴
  --jring: LSWT ring renormalisation applied (JringS⁴ = 2.0433 meV …

Advanced 4 — Diagnostics you should know about

  • The four-state method prints contamination checks by default: how much the ring term (--no-jring to disable) and beyond-cutoff shells out to 2×cutoff (--no-jprime to disable, --jprime R to widen) leak into each extracted J, ending with a corrected-J formula (J2_true = J2_4state − 2·Jring·S² − …). They cost nothing and change nothing in the generated inputs.

  • --check-degeneracy (all methods) emits the time-reversal-degenerate pairs for explicit verification.

  • The WIEN2k back-end reads the global band gap (:GAP (global)), uses symmetric configN.struct files, and mag4-extract rescales per-config cell sizes exactly through DFT extensivity.

  • Old pre-mag4 Ising datasets can be converted with tools/ising2report.py and refitted without redoing any DFT.

Where to go next