coset_table.py 42 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798991001011021031041051061071081091101111121131141151161171181191201211221231241251261271281291301311321331341351361371381391401411421431441451461471481491501511521531541551561571581591601611621631641651661671681691701711721731741751761771781791801811821831841851861871881891901911921931941951961971981992002012022032042052062072082092102112122132142152162172182192202212222232242252262272282292302312322332342352362372382392402412422432442452462472482492502512522532542552562572582592602612622632642652662672682692702712722732742752762772782792802812822832842852862872882892902912922932942952962972982993003013023033043053063073083093103113123133143153163173183193203213223233243253263273283293303313323333343353363373383393403413423433443453463473483493503513523533543553563573583593603613623633643653663673683693703713723733743753763773783793803813823833843853863873883893903913923933943953963973983994004014024034044054064074084094104114124134144154164174184194204214224234244254264274284294304314324334344354364374384394404414424434444454464474484494504514524534544554564574584594604614624634644654664674684694704714724734744754764774784794804814824834844854864874884894904914924934944954964974984995005015025035045055065075085095105115125135145155165175185195205215225235245255265275285295305315325335345355365375385395405415425435445455465475485495505515525535545555565575585595605615625635645655665675685695705715725735745755765775785795805815825835845855865875885895905915925935945955965975985996006016026036046056066076086096106116126136146156166176186196206216226236246256266276286296306316326336346356366376386396406416426436446456466476486496506516526536546556566576586596606616626636646656666676686696706716726736746756766776786796806816826836846856866876886896906916926936946956966976986997007017027037047057067077087097107117127137147157167177187197207217227237247257267277287297307317327337347357367377387397407417427437447457467477487497507517527537547557567577587597607617627637647657667677687697707717727737747757767777787797807817827837847857867877887897907917927937947957967977987998008018028038048058068078088098108118128138148158168178188198208218228238248258268278288298308318328338348358368378388398408418428438448458468478488498508518528538548558568578588598608618628638648658668678688698708718728738748758768778788798808818828838848858868878888898908918928938948958968978988999009019029039049059069079089099109119129139149159169179189199209219229239249259269279289299309319329339349359369379389399409419429439449459469479489499509519529539549559569579589599609619629639649659669679689699709719729739749759769779789799809819829839849859869879889899909919929939949959969979989991000100110021003100410051006100710081009101010111012101310141015101610171018101910201021102210231024102510261027102810291030103110321033103410351036103710381039104010411042104310441045104610471048104910501051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120112111221123112411251126112711281129113011311132113311341135113611371138113911401141114211431144114511461147114811491150115111521153115411551156115711581159116011611162116311641165116611671168116911701171117211731174117511761177117811791180118111821183118411851186118711881189119011911192119311941195119611971198119912001201120212031204120512061207120812091210121112121213121412151216121712181219122012211222122312241225122612271228122912301231123212331234123512361237123812391240124112421243124412451246124712481249125012511252125312541255
  1. from sympy.combinatorics.free_groups import free_group
  2. from sympy.printing.defaults import DefaultPrinting
  3. from itertools import chain, product
  4. from bisect import bisect_left
  5. ###############################################################################
  6. # COSET TABLE #
  7. ###############################################################################
  8. class CosetTable(DefaultPrinting):
  9. # coset_table: Mathematically a coset table
  10. # represented using a list of lists
  11. # alpha: Mathematically a coset (precisely, a live coset)
  12. # represented by an integer between i with 1 <= i <= n
  13. # alpha in c
  14. # x: Mathematically an element of "A" (set of generators and
  15. # their inverses), represented using "FpGroupElement"
  16. # fp_grp: Finitely Presented Group with < X|R > as presentation.
  17. # H: subgroup of fp_grp.
  18. # NOTE: We start with H as being only a list of words in generators
  19. # of "fp_grp". Since `.subgroup` method has not been implemented.
  20. r"""
  21. Properties
  22. ==========
  23. [1] `0 \in \Omega` and `\tau(1) = \epsilon`
  24. [2] `\alpha^x = \beta \Leftrightarrow \beta^{x^{-1}} = \alpha`
  25. [3] If `\alpha^x = \beta`, then `H \tau(\alpha)x = H \tau(\beta)`
  26. [4] `\forall \alpha \in \Omega, 1^{\tau(\alpha)} = \alpha`
  27. References
  28. ==========
  29. .. [1] Holt, D., Eick, B., O'Brien, E.
  30. "Handbook of Computational Group Theory"
  31. .. [2] John J. Cannon; Lucien A. Dimino; George Havas; Jane M. Watson
  32. Mathematics of Computation, Vol. 27, No. 123. (Jul., 1973), pp. 463-490.
  33. "Implementation and Analysis of the Todd-Coxeter Algorithm"
  34. """
  35. # default limit for the number of cosets allowed in a
  36. # coset enumeration.
  37. coset_table_max_limit = 4096000
  38. # limit for the current instance
  39. coset_table_limit = None
  40. # maximum size of deduction stack above or equal to
  41. # which it is emptied
  42. max_stack_size = 100
  43. def __init__(self, fp_grp, subgroup, max_cosets=None):
  44. if not max_cosets:
  45. max_cosets = CosetTable.coset_table_max_limit
  46. self.fp_group = fp_grp
  47. self.subgroup = subgroup
  48. self.coset_table_limit = max_cosets
  49. # "p" is setup independent of Omega and n
  50. self.p = [0]
  51. # a list of the form `[gen_1, gen_1^{-1}, ... , gen_k, gen_k^{-1}]`
  52. self.A = list(chain.from_iterable((gen, gen**-1) \
  53. for gen in self.fp_group.generators))
  54. #P[alpha, x] Only defined when alpha^x is defined.
  55. self.P = [[None]*len(self.A)]
  56. # the mathematical coset table which is a list of lists
  57. self.table = [[None]*len(self.A)]
  58. self.A_dict = {x: self.A.index(x) for x in self.A}
  59. self.A_dict_inv = {}
  60. for x, index in self.A_dict.items():
  61. if index % 2 == 0:
  62. self.A_dict_inv[x] = self.A_dict[x] + 1
  63. else:
  64. self.A_dict_inv[x] = self.A_dict[x] - 1
  65. # used in the coset-table based method of coset enumeration. Each of
  66. # the element is called a "deduction" which is the form (alpha, x) whenever
  67. # a value is assigned to alpha^x during a definition or "deduction process"
  68. self.deduction_stack = []
  69. # Attributes for modified methods.
  70. H = self.subgroup
  71. self._grp = free_group(', ' .join(["a_%d" % i for i in range(len(H))]))[0]
  72. self.P = [[None]*len(self.A)]
  73. self.p_p = {}
  74. @property
  75. def omega(self):
  76. """Set of live cosets. """
  77. return [coset for coset in range(len(self.p)) if self.p[coset] == coset]
  78. def copy(self):
  79. """
  80. Return a shallow copy of Coset Table instance ``self``.
  81. """
  82. self_copy = self.__class__(self.fp_group, self.subgroup)
  83. self_copy.table = [list(perm_rep) for perm_rep in self.table]
  84. self_copy.p = list(self.p)
  85. self_copy.deduction_stack = list(self.deduction_stack)
  86. return self_copy
  87. def __str__(self):
  88. return "Coset Table on %s with %s as subgroup generators" \
  89. % (self.fp_group, self.subgroup)
  90. __repr__ = __str__
  91. @property
  92. def n(self):
  93. """The number `n` represents the length of the sublist containing the
  94. live cosets.
  95. """
  96. if not self.table:
  97. return 0
  98. return max(self.omega) + 1
  99. # Pg. 152 [1]
  100. def is_complete(self):
  101. r"""
  102. The coset table is called complete if it has no undefined entries
  103. on the live cosets; that is, `\alpha^x` is defined for all
  104. `\alpha \in \Omega` and `x \in A`.
  105. """
  106. return not any(None in self.table[coset] for coset in self.omega)
  107. # Pg. 153 [1]
  108. def define(self, alpha, x, modified=False):
  109. r"""
  110. This routine is used in the relator-based strategy of Todd-Coxeter
  111. algorithm if some `\alpha^x` is undefined. We check whether there is
  112. space available for defining a new coset. If there is enough space
  113. then we remedy this by adjoining a new coset `\beta` to `\Omega`
  114. (i.e to set of live cosets) and put that equal to `\alpha^x`, then
  115. make an assignment satisfying Property[1]. If there is not enough space
  116. then we halt the Coset Table creation. The maximum amount of space that
  117. can be used by Coset Table can be manipulated using the class variable
  118. ``CosetTable.coset_table_max_limit``.
  119. See Also
  120. ========
  121. define_c
  122. """
  123. A = self.A
  124. table = self.table
  125. len_table = len(table)
  126. if len_table >= self.coset_table_limit:
  127. # abort the further generation of cosets
  128. raise ValueError("the coset enumeration has defined more than "
  129. "%s cosets. Try with a greater value max number of cosets "
  130. % self.coset_table_limit)
  131. table.append([None]*len(A))
  132. self.P.append([None]*len(self.A))
  133. # beta is the new coset generated
  134. beta = len_table
  135. self.p.append(beta)
  136. table[alpha][self.A_dict[x]] = beta
  137. table[beta][self.A_dict_inv[x]] = alpha
  138. # P[alpha][x] = epsilon, P[beta][x**-1] = epsilon
  139. if modified:
  140. self.P[alpha][self.A_dict[x]] = self._grp.identity
  141. self.P[beta][self.A_dict_inv[x]] = self._grp.identity
  142. self.p_p[beta] = self._grp.identity
  143. def define_c(self, alpha, x):
  144. r"""
  145. A variation of ``define`` routine, described on Pg. 165 [1], used in
  146. the coset table-based strategy of Todd-Coxeter algorithm. It differs
  147. from ``define`` routine in that for each definition it also adds the
  148. tuple `(\alpha, x)` to the deduction stack.
  149. See Also
  150. ========
  151. define
  152. """
  153. A = self.A
  154. table = self.table
  155. len_table = len(table)
  156. if len_table >= self.coset_table_limit:
  157. # abort the further generation of cosets
  158. raise ValueError("the coset enumeration has defined more than "
  159. "%s cosets. Try with a greater value max number of cosets "
  160. % self.coset_table_limit)
  161. table.append([None]*len(A))
  162. # beta is the new coset generated
  163. beta = len_table
  164. self.p.append(beta)
  165. table[alpha][self.A_dict[x]] = beta
  166. table[beta][self.A_dict_inv[x]] = alpha
  167. # append to deduction stack
  168. self.deduction_stack.append((alpha, x))
  169. def scan_c(self, alpha, word):
  170. """
  171. A variation of ``scan`` routine, described on pg. 165 of [1], which
  172. puts at tuple, whenever a deduction occurs, to deduction stack.
  173. See Also
  174. ========
  175. scan, scan_check, scan_and_fill, scan_and_fill_c
  176. """
  177. # alpha is an integer representing a "coset"
  178. # since scanning can be in two cases
  179. # 1. for alpha=0 and w in Y (i.e generating set of H)
  180. # 2. alpha in Omega (set of live cosets), w in R (relators)
  181. A_dict = self.A_dict
  182. A_dict_inv = self.A_dict_inv
  183. table = self.table
  184. f = alpha
  185. i = 0
  186. r = len(word)
  187. b = alpha
  188. j = r - 1
  189. # list of union of generators and their inverses
  190. while i <= j and table[f][A_dict[word[i]]] is not None:
  191. f = table[f][A_dict[word[i]]]
  192. i += 1
  193. if i > j:
  194. if f != b:
  195. self.coincidence_c(f, b)
  196. return
  197. while j >= i and table[b][A_dict_inv[word[j]]] is not None:
  198. b = table[b][A_dict_inv[word[j]]]
  199. j -= 1
  200. if j < i:
  201. # we have an incorrect completed scan with coincidence f ~ b
  202. # run the "coincidence" routine
  203. self.coincidence_c(f, b)
  204. elif j == i:
  205. # deduction process
  206. table[f][A_dict[word[i]]] = b
  207. table[b][A_dict_inv[word[i]]] = f
  208. self.deduction_stack.append((f, word[i]))
  209. # otherwise scan is incomplete and yields no information
  210. # alpha, beta coincide, i.e. alpha, beta represent the pair of cosets where
  211. # coincidence occurs
  212. def coincidence_c(self, alpha, beta):
  213. """
  214. A variation of ``coincidence`` routine used in the coset-table based
  215. method of coset enumeration. The only difference being on addition of
  216. a new coset in coset table(i.e new coset introduction), then it is
  217. appended to ``deduction_stack``.
  218. See Also
  219. ========
  220. coincidence
  221. """
  222. A_dict = self.A_dict
  223. A_dict_inv = self.A_dict_inv
  224. table = self.table
  225. # behaves as a queue
  226. q = []
  227. self.merge(alpha, beta, q)
  228. while len(q) > 0:
  229. gamma = q.pop(0)
  230. for x in A_dict:
  231. delta = table[gamma][A_dict[x]]
  232. if delta is not None:
  233. table[delta][A_dict_inv[x]] = None
  234. # only line of difference from ``coincidence`` routine
  235. self.deduction_stack.append((delta, x**-1))
  236. mu = self.rep(gamma)
  237. nu = self.rep(delta)
  238. if table[mu][A_dict[x]] is not None:
  239. self.merge(nu, table[mu][A_dict[x]], q)
  240. elif table[nu][A_dict_inv[x]] is not None:
  241. self.merge(mu, table[nu][A_dict_inv[x]], q)
  242. else:
  243. table[mu][A_dict[x]] = nu
  244. table[nu][A_dict_inv[x]] = mu
  245. def scan(self, alpha, word, y=None, fill=False, modified=False):
  246. r"""
  247. ``scan`` performs a scanning process on the input ``word``.
  248. It first locates the largest prefix ``s`` of ``word`` for which
  249. `\alpha^s` is defined (i.e is not ``None``), ``s`` may be empty. Let
  250. ``word=sv``, let ``t`` be the longest suffix of ``v`` for which
  251. `\alpha^{t^{-1}}` is defined, and let ``v=ut``. Then three
  252. possibilities are there:
  253. 1. If ``t=v``, then we say that the scan completes, and if, in addition
  254. `\alpha^s = \alpha^{t^{-1}}`, then we say that the scan completes
  255. correctly.
  256. 2. It can also happen that scan does not complete, but `|u|=1`; that
  257. is, the word ``u`` consists of a single generator `x \in A`. In that
  258. case, if `\alpha^s = \beta` and `\alpha^{t^{-1}} = \gamma`, then we can
  259. set `\beta^x = \gamma` and `\gamma^{x^{-1}} = \beta`. These assignments
  260. are known as deductions and enable the scan to complete correctly.
  261. 3. See ``coicidence`` routine for explanation of third condition.
  262. Notes
  263. =====
  264. The code for the procedure of scanning `\alpha \in \Omega`
  265. under `w \in A*` is defined on pg. 155 [1]
  266. See Also
  267. ========
  268. scan_c, scan_check, scan_and_fill, scan_and_fill_c
  269. Scan and Fill
  270. =============
  271. Performed when the default argument fill=True.
  272. Modified Scan
  273. =============
  274. Performed when the default argument modified=True
  275. """
  276. # alpha is an integer representing a "coset"
  277. # since scanning can be in two cases
  278. # 1. for alpha=0 and w in Y (i.e generating set of H)
  279. # 2. alpha in Omega (set of live cosets), w in R (relators)
  280. A_dict = self.A_dict
  281. A_dict_inv = self.A_dict_inv
  282. table = self.table
  283. f = alpha
  284. i = 0
  285. r = len(word)
  286. b = alpha
  287. j = r - 1
  288. b_p = y
  289. if modified:
  290. f_p = self._grp.identity
  291. flag = 0
  292. while fill or flag == 0:
  293. flag = 1
  294. while i <= j and table[f][A_dict[word[i]]] is not None:
  295. if modified:
  296. f_p = f_p*self.P[f][A_dict[word[i]]]
  297. f = table[f][A_dict[word[i]]]
  298. i += 1
  299. if i > j:
  300. if f != b:
  301. if modified:
  302. self.modified_coincidence(f, b, f_p**-1*y)
  303. else:
  304. self.coincidence(f, b)
  305. return
  306. while j >= i and table[b][A_dict_inv[word[j]]] is not None:
  307. if modified:
  308. b_p = b_p*self.P[b][self.A_dict_inv[word[j]]]
  309. b = table[b][A_dict_inv[word[j]]]
  310. j -= 1
  311. if j < i:
  312. # we have an incorrect completed scan with coincidence f ~ b
  313. # run the "coincidence" routine
  314. if modified:
  315. self.modified_coincidence(f, b, f_p**-1*b_p)
  316. else:
  317. self.coincidence(f, b)
  318. elif j == i:
  319. # deduction process
  320. table[f][A_dict[word[i]]] = b
  321. table[b][A_dict_inv[word[i]]] = f
  322. if modified:
  323. self.P[f][self.A_dict[word[i]]] = f_p**-1*b_p
  324. self.P[b][self.A_dict_inv[word[i]]] = b_p**-1*f_p
  325. return
  326. elif fill:
  327. self.define(f, word[i], modified=modified)
  328. # otherwise scan is incomplete and yields no information
  329. # used in the low-index subgroups algorithm
  330. def scan_check(self, alpha, word):
  331. r"""
  332. Another version of ``scan`` routine, described on, it checks whether
  333. `\alpha` scans correctly under `word`, it is a straightforward
  334. modification of ``scan``. ``scan_check`` returns ``False`` (rather than
  335. calling ``coincidence``) if the scan completes incorrectly; otherwise
  336. it returns ``True``.
  337. See Also
  338. ========
  339. scan, scan_c, scan_and_fill, scan_and_fill_c
  340. """
  341. # alpha is an integer representing a "coset"
  342. # since scanning can be in two cases
  343. # 1. for alpha=0 and w in Y (i.e generating set of H)
  344. # 2. alpha in Omega (set of live cosets), w in R (relators)
  345. A_dict = self.A_dict
  346. A_dict_inv = self.A_dict_inv
  347. table = self.table
  348. f = alpha
  349. i = 0
  350. r = len(word)
  351. b = alpha
  352. j = r - 1
  353. while i <= j and table[f][A_dict[word[i]]] is not None:
  354. f = table[f][A_dict[word[i]]]
  355. i += 1
  356. if i > j:
  357. return f == b
  358. while j >= i and table[b][A_dict_inv[word[j]]] is not None:
  359. b = table[b][A_dict_inv[word[j]]]
  360. j -= 1
  361. if j < i:
  362. # we have an incorrect completed scan with coincidence f ~ b
  363. # return False, instead of calling coincidence routine
  364. return False
  365. elif j == i:
  366. # deduction process
  367. table[f][A_dict[word[i]]] = b
  368. table[b][A_dict_inv[word[i]]] = f
  369. return True
  370. def merge(self, k, lamda, q, w=None, modified=False):
  371. """
  372. Merge two classes with representatives ``k`` and ``lamda``, described
  373. on Pg. 157 [1] (for pseudocode), start by putting ``p[k] = lamda``.
  374. It is more efficient to choose the new representative from the larger
  375. of the two classes being merged, i.e larger among ``k`` and ``lamda``.
  376. procedure ``merge`` performs the merging operation, adds the deleted
  377. class representative to the queue ``q``.
  378. Parameters
  379. ==========
  380. 'k', 'lamda' being the two class representatives to be merged.
  381. Notes
  382. =====
  383. Pg. 86-87 [1] contains a description of this method.
  384. See Also
  385. ========
  386. coincidence, rep
  387. """
  388. p = self.p
  389. rep = self.rep
  390. phi = rep(k, modified=modified)
  391. psi = rep(lamda, modified=modified)
  392. if phi != psi:
  393. mu = min(phi, psi)
  394. v = max(phi, psi)
  395. p[v] = mu
  396. if modified:
  397. if v == phi:
  398. self.p_p[phi] = self.p_p[k]**-1*w*self.p_p[lamda]
  399. else:
  400. self.p_p[psi] = self.p_p[lamda]**-1*w**-1*self.p_p[k]
  401. q.append(v)
  402. def rep(self, k, modified=False):
  403. r"""
  404. Parameters
  405. ==========
  406. `k \in [0 \ldots n-1]`, as for ``self`` only array ``p`` is used
  407. Returns
  408. =======
  409. Representative of the class containing ``k``.
  410. Returns the representative of `\sim` class containing ``k``, it also
  411. makes some modification to array ``p`` of ``self`` to ease further
  412. computations, described on Pg. 157 [1].
  413. The information on classes under `\sim` is stored in array `p` of
  414. ``self`` argument, which will always satisfy the property:
  415. `p[\alpha] \sim \alpha` and `p[\alpha]=\alpha \iff \alpha=rep(\alpha)`
  416. `\forall \in [0 \ldots n-1]`.
  417. So, for `\alpha \in [0 \ldots n-1]`, we find `rep(self, \alpha)` by
  418. continually replacing `\alpha` by `p[\alpha]` until it becomes
  419. constant (i.e satisfies `p[\alpha] = \alpha`):w
  420. To increase the efficiency of later ``rep`` calculations, whenever we
  421. find `rep(self, \alpha)=\beta`, we set
  422. `p[\gamma] = \beta \forall \gamma \in p-chain` from `\alpha` to `\beta`
  423. Notes
  424. =====
  425. ``rep`` routine is also described on Pg. 85-87 [1] in Atkinson's
  426. algorithm, this results from the fact that ``coincidence`` routine
  427. introduces functionality similar to that introduced by the
  428. ``minimal_block`` routine on Pg. 85-87 [1].
  429. See Also
  430. ========
  431. coincidence, merge
  432. """
  433. p = self.p
  434. lamda = k
  435. rho = p[lamda]
  436. if modified:
  437. s = p[:]
  438. while rho != lamda:
  439. if modified:
  440. s[rho] = lamda
  441. lamda = rho
  442. rho = p[lamda]
  443. if modified:
  444. rho = s[lamda]
  445. while rho != k:
  446. mu = rho
  447. rho = s[mu]
  448. p[rho] = lamda
  449. self.p_p[rho] = self.p_p[rho]*self.p_p[mu]
  450. else:
  451. mu = k
  452. rho = p[mu]
  453. while rho != lamda:
  454. p[mu] = lamda
  455. mu = rho
  456. rho = p[mu]
  457. return lamda
  458. # alpha, beta coincide, i.e. alpha, beta represent the pair of cosets
  459. # where coincidence occurs
  460. def coincidence(self, alpha, beta, w=None, modified=False):
  461. r"""
  462. The third situation described in ``scan`` routine is handled by this
  463. routine, described on Pg. 156-161 [1].
  464. The unfortunate situation when the scan completes but not correctly,
  465. then ``coincidence`` routine is run. i.e when for some `i` with
  466. `1 \le i \le r+1`, we have `w=st` with `s = x_1 x_2 \dots x_{i-1}`,
  467. `t = x_i x_{i+1} \dots x_r`, and `\beta = \alpha^s` and
  468. `\gamma = \alpha^{t-1}` are defined but unequal. This means that
  469. `\beta` and `\gamma` represent the same coset of `H` in `G`. Described
  470. on Pg. 156 [1]. ``rep``
  471. See Also
  472. ========
  473. scan
  474. """
  475. A_dict = self.A_dict
  476. A_dict_inv = self.A_dict_inv
  477. table = self.table
  478. # behaves as a queue
  479. q = []
  480. if modified:
  481. self.modified_merge(alpha, beta, w, q)
  482. else:
  483. self.merge(alpha, beta, q)
  484. while len(q) > 0:
  485. gamma = q.pop(0)
  486. for x in A_dict:
  487. delta = table[gamma][A_dict[x]]
  488. if delta is not None:
  489. table[delta][A_dict_inv[x]] = None
  490. mu = self.rep(gamma, modified=modified)
  491. nu = self.rep(delta, modified=modified)
  492. if table[mu][A_dict[x]] is not None:
  493. if modified:
  494. v = self.p_p[delta]**-1*self.P[gamma][self.A_dict[x]]**-1
  495. v = v*self.p_p[gamma]*self.P[mu][self.A_dict[x]]
  496. self.modified_merge(nu, table[mu][self.A_dict[x]], v, q)
  497. else:
  498. self.merge(nu, table[mu][A_dict[x]], q)
  499. elif table[nu][A_dict_inv[x]] is not None:
  500. if modified:
  501. v = self.p_p[gamma]**-1*self.P[gamma][self.A_dict[x]]
  502. v = v*self.p_p[delta]*self.P[mu][self.A_dict_inv[x]]
  503. self.modified_merge(mu, table[nu][self.A_dict_inv[x]], v, q)
  504. else:
  505. self.merge(mu, table[nu][A_dict_inv[x]], q)
  506. else:
  507. table[mu][A_dict[x]] = nu
  508. table[nu][A_dict_inv[x]] = mu
  509. if modified:
  510. v = self.p_p[gamma]**-1*self.P[gamma][self.A_dict[x]]*self.p_p[delta]
  511. self.P[mu][self.A_dict[x]] = v
  512. self.P[nu][self.A_dict_inv[x]] = v**-1
  513. # method used in the HLT strategy
  514. def scan_and_fill(self, alpha, word):
  515. """
  516. A modified version of ``scan`` routine used in the relator-based
  517. method of coset enumeration, described on pg. 162-163 [1], which
  518. follows the idea that whenever the procedure is called and the scan
  519. is incomplete then it makes new definitions to enable the scan to
  520. complete; i.e it fills in the gaps in the scan of the relator or
  521. subgroup generator.
  522. """
  523. self.scan(alpha, word, fill=True)
  524. def scan_and_fill_c(self, alpha, word):
  525. """
  526. A modified version of ``scan`` routine, described on Pg. 165 second
  527. para. [1], with modification similar to that of ``scan_anf_fill`` the
  528. only difference being it calls the coincidence procedure used in the
  529. coset-table based method i.e. the routine ``coincidence_c`` is used.
  530. See Also
  531. ========
  532. scan, scan_and_fill
  533. """
  534. A_dict = self.A_dict
  535. A_dict_inv = self.A_dict_inv
  536. table = self.table
  537. r = len(word)
  538. f = alpha
  539. i = 0
  540. b = alpha
  541. j = r - 1
  542. # loop until it has filled the alpha row in the table.
  543. while True:
  544. # do the forward scanning
  545. while i <= j and table[f][A_dict[word[i]]] is not None:
  546. f = table[f][A_dict[word[i]]]
  547. i += 1
  548. if i > j:
  549. if f != b:
  550. self.coincidence_c(f, b)
  551. return
  552. # forward scan was incomplete, scan backwards
  553. while j >= i and table[b][A_dict_inv[word[j]]] is not None:
  554. b = table[b][A_dict_inv[word[j]]]
  555. j -= 1
  556. if j < i:
  557. self.coincidence_c(f, b)
  558. elif j == i:
  559. table[f][A_dict[word[i]]] = b
  560. table[b][A_dict_inv[word[i]]] = f
  561. self.deduction_stack.append((f, word[i]))
  562. else:
  563. self.define_c(f, word[i])
  564. # method used in the HLT strategy
  565. def look_ahead(self):
  566. """
  567. When combined with the HLT method this is known as HLT+Lookahead
  568. method of coset enumeration, described on pg. 164 [1]. Whenever
  569. ``define`` aborts due to lack of space available this procedure is
  570. executed. This routine helps in recovering space resulting from
  571. "coincidence" of cosets.
  572. """
  573. R = self.fp_group.relators
  574. p = self.p
  575. # complete scan all relators under all cosets(obviously live)
  576. # without making new definitions
  577. for beta in self.omega:
  578. for w in R:
  579. self.scan(beta, w)
  580. if p[beta] < beta:
  581. break
  582. # Pg. 166
  583. def process_deductions(self, R_c_x, R_c_x_inv):
  584. """
  585. Processes the deductions that have been pushed onto ``deduction_stack``,
  586. described on Pg. 166 [1] and is used in coset-table based enumeration.
  587. See Also
  588. ========
  589. deduction_stack
  590. """
  591. p = self.p
  592. table = self.table
  593. while len(self.deduction_stack) > 0:
  594. if len(self.deduction_stack) >= CosetTable.max_stack_size:
  595. self.look_ahead()
  596. del self.deduction_stack[:]
  597. continue
  598. else:
  599. alpha, x = self.deduction_stack.pop()
  600. if p[alpha] == alpha:
  601. for w in R_c_x:
  602. self.scan_c(alpha, w)
  603. if p[alpha] < alpha:
  604. break
  605. beta = table[alpha][self.A_dict[x]]
  606. if beta is not None and p[beta] == beta:
  607. for w in R_c_x_inv:
  608. self.scan_c(beta, w)
  609. if p[beta] < beta:
  610. break
  611. def process_deductions_check(self, R_c_x, R_c_x_inv):
  612. """
  613. A variation of ``process_deductions``, this calls ``scan_check``
  614. wherever ``process_deductions`` calls ``scan``, described on Pg. [1].
  615. See Also
  616. ========
  617. process_deductions
  618. """
  619. table = self.table
  620. while len(self.deduction_stack) > 0:
  621. alpha, x = self.deduction_stack.pop()
  622. for w in R_c_x:
  623. if not self.scan_check(alpha, w):
  624. return False
  625. beta = table[alpha][self.A_dict[x]]
  626. if beta is not None:
  627. for w in R_c_x_inv:
  628. if not self.scan_check(beta, w):
  629. return False
  630. return True
  631. def switch(self, beta, gamma):
  632. r"""Switch the elements `\beta, \gamma \in \Omega` of ``self``, used
  633. by the ``standardize`` procedure, described on Pg. 167 [1].
  634. See Also
  635. ========
  636. standardize
  637. """
  638. A = self.A
  639. A_dict = self.A_dict
  640. table = self.table
  641. for x in A:
  642. z = table[gamma][A_dict[x]]
  643. table[gamma][A_dict[x]] = table[beta][A_dict[x]]
  644. table[beta][A_dict[x]] = z
  645. for alpha in range(len(self.p)):
  646. if self.p[alpha] == alpha:
  647. if table[alpha][A_dict[x]] == beta:
  648. table[alpha][A_dict[x]] = gamma
  649. elif table[alpha][A_dict[x]] == gamma:
  650. table[alpha][A_dict[x]] = beta
  651. def standardize(self):
  652. r"""
  653. A coset table is standardized if when running through the cosets and
  654. within each coset through the generator images (ignoring generator
  655. inverses), the cosets appear in order of the integers
  656. `0, 1, \dots, n`. "Standardize" reorders the elements of `\Omega`
  657. such that, if we scan the coset table first by elements of `\Omega`
  658. and then by elements of A, then the cosets occur in ascending order.
  659. ``standardize()`` is used at the end of an enumeration to permute the
  660. cosets so that they occur in some sort of standard order.
  661. Notes
  662. =====
  663. procedure is described on pg. 167-168 [1], it also makes use of the
  664. ``switch`` routine to replace by smaller integer value.
  665. Examples
  666. ========
  667. >>> from sympy.combinatorics import free_group
  668. >>> from sympy.combinatorics.fp_groups import FpGroup, coset_enumeration_r
  669. >>> F, x, y = free_group("x, y")
  670. # Example 5.3 from [1]
  671. >>> f = FpGroup(F, [x**2*y**2, x**3*y**5])
  672. >>> C = coset_enumeration_r(f, [])
  673. >>> C.compress()
  674. >>> C.table
  675. [[1, 3, 1, 3], [2, 0, 2, 0], [3, 1, 3, 1], [0, 2, 0, 2]]
  676. >>> C.standardize()
  677. >>> C.table
  678. [[1, 2, 1, 2], [3, 0, 3, 0], [0, 3, 0, 3], [2, 1, 2, 1]]
  679. """
  680. A = self.A
  681. A_dict = self.A_dict
  682. gamma = 1
  683. for alpha, x in product(range(self.n), A):
  684. beta = self.table[alpha][A_dict[x]]
  685. if beta >= gamma:
  686. if beta > gamma:
  687. self.switch(gamma, beta)
  688. gamma += 1
  689. if gamma == self.n:
  690. return
  691. # Compression of a Coset Table
  692. def compress(self):
  693. """Removes the non-live cosets from the coset table, described on
  694. pg. 167 [1].
  695. """
  696. gamma = -1
  697. A = self.A
  698. A_dict = self.A_dict
  699. A_dict_inv = self.A_dict_inv
  700. table = self.table
  701. chi = tuple([i for i in range(len(self.p)) if self.p[i] != i])
  702. for alpha in self.omega:
  703. gamma += 1
  704. if gamma != alpha:
  705. # replace alpha by gamma in coset table
  706. for x in A:
  707. beta = table[alpha][A_dict[x]]
  708. table[gamma][A_dict[x]] = beta
  709. table[beta][A_dict_inv[x]] == gamma
  710. # all the cosets in the table are live cosets
  711. self.p = list(range(gamma + 1))
  712. # delete the useless columns
  713. del table[len(self.p):]
  714. # re-define values
  715. for row in table:
  716. for j in range(len(self.A)):
  717. row[j] -= bisect_left(chi, row[j])
  718. def conjugates(self, R):
  719. R_c = list(chain.from_iterable((rel.cyclic_conjugates(), \
  720. (rel**-1).cyclic_conjugates()) for rel in R))
  721. R_set = set()
  722. for conjugate in R_c:
  723. R_set = R_set.union(conjugate)
  724. R_c_list = []
  725. for x in self.A:
  726. r = {word for word in R_set if word[0] == x}
  727. R_c_list.append(r)
  728. R_set.difference_update(r)
  729. return R_c_list
  730. def coset_representative(self, coset):
  731. '''
  732. Compute the coset representative of a given coset.
  733. Examples
  734. ========
  735. >>> from sympy.combinatorics import free_group
  736. >>> from sympy.combinatorics.fp_groups import FpGroup, coset_enumeration_r
  737. >>> F, x, y = free_group("x, y")
  738. >>> f = FpGroup(F, [x**3, y**3, x**-1*y**-1*x*y])
  739. >>> C = coset_enumeration_r(f, [x])
  740. >>> C.compress()
  741. >>> C.table
  742. [[0, 0, 1, 2], [1, 1, 2, 0], [2, 2, 0, 1]]
  743. >>> C.coset_representative(0)
  744. <identity>
  745. >>> C.coset_representative(1)
  746. y
  747. >>> C.coset_representative(2)
  748. y**-1
  749. '''
  750. for x in self.A:
  751. gamma = self.table[coset][self.A_dict[x]]
  752. if coset == 0:
  753. return self.fp_group.identity
  754. if gamma < coset:
  755. return self.coset_representative(gamma)*x**-1
  756. ##############################
  757. # Modified Methods #
  758. ##############################
  759. def modified_define(self, alpha, x):
  760. r"""
  761. Define a function p_p from from [1..n] to A* as
  762. an additional component of the modified coset table.
  763. Parameters
  764. ==========
  765. \alpha \in \Omega
  766. x \in A*
  767. See Also
  768. ========
  769. define
  770. """
  771. self.define(alpha, x, modified=True)
  772. def modified_scan(self, alpha, w, y, fill=False):
  773. r"""
  774. Parameters
  775. ==========
  776. \alpha \in \Omega
  777. w \in A*
  778. y \in (YUY^-1)
  779. fill -- `modified_scan_and_fill` when set to True.
  780. See Also
  781. ========
  782. scan
  783. """
  784. self.scan(alpha, w, y=y, fill=fill, modified=True)
  785. def modified_scan_and_fill(self, alpha, w, y):
  786. self.modified_scan(alpha, w, y, fill=True)
  787. def modified_merge(self, k, lamda, w, q):
  788. r"""
  789. Parameters
  790. ==========
  791. 'k', 'lamda' -- the two class representatives to be merged.
  792. q -- queue of length l of elements to be deleted from `\Omega` *.
  793. w -- Word in (YUY^-1)
  794. See Also
  795. ========
  796. merge
  797. """
  798. self.merge(k, lamda, q, w=w, modified=True)
  799. def modified_rep(self, k):
  800. r"""
  801. Parameters
  802. ==========
  803. `k \in [0 \ldots n-1]`
  804. See Also
  805. ========
  806. rep
  807. """
  808. self.rep(k, modified=True)
  809. def modified_coincidence(self, alpha, beta, w):
  810. r"""
  811. Parameters
  812. ==========
  813. A coincident pair `\alpha, \beta \in \Omega, w \in Y \cup Y^{-1}`
  814. See Also
  815. ========
  816. coincidence
  817. """
  818. self.coincidence(alpha, beta, w=w, modified=True)
  819. ###############################################################################
  820. # COSET ENUMERATION #
  821. ###############################################################################
  822. # relator-based method
  823. def coset_enumeration_r(fp_grp, Y, max_cosets=None, draft=None,
  824. incomplete=False, modified=False):
  825. """
  826. This is easier of the two implemented methods of coset enumeration.
  827. and is often called the HLT method, after Hazelgrove, Leech, Trotter
  828. The idea is that we make use of ``scan_and_fill`` makes new definitions
  829. whenever the scan is incomplete to enable the scan to complete; this way
  830. we fill in the gaps in the scan of the relator or subgroup generator,
  831. that's why the name relator-based method.
  832. An instance of `CosetTable` for `fp_grp` can be passed as the keyword
  833. argument `draft` in which case the coset enumeration will start with
  834. that instance and attempt to complete it.
  835. When `incomplete` is `True` and the function is unable to complete for
  836. some reason, the partially complete table will be returned.
  837. # TODO: complete the docstring
  838. See Also
  839. ========
  840. scan_and_fill,
  841. Examples
  842. ========
  843. >>> from sympy.combinatorics.free_groups import free_group
  844. >>> from sympy.combinatorics.fp_groups import FpGroup, coset_enumeration_r
  845. >>> F, x, y = free_group("x, y")
  846. # Example 5.1 from [1]
  847. >>> f = FpGroup(F, [x**3, y**3, x**-1*y**-1*x*y])
  848. >>> C = coset_enumeration_r(f, [x])
  849. >>> for i in range(len(C.p)):
  850. ... if C.p[i] == i:
  851. ... print(C.table[i])
  852. [0, 0, 1, 2]
  853. [1, 1, 2, 0]
  854. [2, 2, 0, 1]
  855. >>> C.p
  856. [0, 1, 2, 1, 1]
  857. # Example from exercises Q2 [1]
  858. >>> f = FpGroup(F, [x**2*y**2, y**-1*x*y*x**-3])
  859. >>> C = coset_enumeration_r(f, [])
  860. >>> C.compress(); C.standardize()
  861. >>> C.table
  862. [[1, 2, 3, 4],
  863. [5, 0, 6, 7],
  864. [0, 5, 7, 6],
  865. [7, 6, 5, 0],
  866. [6, 7, 0, 5],
  867. [2, 1, 4, 3],
  868. [3, 4, 2, 1],
  869. [4, 3, 1, 2]]
  870. # Example 5.2
  871. >>> f = FpGroup(F, [x**2, y**3, (x*y)**3])
  872. >>> Y = [x*y]
  873. >>> C = coset_enumeration_r(f, Y)
  874. >>> for i in range(len(C.p)):
  875. ... if C.p[i] == i:
  876. ... print(C.table[i])
  877. [1, 1, 2, 1]
  878. [0, 0, 0, 2]
  879. [3, 3, 1, 0]
  880. [2, 2, 3, 3]
  881. # Example 5.3
  882. >>> f = FpGroup(F, [x**2*y**2, x**3*y**5])
  883. >>> Y = []
  884. >>> C = coset_enumeration_r(f, Y)
  885. >>> for i in range(len(C.p)):
  886. ... if C.p[i] == i:
  887. ... print(C.table[i])
  888. [1, 3, 1, 3]
  889. [2, 0, 2, 0]
  890. [3, 1, 3, 1]
  891. [0, 2, 0, 2]
  892. # Example 5.4
  893. >>> F, a, b, c, d, e = free_group("a, b, c, d, e")
  894. >>> f = FpGroup(F, [a*b*c**-1, b*c*d**-1, c*d*e**-1, d*e*a**-1, e*a*b**-1])
  895. >>> Y = [a]
  896. >>> C = coset_enumeration_r(f, Y)
  897. >>> for i in range(len(C.p)):
  898. ... if C.p[i] == i:
  899. ... print(C.table[i])
  900. [0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
  901. # example of "compress" method
  902. >>> C.compress()
  903. >>> C.table
  904. [[0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]
  905. # Exercises Pg. 161, Q2.
  906. >>> F, x, y = free_group("x, y")
  907. >>> f = FpGroup(F, [x**2*y**2, y**-1*x*y*x**-3])
  908. >>> Y = []
  909. >>> C = coset_enumeration_r(f, Y)
  910. >>> C.compress()
  911. >>> C.standardize()
  912. >>> C.table
  913. [[1, 2, 3, 4],
  914. [5, 0, 6, 7],
  915. [0, 5, 7, 6],
  916. [7, 6, 5, 0],
  917. [6, 7, 0, 5],
  918. [2, 1, 4, 3],
  919. [3, 4, 2, 1],
  920. [4, 3, 1, 2]]
  921. # John J. Cannon; Lucien A. Dimino; George Havas; Jane M. Watson
  922. # Mathematics of Computation, Vol. 27, No. 123. (Jul., 1973), pp. 463-490
  923. # from 1973chwd.pdf
  924. # Table 1. Ex. 1
  925. >>> F, r, s, t = free_group("r, s, t")
  926. >>> E1 = FpGroup(F, [t**-1*r*t*r**-2, r**-1*s*r*s**-2, s**-1*t*s*t**-2])
  927. >>> C = coset_enumeration_r(E1, [r])
  928. >>> for i in range(len(C.p)):
  929. ... if C.p[i] == i:
  930. ... print(C.table[i])
  931. [0, 0, 0, 0, 0, 0]
  932. Ex. 2
  933. >>> F, a, b = free_group("a, b")
  934. >>> Cox = FpGroup(F, [a**6, b**6, (a*b)**2, (a**2*b**2)**2, (a**3*b**3)**5])
  935. >>> C = coset_enumeration_r(Cox, [a])
  936. >>> index = 0
  937. >>> for i in range(len(C.p)):
  938. ... if C.p[i] == i:
  939. ... index += 1
  940. >>> index
  941. 500
  942. # Ex. 3
  943. >>> F, a, b = free_group("a, b")
  944. >>> B_2_4 = FpGroup(F, [a**4, b**4, (a*b)**4, (a**-1*b)**4, (a**2*b)**4, \
  945. (a*b**2)**4, (a**2*b**2)**4, (a**-1*b*a*b)**4, (a*b**-1*a*b)**4])
  946. >>> C = coset_enumeration_r(B_2_4, [a])
  947. >>> index = 0
  948. >>> for i in range(len(C.p)):
  949. ... if C.p[i] == i:
  950. ... index += 1
  951. >>> index
  952. 1024
  953. References
  954. ==========
  955. .. [1] Holt, D., Eick, B., O'Brien, E.
  956. "Handbook of computational group theory"
  957. """
  958. # 1. Initialize a coset table C for < X|R >
  959. C = CosetTable(fp_grp, Y, max_cosets=max_cosets)
  960. # Define coset table methods.
  961. if modified:
  962. _scan_and_fill = C.modified_scan_and_fill
  963. _define = C.modified_define
  964. else:
  965. _scan_and_fill = C.scan_and_fill
  966. _define = C.define
  967. if draft:
  968. C.table = draft.table[:]
  969. C.p = draft.p[:]
  970. R = fp_grp.relators
  971. A_dict = C.A_dict
  972. p = C.p
  973. for i in range(len(Y)):
  974. if modified:
  975. _scan_and_fill(0, Y[i], C._grp.generators[i])
  976. else:
  977. _scan_and_fill(0, Y[i])
  978. alpha = 0
  979. while alpha < C.n:
  980. if p[alpha] == alpha:
  981. try:
  982. for w in R:
  983. if modified:
  984. _scan_and_fill(alpha, w, C._grp.identity)
  985. else:
  986. _scan_and_fill(alpha, w)
  987. # if alpha was eliminated during the scan then break
  988. if p[alpha] < alpha:
  989. break
  990. if p[alpha] == alpha:
  991. for x in A_dict:
  992. if C.table[alpha][A_dict[x]] is None:
  993. _define(alpha, x)
  994. except ValueError as e:
  995. if incomplete:
  996. return C
  997. raise e
  998. alpha += 1
  999. return C
  1000. def modified_coset_enumeration_r(fp_grp, Y, max_cosets=None, draft=None,
  1001. incomplete=False):
  1002. r"""
  1003. Introduce a new set of symbols y \in Y that correspond to the
  1004. generators of the subgroup. Store the elements of Y as a
  1005. word P[\alpha, x] and compute the coset table similar to that of
  1006. the regular coset enumeration methods.
  1007. Examples
  1008. ========
  1009. >>> from sympy.combinatorics.free_groups import free_group
  1010. >>> from sympy.combinatorics.fp_groups import FpGroup
  1011. >>> from sympy.combinatorics.coset_table import modified_coset_enumeration_r
  1012. >>> F, x, y = free_group("x, y")
  1013. >>> f = FpGroup(F, [x**3, y**3, x**-1*y**-1*x*y])
  1014. >>> C = modified_coset_enumeration_r(f, [x])
  1015. >>> C.table
  1016. [[0, 0, 1, 2], [1, 1, 2, 0], [2, 2, 0, 1], [None, 1, None, None], [1, 3, None, None]]
  1017. See Also
  1018. ========
  1019. coset_enumertation_r
  1020. References
  1021. ==========
  1022. .. [1] Holt, D., Eick, B., O'Brien, E.,
  1023. "Handbook of Computational Group Theory",
  1024. Section 5.3.2
  1025. """
  1026. return coset_enumeration_r(fp_grp, Y, max_cosets=max_cosets, draft=draft,
  1027. incomplete=incomplete, modified=True)
  1028. # Pg. 166
  1029. # coset-table based method
  1030. def coset_enumeration_c(fp_grp, Y, max_cosets=None, draft=None,
  1031. incomplete=False):
  1032. """
  1033. >>> from sympy.combinatorics.free_groups import free_group
  1034. >>> from sympy.combinatorics.fp_groups import FpGroup, coset_enumeration_c
  1035. >>> F, x, y = free_group("x, y")
  1036. >>> f = FpGroup(F, [x**3, y**3, x**-1*y**-1*x*y])
  1037. >>> C = coset_enumeration_c(f, [x])
  1038. >>> C.table
  1039. [[0, 0, 1, 2], [1, 1, 2, 0], [2, 2, 0, 1]]
  1040. """
  1041. # Initialize a coset table C for < X|R >
  1042. X = fp_grp.generators
  1043. R = fp_grp.relators
  1044. C = CosetTable(fp_grp, Y, max_cosets=max_cosets)
  1045. if draft:
  1046. C.table = draft.table[:]
  1047. C.p = draft.p[:]
  1048. C.deduction_stack = draft.deduction_stack
  1049. for alpha, x in product(range(len(C.table)), X):
  1050. if C.table[alpha][C.A_dict[x]] is not None:
  1051. C.deduction_stack.append((alpha, x))
  1052. A = C.A
  1053. # replace all the elements by cyclic reductions
  1054. R_cyc_red = [rel.identity_cyclic_reduction() for rel in R]
  1055. R_c = list(chain.from_iterable((rel.cyclic_conjugates(), (rel**-1).cyclic_conjugates()) \
  1056. for rel in R_cyc_red))
  1057. R_set = set()
  1058. for conjugate in R_c:
  1059. R_set = R_set.union(conjugate)
  1060. # a list of subsets of R_c whose words start with "x".
  1061. R_c_list = []
  1062. for x in C.A:
  1063. r = {word for word in R_set if word[0] == x}
  1064. R_c_list.append(r)
  1065. R_set.difference_update(r)
  1066. for w in Y:
  1067. C.scan_and_fill_c(0, w)
  1068. for x in A:
  1069. C.process_deductions(R_c_list[C.A_dict[x]], R_c_list[C.A_dict_inv[x]])
  1070. alpha = 0
  1071. while alpha < len(C.table):
  1072. if C.p[alpha] == alpha:
  1073. try:
  1074. for x in C.A:
  1075. if C.p[alpha] != alpha:
  1076. break
  1077. if C.table[alpha][C.A_dict[x]] is None:
  1078. C.define_c(alpha, x)
  1079. C.process_deductions(R_c_list[C.A_dict[x]], R_c_list[C.A_dict_inv[x]])
  1080. except ValueError as e:
  1081. if incomplete:
  1082. return C
  1083. raise e
  1084. alpha += 1
  1085. return C