1# -*- tcl -*-
2# interpolate.test --
3#    Test cases for the ::math::interpolate package
4#
5
6# -------------------------------------------------------------------------
7
8source [file join \
9	[file dirname [file dirname [file join [pwd] [info script]]]] \
10	devtools testutilities.tcl]
11
12testsNeedTcl     8.5
13testsNeedTcltest 2.1
14
15support {
16    use      struct/matrix.tcl struct::matrix
17    useLocal math.tcl          math
18}
19testing {
20    useLocal interpolate.tcl math::interpolate
21}
22
23# -------------------------------------------------------------------------
24
25#
26# Compare a list of numbers
27#
28proc matchNumbers {expected actual} {
29    set match 1
30    foreach a $actual e $expected {
31        if {$e != 0.0} {
32            if {abs($a-$e) > 0.5e-4*abs($a+$e)} {
33                set match 0
34                break
35            }
36        } else {
37            if {abs($a-$e) > 1.0e-5} {
38                set match 0
39                break
40            }
41        }
42    }
43    return $match
44}
45
46customMatch numbers matchNumbers
47
48# -------------------------------------------------------------------------
49
50#
51# Test cases: interpolation in tables
52#
53# Add a dummy row to the table - ticket b25b826973edcbb5b3a95f6c284214925a1d5e67
54# This makes it possible to use the same table in both 1D and 2D interpolations
55#
56set t [::math::interpolate::defineTable table1 \
57   {    x     v1    v2    v3 } \
58   {    -      1     2     3
59        0      0    10     1
60        1      1     9     4
61        2      2     8     9
62        5      5     5    25
63        7      7     3    49
64       10     10     0   100 }]
65
66test "Interpolate-1.1" "Interpolate in a one-dimensional table" \
67     -match numbers -body {
68   set result {}
69   foreach x { -1.0 0.0 3.0 5.0 9.9 11.0 } {
70      set result [concat $result \
71                 [::math::interpolate::interp-1d-table $t $x]]
72   }
73   set result
74} -result {
75   -1    0    10     1
76    0    0    10     1
77    3    3     7    14.333333
78    5    5     5    25
79  9.9  9.9     0.1  98.3
80   11   10     0   100 }
81
82
83# value = x+y
84set t2 [::math::interpolate::defineTable table2 \
85   {    x      y1   y2    y3 } \
86   {    -      0     3    10
87        1      1     4    11
88        2      2     5    12
89        5      5     8    15
90        7      7    10    17
91       10     10    13    20 }]
92
93test "Interpolate-1.2" "Interpolate in a two-dimensional table" \
94     -match numbers -body {
95   set result {}
96   foreach y { -1.0 0.0 3.0 5.0 9.9 11.0 } {
97      foreach x { -1.0 0.0 3.0 5.0 9.9 11.0 } {
98         set result [concat $result \
99            $x $y [::math::interpolate::interp-table $t2 $x $y]]
100      }
101   }
102   set result
103} -result {
104    -1.0 -1.0 1.0
105     0.0 -1.0 1.0
106     3.0 -1.0 3.0
107     5.0 -1.0 5.0
108     9.9 -1.0 9.9
109    11.0 -1.0 10.0
110    -1.0 0.0 1.0
111     0.0 0.0 1.0
112     3.0 0.0 3.0
113     5.0 0.0 5.0
114     9.9 0.0 9.9
115    11.0 0.0 10.0
116    -1.0 3.0 4.0
117     0.0 3.0 4.0
118     3.0 3.0 6.0
119     5.0 3.0 8.0
120     9.9 3.0 12.9
121    11.0 3.0 13.0
122    -1.0 5.0 6.0
123     0.0 5.0 6.0
124     3.0 5.0 8.0
125     5.0 5.0 10.0
126     9.9 5.0 14.9
127    11.0 5.0 15.0
128    -1.0 9.9 10.9
129     0.0 9.9 10.9
130     3.0 9.9 12.9
131     5.0 9.9 14.9
132     9.9 9.9 19.8
133    11.0 9.9 19.9
134    -1.0 11.0 11.0
135     0.0 11.0 11.0
136     3.0 11.0 13.0
137     5.0 11.0 15.0
138     9.9 11.0 19.9
139    11.0 11.0 20.0
140}
141
142test "Interpolate-1.3" "Interpolate with integers" \
143     -match numbers -body {
144    set equalResults {}
145
146    set table [::math::interpolate::defineTable table3 \
147                   {A B C D E} \
148                   {0   0.00000   0.00000   0.00000   0.0000
149                   1   0.52000   0.52000   0.52000   0.5200
150                   3   0.69831   0.63142   0.67758   0.68457
151                   5   0.86111   0.71690   0.81118   0.80365
152                   7   1.01367   0.78725   0.92891   0.89851}]
153
154    foreach A {2 4 6} A2 {2.0 4.0 6.0} {
155        set intResults   [::math::interpolate::interp-1d-table $table $A]
156        set floatResults [::math::interpolate::interp-1d-table $table $A2]
157        set equal 1
158        foreach i $intResults f $floatResults {
159            if { $i != $f } {
160                set equal 0
161                break
162            }
163       }
164       lappend equalResults $equal
165   }
166
167   return $equalResults
168} -result {1 1 1}
169
170# linear interpolation: y = x + 1 and y = 2*x, x<5, or 20-2*x, x>5
171
172test "Interpolate-2.1" "Linear interpolation - 1" \
173     -match numbers -body {
174   set result {}
175
176   set xyvalues { 0.0 1.0  10.0 11.0 }
177   foreach x { 0.0 4.0 7.0 10.0 101.0 } {
178      lappend result [::math::interpolate::interp-linear $xyvalues $x]
179   }
180   set result
181} -result { 1.0 5.0 8.0 11.0 11.0 }
182
183test "Interpolate-2.2" "Linear interpolation - 2" \
184     -match numbers -body {
185   set result {}
186   set xyvalues { 0.0 0.0  5.0 10.0 10.0 0.0 }
187   foreach x { 0.0 4.0 7.0 10.0 11.0 } {
188      lappend result [::math::interpolate::interp-linear $xyvalues $x]
189   }
190   set result
191} -result { 0.0 8.0 6.0 0.0 0.0 }
192
193# Lagrange interpolation: y = x + 1
194test "Interpolate-3.1" "Lagrange interpolation - 1" \
195     -match numbers -body {
196   set result {}
197   set xyvalues { 0.0 1.0  10.0 11.0 }
198   foreach x { 0.0 4.0 7.0 10.0 101.0 } {
199      lappend result [::math::interpolate::interp-lagrange $xyvalues $x]
200   }
201   set result
202} -result { 1.0  5.0  8.0  11.0  102.0 }
203
204
205#Lagrange interpolation (2) - expected: y=10-2*(x-5)**2/5
206test "Interpolate-3.2" "Lagrange interpolation - 2" \
207     -match numbers -body {
208   set result {}
209   set xyvalues { 0.0 0.0  5.0 10.0 10.0 0.0 }
210   foreach x { 0.0 4.0 7.0 10.0 11.0 } {
211      lappend result [::math::interpolate::interp-lagrange $xyvalues $x]
212   }
213   set result
214} -result { 0.0 9.6 8.4 0.0 -4.4 }
215
216# Spatial interpolation
217test "Interpolate-4.1" "Spatial interpolation - 1" \
218     -match numbers -body {
219   set result {}
220   set xyzvalues { {-1.0 0.0 -2.0 }
221                   { 1.0 0.0  2.0 } }
222   foreach coord { {0.0 0.0} {0.0 1.0} {3.0 0.0} {100.0 0.0} } {
223      lappend result [::math::interpolate::interp-spatial $xyzvalues $coord]
224   }
225   set result
226} -result { 0.0 0.0 1.2 0.039996 }
227
228test "Interpolate-4.2" "Spatial interpolation - 2" \
229   -match numbers -body {
230   set result {}
231
232   set xyzvalues { {-1.0 0.0 { -2.0  1.0 } }
233                   { 1.0 0.0 {  2.0 -1.0 } } }
234   foreach coord { {0.0 0.0} {0.0 1.0} {3.0 0.0} {100.0 0.0} } {
235      set result [concat $result \
236            [::math::interpolate::interp-spatial $xyzvalues $coord]]
237   }
238   set result
239} -result { 0.0       0.0
240            0.0       0.0
241            1.2      -0.6
242            0.039996 -0.019998 }
243
244test "Interpolate-4.3" "Spatial interpolation - 3 - coincident points" \
245   -match numbers -body {
246   set result {}
247
248   set xyzvalues { {-1.0 0.0 { -2.0  1.0 } }
249                   { 1.0 0.0 {  2.0 -1.0 } } }
250   set coord {-1.0 0.0}
251      set result [::math::interpolate::interp-spatial $xyzvalues $coord]
252
253   set result
254} -result { -2.0 1.0 }
255
256#
257# Test TODO: parameters for spatial interpolation
258#
259
260test interpolate-5.1 "neville algorithm" \
261    -body {
262	set problems {}
263	namespace import ::math::interpolate::neville
264	set xtable [list 0. 30. 45. 60. 90. 120. 135. 150. 180.]
265	set ytable [list 0. 0.5 [expr sqrt(0.5)] [expr sqrt(0.75)] 1. \
266			[expr sqrt(0.75)] [expr sqrt(0.5)] 0.5 0.]
267	for { set x -15 } { $x <= 195 } { incr x } {
268	    foreach { y error } [neville $xtable $ytable $x] break
269	    set diff [expr { abs( $y - sin( $x*3.1415926535897932/180. ) ) }]
270	    if { $error > 3.e-4 || ( $diff > $error && $diff > 1.e-8 ) } {
271		append problems \n "interpolating for sine of " $x " degrees" \
272		    \n "value was " $y " +/- " $error \
273		    \n "actual error was " $diff
274	    }
275	}
276	set problems
277    } \
278    -result {}
279
280proc matchNumbers {expected actual} {
281    set match 1
282    foreach a $actual e $expected {
283        if {abs($a-$e) > 0.1e-6} {
284            set match 0
285            break
286        }
287    }
288    return $match
289}
290
291customMatch numbers matchNumbers
292
293test "cubic-splines-1.0" "Interpolate linear function" \
294   -match numbers -body {
295    set xcoord {1 2 3 4 5}
296    set ycoord {1 2 3 4 5}
297    set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
298    set yvalues {}
299    foreach x {1.5 2.5 3.5 4.5} {
300        lappend yvalues [::math::interpolate::interp-cubic-splines $coeffs $x]
301    }
302    set yvalues
303} -result {1.5 2.5 3.5 4.5}
304
305test "cubic-splines-1.1" "Interpolate quadratic function" \
306   -match numbers -body {
307    set xcoord {1 2 3 4 5}
308    set ycoord {1 4 9 16 25}
309    set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
310    set yvalues {}
311    foreach x $xcoord {
312        lappend yvalues [::math::interpolate::interp-cubic-splines $coeffs $x]
313    }
314    set yvalues
315} -result {1 4 9 16 25}
316
317test "cubic-splines-1.2" "Interpolate arbitrary function" \
318   -match numbers -body {
319    set coeffs [::math::interpolate::prepare-cubic-splines \
320                  {0.1 0.3 0.4 0.8  1.0} \
321                  {1.0 2.1 2.2 4.11 4.12}]
322    set yvalues {}
323    foreach x {0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0} {
324        lappend yvalues [::math::interpolate::interp-cubic-splines $coeffs $x]
325    }
326    set yvalues
327} -result {1.0 1.6804411764705884 2.1 2.2 2.5380974264705882
328           3.1041911764705885 3.695689338235294 4.11 4.2099448529411765 4.12}
329
330test "cubic-splines-2.1" "Too few data" \
331-match glob -body {
332   set xcoord {1 2}
333   set ycoord {1 4}
334   set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
335} -result "At least *" -returnCodes error
336
337test "cubic-splines-2.2" "Unequal lengths" \
338-match glob -body {
339   set xcoord {1 2 4 5}
340   set ycoord {1 4 5 5 6}
341   set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
342} -result "Equal number *" -returnCodes error
343
344test "cubic-splines-2.3" "Not-ascending x-coordinates" \
345-match glob -body {
346   set xcoord {1 2 1.5}
347   set ycoord {1 4 5}
348   set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
349} -result "* ascending" -returnCodes error
350
351test "cubic-splines-2.4" "X too small" \
352-match glob -body {
353   set xcoord {1 2 3}
354   set ycoord {1 4 5}
355   set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
356   set yvalue [::math::interpolate::interp-cubic-splines $coeffs -1]
357} -result "* too small" -returnCodes error
358
359test "cubic-splines-2.5" "X too large" \
360-match glob -body {
361   set xcoord {1 2 3}
362   set ycoord {1 4 5}
363   set coeffs [::math::interpolate::prepare-cubic-splines $xcoord $ycoord]
364   set yvalue [::math::interpolate::interp-cubic-splines $coeffs 6]
365} -result "* too large" -returnCodes error
366
367
368
369# -------------------------------------------------------------------------
370testsuiteCleanup
371
372# Local Variables:
373# mode: tcl
374# End:
375