� >Coef
@� Calculate polynomial approximations to arbitrary functions
C� in range [-1,+1], using Chebyshev Forced Oscillation of Error
(.� method to derive a near minimax solution
2�
<� �
FZ� "Please enter function to be approximated in range [-1,+1],  eg SIN(PI/4*(x+1)) ";f$
P<� "Please enter degree of approximation (1 to 100) ";n%'
Z)� f(101), c(101), t(100, 100), p(100)
d� i%=0 � n%+1
n x = �(i%*�/(n%+1))
x f(i%) = � f$
��
�� j%=0 � n%+1
�= � (j%�1) sum = (f(0)-f(n%+1))/2 � sum = (f(0)+f(n%+1))/2
� � i%=1 � n%
�$  sum += f(i%)*�(i%*j%*�/(n%+1))
� �
� c(j%) = sum/(n%+1)
� � j%<>0 c(j%) = 2*c(j%)
��
�E = �(c(n%+1)/2)
�&� '"Polynomial approximation is:"'
�� i%=0 � n%
�# � "T[";i%;"](x) *",c(i%);�27);

 � i%<n% � "+" � �

�

� '"where T[0](x)   = 1;"

"�  "      T[1](x)   = x;"

,9�  "and   T[n+1](x) = 2xT[n](x)-T[n-1](x), for n>=1"'

6$� Now we need to gen the T polys

@=� note we use t(i,j) to hold the coef of x^j in poly T[i]

Jt(0,0)=1

Tt(1,1)=1

^i%=1

hȕ i%<n%

r t(i%+1,0) = -t(i%-1,0)

| � j% = 1 � i%+1

�*  t(i%+1,j%) = 2*t(i%,j%-1)-t(i%-1,j%)

� �

�
 i%+=1

��

�� i%=0 � n%

� sum = 0

� � j%=i% � n%

�  sum += t(j%, i%)*c(j%)

� �

� p(i%) = sum

��

�$� '"Which may also be written:"'

�� i%=0 � n%
  � "x^";i%;"  *",p(i%);�27);
 � i%<n% � "+" � �
�
&maxerror = 0
0� x=-1 � 1 � 0.001
: v = p(n%)
D � i%=n%-1 � 0 � -1
N  v = v*x+p(i%)
X �
b e = �(v - � f$)
l � e>maxerror maxerror=e
v�
�E� ''"For this approximation an upper bound on error is ";maxerror
�Y�   "Also we now know minimax error p (the error made by the best possible order ";n%
�$�   "polynomial) lies between:"'
�� E;"  and  ";maxerror;"."'
�<� "This indicates how close our approximation is to it."
�V� ''"Finally, if only need say 16 bit accuracy, note that our poly's coefficients"
�>�   "multiplied by 2^20 & written in hex are as follows:"'
�� i%=n% � 0 � -1
� v = p(i%)
� � v<0 s$="-":v=-v � s$="+"
� v2 = v*1048576+0.5
�1 � "x^";i%;"  *",s$;�"00000000"+�~v2,8);�27);
� � i%>0 � "+" � �
�
�

00000000  0d 00 0a 0b f4 20 3e 43  6f 65 66 0d 00 14 40 f4  |..... >Coef...@.|
00000010  20 43 61 6c 63 75 6c 61  74 65 20 70 6f 6c 79 6e  | Calculate polyn|
00000020  6f 6d 69 61 6c 20 61 70  70 72 6f 78 69 6d 61 74  |omial approximat|
00000030  69 6f 6e 73 20 74 6f 20  61 72 62 69 74 72 61 72  |ions to arbitrar|
00000040  79 20 66 75 6e 63 74 69  6f 6e 73 0d 00 1e 43 f4  |y functions...C.|
00000050  20 69 6e 20 72 61 6e 67  65 20 5b 2d 31 2c 2b 31  | in range [-1,+1|
00000060  5d 2c 20 75 73 69 6e 67  20 43 68 65 62 79 73 68  |], using Chebysh|
00000070  65 76 20 46 6f 72 63 65  64 20 4f 73 63 69 6c 6c  |ev Forced Oscill|
00000080  61 74 69 6f 6e 20 6f 66  20 45 72 72 6f 72 0d 00  |ation of Error..|
00000090  28 2e f4 20 6d 65 74 68  6f 64 20 74 6f 20 64 65  |(.. method to de|
000000a0  72 69 76 65 20 61 20 6e  65 61 72 20 6d 69 6e 69  |rive a near mini|
000000b0  6d 61 78 20 73 6f 6c 75  74 69 6f 6e 0d 00 32 05  |max solution..2.|
000000c0  f4 0d 00 3c 07 eb 20 eb  0d 00 46 5a e8 20 22 50  |...<.. ...FZ. "P|
000000d0  6c 65 61 73 65 20 65 6e  74 65 72 20 66 75 6e 63  |lease enter func|
000000e0  74 69 6f 6e 20 74 6f 20  62 65 20 61 70 70 72 6f  |tion to be appro|
000000f0  78 69 6d 61 74 65 64 20  69 6e 20 72 61 6e 67 65  |ximated in range|
00000100  20 5b 2d 31 2c 2b 31 5d  2c 20 20 65 67 20 53 49  | [-1,+1],  eg SI|
00000110  4e 28 50 49 2f 34 2a 28  78 2b 31 29 29 20 22 3b  |N(PI/4*(x+1)) ";|
00000120  66 24 0d 00 50 3c e8 20  22 50 6c 65 61 73 65 20  |f$..P<. "Please |
00000130  65 6e 74 65 72 20 64 65  67 72 65 65 20 6f 66 20  |enter degree of |
00000140  61 70 70 72 6f 78 69 6d  61 74 69 6f 6e 20 28 31  |approximation (1|
00000150  20 74 6f 20 31 30 30 29  20 22 3b 6e 25 27 0d 00  | to 100) ";n%'..|
00000160  5a 29 de 20 66 28 31 30  31 29 2c 20 63 28 31 30  |Z). f(101), c(10|
00000170  31 29 2c 20 74 28 31 30  30 2c 20 31 30 30 29 2c  |1), t(100, 100),|
00000180  20 70 28 31 30 30 29 0d  00 64 11 e3 20 69 25 3d  | p(100)..d.. i%=|
00000190  30 20 b8 20 6e 25 2b 31  0d 00 6e 17 20 78 20 3d  |0 . n%+1..n. x =|
000001a0  20 9b 28 69 25 2a af 2f  28 6e 25 2b 31 29 29 0d  | .(i%*./(n%+1)).|
000001b0  00 78 11 20 66 28 69 25  29 20 3d 20 a0 20 66 24  |.x. f(i%) = . f$|
000001c0  0d 00 82 05 ed 0d 00 8c  11 e3 20 6a 25 3d 30 20  |.......... j%=0 |
000001d0  b8 20 6e 25 2b 31 0d 00  96 3d 20 e7 20 28 6a 25  |. n%+1...= . (j%|
000001e0  80 31 29 20 73 75 6d 20  3d 20 28 66 28 30 29 2d  |.1) sum = (f(0)-|
000001f0  66 28 6e 25 2b 31 29 29  2f 32 20 8b 20 73 75 6d  |f(n%+1))/2 . sum|
00000200  20 3d 20 28 66 28 30 29  2b 66 28 6e 25 2b 31 29  | = (f(0)+f(n%+1)|
00000210  29 2f 32 0d 00 a0 10 20  e3 20 69 25 3d 31 20 b8  |)/2.... . i%=1 .|
00000220  20 6e 25 0d 00 aa 24 20  20 73 75 6d 20 2b 3d 20  | n%...$  sum += |
00000230  66 28 69 25 29 2a 9b 28  69 25 2a 6a 25 2a af 2f  |f(i%)*.(i%*j%*./|
00000240  28 6e 25 2b 31 29 29 0d  00 b4 06 20 ed 0d 00 be  |(n%+1)).... ....|
00000250  17 20 63 28 6a 25 29 20  3d 20 73 75 6d 2f 28 6e  |. c(j%) = sum/(n|
00000260  25 2b 31 29 0d 00 c8 1c  20 e7 20 6a 25 3c 3e 30  |%+1).... . j%<>0|
00000270  20 63 28 6a 25 29 20 3d  20 32 2a 63 28 6a 25 29  | c(j%) = 2*c(j%)|
00000280  0d 00 d2 05 ed 0d 00 dc  14 45 20 3d 20 94 28 63  |.........E = .(c|
00000290  28 6e 25 2b 31 29 2f 32  29 0d 00 e6 26 f1 20 27  |(n%+1)/2)...&. '|
000002a0  22 50 6f 6c 79 6e 6f 6d  69 61 6c 20 61 70 70 72  |"Polynomial appr|
000002b0  6f 78 69 6d 61 74 69 6f  6e 20 69 73 3a 22 27 0d  |oximation is:"'.|
000002c0  00 f0 0f e3 20 69 25 3d  30 20 b8 20 6e 25 0d 00  |.... i%=0 . n%..|
000002d0  fa 23 20 f1 20 22 54 5b  22 3b 69 25 3b 22 5d 28  |.# . "T[";i%;"](|
000002e0  78 29 20 2a 22 2c 63 28  69 25 29 3b 8a 32 37 29  |x) *",c(i%);.27)|
000002f0  3b 0d 01 04 16 20 e7 20  69 25 3c 6e 25 20 f1 20  |;.... . i%<n% . |
00000300  22 2b 22 20 8b 20 f1 0d  01 0e 05 ed 0d 01 18 1d  |"+" . ..........|
00000310  f1 20 27 22 77 68 65 72  65 20 54 5b 30 5d 28 78  |. '"where T[0](x|
00000320  29 20 20 20 3d 20 31 3b  22 0d 01 22 1d f1 20 20  |)   = 1;".."..  |
00000330  22 20 20 20 20 20 20 54  5b 31 5d 28 78 29 20 20  |"      T[1](x)  |
00000340  20 3d 20 78 3b 22 0d 01  2c 39 f1 20 20 22 61 6e  | = x;"..,9.  "an|
00000350  64 20 20 20 54 5b 6e 2b  31 5d 28 78 29 20 3d 20  |d   T[n+1](x) = |
00000360  32 78 54 5b 6e 5d 28 78  29 2d 54 5b 6e 2d 31 5d  |2xT[n](x)-T[n-1]|
00000370  28 78 29 2c 20 66 6f 72  20 6e 3e 3d 31 22 27 0d  |(x), for n>=1"'.|
00000380  01 36 24 f4 20 4e 6f 77  20 77 65 20 6e 65 65 64  |.6$. Now we need|
00000390  20 74 6f 20 67 65 6e 20  74 68 65 20 54 20 70 6f  | to gen the T po|
000003a0  6c 79 73 0d 01 40 3d f4  20 6e 6f 74 65 20 77 65  |lys..@=. note we|
000003b0  20 75 73 65 20 74 28 69  2c 6a 29 20 74 6f 20 68  | use t(i,j) to h|
000003c0  6f 6c 64 20 74 68 65 20  63 6f 65 66 20 6f 66 20  |old the coef of |
000003d0  78 5e 6a 20 69 6e 20 70  6f 6c 79 20 54 5b 69 5d  |x^j in poly T[i]|
000003e0  0d 01 4a 0c 74 28 30 2c  30 29 3d 31 0d 01 54 0c  |..J.t(0,0)=1..T.|
000003f0  74 28 31 2c 31 29 3d 31  0d 01 5e 08 69 25 3d 31  |t(1,1)=1..^.i%=1|
00000400  0d 01 68 0c c8 95 20 69  25 3c 6e 25 0d 01 72 1b  |..h... i%<n%..r.|
00000410  20 74 28 69 25 2b 31 2c  30 29 20 3d 20 2d 74 28  | t(i%+1,0) = -t(|
00000420  69 25 2d 31 2c 30 29 0d  01 7c 14 20 e3 20 6a 25  |i%-1,0)..|. . j%|
00000430  20 3d 20 31 20 b8 20 69  25 2b 31 0d 01 86 2a 20  | = 1 . i%+1...* |
00000440  20 74 28 69 25 2b 31 2c  6a 25 29 20 3d 20 32 2a  | t(i%+1,j%) = 2*|
00000450  74 28 69 25 2c 6a 25 2d  31 29 2d 74 28 69 25 2d  |t(i%,j%-1)-t(i%-|
00000460  31 2c 6a 25 29 0d 01 90  06 20 ed 0d 01 9a 0a 20  |1,j%).... ..... |
00000470  69 25 2b 3d 31 0d 01 a4  05 ce 0d 01 ae 0f e3 20  |i%+=1.......... |
00000480  69 25 3d 30 20 b8 20 6e  25 0d 01 b8 0c 20 73 75  |i%=0 . n%.... su|
00000490  6d 20 3d 20 30 0d 01 c2  11 20 e3 20 6a 25 3d 69  |m = 0.... . j%=i|
000004a0  25 20 b8 20 6e 25 0d 01  cc 1c 20 20 73 75 6d 20  |% . n%....  sum |
000004b0  2b 3d 20 74 28 6a 25 2c  20 69 25 29 2a 63 28 6a  |+= t(j%, i%)*c(j|
000004c0  25 29 0d 01 d6 06 20 ed  0d 01 e0 10 20 70 28 69  |%).... ..... p(i|
000004d0  25 29 20 3d 20 73 75 6d  0d 01 ea 05 ed 0d 01 f4  |%) = sum........|
000004e0  24 f1 20 27 22 57 68 69  63 68 20 6d 61 79 20 61  |$. '"Which may a|
000004f0  6c 73 6f 20 62 65 20 77  72 69 74 74 65 6e 3a 22  |lso be written:"|
00000500  27 0d 01 fe 0f e3 20 69  25 3d 30 20 b8 20 6e 25  |'..... i%=0 . n%|
00000510  0d 02 08 20 20 f1 20 22  78 5e 22 3b 69 25 3b 22  |...  . "x^";i%;"|
00000520  20 20 2a 22 2c 70 28 69  25 29 3b 8a 32 37 29 3b  |  *",p(i%);.27);|
00000530  0d 02 12 16 20 e7 20 69  25 3c 6e 25 20 f1 20 22  |.... . i%<n% . "|
00000540  2b 22 20 8b 20 f1 0d 02  1c 05 ed 0d 02 26 10 6d  |+" . ........&.m|
00000550  61 78 65 72 72 6f 72 20  3d 20 30 0d 02 30 16 e3  |axerror = 0..0..|
00000560  20 78 3d 2d 31 20 b8 20  31 20 88 20 30 2e 30 30  | x=-1 . 1 . 0.00|
00000570  31 0d 02 3a 0e 20 76 20  3d 20 70 28 6e 25 29 0d  |1..:. v = p(n%).|
00000580  02 44 17 20 e3 20 69 25  3d 6e 25 2d 31 20 b8 20  |.D. . i%=n%-1 . |
00000590  30 20 88 20 2d 31 0d 02  4e 13 20 20 76 20 3d 20  |0 . -1..N.  v = |
000005a0  76 2a 78 2b 70 28 69 25  29 0d 02 58 06 20 ed 0d  |v*x+p(i%)..X. ..|
000005b0  02 62 14 20 65 20 3d 20  94 28 76 20 2d 20 a0 20  |.b. e = .(v - . |
000005c0  66 24 29 0d 02 6c 1c 20  e7 20 65 3e 6d 61 78 65  |f$)..l. . e>maxe|
000005d0  72 72 6f 72 20 6d 61 78  65 72 72 6f 72 3d 65 0d  |rror maxerror=e.|
000005e0  02 76 05 ed 0d 02 80 45  f1 20 27 27 22 46 6f 72  |.v.....E. ''"For|
000005f0  20 74 68 69 73 20 61 70  70 72 6f 78 69 6d 61 74  | this approximat|
00000600  69 6f 6e 20 61 6e 20 75  70 70 65 72 20 62 6f 75  |ion an upper bou|
00000610  6e 64 20 6f 6e 20 65 72  72 6f 72 20 69 73 20 22  |nd on error is "|
00000620  3b 6d 61 78 65 72 72 6f  72 0d 02 8a 59 f1 20 20  |;maxerror...Y.  |
00000630  20 22 41 6c 73 6f 20 77  65 20 6e 6f 77 20 6b 6e  | "Also we now kn|
00000640  6f 77 20 6d 69 6e 69 6d  61 78 20 65 72 72 6f 72  |ow minimax error|
00000650  20 70 20 28 74 68 65 20  65 72 72 6f 72 20 6d 61  | p (the error ma|
00000660  64 65 20 62 79 20 74 68  65 20 62 65 73 74 20 70  |de by the best p|
00000670  6f 73 73 69 62 6c 65 20  6f 72 64 65 72 20 22 3b  |ossible order ";|
00000680  6e 25 0d 02 94 24 f1 20  20 20 22 70 6f 6c 79 6e  |n%...$.   "polyn|
00000690  6f 6d 69 61 6c 29 20 6c  69 65 73 20 62 65 74 77  |omial) lies betw|
000006a0  65 65 6e 3a 22 27 0d 02  9e 1f f1 20 45 3b 22 20  |een:"'..... E;" |
000006b0  20 61 6e 64 20 20 22 3b  6d 61 78 65 72 72 6f 72  | and  ";maxerror|
000006c0  3b 22 2e 22 27 0d 02 a8  3c f1 20 22 54 68 69 73  |;"."'...<. "This|
000006d0  20 69 6e 64 69 63 61 74  65 73 20 68 6f 77 20 63  | indicates how c|
000006e0  6c 6f 73 65 20 6f 75 72  20 61 70 70 72 6f 78 69  |lose our approxi|
000006f0  6d 61 74 69 6f 6e 20 69  73 20 74 6f 20 69 74 2e  |mation is to it.|
00000700  22 0d 02 b2 56 f1 20 27  27 22 46 69 6e 61 6c 6c  |"...V. ''"Finall|
00000710  79 2c 20 69 66 20 6f 6e  6c 79 20 6e 65 65 64 20  |y, if only need |
00000720  73 61 79 20 31 36 20 62  69 74 20 61 63 63 75 72  |say 16 bit accur|
00000730  61 63 79 2c 20 6e 6f 74  65 20 74 68 61 74 20 6f  |acy, note that o|
00000740  75 72 20 70 6f 6c 79 27  73 20 63 6f 65 66 66 69  |ur poly's coeffi|
00000750  63 69 65 6e 74 73 22 0d  02 bc 3e f1 20 20 20 22  |cients"...>.   "|
00000760  6d 75 6c 74 69 70 6c 69  65 64 20 62 79 20 32 5e  |multiplied by 2^|
00000770  32 30 20 26 20 77 72 69  74 74 65 6e 20 69 6e 20  |20 & written in |
00000780  68 65 78 20 61 72 65 20  61 73 20 66 6f 6c 6c 6f  |hex are as follo|
00000790  77 73 3a 22 27 0d 02 c6  14 e3 20 69 25 3d 6e 25  |ws:"'..... i%=n%|
000007a0  20 b8 20 30 20 88 20 2d  31 0d 02 d0 0e 20 76 20  | . 0 . -1.... v |
000007b0  3d 20 70 28 69 25 29 0d  02 da 1f 20 e7 20 76 3c  |= p(i%).... . v<|
000007c0  30 20 73 24 3d 22 2d 22  3a 76 3d 2d 76 20 8b 20  |0 s$="-":v=-v . |
000007d0  73 24 3d 22 2b 22 0d 02  e4 17 20 76 32 20 3d 20  |s$="+".... v2 = |
000007e0  76 2a 31 30 34 38 35 37  36 2b 30 2e 35 0d 02 ee  |v*1048576+0.5...|
000007f0  31 20 f1 20 22 78 5e 22  3b 69 25 3b 22 20 20 2a  |1 . "x^";i%;"  *|
00000800  22 2c 73 24 3b c2 22 30  30 30 30 30 30 30 30 22  |",s$;."00000000"|
00000810  2b c3 7e 76 32 2c 38 29  3b 8a 32 37 29 3b 0d 02  |+.~v2,8);.27);..|
00000820  f8 15 20 e7 20 69 25 3e  30 20 f1 20 22 2b 22 20  |.. . i%>0 . "+" |
00000830  8b 20 f1 0d 03 02 05 ed  0d ff                    |. ........|
0000083a