Improvement of Jij calculation forextended basis sets having large
non-orthogonality:
current situation and problem
Asako Terasawa
Research Organization for Information Science andTechnology (RIST)
2023-11-09, The 3rd OpenMX Developers Meeting
1/15

Experimental background: complexity ofmulti-phase structure of permanent magnet
To create strong permanent magnets…
• Main phase: strong ferromagnetism and strong magnetic
anisotropy• Grain boundary phase: magnetic insulation of main phasegrainsIf not, magnetic reversal
easily spread in thepresence ofdemagnetization eld

Magnetic insulation mechanism
Mat. Trans. 57 (2016), 1221-1229
2/15fi•Computational study ofpermanent magnets
Calculation of exchange coupling constant Jij••for systems including rare earth elements, e.g. Nd, Sm, Dy
for suf ciently large systems that can represent different
phases and multi-phase structures of permanent magnets• Minimum ~ 50 atoms• Maximum ~ a few hundreds or a few thousands of atoms
A bit too large,
but anyway this cannot be
attained without development of code

3/15
fi

•

jx: Jij calculation code
for OpenMX post-processing•Post-process calculation ofexchange coupling constant Jij using
ground state density and second
perturbation theory [J. Magn. Magn. Mat 67,65–74 (1987)]

••

•

•

Applicable for isolated systems andperiodic systems
MPI parallelization for periodicsystems
Ef cient algorithm for energy
integration
To use jx, you need to choose
relatively small basis set.

Im (z)
C2

Poles ofFermi
functions

SC1
Re(z)
maxmaxPoles of retarded
Green’s functions[T. Ozaki, PRB,75, 035123 (2007)][AT et al., JPSJ 88, 114706 (2019)]

fi4/15

Jij for transition metals
Fair correpondence with experimental Curie temperature•
181616141410J ij [meV]

J ij [meV]122.866 Å

128644.070 Å

10862.507 Å

422x0-2

23456

rij [Å]708910(a) bcc Fe, k-grid = 27×27×27-2

2

332.5Systembcc Fehcp Co
fcc Ni3.530 Å2456

rij [Å]78910(b) hcp Co, k-grid = 32×32×20

3.5

J ij [meV]•1.5
1TC [K]
calculated experimental
1321
1040
1640
1131
4456270.50-0.5

23456

rij [Å]78910

XC: GGA-PBE
PAO: s2p2d2(c) fcc Ni, k-grid = 24×24×245/15

Jij for rare earth metals,
for various basis sets
410fcc Nd

20-10

Jij [meV]Jij [meV]0-2-20-4-30-6-40-8-50-10-60-12-70
2345678910

hcp Sm2345

rij [A]s2p1d1s2p1d2•s2p2d1s2p2d2678910rij [A]s3p2d2s2p1d1s2p1d2s2p2d1s2p2d2s3p2d2

However, jx seems to work bad in calculations of rare earth
metals…
•

Strong dependence on basis sets in density functional calculation•Diverges when choosing large number of basis!
6/15

Problem
What is the origin of divergence in rare earth metal
calculation?
• How can we eliminate it?
• How can we calculate very accurate and converged Jij
for large number of basis?
•What is clear is that …
Electronic state itself will be converged when choosing
very large number of basis• This problem comes from the calculation method of Jij•

7/15

Non-Orthogonality problem
Implemented ecuation to calculate Jij using Liechtensteinmethod [J. Magn. Magn. Mat 67, 65–74 (1987)][AT et al., JPSJ 88, 114706 (2019)]•Jij<latexit sha1_base64="4qh1EGohIftE+dBKcCAzj39mOrU=">AAADHXicdVJNbxMxEPUuX2X5aApHLhYRTYEqyqaolEOlShxAHFBApK0Upyuv402c7nq3/qAKlv8IF/4KFw4gxIEL4t9gJytIKYy00ui9mfdmx5NWOZOq0/kZhBcuXrp8ZeVqdO36jZurjbVb+7LUgtA+KfNSHKZY0pxx2ldM5fSwEhQXaU4P0uOnnj94S4VkJX+jZhUdFnjMWcYIVg5K1oItlNIx44aecCwEnj2w0YvEsKmN1ncjuI4ygYmJrenaCEIkdZGYaje2R+ZlYlCB1UQUpmetRYrlI2peW8fbRR0q9CbiGjEO2R+o5bEW9KgzQTnNFDJOe4AmWDkpZ26HrtI1FtqbPjsyDz3q5aa+126gUXnqxy1PN2vfd3Pf+2d0prVOy7suKU3nUzDHOCVd/UcHCTaeKGTbEaJ8tLyfKGk0O+3OPOD5JK6TJqijlzS+u4mJLihXJMdSDuJOpYYGC8VITt0WtKQVJsd4TAcu5bigcmjmr2vhPYeMYFYK93EF5+hyh8GFlLMidZX+PeTfnAf/xQ20ynaGhvFKK8rJwijTOVQl9KcCR0xQovKZSzARzM0KyQS7e1DuoBZLeOJj+/cvn0/2u+14q/3oVbe5t1OvYwXcAXfBBojBY7AHnoMe6AMSvA8+Bp+DL+GH8FP4Nfy2KA2Duuc2OBPhj18zuv+m</latexit>=NPX
12 p=1R̃p

