Fix order computation of linear groups

src/sage/groups/matrix_gps/linear.py

@@ -68,6 +68,8 @@
 from sage.misc.latex import latex
 from sage.misc.misc_c import prod
 from sage.rings.infinity import Infinity
+from sage.rings.integer_ring import ZZ
+from sage.rings.finite_rings.integer_mod_ring import Integers
 
 
 ###############################################################################
@@ -311,7 +313,7 @@ def _check_matrix(self, x, *args):
                 raise TypeError('matrix must non-zero determinant')
 
     def order(self):
-        """
+        r"""
         Return the order of ``self``.
 
         EXAMPLES::
@@ -320,6 +322,23 @@ def order(self):
             sage: G.order()
             372000
 
+        The order computation also works over the base rings `\ZZ/n\ZZ`::
+
+            sage: GL(4, Integers(15)).order()
+            2815842631680000000
+
+            sage: SL(4, Integers(15)).order()
+            351980328960000000
+
+        Arbitrary base rings are currently not fully supported::
+
+            sage: R.<x> = PolynomialRing(GF(7))
+            sage: S = R.quotient(x^2 + 5)
+            sage: GL(2, S).order()
+            Traceback (most recent call last):
+            ...
+            NotImplementedError: order computation of linear groups not fully supported for arbitrary base rings
+
         TESTS:
 
         Check if :issue:`36876` is fixed::
@@ -340,18 +359,71 @@ def order(self):
             True
             sage: S.order()
             117600
-        """
-        n = self.degree()
 
-        if self.base_ring().is_finite():
-            q = self.base_ring().order()
+        Check if :issue:`37934` is fixed::
+
+            sage: GL(2, Integers(4)).order()
+            96
+
+            sage: GL(2, Integers(1)).order()
+            1
+
+            sage: GL(1, ZZ).order()
+            2
+        """
+        def order_over_finite_field(q, n, special):
             ord = prod(q**n - q**i for i in range(n))
-            if self._special:
-                return ord / (q-1)
+            if special:
+                ord = ord // (q-1)
             return ord
 
-        if self._special and n == 1:
+        n = self.degree()
+        R = self.base_ring()
+
+        if R.is_finite():
+            q = R.order()
+
+            if q == 1:
+                return 1
+
+            if R.is_field():
+                return order_over_finite_field(q, n, self._special)
+
+            if R == Integers(q):
+                ord = 1
+
+                # By the Chinese remainder theorem we need to build the product
+                # over the orders of GL(n, ZZ/p^e ZZ) (or SL) for all prime
+                # powers in the factorization of q
+                for (p,e) in q.factor():
+                    ord_base = order_over_finite_field(p, n, self._special)
+
+                    if not self._special:
+                        ord *= p**((e-1)*n**2) * ord_base
+
+                    # We apply |SL(n, R)| = |GL(n, R)| / euler_phi(q), but since we
+                    # already iterate over the prime factorization of q, we divide
+                    # out euler_phi(q) iteratively. Noting that (p-1) is already
+                    # handled in the call to order_over_finite_field, we only
+                    # need to remove p^(e-1) compared to the above formula 
+                    else:
+                        ord *= p**((e-1)*(n**2 - 1)) * ord_base
+
+                return ord
+
+            raise NotImplementedError("order computation of linear groups not "
+                                      "fully supported for arbitrary base rings")
+
+        if n > 1 or (R.is_field() and not self._special):
+            return Infinity
+
+        if self._special:
             return 1
-        return Infinity
+
+        if R == ZZ:
+            return 2
+
+        raise NotImplementedError("order computation of linear groups not "
+                                  "fully supported for arbitrary base rings")
 
     cardinality = order