0
$\begingroup$

I try to compute the coefficients b[i] , $i=2,...,m$. I got the following output:

m := 20; 
μ := 0.5
u[x_] := (12*μ^2)/25 + (6/25)*μ^2* Sech[(x*μ)/5]^2 - (12/25)*μ^2*Tanh[(x*μ)/5]; 

T[n_, x_] := 
  Sum[(-1)^i*2^(n - 2*i - 1)*((n*(n - i - 1)!)/(i!*(n - 2*i)!))*x^(n - 2*i), 
       {i, 0, Ceiling[n/2]}];

Table[b[n] = (2/Pi)*NIntegrate[(u[0.5*x + 0.5]*T[n, x])/Sqrt[1 - x^2], 
          {x, -1, 1}, PrecisionGoal -> 12], {n, 2, m}]

Here the message that before the results.

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in x near {x} = {-0.99999999999999667104594803109799410885580269090806374709815354529}. NIntegrate obtained -0.0001047215230431864 and 7.698445916412437`*^-11 for the integral and error estimates. >>

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in x near {x} = {0.99999999999999667104594803109799410885580269090806374709815354529}. NIntegrate obtained 2.330428756121359`*^-6 and 7.096561136045723`*^-11 for the integral and error estimates. >>

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in x near {x} = {-0.99999999999999667104594803109799410885580269090806374709815354529}. NIntegrate obtained 4.308639728030708`*^-8 and 7.540153824713211`*^-11 for the integral and error estimates. >>

General::stop: Further output of NIntegrate::ncvb will be suppressed during this calculation. >>

 {-0.00006666779216173972, 1.483597024240845*10^-6, 
 2.742965242872827*10^-8, -4.361781616250265*10^-10,
-9.55048318138621*10^-12, 1.324805938959673*10^-13, 
 3.053277276241757*10^-15, 
 6.813896661329826*10^-16, -2.257310336427579*10^-15,
-1.279952387436186*10^-15, 4.762549327712297*10^-15, 
 8.651550477286525*10^-15, 6.111965912879058*10^-14, 
 5.913733424919787*10^-14, -2.416647291094786*10^-13,
-5.98964708075444*10^-13, -5.167457130207156*10^-13, 
 8.120386801889474*10^-13, 3.653845513969105*10^-12}
$\endgroup$
1
  • 1
    $\begingroup$ By the way, the function T is already in Mathematica as ChebyshevT. $\endgroup$
    – Michael E2
    Commented Sep 19, 2017 at 3:03

2 Answers 2

1
$\begingroup$

Use exact numbers, and adjust WorkingPrecision as well:

m = 20; 
\[Mu] = 1/2;
u[x_] := (12*\[Mu]^2)/25 + (6/25)*\[Mu]^2* Sech[(x*\[Mu])/5]^2 - (12/25)*\[Mu]^2*Tanh[(x*\[Mu])/5]; 

T[n_, x_] := Sum[(-1)^i*2^(n - 2*i - 1)*((n*(n - i - 1)!)/(i!*(n - 2*i)!))*x^(n - 2*i), 
   {i, 0, Ceiling[n/2]}];

Table[b[n] = (2/Pi)*NIntegrate[(u[1/2*x + 6/5]*T[n, x])/Sqrt[1 - x^2], 
      {x, -1, 1}, WorkingPrecision->60, PrecisionGoal -> 12], {n, 2, m}]

{-0.0000530269922322758700872406031298674120821932821983747589387163, 1.75261542689904550105630367593488021854698740535599027651498*10^-6, 2.03174943882974894986350353242628992979051384185133054522738*10^-8, \ -5.71381464509516823277911551847652467633881431176887724116102*10^-10, \ -6.37328149248641908321589282844941763150469420469053816298615*10^-12, 1.87011744413795397954172376476016957818471258995457314173343*10^-13, 1.75068011813384050090653986140979403776340305968235591804942*10^-15, \ -5.97686457967731296979103077829952153428513260249359486336513*10^-17, \ -4.26446546741900038202004327572870149796469823164759633901685*10^-19, 1.85208901415995012427155570308574586021728606537005677375062*10^-20, 8.94550465566852889943538290071348096418450173316421494793794*10^-23, \ -5.55865347380339884310848480886507901940568338741575442041513*10^-24, \ -1.42039361968678957657596810831771654198253211270604202264774*10^-26, 1.61716665239712597324147883468714582375220102302954973310439*10^-27, 5.12461112860934563767464690082079472604671608476673843091573*10^-31, \ -4.56568185693845647383269288596444770199336794356696282143777*10^-31, 8.78178394761340122535830010605721994209448116711339777120330*10^-34, 1.25183364269025093541820893853862105582025612063513904034290*10^-34, \ -5.29261315659887056938709191085499138650450119587992124610381*10^-37}

$\endgroup$
1
  • $\begingroup$ Thank you very much sir Carl Woll. $\endgroup$
    – Khaled
    Commented Sep 18, 2017 at 17:37
2
$\begingroup$

The Chebyshev interpolation coefficients converge to the Chebyshev series coefficients b[i] quite rapidly, and the FFT is an efficient way to compute the interpolation coefficients. See this answer by J.M. for more, or search the site for "Chebyshev Boyd" or "Chebfun", which is a MATLAB package for using Chebyshev approximations to solve a wide variety of numerical problems.

m = 20; μ = 1/2;
u[x_] := (12*μ^2)/25 + (6/25)*μ^2*Sech[(x*μ)/5]^2 - (12/25)*μ^2*Tanh[(x*μ)/5];

Clear[chebnodes, chebcoeffs];
chebnodes[n_] := chebnodes[n] = Sin[Pi/2*Range[-n, n, 2]/n];
chebcoeffs[f_, n_, wp_] := Module[{xx, cc},
   xx = N[chebnodes[n], wp];
   cc = Sqrt[2/n] FourierDCT[f /@ xx, 1];
   cc[[{1, -1}]] /= 2;
   cc
   ];

prec = 52;  (* this working precision yields a precision of 16+ in the b[i] *)
n = 2^Ceiling[Log2[2. m]];  (* the interpolation of order n >= 2m is sufficient *)
cc = chebcoeffs[u[1/2*# + 6/5] &, n, prec];
cc[[3 ;; m + 1]]   (* the array cc equals {b[0], b[1], ..., b[n]} *)
(*
{-0.00005302699223227587008724060312986741208219328219847, 
 -1.75261542689904550105630367593488021854698740536*10^-6, 
  2.031749438829748949863503532426289929790513832*10^-8, 
  5.7138146450951682327791155184765246763388143*10^-10, 
 -6.37328149248641908321589282844941763150480*10^-12, 
 -1.8701174441379539795417237647601695781848*10^-13, 
  1.75068011813384050090653986140979403765*10^-15, 
  5.976864579677312969791030778299521534*10^-17, 
 -4.2644654674190003820200432757287029*10^-19, 
 -1.852089014159950124271555703085746*10^-20, 
  8.945504655668528899435382900713*10^-23, 
  5.55865347380339884310848480887*10^-24, 
 -1.420393619686789576575968108*10^-26, 
 -1.61716665239712597324147883*10^-27,
  5.1246111286093456376746*10^-31,
  4.5656818569384564738327*10^-31, 
  8.7817839476134012254*10^-34,
 -1.2518336426902509354*10^-34, 
 -5.2926131565988706*10^-37}
*)

The loss of precision is really because the accuracy is roughly constant and the magnitude of the coefficients decreases:

Accuracy /@ cc // MinMax
(*  {53.2505, 53.8266}  *)
$\endgroup$
1
  • $\begingroup$ Wow, we sure are getting a lot of mileage out of Chebyshev... :D $\endgroup$ Commented Sep 19, 2017 at 5:40

Not the answer you're looking for? Browse other questions tagged or ask your own question.