X

X nµ,⌫2i µ0 ,⌫ 0 2jo+00[P̂i ]⌫µ G+(#,z̃)[P̂]
Gpj ⌫ µiµ,j⌫ 0jµ0 ,i⌫ (", z̃p ) .

Contribution from atomic sites i and j•Currently, in Pi and Pj, we justput the difference of diagonal
element of Hamiltonian for spin
up and down.

•For large basis sets, Atomic
sites will be ill-depicted by thisde nition, because of large
overlap of wide basis functions

fi

8/15

Atomic sites by valence state
First, we tried de nition of atomic sites by valence state only.
In this method, we sum up the constribution of valence state
for Liechtenstein calculation, while using extended basis set
for electronic state calculation.

••

30

825

6

fcc Nd4

Jij [meV]Jij [meV]

20
15
10
5

20-20-4-5-6-10-82345678910

hcp Sm234

rij [A]s2p1d1s2p1d2•678910rij [A]s3p2d2s2p1d1s2p1d2This didn’t solve the problem!
Strong dependence on basis set for 1st NNfi•

s2p2d1s2p2d25s2p2d1s2p2d2s3p2d2

9/15

Orthogonalization of atomic orbitals•

Then, we tried orthogonalized atomic orbitals:
X|ii(LO) ⌘|j 0 iS 1/2 |j 0 i .
Löwdin orthogonalization (LO)j0<latexit sha1_base64="bybtTbevSOlunEKH2OXwEyHNTP0=">AAACPnicbVBNT9wwFHQobSH9YCnHXixWFfTQbbKt+LitxKUHJEB0AWkTVo73ZTHYTmo7SCtvfhmX/obeOHLhAKp65YgToqotHcnSaN48ed4kOWfaBMGlN/Nk9umz53Pz/ouXr14vtBbfHOisUBT6NOOZOkqIBs4k9A0zHI5yBUQkHA6Ts61qfngOSrNMfjWTHGJBxpKljBLjpGGrHw38KWaRInLM4diuRoKYEyXs9k75vvQj+Fawcz/ShRja05VyerrSWHFtTFK7Xx7bD+HHbjmtHKzs+FE8bLWDTlADPyZhQ9qowe6w9SMaZbQQIA3lROtBGOQmtkQZRjm4HIWGnNAzMoaBo5II0LGtzy/xO6eMcJop96TBtfrnhiVC64lInLPKrP+dVeL/ZoPCpBuxZTIvDEj68FFacGwyXHWJR0wBNXziCKGKuayYnhBFqHGN+3UJmxXWfp/8mBx0O+Gnzue9brvXa+qYQ2/RMlpFIVpHPfQF7aI+ougCXaEbdOt99669n96vB+uM1+wsob/g3d0D43+vDA==</latexit>

Schematic of nonorthogonal basis
function

bcc FeSchematic ofLöwdin basis
function

dhcp Nd25

415

Jij [meV]Jij [meV]20

1050-20-5-102-4234567891034

rij [A]s2p1d1s2p2d1s2p2d2s3p2d2s3p2d2f1
s3p3d25678910rij [A]s3p3d2f1
s3p3d3s3p3d3f1

Do not converge at all!

s2p1d1s2p1d2s2p2d2s3p2d2s3p2d2f1
s3p2d3s3p2d3f1
s3p3d3s3p3d3f1

10/15

Orthogonalization of atomic orbitals•

We also tried another orthogonalization method:
Single-site orthogonalization (SO)XX

Schematic of
single-site
orthogonalizedbasis function

(SO)
|iii

⌘

|ii

|jii

(SO)

⌘|ji, j 6= i,Si,i1 :Submatrix of overlap matrix for i-site element<latexit sha1_base64="sWfVlWnECFipTdWkVsi/ZMOpCuk=">AAADOHicbVJLb9NAEF6bVzGvFI5cVkQ0BaVRHBAvcajEhRMUlbSVssFabybpJuu1u7uuFG38s7jwM7ghLhxAiCu/gLXrQGg7kqWZb2a+mW+8cSa4Nt3uF8+/cPHS5StrV4Nr12/cvNVYv72n01wx6LNUpOogphoEl9A33Ag4yBTQJBawH89elfn9Y1Cap/K9mWcwTOhE8jFn1DgoWvfekBgmXFo4klQpOn9YBAuOiaJyIiCyvPhgN0lCzaFK7O7b4kERbBA4yvlxgDcWvK4LsLMtovMkstMWkYB5sYzqsCrBeIGnLbzahXHFHo/tbhFZkrplSy1ucHs1cGtshcWiJHSUSzYiKiLH2cL/lnFoDMIpemG4GMEMTBQWASHBYrocfVpUKXRV2LSuaxPsekoB7fNoe44W5GjleFGj2e10K8NnnbB2mqi2najxmYxSlicgDRNU60HYzczQUmU4E+D4cw0ZZTM6gYFzJU1AD2314wt83yEjPE6V+6TBFbraYWmi9TyJXWWpVp/OleB5uUFuxs+GlsssNyDZyaBxLrBJcfmK8IgrYEbMnUOZ4m5XzA6posy4txZUR3he2pO/ks86e71O+Kjz+F2vuf2yPscauovuoU0UoqdoG71GO6iPmPfR++p99374n/xv/k//10mp79U9d9B/5v/+A0baBNg=</latexit>

j 0 6=i j 00 6=i|j 0 iSi,i1 |j 0 j 00 hj 00 |ii<latexit sha1_base64="J9+wWpHiJJ1wFtLnfHNhaaPc24k=">AAACWHicbVHLTsMwEHTCu7xKOXKxKEgcoEoA8TohceEIggJSEyrH3YCFEwfbQa2s/CQSB/gVLjghQrzmNJ7dWa/HUcaZ0p736rhj4xOTU9Mzjdm5+YXF5lLrSolcUuhSwYW8iYgCzlLoaqY53GQSSBJxuI4eTsr69RNIxUR6qUcZhAm5S1nMKNFW6jdF0GsECdH3UWwuir4JhO0uhxlWbH4/FLdmyy+OgsecDKwjEkNzkUfWKdkQixhXrSTDtRILidfY2pZiGjBwSCDVRSMI+8221/Eq4L/Er0kb1TjrN5+DgaB56aecKNXzvUyHhkjNKAc7MleQEfpA7qBnaUoSUKGpginwulUG1S6xSDWu1O8OQxKlRklkO8sM1O9aKf5X6+U6PggNS7NcQ0o/L4pzjrXAZcp4wCRQzUeWECqZ3RXTeyIJ1fYvGlUIhyX2vp78l1xtd/ydzu75dvv4oI5jGq2gVbSBfLSPjtEpOkNdRNELenfGnQnnzUXulDvz2eo6tWcZ/YDb+gAdVbVE</latexit>

bcc Fedhcp Nd

Become slightly smaller as increasing the number of basis functions

10/15

Spin population scaling•

Because it turns out that the SO basis underestimate the spin
population slightly, we tried spin-population scaling for SO results

Spin-population scaling

(SOS)
Jij
=<latexit sha1_base64="juxdwFppekxTNQrEK+4uacB1NM0=">AAACinicjVHdSsMwGE3rz7ROrXrpTXAIeuHopkxFhIFeiDcqc5uw1pFm6ZYtTUuSCqP0YXwl73wb0zlEt114IHByvu8kX078mFGpHOfTMJeWV1YLa+vWRnFza9ve2W3JKBGYNHHEIvHiI0kY5aSpqGLkJRYEhT4jbX90k9fbb0RIGvFnNY6JF6I+pwHFSGmpa7+7Heu+m9Jh9pq6IVIDEaZHjYfGcZZZ124gEE4t95YwhSDXbZmV/d7PuHJTZs25hv9xLRgiL7he1y45ZWcCOE8qU1ICUzx27Q+3F+EkJFxhhqTsVJxYeSkSimJG9JGJJDHCI9QnHU05Con00kmUGTzUSg8GkdCLKzhRfztSFEo5Dn3dmQ8qZ2u5uKjWSVRw4aWUx4kiHH9fFCQMqgjm/wJ7VBCs2FgThAXVs0I8QDpJpX/PmoRwmaP28+R50qqWK6fls6dqqV6fxrEG9sEBOAIVcA7q4A48gibARsE4MWrGuVk0q+alefXdahpTzx74A/P2C02txFc=</latexit>

bcc Feni(SO)nj

nj(SO)
J(SO) ijnj

dhcp Nd

Good convergence, but physical meaning?

12/15

Discussion
jx: fails when taking large number of basis functions• Rede nition of atomic sites by valence state only
• Works bad for rare earth metals, particularly when thenumber of basis function is very large.
• Rede nition of atomic sites by orthogonalized orbitals• Löwdin orthogonalization scheme does not converge Jij atall.
• Single-site orthogonalization scheme underestimates thespin population slightly, and it affects the calculated value ofJij.
••

Although the convergence can be improved by spin
population scaling, it makes the physical meaning ofcalculated values unclear.

fifi

13/15

Discussion•

Physically-meaningful de nition of atomic sites?
−θ/2−θ/2+θ/2+θ/2
i

i
j

i

i
j

j

j

Because Liechtenstein method is based on Heisenberg
model of localized spins, we need to de ne atomic sites
by the spin rotation on the atoms.
• How does it rotate? …
• De nition of atomic sites by its electronic occupation• We should try (or, have tried) Wannier functions…
•fififi

14/15

Future works•

Implementation of closest Wannier functions to a given
set of localized orbitals• https://arxiv.org/abs/2306.15296

We then should make our way by just, simply coding.

15/15