HCP tracer calculation

What follows is the calculation of the tracer correlation for an HCP crystal (ideal $c/a = \sqrt{8/3}$, $a_0 = 1$, and $\nu = 1$, all for convenience).


In [1]:
import numpy as np
from onsager import OnsagerCalc
from onsager import crystal
import matplotlib.pyplot as plt
%matplotlib inline

In [2]:
HCP = crystal.Crystal.HCP(1., chemistry="ideal HCP")
print(HCP)


#Lattice:
  a1 = [ 0.5       -0.8660254  0.       ]
  a2 = [ 0.5        0.8660254  0.       ]
  a3 = [ 0.          0.          1.63299316]
#Basis:
  (ideal HCP) 0.0 = [ 0.33333333  0.66666667  0.25      ]
  (ideal HCP) 0.1 = [ 0.66666667  0.33333333  0.75      ]

In [3]:
sitelist = HCP.sitelist(0)
vacancyjumps = HCP.jumpnetwork(0, 1.01)
for jlist in vacancyjumps:
    print('---')
    for (i,j), dx in jlist:
        print(i, '-', j, dx)


---
1 - 1 [ -1.00000000e+00  -5.55111512e-17   0.00000000e+00]
1 - 1 [  1.00000000e+00   5.55111512e-17  -0.00000000e+00]
1 - 1 [-0.5       -0.8660254  0.       ]
1 - 1 [ 0.5        0.8660254 -0.       ]
0 - 0 [ 0.5        0.8660254  0.       ]
0 - 0 [-0.5       -0.8660254 -0.       ]
0 - 0 [ 0.5       -0.8660254  0.       ]
0 - 0 [-0.5        0.8660254 -0.       ]
1 - 1 [ 0.5       -0.8660254  0.       ]
1 - 1 [-0.5        0.8660254 -0.       ]
0 - 0 [ -1.00000000e+00   5.55111512e-17   0.00000000e+00]
0 - 0 [  1.00000000e+00  -5.55111512e-17  -0.00000000e+00]
---
1 - 0 [-0.5        -0.28867513 -0.81649658]
0 - 1 [ 0.5         0.28867513  0.81649658]
1 - 0 [ 0.5        -0.28867513 -0.81649658]
0 - 1 [-0.5         0.28867513  0.81649658]
0 - 1 [ 0.5         0.28867513 -0.81649658]
1 - 0 [-0.5        -0.28867513  0.81649658]
0 - 1 [ 0.         -0.57735027 -0.81649658]
1 - 0 [-0.          0.57735027  0.81649658]
0 - 1 [ 0.         -0.57735027  0.81649658]
1 - 0 [-0.          0.57735027 -0.81649658]
0 - 1 [-0.5         0.28867513 -0.81649658]
1 - 0 [ 0.5        -0.28867513  0.81649658]

In [4]:
HCPdiffuser = OnsagerCalc.VacancyMediated(HCP, 0, sitelist, vacancyjumps, 1)

In [5]:
for state in HCPdiffuser.interactlist():
    print(state)


1.[0,0,0]:0.[0,-1,0] (dx=[-0.5,-0.28867513459481287,-0.8164965809277259])
1.[0,0,0]:1.[0,-1,0] (dx=[-0.5,-0.8660254037844386,0.0])

In [6]:
nu0 = 1
dE0 = 1
HCPtracer = {'preV': np.ones(1), 'eneV': np.zeros(1), 
             'preT0': np.array([nu0, nu0]), 'eneT0': np.array([dE0, dE0])}
HCPtracer.update(HCPdiffuser.maketracerpreene(**HCPtracer))
for k,v in zip(HCPtracer.keys(), HCPtracer.values()): print(k,v)


eneV [ 0.]
preS [ 1.]
preV [ 1.]
preT2 [ 1.  1.]
eneSV [ 0.  0.]
eneT0 [1 1]
preT0 [1 1]
eneS [ 0.]
eneT1 [ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.]
preT1 [ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.]
preSV [ 1.  1.]
eneT2 [ 1.  1.]

In [7]:
Lvv, Lss, Lsv, L1vv = HCPdiffuser.Lij(*HCPdiffuser.preene2betafree(1, **HCPtracer))

Correlation coefficient = $L_\text{ss} / L_\text{sv}$ (as $L_\text{sv} = L_\text{vv}$). Should be very close to the FCC correlation coefficient of 0.78145, for this purely isotropic diffusion case.


In [8]:
np.dot(Lss, np.linalg.inv(Lsv))


Out[8]:
array([[ -7.81205838e-01,   9.34788042e-18,  -3.00561526e-34],
       [  9.81269863e-18,  -7.81205838e-01,  -2.06843412e-18],
       [ -5.00069150e-21,  -2.05361853e-18,  -7.81451187e-01]])

In [9]:
print(HCPdiffuser.GFvalues)


