File: rtest_interpol.mac

package info (click to toggle)
maxima 5.47.0-9
  • links: PTS
  • area: main
  • in suites: forky, sid
  • size: 193,104 kB
  • sloc: lisp: 434,678; fortran: 14,665; tcl: 10,990; sh: 4,577; makefile: 2,763; ansic: 447; java: 328; python: 262; perl: 201; xml: 60; awk: 28; sed: 15; javascript: 2
file content (154 lines) | stat: -rw-r--r-- 5,721 bytes parent folder | download | duplicates (6)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
(kill (all),
 load (interpol),
 float_approx_equal_tolerance : 1e-8,
 dat1: matrix([6.5,3],[1.8,-8],[9,2],[-2,4],[1/8,1],[9/12,3],[10,%pi]),
 dat2: [[2,4],[7,1],[-7,2/7],[5,9],[8,4/9],[9,4]],
 dat3: [8,5.9,4.1,8.6,9.7,3,99.8],
 0);
0;


/* Lagrange interpolation */

lagrange(dat1);
3.178762054958412E-5*%pi*(x-9)*(x-6.5)*(x-1.8)*(x-3/4)*(x-1/8)*(x+2)
-1.3795660402574614E-4*(x-10)*(x-6.5)*(x-1.8)*(x-3/4)*(x-1/8)*(x+2)
+2.3412532015037063E-4*(x-10)*(x-9)*(x-1.8)*(x-3/4)*(x-1/8)*(x+2)
+.004313780800478495*(x-10)*(x-9)*(x-6.5)*(x-3/4)*(x-1/8)*(x+2)
+.003788399045316041*(x-10)*(x-9)*(x-6.5)*(x-1.8)*(x-1/8)*(x+2)
-8.045639731099019E-4*(x-10)*(x-9)*(x-6.5)*(x-1.8)*(x-3/4)*(x+2)
+1.605431979101289E-4*(x-10)*(x-9)*(x-6.5)*(x-1.8)*(x-3/4)*(x-1/8);

lagrange(dat2);
((x-8)*(x-7)*(x-5)*(x-2)*(x+7))/224
-2*(x-9)*(x-7)*(x-5)*(x-2)*(x+7)/1215
+(x-9)*(x-8)*(x-5)*(x-2)*(x+7)/280
+(x-9)*(x-8)*(x-7)*(x-2)*(x+7)/-96
+2*(x-9)*(x-8)*(x-7)*(x-5)*(x+7)/2835
+(x-9)*(x-8)*(x-7)*(x-5)*(x-2)/-1270080;

lagrange(dat3, varname=z);
.1386111111111111*(z-6)*(z-5)*(z-4)*(z-3)*(z-2)*(z-1)
-((z-7)*(z-5)*(z-4)*(z-3)*(z-2)*(z-1))/40
+.2020833333333333*(z-7)*(z-6)*(z-4)*(z-3)*(z-2)*(z-1)
-.2388888888888889*(z-7)*(z-6)*(z-5)*(z-3)*(z-2)*(z-1)
+.08541666666666665*(z-7)*(z-6)*(z-5)*(z-4)*(z-2)*(z-1)
-.04916666666666667*(z-7)*(z-6)*(z-5)*(z-4)*(z-3)*(z-1)
-((2-z)*(z-7)*(z-6)*(z-5)*(z-4)*(z-3))/90;

/* Linear interpolation */

linearinterpol(dat1);
(20/17-24*x/17)*charfun(x<1/8)
+(%pi*x-2*x-9*%pi+20)*charfun(9<=x)+(5.6-0.4*x)*charfun(6.5<=x and x<9)
+(2.340425531914894*x-12.21276595744681)*charfun(1.8<=x and x<6.5)
+(10.85714285714286-10.47619047619048*x)*charfun(3/4<=x and x<1.8)
+(16*x/5+3/5)*charfun(1/8<=x and x<3/4);

linearinterpol(dat2);
(26*x/63+200/63)*charfun(x<2)+((32*x)/9-28)*charfun(8<=x)
+(44/9-5*x/9)*charfun(7<=x and x<8)
+(29-4*x)*charfun(5<=x and x<7)
+(5*x/3+2/3)*charfun(2<=x and x<5);

linearinterpol(dat3,varname=w);
(10.1-2.1*w)*charfun(w<2)+(96.8*w-577.8)
*charfun(6<=w)+(43.2-6.7*w)*charfun(5<=w and w<6)
+(1.1*w+4.2)*charfun(4<=w and w<5)
+(4.5*w-9.4)*charfun(3<=w and w<4)+(9.5-1.8*w)*
charfun(2<=w and w<3);

/* Cubic splines interpolation */

