Provide option to set an explicit ibrav value

[sphuber]

Currently this parameter from the CONTROL list is explicitly forbidden from being set directly in the input parameters, as it is always set to 0 by the calculation plugin. However, for some purposes people want to enforce a specific ibrav if for example they know their structure is FCC and want to keep that symmetry even during volume relaxation for example.

@giovannipizzi are there any objections to providing this as an option through, for example, a key force_ibrav in the settings input?

[sphuber]

CC'ing @lekah

[lekah]

CC'ing @lorisercole

[lorisercole]

It looks like we will need to implement this option, as Giannozzi deprecated the use of ibrav=0 when using symmetries. Read this issue:
https://gitlab.com/QEF/q-e/issues/131

Also reported in the release notes of QE 6.4.1:

New in 6.4.1 branch :

• A warning is issued if the lattice parameter seems to be a conversion
  factor instead of a true lattice parameter. Conversion should be achieved
  with the appropriate options, not with dirty tricks. In the future this
  will no longer be allowed
• A warning is issued if ibrav=0 is used for systems having symmetry. If not
  properly done this may lead to strange problems with symmetry detection
  and symmetrization. Lattice information should be used if available.

[giovannipizzi]

Indeed I have discussed with him and I don't agree with the term 'deprecated'. Deprecated means it's going to be removed, while here I think the point is that he wants to warn people of potential problems if not used properly. The tricky part is that it's not trivial to convert from ibrav=0 to another ibrav if you have numerical inaccuracies.

I'm ok in adding a 'force_ibrav' but how would this exactly work? is there an error if I say "force ibrav" to be cubic, but the system isn't? And for e.g. an orthorhombic system, who decides the order of the three a, b, c axes?

[lekah]

@giovannipizzi Yes, the idea is to have force_ibrav, and then check with a tolerance (i.e. force_ibrav_tol, by default 1e-5 or similar) that the cell vectors of the structure can be recreated. I.e. if someone wants the system to be cubic, i.e. take the mean cell vector length as alat. Then I can check if the \sum_i (v_i - alat)^2 < tol, where v is the cell vector before any conversion using the forced ibrav. Then of course, certain standard must be enforced, i.e. that in a tetragonal system a=b \neq c, and a < b < c for orthorhombic. Can we take this as a reference for this, you think: https://arxiv.org/abs/1506.01455 ?

[giovannipizzi]

I see. I would strongly advise not to reimplement the logic in that paper. This is now part of spglib + seekpath; this however means an stronger dependency on these two libraries. There might still be the need to do some patching between the plugin and seekpath (e.g. seekpath also standardises with respect to the Cartesian axes, while with this flag you need to accept the orientation of axes decided by QE - but if a user wants this explicitly with a force_ibrav flag, it's fine to do it probably.

Any takers for an implementation?

[greschd]

Is the goal here to determine ibrav automatically and have force_ibrav just be an on / off switch, or just to allow specifying ibrav manually?

For the latter, I think all that is needed is to convert between matrix form and the celldm QE input for the structure. If the given structure doesn't match the QE convention, we can just reject it.
@giovannipizzi does seekpath or spglib have tools for this conversion, or did you mean something else?

[giovannipizzi]

I think the idea is to have the user specify the ibrav. Indeed, maybe it's OK as you mention to just have an hardcoded list of ibrav cells, and check that it matches? This would go in the settings (maybe together with the celldm or a,b,c values)? Even all in a single dict:
force_ibrav = {'ibrav': 2, 'celldm': [2,0,0,0,0,0]} or something like this.

Note: the logic for this is already in qe-tools

[greschd]

Even all in a single dict:
force_ibrav = {'ibrav': 2, 'celldm': [2,0,0,0,0,0]} or something like this.

I think we should put it in the SYSTEM key of the parameters dict, to be consistent with all the other inputs.

We could allow any of ibrav, A, B, C, cosAB, cosBC, cosAC and celldm, with some validation (e.g., A, B, C, cosAB, cosBC, cosAC is mutually exclusive with celldm). Or would you prefer only one of the two being available? @sphuber, what's your opinion?

Note: the logic for this is already in qe-tools

Nice!

[greschd]

Any takers for an implementation?

✋

Once we're agreeing on the interface I'll have a crack at this.

[greschd]

Actually, regarding the interface: Since we anyway want a StructureData to be passed, we could also allow only ibrav to be set, and then get the other values from the cell.

Of course, we would then need to validate that the cell matches what QE will produce.

Does qe-tools also provide a "get parameters from cell", or would we need to implement that ourselves?

[giovannipizzi]

Of course, we would then need to validate that the cell matches what QE will produce.
Does qe-tools also provide a "get parameters from cell", or would we need to implement that ourselves?

That was the part that worried me most ;-) I'm not sure qe-tools has it (@lekah correct me if I'm wrong). If you feel confident trying to implement it I'm ok, but I fear we will end up in a lot of tiny checks and thresholds... but OK to try, maybe I'm too pessimistic :-D

[sphuber]

I think the real problem is exactly going to be to define what it means for the StructureData to match how QE interprets the ibrav parameters. If we are going to add validation, then we will have (multiple) thresholds like giovanni said. These most likely should then also be configurable by the user, maybe through the settings node, but this all adds a lot of complexity of course.

[greschd]

I would just extract the parameters without checks (e.g. for hexagonal, celldm(1) is the norm of cell[0], and celldm(2) = (norm of cell[2]) / celldm(1)).

Then add one check at the end, whether qe-tools's get_cell_from_parameters produces the initial cell, of course with some tolerance on the unit cell vectors. But that is not that complex..

[greschd]

I think 99% of use cases would be covered by just allowing a (single) absolute tolerance on the individual unit cell entries to be set by the user.

[giovannipizzi]

I guess the only way to know is to see an implementation :-D

[greschd]

@sphuber could you assign me, so I don't forget about it?

[sphuber]

With pleasure :)