Repository: sympy/sympy

Issue:
fresnels: write as integral missing pi
@arnabnandikgp I noticed something about the fresnels rewrite:
```
In [12]: S = fresnels(x)

In [13]: S.diff(x)
Out[13]: 
   ⎛   2⎞
   ⎜π⋅x ⎟
sin⎜────⎟
   ⎝ 2  ⎠

In [14]: S.rewrite(Integral)
Out[14]: 
x           
⌠           
⎮    ⎛ 2⎞   
⎮ sin⎝t ⎠ dt
⌡           
0           

In [15]: _.diff(x)
Out[15]: 
   ⎛ 2⎞
sin⎝x ⎠
```


Discussion:

Patch:
diff --git a/.mailmap b/.mailmap
index 909f44816b..20b4d5efc5 100644
--- a/.mailmap
+++ b/.mailmap
@@ -984,7 +984,6 @@ Mikhail Remnev <maremnev@gmail.com> Mikhail Remnev <141655736+maremnev@users.nor
 Milan Jolly <milan.cs16@iitp.ac.in> mijo2 <milan.cs16@iitp.ac.in>
 Min Ragan-Kelley <benjaminrk@gmail.com>
 Miro Hrončok <miro@hroncok.cz>
-Mohamed Rezk <mohrizq895@gmail.com>
 Mohammad Sadeq Dousti <msdousti@gmail.com>
 Mohammed Bilal <r.mohammedbilal@gmail.com>
 Mohit Balwani <mohitbalwani.ict17@gmail.com> Mohit Balwani <44258119+Mohitbalwani26@users.noreply.github.com>
diff --git a/doc/src/modules/physics/vector/advanced.rst b/doc/src/modules/physics/vector/advanced.rst
index 316487f64a..4f0f922d62 100644
--- a/doc/src/modules/physics/vector/advanced.rst
+++ b/doc/src/modules/physics/vector/advanced.rst
@@ -198,34 +198,3 @@ so dynamic symbols created before or after will print the same way.
 Also note that ``Vector``'s ``.dt`` method uses the ``._t`` attribute of
 ``dynamicsymbols``, along with a number of other important functions and
 methods. Don't mix and match symbols representing time.
-
-Solving Vector Equations
-========================
-
-To solve equations involving vectors, you cannot directly use the solve
-functions on a vector. Instead, you must convert the vector to a set of scalar
-equations.
-
-Suppose that we have two frames ``N`` and ``A``, where ``A`` is rotated 30
-degrees about the z-axis with respect to ``N``. ::
-
-  >>> from sympy import pi, symbols, solve
-  >>> from sympy.physics.vector import ReferenceFrame
-  >>> N = ReferenceFrame("N")
-  >>> A = ReferenceFrame("A")
-  >>> A.orient_axis(N, pi / 6, N.z)
-
-Suppose that we have two vectors ``v1`` and ``v2``, which represent the same
-vector using different symbols. ::
-
-  >>> v1x, v1y, v1z = symbols("v1x v1y v1z")
-  >>> v2x, v2y, v2z = symbols("v2x v2y v2z")
-  >>> v1 = v1x * N.x + v1y * N.y + v1z * N.z
-  >>> v2 = v2x * A.x + v2y * A.y + v2z * A.z
-
-Our goal is to find the relationship between the symbols used in ``v2`` and the
-symbols used in ``v1``. We can achieve this by converting the vector to a matrix
-and then solving the matrix using :meth:`sympy.solvers.solvers.solve`. ::
-
-  >>> solve((v1 - v2).to_matrix(N), [v2x, v2y, v2z])
-  {v2x: sqrt(3)*v1x/2 + v1y/2, v2y: -v1x/2 + sqrt(3)*v1y/2, v2z: v1z}
diff --git a/sympy/functions/special/error_functions.py b/sympy/functions/special/error_functions.py
index 292ca7089b..eeaa3ddb5e 100644
--- a/sympy/functions/special/error_functions.py
+++ b/sympy/functions/special/error_functions.py
@@ -1795,6 +1795,7 @@ def _eval_rewrite_as_uppergamma(self, z, **kwargs):
 
     def _eval_nseries(self, x, n, logx, cdir=0):
         # NOTE this is fairly inefficient
+        n += 1
         if self.args[0].subs(x, 0) != 0:
             return super()._eval_nseries(x, n, logx)
         baseseries = self._trigfunc(x)._eval_nseries(x, n, logx)
@@ -1927,11 +1928,11 @@ def _eval_aseries(self, n, args0, x, logx):
         # Expansion at oo
         if point is S.Infinity:
             z = self.args[0]
