7+͈C:)=_gFLUF10 CRCFLUF10 DOC!FLUF10 BAS%FLUF10/TPAS?3FLUF10/TCOMrP!xv owOo!v4o<552vx!v2vviKgKqo.og!vCRC Ver 5.0 CTL-S pauses, CTL-C aborts --> FILE: FLUF10 .DOC CRC = CA 49 --> FILE: FLUF10 .BAS CRC = E1 F5 --> FILE: FLUF10/T.PAS CRC = FC 97 --> FILE: FLUF10/T.COM CRC = EC 7E DONE  FLUFF Version 1.0 MBASIC (tm) and TURBO Pascal (tm) July 7, 1984 D.M. Fritz-Rohner Post Office Box 9080 Akron, Ohio 44305 Thes program deriv fro th BASI languag versio pre sente i th articl ' Simpl Minima Algorithm b Steve A Ruzinsk foun i Dr Dobb' Journal Numbe 93 July 1984. Th fil FLUFF10.BA i carefull copie fro th articl an i MBASI compatible. Th fil FLUFF10.PAӠ i transliteratio o th BASI versio t TURB Pascal Th syntactica an structura change ar obvious Th variabl i substitute fo intermediat variabl usag o q an xk Exponentiatio operato usag i circumvented Th covarianc matri invers an coefficien vecto ar completel initialized. Th articl propose a algorith tha th autho call 'Erro Fluffing t perfor minima fit Th articl present bot BASI an FORT derive languag application o th algo rith t polynomia fi o elementar transcendenta functions. Th堠 progra initialize tabl wit th堠 normalize cofactor o th term o th polynomial Normalizatio use th fitted functio valu t produc relativ erro minima fit. Th diagona element o th covarianc matri invers ar initialize t somethin larg t reflec uncertaint o initia coefficien vecto value whic ar arbitraril initialize t zero. Th progra applie th matri inversio lemm t updat an estimat th covarianc matri invers an coefficien vecto fo eac o th compute discret SI functio values Thi i th sequentia leas square estimat o th polynomia coefficients Thi portio o th algorith migh equall wel b don usin batc leas square wit direc computatio o th covarianc matrix it inverse an a initia coefficien vecto fo th 'erro fluffing process. Th novelt o th algorith consist i wha i no does Th progra searche th se o fitte functio value t fin th on a whic th relativ erro o th leas square fi i poorest Thi valu i the applie t th sequentia leas square estimato t forc adjustmen o th coefficient i a attemp t improv th poores fit Presto Minimax. Th BASIà progra illustrate th algorith b estimatin th coefficient o a od polynomia o fift orde t fi th BASIà SIΠ functio representatio o th mathematica sin functio o th interva zer t 9 degrees Th SI functio i sample a fift uniforml space point o th interval Th progra specifie on thousan iterations. Fo Mh Z-8 base system th MBASI versio require roughl si second pe iteration Afte fift iteration th progra report th followin values: Coefficients: 1.570627592139162 -.6432238295815594 .07270234361204301 Iterations = 51 Maximum Absolute Error = 1.088142336302916D-04 wit Figur i th referenc articl assume t b th fina resul fo thi case. Th TURBϠ Pasca versio take abou on hal secon pe iteratio an report th followin value afte on thousan iterations: Coefficients: 1.570626E+00 -6.43222E-01 7.270485E-02 Iteration = 1000 Maximum Absolute Error = 1.080364E-04 Not tha th lis produce b th progra say 'Maximu Absolut Erro ' Th valu is however th absolut valu o th relativ error Thi progra illustrate th algorith fo particula case Th normalizatio o BASI progra lin 21 wil necessaril caus problem fo a independen variabl value of zero. MBASIC (tm) MicroSoft TURBO Pascal (tm) Borland International ad: MINIMAX VIA ERROR FLUFFING ALGORITHM by Steven A. Ruzinsky This program demonstrates the determination of polynomial coefficients to a mimimum maximum absolute error criterion. The resultbn: is better than that of Hastings, Approximations for Digital Computers, 1955, Princeton University Press, p 138. This program is in IBM Basic.bx: -------------------------------------------------------------------------b AZ : INcN  : M 2 : ITERATIONS [c X(M,N) , Y(M) , XK(N) , Q(N,N) , QX(N) , AK(N) , A(N)c: -------------------------------------------------------------------------xd: The following fills matrices Y(i) and X(i,j) with data. The function used is SIN(PI*X/2), however, the data is modified so that the minimax criterion is applied to the relative error :d I  MdY(I) dE IMdD (*IE)d J  NdX(I,J) DE(JJ)d J , ICe: -------------------------------------------------------------------------e: The following initiates the Q matrix. It may be necessary to to adjust the number in line 280 for best results.e I  NeQ(I,I) $te" I>f,: -------------------------------------------------------------------------.g6: The following loop with index, K, reiterates the sequential least squares algorithm. Up to limit M, each data point is incorporated once into Q and AK. This results in a least squares fit to the data. Afterwards, the datag@: corresponding to the maximum absolute error are reincorporated back into Q and AK.gJEBEST gT K  ITERATIONS Mg^ K M  : ghD gr J  Ng|QX  h I  N'hQX QX XK(I)Q(J,I)/h I>hQX(J) QXShD D XK(J)QX[h Jkh J  N|hQX QX(J)Dh I  NhQ(I,J) Q(I,J) QX(I)QXh I , Jh J  NhQX h I  NhQX QX XK(I)Q(J,I)i Ii&AK(J) AK(J) QXE*i0 J , Kzi:: -----------------------------------------------------------------------iD: The following prints the results :iN "Coefficients:" : I  NiXAK(I) A(I)ib AK(I)il I :  : Ijv: -----------------------------------------------------------------------j: Subroutine for incorporating each data point once :jE Y(K)j I  NjXK X(K,I)jXK(I) XKjE E AK(I)XKj Ij4k: -----------------------------------------------------------------------k: Subroutine for finding maximum absolute error and correspnding data point :kEMAX  : JMAX k J  MkE Y(J)k I  Nk E E X(J,I)AK(I)k Il E (E)&l* E EMAX EMAX E : JMAX J.l4 Jjl> "Iterations =", KM, "Maximum Absolute Error =", EMAXzlHE Y(JMAX)lR I  Nl\XK X(JMAX,I)lfXK(I) XKlpE E AK(I)XKlz Il EMAX EBEST EBEST EMAX : lHm: ----------------------------------------------------------------------}m: Subroutine for saving best coefficients, A :m I NmA(I) AK(I)m Imm, PRINCETON UNIVERSITY PRESS, P 138. THIS PROGRAM IS IN IBM BASIC.s, A :m I NmA(I) AK(I)m Im { Copyright (C) 1984 D.M. Fritz-Rohner } program main ; { Name: } { FLUFF - Error FLUFFing Algorithm. } { } { Version: Date: } { 1.0/TURBO (tm) 1984 July 7 } { } { Purpose: } { Demonstrate PASCAL application of error fluff- } { ing algorithm proposed by the reference. } { } { Author: } { D.M. Fritz-Rohner } { Post Office Box 9080 } { Akron, Ohio 44305 } { } { Reference: } { 'A Simple Minimax Algorithm' } { Steven A. Ruzinsky } { Dr. Dobbs Journal } { Number 93, July, 1984 } { pg. 84 et seq } { } { Description: } { This program estimates polynomial coefficients } { for a minimax fit to an elementary transcenden- } { tal function. See FLUFF10.DOC. } { } { This version is written in TURBO Pascal and is } { transliterated from the Basic version pre- } { sented in the reference. See FLUFF10.BAS. } { } { The variable t is substituted for intermediate } { variable usage of the symbols qx and xk. Ex- } { ponentiating operator usage is circumvented. } { Covariance matrix inverse and coefficient } { vector are fully initialized. } { } { TURBO Pascal (tm) Borland International } const N = 3 ; M = 50 ; Iterations = 1000 ; var x : array[1..M,1..N] of real ; y : array[1..M] of real ; xk : array[1..N] of real ; q : array[1..N,1..N] of real ; qx : array[1..N] of real ; ak : array[1..N] of real ; a : array[1..N] of real ; d : real ; e, e2 : real ; eMax : real ; eBest : real ; t : real ; i, j, k : integer ; jMax : integer ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure IncludeData ; { Incorporate each data point once. } begin e := y[k] ; for i := 1 to N do begin t := x[k,i] ; xk[i] := t ; e := e - ak[i]*t end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure SaveBest ; { Save coefficients of best fit. } begin for i := 1 to N do a[i] := ak[i] end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} procedure FindMax ; { Find maximum absolute error and corresponding index. } begin eMax := 0.0 ; jMax := 1 ; for j := 1 to M do begin e := y[j] ; for i := 1 to N do e := e - x[j,i]*ak[i] ; e := ABS(e) ; if e > eMax then begin eMax := e ; jMax := j end end ; writeln('Iterations = ',k-M:4,' Maximum Absolute Error = ',eMax:12) ; e := y[jMax] ; for i := 1 to N do begin t := x[jMax,i] ; xk[i] := t ; e := e - ak[i]*t end ; if eMax < eBest then begin eBest := eMax ; SaveBest end end ; {~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~} begin { Fill y and x with data. The function used } { is SIN(Pi*x/2), however, the data is modified so } { that the minimax criterion is applied to the rela- } { tive error. } for i := 1 to M do begin y[i] := 1.0 ; e := i/M ; e2 := e*e ; d := 1.0/SIN(1.570796327*e) ; for j := 1 to N do begin x[i,j] := d*e ; e := e*e2 end end ; { The following initiates the q matrix. It may be } { necessary to to adjust the number below for } { best results. } for i := 1 to N do begin q[i,i] := 1.0e+06 ; if i < N then begin for j := i+1 to N do begin q[i,j] := 0.0 ; q[j,i] := 0.0 end end ; xk[i] := 0.0 ; ak[i] := 0.0 end ; { The following loop with index, k, reiterates the } { sequential least squares algorithm. Up to limit M, } { each data point is incorporated once into q and ak. } { This results in a least squares fit to the data. } { Afterwards, the data corresponding to the maximum } { absolute error are incorporated back into q and ak. } eBest := 1.0 ; for k := 1 to Iterations+M do begin if k > M then FindMax else IncludeData ; d := 1.0 ; for j := 1 to N do begin t := 0.0 ; for i := 1 to N do t := t + xk[i]*q[j,i] ; qx[j] := t ; d := d + xk[j]*t ; end ; for j := 1 to N do begin t := qx[j]/d ; for i := 1 to N do q[i,j] := q[i,j] - qx[i]*t ; end ; for j := 1 to N do begin t := 0.0 ; for i := 1 to N do t := t + xk[i]*q[j,i] ; ak[j] := ak[j] + t*e end end ; { Print the results. } writeln('Coefficients:') ; for i := 1 to N do begin ak[i] := a[i] ; writeln(ak[i]:12) end ; FindMax end. ͫCopyright (C) 1984 BORLAND IncA SDS VDB-8024selectedP=  ~7#~=% o&ͦoͦܐԩͣ}!!"8~#(}:$= +*!Z!*B!!:(=2!Z: <2!!!:O::O:!*B! !45(!.+/ 0y0( d!kZ!{Z͈͈o&  :(y ͠|( *"x2y( >28!?"9!!>2 :D]SXN]D [ (!e}̈́A8Q0G: x@!\w# (   yV. V!h6# (*(.(!8}(*(̈́w#>?> w#a{ |͒}͛Ɛ'@'7||}>"C"6# ""͐ͩ*B"[R5*"^#V#^#V#N#FO/o&9O/o&9!9(> (G!9 w#E͊w}8uRB0 >R@RR!+ͨ z R!+ͨ z <!+ͨ z <!+ͨ z <!#ͨ z <!+ͨ z T]KB!z> S>))0 = |JJDMgo>jB0 7?= H\<z5+)+<z {0Gɯgo||H}||/g}/o#}o&K[xAJSJDM!b"!6J"DM'ͬͬdͬ ͬ} wͦWͧ _}8(8J`9{T]=o`9y w >uJ u` }>(; xQ }} ˸T}ٕ(0D=C ,= ( [ 0%D , 7 ͏ ?(8u x O - ; 8˸x X ,-xG}; }م 9; .>#n0[ D = - nx P ,-(-˸G,-; }ٕ? 9.>͏ 8u ?= u+-(>O 0u O 8͏ ?x P , 78ƀ8ƀ8ox٨!دoGOW_gɷɷ|لg{ً_zيWyىOxوG|ٔg{ٛ_zٚWyٙOx٘Gxٸyٹzٺ{ٻ|ټx٨ xx(ͼ ?}ٽÏ }ց; <(; 7D = |٤g{٣_z٢Wy١Ox٠GD u J }x>uu}ƀ/ƀo; -J }0W-J W,}l˸ͨ 8 ; ` x( -ͨ 8J -ͨ 8,J }l8;*!` ! >u` ` u--- J ,,,-xGg?+2n*8t z~,->uxua}.; OJ , ; !U >,k- o&0%,` }g; }؉}颋.:}8c~I$I~L*kٷx˸; }0G,͙<},-(-J ! >0 a` o8 Oþ >um.`1pF,t6|!wS<.z}[|%FXc~ur1}Oٯx(<˸ͨ 8; !~Jͨ 0O!><ͨ 8 =  7 <` O ; 7 0 W-J OT0 j oD,:j !I}袋.}8c~I$I~L!>u` ` 77 ` = O nf^VNF!DLT\I!!53!r1!\!> x #-= o˸xO(- }(x>8(C ,C `iM!>u|; |J>| )=|(DMbo˸ͦ88ͦx(0 8> Mx(>-Ͳ{(ay(Ͱͦ \z(>.Ͳ (Ͱ ~ͦ{>EͲ>+|(|Dg>-Ͳ|/ 0:p# ~# +>0w#,-  60#J˸}րogM| .(C = ~> x0w#xG%P %P ZJDM%P = _~65i+~hìx-Sx9?+{Η@}|C C gZJDM0D ,7}o˸  #yO!@9i&   # w# /w# w#!9! E9!!9~(+F͊!"9!(#>2*Ͳ"|>" :( ͆ *6#w*6#6 !\$![ (̈́( #:~CONTRMKBDLSTAUXUSR>2$*#~ Ͷ$*:> >w###6  #6++p>2S-$Ͷ:*6###ww#w$w#w: ##N#F*B> w#w#[s#r>2S$Ͷ$*6 #-Nw#Fwq#p#6#w#w#w* :( ͒: *^ F* < >26"~͟*-w#ww#͟"~ <@*Ͳ!\  <ʮ!\$> >2*|>! * \$\<(!: [1Á\!(f"> 2:!<"F( #~#6e>!["N>!~8>O6*"w (=(&("( :(N 8y(~#x+% (6*#~[*#~ *~(h#"b=  8 J= B== ͯ}8= ͵}/ͭ !*###~-_~(4Q6*>2>*##w:>*##~*#~(E[ ( ( ( !][ ( ( ((w#(6!]~-#8~>7  [>OkͼMs #rkͼpX á[ [ (( #w(q*#~[ (  *##~6͜O$*#~(08ʦ=ʦ==ʩ=ʬò+###~-_q46͡> *:4^q}Ò*|(M|( M6-#͐ͦ[R8 (G> ͒C~͒#*ͦC!h !lTRUEFALSEͦ!9^#(~#(G~͒#> ͒> Ò "F![(#RR0*4#4> RR *4 #4(>>2$*V(/˖:(#~+ x y2!͵( =( X:(R*:(###~-_-͌X> :("͟"*^˞*V˖0 SRѷR8A* N#F#s#r$ 0})jS\*###w* N#FB ͟r+s> !T]>)j)0 0= UR!#U*^#V#N#F#^#V>">!2DM"~x(L* :O(o:" C}=( ?*-N#Fp+qq#p! * F+N+++V+^Bq#p>>> SRѷR* s#r$ s#r"S"! N#FB(^x * 6#[<(H*! Kq#p##K[! *! 4 #4! x *$ *>w""{_!"*nf}(HR0nf" ^VMDnfutqp*s#r*s#r"* 5KB!>u~#fo{_"*R0RnfR0KqputsrNF( ^VNF^V*SutKqp R*R(~w~wnf ut"6# * *!""*NFy(* "*B0Cnf* [R*"*RS[s#r^#VS>O"w2x2!"" @*>2"!"""!\Ͳ*: !~6go(\R*s#r_2x( s x(T]DMR0 -a%}̈́o*!~6o&͠|ͣ}%^C User break1:% I/O% Run-time% error ͒%, PC=[R"͍% Program aborted*1!͍!(VͲ#! *d+)]T)!͡!!5zʘ "h!8*d+))))*h+)]T)!j͡!*h+)]T)!j͡!!*h+)]T)!j ͼ !͡*h# !!5z "h!*h+)]T)!*h+)]T)͡*h#ä !!v͡!"b!!25z!"f! *f+)]T)!͡!!5zʑ!"h!!8*f+))))*h+)]T)!*h+)]T) ͼ !͡*h#,!!˸!͡!!vͥE!!!v͡*f"b*f# Ŕ Iterations = *d!2R!́ Maximum Absolute Error = !v! !@͐b! *b+)]T)!͡!!5z""h!8*b+))))*h+)]T)!j͡!*h+)]T)!j͡!!*h+)]T)!j ͼ !͡*h#["!v!p͸E#!v!p͙͡ !!25z.$"h! *h+)]T)!͡*h!2 !͡!! !|͡!!I! P !͡!!5z%$"f!8*h+))))*f+)]T)!! ͡!!| !͡*f#ý#*h##!!5z>%"h!*h+))))*h+)]T)!$t͡*h!ͯE$*h!!5z$"f!*h+))))*f+)]T)!͡!*f+))))*h+)]T)!͡*f#Ï$!*h+)]T)!͡!*h+)]T)!͡*h#9$!!p͡!!!25zA("d*d!2͛E~% Á%!!͡!!5zʂ&"f!!j͡!!5z%&"h!j!*h+)]T)!*f+))))*h+)]T) ͳ !j͡*h#%!*f+)]T)!j͡!!*f+)]T)!j ͳ !͡*f#Ü%!!5zT'"f!*f+)]T)! !j͡!!5zK'"h!*h+))))*f+)]T)!*h+))))*f+)]T)!*h+)]T)!j ͼ ͡*h#&*f#Í&!!5z8("f!!j͡!!5z'"h!j!*h+)]T)!*f+))))*h+)]T) ͳ !j͡*h#Ã'!*f+)]T)!*f+)]T)!j! ͳ ͡*f#_'*d#_%Ŕ Coefficients:͐b!!5z("h!*h+)]T)!*h+)]T)͡!*h+)]T)R! !@͐b*h#f( icients:͐b!!5z"h!*h+)]T)!