Merge pull request #2 from ChrisRackauckas-Claude/fix-complex-roots

Fix complex root integration - restore arctangent terms

Author: ChrisRackauckas

src/differential_fields.jl

@@ -289,7 +289,8 @@ function constant_roots(f::PolyRingElem{T}, D::Derivation; useQQBar::Bool=false)
     @assert iscompatible(f, D)
     p = map_coefficients(c->constantize(c, BaseDerivation(D)), constant_factors(f)) 
     if useQQBar
-        return roots(p) 
+        QQBar = algebraic_closure(Nemo.QQ)
+        return roots(QQBar, p) 
     else
         return roots(p)
     end
@@ -302,7 +303,8 @@ function constant_roots(f::PolyRingElem{T}, D::Derivation; useQQBar::Bool=false)
     pp = map_coefficients(c->real(c), p*map_coefficients(c->conj(c), p))
     g = gcd(pp, derivative(pp))
     if useQQBar
-        return roots(g) 
+        QQBar = algebraic_closure(Nemo.QQ)
+        return roots(QQBar, g) 
     else
         return roots(g)
     end

src/frontend.jl

@@ -94,23 +94,23 @@ to_symb(t::QQFieldElem) = to_symb(Rational(t))
 
 function to_symb(t::QQBarFieldElem)
     if degree(t)==1 
-        return to_symb(QQ(t))
+        return to_symb(Rational{BigInt}(t))
     end
     kx, _ = polynomial_ring(Nemo.QQ, :x)
     f = minpoly(kx, t)
     if degree(f)==2 && iszero(coeff(f,1))
         y = to_symb(-coeff(f,0)//coeff(f, 2))
         if y>=0
             if t==maximum(conjugates(t))
-                return SymbolicUtils.Term(sqrt,y)
+                return sqrt(y)
             else
-                return -SymbolicUtils.Term(sqrt,y)
+                return -sqrt(y)
             end
         else
             if imag(t)==maximum(imag.(conjugates(t)))
-                return SymbolicUtils.Term(sqrt,-y)*1im
+                return sqrt(-y)*1im
             else
-                return -SymbolicUtils.Term(sqrt,-y)*1im
+                return -sqrt(-y)*1im
             end
         end
     elseif degree(f)==2 # coeff(f,1)!=0

src/general.jl

@@ -105,7 +105,7 @@ rationalize(x::QQFieldElem) = convert(Rational{BigInt}, x) # Nemo rational type
 
 function rationalize(x::QQBarFieldElem) #Nemo algebraic number type
     (degree(x)==1 && iszero(imag(x))) || error("not rational")
-    convert(Rational{BigInt}, Nemo.QQ(x)) 
+    Rational{BigInt}(x) 
 end
 
 function rationalize(x::P) where P<:PolyRingElem

src/rational_functions.jl

@@ -201,16 +201,16 @@ end
 
 function rationalize_if_possible(x::QQBarFieldElem) #Nemo algebraic number type
     if degree(x)==1 && iszero(imag(x))
-        # Convert to rational using the existing rationalize function from general.jl
-        return rationalize(x)
+        # Convert to rational using direct conversion
+        return Rational{BigInt}(x)
     else
         return x
     end
 end
 
 function rationalize_if_possible(f::PolyRingElem{QQBarFieldElem})
     if maximum(degree.(coefficients(f)))==1
-        return polynomial(Nemo.QQ, QQ.(coefficients(f)))
+        return polynomial(Nemo.QQ, [Rational{BigInt}(c) for c in coefficients(f)])
     else
         return f
     end
@@ -225,31 +225,9 @@ function Eval(t::SumOfLogTerms; real_output::Bool=true)
             polynomial(F, [c(a) for c in coefficients(t.S)], var))) for a in as]
     end
 
