Re: integration using PSQL algorithm
- To: mathgroup at smc.vnet.net
- Subject: [mg52368] Re: integration using PSQL algorithm
- From: Paul Abbott <paul at physics.uwa.edu.au>
- Date: Thu, 25 Nov 2004 05:49:38 -0500 (EST)
- Organization: The University of Western Australia
- References: <co1dl9$sb8$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
In article <co1dl9$sb8$1 at smc.vnet.net>, Arturas Acus <acus at itpa.lt> wrote: > I have some reasons to suspect that some integral may have closed form > in terms of powers of Pi: a+b*Pi+c*Pi^2+d*Pi^3+e*Pi^4 (a,b,c,d,e being > algebraic numbers, probably with a=d=0). If the coefficients were positive integers less than the base Pi you could just use RealDigits. For example, num = 1 + Pi + 2 Pi^2 + Pi^3 + 2 Pi^4 coefs = RealDigits[N[num, 25], Pi] num == Expand[FromDigits[coefs, Pi]] > Assuming I can calculate the > integral numerically to desired precision, is it possible to verify this > hypothesis. As far as I know, the PSLQ algorithm, which helps to find > algebraic relations between real numbers, can be useful here. The only > (Mathematica related) reference I know is > http://library.wolfram.com/infocenter/MathSource/4263/ > , but the example here is too simple. In The Mathematica Journal 6(2): 29-30, you will find code for TranscendentalRecognize[] that rationalizes the given transcendental basis and then uses rational arithmetic and LatticeReduce to find the "simplest" (rational) representation for the floating point number in that basis: TranscendentalRecognize[n_, basis_] := Module[{c, d, digs, e, id, lat, powerten, r, s, vals}, {d, e} = RealDigits[n]; s = Sign[n]; c = FromDigits[d]; powerten = 10^(Length[d] - e); digs = (RealDigits[N[#1, -e + Length[d] + 5]] & ) /@ basis; r = (FromDigits[Take[First[#1], -e + Last[#1] + Length[d]]] & ) /@ digs; lat = Transpose[Append[IdentityMatrix[Length[basis] + 2], Flatten[{powerten, r, c}]]]; vals = Take[First[LatticeReduce[lat]], Length[basis] + 2]; Expand[-((s*(Take[vals, {2, -2}] . basis + First[vals]))/Last[vals])] ] This code works fine on the example above TranscendentalRecognize[N[num], Pi^Range[4]] but it also works with more general examples such as num = 1/4 - (2/3) Pi + (40/11) Pi^2 + Pi^4; TranscendentalRecognize[N[num, 30], Pi^Range[4]] > As for details, the integral in consideration is > > integral = (-6*(-1 + a)*(4*Sqrt[7 + a] + Sqrt[2]*(5 + a)*Pi))/((5 + > a)*(7 + a)^(3/2)) - > (24*(-1 + a)*(6 + a)*Log[6 + a - Sqrt[5 + a]*Sqrt[7 + a]])/((5 + a)*(7 > + a))^(3/2) > > with a= Cos[4* phi], integrated from 0 to Pi/2 (or to 2 Pi, no matter) > The problems arises from second term. Mathematica can integrate the first term in closed form. The result involves EllipticE[4/3] and EllipticE[Pi/3, 4/3]. The first of these can be written as Hypergeometric2F1[-1/2, 1/2, 1, 4/3] times Pi/2 but there is no simple form for the second. Unless there is cancellation between the two terms, I don't see how algebraic coefficients can result. Working to high precision, Clear[num] num[n_] := num[n] = NIntegrate[Evaluate[integral /. a -> Cos[4 phi]], {phi, 0, Pi/2}, WorkingPrecision -> n] I find that tr = TranscendentalRecognize[num[150], Pi^Range[4]] yields -(402550287036105137464991/240420095219532969056072) - (398516724123798207595259*Pi)/240420095219532969056072 - (177542572056568175522549*Pi^2)/240420095219532969056072 - (502586551534098939521059*Pi^3)/240420095219532969056072 + (196649748041515569157439*Pi^4)/240420095219532969056072 and this agrees with the numerical value to better than 1 part in 10^(142): num[150] - tr The rational coefficients in this basis are coefs = CoefficientList[tr, Pi] However, changing the precision causes the coefficients to vary in such a way that makes me quite confident your suspicion is wrong. Alternatively, you can use << NumberTheory`Recognize` (which itself uses LatticeReduce) to attempt to find simple algebraic numbers for each of the coefficients. The resulting low order polynomials have large coefficients indicating that, again, your suspicion is likely unfounded. Cheers, Paul -- Paul Abbott Phone: +61 8 6488 2734 School of Physics, M013 Fax: +61 8 6488 1014 The University of Western Australia (CRICOS Provider No 00126G) 35 Stirling Highway Crawley WA 6009 mailto:paul at physics.uwa.edu.au AUSTRALIA http://physics.uwa.edu.au/~paul