Re: Rule-based integration
clicl...@freenet.de
Albert schrieb:
>
> I would like to announce the launch of the website
>
> www.apmaths.uwo.ca/RuleBasedMathematics
>
> It is dedicated to dedicated to demonstrating the numerous advantages
> of the rule-based approach to automating mathematics. In systems
> implemented using this approach, rules are expressed as elegant
> mathematical formulas, rather than embedded in conventional program
> code.
>
> As proof-of-concept, I have implemented an efficient and robust Rule-
> based Integrator, nicknamed Rubi. Not only can Rubi compute the
> antiderivative for a broad class of integrands, but the results are
> often significantly superior to those produced by the commercial
> computer algebra systems.
>
> The 1500 or so integration rules Rubi uses are freely available on the
> website in both human and machine readable form. Also available is a
> test suite of over 9400 integration problems developed in conjunction
> with the rules. After reviewing the homepage, I recommend clicking on
> "Highlights of the Indefinite Integration Test Results" for an eye-
> opening comparison of the rule-based integrator (Rubi) with
> Mathematica's and Maple's built-in integrators.
>
> Also if you have access to Mathematica 6 or better, there is a link
> near the bottom of the homepage to download Rubi so you can verify the
> results for yourself. Rubi also provides the option to show the rules
> required to integrate expressions, along with the intermediate
> results. I think this show-step ability has great potential
> pedagogical and research value.
>
> Currently the website is pretty Spartan in format, and limited to
> indefinite integration. However with the help of the computer algebra
> community, I hope it evolves into a true repository of mathematical
> knowledge.
>
Thanks for your announcement. I have a few questions and comments:
1500 rules for indefinite integration amount to 2.5 times the number
implemented in Derive 6.10. A large part of the new rules probably
involve special functions unknown to the Derive kernel, like exponential
integrals, sine and cosine integrals, Fresnel integrals, the Lambert W
function, elliptic integrals, polylogarithms, incomplete gamma
functions, and perhaps even general hypergeometric functions. But there
seem to be additions in the original Derive area of elementary functions
as well (Derive 6.10 knows about ERF, DIGAMMA, DILOG, and Hurwitz ZETA
only).
Does this imply that Rubi can do more elementary antiderivatives than
Derive 6.10 can? That would require that the Rubi tries about as many
integrand transformations, including various substitutions, partial
fraction decomposition and integration by parts, as Derive 6.10 does.
For example, Derive succeeds on INT(x*(x^2+3)/(2*a^2+b^2*(x^2+1))^
(5/2)*LN((SQRT(2)*x*SQRT(2*a^2+b^2*(x^2+1))-2*x*a+b*(x^2+1))/x),x), but
the task is obviously hard - and Derive's answer is characteristically
incomplete. Does Rubi do better on this one?
Derive 6.10 can do INT(1/SQRT(1-a*x),x), of course, but produces an
incorrect result for the equivalent INT((LN(a*x-1)-2*LN(-SQRT(a*x-1)))/
(2*pi*SQRT(a*x-1)),x), as a graphic comparison of the antiderivatives
will show. Has this bug been repaired in Rubi?
In my feeling, the choice to base Rubi on Mathematica may seriously
impede a pedagogical mission. What capabilities of MMA are actually used
in the implementation of Rubi? How many man-years would it take to move
Rubi to a freely available programming language like LISP or Python?
Martin.
clicl...@freenet.de
Albert schrieb:
>
> [...]
>
> Rubi is unable to integrate
>
> x*(x^2+3)/
> (2*a^2+b^2*(x^2+1))^(5/2)*log((sqrt(2)*x*sqrt(2*a^2+b^2*(x^2+1))-2*x*a
> +b*(x^2+1))/x)
>
> and Mathematica returns a huge, multipage result...
>
I wonder why Rubi fails here. Integration by parts produces an algebraic
integrand containing SQRT(b^2*x^2 + 2*a^2 + b^2) and no other root.
Derive 6.10 discovers this route, but then does not fully succeed on the
algebraic problem. Yet this algebraic integral would obviously yield to
one of Euler's substitutions, for example. In particular, Derive's
result contains the following three unevaluated integrals:
INT(SQRT(b^2*x^2+2*a^2+b^2)/(b^3*x^4+4*a*b^2*x^3+2*a^2*b*x^2+4*a~
*x*(2*a^2+b^2)-b*(2*a^2+b^2)),x)
INT(SQRT(b^2*x^2+2*a^2+b^2)/(b^3*x^6+4*a*b^2*x^5+b*x^4*(2*a^2+b^~
2)+8*a*x^3*(a^2+b^2)-b^3*x^2+4*a*x*(2*a^2+b^2)-b*(2*a^2+b^2)),x)
INT(x/((b^2*x^2+2*a^2+b^2)*((b*x^2-2*a*x+b)*SQRT(b^2*x^2+2*a^2+b~
^2)+SQRT(2)*b^2*x^3+SQRT(2)*x*(2*a^2+b^2))),x)
However, Derive can evaluate these too, once they are rewritten as:
INT(SQRT(b^2*x^2+2*a^2+b^2)*(-a*b^2*x/((9*a^4+6*a^2*b^2+b^4)*(b^~
2*x^2+2*a^2+b^2))-b*(a^2+b^2)/(2*(9*a^4+6*a^2*b^2+b^4)*(b^2*x^2+~
2*a^2+b^2))+a*b*x/((9*a^4+6*a^2*b^2+b^4)*(b*x^2+4*a*x-b))+(9*a^2~
+b^2)/(2*(9*a^4+6*a^2*b^2+b^4)*(b*x^2+4*a*x-b))),x)
INT(SQRT(b^2*x^2+2*a^2+b^2)*(b^4*x/(2*a*(9*a^4+6*a^2*b^2+b^4)*(b~
^2*x^2+2*a^2+b^2))+b^3*(a^2+b^2)/(4*a^2*(9*a^4+6*a^2*b^2+b^4)*(b~
^2*x^2+2*a^2+b^2))+a*b*x*(9*a^2+2*b^2)/(2*(4*a^2+b^2)*(9*a^4+6*a~
^2*b^2+b^4)*(b*x^2+4*a*x-b))+(72*a^4+21*a^2*b^2+b^4)/(4*(4*a^2+b~
^2)*(9*a^4+6*a^2*b^2+b^4)*(b*x^2+4*a*x-b))-x/(2*a*(x^2+1)*(4*a^2~
+b^2))-b/(4*a^2*(x^2+1)*(4*a^2+b^2))),x)
INT((SQRT(b^2*x^2+2*a^2+b^2)-SQRT(2)*x)*(b*x*(a^2+1)*(b^2-2)^2/(~
4*(2*a^2+b^2)*(a^4*(3*b^2-4)+a^2*b^4+b^2)*(x^2*(b^2-2)+2*a^2+b^2~
))+a*(b^2-2)^2/(4*(a^4*(3*b^2-4)+a^2*b^4+b^2)*(x^2*(b^2-2)+2*a^2~
+b^2))-b^3*x/(4*(2*a^2+b^2)*(3*a^2+b^2)*(b^2*x^2+2*a^2+b^2))-b^2~
/(4*a*(3*a^2+b^2)*(b^2*x^2+2*a^2+b^2))),x)+INT((SQRT(b^2*x^2+2*a~
^2+b^2)-SQRT(2)*x)*(b^2*x*(a^2-b^2+3)/(4*(3*a^2+b^2)*(a^4*(3*b^2~
-4)+a^2*b^4+b^2)*(b*x^2-2*a*x+b))+b*(2*a^2*(2*b^2-3)+b^2)/(4*a*(~
3*a^2+b^2)*(a^4*(3*b^2-4)+a^2*b^4+b^2)*(b*x^2-2*a*x+b))),x)
In the first two cases, the rational part has simply been expanded into
partial fractions. In the third case the denominator has been
rationalized, the rational part subsequently been expanded into partial
fractions, and the single integral then split into two. The latter step
seems to be necessary if you want Derive to evaluate the third integral
in reasonable time.
I conclude that, on Rubi too, work remains to be done in this area.
Martin.
Vladimir Bondarenko
One tiny example.
The VM machine says:
Mathematica 7.0.1.0 returns
Integrate[1/(Sqrt[2] + Sin[z] + Cos[z]), z]
unevaluated (!).
Rubi returns
-(2/(1 + (-1 + Sqrt[2]) Tan[z/2]))
which is a correct answer.
Best wishes,
Vladimir Bondarenko
VM and GEMM architect
Co-founder, CEO, Mathematical Director
http://www.cybertester.com/ Cyber Tester, LLC
http://maple.bug-list.org/ Maple Bugs Encyclopaedia
http://www.CAS-testing.org/ CAS Testing
-----------------------------------------------------
"We must understand that technologies
like these are the way of the future."
-----------------------------------------------------
> Aloha from Hawaii,
> Albert D. Rich
David Bernier
>
> For most integrands of the form 1/(a+b*Sin[x]+c*Cos[x]) Mathematica
> can use the rule
>
> Integrate[1/(a+b*Sin[x]+c*Cos[x]), x] -->
> 2 ArcTan[(b+(a-c)*Tan[x/2])/Sqrt[a^2-b^2-c^2]] / Sqrt[a^2-b^2-c^2]
>
> Obviously this rule cannot be used when a^2-b^2-c^2 equals 0, which is
> probably the reason Mathematica 7 cannot integrate 1/(Sqrt[2]+Sin[x]
> +Cos[x]). But, interestingly Mathematica is able to integrate 1/
> (5+4*Sin[x]+3*Cos[x])...
>
> When a^2-b^2-c^2 equal 0, Rubi uses the rule
>
> Integrate[1/(a+b*Sin[x]+c*Cos[x]), x] --> -2/(b+(a-c)*Tan[x/2])
>
> which gives the simple answer Vladimir reports.
This is quite interesting.
Suppose we add real coefficients h and k for sin(2x)
and cos(2x) respectively.
Then, in the Mathematica notation, I arrive at this:
Integrate[1/(a+b*Sin[x]+c*Cos[x]+h*Sin[2x]+k*Cos[2x]), x]
Of course, there could be quite a few special cases here.
David Bernier
Peter Pein
Hi Vladimir,
nice to read you again :-))
my solution involves "PetSie.nb" which can be found at
http://dl.dropbox.com/u/3030567/Mathematica/PetSie.nb
In[3]:= <<PetSie`
In[6]:=
Assuming[Pi>x>0,ChangeVar[Tan[z/2]==t,Integrate[1/(Sqrt[2]+Sin[z]+Cos[z]),{z,0,x}],t]]
Out[6]= (2 (-1+Sqrt[2])Sin[x/2])/(Cos[x/2]+(-1+Sqrt[2]) Sin[x/2])
there is a difference between the usage of a cas and simply asking a
cas.
Greetings, Peter
Am Wed, 9 Jun 2010 05:14:21
-0700 (PDT) schrieb Vladimir Bondarenko <v...@cybertester.com>:
jakubi
> I would like to announce the launch of the website
>
> www.apmaths.uwo.ca/RuleBasedMathematics
>
> It is dedicated to dedicated to demonstrating the numerous advantages
> of the rule-based approach to automating mathematics. In systems
> implemented using this approach, rules are expressed as elegant
> mathematical formulas, rather than embedded in conventional program
> code.
>
> The 1500 or so integration rules Rubi uses are freely available on the
> website in both human and machine readable form.
Very interesting. I wonder how have you checked that these O(10^3)
rules are mutually exclusive. And how would newly contributed rules be
checked against previous ones.
Regards, Alejandro Jakubi
clicl...@freenet.de
Albert schrieb:
>
> [...]
>
> Mutual exclusivity of rules (or conversely missing rules) for
> integrands of a particular form can be determined by examining the
> application conditions on those rules. For example, the following is
> an excerpt of the rational function integration rules posted on the
> website:
>
> http://www.apmaths.uwo.ca/~arich/RuleRepository/MathematicalNotation/x^m(a+bx)^nRules.pdf
>
> It shows 7 rules required to integrate expressions of the form x^m (a
> +b x)^n, where m and n are integers. Their mutual exclusivity can be
> verified by examining the application conditions given before each
> rule. Polynomial expansion or partial fraction expansion is used to
> expand those integrands for which none of these 7 rules apply. For
> example, x^4/(a+b x)^10 is expanded using pfe.
>
So there are two kinds of rules: "ordinary" rules where mutual
exclusivity is desirable (otherwise the result would depend on the order
of application by the pattern matcher), and "super" rules that recast
the problem when the pattern matcher has failed to find a match. Then
the ordinary rules are traversed again, and perhaps this cycle is
iterated with the same or a different super rule. I suppose
trigonometric transformations and (very importantly) integration by
parts are also included among the super rules.
If not included already, I propose to add the rationalization of
algebraic denominators to this group, where integrands of the type
C(x)/(A(x)*SQRT(F(x))+B(x)*SQRT(G(x))) are replaced by
C(x)*(A(x)*SQRT(F(x))-B(x)*SQRT(G(x))) / (A(x)^2*F(x)-B(x)^2*G(x)). The
fairly common functions C(x)/P(x,SQRT(F(x))), where P(x,y) is a
polynomial in two variables, represent just a special case of this (with
G(x)=1).
Here is a small gallery of fairly simple algebraic integrands on which I
have found Derive 6.10 to fail (as in an earlier example, more or less
obvious transformations on part of the operator help it succeed though):
INT(1/(SQRT(x^2+1)+2*x)^2,x)
INT(1/(SQRT(x^2-1)*(3*x^2-4)^2),x)
INT(1/(2*SQRT(x)+SQRT(x+1))^2,x)
INT(SQRT(x^2-1)/(x-#i)^2,x)
INT(1/(SQRT(x^2-1)*(x^2+1)^2),x)
INT(1/(SQRT(x-1)*(SQRT(x-1)+2*SQRT(x))^2),x)
INT(1/(SQRT(x^2-1)*(SQRT(x^2-1)+SQRT(x))^2),x)
INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
INT(((x-1)^(3/2)+(x+1)^(3/2))/((x+1)^(3/2)*(x-1)^(3/2)),x)
Martin.
Richard Fateman
If you use a program that looks like this:
toplevel(integrand, var) :=
repeat until no more changes {
if match(super_pattern1,integrand,var)
then
integrand := apply(super-rhs1, run_some_algorithm(integrand), var);
else
if match(.....)
}
If the run_some_algorithm includes much of the traditional computer
algebra repertoire, it is not surprising that the integration
program would be at least as successful. An important but tricky
optimization is ordering the matching so that no time is wasted
on rules that cannot possibly match.
An advantage of the rule-based approach is that it is possible,
though not necessarily as easy as it might seem, to insert an
"optimal patch". This would be a rule that adds a new capability, but
does it in such a way as to interfere with nothing else. And to not
slow down anything else. Finding a place to put an optimal patch
in a conventional program requires tracking down where some critical
decision is made and putting in an extra test.
Ordering the rules automatically so that all possibly applicable rules
are tried, and the ruleset converges, are challenges.
There was period during which people wrote rule-based expert
systems (and meta-systems), and tried various things to keep correctness
and cost under control as the systems were scaled up in size.
(1980's?) I think people had enough bad experiences to sour them
on "Artificial Intelligence" for a few decades.
Nevertheless, if Albert's stuff works, there is no arguing with
success!
RJF
Axel Vogt
> have found Derive 6.10 to fail (as in an earlier example, more or less
> obvious transformations on part of the operator help it succeed though):
>
> INT(1/(SQRT(x^2+1)+2*x)^2,x)
>
> INT(1/(SQRT(x^2-1)*(3*x^2-4)^2),x)
>
> INT(1/(2*SQRT(x)+SQRT(x+1))^2,x)
>
> INT(SQRT(x^2-1)/(x-#i)^2,x)
>
> INT(1/(SQRT(x^2-1)*(x^2+1)^2),x)
>
> INT(1/(SQRT(x-1)*(SQRT(x-1)+2*SQRT(x))^2),x)
>
> INT(1/(SQRT(x^2-1)*(SQRT(x^2-1)+SQRT(x))^2),x)
>
> INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
>
> INT(((x-1)^(3/2)+(x+1)^(3/2))/((x+1)^(3/2)*(x-1)^(3/2)),x)
>
> Martin.
Side remark: Maple 12 finds them directly, except the one
before the last (though I neither checked the answers nor
want to / can to judge the 'simplicity' of the results).
Nasser M. Abbasi
That is good.
I just did a quick test on the above 9 integrations posted by Martin,
using Mathematica 7 and Rubi, and numbers 7 and 8 in the list did not
evaluate.
To compare the result from Matematica and Rubi, I cheated, I just
plotted the difference over some range, and in all cases I get a zero
plot. So numerically the result is the same. Not the best way to check
equivalence of the 2 analytical results, but a quick and dirty way for now.
If someone wants to see the output, here is a link
http://12000.org/my_notes/rubiTest/rubi_1/rubi.html
and the notebook is here
http://12000.org/my_notes/rubiTest/rubi_1/rubi.nb
--Nasser
clicl...@freenet.de
"Nasser M. Abbasi" schrieb:
> On 6/20/2010 11:36 AM, Axel Vogt wrote:
> > clicl...@freenet.de wrote:
> >>
> >> Here is a small gallery of fairly simple algebraic integrands on which I
> >> have found Derive 6.10 to fail (as in an earlier example, more or less
> >> obvious transformations on part of the operator help it succeed though):
> >>
> >> INT(1/(SQRT(x^2+1)+2*x)^2,x)
> >>
> >> INT(1/(SQRT(x^2-1)*(3*x^2-4)^2),x)
> >>
> >> INT(1/(2*SQRT(x)+SQRT(x+1))^2,x)
> >>
> >> INT(SQRT(x^2-1)/(x-#i)^2,x)
> >>
> >> INT(1/(SQRT(x^2-1)*(x^2+1)^2),x)
> >>
> >> INT(1/(SQRT(x-1)*(SQRT(x-1)+2*SQRT(x))^2),x)
> >>
> >> INT(1/(SQRT(x^2-1)*(SQRT(x^2-1)+SQRT(x))^2),x)
> >>
> >> INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
> >>
> >> INT(((x-1)^(3/2)+(x+1)^(3/2))/((x+1)^(3/2)*(x-1)^(3/2)),x)
> >>
> >
> > before the last (though I neither checked the answers nor
> > want to / can to judge the 'simplicity' of the results).
>
> That is good.
>
> I just did a quick test on the above 9 integrations posted by Martin,
> using Mathematica 7 and Rubi, and numbers 7 and 8 in the list did not
> evaluate.
>
> To compare the result from Matematica and Rubi, I cheated, I just
> plotted the difference over some range, and in all cases I get a zero
> plot. So numerically the result is the same. Not the best way to check
> equivalence of the 2 analytical results, but a quick and dirty way for now.
>
> If someone wants to see the output, here is a link
>
> http://12000.org/my_notes/rubiTest/rubi_1/rubi.html
>
> and the notebook is here
>
> http://12000.org/my_notes/rubiTest/rubi_1/rubi.nb
>
In the case of MMA the effect is called "Risch integration (not fully
implemented)". Could someone perhaps try them on Maxima and on FriCAS
too?
Martin.
clicl...@freenet.de
Albert schrieb:
>
> Maple 13 is able to simplify to 0 the integrand minus the derivative
> of the huge antiderivative it returns for problem 7 in Martin's list.
> So Maple beats Mathematica 7 and Rubi on this problem! However, I
> would be willing to bet there is a much simpler antiderivative...
>
> However, none of these systems can integrate problem 8. Therefore it
> seems to me the challenger (i.e. Martin) is now obligated to provide
> us with a closed-form antiderivative...
>
I do accept this obligation. But note that the performance statistics
for the venerable Lady Maxima and for the renowned Sir Axiom, also
performing under the alias of FriCAS, are not yet in. In fact, I expect
the latter to be able to relieve me from my obligation. Hence only this
much for now: compact (though not optimized) forms of the
antiderivatives 7 and 8 are roughly 350 and 500 Bytes in size.
I realize that the advance that Rubi represents over Derive 6.10 is
truly substantial!
Martin.
Waldek Hebisch
>
> "Nasser M. Abbasi" schrieb:
> > On 6/20/2010 11:36 AM, Axel Vogt wrote:
> > > clicl...@freenet.de wrote:
> > >>
> > >> Here is a small gallery of fairly simple algebraic integrands on which I
> > >> have found Derive 6.10 to fail (as in an earlier example, more or less
> > >> obvious transformations on part of the operator help it succeed though):
> > >>
> > >> INT(1/(SQRT(x^2+1)+2*x)^2,x)
> > >>
> > >> INT(1/(SQRT(x^2-1)*(3*x^2-4)^2),x)
> > >>
> > >> INT(1/(2*SQRT(x)+SQRT(x+1))^2,x)
> > >>
> > >> INT(SQRT(x^2-1)/(x-#i)^2,x)
> > >>
> > >> INT(1/(SQRT(x^2-1)*(x^2+1)^2),x)
> > >>
> > >> INT(1/(SQRT(x-1)*(SQRT(x-1)+2*SQRT(x))^2),x)
> > >>
> > >> INT(1/(SQRT(x^2-1)*(SQRT(x^2-1)+SQRT(x))^2),x)
> > >>
> > >> INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
> > >>
> > >> INT(((x-1)^(3/2)+(x+1)^(3/2))/((x+1)^(3/2)*(x-1)^(3/2)),x)
> > >>
> > >
> > > Side remark: Maple 12 finds them directly, except the one
> > > before the last (though I neither checked the answers nor
> > > want to / can to judge the 'simplicity' of the results).
> >
> > That is good.
> >
> > I just did a quick test on the above 9 integrations posted by Martin,
> > using Mathematica 7 and Rubi, and numbers 7 and 8 in the list did not
> > evaluate.
> >
>
> implemented)". Could someone perhaps try them on Maxima and on FriCAS
> too?
>
1-5 and 9 evaluate quickly (I did not check correctness). 6, 7, 8
are slow. In my experience integrals like 6 and 7 worked but
needed a lot of time (days or weeks). I left nr 8 running, will
report about result.
--
Waldek Hebisch
heb...@math.uni.wroc.pl
Axel Vogt
((x^4+1)^(1/2)+x^2)^(1/2)/(x+1)^2/(x^4+1)^(1/2) =
1/2*2^(1/2)/(x+1)^2/(x^2+I)^(1/2) +
1/2*2^(1/2)/(x+1)^2/(x^2-I)^(1/2)
and for that Maple finds anti-derivatives (though they
lengthy and involve ln and sqrt, but they are real and
obey the problem in x=-1 and can be simplified):
((-1-I)*(2*x^2-2*I)^(1/2)+(-1+I)*(2*x^2+2*I)^(1/2)+
(-I*(1-I)^(1/2)*ln((-2*I-2*x+2*(1-I)^(1/2)*(x^2-I)^(1/2))/
(x+1))+I*(1+I)^(1/2)*ln((2*I-2*x+2*(1+I)^(1/2)*
(x^2+I)^(1/2))/(x+1)))*(x+1)*2^(1/2))/(4*x+4)
I did cross checking by differentiation against the task
only by plotting or a (finite) Taylor series (seems some
work is needed to convince Maple to proof it).
Axel Vogt
> On Jun 21, 8:26 am, Axel Vogt <&nore...@axelvogt.de> wrote:
>
>> problem 8:
>>
>> ((x^4+1)^(1/2)+x^2)^(1/2)/(x+1)^2/(x^4+1)^(1/2) =
>>
>> 1/2*2^(1/2)/(x+1)^2/(x^2+I)^(1/2) +
>> 1/2*2^(1/2)/(x+1)^2/(x^2-I)^(1/2)
>
> terms. Unfortunately, however, the above identity is not valid for
> all x. For example, when x=1+2*I, the left and rights sides evaluate
> to 0.00779279 - 0.00161297 I and -0.0714979 - 0.0344563 I,
> respectively.
Agreed, there is a brunch cut problem (term1 - term2 would be
correct, but I do not want to correct it by paper and pencil).
Fine, if your system finds it without pain!
Vladimir Bondarenko
z)
while Rubi fails to handle Tan[z]/(Sin[2 z] + Cos[2 z]^2)...
I wonder how the Derive 6.10 and Rubi indefinite integrators compare
on
the elementary function integrands?
Do you expect with Rubi to surpass the Derive 6.10 indefinite
integration
facility?
Do you expect to reinforce the Rubi so it would be stronger in
indefinite
integration than the Mathematica 7 integrator?
On Jun 4, 6:18 am, Albert <Albert_R...@msn.com> wrote:
> I would like to announce the launch of the website
>
> www.apmaths.uwo.ca/RuleBasedMathematics
>
> It is dedicated to dedicated to demonstrating the numerous advantages
> of the rule-based approach to automating mathematics. In systems
> implemented using this approach, rules are expressed as elegant
> mathematical formulas, rather than embedded in conventional program
> code.
>
> As proof-of-concept, I have implemented an efficient and robust Rule-
> based Integrator, nicknamed Rubi. Not only can Rubi compute the
> antiderivative for a broad class of integrands, but the results are
> often significantly superior to those produced by the commercial
> computer algebra systems.
>
> The 1500 or so integration rules Rubi uses are freely available on the
Daniel Lichtblau
In[12]:= expr = y(y^2+1)/(y^4+2y^3-2y^2+2y+1);
In[13]:= solns = y /. Solve[Denominator[expr]==0, y];
In[14]:= quadratics = Union[Cases[solns, Power[_Integer,1/2], Infinity]];
In[15]:= InputForm[Apart[Factor[expr, Extension->quadratics]]]
Out[15]//InputForm=
((-1 + Sqrt[5])*y)/(2*Sqrt[5]*(1 + y - Sqrt[5]*y + y^2)) +
((1 + Sqrt[5])*y)/(2*Sqrt[5]*(1 + y + Sqrt[5]*y + y^2))
Daniel Lichtblau
Wolfram Research
clicl...@freenet.de
Axel Vogt schrieb:
Your transformation is valid for Re(x^2) > 0, and so your result will be
valid in this region (which includes the positive and negative real
axis); but because the transformation is invalid outside this region, so
should be your result.
However, there is a indeed a clean way of splitting the factor
SQRT(SQRT(x^4+1) + x^2) / (SQRT(x^4+1)) in the integrand, namely as
1/SQRT(2) (SQRT(#i)/SQRT(1 + #i*x^2) + SQRT(-#i)/SQRT(1 - #i*x^2)).
Subsequent integration by Rubi, Maple, or Mathematica should be no
problem. But can the systems be instructed to remove the imaginary units
from the antiderivative? Can Maple and Mathematica verify the result for
the entire complex plane?
This complex splitting method of solving integral no 8 is simple and
elegant; I therefore propose to include the most general version of the
transformation that can be worked out and all possible variants of it
(such as the case resulting for x = #i*y) with Rubi's set of integration
rules.
Martin.
Axel Vogt
>
Aha, one has to care only for the radicals and uses 1/(x+1)^2 later
(and then Maple finds an antiderivative).
I could not convince Maple to proof your identity A=B directly
(where I assume there is a multiplication after 1/SQRT(2), yes?).
# translate Derive to Maple notation, saves typing ...
"SQRT(SQRT(x^4+1) + x^2) / (SQRT(x^4+1)) = 1/SQRT(2) *
(SQRT(#i)/SQRT(1 + #i*x^2) + SQRT(-#i)/SQRT(1 - #i*x^2))";
StringTools[SubstituteAll]( %, "INT", "Int" ):
StringTools[SubstituteAll]( %, "SQRT", "sqrt" ):
StringTools[SubstituteAll]( %, "#i", "I" ):
A,B:=lhs(%), rhs(%);
Now
(A/B); # this should be 1, we only take the square
simplify(%): simplify(%, size); # simplify it
%^2;
numer(%) - denom(%); # should be 0
expand(%):
C:=combine(%): # simplified, show C=0
RootOf(%, x); # solve this algebraic equation
[allvalues(%)]; # and show all solutions
(1/4)
[_X, -I (-1) ]
This says: _any_ x is a solution. And the second value as well.
(if one trusts the system ...).
Both MultiSeries:-limit(A/B, x=-I*(-1)^(1/4), complex)
and MultiSeries:-limit(C, x=-I*(-1)^(1/4), complex) give
what the should give (1 and 0).
Now I am too lazy to care for +-1 = sqrt(1) ...
Axel Vogt
> Your transformation is valid for Re(x^2) > 0, and so your result will be
> valid in this region (which includes the positive and negative real
> axis); but because the transformation is invalid outside this region, so
> should be your result.
>
> However, there is a indeed a clean way of splitting the factor
> This complex splitting method of solving integral no 8 is simple and
> elegant; I therefore propose to include the most general version of the
> transformation that can be worked out and all possible variants of it
> (such as the case resulting for x = #i*y) with Rubi's set of integration
> rules.
I think I should sketch (in a _polished_ way) how one can come across
it (since just looking at the (corrected) result does not give it),
may be it is a kind of recipe as well.
The integrand is algebraic and univariate, hence it satisfies a DE
with polynomial coefficients (have no reference).
In that case (a standard attack with 'dpolyform' fails) one can use
the package 'gfun' to find first a polynomial equation, which is
satisfied by the integrand by 'algfuntoalgeq').
In a second step the package has 'algeqtodiffeq' to give the DE for
the problem and 'dsolve' two linear independent solutions; they do
involve integrand(0) and the 2nd derivative of it in 0.
I got my coefficients by trusting that and using the known values
for that (hence there might be a bug in the used routines).
However it is *much* more simple to try a splitting method directly,
if having that in mind. Just thought it might be worth to sketch a
way, how to get a different view on such integrands.
clicl...@freenet.de
Albert schrieb:
>
> [...]
>
> The Rubi result differentiates back to the expanded form of the
> original integrand. So it is a valid antiderivative for problem 7,
> assuming your expansion is valid for all real and complex x. It
> satisfies all numerical tests I have run, but I am unable to get
> Mathematica to symbolically confirm it. You wrote
>
> > This complex splitting method of solving integral no 8 is simple and
> > elegant; I therefore propose to include the most general version of the
> > transformation that can be worked out and all possible variants of it
> > (such as the case resulting for x = #i*y) with Rubi's set of integration
> > rules.
>
> expansion identities, preferably in the form of transformation rules
> like those posted on the rule-based integration website. Then they
> can be directly incorporated into Rubi.
>
Some clarifications first: I was not aware that this kind of splitting
is possible until I saw the post by Axel Vogt. As Axel explained, he
discovered (a somewhat incorrect version of) the formula by deriving and
solving a differential equation for the entire integrand using Maple.
When you pointed out the branch-cut problem with his transformation, my
intuition told me that a correct version should exist; the problem thus
remained on my mind until three days later (while awaking in the
morning) I "saw" the necessary correction, and later that day confirmed
it graphically and posted it here. I also used it to integrate problem 8
on Derive 6.10, and (with some massaging) obtained the simple
antiderivative:
INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
(LN(-SQRT(#i-1)*SQRT(1-#i*x^2)-x+#i)/(2*SQRT(#i-1))+LN(-SQRT(-#i~
-1)*SQRT(1+#i*x^2)-x-#i)/(2*SQRT(-#i-1)))-SQRT(SQRT(2)-1)*LN(x+1~
)/2-SQRT(SQRT(x^4+1)+1)/(SQRT(2)*(x+1))
As regards a manual derivation of the splitting relation, I can only
offer what you are probably aware of already:
SQRT(SQRT(1+x^4) + x^2) / SQRT(1+x^4)
SQRT(SQRT((1 - #i*x^2)*(1 + #i*x^2)) + x^2) /
SQRT((1 + #i*x^2)*(1 - #i*x^2))
SQRT(SQRT(1 - #i*x^2)*SQRT(1 + #i*x^2) + x^2) /
(SQRT(1 + #i*x^2) * SQRT(1 - #i*x^2))
SQRT(#i*(1 - #i*x^2) + 2*#i*SQRT(1 - #i*x^2)*(-#i)*SQRT(1 + #i*x^2) +
(-#i)*(1 + #i*x^2)) / (SQRT(2) * SQRT(1 + #i*x^2) * SQRT(1 - #i*x^2))
(SQRT(#i)*SQRT(1 - #i*x^2) + SQRT(-#i)*SQRT(1 + #i*x^2)) /
(SQRT(2) * SQRT(1 + #i*x^2) * SQRT(1 - #i*x^2))
1/SQRT(2) * (SQRT(#i)/SQRT(1 + #i*x^2) + SQRT(-#i)/SQRT(1 - #i*x^2))
I am in doubt whether this derivation should be called "simple and
elegant"; the validity of the numerator transformation may not be
obvious to many.
As regards possible generalizations needed for a good Rubi rule, I
suspect that SQRT(1+x^4) and SQRT(1 +- #i*x^2) in the formula can be
replaced by SQRT(a+x^4) and SQRT(SQRT(a) +- #i*x^2), respectively, where
RE(a) must probably be > 0 or perhaps only >= 0 (as I haven't looked
into this, the range of values actually allowed in the complex plane may
also be a smaller "wedge" only). For pattern matching, x should
presumably be generalized to b*x+c. I don't know if a separate rule
based on SQRT(SQRT(1+x^4) - x^2) / SQRT(1+x^4) = 1/SQRT(2) *
(SQRT(-#i)/SQRT(1 + #i*x^2) + SQRT(#i)/SQRT(1 - #i*x^2)) would still be
needed then; this variant is obtained upon replacing x by #i*x.
Martin.
Axel Vogt
> On Jun 26, 6:31 am, cliclic...@freenet.de wrote:
>
>> I am in doubt whether this derivation should be called "simple and
>> elegant"; the validity of the numerator transformation may not be
>> obvious to many.
>
> suggested, I generalized your identity to
>
> sqrt(z+sqrt(a+z^2))/sqrt(a+z^2) =
> sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*z)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*z))
>
> which appears to be valid for all z provided a>0. Its derivation
> tracks yours as follows
...
But that follows from special case a=1 if one uses x=z/sqrt(a)
and then multiplies the identity by (sqrt(a))^2, since for 0<a
one has sqrt(a)*sqrt(w) = sqrt(a*w), no?
> sqrt(a+z^2) = sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z)
> which is valid for all z provided a>0.
...
Just for the records (as you certainly do not have the 'obvious
and false proof' in mind): sqrt(1-x)*sqrt(1+x)=sqrt(1-x^2) is
shown in Davenport & Corless, Reasoning about the Elementary
Functions of Complex Analysis (2000). Now using z=x/sqrt(a)*I
and multiplying by a gives it (Maple does not know that, it
even states it as false for a=1).
...
> Aloha,
> Albert
Another way to look at these identities is 'denesting radicals'
(as you did with Jeffrey in 'Simplifying square roots of square
roots by denesting').
Is there a way to make use of Rubi's abilities calling it from
Maple?
Nasser M. Abbasi
>
> If someone wants to see the output, here is a link
>
> http://12000.org/my_notes/rubiTest/rubi_1/rubi.html
>
fyi,
I added Maple 14 output to the above. The link to Maple output is
contained in the above, near the top.
Maple 14 did them all, except for Martin's number 8 problem.
--Nasser
Vladimir Bondarenko
The Cyber Tester's team is interested in improvements along
the direction you offer, rule-based approach to automating
mathematics.
As you know, we do not test software ourselves; instead, this
work is done by our VM machine.
It is important to learn how much we could help you and the
community in reinforcing Rubi as the very first step of your
project.
To know this, it would be nice and fair to learn, before we
start presenting the results by the VM machine,
1) how much time have you spent writing the Rubi?
2) how much time you and your beta testers spent for QA?
3) is there a list of known bugs in the Rubi?
About the time, ANY estimation however approximate is good.
Best wishes from algorithmic Simferopol,
Vladimir
clicl...@freenet.de
Albert schrieb:
> > elegant"; the validity of the numerator transformation may not be
> > obvious to many.
>
> suggested, I generalized your identity to
>
> sqrt(z+sqrt(a+z^2))/sqrt(a+z^2) =
> sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*z)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*z))
>
> which appears to be valid for all z provided a>0. Its derivation
> tracks yours as follows
>
>
> 2) sqrt(z+sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))/
> (sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 3) sqrt(2*z+2*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 4) sqrt(-I*(sqrt(a)+I*z) +
> 2*sqrt(-I)*sqrt(sqrt(a)+I*z)*sqrt(I)*sqrt(sqrt(a)-I*z) +
> I*(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 5) sqrt((sqrt(-I)*sqrt(sqrt(a)+I*z) + sqrt(I)*sqrt(sqrt(a)-I*z))^2)/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 6) (sqrt(-I)*sqrt(sqrt(a)+I*z) + sqrt(I)*sqrt(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 7) sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*z)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*z))
>
> Expression 2) above is based on the identity
>
> sqrt(a+z^2) = sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z)
>
> which is valid for all z provided a>0. The only questionable step
> above is from 5) to 6). It assumes sqrt(u^2) equals u, which is only
> true if the real part of u is positive. So without loss of
> generality, all that remains to be proved is that the real part of
>
> sqrt(-I)*sqrt(1+I*z) + sqrt(I)*sqrt(1-I*z)
>
> is positive for all z. If any reader can prove it, please let me
> know.
>
> Substituting sqrt(b)*x^2 for z in the generalized identity above
> yields the special case
>
> sqrt(sqrt(b)*x^2+sqrt(a+b*x^4))/sqrt(a+b*x^4) =
> sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*sqrt(b)*x^2)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*sqrt(b)*x^2))
>
> which is valid provided a>0. Integrating both sides and translating
> it into Mathematica syntax yields the nice new rule
>
> Int[u_.*Sqrt[c_.*x_^2+Sqrt[a_.+b_.*x_^4]]/Sqrt[a_.+b_.*x_^4],x_] :=
> Sqrt[-I]/Sqrt[2]*Int[u/Sqrt[Sqrt[a]-I*c*x^2],x]] +
> Sqrt[ I]/Sqrt[2]*Int[u/Sqrt[Sqrt[a]+I*c*x^2],x]] /;
> FreeQ[{a,b,c},x] && Sqrt[b]===c && a>0
>
> So now Rubi is able to integrate Martin’s problem 8, and a lot more.
> For example,
>
> sqrt(2*x^2+sqrt(3+4*x^4))/((c+d*x)*sqrt(3+4*x^4))
>
> and
>
> sqrt(2*x^2+sqrt(3+4*x^4))/((c+d*x)^2*sqrt(3+4*x^4))
>
> This whole exercise has been typical of the way knowledge is added to
> a rule-based system, like Rubi: Generalizing the procedure required to
> solve a specific problem, yields a general rule able to handle a whole
> class of problems. I hope you have found it instructive.
>
Experiments of mine confirm that SQRT(a^2+z^2) = SQRT(a + #i*z) * SQRT(a
- #i*z) requires a > 0, and doesn't permit a larger domain. I think this
particular relation (in suitably generalized form) should be available
both ways, to help matching a given integrand to rules in Rubi's
collection. The same applies to SQRT(a^2-z^2) = SQRT(a+z) * SQRT(a-z),
provided this was indeed needed as a separate rule (where again a > 0).
It appears that Rubi's newly crafted rule involving the fourth power x^4
covers the case resulting for x <- #i*x automatically. But does it also
cover the simpler quadratic case for x^2 <- x and x^2 <- -x (unless
never needed), or that of arbitrary functions x^2 <- f(x) (unless
completely useless)? How about related root splitting or denesting
rules? Some candidate relations:
SQRT(SQRT(1+z^2) + 1) = (SQRT(1 + #i*z) + SQRT(1 - #i*z)) / SQRT(2)
1 / SQRT(SQRT(1+z^2) + 1) = (SQRT(1 + #i*z) - SQRT(1 - #i*z)) /
(SQRT(2)*#i*z)
SQRT(SQRT(1+z^2) - 1) = SQRT(z^2) (SQRT(1 + #i*z) - SQRT(1 - #i*z)) /
(SQRT(2)*#i*z)
1 / SQRT(SQRT(1+z^2) - 1) = (SQRT(1 + #i*z) + SQRT(1 - #i*z)) /
SQRT(2*z^2)
Needless to say, these can be generalized in similar ways (z <- x^2,
z <- #i*x, z <- #i*x^2, z <- SQRT(#i)*x, z <- SQRT(#i)*x^2, etc.).
Finally, wouldn't it be more profitable to concentrate on the (suitably
generalized) true root splitting (or denesting) step
SQRT(SQRT(1+z^2) + z) = 1 / SQRT(SQRT(1+z^2) - z) =
(SQRT(-#i)*SQRT(1 + #i*z) + SQRT(#i)*SQRT(1 - #i*z)) / SQRT(2),
as a new Rubi rule, and let everything else follow automatically via
other rules? But I haven't looked at your code.
Martin.
PS: If you want sci.math.symbolic readers to suggest new rules, it would
be helpful if your site could present the tests with unsatisfactory or
missing Rubi answers in a form readable without access to MMA (or
some MMA "player"). Or perhaps you could post a selection of the more
"promising" failures here?
PPS: I didn't know that the Risch integrator of Axiom (or FriCAS) was
this abysmally slow (it's been running for _ten_ days now). In 1984 (the
year of Trager's thesis), the code must have been near impossible to
debug even on a Cray!
Waldek Hebisch
>
> PPS: I didn't know that the Risch integrator of Axiom (or FriCAS) was
> this abysmally slow (it's been running for _ten_ days now). In 1984 (the
> year of Trager's thesis), the code must have been near impossible to
> debug even on a Cray!
Axiom code is mostly due to Manuel Bronstein in 1990-1994 (I do not
think is make much difference, in 1984 Axiom team should have access
to 30 MIPS mainframe (maybe better) and in 1994 Bronstein probably
had 30-100 MIPS workstation).
Actually, the Risch integrator itself runs reasonably fast. The
problem is arithmetic with expressions containing two square
roots. Risch integrator creates intermediate expressions which are
significantly larger than input expression -- not very large but
enough to make arithmetic very slow. I think that all routines
were debugged using single algebraics (in such case arithmetic is
reasonably fast).
BTW1: Arithmetic on algebraic expression can be done much faster
then currently in FriCAS. New routines are under developement.
BTW2: Currently Risch integrator makes some shortcuts for case
of roots (as opposed to general algebraics). But I think that
in case of two roots it is possible to make shortcut doing
only slightly more computatins as for single root, while
current code treats roots in sequence doing much more
expensive operations.
--
Waldek Hebisch
heb...@math.uni.wroc.pl
clicl...@freenet.de
Albert schrieb:
>
> Is it your desire to have the Mathematica and Rubi test suite results
> expressed in a form readable without Mathematica? If so, I can also
> make them available as pdf files.
>
This is what I was getting at. These files are awfully hard to read
without access to Mathematica; even Euler would have despaired. Wide
dissemination of Rubi's failures could help increase feedback in the
form of suggestions, which might lead to improved rules.
Martin.
Nasser M. Abbasi
Fyi;
I've added Sage 4.4.4 (downloaded on july 2,2010) result to the above
page (the link to Sage worksheet output can be found near the top of the
above page) it is an HTML page.
Sage did only integral #9.
I am not a Sage user, so may be someone with more Sage experience can
try some other tricks or commands to get a better result than what I
got. I used sage integral() command with no assumptions. Run sage from
inside VM machine on windows.
--Nasser
Nasser M. Abbasi
> integrator...
Sage itself is open source, so adding open source code to it should not
be a problem.
For integration, using integral(), it seems to use Maxima interface
based on my reading of the top of this page:
http://www.sagemath.org/doc/constructions/calculus.html
>Does it support a pattern-matching programming language?
>
> Albert
Some pattern commands on expressions can be seen here:
http://www.sagemath.org/doc/reference/sage/symbolic/expression.html
Since sage is phython based, I also assume that any pattern matching
functionality in python or sympy can be used from within sage as well.
I know very little about sage, every once in a while I download the
latest version and play with it a little and that is all.
--Nasser
clicl...@freenet.de
"Nasser M. Abbasi" schrieb:
>
> > On Jul 2, 9:18 pm, "Nasser M. Abbasi"<n...@12000.org> wrote:
> >
> >> Sage did only integral #9.
> >
>
>
> For integration, using integral(), it seems to use Maxima interface
> based on my reading of the top of this page:
>
> http://www.sagemath.org/doc/constructions/calculus.html
>
That Maxima is used seems to follow from the output for problem #4 as
well:
#4 f=sqrt(x^2-1)/(x-I)^2 f.integral(x)
Traceback (click to the left of this block for traceback)
...
Is (I-1)*(I+1) zero or nonzero?
This kind of question to the operator is typical of the Maxima
integrator. Nobody will be surprised to learn that Derive 6.10
simplifies IDENTICAL?((#i-1)*(#i+1), 0) to FALSE automatically, no
questions asked. Perhaps your "I" was not really #i, the imaginary unit
- in spite of the definition at the top of your page?
Martin.
Richard Fateman
I too understand that Sage uses Maxima for integration, if so, if you
want to test integrals, you can just download and use Maxima.
Maxima uses %i for sqrt(-1), not I.
and (%i-1)*(%i+1) expands to -2, as expected.
And if you want to add rules to the integration procedure in Maxima,
you can do so using procedures or patterns in Maxima. There are
several pattern matchers, but the one that is only incompletely
accessible in Maxima is the one probably of most interest, which is the
one from MockMMA -- it looks like Mathematica's.
I think that Sage and Python are, in this case, merely a distraction.
You might try one of the Axiom variants, too.
As for the possibility of expressing integrals without algebraic
extensions, I think that Barry Trager's MIT thesis describes the
conditions under which this works. It looks like one direction
to head is to encode the Risch algorithm as a side-light to certain
rules.
RJF
Nasser M. Abbasi
> Perhaps your "I" was not really #i, the imaginary unit
> - in spite of the definition at the top of your page?
>
> Martin.
Thanks. It looks like you are correct. I corrected #4 and just used
sqrt(-1), and now sage does this problem. So sage does #4 and #9 now.
I updated the page
http://12000.org/my_notes/rubiTest/rubi_1/rubi.html
#4
I=sqrt(-1)
x=var('x')
f=sqrt(x^2-1)/(x-I)^2
f.integral(x)
(answer below)
I also plotted difference of sage answer against Mathematica's answer
Zero plot. good. at least numerically same answer.
Btw, for this problem, sage (or maxima really) has the smallest leaf
count. Using Mathematica LeafCount on sage answer, Rubi and Mathematica
own answers, I get
LeafCount[sageAnswer]
69
LeafCount[mmaAnswer]
154
LeafCount[rubiAnswer]
136
sageAnswer = (1/2)*I*Sqrt[2]*ArcSin[I*(x/Sqrt[x^2 + 1]) - 1/Sqrt[x^2 +
1]] - Sqrt[x^2 - 1]/(x - I) + Log[2*x + 2*Sqrt[x^2 - 1]];
mmaAnswer = (1/4)*(-((4*(-1 + x^2)^(1/2))/(-I + x)) -
2*I*2^(1/2)*ArcTan[(1/2)*(-I + x - 2^(1/2) - 1 + x^2^(1/2))] +
4*ArcTanh[(2*x)/(I - x + (-1 + x^2)^(1/2))] - 2^(1/2)*Log[-I + x]
+ 2^(1/2)*Log[2*I + 6*x - 4*2^(1/2) - 1 + x^2^(1/2)] +
2*Log[-4 + 8*x^2 - 8*I - 1 + x^2^(1/2) + 8*x*(-I - 1 + x^2^(1/2))]);
rubiAnswer = 4/(1 - 2*I*(x + Sqrt[-1 + x^2]) + (x + Sqrt[-1 + x^2])^2) -
(I*(-2*I + 2*(x + Sqrt[-1 + x^2])))/(1 - 2*I*(x + Sqrt[-1 + x^2]) +
(x + Sqrt[-1 + x^2])^2) + I*Sqrt[2]*ArcTan[(-2*I + 2*(x + Sqrt[-1 +
x^2]))/ (2*Sqrt[2])] + Log[x + Sqrt[-1 + x^2]];
--Nasser
clicl...@freenet.de
Albert schrieb:
>
> On Jul 1, 2:12 am, cliclic...@freenet.de wrote:
>
> >> Expression 2) above is based on the identity
> >> sqrt(a+z^2) = sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z)
> >> which is valid for all z provided a>0.
> >
> > SQRT(a^2+z^2) = SQRT(a + #i*z) * SQRT(a - #i*z)
> > requires a > 0, and doesn't permit a larger domain. I think this
> > particular relation (in suitably generalized form) should be available
> > both ways, to help matching a given integrand to rules in Rubi's
> > collection. The same applies to
> > SQRT(a^2-z^2) = SQRT(a+z) * SQRT(a-z),
> > provided this was indeed needed as a separate rule (where again a > 0).
> >
> > It appears that Rubi's newly crafted rule involving the fourth power x^4
> > covers the case resulting for x <- #i*x automatically. But does it also
> > cover the simpler quadratic case for x^2 <- x and x^2 <- -x (unless
> > never needed), or that of arbitrary functions x^2 <- f(x) (unless
> > completely useless)?
>
>
> Int[u_.*Sqrt[v_+Sqrt[a_+w_]]/Sqrt[a_+w_],x_Symbol] :=
> (1-I)/2*Int[u/Sqrt[Sqrt[a]-I*v],x]] +
> (1+I)/2*Int[u/Sqrt[Sqrt[a]+I*v],x]] /;
> FreeQ[a,x] && w==v^2 && a>0
>
> It is based on the previously discussed identity
>
> sqrt(z+sqrt(a+z^2))/sqrt(a+z^2) =
> (1-I)/(2*sqrt(sqrt(a)-I*z)) + (1+I)/(2*sqrt(sqrt(a)+I*z))
>
> modified to use the fact that sqrt(I/2)=(1+I)/2. Since replacing z
> with I*z yields
>
> sqrt(I*z+sqrt(a-z^2))/sqrt(a-z^2) =
> (1-I)/(2*sqrt(sqrt(a)+z)) + (1+I)/(2*sqrt(sqrt(a)-z)),
>
> there is no need to have a separate rule for the sqrt(a^2-z^2) case.
>
> This rule also covers "the simpler quadratic case" like you
> requested. However, Rubi already uses Euler's substitution for this
> case. Unlike the above algebraic expansion, Euler's substitution does
> not introduce the imaginary unit, but it does make for a more
> difficult derivation and a bigger antiderivative. The pdf file
>
> http://www.apmaths.uwo.ca/~arich/SubstitutionVsExpansion.pdf
>
> provides examples of this by showing the steps Rubi uses to integrate
>
> sqrt(2*x+sqrt(3+4*x^2))/sqrt(3+4*x^2)
> and
> sqrt(2*x+sqrt(3+4*x^2))/((c+d*x)*sqrt(3+4*x^2))
>
> first using Euler's substitution, and then using algebraic expansion.
I have difficulties matching the "Euler substitution" displayed in the
pdf file with what seems to be actually going on. Could the wrong rule
be displayed for both integrands? Or am I having a bad day? I think
there is only one way of matching SQRT(a+b*x+c*x^2) with the integrand,
which should result in a=3, b=0, c=4. The actual transformation,
however, apparently just takes the argument of the outer root as the new
integration variable. This procedure works with problem #8 as well: it
is the solution I was originally aware of (I believe it was once shown
on this group by Oleksandr at WRI), where the resulting antiderivative
is larger but doesn't involve the imaginary unit. For the general user,
this would often be preferable, I think.
>
> Euler's ability to avoid introducing the imaginary unit makes me
> wonder whether there is an analogous substitution for integrands
> having subexpressions of the form sqrt(a+b*x^4), thus obviating the
> need for our algebraic expansion rule. If anyone knows of such a
> substitution, please let me know.
>
All that seems to be needed would be to extend the range of
applicability of your secret substitution rule to a larger class of
integrands, namely those involving the denestable radicands
SQRT(1+x^4)+x^2, SQRT(1+x^4)+1, SQRT(1+x^4)-1 (all of them suitably
generalized), and perhaps even many many more, and to see if Rubi can
manage then. Right?
Martin.
Valeri Astanoff
> On Jun 26, 6:31 am, cliclic...@freenet.de wrote:
>
> > elegant"; the validity of the numerator transformation may not be
> > obvious to many.
>
> suggested, I generalized your identity to
>
> sqrt(z+sqrt(a+z^2))/sqrt(a+z^2) =
> sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*z)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*z))
>
> which appears to be valid for all z provided a>0. Its derivation
> tracks yours as follows
>
> 1) sqrt(z+sqrt(a+z^2))/sqrt(a+z^2)
>
> 2) sqrt(z+sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))/
> (sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 3) sqrt(2*z+2*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 4) sqrt(-I*(sqrt(a)+I*z) +
> 2*sqrt(-I)*sqrt(sqrt(a)+I*z)*sqrt(I)*sqrt(sqrt(a)-I*z) +
> I*(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 5) sqrt((sqrt(-I)*sqrt(sqrt(a)+I*z) + sqrt(I)*sqrt(sqrt(a)-I*z))^2)/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 6) (sqrt(-I)*sqrt(sqrt(a)+I*z) + sqrt(I)*sqrt(sqrt(a)-I*z))/
> (sqrt(2)*sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z))
>
> 7) sqrt(-I)/(sqrt(2)*sqrt(sqrt(a)-I*z)) +
> sqrt(I)/(sqrt(2)*sqrt(sqrt(a)+I*z))
>
>
> sqrt(a+z^2) = sqrt(sqrt(a)+I*z)*sqrt(sqrt(a)-I*z)
>
> above is from 5) to 6). It assumes sqrt(u^2) equals u, which is only
> generality, all that remains to be proved is that therealpartof
>
> sqrt(-I)*sqrt(1+I*z) + sqrt(I)*sqrt(1-I*z)
>
> is positive for all z. If any reader can prove it, please let me
> know.
Good day,
I tried and show (with some help from Mathematica)
that the real part of
f(z) = sqrt(-I)*sqrt(1+I*z) + sqrt(I)*sqrt(1-I*z)
was positive :
In[1]:= f[z_] = Sqrt@(-I)*Sqrt@(1 + I*z) + Sqrt@(I)*Sqrt@(1 - I*z);
Let r(a,b) be the real part of f(a + i b)
In[2]:= r[0, 0] = Sqrt[2]; r[0, -1] = 1; r[0, 1] = 1;
In[3]:= r[a_, 0] = Sqrt[2*(a + Sqrt[a^2 + 1])];
In[4]:= r[0, b_ /; b < -1] := (-Sqrt[-1 - b] + Sqrt[1 - b])/Sqrt[2];
Reduce[(-Sqrt[-1 - b] + Sqrt[1 - b])/Sqrt[2] > 0, {a, b}, Reals] //
Simplify[#, b < -1] &
Out[5]= True
In[6]:= r[0, b_ /; -1 < b < 1] := (Sqrt[1 - b] + Sqrt[1 + b])/
Sqrt[2];
Reduce[-1 < b < 1 && (Sqrt[1 - b] + Sqrt[1 + b])/Sqrt[2] > 0, {a, b},
Reals] //
Simplify[#, -1 < b < 1] &
Out[7]= True
In[8]:= r[0, b_ /; b > 1] := (Sqrt[-1 + b] + Sqrt[1 + b])/Sqrt[2];
Reduce[b > 1 && (Sqrt[-1 + b] + Sqrt[1 + b])/Sqrt[2] > 0, {a, b},
Reals] //
Simplify[#, b > 1] &
Out[9]= True
In[10]:= r[a_, b_ /; -1 < b < 1] := Sqrt[a/2 + (1/2)*Sqrt[1 + a^2 -
2*b + b^2]] +
Sqrt[a/2 + (1/2)*Sqrt[1 + a^2 + 2*b + b^2]];
Reduce[-1 < b < 1 && Sqrt[a/2 + (1/2)*Sqrt[1 + a^2 - 2*b + b^2]] +
Sqrt[a/2 + (1/2)*Sqrt[1 + a^2 + 2*b + b^2]] > 0,
{a, b}, Reals] // Simplify[#, -1 < b < 1] &
Out[11]= True
In[12]:= r[a_, b_ /;
Not[-1 < b < 1]] := (If[a < 0 && b > 0, -1, 1]*
Sqrt[a + Sqrt[a^2 + (b - 1)^2]] +
If[a < 0 && b < 0, -1, 1]*Sqrt[a + Sqrt[a^2 + (b + 1)^2]])/
Sqrt[2];
Reduce[(If[a < 0 && b > 0, -1, 1]*Sqrt[a + Sqrt[a^2 + (b - 1)^2]] +
If[a < 0 && b < 0, -1, 1]*Sqrt[a + Sqrt[a^2 + (b + 1)^2]])/
Sqrt[2] > 0,
{a, b}, Reals] // Simplify
Out[13]= True
Testing the 'r' formulas with a few numerical values :
In[14]:= t = Table[{a, b, r[a, b] // N, Re[N[f[a + b*I]]]},
{a, -2, 2, 1/2}, {b, -2, 2, 1/2}] // Flatten[#, 1] & ;
In[15]:= t[[All, -1]] - t[[All, -2]] // Chop
Out[15]= {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0,
0, 0, 0, 0}
Seems ok
--
V.Astanoff
clicl...@freenet.de
Albert schrieb:
>
> On Jul 4, 2:31 am, cliclic...@freenet.de wrote:
>
> > I have difficulties matching the "Euler substitution" displayed in the
> > pdf file with what seems to be actually going on. Could the wrong rule
> > be displayed for both integrands? Or am I having a bad day? I think
> > there is only one way of matching SQRT(a+b*x+c*x^2) with the integrand,
> > which should result in a=3, b=0, c=4. The actual transformation,
> > however, apparently just takes the argument of the outer root as the new
> > integration variable. This procedure works with problem #8 as well: it
> > is the solution I was originally aware of (I believe it was once shown
> > on this group by Oleksandr at WRI), where the resulting antiderivative
> > is larger but doesn't involve the imaginary unit. For the general user,
> > this would often be preferable, I think.
>
> displayed the wrong variant of Euler's substitution, although Rubi
> itself was using the correct variant. The revised pdf file
>
> http://www.apmaths.uwo.ca/~arich/SubstitutionVsExpansion.pdf
>
> now displays the correct substitution. The "Rubi Package zip file" on
> the website has been updated as well.
Oh. I am relieved. My brain is still in working order.
>
> >> Euler's ability to avoid introducing the imaginary unit makes me
> >> wonder whether there is an analogous substitution for integrands
> >> having subexpressions of the form sqrt(a+b*x^4), thus obviating the
> >> need for our algebraic expansion rule. If anyone knows of such a
> >> substitution, please let me know.
> >
> > All that seems to be needed would be to extend the range of
> > applicability of your secret substitution rule to a larger class of
> > integrands, namely those involving the denestable radicands
> > SQRT(1+x^4)+x^2, SQRT(1+x^4)+1, SQRT(1+x^4)-1 (all of them suitably
> > generalized), and perhaps even many many more, and to see if Rubi
> > can manage then. Right?
>
> As the above pdf file illustrates, Euler's substitution transforms
> subexpressions of the form sqrt(a+b*x+c*x^2) to rational functions of
> the form (d+e*x+f*x^2)/(g+h*x) free of squareroots.
>
> Can you propose an analogous substitution that transforms
> subexpressions of the form sqrt(a+b*x^4) to rational functions free of
> squareroots?
>
I think this is not possible. As the algebraic denesting showed,
integrands containing any combination of SQRT(SQRT(1+x^4)+x^2),
SQRT(SQRT(1+x^4)-x^2), SQRT(SQRT(1+x^4)+1), SQRT(SQRT(1+x^4)-1), and
SQRT(1+x^4) are equivalent to integrands containing two square roots of
quadratics only, namely SQRT(1+#i*x^2) and SQRT(1-#i*x^2). Either one of
them can be "transformed away", but the other becomes worse in the
process. Problem #8 can be integrated because the denested integrand
separates into two addditive terms, each of which contains just one kind
of these two roots.
For nested roots involving an inner SQRT(1+x^2), Rubi's Euler
substitution means that the outer radicand became the new integration
variable, t, and the transformed integrand involves a single square root
only, SQRT(t), which is easily dealt with. For the above nested roots
which involve an inner SQRT(1+x^4), the equivalent transformation would
be to again make the outer radicand the new variable. The transformed
integrand then contains two square roots, SQRT(t) and SQRT(t^2-1), and
its integration is straightforward if it separates into two terms where
each contains one kind of root only. This behavior closely parallels
that for the denesting approach.
So, if we want to avoid the imaginary unit in the antiderivative, we
should take t = SQRT(1+x^4)+x^2 for the new integration variiable, and
thus have to put x <- SQRT(t^2-1)/SQRT(2*t) and dx <- (t^2+1)/
(2*t*SQRT(2*t)*SQRT(t^2-1)) dt. Accordingly, we put SQRT(1+#i*x^2) <-
SQRT(-#i)*(1+#i*t)/SQRT(2*t), SQRT(1-#i*x^2) <- SQRT(#i)*(1-#i*t)/
SQRT(2*t), and therefore SQRT(x^4+1) <- (t^2+1)/(2*t),
SQRT(SQRT(x^4+1)+x^2) <- SQRT(t), SQRT(SQRT(x^4+1)-x^2) <- 1/SQRT(t),
SQRT(SQRT(x^4+1)+1) <- (t+1)/SQRT(2*t), and so long as nobody really
notices, SQRT(SQRT(x^4+1)-1) <- (t-1)/SQRT(2*t). Thus, a single
transformation is good for any combination of our four kinds of nested
roots.
Like so many successful 18th century mathematicians, we have here
blithely ignored all branch cuts, in order to expose the two square
roots SQRT(t) and SQRT(t^2-1); we may therefore expect our
transformation to be valid in some part of the complex t plane only,
corresponding to some part of the complex x plane only. Once we have
separated the integrand into two parts containing only one kind of root,
and then laboriously performed the integration, we replace t <-
SQRT(1+x^4)+x^2, not forgetting SQRT(t) <- SQRT(SQRT(1+x^4)+x^2) and
SQRT(t^2-1) <- SQRT(2)*x*SQRT(SQRT(1+x^4)+x^2), and finally try to
bristle with confidence that our hard-won result simply cannot but
provide the correct analytic continuation from the restricted domain of
validity to the entire complex x plane automatically.
Note in this context that SQRT(SQRT(x^4+1)+x^2), SQRT(SQRT(x^4+1)-x^2),
and SQRT(SQRT(x^4+1)+1) all possess branch cuts radiating outward from
the four roots at x = SQRT(#i), x = -SQRT(#i), x = SQRT(-#i), and x =
-SQRT(-#i), which leave the complex plane connected around the origin.
By contrast, SQRT(SQRT(x^4+1)-1) cuts the plane into four sectors, on
whose the mutual boundaries the function appears to be non-analytic.
Martin.
Vladimir Bondarenko
> On Jul 2, 9:18 pm, "Nasser M. Abbasi" <n...@12000.org> wrote:
>
>
>
> Albert
Hello Albert,
The VM machine says the most powerful open source indefinite
integrator belongs to Axiom, the second place keeps Maxima.
I am not quite sure if the Axiom group would consider non-
algorithmic approach compatible with their long-term goals.
Best wishes,
Vladimir
clicl...@freenet.de
Albert schrieb:
>
> On Jul 5, 1:25 pm, cliclic...@freenet.de wrote:
>
> > So, if we want to avoid the imaginary unit in the antiderivative, we
> > should take t = SQRT(1+x^4)+x^2 for the new integration variiable, and
> > thus have to put x <- SQRT(t^2-1)/SQRT(2*t) and dx <- (t^2+1)/
> > (2*t*SQRT(2*t)*SQRT(t^2-1)) dt. Accordingly, we put SQRT(1+#i*x^2) <-
> > SQRT(-#i)*(1+#i*t)/SQRT(2*t), SQRT(1-#i*x^2) <- SQRT(#i)*(1-#i*t)/
> > SQRT(2*t), and therefore SQRT(x^4+1) <- (t^2+1)/(2*t),
> > SQRT(SQRT(x^4+1)+x^2) <- SQRT(t), SQRT(SQRT(x^4+1)-x^2) <- 1/SQRT(t),
> > SQRT(SQRT(x^4+1)+1) <- (t+1)/SQRT(2*t), and so long as nobody really
> > notices, SQRT(SQRT(x^4+1)-1) <- (t-1)/SQRT(2*t). Thus, a single
> > transformation is good for any combination of our four kinds of nested
> > roots.
>
> integrate
>
> sqrt(x^2+sqrt(1+x^4))/sqrt(1+x^4)
>
> yields the integrand 1/(sqrt(2)*sqrt(-1+t^2)). Integrating wrt t
> yields
>
> arctanh(sqrt(-1+t^2)/t)/sqrt(2)
>
> Back substitution and simplification yields
>
> arctanh(sqrt(2)*x^2/sqrt(x^2*(x^2+sqrt(1+x^4))))/sqrt(2)
>
> This differentiates to
>
> sqrt(x^2*(x^2+sqrt(1+x^4)))/(x*sqrt(1+x^4))
>
> which equals the original integrand only for x>0.
That's why I explicitely prescribed (in the paragraph following the one
quoted above) that SQRT(t^2-1) be replaced by SQRT(2)*x*
SQRT(SQRT(1+x^4)+x^2). In your case too, this recipe immediately
provides a valid antiderivative for the entire complex plane,
ATANH(SQRT(2)*x*SQRT(SQRT(1+x^4)+x^2)/(SQRT(1+x^4)+x^2))/SQRT(2). Can
you come up with a counterexample where it doesn't? Perhaps somebody
could even prove its validity in general?
This kind of difficulty is the price one pays for playing fast and loose
with square roots. And this difficulty was on my mind when I labelled
the complex splitting approach "simple and elegant".
>
> However, this is easily remedied by multiplying it by the piecewise
> constant
>
> sqrt(x^2*(x^2+sqrt(1+x^4)))/(x*sqrt(x^2+sqrt(1+x^4)))
>
> whose derivative is 0. Mathematica's FullSimplify command miraculously
> simplifies this product to the elegant antiderivative
>
> arctanh(sqrt(2)*x/sqrt(x^2+sqrt(1+x^4)))/sqrt(2)
>
> which easily differentiates back to the original integrand!
The simplified MMA result is trivially related to the correct
antiderivative given above. In the simplification, MMA's simplifier has
presumably made use of ATANH(-z) = -ATANH(z).
>
> The algebraic expansion we discussed before produces the
> antiderivative
>
> -I*arcsin((1+I)*x/sqrt(2))/sqrt(2) + arcsinh((1+I)*x/sqrt(2))/sqrt(2)
>
> involving the imaginary unit. Mathematica produces a truly scary
> result for this integrand. Although both are valid antiderivatives,
> Mathematica is unable to differentiate either of them back to the
> original integrand.
>
> Generalizing on the elegant result above, the antiderivative of
>
> sqrt(b*x^2 + sqrt(a+b^2*x^4))/sqrt(a+b^2*x^4)
> is
> arctanh(sqrt(2*b)*x/sqrt(b*x^2 + sqrt(a+b^2*x^4)))/sqrt(2*b)
>
> Also the antiderivative of
>
> sqrt(-b*x^2 + sqrt(a+b^2*x^4))/sqrt(a+b^2*x^4)
> is
> arctan(sqrt(2*b)*x/sqrt(-b*x^2 + sqrt(a+b^2*x^4)))/sqrt(2*b)
>
> Note that Martin's algebraic expansion technique cannot integrate
> either of the above integrands, since it requires knowing a>0. Also
> Mathematica and Maple cannot integrate either of them, despite their
> use of the Risch algorithm...
>
> I cannot find either of these antiderivatives in published tables. Can
> you?
>
So, you are now going to teach Rubi a correspondingly generalized
version of the audacious substitution x <- SQRT(t^2-1)/SQRT(2*t), dx <-
(t^2+1)/(2*t*SQRT(2*t)*SQRT(t^2-1)) dt, one that can deal with
SQRT(SQRT(a)+#i*b*x^2), SQRT(SQRT(a)-#i*b*x^2), SQRT(a+b^2*x^4),
SQRT(SQRT(a+b^2*x^4)+b*x^2), SQRT(SQRT(a+b^2*x^4)-b*x^2),
SQRT(SQRT(a+b^2*x^4)+SQRT(a)), and perhaps even
SQRT(SQRT(a+b^2*x^4)-SQRT(a)) all in one swat - and without requiring
a>0 anymore?
Can Rubi's Euler substitutions already handle integrands of the type
SQRT(SQRT(a+z^2)+z) with z = b*x+c? Should the new substitution be
enabled to handle z = b*x^2+c*x+d as well?
Martin.
Waldek Hebisch
> On Jul 3, 10:59?am, Albert <Albert_R...@msn.com> wrote:
> > On Jul 2, 9:18?pm, "Nasser M. Abbasi" <n...@12000.org> wrote:
> >
> > > Sage did only integral #9.
> >
> > Perhaps Sage could take advantage of an open-source rule-based
> > integrator... Does it support a pattern-matching programming language?
> >
> > Albert
>
> Hello Albert,
>
> The VM machine says the most powerful open source indefinite
> integrator belongs to Axiom, the second place keeps Maxima.
>
> I am not quite sure if the Axiom group would consider non-
> algorithmic approach compatible with their long-term goals.
>
I believe that implementation of Risch (and Trager, Bronstein, etc)
algorithm can be significantly improved (some current weak
point are outside of algorithm and there are few interesting
twists waiting for implementation). However, I am also interested
in rule-based approach -- at least having implementation in
FriCAS would allow easier comparison. And if rule-based
method is better for some integrals I would use it.
AFAICS the main issue is having Mathematica-like pattern
matcher and porting about 4000 lines of Mathematica
code in "utility routines" that Rubi uses. This may
be problematic without asscess to Mathematica and ATM
my acccess to Mathematica is very limited so currently
I will work on other things.
--
Waldek Hebisch
heb...@math.uni.wroc.pl
Nasser M. Abbasi
> This may
> be problematic without asscess to Mathematica and ATM
> my acccess to Mathematica is very limited so currently
> I will work on other things.
>
fyi,
Assuming you meant by access, as the cost of Mathematica, there is a
"home edition" for $295
http://www.wolfram.com/products/mathematicahomeedition/
It has everything the commercial edition has, but only for home use:
"It is not licensed for commercial, nonprofit, academic, or government use"
So, if that applies to you, that might be an option to get Mathematica
for home use. For students, the price is about $140.
--Nasser
clicl...@freenet.de
Waldek Hebisch schrieb:
FriCAS has been working on
INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
for twenty days now!
Martin.
TimDaly
Axiom is currently working with Albert to import the Rubi
pattern functionality. The latest version of Axiom has
a test suite of several thousand integrands (see the
src/input/rich*.input.pamphlet files) and their current
results from Axiom's integration using Albert's test suite.
Dave Parnas has published a paper proposing that table forms
be used to specify the input/output relationship of a function.
I am trying to classify Alberts rules into table form. With
Parnas's table idea it should be possible to specify the exact
input/output of portions of Axiom's integration routines. This
will greatly enhance both the documentation and testing. (Vlad
will be happy :-) )
I see this effort as not only compatible with Axiom but
actually enhancing the long term goals of documentation,
verification, and validation as well as having a stronger
integration facility.
Of course, it is a large effort and will take several releases
to accomplish. Just converting the integrals into Axiom regression
test suites took a month.
In fact, Axiom has a Computer Algebra Test Suite (CATS) subgoal.
There are two full test suites, one on Schaums integral set
and one on the Kamke differential equations set.
(http://axiom-developer.org/axiom-website/CATS/index.html)
Eventually Albert's test suite will be added.
Tim Daly
Richard Fateman
not be obvious to the casual reader (or programmer).
For example, you may have a rule that works to integrate
exp(n*x) * x^m wrt x.
Mathematica 7 gives
-E^(-n x) (E^x)^n x^(1 + m) (-n x)^(-1 - m) Gamma[1 + m, -n x]
now if we substitute n->0, we should get something like
If m=-1 then log(x) else x^(m+1)/(m+1)
but
mathematica gives
-0^(-1 - m) x^(1 + m) Gamma[1 + m, 0]
and actually, for integral of x^m, it still does not include the log
possibilities.
But back to patterns. If we try to classify the pattern for
exp(n*x)*x^m, we might say that this rule applies to
products of exponentials and powers.
but perhaps it is also a rule that applies to
products of exponentials and symbols like x, where of course x^1 is a
power, but simplified to x.
or maybe it is a rule that applies just to exponentials, since x^0 is a
power, but but simplified to 1, after which the "product" disappears
entirely.
This is all to provide a specific example of a kind of rule that may
crop up fairly often..
Consider a rule that for particular functions f,g,h like (say) exp,
log, sin, power, ...erf, gamma, ... says how to integrate f^n*g^m*h^k ..
given a particular expression to integrate, it might be missing some
power, or even missing some function entirely. This makes finding the
right rules more difficult.
The trick in assembling large numbers of rules is to make sure you don't
have to run them all, but only that (one hopes) very small subset that
is likely or even possible to match.
I think it is great that Rubi may be adapted for free software systems.
RJF
TimDaly
There has been some initial discussions between myself, Albert Rich,
and David Parnas about how to think about this problem.
My current thinking involves a tree-structuring of the rule patterns
where each node in the tree represents one or more pattern
tables and each edge is a function (e.g. division operation, etc).
Walking the input expression walks the tree to find applicable
patterns.
Because a pattern table can be reached from many different paths the
"tree" is actually a graph. Each table specifies an input pattern and
an output pattern and the output pattern points to another node.
So far all I have are "paper studies" using index cards as the tables.
I looked at constructing a Rete pattern match (Charles Forgy, OPS5)
but there are only a small number of rules that can be applied so
there is not much depth. I also don't know quite how to factor the
"meta rules" into a Rete network.
There are a lot of interesting questions that arise, such as the
issue of branch cuts and choice of trigonometric identities that
get applied (which can easily lead to cycles).
The problem also seems ideal for a parallel decomposition because
I can see multiple patterns applied in any situation (e.g. making
various choices for integration by parts).
Special functions are also a question... is it better to leave the
result as sin(x)/x or to choose sinc(x). The first is easier to
combine with other results but the second is shorter to represent.
Albert measures the number of leaf nodes as a measure of simplicity.
Unfortunately, "simplicity" is in the eye of the beholder. Axiom
generally deals with the simplicity issue by making different domains.
Thus in the domain of Polynomials with fractional coefficients the
"simple" form is:
(1/3)*x^2+7
but in the domain of fractions of polynomials with integer
coefficients
we see the same answer as:
(x^2+21)/3
The "simple" form depends on the domain. This raises interesting
questions for some intermediate domains to express patterns rather
than using the general "Expression(Integer)" domain. The same
"pattern table" advocated by Parnas might generate multiple results
depending on the target domain. Axiom is ideal for this but it
requires a lot of thought.
clicl...@freenet.de
Albert schrieb:
>
< [...] I wonder if anybody could help me out with the following:
>
> I would like Rubi to expand integrands of the form x^m*f(x)/(a+b*x^n)
> into a sum whose terms are expressions of the form f(x)/(c+d*x) where
> m and n are integers and 0<=m<n. Currently Rubi can expand such
> expressions for m=0 and m=n-1 using the identities
>
> 1/(a+b*z^n) = Sum[1/(1 - z/((-a/b)^(1/n)*(-1)^(2*k/n))),{k,1,n}] /
> (a*n)
>
> and
>
> z^(n-1)/(a+b*z^n) = Sum[1/(z - (-a/b)^(1/n)*(-1)^(2*k/n)), {k, 1,
> n}] / (b*n)
>
> So what I want is a generalization of these identities for expanding
> z^m/(a+b*z^n) where m and n are integers and 0<=m<n. I know it must
> be easy, but I'm getting lazy in my old age... Of course, the first
> person to provide it will be credited as the author of the
> corresponding rule in the Rubi source.
>
I suppose this (or a trivial variant of it) will do:
z^m/(a+b*z^n) = (-a/b)^((m+1)/n)/(a*n)*
SUM((-1)^(2*m*k/n)/((-a/b)^(1/n)-(-1)^(-2*k/n)*z),k,1,n)
Holmes (after two opium pipes).
Nasser M. Abbasi
May be there should be restrictions on a and b also?
For example, for a=0, there is a divide by zero in the RHS side above.
Also for b=0. May be these can be special cased?
--Nasser
Nasser M. Abbasi
>>
>
>
>> I suppose this (or a trivial variant of it) will do:
>>
>> z^m/(a+b*z^n) = (-a/b)^((m+1)/n)/(a*n)*
>> SUM((-1)^(2*m*k/n)/((-a/b)^(1/n)-(-1)^(-2*k/n)*z),k,1,n)
>>
>> Holmes (after two opium pipes).
>
> May be there should be restrictions on a and b also?
>
> For example, for a=0, there is a divide by zero in the RHS side above.
> Also for b=0. May be these can be special cased?
>
May be with the way the pattern matching work inside Rubi, I think now
that if a was zero or b was zero, then this pattern will not even 'fire'
up to start with, since a different pattern will be used in those case.
So, I do not think one needs to worry about a=0 or b=0 ?
--Nasser
clicl...@freenet.de
Albert schrieb:
>
> Thanks Martin. That altered state of consciousness must have done the
> trick! I ended up using your identity in the form
>
> z^m/(a+b*z^n) =
> q^(m+1)*Sum[(-1)^(2*k*(m+1)/n)/(q*(-1)^(2*k/n) - z), {k,1,n}]/(a*n)
>
> where q=(-a/b)^(1/n). And yes, the host CAS's pattern matcher along
> with the application conditions on the new rule must ensure that
> neither a nor b is 0, that m and n are integers, and 0<=m<n.
>
> For integrands of the form x^m*log(a+b*x)/(c+d*x^n) when n is odd, the
> new rule reduces the number of terms in the antiderivative from 4 n to
> 2 n. Thus its size was also reduced by a factor of about 2. (Note that
> when n is even another rule is used.)
>
> This new rule represents a clear increase in Rubi's knowledge density
> since two rules have been replaced by one, no more complicated than
> either of the originals.
>
I am a bit worried about the imaginary unit appearing in results. This
was explicitly avoided in Derive 6.10:
If n>0 is even and 0<=m<n is an integer, INT(x^m/(a+b*x^n),x) ->
1/(n*(a^(1/n))^(n-(m+1))*(b^(1/n))^(m+1))*SUM(2*ATAN((x*b^(1/n)/a^(1/n
)-COS((2*k-1)/n*pi))/SIN((2*k-1)/n*pi))*SIN((2*k-1)/n*(m+1)*pi)-LN(x^2
*b^(2/n)-2*x*a^(1/n)*b^(1/n)*COS((2*k-1)/n*pi)+a^(2/n))*COS((2*k-1)/n*
(m+1)*pi),k,1,n/2)
If n>0 is odd and 0<=m<n is an integer, INT(x^m/(a+b*x^n),x) ->
1/(n*(a^(1/n))^(n-(m+1))*(b^(1/n))^(m+1))*((-1)^m*LN(a^(1/n)+x*b^(1/n)
)+SUM(2*(ATAN((x*b^(1/n)/a^(1/n)-COS((2*k-1)/n*pi))/SIN((2*k-1)/n*pi))
*SIN((2*k-1)/n*(m+1)*pi))-LN(x^2*b^(2/n)-2*x*a^(1/n)*b^(1/n)*COS((2*k-
1)/n*pi)+a^(2/n))*COS((2*k-1)/n*(m+1)*pi),k,1,(n-1)/2))
If n>0 is even and 0<=m<n is even, INT(x^m/(b*x^n-a),x) ->
1/(n*(a^(1/n))^(n-(m+1))*(b^(1/n))^(m+1))*(LN((x*b^(1/n)-a^(1/n))/(x*b
^(1/n)+a^(1/n)))-SUM(2*ATAN((x*b^(1/n)/a^(1/n)-COS(2*k/n*pi))/SIN(2*k/
n*pi))*SIN(2*k/n*(m+1)*pi)-LN(x^2*b^(2/n)-2*x*a^(1/n)*b^(1/n)*COS(2*k/
n*pi)+a^(2/n))*COS(2*k/n*(m+1)*pi),k,1,n/2-1))
If n>0 is even and 0<=m<n is odd, INT(x^m/(b*x^n-a),x) ->
1/(n*(a^(1/n))^(n-(m+1))*(b^(1/n))^(m+1))*(LN(x^2*b^(2/n)-a^(2/n))-SUM
(2*ATAN((x*b^(1/n)/a^(1/n)-COS(2*k/n*pi))/SIN(2*k/n*pi))*SIN(2*k/n*(m+
1)*pi)-LN(x^2*b^(2/n)-2*x*a^(1/n)*b^(1/n)*COS(2*k/n*pi)+a^(2/n))*COS(2
*k/n*(m+1)*pi),k,1,n/2-1))
If n>0 is odd and 0<=m<n is an integer, INT(x^m/(b*x^n-a),x) ->
1/(n*(a^(1/n))^(n-(m+1))*(b^(1/n))^(m+1))*(LN(x*b^(1/n)-a^(1/n))-(-1)^
m*SUM(2*(ATAN((x*b^(1/n)/a^(1/n)+COS((2*k-1)/n*pi))/SIN((2*k-1)/n*pi))
*SIN((2*k-1)/n*(m+1)*pi))+LN(x^2*b^(2/n)+2*x*a^(1/n)*b^(1/n)*COS((2*k-
1)/n*pi)+a^(2/n))*COS((2*k-1)/n*(m+1)*pi),k,1,(n-1)/2))
Differentiating these replacement rules will give real versions of the
partial-fraction expansion of x^m/(a+b*x^n) and x^m/(b*x^n-a).
Martin.
clicl...@freenet.de
clicl...@freenet.de schrieb:
>
> FriCAS has been working on
>
> INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
>
> for twenty days now!
>
FriCAS has been working on
INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
for thirty days now!
Martin.
Vladimir Bondarenko
Guinness record book application :)
As for Rubi, we usually kill the tasks
after 100,000 seconds maximally. But
often Rubi kills Mathematica quicker.
Cheers,
Vladimir
On Jul 21, 9:08 pm, cliclic...@freenet.de wrote:
> cliclic...@freenet.de schrieb:
Nasser M. Abbasi
>
> FriCAS has been working on
>
> INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
>
> for thirty days now!
>
> Martin.
Martin,
How do you know it is actually doing something useful and not stuck in a
loop?
Does FriCAS has an indicator saying it is actually doing different
things? Do you have tracing on?
What possibly it can be doing in these 30 days? trying what? is risch
algorithm used here?
Let us know when and if it completes ;)
--Nasser
Rob Johnson
Vladimir Bondarenko <v...@cybertester.com> wrote:
>On Jun 4, 6:18 am, Albert <Albert_R...@msn.com> wrote:
>> I would like to announce the launch of the website
>>
>> www.apmaths.uwo.ca/RuleBasedMathematics
>>
>> It is dedicated to dedicated to demonstrating the numerous advantages
>> of the rule-based approach to automating mathematics. In systems
>> implemented using this approach, rules are expressed as elegant
>> mathematical formulas, rather than embedded in conventional program
>> code.
>>
>> As proof-of-concept, I have implemented an efficient and robust Rule-
>> based Integrator, nicknamed Rubi. Not only can Rubi compute the
>> antiderivative for a broad class of integrands, but the results are
>> often significantly superior to those produced by the commercial
>> computer algebra systems.
>>
>> The 1500 or so integration rules Rubi uses are freely available on the
>> website in both human and machine readable form. Also available is a
>> test suite of over 9400 integration problems developed in conjunction
>> with the rules. After reviewing the homepage, I recommend clicking on
>> "Highlights of the Indefinite Integration Test Results" for an eye-
>> opening comparison of the rule-based integrator (Rubi) with
>> Mathematica's and Maple's built-in integrators.
>>
>> Also if you have access to Mathematica 6 or better, there is a link
>> near the bottom of the homepage to download Rubi so you can verify the
>> results for yourself. Rubi also provides the option to show the rules
>> required to integrate expressions, along with the intermediate
>> results. I think this show-step ability has great potential
>> pedagogical and research value.
>>
>> Currently the website is pretty Spartan in format, and limited to
>> indefinite integration. However with the help of the computer algebra
>> community, I hope it evolves into a true repository of mathematical
>> knowledge.
>>
>> Aloha from Hawaii,
>> Albert D. Rich
>
>There is much to say about the Rubi.
>
>One tiny example.
>
>The VM machine says:
>
>Mathematica 7.0.1.0 returns
>
> Integrate[1/(Sqrt[2] + Sin[z] + Cos[z]), z]
>
>unevaluated (!).
>
>Rubi returns
>
> -(2/(1 + (-1 + Sqrt[2]) Tan[z/2]))
>
>which is a correct answer.
Integrate[TrigFactor[1/(Sqrt[2] + Sin[z] + Cos[z])], z]
gives
Pi + 4 z
Sqrt[2] Sin[--------]
8
-----------------------------
Pi + 4 z Pi + 4 z
Cos[--------] + Sin[--------]
8 8
which is sqrt[2] + 1 greater than the Rubi answer, but still a valid
integral.
Mathematica has the ability to finish off this integral, but it
simply stops too soon.
Rob Johnson <r...@trash.whim.org>
take out the trash before replying
to view any ASCII art, display article in a monospaced font
clicl...@freenet.de
"Nasser M. Abbasi" schrieb:
>
> On 7/21/2010 11:08 AM, clicl...@freenet.de wrote:
> >
> > FriCAS has been working on
> >
> > INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
> >
> > for thirty days now!
> >
>
>
> Does FriCAS has an indicator saying it is actually doing different
> things? Do you have tracing on?
>
I have to pass these questions on to Waldek Hebisch, who announced this
calculation one month ago - I just provided the integral.
> What possibly it can be doing in these 30 days? trying what? is risch
> algorithm used here?
>
> Let us know when and if it completes ;)
>
Waldek seemed to be confident that FriCAS would ultimately be
successful, but he has warned already about running times reaching weeks
or months. The latter prediction is already borne out!
Are you and Vladimir advocating to show mercy by terminating this
calculation? There are no objections from my side, but it is for Wakdek
to decide.
Martin.
Waldek Hebisch
>
> "Nasser M. Abbasi" schrieb:
> >
> > On 7/21/2010 11:08 AM, clicl...@freenet.de wrote:
> > >
> > > FriCAS has been working on
> > >
> > > INT(SQRT(SQRT(x^4+1)+x^2)/((x+1)^2*SQRT(x^4+1)),x)
> > >
> > > for thirty days now!
> > >
> >
> > How do you know it is actually doing something useful and not stuck in > a loop?
> >
> > Does FriCAS has an indicator saying it is actually doing different
> > things? Do you have tracing on?
> >
>
> I have to pass these questions on to Waldek Hebisch, who announced this
> calculation one month ago - I just provided the integral.
>
> > What possibly it can be doing in these 30 days? trying what? is risch
> > algorithm used here?
> >
The algorithm is mostly from Trager thesis, with some modifications
due to Bronstein. AFAICS the first time sink is computing
determinant of 4x4 matrix with entries beeing polynomials
with coefficients beeing rational functions with algebraic
coefficients. The matrix when printed takes about 16 kB,
so it is not really large. But the problem is in algebraic
coefficients. Normally FriCAS simplifies fractions, computing
GCD of numerator and denomimator. But in case of algebraic
coefficients we (may) loose unique factorization, so there
is no GCD. So FriCAS keeps algeraic fraction only partially
simplified. For example:
1 - sqrt(2) = ((1-sqrt(2))*(1+sqrt(2)))*(1+sqrt(2))^(n-1)/(1+sqrt(2))^n
-(1+sqrt(2))^(n-1)
= --------------------
(1+sqrt(2))^n
FriCAS simplifies here numerator and denominator (separately),
but missess the fact that there is cancellation between
numerator and denominator. Such spurious common factors
tends to increase exponentially with number of operations.
The result is that FriCAS has to compute with really large
numbers (from hundreds to millions of digits).
The problem is made worse by fact that FriCAS uses Bareiss
method to compute determiant. Normally Bareiss is expected
to cancel common factors without computing GCD, but is
this case (algebraics wrappend into rational functions wrapped
into polynomials) it make matters much worse. In fact,
when I replaced use of Bareiss method by expansion into
minors FriCAS can compute the determinant in question in
few seconds.
After that FriCAS has to compute a few GCD-s of polynomials
coming from coefficients of determinant above. Then a few
more GCD-s... GCD-s have similar problems with coefficient
growth (like determinant) but this time there is no simple
workaround (there is a better method, I am working on
implementing it but this is not finished).
> > Let us know when and if it completes ;)
> >
>
> Waldek seemed to be confident that FriCAS would ultimately be
> successful, but he has warned already about running times reaching weeks
> or months. The latter prediction is already borne out!
>
> Are you and Vladimir advocating to show mercy by terminating this
> calculation? There are no objections from my side, but it is for Wakdek
> to decide.
>
I am still think that given sufficient time and memory FriCAS would
be successful. But time is clearly excessive and it seems that
before finishing FriCAS will run out of memory. Namely, in
developement version after using different determinant routine
and with significant speedup to bignum routines after few
minutes I got the first bunch of GCD-s, but then after few
hours FriCAS run out of memory computing next GCD. So
it make sense to kill the calculation.
--
Waldek Hebisch
heb...@math.uni.wroc.pl
Vladimir Bondarenko
the VM machine to test the Albert Rich's Rubi.
http://www.apmaths.uwo.ca/~arich/
Our goal is to see how the human beings compare
against the VM machine, so we urge all the
interested person to publish the Rubi bugs.
Out comment is that the Rubi bug list Albert puts
at his site is, to put it mildly, too short.
Best wishes,
Vladimir Bondarenko
VM and GEMM architect
Co-founder, CEO, Mathematical Director
http://www.cybertester.com/ Cyber Tester, LLC
http://maple.bug-list.org/ Maple Bugs Encyclopaedia
http://www.CAS-testing.org/ CAS Testing
-----------------------------------------------------
"We must understand that technologies
like these are the way of the future."
-----------------------------------------------------
Richard Fateman
...
>
> There has been some initial discussions between myself, Albert Rich,
> and David Parnas about how to think about this problem.
>
> My current thinking involves a tree-structuring of the rule patterns
> where each node in the tree represents one or more pattern
> tables and each edge is a function (e.g. division operation, etc).
> Walking the input expression walks the tree to find applicable
> patterns.
> Because a pattern table can be reached from many different paths the
> "tree" is actually a graph. Each table specifies an input pattern and
> an output pattern and the output pattern points to another node.
>
...
Tilu would take a rule and place it several places in a tree. As an
overly-simplified example, consider
a*cos(x)+b. It appears to have a head operator of "+" and so that is
where the rule would seem to fit in a tree of rules.
However, if b=0, the head operator is "*", and so the pattern a*cos(x)
--> .. would also appear in the tree at "*".
But if a=1, then the rule would appear with head operator "cos".
As I recall, the head operators are not really in a tree but a hash table.
The entries in the hash table may themselves be hash tables, where the
presence or absence of some sub-operators is allowed, in order to narrow
down the number of possibly-appropriate rules.
Note that one may end up in Tilu with some number of matching rules --
with possibly different sets of side conditions, or maybe not. Tilu does
not exclude the possibility of any number of applicable rules giving
perhaps different versions of the answer.
RJF
Sam Blake
In[5296]:= $Version
Out[5296]= "13.2.0 for Mac OS X x86 (64-bit) (November 18, 2022)"
In[5297]:= Integrate[Sqrt[x^2 + Sqrt[1 + x^4]]/((1 + x)^2 Sqrt[1 + x^4]), x]
Out[5297]= 1/2 ((-1 - 2 x^4 - Sqrt[1 + x^4] - x^2 (1 + 2 Sqrt[1 + x^4]))/((1 + x) (x^2 + Sqrt[1 + x^4])^(3/2)) +
ArcTan[Sqrt[1 + Sqrt[2]] Sqrt[x^2 + Sqrt[1 + x^4]]]/Sqrt[-1 + Sqrt[2]] -
Sqrt[1 + Sqrt[2]]ArcTan[(Sqrt[2 (-1 + Sqrt[2])] x Sqrt[x^2 + Sqrt[1 + x^4]])/(
1 + x^2 + Sqrt[1 + x^4])] - ArcTanh[Sqrt[-1 + Sqrt[2]] Sqrt[x^2 + Sqrt[1 + x^4]]]/Sqrt[
1 + Sqrt[2]] + Sqrt[-1 + Sqrt[2]]ArcTanh[(Sqrt[2 (1 + Sqrt[2])] x Sqrt[x^2 + Sqrt[1 + x^4]])/(1 + x^2 + Sqrt[1 + x^4])])
If the time constraint to IntegrateAlgebraic in increased then it can also compute
Integrate[Sqrt[x^2 + Sqrt[1 + x^4]]/((a x + b)^2 Sqrt[1 + x^4]), x]