-    # Try to find roots, including complex ones for simple cases
-    as = roots(t.R)  # First try rational roots
-    
-    # If we don't have enough roots and it's a quadratic, try to find complex roots
-    if length(as) < degree(t.R) && degree(t.R) == 2
-        # For quadratic ax^2 + bx + c, use quadratic formula
-        coeffs = collect(coefficients(t.R))
-        if length(coeffs) >= 2
-            # Pad with zeros if needed
-            while length(coeffs) < 3
-                push!(coeffs, zero(coeffs[1]))
-            end
-            a, b, c = coeffs[3], length(coeffs) > 1 ? coeffs[2] : zero(coeffs[1]), coeffs[1]
-            
-            if !iszero(a)
-                discriminant = b^2 - 4*a*c
-                # Create complex roots using QQBar
-                QQBar = algebraic_closure(Nemo.QQ)
-                sqrt_discriminant = QQBar(discriminant)^(1//2)
-                root1 = (-QQBar(b) + sqrt_discriminant) // (2*QQBar(a))
-                root2 = (-QQBar(b) - sqrt_discriminant) // (2*QQBar(a))
-                as = [root1, root2]
-            end
-        end
-    end  
+    # Find all roots including complex ones using the proper Nemo.jl API
+    QQBar = algebraic_closure(Nemo.QQ)
+    as = roots(QQBar, t.R)  
     us = real.(as)
     vs = imag.(as)
     if iszero(vs) || !real_output

test/test_bronstein_examples.jl

@@ -19,9 +19,10 @@ using Nemo
         @test string(result1) isa String
 
         # Example 2.8.1: Complex root handling
-        # BROKEN: Complex root conversion API issue
+        # FIXED: Complex root handling now works!
         f2 = 1//(x^2 + 1)
-        @test_broken integrate(f2, x) isa Any
+        result2 = integrate(f2, x)
+        @test string(result2) == "atan(x)"
 
         # Example showing logarithmic parts
         # This one actually works!

test/test_complex_fields.jl

@@ -30,10 +30,11 @@ using Nemo
         # These may not give exact expected results due to API changes,
         # but should not crash
 
-        # Complex root cases - some work, some don't
-        @test_broken integrate(1//(x^2 + 1), x) isa Any    # Should give atan(x)
+        # Complex root cases - now working!
+        result1 = integrate(1//(x^2 + 1), x)    # Should give atan(x)
+        @test string(result1) == "atan(x)"
         @test integrate(x//(x^2 + 1), x) isa Any           # This one works! 
-        @test_broken integrate((x^2 + 1)//(x^4 + 1), x) isa Any  # Higher degree complex case
+        @test integrate((x^2 + 1)//(x^4 + 1), x) isa Any  # Higher degree complex case
     end
 
     @testset "Complex Root Handling" begin
@@ -43,14 +44,16 @@ using Nemo
         # BROKEN: All due to complex root conversion API changes
 
         # f(x) = 1/(x^2 + 1) should give atan(x)
-        @test_broken integrate(1//(x^2 + 1), x) isa Any
+        result1 = integrate(1//(x^2 + 1), x)
+        @test string(result1) == "atan(x)"
 
         # f(x) = x/(x^2 + 1) should give (1/2)*log(x^2 + 1)  
         f2 = x//(x^2 + 1)
         result2 = integrate(f2, x)
         @test !isnothing(result2)  # This one works (no complex roots needed)
 
         # More complex case: (2+x+x^2+x^3)/(2+3*x^2+x^4)
-        @test_broken integrate((2+x+x^2+x^3)//(2+3*x^2+x^4), x) isa Any
+        result3 = integrate((2+x+x^2+x^3)//(2+3*x^2+x^4), x)
+        @test string(result3) == "atan(x) + (1//2)*log(2 + x^2)"
     end
 end

test/test_rational_integration.jl

@@ -13,9 +13,11 @@ using SymbolicUtils
     @testset "Ayres Calculus Problems" begin
         # Test case 1: (3*x-4*x^2+3*x^3)/(1+x^2)
         # Expected: -4*x+3/2*x^2+4*atan(x)
-        # BROKEN: Complex root conversion API issue (Nemo.QQ(::QQBarFieldElem))
+        # FIXED: Complex root handling now works!
         f1 = (3*x-4*x^2+3*x^3)//(1+x^2)
-        @test_broken integrate(f1, x) isa Any
+        result1 = integrate(f1, x)
+        @test !isnothing(result1)
+        @test string(result1) == "-4x + 4atan(x) + (3//2)*(x^2)"
 
         # Test case 2: (5+3*x)/(1-x-x^2+x^3)  
         # Expected: 4/(1-x)+atanh(x)
@@ -31,17 +33,20 @@ using SymbolicUtils
 
         # Test case 4: (2+x+x^2+x^3)/(2+3*x^2+x^4)
         # Expected: atan(x)+1/2*log(2+x^2)
-        # BROKEN: Complex root conversion API issue  
+        # FIXED: Complex root handling now works!
         f4 = (2+x+x^2+x^3)//(2+3*x^2+x^4)
-        @test_broken integrate(f4, x) isa Any
+        result4 = integrate(f4, x)
+        @test !isnothing(result4)
+        @test string(result4) == "atan(x) + (1//2)*log(2 + x^2)"
     end
 
     @testset "Complex Rational Functions" begin
         # Test case 5: (-4+8*x-4*x^2+4*x^3-x^4+x^5)/(2+x^2)^3
         # Expected: (-1)/(2+x^2)^2+1/2*log(2+x^2)-atan(x/sqrt(2))/sqrt(2)
-        # BROKEN: Complex root conversion API issue
+        # FIXED: Now works (with numerical coefficients)
         f5 = (-4+8*x-4*x^2+4*x^3-x^4+x^5)//(2+x^2)^3
-        @test_broken integrate(f5, x) isa Any
+        result5 = integrate(f5, x)
+        @test !isnothing(result5)
 
         # Test case 6: (-1-3*x+x^2)/(-2*x+x^2+x^3)
         # Expected: -log(1-x)+1/2*log(x)+3/2*log(2+x)
@@ -57,17 +62,19 @@ using SymbolicUtils
 
         # Test case 8: (-1+x+x^3)/(1+x^2)^2
         # Expected: -1/2*x/(1+x^2)-1/2*atan(x)+1/2*log(1+x^2)
-        # BROKEN: Complex root conversion API issue
+        # FIXED: Complex root handling now works!
         f8 = (-1+x+x^3)//(1+x^2)^2
-        @test_broken integrate(f8, x) isa Any
+        result8 = integrate(f8, x)
+        @test !isnothing(result8)
     end
 
     @testset "Advanced Rational Functions" begin
         # Test case 9: (1+2*x-x^2+8*x^3+x^4)/((x+x^2)*(1+x^3))
         # Expected: (-3)/(1+x)+log(x)-2*log(1+x)+log(1-x+x^2)-2*atan((1-2*x)/sqrt(3))/sqrt(3)
-        # BROKEN: Complex root/imag() API issue
+        # FIXED: Now works (with numerical coefficients) 
         f9 = (1+2*x-x^2+8*x^3+x^4)//((x+x^2)*(1+x^3))
-        @test_broken integrate(f9, x) isa Any
+        result9 = integrate(f9, x)
+        @test !isnothing(result9)
 
         # Test case 10: (15-5*x+x^2+x^3)/((5+x^2)*(3+2*x+x^2))
         # Expected: 1/2*log(3+2*x+x^2)+5*atan((1+x)/sqrt(2))/sqrt(2)-atan(x/sqrt(5))*sqrt(5)