{vacancyThermoKinetics(pre=[ 1.], betaene=[ 0.], preT=[ 1.  1.], betaeneT=[ 1.  1.]): array([-0.30459844, -0.07807519, -0.0780751 , -0.05208747, -0.04717321,
       -0.04427962, -0.04427942, -0.04427942, -0.04004802, -0.0387123 ,
       -0.03390849, -0.03199761, -0.03199761, -0.0313453 , -0.03052232,
       -0.02955949, -0.0289876 , -0.02898708, -0.02898708, -0.02827   ,
       -0.02652049, -0.02539412, -0.02559413, -0.02455041, -0.02455041,
       -0.02415871, -0.02374272, -0.02349651, -0.02312866, -0.02271453,
       -0.02243828, -0.02237528, -0.02210814, -0.02210814, -0.02176107,
       -0.02125203, -0.02125072, -0.02125072, -0.02069891, -0.02069891,
       -0.01997263, -0.01915803])}

In [10]:
print(crystal.yaml.dump(HCPdiffuser.GFvalues))


? !VacancyThermoKinetics
  betaene: !numpy.ndarray [0.0]
  betaeneT: !numpy.ndarray [1.0, 1.0]
  pre: !numpy.ndarray [1.0]
  preT: !numpy.ndarray [1.0, 1.0]
: !numpy.ndarray [-0.3045984378441563, -0.07807519247549728, -0.07807510314368224,
  -0.05208747289082072, -0.047173209553052486, -0.04427961918863266, -0.044279420733164286,
  -0.044279420733164286, -0.04004801870040598, -0.038712304675323134, -0.03390849411746859,
  -0.031997611830514014, -0.031997611830514014, -0.03134530105697171, -0.030522319402592056,
  -0.029559490801805685, -0.028987599285963193, -0.02898708137907353, -0.028987081379073548,
  -0.0282700014920762, -0.02652049347674716, -0.025394123802558923, -0.02559413053894741,
  -0.024550412041830297, -0.024550412041830297, -0.024158713790440225, -0.023742720504815785,
  -0.02349650870893697, -0.02312866339151506, -0.022714525195403047, -0.022438283372584703,
  -0.022375279438156537, -0.02210814352729643, -0.02210814352729643, -0.02176107403972803,
  -0.021252027188492804, -0.021250716859829214, -0.021250716859829214, -0.020698905201038707,
  -0.02069890520103871, -0.019972630348539536, -0.019158030330971115]

/opt/local/Library/Frameworks/Python.framework/Versions/3.4/lib/python3.4/site-packages/yaml/representer.py:135: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if data in [None, ()]:

In [11]:
diffcopy = crystal.yaml.load(crystal.yaml.dump(HCPdiffuser))


/opt/local/Library/Frameworks/Python.framework/Versions/3.4/lib/python3.4/site-packages/yaml/representer.py:135: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if data in [None, ()]:

In [12]:
diffcopy.Lij(*diffcopy.preene2betafree(1, **HCPtracer))


Out[12]:
(array([[ 0.73575888, -0.        , -0.        ],
        [-0.        ,  0.73575888, -0.        ],
        [-0.        , -0.        ,  0.73575888]]),
 array([[  5.74779134e-01,  -1.77191949e-17,   1.60815019e-34],
        [ -1.80611890e-17,   5.74779134e-01,   6.94257321e-18],
        [  3.67930319e-21,   6.93337495e-18,   5.74959652e-01]]),
 array([[ -7.35758882e-01,   1.38777878e-17,  -5.80917035e-36],
        [  1.38777878e-17,  -7.35758882e-01,  -6.93889390e-18],
        [ -1.28215242e-36,  -6.93889390e-18,  -7.35758882e-01]]),
 array([[ -1.65173731e-32,  -1.30519932e-48,   1.10387665e-51],
        [ -1.19747848e-48,  -1.65173731e-32,  -1.43166826e-50],
        [ -1.63586329e-50,   2.45379493e-50,  -1.95268471e-32]]))

In [13]:
print(crystal.yaml.dump(HCPtracer))


eneS: !numpy.ndarray [0.0]
eneSV: !numpy.ndarray [0.0, 0.0]
eneT0: !numpy.ndarray [1, 1]
eneT1: !numpy.ndarray [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0,
  1.0]
eneT2: !numpy.ndarray [1.0, 1.0]
eneV: !numpy.ndarray [0.0]
preS: !numpy.ndarray [1.0]
preSV: !numpy.ndarray [1.0, 1.0]
preT0: !numpy.ndarray [1, 1]
preT1: !numpy.ndarray [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0,
  1.0]
preT2: !numpy.ndarray [1.0, 1.0]
preV: !numpy.ndarray [1.0]

/opt/local/Library/Frameworks/Python.framework/Versions/3.4/lib/python3.4/site-packages/yaml/representer.py:135: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.
  if data in [None, ()]:

In [ ]: