Make reverse-supercell-matrix determination for `generate_supercell=F…

…alse` more robust, and add test for this

doped/generation.py

@@ -1285,7 +1285,6 @@ class (such as ``clustering_tol``, ``stol``, ``min_dist`` etc), or to
  f"using input structure as defect & bulk supercells. Caution advised!"
  )
 
- # else input structure is greater than ``min_image_distance`` Å in each direction, and
  # ``generate_supercell=False`` or input structure has fewer or same number of atoms as
  # doped supercell, so use input structure:
 

doped/utils/symmetry.py

@@ -301,8 +301,18 @@ def _rotate_and_get_supercell_matrix(prim_struct, target_struct):
  returns the supercell matrix to convert from the rotated prim_struct to the
  target_struct.
  """
- # first rotate primitive structure to match target structure:
- mapping = prim_struct.lattice.find_mapping(target_struct.lattice)
+ possible_mappings = list(prim_struct.lattice.find_all_mappings(target_struct.lattice))
+ mapping = next( # get possible mappings, then sort by R*S, S, R, then return first
+ sorted(
+ possible_mappings,
+ key=lambda x: (
+ _lattice_matrix_sorting_func(x[1].T @ x[2]),
+ _lattice_matrix_sorting_func(x[2]),
+ _lattice_matrix_sorting_func(x[1]),
+ ),
+ )
+ )
+
  rotation_matrix = mapping[1]
  if np.allclose(rotation_matrix, -1 * np.eye(3)):
  # pymatgen sometimes gives a rotation matrix of -1 * identity matrix, which is
@@ -311,6 +321,7 @@ def _rotate_and_get_supercell_matrix(prim_struct, target_struct):
  supercell_matrix = -1 * mapping[2]
  else:
  supercell_matrix = mapping[2]
+
  rotation_symm_op = SymmOp.from_rotation_and_translation(
  rotation_matrix=rotation_matrix.T
  ) # Transpose = inverse of rotation matrices (orthogonal matrices), better numerical stability
@@ -343,21 +354,18 @@ def _get_supercell_matrix_and_possibly_rotate_prim(prim_struct, target_struct):
  """
  try:
  # supercell transform matrix is T in `T*P = S` (P = prim, S = super), so `T = S*P^-1`:
- transformation_matrix = np.rint(
- target_struct.lattice.matrix @ np.linalg.inv(prim_struct.lattice.matrix)
- )
- if not np.allclose(
- (prim_struct * transformation_matrix).lattice.matrix,
- target_struct.lattice.matrix,
- rtol=5e-3,
- ):
+ transformation_matrix = target_struct.lattice.matrix @ np.linalg.inv(prim_struct.lattice.matrix)
+ if not np.allclose(np.rint(transformation_matrix), transformation_matrix, atol=1e-3):
  raise ValueError # if non-integer transformation matrix
 
- return prim_struct, transformation_matrix
+ return prim_struct, np.rint(transformation_matrix)
 
  except ValueError: # if non-integer transformation matrix
  prim_struct, transformation_matrix = _rotate_and_get_supercell_matrix(prim_struct, target_struct)
 
+ if np.allclose(np.rint(transformation_matrix), transformation_matrix, atol=1e-3):
+ return prim_struct, np.rint(transformation_matrix)
+
  return prim_struct, transformation_matrix
 
 

tests/data/CsPbCl3_supercell_POSCAR

