Skip to content

[ERROR] Incorrect lattice vectors for BCT structures, issues with polar correction #177

Description

@hurricane642

Dear TDEP developers,
I work with the last TDEP release (25.03).
I followed advices on the compilation of the code from the official webpage.
Then I started running the calculations for the BPO_4, material, which is identified by the TDEP as 'BCT2', which, I believe, correct.

However, if i run
extract_forceconstants -rc2 10.0
For any value of rc2, code just freezes. The reasoning behind it is, I believe, this block of code:

    case('BCT')
    ! BCT  a=b=c, al=be
        call get_al_eq_be_order(latticevectors,m0)
        call lo_get_axis_angles(m0,ap,bp,cp,ra,rb,rg)
        a = sqrt(2.0_flyt)*ap*Sqrt(-Cos(ra) - Cos(rg))
        b = sqrt(2.0_flyt)*bp*Sqrt(-Cos(ra) - Cos(rg))
        c = 2*cp*Sqrt(-Cos(rb))
        pri(:,1)=[-a/2,a/2,c/2]
        pri(:,2)=[a/2,-a/2,c/2]
        pri(:,3)=[a/2,a/2,-c/2]

        con(:,1)=[a,0.0_flyt,0.0_flyt]
        con(:,2)=[0.0_flyt,a,0.0_flyt]
        con(:,3)=[0.0_flyt,0.0_flyt,c]

And, to be more precise, lines with sqrt of cosines. The assumes that cos values will be positive, which is not true. As a result, I obtain Nan for the values of lattice.
I changed the lattice value determination with this code:

    case('BCT', 'BCT1', 'BCT2')
        call get_al_eq_be_order(latticevectors,m0)

        ap = norm2(m0(:,1))
        d12 = dot_product(m0(:,1), m0(:,2))
        d13 = dot_product(m0(:,1),m0(:,3))
        d23 = dot_product(m0(:,2),m0(:,3))
        dmin = min(d12, d13, d23)

        a = sqrt(max(0.0_flyt, ap**2 - dmin))
        b = a
        c = sqrt(max(0.0_flyt, 2.0_flyt*(ap**2 + dmin)))

        pri(:,1) = [-a/2,  a/2,  c/2]
        pri(:,2) = [ a/2, -a/2,  c/2]
        pri(:,3) = [ a/2,  a/2, -c/2]

        con(:,1) = [a, 0.0_flyt, 0.0_flyt]
        con(:,2) = [0.0_flyt, a, 0.0_flyt]
        con(:,3) = [0.0_flyt, 0.0_flyt, c]

But high-symmetry path labeling doesn't work in this case correctly, which is make sense.
I've just manually turned off the error in this subroutine, and calculation passed.

Finally, after it I run
extract_forceconstants -rc2 10.0 --polar

But I obtain negative phonon frequencies. I did the calculations for 150 structures, which i suppose should be enough. Before it I also had calculations with 250 different structures, and result still was the same. So I don't think this is the problem related to the number of the structures.
Calculation without '--polar' options works perfectly fine

I provide all supplementary files, lmk if you'll need something else. I'll be glad to obtain any recommendations and feedback from you

Best,
Sergei

bpo.zip

Activity

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions