lambdaway :: simpson

6 | list | login |load

| |

_h1 [[integral]] | adaptative simpson's method

_p We compute {b ∫{sub 0}{sup 1}sin(x)} using the Adaptive Simpson's method, as it can be seen in  [[rosettacode.org/Adaptive Simpson's method|http://rosettacode.org/wiki/Numerical_integration/Adaptive_Simpson%27s_method]].

_h3 using the {b sin}'s primitive {b -cos(x)}

{pre
∫{sub 0}{sup 1} sin(x) = '{- {cos 0} {cos 1}} 
-> {- {cos 0} {cos 1}} 
}

_h3 using the adaptative simpson method

{pre 
'{def adasimp

{def adasimp.calc
  {lambda {:fun :a :fa :b :fb}
   {let { {:delta {- :b :a}}
          {:m {/ {+ :a :b} 2}}
          {:fa :fa} {:fb :fb}
          {:fm {:fun {/ {+ :a :b} 2}}}
        } {A.new :m :fm {/ {* {abs :delta} {+ :fa {* 4 :fm} :fb}} 6}} }}}

{def adasimp.rec
  {lambda {:fun :a :fa :b :fb :eps :mfs}
   {let { {:fun :fun}
          {:a :a} {:fa :fa} {:b :b} {:fb :fb}
          {:eps :eps} {:mfs :mfs}
          {:left  {adasimp.calc :fun :a :fa {A.get 0 :mfs} {A.get 1 :mfs}}}
          {:right {adasimp.calc :fun {A.get 0 :mfs} {A.get 1 :mfs} :b :fb}}
        } {let { {:fun :fun} {:a :a} {:fa :fa} {:b :b} {:fb :fb}
                 {:eps :eps} {:mfs :mfs} {:left :left} {:right :right}
                 {:delta {+ {A.get 2 :left}
                            {A.get 2 :right}
                         {- {A.get 2 :mfs}}}}
               } {if {<= {abs :delta} {* :eps 15}}
                  then {+ {A.get 2 :left}
                          {A.get 2 :right}
                          {/ :delta 15}}
                  else {+ {adasimp.rec :fun :a :fa 
                                       {A.get 0 :mfs} {A.get 1 :mfs}
                                       {/ :eps 2} :left }
                          {adasimp.rec :fun 
                                       {A.get 0 :mfs} {A.get 1 :mfs}
                                       :b :fb {/ :eps 2} :right }}} }}}}
 
 {lambda {:fun :a :b :eps}
  {adasimp.rec :fun :a {:fun :a} :b {:fun :b} :eps 
               {adasimp.calc :fun :a {:fun :a} :b {:fun :b}}} }}
-> {def adasimp

{def adasimp.rec
  {lambda {:fun :a :fa :b :fb :eps :mfs}
   {let { {:fun :fun} {:a :a} {:fa :fa} {:b :b} {:fb :fb}
          {:eps :eps} {:mfs :mfs}
          {:left  {adasimp.calc :fun :a :fa {A.get 0 :mfs} {A.get 1 :mfs}}}
          {:right {adasimp.calc :fun {A.get 0 :mfs} {A.get 1 :mfs} :b :fb}}
        } {let { {:fun :fun} {:a :a} {:fa :fa} {:b :b} {:fb :fb}
                 {:eps :eps} {:mfs :mfs} {:left :left} {:right :right}
                 {:delta {+ {A.get 2 :left} {A.get 2 :right} {- {A.get 2 :mfs}}}}
               } {if {<= {abs :delta} {* :eps 15}}
                  then {+ {A.get 2 :left} {A.get 2 :right} {/ :delta 15}}
                  else {+ {adasimp.rec :fun :a :fa {A.get 0 :mfs} {A.get 1 :mfs} {/ :eps 2} :left }
                          {adasimp.rec :fun {A.get 0 :mfs} {A.get 1 :mfs} :b :fb {/ :eps 2} :right }}} }}}}
 
 {lambda {:fun :a :b :eps}
  {adasimp.rec :fun :a {:fun :a} :b {:fun :b} :eps 
               {adasimp.calc :fun :a {:fun :a} :b {:fun :b}}} }}

'{adasimp sin 0 1 1.0e-11}
-> {adasimp sin 0 1 1.0e-11}  
  (0.45969769413186023)  // reference value
}

_p {i alain marty | 2022/06/05}