@@ -0,0 +1,143 @@
+Cl Pb Cs 
+ 1.0000000000000000
+ 17.0793930733230077 0.0000002765109006 0.0000002770046049
+ 0.0000002770046049 17.0793930733230077 0.0000002765109006
+ 0.0000002765109006 0.0000002770046049 17.0793930733230077
+ Cl Pb Cs 
+ 81 27 27
+Cartesian
+ 0.0000000922525842 2.8465655583879355 2.8465655583056511
+ 2.8465655583056511 0.0000000922525842 2.8465655583879359
+ 2.8465655583879355 2.8465655583056506 0.0000000922525842
+ 0.0000001844228844 2.8465656507228032 8.5396965827466538
+ 2.8465656504759513 0.0000001845874525 8.5396965828289400
+ 2.8465656505582357 2.8465656506405188 5.6931311166935874
+ 0.0000002765931846 2.8465657430576714 14.2328276071876569
+ 2.8465657426462516 0.0000002769223208 14.2328276072699413
+ 2.8465657427285360 2.8465657429753866 11.3862621411345906
+ 0.0000001845874525 8.5396965828289364 2.8465656504759513
+ 2.8465656506405193 5.6931311166935856 2.8465656505582357
+ 2.8465656507228037 8.5396965827466538 0.0000001844228844
+ 0.0000002767577527 8.5396966751638050 8.5396966749169536
+ 2.8465657428108195 5.6931312090284543 8.5396966749992380
+ 2.8465657428931039 8.5396966750815224 5.6931312088638864
+ 0.0000003689280529 8.5396967674986737 14.2328276993579550
+ 2.8465658349811198 5.6931313013633220 14.2328276994402412
+ 2.8465658350634042 8.5396967674163893 11.3862622333048886
+ 0.0000002769223208 14.2328276072699396 2.8465657426462521
+ 2.8465657429753874 11.3862621411345888 2.8465657427285360
+ 2.8465657430576718 14.2328276071876552 0.0000002765931846
+ 0.0000003690926210 14.2328276996048100 8.5396967670872534
+ 2.8465658351456877 11.3862622334694574 8.5396967671695378
+ 2.8465658352279721 14.2328276995225238 5.6931313010341880
+ 0.0000004612629212 14.2328277919396786 14.2328277915282548
+ 2.8465659273159880 11.3862623258043261 14.2328277916105392
+ 2.8465659273982724 14.2328277918573924 11.3862623254751885
+ 5.6931311166935865 2.8465656505582357 2.8465656506405197
+ 8.5396965827466520 0.0000001844228844 2.8465656507228041
+ 8.5396965828289364 2.8465656504759509 0.0000001845874525
+ 5.6931312088638864 2.8465657428931035 8.5396966750815242
+ 8.5396966749169518 0.0000002767577527 8.5396966751638086
+ 8.5396966749992362 2.8465657428108191 5.6931312090284560
+ 5.6931313010341871 2.8465658352279717 14.2328276995225256
+ 8.5396967670872517 0.0000003690926210 14.2328276996048100
+ 8.5396967671695361 2.8465658351456868 11.3862622334694592
+ 5.6931312090284552 8.5396966749992380 2.8465657428108195
+ 8.5396966750815206 5.6931312088638864 2.8465657428931044
+ 8.5396966751638050 8.5396966749169518 0.0000002767577527
+ 5.6931313011987550 8.5396967673341067 8.5396967672518223
+ 8.5396967672518205 5.6931313011987541 8.5396967673341084
+ 8.5396967673341049 8.5396967672518223 5.6931313011987559
+ 5.6931313933690557 8.5396968596689735 14.2328277916928254
+ 8.5396968594221203 5.6931313935336227 14.2328277917751098
+ 8.5396968595044047 8.5396968595866891 11.3862623256397590
+ 5.6931313013633229 14.2328276994402394 2.8465658349811198
+ 8.5396967674163875 11.3862622333048886 2.8465658350634047
+ 8.5396967674986719 14.2328276993579550 0.0000003689280529
+ 5.6931313935336227 14.2328277917751080 8.5396968594221221
+ 8.5396968595866873 11.3862623256397573 8.5396968595044065
+ 8.5396968596689717 14.2328277916928236 5.6931313933690548
+ 5.6931314857039235 14.2328278841099785 14.2328278838631235
+ 8.5396969517569890 11.3862624179746259 14.2328278839454097
+ 8.5396969518392734 14.2328278840276923 11.3862624178100571
+ 11.3862621411345888 2.8465657427285360 2.8465657429753879
+ 14.2328276071876552 0.0000002765931846 2.8465657430576723
+ 14.2328276072699396 2.8465657426462512 0.0000002769223208
+ 11.3862622333048886 2.8465658350634042 8.5396967674163893
+ 14.2328276993579550 0.0000003689280529 8.5396967674986737
+ 14.2328276994402394 2.8465658349811194 5.6931313013633229
+ 11.3862623254751885 2.8465659273982720 14.2328277918573942
+ 14.2328277915282548 0.0000004612629212 14.2328277919396786
+ 14.2328277916105392 2.8465659273159876 11.3862623258043278
+ 11.3862622334694574 8.5396967671695378 2.8465658351456882
+ 14.2328276995225238 5.6931313010341862 2.8465658352279726
+ 14.2328276996048082 8.5396967670872517 0.0000003690926210
+ 11.3862623256397573 8.5396968595044065 8.5396968595866927
+ 14.2328277916928236 5.6931313933690548 8.5396968596689771
+ 14.2328277917751080 8.5396968594221203 5.6931313935336245
+ 11.3862624178100571 8.5396969518392751 14.2328278840276941
+ 14.2328278838631235 5.6931314857039226 14.2328278841099785
+ 14.2328278839454079 8.5396969517569907 11.3862624179746277
+ 11.3862623258043261 14.2328277916105392 2.8465659273159880
+ 14.2328277918573924 11.3862623254751867 2.8465659273982729
+ 14.2328277919396768 14.2328277915282531 0.0000004612629212
+ 11.3862624179746259 14.2328278839454079 8.5396969517569907
+ 14.2328278840276923 11.3862624178100571 8.5396969518392769
+ 14.2328278841099767 14.2328278838631235 5.6931314857039244
+ 11.3862625101449257 14.2328279762802765 14.2328279761979939
+ 14.2328279761979921 11.3862625101449257 14.2328279762802783
+ 14.2328279762802765 14.2328279761979921 11.3862625101449275
+ 2.8465656044730854 2.8465656044730849 2.8465656044730858
+ 2.8465656966433857 2.8465656968079531 8.5396966289140863
+ 2.8465657888136859 2.8465657891428213 14.2328276533550913
+ 2.8465656968079536 8.5396966289140881 2.8465656966433861
+ 2.8465657889782539 8.5396967212489567 8.5396967210843879
+ 2.8465658811485541 8.5396968135838254 14.2328277455253911
+ 2.8465657891428218 14.2328276533550895 2.8465657888136864
+ 2.8465658813131220 14.2328277456899599 8.5396968132546878
+ 2.8465659734834223 14.2328278380248285 14.2328278376956892
+ 8.5396966289140863 2.8465656966433852 2.8465656968079536
+ 8.5396967210843862 2.8465657889782534 8.5396967212489585
+ 8.5396968132546860 2.8465658813131216 14.2328277456899599
+ 8.5396967212489550 8.5396967210843879 2.8465657889782543
+ 8.5396968134192548 8.5396968134192566 8.5396968134192548
+ 8.5396969055895546 8.5396969057541234 14.2328278378602597
+ 8.5396968135838218 14.2328277455253893 2.8465658811485546
+ 8.5396969057541217 14.2328278378602580 8.5396969055895564
+ 8.5396969979244233 14.2328279301951266 14.2328279300305596
+ 14.2328276533550895 2.8465657888136859 2.8465657891428218
+ 14.2328277455253893 2.8465658811485537 8.5396968135838236
+ 14.2328278376956892 2.8465659734834219 14.2328278380248285
+ 14.2328277456899581 8.5396968132546878 2.8465658813131220
+ 14.2328278378602580 8.5396969055895564 8.5396969057541270
+ 14.2328279300305578 8.5396969979244250 14.2328279301951284
+ 14.2328278380248268 14.2328278376956892 2.8465659734834228
+ 14.2328279301951266 14.2328279300305578 8.5396969979244233
+ 14.2328280223654264 14.2328280223654264 14.2328280223654282
+ 0.0000000000000000 0.0000000000000000 0.0000000000000000
+ 0.0000000921703002 0.0000000923348683 5.6931310244410032
+ 0.0000001843406004 0.0000001846697366 11.3862620488820063
+ 0.0000000923348683 5.6931310244410014 0.0000000921703002
+ 0.0000001845051685 5.6931311167758709 5.6931311166113021
+ 0.0000002766754687 5.6931312091107378 11.3862621410523044
+ 0.0000001846697366 11.3862620488820028 0.0000001843406004
+ 0.0000002768400368 11.3862621412168714 5.6931312087816037
+ 0.0000003690103370 11.3862622335517418 11.3862622332226042
+ 5.6931310244410023 0.0000000921703002 0.0000000923348683
+ 5.6931311166113021 0.0000001845051685 5.6931311167758718
+ 5.6931312087816028 0.0000002768400368 11.3862621412168750
+ 5.6931311167758709 5.6931311166113012 0.0000001845051685
+ 5.6931312089461708 5.6931312089461699 5.6931312089461716
+ 5.6931313011164715 5.6931313012810376 11.3862622333871748
+ 5.6931312091107387 11.3862621410523026 0.0000002766754687
+ 5.6931313012810385 11.3862622333871712 5.6931313011164706
+ 5.6931313934513392 11.3862623257220399 11.3862623255574729
+ 11.3862620488820045 0.0000001843406004 0.0000001846697366
+ 11.3862621410523044 0.0000002766754687 5.6931312091107387
+ 11.3862622332226042 0.0000003690103370 11.3862622335517436
+ 11.3862621412168732 5.6931312087816019 0.0000002768400368
+ 11.3862622333871730 5.6931313011164715 5.6931313012810403
+ 11.3862623255574729 5.6931313934513383 11.3862623257220434
+ 11.3862622335517418 11.3862622332226024 0.0000003690103370
+ 11.3862623257220417 11.3862623255574711 5.6931313934513401
+ 11.3862624178923415 11.3862624178923397 11.3862624178923433

