Lattice custom orient: confusion/bug

Hello,

I’d like to do a simple test to calculate the cohesive energy of HCP-Mg for different directions. Using this command to create lattice

lattice custom 4 a1 1 0 0 a2 0 sqrt(3) 0 a3 0 0 sqrt(8/3) basis 0.0 0.0 0.0 basis 0.5 0.5 0.0 basis 0.5 0.83333333 0.5 basis 0.0 0.33333333 0.5 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1

I get: Cohesive energy per atom=-1.510114302, a= 3.202795516, and c/a=1.619608977 which seem to be right for me.

The problem rises when I want to find orguments of x, y, and z for orient command above for one of pyramidal slip systems; e.g. {1-101}<11-20>. Here is what I figured to do to find the directions for lammps:

  1. I found another direction in {1-101} plane which is <-1102>.

  2. I converted these two directions from [uvtw] format to [UVw] and then to cartesian coordinate system as a1=<100> and a2=<012> keeping in mind that unit vectors in cartesian coordinate system aren’t really unit from the definitions of a1, a2, and a3.

  3. I calcuate a3=a1 cross a2=-2(c/4)/sqrt(3)j+sqrt(3)/(c/a). and found the closest integer number for this as <0 -7 4>

  4. recalculated a2=a3 cross a1= <0 4 7>

Then I used the following commands:

variable facx equal sqrt((v_x1)^2+(v_x2v_ay)^2+(v_x3v_ca)^2)

variable facy equal sqrt((v_y1)^2+(v_y2v_ay)^2+(v_y3v_ca)^2)

variable facz equal sqrt((v_z1)^2+(v_z2v_ay)^2+(v_z3v_ca)^2)

variable latdisx equal “1v_latparamv_facx-0.05”

variable latdisy equal “1v_latparamv_facy-0.05”

variable latdisz equal “1v_latparamv_facz-0.05”

lattice custom 4 a1 1 0 0 a2 0 sqrt(3) 0 a3 0 0 sqrt(8/3) basis 0.0 0.0 0.0 basis 0.5 0.5 0.0 basis 0.5 0.83333333 0.5 basis 0.0 0.33333333 0.5 orient x 1 0 0 orient y 0 4 7 orient z 0 -7 4

region box block 0.05 {latdisx} 0.05 {latdisy} 0.05 ${latdisz} units box

and got Cohesive energy=-1.502639033

I am pretty much sure that I don’t have the periodic simulation box. I guess the problem is that the unit vectors are not actually unite in the lattice command (a1, a2, and a3).

Is there anybody here who worked on HCP structures before and see where I went wrong?

Thanks,

Ebrahim

[…]

The problem rises when I want to find orguments of x, y, and z for orient command above for one of pyramidal slip systems; e.g. {1-101}<11-20>. Here is what I figured to do to find the directions for lammps:

  1. I found another direction in {1-101} plane which is <-1102>.

  2. I converted these two directions from [uvtw] format to [UVw] and then to cartesian coordinate system as a1=<100> and a2=<012> keeping in mind that unit vectors in cartesian coordinate system aren’t really unit from the definitions of a1, a2, and a3.

I don’t think these conversions were correct.

Converting [UVW] to [uvtw]:
3u = n*(2U-V)
3v = n*(2V-U)
t = -(u+v)
w = n*W

If opposite (converting [uvtw] to [UVW]):
nU = 2u-v
n
V = 2v-u
n*W = w

where n is the smallest integer so that u, v, t, w are integers.

Neither [11-20] nor [-1102] convert correctly to [100] or [012].

Ray

Ray,
You are right but I went one step more to convert to prependicular coordinate system because the angle between U and V is 120. So, I used x=(U+V)/2, y=(V-U)/2, and z=w. Also, I think the correct opposite relations are:
nU=2u+v
n
V=2v+u

Ebrahim

Ebrahim,

Yes, you are right about the opposite relations.

Try with this “region box block 0.05 1.05 0.05 7.05 0.05 13.05 units lattice”

If it does not work, then something is still wrong with the indices.

Cheers,
Ray

Ray,

It didn’t work. Would you please let me know how did you get those numbers that I can double check them.

Ebrahim

Yes, as I said before, if the indices were not correct so that the calculated spacings were not correct, then the number I gave you will not work.

As to how I came up with the numbers:

1). Lattice spacing in x,y,z from LAMMPS before orient: A1, B1, C1
2). Lattice spacing in x,y,z from LAMMPS after orient: A2, B2, C2.
3). Index used in orient for, say, C: h, k, l

Then mfactor_C = C2/{C1sqrt[ (h/n)^2+(k/n)^2+(l/n)^2] },

where n is an integer so that the max of (h/n), (k/n) or (l/n) is 1, and m is an integer so that m*factor_C is an integer. m is then the number to replicate.

Ray