|
|
本帖最后由 tao11 于 2025-11-25 21:50 编辑
我在使用cp2k计算弹性性质时,发现没有现成的计算程序,
就考虑cp2k是否能像vasp那样,通过计算应力(σ)和应变(ε)的比例系数,即弹性张量Cᵢⱼ
我的具体操作如下(在ai的帮助下完成,使用有限差分法进行计算):1、以Cu晶体为例,通过改变晶胞的晶格系数施加应变,在保持晶胞不变的前提下进行结构优化 GEO_OPT,并分别计算拉伸和压缩应力 &STRESS_TENSOR
strain_x_neg表示压缩,strain_x_pos表示拉伸
1、2、3分别表示在A、B、C方向的压缩和拉伸
4、5、6分别表示在BC、AC、AB方向的压缩和拉伸
2、cp2k计算后,只列举strain_1_neg和strain_1_pos的计算结果:
(1)strain_1_neg:
STRESS| Analytical stress tensor [GPa]
STRESS| x y z
STRESS| x -2.63765793513E+01 1.58809768258E-02 1.30689040856E+01
STRESS| y 1.58809768258E-02 1.61588502990E+02 9.68472654766E-03
STRESS| z 1.30689040856E+01 9.68472654766E-03 -2.50782827150E+01
(2)strain_1_pos:
STRESS| Analytical stress tensor [GPa]
STRESS| x y z
STRESS| x 1.54834198652E+02 1.13067181226E-03 -1.36712827270E-03
STRESS| y 1.13067181226E-03 -2.67613582510E+01 1.19137343314E+01
STRESS| z -1.36712827270E-03 1.19137343314E+01 -2.67683212736E+01
3、通过python提取out文件计算得到的Cij如下,可以看出数值特别不正常,计算得到的弹性性质要么过大、要么过小、要么是负值。尝试了很多次,也尝试过输入文件的计算类型为能量+力(ENERGY_FORCE),但计算得到的结果都不对
18121.078 -9503.612 -170.128 682.342 -645.744 -2272.390
-9503.612 18120.729 -18835.472 7.007 -509.598 -400.299
-170.128 -18835.472 18116.877 7.641 -639.552 2518.376
682.342 7.007 7.641 -135.444 1.069 1.357
-645.744 -509.598 -639.552 1.069 -138.951 -7.334
-2272.390 -400.299 2518.376 1.357 -7.334 -1279.485
我不确定这样计算思路是否正确,请大佬们帮忙指点下。
|
|