tests/test_generation.py

@@ -691,6 +691,8 @@ def setUp(self):
 
  self.conv_si = Structure.from_file(f"{self.data_dir}/Si_MP_conv_POSCAR")
 
+ self.cspbcl3_supercell = Structure.from_file(f"{self.data_dir}/CsPbCl3_supercell_POSCAR")
+
  def _save_defect_gen_jsons(self, defect_gen):
  defect_gen.to_json("test.json")
  dumpfn(defect_gen, "test_defect_gen.json")
@@ -3323,3 +3325,28 @@ def test_P1_sanity_check(self):
  rattled_struc = rattle(self.prim_cdte, 0.25)
  defect_gen, output = self._generate_and_test_no_warnings(rattled_struc)
  assert "Note that the detected symmetry of the input structure is P1" in output
+
+ def test_cspbcl3_no_generate_supercell(self):
+ """
+ Test the created supercell, primitive structure and rotation matrix for
+ a CsPbCl3 supercell input with ``generate_supercell=False``.
+
+ Previously gave a rotated supercell (but same lattice definition, just
+ rotated octahedra) due to falsely mis-matched integer transformation
+ matrices (and using ``find_mapping`` rather than ``find_all_mappings``),
+ but now fixed.
+ """
+ defect_gen, output = self._generate_and_test_no_warnings(
+ self.cspbcl3_supercell, generate_supercell=False
+ )
+ self._general_defect_gen_check(defect_gen)
+
+ assert np.allclose(defect_gen.supercell_matrix, np.eye(3) * 3)
+ assert np.allclose(
+ defect_gen.bulk_supercell.lattice.matrix, self.cspbcl3_supercell.lattice.matrix, atol=1e-4
+ )
+ assert {tuple(i) for i in np.round(defect_gen.bulk_supercell.frac_coords, 4)} == {
+ tuple(i) for i in np.round(self.cspbcl3_supercell.frac_coords, 4)
+ }
+ assert np.allclose(defect_gen._T, np.eye(3))
+ assert np.allclose(defect_gen.primitive_structure, np.eye(3) * 5.69313, atol=1e-4)