-            p = [S.NegativeOne**k * factorial(2*k) / z**(2*k + 1)
-                    for k in range(n//2 + 1)] + [Order(1/z**n, x)]
-            q = [S.NegativeOne**k * factorial(2*k + 1) / z**(2*(k + 1))
-                    for k in range(n//2)] + [Order(1/z**n, x)]
-            return pi/2 - cos(z)*Add(*p) - sin(z)*Add(*q)
+            p = [S.NegativeOne**k * factorial(2*k) / z**(2*k)
+                    for k in range(int((n - 1)/2))] + [Order(1/z**n, x)]
+            q = [S.NegativeOne**k * factorial(2*k + 1) / z**(2*k + 1)
+                    for k in range(int(n/2) - 1)] + [Order(1/z**n, x)]
+            return pi/2 - (cos(z)/z)*Add(*p) - (sin(z)/z)*Add(*q)
 
         # All other points are not handled
         return super(Si, self)._eval_aseries(n, args0, x, logx)
@@ -2070,11 +2071,11 @@ def _eval_aseries(self, n, args0, x, logx):
 
         if point in (S.Infinity, S.NegativeInfinity):
             z = self.args[0]
-            p = [S.NegativeOne**k * factorial(2*k) / z**(2*k + 1)
-                    for k in range(n//2 + 1)] + [Order(1/z**n, x)]
-            q = [S.NegativeOne**k * factorial(2*k + 1) / z**(2*(k + 1))
-                    for k in range(n//2)] + [Order(1/z**n, x)]
-            result = sin(z)*(Add(*p)) - cos(z)*(Add(*q))
+            p = [S.NegativeOne**k * factorial(2*k) / z**(2*k)
+                    for k in range(int((n - 1)/2))] + [Order(1/z**n, x)]
+            q = [S.NegativeOne**k * factorial(2*k + 1) / z**(2*k + 1)
+                    for k in range(int(n/2) - 1)] + [Order(1/z**n, x)]
+            result = (sin(z)/z)*Add(*p) - (cos(z)/z)*Add(*q)
 
             if point is S.NegativeInfinity:
                 result += I*pi
@@ -2475,7 +2476,7 @@ def _eval_rewrite_as_meijerg(self, z, **kwargs):
     def _eval_rewrite_as_Integral(self, z, **kwargs):
         from sympy.integrals.integrals import Integral
         t = Dummy(uniquely_named_symbol('t', [z]).name)
-        return Integral(sin(t**2), (t, 0, z))
+        return Integral(sin(pi*t**2/2), (t, 0, z))
 
     def _eval_as_leading_term(self, x, logx=None, cdir=0):
         from sympy.series.order import Order
@@ -2636,7 +2637,7 @@ def _eval_rewrite_as_meijerg(self, z, **kwargs):
     def _eval_rewrite_as_Integral(self, z, **kwargs):
         from sympy.integrals.integrals import Integral
         t = Dummy(uniquely_named_symbol('t', [z]).name)
-        return Integral(cos(t**2), (t, 0, z))
+        return Integral(cos(pi*t**2/2), (t, 0, z))
 
     def _eval_as_leading_term(self, x, logx=None, cdir=0):
         from sympy.series.order import Order
diff --git a/sympy/physics/control/lti.py b/sympy/physics/control/lti.py
index 0c14574138..2524da6d21 100644
--- a/sympy/physics/control/lti.py
+++ b/sympy/physics/control/lti.py
@@ -618,8 +618,8 @@ def __new__(cls, num, den, var):
         if den == 0:
             raise ValueError("TransferFunction cannot have a zero denominator.")
 
-        if (((isinstance(num, (Expr, TransferFunction, Series, Parallel)) and num.has(Symbol)) or num.is_number) and
-            ((isinstance(den, (Expr, TransferFunction, Series, Parallel)) and den.has(Symbol)) or den.is_number)):
+        if (((isinstance(num, Expr) and num.has(Symbol)) or num.is_number) and
+            ((isinstance(den, Expr) and den.has(Symbol)) or den.is_number)):
             return super(TransferFunction, cls).__new__(cls, num, den, var)
 
         else:
@@ -2205,7 +2205,7 @@ def __neg__(self):
         return MIMOParallel(*arg_list)
 
 
-class Feedback(TransferFunction):
+class Feedback(SISOLinearTimeInvariant):
     r"""
     A class for representing closed-loop feedback interconnection between two
     SISO input/output systems.
@@ -2296,8 +2296,8 @@ def __new__(cls, sys1, sys2=None, sign=-1):
         if not sys2:
             sys2 = TransferFunction(1, 1, sys1.var)
 
-        if not (isinstance(sys1, (TransferFunction, Series, Feedback))
-            and isinstance(sys2, (TransferFunction, Series, Feedback))):
+        if not (isinstance(sys1, (TransferFunction, Series))
+            and isinstance(sys2, (TransferFunction, Series))):
             raise TypeError("Unsupported type for `sys1` or `sys2` of Feedback.")
 
         if sign not in [-1, 1]:
@@ -2314,7 +2314,7 @@ def __new__(cls, sys1, sys2=None, sign=-1):
                 Both `sys1` and `sys2` should be using the
                 same complex variable."""))
 
-        return super(TransferFunction, cls).__new__(cls, sys1, sys2, _sympify(sign))
+        return super().__new__(cls, sys1, sys2, _sympify(sign))
 
     @property
     def sys1(self):
@@ -2400,24 +2400,6 @@ def sign(self):
         """
         return self.args[2]
 
-    @property
-    def num(self):
-        """
-        Returns the numerator of the closed loop feedback system.
-        """
-        return self.sys1
-
-    @property
-    def den(self):
-        """
-        Returns the denominator of the closed loop feedback model.
-        """
-        unit = TransferFunction(1, 1, self.var)
-        arg_list = list(self.sys1.args) if isinstance(self.sys1, Series) else [self.sys1]
-        if self.sign == 1:
-            return Parallel(unit, -Series(self.sys2, *arg_list))
-        return Parallel(unit, Series(self.sys2, *arg_list))
-
     @property
     def sensitivity(self):
         """
@@ -2497,27 +2479,6 @@ def doit(self, cancel=False, expand=False, **hints):
     def _eval_rewrite_as_TransferFunction(self, num, den, sign, **kwargs):
         return self.doit()
 
-    def to_expr(self):
-        """
-        Converts a ``Feedback`` object to SymPy Expr.
-
-        Examples
-        ========
-
-        >>> from sympy.abc import s, a, b
-        >>> from sympy.physics.control.lti import TransferFunction, Feedback
-        >>> from sympy import Expr
-        >>> tf1 = TransferFunction(a+s, 1, s)
-        >>> tf2 = TransferFunction(b+s, 1, s)
-        >>> fd1 = Feedback(tf1, tf2)
-        >>> fd1.to_expr()
-        (a + s)/((a + s)*(b + s) + 1)
-        >>> isinstance(_, Expr)
-        True
-        """
-
-        return self.doit().to_expr()
-
     def __neg__(self):
         return Feedback(-self.sys1, -self.sys2, self.sign)
 
diff --git a/sympy/physics/optics/utils.py b/sympy/physics/optics/utils.py
index 881f1cdfc0..4f1f612b21 100644
--- a/sympy/physics/optics/utils.py
+++ b/sympy/physics/optics/utils.py
@@ -671,8 +671,8 @@ def hyperfocal_distance(f, N, c):
 def transverse_magnification(si, so):
     """
 
-    Calculates the transverse magnification upon reflection in a mirror,
-    which is the ratio of the image size to the object size.
+    Calculates the transverse magnification, which is the ratio of the
+    image size to the object size.
 
     Parameters
     ==========
diff --git a/sympy/polys/matrices/sdm.py b/sympy/polys/matrices/sdm.py
index 3fc3643cde..a465fe66b7 100644
--- a/sympy/polys/matrices/sdm.py
+++ b/sympy/polys/matrices/sdm.py
@@ -1805,13 +1805,6 @@ def sdm_rref_den(A, K):
         Aij = Ai[j]
         return ({0: Ai.copy()}, Aij, [j])
 
-    # For inexact domains like RR[x] we use quo and discard the remainder.
-    # Maybe it would be better for K.exquo to do this automatically.
-    if K.is_Exact:
-        exquo = K.exquo
-    else:
-        exquo = K.quo
-
     # Make sure we have the rows in order to make this deterministic from the
     # outset.
     _, rows_in_order = zip(*sorted(A.items()))
@@ -1908,7 +1901,7 @@ def sdm_rref_den(A, K):
                 for l, Akl in Ak.items():
                     Akl = Akl * Aij
                     if divisor is not None:
-                        Akl = exquo(Akl, divisor)
+                        Akl = K.exquo(Akl, divisor)
                     Ak[l] = Akl
                 continue
 
@@ -1919,19 +1912,19 @@ def sdm_rref_den(A, K):
             for l in Ai_nz - Ak_nz:
                 Ak[l] = - Akj * Ai[l]
                 if divisor is not None:
-                    Ak[l] = exquo(Ak[l], divisor)
+                    Ak[l] = K.exquo(Ak[l], divisor)
 
             # This loop also not needed in sdm_irref.
             for l in Ak_nz - Ai_nz:
                 Ak[l] = Aij * Ak[l]
                 if divisor is not None:
-                    Ak[l] = exquo(Ak[l], divisor)
+                    Ak[l] = K.exquo(Ak[l], divisor)
 
             for l in Ai_nz & Ak_nz:
                 Akl = Aij * Ak[l] - Akj * Ai[l]
                 if Akl:
                     if divisor is not None:
-                        Akl = exquo(Akl, divisor)
+                        Akl = K.exquo(Akl, divisor)
                     Ak[l] = Akl
                 else:
                     Ak.pop(l)
@@ -1955,7 +1948,7 @@ def sdm_rref_den(A, K):
                 denom *= Aij
 
         if divisor is not None:
-            denom = exquo(denom, divisor)
+            denom = K.exquo(denom, divisor)
 
         # Update the divisor.
         divisor = denom
diff --git a/sympy/printing/latex.py b/sympy/printing/latex.py
index 7b70e294e6..3e7f779120 100644
--- a/sympy/printing/latex.py
+++ b/sympy/printing/latex.py
@@ -695,8 +695,6 @@ def _helper_print_standard_power(self, expr, template: str) -> str:
         base = self.parenthesize(expr.base, PRECEDENCE['Pow'])
         if expr.base.is_Symbol:
             base = self.parenthesize_super(base)
-        elif expr.base.is_Float:
-            base = r"{%s}" % base
         elif (isinstance(expr.base, Derivative)
             and base.startswith(r'\left(')
             and re.match(r'\\left\(\\d?d?dot', base)
diff --git a/sympy/solvers/diophantine/diophantine.py b/sympy/solvers/diophantine/diophantine.py
index a1198eec84..cba48dabc7 100644
--- a/sympy/solvers/diophantine/diophantine.py
+++ b/sympy/solvers/diophantine/diophantine.py
@@ -1,5 +1,3 @@
-from __future__ import annotations
-
 from sympy.core.add import Add
 from sympy.core.assumptions import check_assumptions
 from sympy.core.containers import Tuple
@@ -2176,26 +2174,44 @@ def _special_diop_DN(D, N):
     #
     # assert (1 < N**2 < D) and (not integer_nthroot(D, 2)[1])
 
-    sqrt_D = isqrt(D)
-    F = {N // f**2: f for f in divisors(square_factor(abs(N)), generator=True)}
+    sqrt_D = sqrt(D)
+    F = [(N, 1)]
+    f = 2
+    while True:
+        f2 = f**2
+        if f2 > abs(N):
+            break
+        n, r = divmod(N, f2)
+        if r == 0:
+            F.append((n, f))
+        f += 1
+
     P = 0
     Q = 1
     G0, G1 = 0, 1
     B0, B1 = 1, 0
 
     solutions = []
+
+    i = 0
     while True:
-        for _ in range(2):
-            a = (P + sqrt_D) // Q
-            P = a*Q - P
-            Q = (D - P**2) // Q
-            G0, G1 = G1, a*G1 + G0
-            B0, B1 = B1, a*B1 + B0
-            if (s := G1**2 - D*B1**2) in F:
-                f = F[s]
-                solutions.append((f*G1, f*B1))
-        if Q == 1:
+        a = floor((P + sqrt_D) / Q)
+        P = a*Q - P
+        Q = (D - P**2) // Q
+        G2 = a*G1 + G0
+        B2 = a*B1 + B0
+
+        for n, f in F:
+            if G2**2 - D*B2**2 == n:
+                solutions.append((f*G2, f*B2))
+
+        i += 1
+        if Q == 1 and i % 2 == 0:
             break
+
+        G0, G1 = G1, G2
+        B0, B1 = B1, B2
+
     return solutions
 
 
@@ -3247,45 +3263,15 @@ def descent(A, B):
     return _remove_gcd(x_0*x_1 + A*z_0*z_1, z_0*x_1 + x_0*z_1, t_1*t_2*y_1)
 
 
-def gaussian_reduce(w:int, a:int, b:int) -> tuple[int, int]:
+def gaussian_reduce(w, a, b):
     r"""
     Returns a reduced solution `(x, z)` to the congruence
-    `X^2 - aZ^2 \equiv 0 \pmod{b}` so that `x^2 + |a|z^2` is minimal.
-    Here ``w`` is a solution of the congruence `x^2 \equiv a \pmod{b}`.
-
-    Explanation
-    ===========
+    `X^2 - aZ^2 \equiv 0 \ (mod \ b)` so that `x^2 + |a|z^2` is minimal.
 
-    We define a special dot product of the vectors `u = (u_{1}, u_{2})` and
-    `v = (v_{1}, v_{2})` which is defined in order to reduce solution of
-    the congruence equation `X^2 - aZ^2 \equiv 0 \pmod{b}`.
-
-    .. math ::
-        u \cdot v := (wu_1 + bu_2)(wv_1 + bv_2) + |a|u_1v_1
-
-    Here, if `u_3` is defined as `wu_1 + bu_2`, `u \cdot v` can be computed in `u_3v_3 + |a|u_1v_1`.
-    Therefore, it is kept in pairs of ``(u_3, u_1)``.
-
-    Parameters
-    ==========
-
-    w : int
-        ``w`` s.t. `w^2 \equiv a \pmod{b}`
-    a : int
-    b : int
-
-    Examples
-    ========
+    Details
+    =======
 
-    >>> from sympy.solvers.diophantine.diophantine import gaussian_reduce
-    >>> from sympy.ntheory.residue_ntheory import sqrt_mod
-    >>> a, b = 19, 101
-    >>> gaussian_reduce(sqrt_mod(a, b), a, b) # 1**2 - 19*(-4)**2 = -303
-    (1, -4)
-    >>> a, b = 11, 14
-    >>> x, z = gaussian_reduce(sqrt_mod(a, b), a, b)
-    >>> (x**2 - a*z**2) % b == 0
-    True
+    Here ``w`` is a solution of the congruence `x^2 \equiv a \ (mod \ b)`
 
     References
     ==========
@@ -3296,25 +3282,48 @@ def gaussian_reduce(w:int, a:int, b:int) -> tuple[int, int]:
            Mathematics of Computation, 72(243), 1417-1441.
            https://doi.org/10.1090/S0025-5718-02-01480-1
     """
-    a = abs(a)
-    def _dot(u, v):
-        return u[0]*v[0] + a*u[1]*v[1]
+    u = (0, 1)
+    v = (1, 0)
 
-    u = (b, 0)
-    v = (w, 1) if b*w >= 0 else (-w, -1)
-    # i.e., _dot(u, v) >= 0
+    if dot(u, v, w, a, b) < 0:
+        v = (-v[0], -v[1])
 
-    if b**2 < w**2 + a:
+    if norm(u, w, a, b) < norm(v, w, a, b):
         u, v = v, u
-    # i.e., norm(u) >= norm(v), where norm(u) := sqrt(_dot(u, u))
-
-    while _dot(u, u) > (dv := _dot(v, v)):
-        k = _dot(u, v) // dv
-        u, v = v, (u[0] - k*v[0], u[1] - k*v[1])
-    c = (v[0] - u[0], v[1] - u[1])
-    if _dot(c, c) <= _dot(u, u) <= 2*_dot(u, v):
-        return c
-    return u
+
+    while norm(u, w, a, b) > norm(v, w, a, b):
+        k = dot(u, v, w, a, b) // dot(v, v, w, a, b)
+        u, v = v, (u[0]- k*v[0], u[1]- k*v[1])
+
+    u, v = v, u
+
+    if dot(u, v, w, a, b) < dot(v, v, w, a, b)/2 or norm((u[0]-v[0], u[1]-v[1]), w, a, b) > norm(v, w, a, b):
+        c = v
+    else:
+        c = (u[0] - v[0], u[1] - v[1])
+
+    return c[0]*w + b*c[1], c[0]
+
+
+def dot(u, v, w, a, b):
+    r"""
+    Returns a special dot product of the vectors `u = (u_{1}, u_{2})` and
+    `v = (v_{1}, v_{2})` which is defined in order to reduce solution of
+    the congruence equation `X^2 - aZ^2 \equiv 0 \ (mod \ b)`.
+    """
+    u_1, u_2 = u
+    v_1, v_2 = v
+    return (w*u_1 + b*u_2)*(w*v_1 + b*v_2) + abs(a)*u_1*v_1
+
+
+def norm(u, w, a, b):
+    r"""
+    Returns the norm of the vector `u = (u_{1}, u_{2})` under the dot product
+    defined by `u \cdot v = (wu_{1} + bu_{2})(w*v_{1} + bv_{2}) + |a|*u_{1}*v_{1}`
+    where `u = (u_{1}, u_{2})` and `v = (v_{1}, v_{2})`.
+    """
+    u_1, u_2 = u
+    return sqrt(dot((u_1, u_2), (u_1, u_2), w, a, b))
 
 
 def holzer(x, y, z, a, b, c):


Test patch:
diff --git a/.github/workflows/runtests.yml b/.github/workflows/runtests.yml
index 1552e84a60..fe46ad79c5 100644
--- a/.github/workflows/runtests.yml
+++ b/.github/workflows/runtests.yml
@@ -29,12 +29,12 @@ jobs:
       - name: Basic code quality tests
         run: bin/test quality
 
-      - name: Run Ruff on the sympy package
-        run: ruff check .
-
       - name: Run flake8 on the sympy package
         run: flake8 sympy
 
+      - name: Run Ruff on the sympy package
+        run: ruff check .
+
       - name: Detect invalid escapes like '\e'
         run: python -We:invalid -We::SyntaxWarning -m compileall -f -q sympy/
 
diff --git a/sympy/functions/special/tests/test_error_functions.py b/sympy/functions/special/tests/test_error_functions.py
index 8fb6645761..774ee183cb 100644
--- a/sympy/functions/special/tests/test_error_functions.py
+++ b/sympy/functions/special/tests/test_error_functions.py
@@ -621,15 +621,16 @@ def test_si():
         Si(1/x)._eval_as_leading_term(x, cdir=-1) == Si(1/x)
 
     assert Si(x).nseries(x, n=8) == \
-        x - x**3/18 + x**5/600 - x**7/35280 + O(x**8)
+        x - x**3/18 + x**5/600 - x**7/35280 + O(x**9)
     assert Shi(x).nseries(x, n=8) == \
-        x + x**3/18 + x**5/600 + x**7/35280 + O(x**8)
-    assert Si(sin(x)).nseries(x, n=5) == x - 2*x**3/9 + O(x**5)
+        x + x**3/18 + x**5/600 + x**7/35280 + O(x**9)
+    assert Si(sin(x)).nseries(x, n=5) == x - 2*x**3/9 + 17*x**5/450 + O(x**6)
     assert Si(x).nseries(x, 1, n=3) == \
         Si(1) + (x - 1)*sin(1) + (x - 1)**2*(-sin(1)/2 + cos(1)/2) + O((x - 1)**3, (x, 1))
 
-    assert Si(x).series(x, oo) == -sin(x)*(-6/x**4 + x**(-2) + O(x**(-6), (x, oo))) - \
-        cos(x)*(24/x**5 - 2/x**3 + 1/x + O(x**(-6), (x, oo))) + pi/2
+    assert Si(x).series(x, oo) == pi/2 - (- 6/x**3 + 1/x \
+        + O(x**(-7), (x, oo)))*sin(x)/x - (24/x**4 - 2/x**2 + 1 \
+        + O(x**(-7), (x, oo)))*cos(x)/x
 
     t = Symbol('t', Dummy=True)
     assert Si(x).rewrite(sinc).dummy_eq(Integral(sinc(t), (t, 0, x)))
@@ -674,15 +675,16 @@ def test_ci():
     assert tn_arg(Chi)
 
     assert Ci(x).nseries(x, n=4) == \
-        EulerGamma + log(x) - x**2/4 + O(x**4)
+        EulerGamma + log(x) - x**2/4 + x**4/96 + O(x**5)
     assert Chi(x).nseries(x, n=4) == \
-        EulerGamma + log(x) + x**2/4 + O(x**4)
+        EulerGamma + log(x) + x**2/4 + x**4/96 + O(x**5)
 
-    assert Ci(x).series(x, oo) == -cos(x)*(-6/x**4 + x**(-2) + O(x**(-6), (x, oo))) + \
-        sin(x)*(24/x**5 - 2/x**3 + 1/x + O(x**(-6), (x, oo)))
+    assert Ci(x).series(x, oo) == -cos(x)*(-6/x**3 + 1/x \
+        + O(x**(-7), (x, oo)))/x + (24/x**4 - 2/x**2 + 1 \
+        + O(x**(-7), (x, oo)))*sin(x)/x
 
-    assert Ci(x).series(x, -oo) == -cos(x)*(-6/x**4 + x**(-2) + O(x**(-6), (x, -oo))) + \
-        sin(x)*(24/x**5 - 2/x**3 + 1/x + O(x**(-6), (x, -oo))) + I*pi
+    assert Ci(x).series(x, -oo) == -cos(x)*(-6/x**3 + 1/x + O(-1/x**7, (x, -oo)))/x + \
+            (24/x**4 - 2/x**2 + 1 + O(-1/x**7, (x, -oo)))*sin(x)/x + I*pi
 
     assert limit(log(x) - Ci(2*x), x, 0) == -log(2) - EulerGamma
     assert Ci(x).rewrite(uppergamma) == -expint(1, x*exp_polar(-I*pi/2))/2 -\
@@ -826,10 +828,12 @@ def test_fresnel_series():
     assert ((3*fresnelc(2*z)).series(z, oo) - 3*fc.subs(z, 2*z)).expand().is_Order
 
 
-def test_integral_rewrites(): #issues 26134, 26144
+def test_integral_rewrites(): #issues 26134, 26144, 26306
     assert expint(n, x).rewrite(Integral).dummy_eq(Integral(t**-n * exp(-t*x), (t, 1, oo)))
     assert Si(x).rewrite(Integral).dummy_eq(Integral(sinc(t), (t, 0, x)))
     assert Ci(x).rewrite(Integral).dummy_eq(log(x) - Integral((1 - cos(t))/t, (t, 0, x)) + EulerGamma)
-    assert fresnels(x).rewrite(Integral).dummy_eq(Integral(sin(t**2), (t, 0, x)))
-    assert fresnelc(x).rewrite(Integral).dummy_eq(Integral(cos(t**2), (t, 0, x)))
+    assert fresnels(x).rewrite(Integral).dummy_eq(Integral(sin(pi*t**2/2), (t, 0, x)))
+    assert fresnelc(x).rewrite(Integral).dummy_eq(Integral(cos(pi*t**2/2), (t, 0, x)))
     assert Ei(x).rewrite(Integral).dummy_eq(Integral(exp(t)/t, (t, -oo, x)))
+    assert fresnels(x).diff(x) == fresnels(x).rewrite(Integral).diff(x)
+    assert fresnelc(x).diff(x) == fresnelc(x).rewrite(Integral).diff(x)
diff --git a/sympy/physics/control/tests/test_lti.py b/sympy/physics/control/tests/test_lti.py
index 9a3f599cbe..a13fe4d20e 100644
--- a/sympy/physics/control/tests/test_lti.py
+++ b/sympy/physics/control/tests/test_lti.py
@@ -19,7 +19,7 @@
     StateSpace, gbt, bilinear, forward_diff, backward_diff, phase_margin, gain_margin)
 from sympy.testing.pytest import raises
 
-a, x, b, c, s, g, d, p, k, tau, zeta, wn, T = symbols('a, x, b, c, s, g, d, p, k,\
+a, x, b, s, g, d, p, k, tau, zeta, wn, T = symbols('a, x, b, s, g, d, p, k,\
     tau, zeta, wn, T')
 a0, a1, a2, a3, b0, b1, b2, b3, c0, c1, c2, c3, d0, d1, d2, d3 = symbols('a0:4,\
     b0:4, c0:4, d0:4')
@@ -1006,75 +1006,6 @@ def test_Feedback_functions():
         TransferFunction(p, a0*p + p + p**a1 - s, p)
 
 
-def test_Feedback_as_TransferFunction():
-    # Solves issue https://github.com/sympy/sympy/issues/26161
-    tf1 = TransferFunction(s+1, 1, s)
-    tf2 = TransferFunction(s+2, 1, s)
-    fd1 = Feedback(tf1, tf2, -1) # Negative Feedback system
-    fd2 = Feedback(tf1, tf2, 1) # Positive Feedback system
-    unit = TransferFunction(1, 1, s)
-
-    # Checking the type
-    assert isinstance(fd1, TransferFunction)
-    assert isinstance(fd1, Feedback)
-
-    # Testing the numerator and denominator
-    assert fd1.num == tf1
-    assert fd2.num == tf1
-    assert fd1.den == Parallel(unit, Series(tf2, tf1))
-    assert fd2.den == Parallel(unit, -Series(tf2, tf1))
-
-    # Testing the Series and Parallel Combination with Feedback and TransferFunction
-    s1 = Series(tf1, fd1)
-    p1 = Parallel(tf1, fd1)
-    assert tf1 * fd1 == s1
-    assert tf1 + fd1 == p1
-    assert s1.doit() == TransferFunction((s + 1)**2, (s + 1)*(s + 2) + 1, s)
-    assert p1.doit() == TransferFunction(s + (s + 1)*((s + 1)*(s + 2) + 1) + 1, (s + 1)*(s + 2) + 1, s)
-
-    # Testing the use of Feedback and TransferFunction with Feedback
-    fd3 = Feedback(tf1*fd1, tf2, -1)
-    assert fd3 == Feedback(Series(tf1, fd1), tf2)
-    assert fd3.num == tf1 * fd1
-    assert fd3.den == Parallel(unit, Series(tf2, Series(tf1, fd1)))
-
-    # Testing the use of Feedback and TransferFunction with TransferFunction
-    tf3 = TransferFunction(tf1*fd1, tf2, s)
-    assert tf3 == TransferFunction(Series(tf1, fd1), tf2, s)
-    assert tf3.num == tf1*fd1
-
-def test_issue_26161():
-    # Issue https://github.com/sympy/sympy/issues/26161
-    Ib, Is, m, h, l2, l1 = symbols('I_b, I_s, m, h, l2, l1',
-                                            real=True, nonnegative=True)
-    KD, KP, v = symbols('K_D, K_P, v', real=True)
-
-    tau1_sq = (Ib + m * h ** 2) / m / g / h
-    tau2 = l2 / v
-    tau3 = v / (l1 + l2)
-    K = v ** 2 / g / (l1 + l2)
-
-    Gtheta = TransferFunction(-K * (tau2 * s + 1), tau1_sq * s ** 2 - 1, s)
-    Gdelta = TransferFunction(1, Is * s ** 2 + c * s, s)
-    Gpsi = TransferFunction(1, tau3 * s, s)
-    Dcont = TransferFunction(KD * s, 1, s)
-    PIcont = TransferFunction(KP, s, s)
-    Gunity = TransferFunction(1, 1, s)
-
-    Ginner = Feedback(Dcont * Gdelta, Gtheta)
-    Gouter = Feedback(PIcont * Ginner * Gpsi, Gunity)
-    assert Gouter == Feedback(Series(PIcont, Series(Ginner, Gpsi)), Gunity)
-    assert Gouter.num == Series(PIcont, Series(Ginner, Gpsi))
-    assert Gouter.den == Parallel(Gunity, Series(Gunity, Series(PIcont, Series(Ginner, Gpsi))))
-    expr = (KD*KP*g*s**3*v**2*(l1 + l2)*(Is*s**2 + c*s)**2*(-g*h*m + s**2*(Ib + h**2*m))*(-KD*g*h*m*s*v**2*(l2*s + v) + \
-            g*v*(l1 + l2)*(Is*s**2 + c*s)*(-g*h*m + s**2*(Ib + h**2*m))))/((s**2*v*(Is*s**2 + c*s)*(-KD*g*h*m*s*v**2* \
-            (l2*s + v) + g*v*(l1 + l2)*(Is*s**2 + c*s)*(-g*h*m + s**2*(Ib + h**2*m)))*(KD*KP*g*s*v*(l1 + l2)**2* \
-            (Is*s**2 + c*s)*(-g*h*m + s**2*(Ib + h**2*m)) + s**2*v*(Is*s**2 + c*s)*(-KD*g*h*m*s*v**2*(l2*s + v) + \
-            g*v*(l1 + l2)*(Is*s**2 + c*s)*(-g*h*m + s**2*(Ib + h**2*m))))/(l1 + l2)))
-
-    assert (Gouter.to_expr() - expr).simplify() == 0
-
-
 def test_MIMOFeedback_construction():
     tf1 = TransferFunction(1, s, s)
     tf2 = TransferFunction(s, s**3 - 1, s)
diff --git a/sympy/polys/matrices/tests/test_inverse.py b/sympy/polys/matrices/tests/test_inverse.py
index 47c8279932..9b3415469e 100644
--- a/sympy/polys/matrices/tests/test_inverse.py
+++ b/sympy/polys/matrices/tests/test_inverse.py
@@ -1,4 +1,4 @@
-from sympy import ZZ, Matrix
+from sympy import ZZ
 from sympy.polys.matrices import DM, DomainMatrix
 from sympy.polys.matrices.dense import ddm_iinv
 from sympy.polys.matrices.exceptions import DMNonInvertibleMatrixError
@@ -6,9 +6,6 @@
 
 import pytest
 from sympy.testing.pytest import raises
-from sympy.core.numbers import all_close
-
-from sympy.abc import x
 
 
 # Examples are given as adjugate matrix and determinant adj_det should match
@@ -154,40 +151,3 @@ def test_Matrix_adjugate(name, A, A_inv, den):
 @pytest.mark.parametrize('name, A, A_inv, den', INVERSE_EXAMPLES)
 def test_dm_adj_det(name, A, A_inv, den):
     assert A.adj_det() == (A_inv, den)
-
-
-def test_inverse_inexact():
-
-    M = Matrix([[x-0.3, -0.06, -0.22],
-                [-0.46, x-0.48, -0.41],
-                [-0.14, -0.39, x-0.64]])
-
-    Mn = Matrix([[1.0*x**2 - 1.12*x + 0.1473, 0.06*x + 0.0474, 0.22*x - 0.081],
-                 [0.46*x - 0.237, 1.0*x**2 - 0.94*x + 0.1612, 0.41*x - 0.0218],
-                 [0.14*x + 0.1122, 0.39*x - 0.1086, 1.0*x**2 - 0.78*x + 0.1164]])
-
-    d = 1.0*x**3 - 1.42*x**2 + 0.4249*x - 0.0546540000000002
-
-    Mi = Mn / d
-
-    M_dm = M.to_DM()
-    M_dmd = M_dm.to_dense()
-    M_dm_num, M_dm_den = M_dm.inv_den()
-    M_dmd_num, M_dmd_den = M_dmd.inv_den()
-
-    # XXX: We don't check M_dm().to_field().inv() which currently uses division
-    # and produces a more complicate result from gcd cancellation failing.
-    # DomainMatrix.inv() over RR(x) should be changed to clear denominators and
-    # use DomainMatrix.inv_den().
-
-    Minvs = [
-        M.inv(),
-        (M_dm_num.to_field() / M_dm_den).to_Matrix(),
-        (M_dmd_num.to_field() / M_dmd_den).to_Matrix(),
-        M_dm_num.to_Matrix() / M_dm_den.as_expr(),
-        M_dmd_num.to_Matrix() / M_dmd_den.as_expr(),
-    ]
-
-    for Minv in Minvs:
-        for Mi1, Mi2 in zip(Minv.flat(), Mi.flat()):
-            assert all_close(Mi2, Mi1)
diff --git a/sympy/printing/tests/test_latex.py b/sympy/printing/tests/test_latex.py
index c6e771e9f4..eaeef06816 100644
--- a/sympy/printing/tests/test_latex.py
+++ b/sympy/printing/tests/test_latex.py
@@ -2453,8 +2453,6 @@ def test_Pow():
     assert latex(x**(Rational(-1, 3))) == r'\frac{1}{\sqrt[3]{x}}'
     x2 = Symbol(r'x^2')
     assert latex(x2**2) == r'\left(x^{2}\right)^{2}'
-    # Issue 11011
-    assert latex(S('1.453e4500')**x) == r'{1.453 \cdot 10^{4500}}^{x}'
 
 
 def test_issue_7180():