Open up multiple angles, and integrate the squares of sine and cosine

Sources/AngouriMath/Functions/Continuous/Integration/IntegralPatterns.cs

@@ -19,6 +19,18 @@ internal static class IntegralPatterns
                 TreeAnalyzer.TryGetPolyLinear(arg, x, out var a, out _) =>
                     MathS.Sin(arg) / a,
 
+            // By power reduction: sin(u)^2 = (1 - cos(2u)) / 2, so the integral is
+            // x/2 - sin(2u)/(4a). Without this, integrating sin(x)^2 fell through to
+            // integration by parts and cycled there.
+            Entity.Powf(Entity.Sinf(var arg), Entity.Number.Integer(2)) when
+                TreeAnalyzer.TryGetPolyLinear(arg, x, out var a, out _) =>
+                    x / 2 - MathS.Sin(2 * arg) / (4 * a),
+
+            // cos(u)^2 = (1 + cos(2u)) / 2
+            Entity.Powf(Entity.Cosf(var arg), Entity.Number.Integer(2)) when
+                TreeAnalyzer.TryGetPolyLinear(arg, x, out var a, out _) =>
+                    x / 2 + MathS.Sin(2 * arg) / (4 * a),
+
             Entity.Secantf(var arg) when
                 TreeAnalyzer.TryGetPolyLinear(arg, x, out var a, out _) =>
                     MathS.Hyperbolic.Artanh(MathS.Sin(arg)) / a,

Sources/AngouriMath/Functions/Simplification/Patterns/Patterns.Trigonometry.cs

@@ -79,6 +79,39 @@ internal static partial class Patterns
             _ => x
         };
 
+        /// <summary>
+        /// The largest whole multiplier worth opening up. The expansion of sin(n x) grows
+        /// with n, and past a handful the result is longer than anything it buys.
+        /// </summary>
+        private const int MaxAngleMultiplier = 8;
+
+        /// <summary>
+        /// <c>sin(n x)</c> and <c>cos(n x)</c> written out in <c>sin(x)</c> and
+        /// <c>cos(x)</c>, for a whole n.
+        /// </summary>
+        /// <remarks>
+        /// The expansions are <see cref="TrigonometricAngleExpansion"/>'s, which were only
+        /// ever reached for numeric angles that are multiples of pi. Nothing opened a
+        /// symbolic <c>sin(2x)</c>, so <c>cos(2x) - (1 - 2sin(x)^2)</c> did not reduce to
+        /// zero and neither did <c>(sin(2t)csc(t))^2/4 - cos(2t) - sin(t)^2</c>, which is
+        /// #557.
+        /// </remarks>
+        internal static Entity ExpandMultipleAngleRules(Entity x) => x switch
+        {
+            Sinf(Mulf(Integer n, var inner)) when IsWorthExpanding(n) =>
+                TrigonometricAngleExpansion.ExpandSineArgumentMultiplied(
+                    new Sinf(inner), new Cosf(inner), n.EInteger.ToInt32Checked()),
+
+            Cosf(Mulf(Integer n, var inner)) when IsWorthExpanding(n) =>
+                TrigonometricAngleExpansion.ExpandCosineArgumentMultiplied(
+                    new Sinf(inner), new Cosf(inner), n.EInteger.ToInt32Checked()),
+
+            _ => x
+        };
+
+        private static bool IsWorthExpanding(Integer n)
+            => n.EInteger.Abs() >= 2 && n.EInteger.Abs() <= MaxAngleMultiplier;
+
         internal static Entity CollapseTrigonometricFunctions(Entity x) => x switch
         {
             // sin / cos = tan

Sources/AngouriMath/Functions/Simplification/Simplificator.cs

@@ -109,6 +109,23 @@ void __IterAddHistory(Entity expr)
                     AddHistory(res1);
                     res = PickSimplest(res, res1);
                     AddHistory(res = res.Replace(Patterns.CollapseTrigonometricFunctions).Replace(Patterns.TrigonometricRules));
+
+                    // Multiple angles opened up, then gathered again by the ordinary rules.
+                    // Offered as a candidate rather than taken: written out, sin(4x) is far
+                    // longer than it started, and only worth it where the pieces cancel --
+                    // which is what the complexity metric is for.
+                    var expandedAngles = res
+                        .Replace(Patterns.ExpandMultipleAngleRules)
+                        .Replace(Patterns.NormalTrigonometricForm)
+                        .InnerSimplified;
+                    if (expandedAngles != res)
+                    {
+                        AddHistory(expandedAngles);
+                        AddHistory(SimplifyChildren(expandedAngles)
+                            .Replace(Patterns.TrigonometricRules)
+                            .Replace(Patterns.CommonRules)
+                            .InnerSimplified);
+                    }
                 }
 
 

