diff --git a/sympy/physics/continuum_mechanics/beam.py b/sympy/physics/continuum_mechanics/beam.py
index edc3b68f5b09..1101ba8d4624 100644
--- a/sympy/physics/continuum_mechanics/beam.py
+++ b/sympy/physics/continuum_mechanics/beam.py
@@ -2001,11 +2001,10 @@ class Beam3D(Beam):
     >>> dy = collect(simplify(dy), x)
     >>> dx == dz == 0
     True
-    >>> dy == (x*(12*A*E*G*I*l**3*q - 24*A*E*G*I*l**2*m + 144*E**2*I**2*l*q +
-    ...           x**3*(A**2*G**2*l**2*q + 12*A*E*G*I*q) +
-    ...           x**2*(-2*A**2*G**2*l**3*q - 24*A*E*G*I*l*q - 48*A*E*G*I*m) +
-    ...           x*(A**2*G**2*l**4*q + 72*A*E*G*I*l*m - 144*E**2*I**2*q)
-    ...           )/(24*A*E*G*I*(A*G*l**2 + 12*E*I)))
+    >>> dy == (x*(12*E*I*l*(A*G*l**2*q - 2*A*G*l*m + 12*E*I*q)
+    ... + x*(A*G*l*(3*l*(A*G*l**2*q - 2*A*G*l*m + 12*E*I*q) + x*(-2*A*G*l**2*q + 4*A*G*l*m - 24*E*I*q))
+    ... + A*G*(A*G*l**2 + 12*E*I)*(-2*l**2*q + 6*l*m - 4*m*x + q*x**2)
+    ... - 12*E*I*q*(A*G*l**2 + 12*E*I)))/(24*A*E*G*I*(A*G*l**2 + 12*E*I)))
     True
 
     References
@@ -2875,9 +2874,9 @@ def plot_loading_results(self, dir='x', subs=None):
             Plot[1]:Plot object containing:
             [0]: cartesian line: -15*x**2/2 for x over (0.0, 20.0)
             Plot[2]:Plot object containing:
-            [0]: cartesian line: x**2*(150 - 5*x)/8000 - x/8 for x over (0.0, 20.0)
+            [0]: cartesian line: -x**3/1600 + 3*x**2/160 - x/8 for x over (0.0, 20.0)
             Plot[3]:Plot object containing:
-            [0]: cartesian line: x*(105*x**4 - 8026000*x**2/43 + 109200000*x/43 + 304000000/43)/4200000 for x over (0.0, 20.0)
+            [0]: cartesian line: x**5/40000 - 4013*x**3/90300 + 26*x**2/43 + 1520*x/903 for x over (0.0, 20.0)
 
         """
 
diff --git a/sympy/physics/continuum_mechanics/tests/test_beam.py b/sympy/physics/continuum_mechanics/tests/test_beam.py
index 748ce5b1f6c7..c8528de23cf4 100644
--- a/sympy/physics/continuum_mechanics/tests/test_beam.py
+++ b/sympy/physics/continuum_mechanics/tests/test_beam.py
@@ -582,9 +582,9 @@ def test_Beam3D():
     assert b2.reaction_loads == {R1: -750, R2: -750}
 
     b2.solve_slope_deflection()
-    assert b2.slope() == [0, 0, x**2*(50*x - 2250)/(6*E*I) + 3750*x/(E*I)]
-    expected_deflection = (x*(25*A*G*x**3/2 - 750*A*G*x**2 + 4500*E*I +
-        15*x*(750*A*G - 10*E*I))/(6*A*E*G*I))
+    assert b2.slope() == [0, 0, 25*x**3/(3*E*I) - 375*x**2/(E*I) + 3750*x/(E*I)]
+    expected_deflection = 25*x**4/(12*E*I) - 125*x**3/(E*I) + 1875*x**2/(E*I) - \
+        25*x**2/(A*G) + 750*x/(A*G)
     dx, dy, dz = b2.deflection()
     assert dx == dz == 0
     assert dy == expected_deflection
diff --git a/sympy/solvers/ode/single.py b/sympy/solvers/ode/single.py
index 13efc95acba9..f42521ab9f40 100644
--- a/sympy/solvers/ode/single.py
+++ b/sympy/solvers/ode/single.py
@@ -2971,7 +2971,10 @@ def _get_general_solution(self, *, simplify_flag: bool = True):
 
 solver_map = {
     'factorable': Factorable,
-    'nth_algebraic': NthAlgebraic,
+    'nth_linear_constant_coeff_homogeneous': NthLinearConstantCoeffHomogeneous,
+    'nth_linear_euler_eq_homogeneous': NthLinearEulerEqHomogeneous,
+    'nth_linear_constant_coeff_undetermined_coefficients': NthLinearConstantCoeffUndeterminedCoefficients,
+    'nth_linear_euler_eq_nonhomogeneous_undetermined_coefficients': NthLinearEulerEqNonhomogeneousUndeterminedCoefficients,
     'separable': Separable,
     '1st_exact': FirstExact,
     '1st_linear': FirstLinear,
@@ -2984,11 +2987,6 @@ def _get_general_solution(self, *, simplify_flag: bool = True):
     'almost_linear': AlmostLinear,
     'linear_coefficients': LinearCoefficients,
     'separable_reduced': SeparableReduced,
-    'lie_group': LieGroup,
-    'nth_linear_constant_coeff_homogeneous': NthLinearConstantCoeffHomogeneous,
-    'nth_linear_euler_eq_homogeneous': NthLinearEulerEqHomogeneous,
-    'nth_linear_constant_coeff_undetermined_coefficients': NthLinearConstantCoeffUndeterminedCoefficients,
-    'nth_linear_euler_eq_nonhomogeneous_undetermined_coefficients': NthLinearEulerEqNonhomogeneousUndeterminedCoefficients,
     'nth_linear_constant_coeff_variation_of_parameters': NthLinearConstantCoeffVariationOfParameters,
     'nth_linear_euler_eq_nonhomogeneous_variation_of_parameters': NthLinearEulerEqNonhomogeneousVariationOfParameters,
     'Liouville': Liouville,
@@ -2997,6 +2995,8 @@ def _get_general_solution(self, *, simplify_flag: bool = True):
     '2nd_hypergeometric': SecondHypergeometric,
     'nth_order_reducible': NthOrderReducible,
     '2nd_nonlinear_autonomous_conserved': SecondNonlinearAutonomousConserved,
+    'nth_algebraic': NthAlgebraic,
+    'lie_group': LieGroup,
     }
 
 # Avoid circular import:
diff --git a/sympy/solvers/ode/tests/test_ode.py b/sympy/solvers/ode/tests/test_ode.py
index 9224995d1038..1a92901bf8ca 100644
--- a/sympy/solvers/ode/tests/test_ode.py
+++ b/sympy/solvers/ode/tests/test_ode.py
@@ -301,7 +301,7 @@ def test_classify_ode():
     #This is for new behavior of classify_ode when called internally with default, It should
     # return the first hint which matches therefore, 'ordered_hints' key will not be there.
     assert sorted(classify_ode(Eq(f(x).diff(x), 0), f(x), dict=True).keys()) == \
-        ['default', 'nth_algebraic', 'nth_algebraic_Integral', 'order']
+        ['default', 'nth_linear_constant_coeff_homogeneous', 'order']
     a = classify_ode(2*x*f(x)*f(x).diff(x) + (1 + x)*f(x)**2 - exp(x), f(x), dict=True, hint='Bernoulli')
     assert sorted(a.keys()) == ['Bernoulli', 'Bernoulli_Integral', 'default', 'order', 'ordered_hints']
 
diff --git a/sympy/solvers/ode/tests/test_single.py b/sympy/solvers/ode/tests/test_single.py
index 7ef52d08b290..1020f8475d19 100644
--- a/sympy/solvers/ode/tests/test_single.py
+++ b/sympy/solvers/ode/tests/test_single.py
@@ -1014,14 +1014,14 @@ def _get_examples_ode_sol_factorable():
     # This is from issue: https://github.com/sympy/sympy/issues/9446
     'fact_18':{
         'eq': Eq(f(2 * x), sin(Derivative(f(x)))),
-        'sol': [Eq(f(x), C1 + pi*x - Integral(asin(f(2*x)), x)), Eq(f(x), C1 + Integral(asin(f(2*x)), x))],
+        'sol': [Eq(f(x), C1 + Integral(pi - asin(f(2*x)), x)), Eq(f(x), C1 + Integral(asin(f(2*x)), x))],
         'checkodesol_XFAIL':True
     },
 
     # This is from issue: https://github.com/sympy/sympy/issues/7093
     'fact_19': {
         'eq': Derivative(f(x), x)**2 - x**3,
-        'sol': [Eq(f(x), C1 - 2*x*sqrt(x**3)/5), Eq(f(x), C1 + 2*x*sqrt(x**3)/5)],
+        'sol': [Eq(f(x), C1 - 2*x**Rational(5,2)/5), Eq(f(x), C1 + 2*x**Rational(5,2)/5)],
     },
 
     'fact_20': {
@@ -1343,7 +1343,7 @@ def _get_examples_ode_sol_nth_order_reducible():
 
     'reducible_10': {
         'eq': f(x).diff(x, 5) + 2*f(x).diff(x, 3) + f(x).diff(x),
-        'sol': [Eq(f(x), C1 + C2*(x*sin(x) + cos(x)) + C3*(-x*cos(x) + sin(x)) + C4*sin(x) + C5*cos(x))],
+        'sol': [Eq(f(x), C1 + C2*x*sin(x) + C2*cos(x) - C3*x*cos(x) + C3*sin(x) + C4*sin(x) + C5*cos(x))],
         'slow': True,
     },
 
@@ -1357,7 +1357,8 @@ def _get_examples_ode_sol_nth_order_reducible():
     # Needs to be a way to know how to combine derivatives in the expression
     'reducible_12': {
         'eq': Derivative(x*f(x), x, x, x) + Derivative(f(x), x, x, x),
-        'sol': [Eq(f(x), C1 + C2*x + C3/Mul(2, (x + 1), evaluate=False))], # 2-arg Mul!
+        'sol': [Eq(f(x), C1 + C3/Mul(2, (x**2 + 2*x + 1), evaluate=False) +
+        x*(C2 + C3/Mul(2, (x**2 + 2*x + 1), evaluate=False)))], # 2-arg Mul!
         'slow': True,
     },
     }