function set3() %Math505-Set3 Convergence for trap rule % smooth, nonperiodic functions (cases 2, 6) get nominal % improvement from Rhomberg integration % case 1 is not diff/ble at x=0 % case 4 is periodic (alredy trap rule very well behaved) % cases 3, 5 have 0 coeffs in leading part of error due to % spurious cancellation of odd derivatives at x=0. n = 14; N = 2^(n-1); h = 1/N; f(1,1) = 0;%1;%0; f(1,2) = 1; f(1,3) = 1; f(1,4) = 1; f(1,5) = 1; f(1,6) = 0; for i = 1:N ih = i*h; f(i+1,1) = sqrt(ih);%1-ih^2;%sqrt(ih); f(i+1,2) = 1/(1+ih^2); f(i+1,3) = cos(8*sin(pi*ih)-pi*ih); f(i+1,4) = 1/(1+.5*sin(2*pi*ih)); f(i+1,5) = exp(ih^2*(1-ih)^2); f(i+1,6) = ih^2*exp(-ih); end format long g(1) = 2/3; g(2) = pi/4; g(3) = 0.234636347; g(4) = 1.15470054; g(5) = 1.03414105; g(6) = -5*exp(-1) + 2; for j = 1:6 Ni = N/2; L = 2; T(1) = .5*(f(1,j) + f(N+1,j)); T(2) = T(1) + f(Ni + 1,j); for i = 3:n Ni = Ni/2; L = 2*L; T(i) = T(i-1); for k = 1:2:L-1 T(i) = T(i) + f(Ni*k+1,j); end end a = 1; for i = 1:n T(i) = a*T(i); a = a/2; end if (j == 1) for i = n:-1:2 T(i) = 2*T(i) - T(i-1); end end format long R=rhomberg(n,T); if (j == 3 || j == 4 || j == 5) % use final Rhomberg estimate in lieu of full-precision answer for 3-5 g(j) = R(n); end fprintf('\n n Trapezoidal Rhomberg Exact \n %1.1i',j) for i = 1:n fprintf('\n %1.2i %1.14f %1.14f %1.14f',i,T(i),R(i), g(j)) end fprintf('\n n T-error R-error \n %1.1i',j) for i = 1:n fprintf('\n %1.2i %1.14f %1.14f ',i,abs(T(i)-g(j)),abs(R(i)-g(j))) end end end function [ R ] = rhomberg( n, T ) %UNTITLED performs Richardson extrapolation on a vecrtor % of results for trapezoidal integration with % step halved for each successive entry % R(n+1;k) = [2^(2n) R(n;k) - R(n;k-1)]/[2^(2n)-1] R(1:n) = T(1:n); a = 1; for k = 2:n a = a*4; b = a -1; for i = n:-1:k R(i) = (a*R(i) - R(i-1))/b; end end end OUTPUT; INCLUDES ADDITIONAL CASE (#6) AS A SECOND EXAMPLE OF A SMOOTH FUNCTION FOR WHICH RHOMBERG INTEGRATION PERFORMS NOMINALLY. n Trapezoidal Rhomberg Exact 1 01 0.50000000000000 0.50000000000000 0.66666666666667 02 0.70710678118655 0.77614237491540 0.66666666666667 03 0.68301270189222 0.66823727327491 0.66666666666667 04 0.67297739700616 0.66929217510186 0.66666666666667 05 0.66903217213000 0.66755586857413 0.66666666666667 06 0.66753667568066 0.66698061406877 0.66666666666667 07 0.66698268647807 0.66677761415593 0.66666666666667 08 0.66678050321514 0.66670588814990 0.66666666666667 09 0.66670744065623 0.66668053316816 0.66666666666667 10 0.66668121412338 0.66667156918109 0.66666666666667 11 0.66667184288802 0.66666839996424 0.66666666666667 12 0.66666850496696 0.66666727947963 0.66666666666667 13 0.66666731866139 0.66666688332875 0.66666666666667 14 0.66666689769597 0.66666674326828 0.66666666666667 n T-error R-error 1 01 0.16666666666667 0.16666666666667 02 0.04044011451988 0.10947570824873 03 0.01634603522555 0.00157060660825 04 0.00631073033950 0.00262550843519 05 0.00236550546334 0.00088920190746 06 0.00087000901399 0.00031394740210 07 0.00031601981141 0.00011094748926 08 0.00011383654848 0.00003922148324 09 0.00004077398956 0.00001386650150 10 0.00001454745671 0.00000490251443 11 0.00000517622135 0.00000173329758 12 0.00000183830029 0.00000061281297 13 0.00000065199472 0.00000021666208 14 0.00000023102930 0.00000007660161 n Trapezoidal Rhomberg Exact 2 01 0.75000000000000 0.75000000000000 0.78539816339745 02 0.77500000000000 0.78333333333333 0.78539816339745 03 0.78279411764706 0.78552941176471 0.78539816339745 04 0.78474712362277 0.78539644594047 0.78539816339745 05 0.78523540301035 0.78539816631943 0.78539816339745 06 0.78535747329374 0.78539816340956 0.78539816339745 07 0.78538799087141 0.78539816339743 0.78539816339745 08 0.78539562026594 0.78539816339745 0.78539816339745 09 0.78539752761457 0.78539816339745 0.78539816339745 10 0.78539800445173 0.78539816339745 0.78539816339745 11 0.78539812366102 0.78539816339745 0.78539816339745 12 0.78539815346334 0.78539816339745 0.78539816339745 13 0.78539816091392 0.78539816339745 0.78539816339745 14 0.78539816277657 0.78539816339745 0.78539816339745 n T-error R-error 2 01 0.03539816339745 0.03539816339745 02 0.01039816339745 0.00206483006411 03 0.00260404575039 0.00013124836726 04 0.00065103977468 0.00000171745698 05 0.00016276038710 0.00000000292198 06 0.00004069010370 0.00000000001211 07 0.00001017252603 0.00000000000002 08 0.00000254313151 0.00000000000000 09 0.00000063578288 0.00000000000000 10 0.00000015894572 0.00000000000000 11 0.00000003973643 0.00000000000000 12 0.00000000993411 0.00000000000000 13 0.00000000248353 0.00000000000000 14 0.00000000062088 0.00000000000000 n Trapezoidal Rhomberg Exact 3 01 0.00000000000000 0.00000000000000 0.23463634685391 02 0.49467912331169 0.65957216441559 0.23463634685391 03 0.04009498113445 -0.16283341496954 0.23463634685391 04 0.23436316573765 0.33425640866840 0.23463634685391 05 0.23463634685391 0.22849681707301 0.23463634685391 06 0.23463634685391 0.23472565059340 0.23463634685391 07 0.23463634685391 0.23463610288150 0.23463634685391 08 0.23463634685391 0.23463634668218 0.23463634685391 09 0.23463634685391 0.23463634685436 0.23463634685391 10 0.23463634685391 0.23463634685391 0.23463634685391 11 0.23463634685391 0.23463634685391 0.23463634685391 12 0.23463634685391 0.23463634685391 0.23463634685391 13 0.23463634685391 0.23463634685391 0.23463634685391 14 0.23463634685391 0.23463634685391 0.23463634685391 n T-error R-error 3 01 0.23463634685391 0.23463634685391 02 0.26004277645778 0.42493581756167 03 0.19454136571947 0.39746976182345 04 0.00027318111626 0.09962006181448 05 0.00000000000000 0.00613952978091 06 0.00000000000000 0.00008930373949 07 0.00000000000000 0.00000024397241 08 0.00000000000000 0.00000000017173 09 0.00000000000000 0.00000000000044 10 0.00000000000000 0.00000000000000 11 0.00000000000000 0.00000000000000 12 0.00000000000000 0.00000000000000 13 0.00000000000000 0.00000000000000 14 0.00000000000000 0.00000000000000 n Trapezoidal Rhomberg Exact 4 01 1.00000000000000 1.00000000000000 1.15470053837925 02 1.00000000000000 1.00000000000000 1.15470053837925 03 1.16666666666667 1.23703703703704 1.15470053837925 04 1.15476190476190 1.14458721760309 1.15470053837925 05 1.15470054000982 1.15512171896038 1.15470053837925 06 1.15470053837925 1.15469621974690 1.15470053837925 07 1.15470053837925 1.15470053097776 1.15470053837925 08 1.15470053837925 1.15470053847798 1.15470053837925 09 1.15470053837925 1.15470053837914 1.15470053837925 10 1.15470053837925 1.15470053837925 1.15470053837925 11 1.15470053837925 1.15470053837925 1.15470053837925 12 1.15470053837925 1.15470053837925 1.15470053837925 13 1.15470053837925 1.15470053837925 1.15470053837925 14 1.15470053837925 1.15470053837925 1.15470053837925 n T-error R-error 4 01 0.15470053837925 0.15470053837925 02 0.15470053837925 0.15470053837925 03 0.01196612828742 0.08233649865779 04 0.00006136638266 0.01011332077616 05 0.00000000163057 0.00042118058113 06 0.00000000000000 0.00000431863235 07 0.00000000000000 0.00000000740149 08 0.00000000000000 0.00000000009873 09 0.00000000000000 0.00000000000011 10 0.00000000000000 0.00000000000000 11 0.00000000000000 0.00000000000000 12 0.00000000000000 0.00000000000000 13 0.00000000000000 0.00000000000000 14 0.00000000000000 0.00000000000000 n Trapezoidal Rhomberg Exact 5 01 1.00000000000000 1.00000000000000 1.03414105175077 02 1.03224722945893 1.04299630594524 1.03414105175077 03 1.03401438324028 1.03404390973776 1.03414105175077 04 1.03413297284034 1.03414535906654 1.03414105175077 05 1.03414054406481 1.03414104627654 1.03414105175077 06 1.03414101997639 1.03414105179021 1.03414105175077 07 1.03414104976419 1.03414105175078 1.03414105175077 08 1.03414105162660 1.03414105175078 1.03414105175077 09 1.03414105174302 1.03414105175078 1.03414105175077 10 1.03414105175029 1.03414105175078 1.03414105175077 11 1.03414105175075 1.03414105175078 1.03414105175077 12 1.03414105175078 1.03414105175078 1.03414105175077 13 1.03414105175078 1.03414105175077 1.03414105175077 14 1.03414105175077 1.03414105175077 1.03414105175077 n T-error R-error 5 01 0.03414105175077 0.03414105175077 02 0.00189382229184 0.00885525419447 03 0.00012666851050 0.00009714201301 04 0.00000807891043 0.00000430731577 05 0.00000050768596 0.00000000547423 06 0.00000003177438 0.00000000003944 07 0.00000000198659 0.00000000000000 08 0.00000000012417 0.00000000000000 09 0.00000000000776 0.00000000000000 10 0.00000000000048 0.00000000000001 11 0.00000000000002 0.00000000000001 12 0.00000000000001 0.00000000000001 13 0.00000000000000 0.00000000000000 14 0.00000000000000 0.00000000000000 n Trapezoidal Rhomberg Exact 6 01 0.18393972058572 0.18393972058572 0.16060279414279 02 0.16778619275694 0.16240168348068 0.16060279414279 03 0.16248840509317 0.16061052869799 0.16060279414279 04 0.16107989607964 0.16060280012980 0.16060279414279 05 0.16072242723629 0.16060279414377 0.16060279414279 06 0.16063272478883 0.16060279414279 0.16060279414279 07 0.16061027820294 0.16060279414279 0.16060279414279 08 0.16060466524525 0.16060279414279 0.16060279414279 09 0.16060326192387 0.16060279414279 0.16060279414279 10 0.16060291108840 0.16060279414279 0.16060279414279 11 0.16060282337921 0.16060279414279 0.16060279414279 12 0.16060280145190 0.16060279414279 0.16060279414279 13 0.16060279597007 0.16060279414279 0.16060279414279 14 0.16060279459961 0.16060279414279 0.16060279414279 n T-error R-error 6 01 0.02333692644293 0.02333692644293 02 0.00718339861415 0.00179888933789 03 0.00188561095038 0.00000773455520 04 0.00047710193685 0.00000000598701 05 0.00011963309350 0.00000000000098 06 0.00002993064604 0.00000000000000 07 0.00000748406015 0.00000000000000 08 0.00000187110246 0.00000000000000 09 0.00000046778108 0.00000000000000 10 0.00000011694561 0.00000000000000 11 0.00000002923642 0.00000000000000 12 0.00000000730911 0.00000000000000 13 0.00000000182728 0.00000000000000 14 0.00000000045682 0.00000000000000