66from ..linear_model import lineardiff , polydiff , savgoldiff , spectraldiff
77from ..total_variation_regularization import velocity , acceleration , jerk , iterative_velocity
88from ..kalman_smooth import * # constant_velocity, constant_acceleration, constant_jerk, known_dynamics
9- from ..smooth_finite_difference import mediandiff , meandiff , gaussiandiff , friedrichsdiff , butterdiff , splinediff
9+ from ..smooth_finite_difference import * # mediandiff, meandiff, gaussiandiff, friedrichsdiff, butterdiff, splinediff
1010from ..finite_difference import first_order , second_order
1111# Function aliases for testing cases where parameters change the behavior in a big way
12- def iterated_first_order ( * args , ** kwargs ): return first_order (* args , ** kwargs )
12+ iterated_first_order = lambda * args , ** kwargs : first_order (* args , ** kwargs )
1313
1414dt = 0.1
15- t = np .arange (0 , 3 + dt , dt ) # sample locations, including the endpoint
15+ t = np .arange (0 , 3 , dt ) # sample locations
1616tt = np .linspace (0 , 3 ) # full domain, for visualizing denser plots
1717np .random .seed (7 ) # for repeatability of the test, so we don't get random failures
1818noise = 0.05 * np .random .randn (* t .shape )
@@ -22,117 +22,70 @@ def iterated_first_order(*args, **kwargs): return first_order(*args, **kwargs)
2222 (0 , r"$x(t)=1$" , lambda t : np .ones (t .shape ), lambda t : np .zeros (t .shape )), # constant
2323 (1 , r"$x(t)=2t+1$" , lambda t : 2 * t + 1 , lambda t : 2 * np .ones (t .shape )), # affine
2424 (2 , r"$x(t)=t^2-t+1$" , lambda t : t ** 2 - t + 1 , lambda t : 2 * t - 1 ), # quadratic
25- (3 , r"$x(t)=\sin(3t)+1/2$" , lambda t : np .sin (3 * t ) + 1 / 2 , lambda t : 3 * np .cos (3 * t )), # sinuoidal
25+ (3 , r"$x(t)=\sin(3t)+1/2$" , lambda t : np .sin (3 * t ) + 1 / 2 , lambda t : 3 * np .cos (3 * t )), # sinuoidal
2626 (4 , r"$x(t)=e^t\sin(5t)$" , lambda t : np .exp (t )* np .sin (5 * t ), # growing sinusoidal
2727 lambda t : np .exp (t )* (5 * np .cos (5 * t ) + np .sin (5 * t ))),
2828 (5 , r"$x(t)=\frac{\sin(8t)}{(t+0.1)^{3/2}}$" , lambda t : np .sin (8 * t )/ ((t + 0.1 )** (3 / 2 )), # steep challenger
2929 lambda t : ((0.8 + 8 * t )* np .cos (8 * t ) - 1.5 * np .sin (8 * t ))/ (0.1 + t )** (5 / 2 ))]
3030
3131# Call both ways, with kwargs (new) and with params list with default options dict (legacy), to ensure both work
3232diff_methods_and_params = [
33- (first_order , {}), # empty dictionary for the case of no parameters. no params -> no diff in new vs old
34- (iterated_first_order , {'num_iterations' :5 }), (iterated_first_order , [5 ], {'iterate' :True }),
35- (second_order , {}),
33+ (first_order , None ), (iterated_first_order , {'num_iterations' :5 }),
34+ (second_order , None ),
3635 #(lineardiff, {'order':3, 'gamma':5, 'window_size':10, 'solver':'CVXOPT'}),
3736 (polydiff , {'polynomial_order' :2 , 'window_size' :3 }), (polydiff , [2 , 3 ]),
3837 (savgoldiff , {'polynomial_order' :2 , 'window_size' :4 , 'smoothing_win' :4 }), (savgoldiff , [2 , 4 , 4 ]),
39- (spectraldiff , {'high_freq_cutoff' :0.1 }), (spectraldiff , [0.1 ]),
40- (mediandiff , {'window_size' :3 , 'num_iterations' :2 }), (mediandiff , [3 , 2 ], {'iterate' :True }),
41- (meandiff , {'window_size' :3 , 'num_iterations' :2 }), (meandiff , [3 , 2 ], {'iterate' :True }),
42- (gaussiandiff , {'window_size' :5 }), (gaussiandiff , [5 ]),
43- (friedrichsdiff , {'window_size' :5 }), (friedrichsdiff , [5 ]),
44- (butterdiff , {'filter_order' :3 , 'cutoff_freq' :0.074 }), (butterdiff , [3 , 0.074 ]),
45- (splinediff , {'order' :5 , 's' :2 }), (splinediff , [5 , 2 ])
38+ (spectraldiff , {'high_freq_cutoff' :0.1 }), (spectraldiff , [0.1 ])
4639 ]
47- diff_methods_and_params = [(velocity , [0.5 ], {'solver' : 'CVXOPT' })]
4840
4941# All the testing methodology follows the exact same pattern; the only thing that changes is the
5042# closeness to the right answer various methods achieve with the given parameterizations. So index a
5143# big ol' table by the method, then the test function, then the pair of quantities we're comparing.
5244error_bounds = {
5345 first_order : [[(- 25 , - 25 ), (- 25 , - 25 ), (0 , 0 ), (1 , 1 )],
54- [(- 25 , - 25 ), (- 13 , - 14 ), (0 , 0 ), (1 , 1 )],
46+ [(- 25 , - 25 ), (- 14 , - 14 ), (0 , 0 ), (1 , 1 )],
5547 [(- 25 , - 25 ), (0 , 0 ), (0 , 0 ), (1 , 0 )],
5648 [(- 25 , - 25 ), (0 , 0 ), (0 , 0 ), (1 , 1 )],
5749 [(- 25 , - 25 ), (2 , 2 ), (0 , 0 ), (2 , 2 )],
5850 [(- 25 , - 25 ), (3 , 3 ), (0 , 0 ), (3 , 3 )]],
59- iterated_first_order : [[(- 8 , - 9 ), (- 11 , - 11 ), (0 , - 1 ), (0 , 0 )],
60- [(- 6 , - 6 ), (- 6 , - 7 ), (0 , - 1 ), (0 , 0 )],
51+ iterated_first_order : [[(- 7 , - 7 ), (- 9 , - 10 ), (0 , - 1 ), (0 , 0 )],
52+ [(- 5 , - 5 ), (- 5 , - 6 ), (0 , - 1 ), (0 , 0 )],
6153 [(- 1 , - 1 ), (0 , 0 ), (0 , - 1 ), (0 , 0 )],
62- [(0 , 0 ), (1 , 0 ), (0 , 0 ), (1 , 0 )],
54+ [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
6355 [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
6456 [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
6557 second_order : [[(- 25 , - 25 ), (- 25 , - 25 ), (0 , 0 ), (1 , 1 )],
66- [(- 25 , - 25 ), (- 13 , - 13 ), (0 , 0 ), (1 , 1 )],
67- [(- 25 , - 25 ), (- 13 , - 13 ), (0 , 0 ), (1 , 1 )],
58+ [(- 25 , - 25 ), (- 14 , - 14 ), (0 , 0 ), (1 , 1 )],
59+ [(- 25 , - 25 ), (- 13 , - 14 ), (0 , 0 ), (1 , 1 )],
6860 [(- 25 , - 25 ), (0 , - 1 ), (0 , 0 ), (1 , 1 )],
6961 [(- 25 , - 25 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
7062 [(- 25 , - 25 ), (3 , 3 ), (0 , 0 ), (3 , 3 )]],
7163 #lineardiff: [TBD when #91 is solved],
7264 polydiff : [[(- 14 , - 15 ), (- 14 , - 14 ), (0 , - 1 ), (1 , 1 )],
7365 [(- 14 , - 14 ), (- 13 , - 13 ), (0 , - 1 ), (1 , 1 )],
74- [(- 14 , - 14 ), (- 13 , - 13 ), (0 , - 1 ), (1 , 1 )],
66+ [(- 14 , - 15 ), (- 13 , - 14 ), (0 , - 1 ), (1 , 1 )],
7567 [(- 2 , - 2 ), (0 , 0 ), (0 , - 1 ), (1 , 1 )],
76- [(0 , 0 ), (1 , 1 ), (0 , - 1 ), (1 , 1 )],
68+ [(0 , 0 ), (2 , 1 ), (0 , 0 ), (2 , 1 )],
7769 [(0 , 0 ), (3 , 3 ), (0 , 0 ), (3 , 3 )]],
78- savgoldiff : [[(- 9 , - 10 ), (- 13 , - 14 ), (0 , - 1 ), (0 , 0 )],
79- [(- 9 , - 10 ), (- 13 , - 13 ), (0 , - 1 ), (0 , 0 )],
70+ savgoldiff : [[(- 7 , - 8 ), (- 13 , - 14 ), (0 , - 1 ), (0 , 0 )],
71+ [(- 7 , - 8 ), (- 13 , - 13 ), (0 , - 1 ), (0 , 0 )],
8072 [(- 1 , - 1 ), (0 , - 1 ), (0 , - 1 ), (0 , 0 )],
8173 [(0 , - 1 ), (0 , 0 ), (0 , - 1 ), (1 , 0 )],
8274 [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
8375 [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
84- spectraldiff : [[(- 9 , - 10 ), (- 14 , - 15 ), (- 1 , - 1 ), (0 , 0 )],
85- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
86- [(1 , 0 ), (1 , 1 ), (1 , 1 ), (1 , 1 )],
76+ spectraldiff : [[(- 7 , - 8 ), (- 14 , - 15 ), (- 1 , - 1 ), (0 , 0 )],
8777 [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
88- [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
89- [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
90- mediandiff : [[(- 25 , - 25 ), (- 25 , - 25 ), (- 1 , - 1 ), (0 , 0 )],
91- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
92- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
93- [(- 1 , - 1 ), (0 , 0 ), (0 , 0 ), (1 , 1 )],
94- [(0 , 0 ), (2 , 2 ), (0 , 0 ), (2 , 2 )],
95- [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
96- meandiff : [[(- 25 , - 25 ), (- 25 , - 25 ), (0 , - 1 ), (0 , 0 )],
97- [(0 , 0 ), (1 , 0 ), (0 , 0 ), (1 , 1 )],
98- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
99- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
100- [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
101- [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
102- gaussiandiff : [[(- 14 , - 15 ), (- 14 , - 14 ), (0 , - 1 ), (1 , 0 )],
103- [(- 1 , - 1 ), (0 , 0 ), (0 , 0 ), (1 , 0 )],
78+ [(1 , 0 ), (1 , 1 ), (1 , 0 ), (1 , 1 )],
10479 [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
105- [(0 , - 1 ), (1 , 0 ), (0 , 0 ), (1 , 1 )],
10680 [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
107- [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
108- friedrichsdiff : [[(- 25 , - 25 ), (- 25 , - 25 ), (0 , - 1 ), (0 , 0 )],
109- [(- 1 , - 1 ), (0 , 0 ), (0 , 0 ), (1 , 0 )],
110- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
111- [(0 , - 1 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
112- [(1 , 1 ), (2 , 2 ), (1 , 1 ), (2 , 2 )],
113- [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]],
114- butterdiff : [[(- 13 , - 14 ), (- 13 , - 13 ), (0 , - 1 ), (0 , 0 )],
115- [(0 , - 1 ), (0 , 0 ), (0 , - 1 ), (0 , 0 )],
116- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
117- [(1 , 0 ), (1 , 1 ), (1 , 0 ), (1 , 1 )],
118- [(2 , 2 ), (3 , 2 ), (2 , 2 ), (3 , 2 )],
119- [(2 , 1 ), (3 , 3 ), (2 , 1 ), (3 , 3 )]],
120- splinediff : [[(- 14 , - 15 ), (- 14 , - 14 ), (- 1 , - 1 ), (0 , 0 )],
121- [(- 14 , - 14 ), (- 13 , - 14 ), (- 1 , - 1 ), (0 , 0 )],
122- [(- 14 , - 14 ), (0 , 0 ), (- 1 , - 1 ), (0 , 0 )],
123- [(0 , 0 ), (1 , 1 ), (0 , 0 ), (1 , 1 )],
124- [(1 , 0 ), (2 , 2 ), (1 , 0 ), (2 , 2 )],
125- [(1 , 0 ), (3 , 3 ), (1 , 0 ), (3 , 3 )]]
81+ [(1 , 1 ), (3 , 3 ), (1 , 1 ), (3 , 3 )]]
12682}
12783
128- # Essentially run the cartesian product of [diff methods] x [test functions] through this one test
12984@mark .filterwarnings ("ignore::DeprecationWarning" ) # I want to test the old and new functionality intentionally
13085@mark .parametrize ("diff_method_and_params" , diff_methods_and_params )
13186@mark .parametrize ("test_func_and_deriv" , test_funcs_and_derivs )
13287def test_diff_method (diff_method_and_params , test_func_and_deriv , request ): # request gives access to context
133- # unpack
134- diff_method , params = diff_method_and_params [:2 ]
135- if len (diff_method_and_params ) == 3 : options = diff_method_and_params [2 ] # optionally pass old-style `options` dict
88+ diff_method , params = diff_method_and_params # unpack
13689 i , latex_name , f , df = test_func_and_deriv
13790
13891 # some methods rely on cvxpy, and we'd like to allow use of pynumdiff without convex optimization
@@ -144,30 +97,27 @@ def test_diff_method(diff_method_and_params, test_func_and_deriv, request): # re
14497 x_noisy = x + noise # add a little noise
14598 dxdt = df (t ) # true values of the derivative at samples
14699
147- # differentiate without and with noise, accounting for new and old styles of calling functions
148- x_hat , dxdt_hat = diff_method (x , dt , ** params ) if isinstance (params , dict ) \
149- else diff_method (x , dt , params ) if (isinstance (params , list ) and len (diff_method_and_params ) < 3 ) \
150- else diff_method (x , dt , params , options )
100+ # differentiate without and with noise
101+ x_hat , dxdt_hat = diff_method (x , dt , ** params ) if isinstance (params , dict ) else diff_method (x , dt , params ) \
102+ if isinstance (params , list ) else diff_method (x , dt )
151103 x_hat_noisy , dxdt_hat_noisy = diff_method (x_noisy , dt , ** params ) if isinstance (params , dict ) \
152- else diff_method (x_noisy , dt , params ) if (isinstance (params , list ) and len (diff_method_and_params ) < 3 ) \
153- else diff_method (x_noisy , dt , params , options )
104+ else diff_method (x_noisy , dt , params ) if isinstance (params , list ) else diff_method (x_noisy , dt )
154105
155106 # check x_hat and x_hat_noisy are close to x and that dxdt_hat and dxdt_hat_noisy are close to dxdt
156- print ("]\n [" , end = "" )
107+ # print("]\n[", end="")
157108 for j ,(a ,b ) in enumerate ([(x ,x_hat ), (dxdt ,dxdt_hat ), (x ,x_hat_noisy ), (dxdt ,dxdt_hat_noisy )]):
158109 l2_error = np .linalg .norm (a - b )
159110 linf_error = np .max (np .abs (a - b ))
160- #print(f"({l2_error},{linf_error})", end=", ")
161- print (f"({ int (np .ceil (np .log10 (l2_error ))) if l2_error > 0 else - 25 } , { int (np .ceil (np .log10 (linf_error ))) if linf_error > 0 else - 25 } )" , end = ", " )
162111
163- # log_l2_bound, log_linf_bound = error_bounds[diff_method][i][j]
164- # assert np.linalg.norm(a - b) < 10**log_l2_bound
165- # assert np.max(np.abs(a - b)) < 10**log_linf_bound
166- # if 0 < np.linalg.norm(a - b) < 10**(log_l2_bound - 1) or 0 < np.max(np.abs(a - b)) < 10**(log_linf_bound - 1):
167- # print(f"Improvement detected for method {diff_method.__name__}")
112+ #print(f"({int(np.ceil(np.log10(l2_error))) if l2_error> 0 else -25}, {int(np.ceil(np.log10(linf_error))) if linf_error > 0 else -25})", end=", ")
113+ log_l2_bound , log_linf_bound = error_bounds [diff_method ][i ][j ]
114+ assert np .linalg .norm (a - b ) < 10 ** log_l2_bound
115+ assert np .max (np .abs (a - b )) < 10 ** log_linf_bound
116+ if 0 < np .linalg .norm (a - b ) < 10 ** (log_l2_bound - 1 ) or 0 < np .max (np .abs (a - b )) < 10 ** (log_linf_bound - 1 ):
117+ print (f"Improvement detected for method { str (diff_method )} " )
168118
169- if request .config .getoption ("--plot" ) and not isinstance ( params , list ) : # Get the plot flag from pytest configuration
170- fig , axes = request .config .plots [diff_method ] # get the appropriate plot, set up by the store_plots fixture in conftest.py
119+ if request .config .getoption ("--plot" ): # Get the plot flag from pytest configuration
120+ fig , axes = request .config .plots [diff_method ]
171121 axes [i , 0 ].plot (t , f (t ), label = r"$x(t)$" )
172122 axes [i , 0 ].plot (t , x , 'C0+' )
173123 axes [i , 0 ].plot (tt , df (tt ), label = r"$\frac{dx(t)}{dt}$" )
0 commit comments