Sources/Tests/UnitTests/PatternsTest/MultipleAngleTest.cs

@@ -0,0 +1,92 @@
+//
+// Copyright (c) 2019-2022 Angouri.
+// AngouriMath is licensed under MIT.
+// Details: https://github.com/asc-community/AngouriMath/blob/master/LICENSE.md.
+// Website: https://am.angouri.org.
+//
+
+using AngouriMath;
+using AngouriMath.Extensions;
+using Xunit;
+
+namespace AngouriMath.Tests.PatternsTest
+{
+    /// <summary>
+    /// sin(n x) and cos(n x) written out in sin(x) and cos(x). The expanded form is
+    /// offered to the simplifier as a candidate, so what these pin is both that it is
+    /// produced and that it is only preferred where the pieces cancel.
+    /// </summary>
+    public sealed class MultipleAngleTest
+    {
+        // Identities that only close once the multiple angle is opened up.
+        [Theory]
+        [InlineData("sin(2 * x) - 2 * sin(x) * cos(x)")]
+        [InlineData("cos(2 * x) - (1 - 2 * sin(x) ^ 2)")]
+        [InlineData("cos(2 * x) - (cos(x) ^ 2 - sin(x) ^ 2)")]
+        [InlineData("sin(4 * x) - 2 * sin(2 * x) * cos(2 * x)")]
+        public void IdentitiesReduceToZero(string input) =>
+            Assert.Equal(Entity.Number.Integer.Create(0), input.ToEntity().Simplify());
+
+        // The same identity written as cos(2x) - (2cos(x)^2 - 1) gets as far as
+        // 1 - cos(x)^2 - sin(x)^2 and stops: sin^2 + cos^2 = 1 is matched as a pair, and
+        // there the pair sits either side of a subtraction. That is the same
+        // pairwise-versus-flattened gap as #531, in a product of trigonometric terms
+        // rather than a sum, and is not something opening the angle can reach.
+        [Fact]
+        public void PythagoreanPairAcrossASubtractionIsStillMissed() =>
+            Assert.Equal("1 - cos(x) ^ 2 - sin(x) ^ 2".ToEntity(),
+                "cos(2 * x) - (2 * cos(x) ^ 2 - 1)".ToEntity().Simplify());
+
+        [Fact]
+        public void DoubleAngleAndSquareCombine() =>
+            Assert.Equal("cos(x) ^ 2".ToEntity(), "cos(2 * x) + sin(x) ^ 2".ToEntity().Simplify());
+
+        // Left alone where opening it up only makes the expression longer.
+        [Theory]
+        [InlineData("sin(2 * x)")]
+        [InlineData("cos(2 * x)")]
+        [InlineData("sin(3 * x)")]
+        public void CompactFormsAreKept(string input) =>
+            Assert.Equal(input.ToEntity(), input.ToEntity().Simplify());
+
+        // Whatever form comes out, it has to be the same function.
+        [Theory]
+        [InlineData("sin(2 * x)")]
+        [InlineData("cos(2 * x)")]
+        [InlineData("sin(3 * x)")]
+        [InlineData("cos(4 * x)")]
+        [InlineData("cos(2 * x) + sin(x) ^ 2")]
+        public void ExpansionPreservesValue(string input)
+        {
+            var expr = input.ToEntity();
+            var simplified = expr.Simplify();
+            foreach (var point in new[] { 0.37, 1.41, 2.71 })
+            {
+                var before = expr.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble();
+                var after = simplified.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble();
+                Assert.Equal(before, after, 9);
+            }
+        }
+
+        // Power reduction, sin(u)^2 = (1 - cos(2u)) / 2, in the integral table. Before it
+        // was there, integrating sin(x)^2 fell through to integration by parts and cycled.
+        [Theory]
+        [InlineData("sin(x) ^ 2")]
+        [InlineData("cos(x) ^ 2")]
+        [InlineData("sin(2 * x + 1) ^ 2")]
+        [InlineData("cos(3 * x) ^ 2")]
+        public void SquaresOfSineAndCosineIntegrate(string integrand)
+        {
+            var f = integrand.ToEntity();
+            var antiderivative = f.Integrate("x");
+            Assert.DoesNotContain("integral(", antiderivative.Stringize());
+            var derivative = antiderivative.Substitute("C", 0).Differentiate("x");
+            foreach (var point in new[] { 0.37, 1.41, 2.71 })
+            {
+                var expected = f.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble();
+                var actual = derivative.Substitute("x", point).EvalNumerical().RealPart.EDecimal.ToDouble();
+                Assert.Equal(expected, actual, 8);
+            }
+        }
+    }
+}