Attached is at x^4 + a case handling both negative and positive a.
I added separate function to compute root of degree 4, but this
still needs improvement. Positive a catches several examples
in mapleok.input.
Example:
(1) -> integrate(1/(x^4 + a^4), x)
(1)
+-+
+-+ 2 2 +-+ 2 2 x\|2 + a
log(a x\|2 + x + a ) - log(- a x\|2 + x + a ) + 2 atan(---------)
a
+
+-+
x\|2 - a
2 atan(---------)
a
/
3 +-+
4 a \|2
Type: Union(Expression(Integer),...)
Unfortunately,
integrate(1/(x^4 + a^2), x)
produces (4*a^2)^(1/4) and consequently answer is more complicated
then it should be.
--
Waldek Hebisch
--
You received this message because you are subscribed to the Google Groups
"FriCAS - computer algebra system" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to fricas-devel+unsubscr...@googlegroups.com.
To view this discussion on the web visit
https://groups.google.com/d/msgid/fricas-devel/ZWx9Al-Qg0JXsxoX%40fricas.org.
diff --git a/src/algebra/irexpand.spad b/src/algebra/irexpand.spad
index ead11d7f..19f05f65 100644
--- a/src/algebra/irexpand.spad
+++ b/src/algebra/irexpand.spad
@@ -150,10 +150,41 @@ IntegrationResultToFunction(R, F) : Exports == Implementation where
bb := -(r.coef1) k
tantrick(aa * a + bb * b, g) + ilog0(aa, bb, r.coef2, -r.coef1, k)
+ -- Computes (a + b*%i)log(lg(a + b*%i)) + (a - b*%i)log(lg(a - b*%i))
+ -- using Rioboo transformation.
+ root_pair(lg : UP, a : F, b : F, x : Symbol) : F ==
+ lge := quadeval(lg, a, b, -1)
+ f := lge.ans1
+ g := lge.ans2
+ a*log(f*f + g*g) + b*ilog(f, g, x)
+
+ lg2cfunc2 : (UP, UP) -> F
+
+ root4(a : F) : F ==
+ rec := froot(a, 4)$PolynomialRoots(IndexedExponents(K), K, R, P, F)
+ p1 := monomial(1, rec.exponent)$UP - rec.radicand::UP
+ rec.coef*zeroOf(p1)
+
+ quartic(p : UP, lg : UP, x : Symbol) : List(F) ==
+ ground?(rp := reductum(p)) =>
+ a := ground(rp)
+ s := sign(a)
+ s case "failed" => [lg2cfunc2(p, lg)]
+ si := s@Integer
+ si = 1 =>
+ r1 := root4(a/(4::F*leadingCoefficient(p)))
+ [root_pair(lg, r1, r1, x) + root_pair(lg, -r1, r1, x)]
+ si = -1 =>
+ r1 := root4(a)
+ [cmplex(r1, lg) + cmplex(-r1, lg) + root_pair(lg, 0, r1, x)]
+ error "impossible"
+ [lg2cfunc2(p, lg)]
+
lg2func(lg, x) ==
zero?(d := degree(p := lg.coeff)) => error "poly has degree 0"
(d = 1) => [linear(p, lg.logand)]
d = 2 => quadratic(p, lg.logand, x)
+ d = 4 => quartic(p, lg.logand, x)
odd? d and
((r := retractIfCan(reductum p)@Union(F, "failed")) case F) =>
pairsum([cmplex(alpha := rootSimp zeroOf p, lg.logand)],
@@ -162,8 +193,10 @@ IntegrationResultToFunction(R, F) : Exports == Implementation where
lg.logand], x))
[lg2cfunc lg]
- lg2cfunc lg ==
- +/[cmplex(alpha, lg.logand) for alpha in zerosOf(lg.coeff)]
+ lg2cfunc(lg) == lg2cfunc2(lg.coeff, lg.logand)
+
+ lg2cfunc2(p : UP, lg : UP) ==
+ +/[cmplex(alpha, lg) for alpha in zerosOf(p)]
mkRealFunc(l, x) ==
ans := empty()$List(F)