cspline(dat1);
(2.233410484248884E-4*%pi*x^3+0.658739927100864*x^3
+0.00134004629054933*%pi*x^2+3.952439562605184*x^2
+.001671568159305024*%pi*x+3.518491936013175*x
-2.3032045618816613E-4*%pi+.4971450384125281)
*charfun(x<1/8)
+(-0.153921193200841*%pi*x^3-.06508981912631977*x^3
+4.61763579602523*%pi*x^2+1.952694573789593*x^2
-45.02243676705146*%pi*x-21.46185591876961*x+143.3819812688326*%pi
+84.43892093505657)*charfun(9<=x)
+(.07396021366527837*%pi*x^3+.3549365250720062*x^3
-1.535162189359993*%pi*x^2-9.388016719565206*x^2
+10.35274510141554*%pi*x+80.60454572142358*x-22.74356433656842*%pi
-221.760283985523)*charfun(6.5<=x and x<9)
+(-0.00936633903933342*%pi*x^3-.5795316182501081*x^3
+.08970558837993714*%pi*x^2+8.834112075216023*x^2
-.2088954538939998*%pi*x-37.83929144465441*x
+.1399901999355959*%pi+34.86802987431265)*charfun(1.8<=x and x<6.5)
+(0.0163989810743406*%pi*x^3+6.51539068532972*x^3
-.04942714023390259*%pi*x^2-29.47846836411505*x^2
+.04154345761091165*%pi*x+31.12335334614151*x
-.01027314696735098*%pi-6.509557000164901)*charfun(3/4<=x and x<1.8)
+(-.007441723733517281*%pi*x^3-10.14309672394197*x^3
+.004214445583777643*%pi*x^2+8.003128306746245*x^2
+.001312268247651485*%pi*x+3.012155842995544*x
-2.153496265359353E-4*%pi+0.518242375621596)*charfun(1/8<=x and x<3/4);

cspline(dat2,d1=2,varname=r);
(635237*r^3/13492332+436937*r^2/1124361+333821*r/642492+6946987/6746166)*charfun(r<2)
+(-71261*r^3/83286+213783*r^2/9254-941613*r/4627+2720884/4627)
*charfun(8<=r)
+(13907*r^3/83286-39995*r^2/27762-14209*r/1983+383564/5949)
*charfun(7<=r and r<8)
+(81073*r^3/111048-491039*r^2/37016+1199465*r/15864-680383/5288)
*charfun(5<=r and r<7)
+(-165773*r^3/499716+110833*r^2/41643-670655*r/166572+1014497/249858)
*charfun(2<=r and r<5);

cspline(dat3,d1=8/5,dn=1/3,varname=k);
(2.410085470085469*k^3-13.34034188034188*k^2+21.05042735042734*k
-2.120170940170938)*charfun(k<2)
+(-119.4588034188034*k^3+2292.709401709402*k^2-14537.15418803419*k
+30491.4882051282)*charfun(6<=k)
+(61.94307692307693*k^3-972.5244444444445*k^2+5054.24888888889*k
-8691.31794871795)*charfun(5<=k and k<6)
+(-17.01350427350427*k^3+211.8242735042735*k^2-867.4947008547009*k
+1178.254700854701)*charfun(4<=k and k<5)
+(1.710940170940171*k^3-12.86905982905983*k^2+31.27863247863247*k
-20.10974358974358)*charfun(3<=k and k<4)
+(.4697435897435902*k^3-1.698290598290602*k^2-2.233675213675204*k
+13.4025641025641)*charfun(2<=k and k<3);

/* Rational interpolation */

float(ratinterpol(dat1,4));
(.1620328314185472*x^4-2.534071721796399*x^3+11.2288985461703*x^2
-7.58855393365858*x+5.924974869319725)/(x^2-6.159704297522396*x+5.901285404551714);

ratinterpol(dat2,3,varname=t);
(6519*t^3/19843-45257*t^2/19843-399142*t/19843+2751280/19843)
/(t^2-274651*t/19843+925960/19843);

ratinterpol(dat3,5,varname=t);
(-38471*t^5/35760+692581*t^4/35760-1525837*t^3/11920+13508819*t^2/35760
-8456297*t/17880+216979/1490)/(t-1253/149);

(kill (h),
 h:[[1,(-(13+2*b)*1)/12],
    [2,(-(15+2*b)*1)/7],
    [3,(-(17+2*b)*3)/16],
    [4,(-(19+2*b)*2)/9],
    [5,(-(21+2*b)*1)/4],
    [6,(-(23+2*b)*3)/11]],
 0);
0;

is (equal (ratinterpol(h,1),
	   ((32*b^5+1552*b^4+30064*b^3+291128*b^2+1412514*b+2766033)*x
	     /(32*b+16)+144)
	     /(x^4-(2*b+53)*x^3/2+(4*b^2+128*b+1283)*x^2/4
	       -(8*b^3+300*b^2+3974*b+19993)*x/8
	       -(80*b^4+3440*b^3+56240*b^2+418500*b+1260063)/(16*b+8))));
true;

is (equal (ratinterpol(h,2), (-x^2-(2*b+11)*x/2)/(x+5)));
true;

is (equal (ratinterpol (subst (b=5/7, h), 2), subst (b=5/7, ratinterpol (h, 2))));
true;

(reset(float_approx_equal_tolerance),0);
0$