bugs in cell-relax with symmetry=1

[QG-phy]

Describe the bug

I am using abacus to run the cell-relax job for a bulk MoS2 system.
The initial structure is adopted from materials project and give honeycomb pattern viewed in c axis.

1. if i run with symmetry=1, the final relaxed configuration is shown in the following figure, where the symmetry is broken.

2. if i set symmetry =0, the final relaxed configuration give the correct symmetry as is shown in the following

Expected behavior

No response

To Reproduce

abacus version : 3.6.0
image: registry.dp.tech/dptech/abacus:v3.6.0
mos2.bulk.tar.gz

Environment

No response

Additional Context

No response

Task list for Issue attackers (only for developers)

• Verify the issue is not a duplicate.
• Describe the bug.
• Steps to reproduce.
• Expected behavior.
• Error message.
• Environment details.
• Additional context.
• Assign a priority level (low, medium, high, urgent).
• Assign the issue to a team member.
• Label the issue with relevant tags.
• Identify possible related issues.
• Create a unit test or automated test to reproduce the bug (if applicable).
• Fix the bug.
• Test the fix.
• Update documentation (if necessary).
• Close the issue and inform the reporter (if applicable).

[maki49]

This is because the default value of symmetry_prec (1e-6) is too small for the cell-relaxed structure.
Setting symmetry_prec 1e-5 can this case converge with the original symmetry in 4 steps.

[QG-phy]

ok, then the question is how or when would someone know he/she should increase the value of symmetry_prec? I mean this time the structure is too different and the discrepancy is too obvious. Besides, I get this message after the relaxation. Is there any signal that can be detected in advance?

[maki49]

how or when would someone know he/she should increase the value of symmetry_prec?
Is there any signal that can be detected in advance?

There are some fatal failures result from an improper symmetry_prec that produce warnings and suggest changing symmetry_prec. For example:

abacus-develop/source/module_cell/klist.cpp

Lines 114 to 118 in 24a27b9

	ModuleBase::WARNING_QUIT("K_Vectors::ibz_kpoint", "Possible solutions: \n \
	1. Refine the lattice parameters in STRU;\n \
	2. Use a different`symmetry_prec`. \n \
	3. Close symemtry: set `symmetry` to 0 in INPUT. \n \
	4. Set `symmetry_autoclose` to 1 in INPUT to automatically close symmetry when this error occurs.");

abacus-develop/source/module_cell/module_symmetry/symmetry_basic.cpp

Lines 1036 to 1037 in 24a27b9

	GlobalV::ofs_running << "\n WARNING: Symmetry operations cannot completely constitute a point group.\n\
	It'll be better to try another `symmetry_prec`.\n Now search the subgroups ..." << std::endl;

abacus-develop/source/module_cell/module_symmetry/symmetry.cpp

Lines 1364 to 1365 in 24a27b9

	GlobalV::ofs_running << "WARNING: Original cell may have more than one primitive cells, \
	but we have to treat it as a primitive cell. Use a larger `symmetry_prec`to avoid this warning." << std::endl;

Symmetry breaking during (cell-)relax is not regarded as a fatal failure, but it still produces warning in OUT.${suffix}/running_${calculation}.log. In this case:

WARNING: Symmetry cannot be kept due to the lost of accuracy with atom position during cell-relax.
Continue cell-relax with a lower symmetry. 

One should adjust (often increase) symmetry_prec if the warning above is produced while symmetry-keeping is expected.

[QG-phy]

其实我还有一点问题。
如果是我设置的分别率要求太高。那么，对于一个结构其真实的空间群为A，对称性精度要求过高，应该是找到了A的 一个子群记作B。

我们已知如下两个事实：

1. 在正确的考虑了群A的约束下优化的结果是正确的，符合预期的。
2. 在完全不考虑对称性（等价为考虑 C1或 P1 的对称性，C1 /P1 自然也是 A的一个子群）优化也是对的。

那么问题来了，为什么在 B 子群下的优化结果是不对的？或者说不符合预期的？

[maki49]

那么问题来了，为什么在 B 子群下的优化结果是不对的？或者说不符合预期的？

本质上是k点在布里渊区采样不均匀造成的误差。

• 关对称性：k点在完整BZ中均匀采样。
• 开对称性且分析出完整空间群A：只计算不可约BZ（IBZ）中的k点的波函数，但之后会对电荷密度做对称化（遍历A的对称操作，旋转得到等价k点的密度，并对所有等价k点取平均），效果相当于在完整BZ中均匀采样。
• 在这个例子中，初始离子步分析出完整空间群A并据此约化k点，第二个离子步分析出子群B，但（由于程序框架所限）k点数没有相应增加（只计算A约化后的IBZ），而在对称化密度时只遍历了B的对称操作，效果相当于k点在BZ中的不均匀采样。

[QG-phy]

本质上是k点在布里渊区采样不均匀造成的误差。

OK，明白了。就是说第一步分析出的对称性是完整空间群 A， 一步离子步优化之后，第二步的对称分析为子群B, 但是不可约的K点不会更新，后面进行k点对称性还原的时候，用的群A的不可约K点，但是却用子群B的对称性操作。

那么这个难道不是个bug吗？ 既然每一个离子步都做对称性分析，如果检查到和上一步对称性不同，就要进行约化 k点的更新，不更新 k 点，却只更新对称性，这个肯定就有问题。

现在的逻辑下，只要在开了对称性结构优化，后面的离子步结构与初始的结构对称性不一致，就会出现问题啊。

[maki49]

现在的逻辑下，只要在开了对称性结构优化，后面的离子步结构与初始的结构对称性不一致，就会出现问题啊。

也不算bug，算是误差或者近似处理。
支持更新k点是一个new feature request，需要对程序框架做不少改动，已经超出了对称性模块的范围，如有需求可以另提一个issue。但也不能为了逻辑一致让cell-relax的过程中不更新对称性，否则无法处理部分方向晶格固定的情况。

[QG-phy]

对称性和不可约k点应该是要保持一致的。 如果不一致了，我感觉怎么看都是错的。我反复思考了一下，如果非要说是近似的话，在什么样的近似下，可以让对称性和不可约k点不一致呢？ 或许尺寸比较大的情况？k 点多一个少一个影响很小？感觉有点牵强。

我个人想法，可能最好的/最正确的做法，应该是更新对称性同时更新不可约k点。

至于这个是属于"new feature request" 还是bug-fix，我想具体看法可能会因人而异。至于说要不要重新提一个issue，在我看来讨论历史，算例复现等都已经在这个issue里面，应该没必要再去提一个 issue？

@mohanchen @WHUweiqingzhou 陈老师，巍青 fyi, 你们的看法呢？

[maki49]

我个人想法，可能最好的/最正确的做法，应该是更新对称性同时更新不可约k点

我同意。

[WHUweiqingzhou]

We will discuss this issue at Issue meeting this Thursday.

[chenjunh]

I rerun the mos2 cell-relax calculation and I do confirm @maki49's statement that ABACUS doesn't increase k points sampling when the symmetry of mos2 changes from high symmetry D_6h to low symmetry C_1h starting from ION STEP 2, where the irreducible k points (nkstot_ibz) should be 254 and 42 for D_6h and C_1h (for the k-spacing provided by @QG-phy), respectively.

[WHUweiqingzhou]

I will put a PR to force close symmetry for cell-relax calculations until this bug solved.

[pxlxingliang]

I use ABACUS-ASE to do the cell-relax with symmetry=1, and the final structure seems no abnormal. This indicates that it is feasible to initial the symmetry in relax step.
ase.zip