Improve numerical accuracy

docs/profiling/results/accuracy_gcc_15.2.0.csv

@@ -3,29 +3,29 @@
 # os: Alpine Linux
 library,filter_type,order,method,max_b_err,max_a_err,max_step_err
 constfilt,butterworth,1,zoh,0.000000e+00,0.000000e+00,0.000000e+00
-constfilt,butterworth,2,zoh,3.608225e-16,2.220446e-16,1.110223e-15
-constfilt,butterworth,3,zoh,9.020562e-16,1.776357e-15,7.771561e-15
-constfilt,butterworth,4,zoh,8.971990e-15,8.881784e-16,1.587619e-14
-constfilt,butterworth,5,zoh,6.945833e-15,3.552714e-15,1.809664e-14
-constfilt,butterworth,6,zoh,3.284352e-14,1.598721e-14,6.572520e-14
-constfilt,butterworth,7,zoh,1.100200e-11,4.796163e-14,1.254094e-12
-constfilt,butterworth,8,zoh,2.955275e-08,5.062617e-13,2.134253e-09
-constfilt,butterworth,9,zoh,7.577398e-05,8.526513e-13,3.692378e-06
-constfilt,butterworth,10,zoh,2.065719e-01,2.030731e-11,6.461003e-03
-constfilt,butterworth,11,zoh,2.543305e+02,3.566485e-08,4.685698e+00
-constfilt,butterworth,12,zoh,1.799631e+03,3.915801e-10,2.955113e+02
+constfilt,butterworth,2,zoh,4.163336e-17,2.220446e-16,1.221245e-15
+constfilt,butterworth,3,zoh,9.714451e-17,1.776357e-15,8.548717e-15
+constfilt,butterworth,4,zoh,1.110223e-15,5.329071e-15,1.487699e-14
+constfilt,butterworth,5,zoh,7.826855e-16,2.220446e-15,1.565414e-14
+constfilt,butterworth,6,zoh,1.077333e-13,8.260059e-14,1.327827e-13
+constfilt,butterworth,7,zoh,4.350617e-13,2.184919e-13,4.956036e-13
+constfilt,butterworth,8,zoh,1.781413e-12,7.887024e-13,1.174533e-12
+constfilt,butterworth,9,zoh,2.683072e-11,2.824407e-12,1.341810e-11
+constfilt,butterworth,10,zoh,2.428169e-10,4.609291e-11,8.981571e-11
+constfilt,butterworth,11,zoh,1.938709e-10,6.137313e-10,8.004775e-11
+constfilt,butterworth,12,zoh,5.023309e-09,3.361052e-09,9.453389e-10
 constfilt,butterworth,1,matchedz,5.551115e-17,0.000000e+00,2.220446e-16
-constfilt,butterworth,2,matchedz,2.775558e-17,2.220446e-16,1.110223e-15
-constfilt,butterworth,3,matchedz,1.387779e-17,1.776357e-15,1.998401e-15
-constfilt,butterworth,4,matchedz,3.122502e-17,1.776357e-15,9.103829e-15
-constfilt,butterworth,5,matchedz,4.163336e-17,4.884981e-15,1.709743e-14
-constfilt,butterworth,6,matchedz,2.168404e-17,4.884981e-14,1.030287e-13
-constfilt,butterworth,7,matchedz,1.170938e-17,9.414691e-14,3.456124e-13
-constfilt,butterworth,8,matchedz,6.505213e-19,2.664535e-14,9.876544e-13
-constfilt,butterworth,9,matchedz,6.179952e-18,6.927792e-13,1.509015e-12
-constfilt,butterworth,10,matchedz,1.523304e-16,5.579892e-11,1.363021e-11
-constfilt,butterworth,11,matchedz,2.168404e-19,1.847411e-13,4.028200e-11
-constfilt,butterworth,12,matchedz,9.662952e-17,3.253149e-10,1.482581e-10
+constfilt,butterworth,2,matchedz,0.000000e+00,2.220446e-16,1.332268e-15
+constfilt,butterworth,3,matchedz,1.387779e-17,4.440892e-16,1.998401e-15
+constfilt,butterworth,4,matchedz,3.816392e-17,1.776357e-15,4.440892e-15
+constfilt,butterworth,5,matchedz,1.040834e-17,3.552714e-15,1.487699e-14
+constfilt,butterworth,6,matchedz,3.469447e-18,1.776357e-15,5.384582e-14
+constfilt,butterworth,7,matchedz,2.602085e-18,7.105427e-15,3.296252e-13
+constfilt,butterworth,8,matchedz,2.168404e-19,5.329071e-15,9.885426e-13
+constfilt,butterworth,9,matchedz,7.589415e-19,4.440892e-15,1.743494e-12
+constfilt,butterworth,10,matchedz,1.626303e-19,7.105427e-15,7.760681e-12
+constfilt,butterworth,11,matchedz,5.421011e-20,1.421085e-14,4.277068e-11
+constfilt,butterworth,12,matchedz,1.490778e-19,1.847411e-13,9.285439e-11
 constfilt,butterworth,1,prewarp,8.326673e-17,2.220446e-16,2.220446e-16
 constfilt,butterworth,2,prewarp,1.249001e-16,0.000000e+00,2.331468e-15
 constfilt,butterworth,3,prewarp,1.283695e-16,2.442491e-15,1.776357e-15
@@ -38,26 +38,26 @@ constfilt,butterworth,9,prewarp,4.832289e-16,2.842171e-14,4.035883e-12
 constfilt,butterworth,10,prewarp,3.041296e-15,2.486900e-14,6.625700e-12
 constfilt,butterworth,11,prewarp,1.426024e-15,3.197442e-13,1.841038e-11
 constfilt,butterworth,12,prewarp,1.350920e-14,1.563194e-13,5.809664e-11
-constfilt,elliptic,2,zoh,2.395861e-13,2.842171e-14,1.409099e-13
-constfilt,elliptic,3,zoh,5.400819e-13,9.858780e-14,3.640976e-13
-constfilt,elliptic,4,zoh,5.393082e-13,1.074696e-13,8.994472e-13
-constfilt,elliptic,5,zoh,1.066647e-12,4.938272e-13,1.474598e-12
-constfilt,elliptic,6,zoh,1.760619e-12,1.481482e-12,1.819322e-12
-constfilt,elliptic,7,zoh,1.055545e-11,1.796252e-11,9.916956e-12
-constfilt,elliptic,8,zoh,2.154665e-09,2.641301e-09,7.933327e-10
-constfilt,elliptic,9,zoh,3.582365e-07,1.304602e-07,9.765907e-07
-constfilt,elliptic,10,zoh,3.461312e-03,3.213908e-06,6.243202e-03
-constfilt,elliptic,11,zoh,5.384745e-01,4.802802e-05,7.544774e-02
-constfilt,elliptic,12,zoh,5.474468e+02,5.003185e-04,2.311049e+01
+constfilt,elliptic,2,zoh,2.395861e-13,2.831069e-14,1.409099e-13
+constfilt,elliptic,3,zoh,5.401790e-13,9.658940e-14,3.639866e-13
+constfilt,elliptic,4,zoh,5.369767e-13,1.052491e-13,9.001133e-13
+constfilt,elliptic,5,zoh,1.077680e-12,4.902745e-13,1.481482e-12
+constfilt,elliptic,6,zoh,1.853295e-12,1.490363e-12,1.784017e-12
+constfilt,elliptic,7,zoh,1.118322e-11,1.780975e-11,1.006395e-11
+constfilt,elliptic,8,zoh,1.626981e-09,2.650339e-09,7.395398e-10
+constfilt,elliptic,9,zoh,7.139073e-08,1.304869e-07,1.850400e-08
+constfilt,elliptic,10,zoh,1.794838e-06,3.207371e-06,2.383539e-07
+constfilt,elliptic,11,zoh,2.479973e-05,4.805415e-05,1.939256e-06
+constfilt,elliptic,12,zoh,2.788489e-04,4.995743e-04,3.710091e-03
 constfilt,elliptic,2,matchedz,3.854481e-11,2.786660e-14,1.928402e-11
-constfilt,elliptic,3,matchedz,5.143178e-13,9.592327e-14,3.533007e-13
-constfilt,elliptic,4,matchedz,1.324132e-13,1.101341e-13,8.861800e-13
-constfilt,elliptic,5,matchedz,1.020281e-12,4.991563e-13,1.489919e-12
-constfilt,elliptic,6,matchedz,6.720215e-13,1.428191e-12,1.833089e-12
-constfilt,elliptic,7,matchedz,1.035810e-11,1.789147e-11,1.008205e-11
-constfilt,elliptic,8,matchedz,6.505636e-10,2.641826e-09,7.457389e-10
-constfilt,elliptic,9,matchedz,6.787615e-08,1.304656e-07,1.862001e-08
-constfilt,elliptic,10,matchedz,7.442555e-07,3.213346e-06,2.409790e-07
+constfilt,elliptic,3,matchedz,5.143108e-13,9.614531e-14,3.533285e-13
+constfilt,elliptic,4,matchedz,1.324114e-13,1.088019e-13,8.872902e-13
+constfilt,elliptic,5,matchedz,1.020448e-12,4.813927e-13,1.496137e-12
+constfilt,elliptic,6,matchedz,6.720492e-13,1.449507e-12,1.837530e-12
+constfilt,elliptic,7,matchedz,1.035799e-11,1.789147e-11,1.007983e-11
+constfilt,elliptic,8,matchedz,6.505633e-10,2.641826e-09,7.456755e-10
+constfilt,elliptic,9,matchedz,6.787615e-08,1.304655e-07,1.861996e-08
+constfilt,elliptic,10,matchedz,7.442555e-07,3.213341e-06,2.409786e-07
 constfilt,elliptic,11,matchedz,2.360767e-05,4.802825e-05,1.950004e-06
 constfilt,elliptic,12,matchedz,1.156955e-04,5.002874e-04,1.090689e-05
 constfilt,elliptic,2,prewarp,1.637024e-13,3.286260e-14,9.213463e-14

docs/profiling/results/compile_times_gcc_15.2.0.csv

@@ -2,54 +2,54 @@
 # compiler: g++ (Alpine 15.2.0) 15.2.0
 # os: Alpine Linux
 filter_type,order,method,compile_time_sec,max_rss_kb,exit_code
-butterworth,10,matchedz,0.37,73000,0
-butterworth,10,tustin,0.45,83112,0
-butterworth,10,zoh,1.07,158552,0
-butterworth,11,matchedz,0.50,88108,0
-butterworth,11,tustin,0.62,102420,0
-butterworth,11,zoh,1.42,203688,0
-butterworth,12,matchedz,0.63,105460,0
-butterworth,12,tustin,0.81,126272,0
-butterworth,12,zoh,1.93,264624,0
-butterworth,1,matchedz,0.07,36764,0
-butterworth,1,tustin,0.06,36776,0
-butterworth,1,zoh,0.07,37024,0
-butterworth,2,matchedz,0.07,37096,0
-butterworth,2,tustin,0.06,37268,0
-butterworth,2,zoh,0.07,38084,0
-butterworth,4,matchedz,0.09,39528,0
-butterworth,4,tustin,0.09,40116,0
-butterworth,4,zoh,0.14,44948,0
-butterworth,6,matchedz,0.12,44468,0
-butterworth,6,tustin,0.14,45600,0
-butterworth,6,zoh,0.27,60956,0
-butterworth,8,matchedz,0.21,54360,0
-butterworth,8,tustin,0.25,58428,0
-butterworth,8,zoh,0.54,94996,0
-butterworth,9,matchedz,0.28,62652,0
-butterworth,9,tustin,0.34,69284,0
-butterworth,9,zoh,0.77,123156,0
-elliptic,10,matchedz,0.72,112964,0
-elliptic,10,tustin,0.53,89548,0
-elliptic,10,zoh,1.13,163684,0
-elliptic,11,matchedz,0.95,140320,0
-elliptic,11,tustin,0.69,109080,0
-elliptic,11,zoh,1.52,210084,0
-elliptic,12,matchedz,1.30,181428,0
-elliptic,12,tustin,0.92,134184,0
-elliptic,12,zoh,2.03,271788,0
-elliptic,2,matchedz,0.11,40780,0
-elliptic,2,tustin,0.10,40608,0
-elliptic,2,zoh,0.11,41888,0
-elliptic,4,matchedz,0.14,45328,0
-elliptic,4,tustin,0.13,43936,0
-elliptic,4,zoh,0.18,49076,0
-elliptic,6,matchedz,0.22,53920,0
-elliptic,6,tustin,0.20,50064,0
-elliptic,6,zoh,0.31,64840,0
-elliptic,8,matchedz,0.40,74308,0
-elliptic,8,tustin,0.31,63908,0
-elliptic,8,zoh,0.61,100160,0
-elliptic,9,matchedz,0.53,89588,0
-elliptic,9,tustin,0.40,74688,0
-elliptic,9,zoh,0.84,128436,0
+butterworth,10,matchedz,0.07,38664,0
+butterworth,10,tustin,0.45,84572,0
+butterworth,10,zoh,0.54,96292,0
+butterworth,11,matchedz,0.07,38672,0
+butterworth,11,tustin,0.61,104092,0
+butterworth,11,zoh,0.72,119492,0
+butterworth,12,matchedz,0.07,38888,0
+butterworth,12,tustin,0.80,127832,0
+butterworth,12,zoh,0.95,146800,0
+butterworth,1,matchedz,0.06,37560,0
+butterworth,1,tustin,0.06,38080,0
+butterworth,1,zoh,0.07,38176,0
+butterworth,2,matchedz,0.06,37732,0
+butterworth,2,tustin,0.08,38680,0
+butterworth,2,zoh,0.07,39148,0
+butterworth,4,matchedz,0.06,38044,0
+butterworth,4,tustin,0.09,41548,0
+butterworth,4,zoh,0.10,42804,0
+butterworth,6,matchedz,0.06,38096,0
+butterworth,6,tustin,0.14,47264,0
+butterworth,6,zoh,0.16,50184,0
+butterworth,8,matchedz,0.07,38272,0
+butterworth,8,tustin,0.24,60016,0
+butterworth,8,zoh,0.30,66640,0
+butterworth,9,matchedz,0.07,38508,0
+butterworth,9,tustin,0.34,71000,0
+butterworth,9,zoh,0.40,80012,0
+elliptic,10,matchedz,0.28,56052,0
+elliptic,10,tustin,0.66,101920,0
+elliptic,10,zoh,0.76,114416,0
+elliptic,11,matchedz,0.27,54756,0
+elliptic,11,tustin,0.82,120820,0
+elliptic,11,zoh,0.94,135896,0
+elliptic,12,matchedz,0.29,57064,0
+elliptic,12,tustin,1.05,147304,0
+elliptic,12,zoh,1.20,166836,0
+elliptic,2,matchedz,0.16,46144,0
+elliptic,2,tustin,0.17,47260,0
+elliptic,2,zoh,0.17,47724,0
+elliptic,4,matchedz,0.19,49120,0
+elliptic,4,tustin,0.22,52524,0
+elliptic,4,zoh,0.23,53976,0
+elliptic,6,matchedz,0.21,50788,0
+elliptic,6,tustin,0.29,59564,0
+elliptic,6,zoh,0.31,62484,0
+elliptic,8,matchedz,0.25,53068,0
+elliptic,8,tustin,0.43,74600,0
+elliptic,8,zoh,0.48,81180,0
+elliptic,9,matchedz,0.26,54128,0
+elliptic,9,tustin,0.52,86224,0
+elliptic,9,zoh,0.59,95288,0

docs/profiling/results/runtime_gcc_15.2.0.csv

@@ -2,42 +2,42 @@
 # compiler: g++ (Alpine 15.2.0) 15.2.0
 # os: Alpine Linux
 library,filter_type,order,method,ns_per_sample,msa_per_s,dc_gain
-constfilt,Butterworth,1,tustin,6.415,155.9,1.00000000
-constfilt,Butterworth,2,tustin,10.876,91.9,1.00000000
-constfilt,Butterworth,4,tustin,18.954,52.8,1.00000000
-constfilt,Butterworth,6,tustin,26.164,38.2,1.00000000
-constfilt,Butterworth,8,tustin,34.550,28.9,1.00000000
-constfilt,Butterworth,1,zoh,6.420,155.8,1.00000000
-constfilt,Butterworth,2,zoh,10.826,92.4,1.00000000
-constfilt,Butterworth,4,zoh,18.535,54.0,1.00000000
-constfilt,Butterworth,6,zoh,26.066,38.4,1.00000000
-constfilt,Butterworth,8,zoh,34.415,29.1,1.00000000
-constfilt,Butterworth,1,matchedz,6.423,155.7,1.00000000
-constfilt,Butterworth,2,matchedz,10.839,92.3,1.00000000
-constfilt,Butterworth,4,matchedz,18.858,53.0,1.00000000
-constfilt,Butterworth,6,matchedz,26.130,38.3,1.00000000
-constfilt,Butterworth,8,matchedz,34.357,29.1,1.00000000
-constfilt,Elliptic,2,tustin,10.844,92.2,0.94406088
-constfilt,Elliptic,4,tustin,18.373,54.4,0.94406088
-constfilt,Elliptic,6,tustin,26.217,38.1,0.94406088
-constfilt,Elliptic,8,tustin,34.485,29.0,0.94406088
-constfilt,Elliptic,2,zoh,10.828,92.4,0.94406088
-constfilt,Elliptic,4,zoh,18.359,54.5,0.94406088
-constfilt,Elliptic,6,zoh,26.095,38.3,0.94406088
-constfilt,Elliptic,8,zoh,34.249,29.2,0.94406088
-constfilt,Elliptic,2,matchedz,11.091,90.2,0.94406088
-constfilt,Elliptic,4,matchedz,18.904,52.9,0.94406088
-constfilt,Elliptic,6,matchedz,26.024,38.4,0.94406088
-constfilt,Elliptic,8,matchedz,34.507,29.0,0.94406088
-iir1,Butterworth,2,runtime,10.732,93.2,1.00000000
-iir1,Butterworth,4,runtime,18.155,55.1,1.00000000
-iir1,Butterworth,6,runtime,30.092,33.2,1.00000000
-iir1,Butterworth,8,runtime,40.322,24.8,1.00000000
-kfr,Butterworth,2,runtime+simd,607.261,1.6,1.00000000
-kfr,Butterworth,4,runtime+simd,402.769,2.5,1.00000000
-kfr,Butterworth,6,runtime+simd,307.412,3.3,1.00000000
-kfr,Butterworth,8,runtime+simd,307.261,3.3,1.00000000
-kfr,Elliptic,2,runtime+simd,605.280,1.7,0.94406088
-kfr,Elliptic,4,runtime+simd,401.691,2.5,0.94406088
-kfr,Elliptic,6,runtime+simd,306.657,3.3,0.94406088
-kfr,Elliptic,8,runtime+simd,308.208,3.2,0.94406088
+constfilt,Butterworth,1,tustin,6.406,156.1,1.00000000
+constfilt,Butterworth,2,tustin,10.824,92.4,1.00000000
+constfilt,Butterworth,4,tustin,18.803,53.2,1.00000000
+constfilt,Butterworth,6,tustin,26.040,38.4,1.00000000
+constfilt,Butterworth,8,tustin,34.177,29.3,1.00000000
+constfilt,Butterworth,1,zoh,6.405,156.1,1.00000000
+constfilt,Butterworth,2,zoh,10.823,92.4,1.00000000
+constfilt,Butterworth,4,zoh,18.407,54.3,1.00000000
+constfilt,Butterworth,6,zoh,25.916,38.6,1.00000000
+constfilt,Butterworth,8,zoh,34.118,29.3,1.00000000
+constfilt,Butterworth,1,matchedz,6.404,156.1,1.00000000
+constfilt,Butterworth,2,matchedz,10.821,92.4,1.00000000
+constfilt,Butterworth,4,matchedz,18.703,53.5,1.00000000
+constfilt,Butterworth,6,matchedz,26.001,38.5,1.00000000
+constfilt,Butterworth,8,matchedz,34.263,29.2,1.00000000
+constfilt,Elliptic,2,tustin,10.847,92.2,0.94406088
+constfilt,Elliptic,4,tustin,18.317,54.6,0.94406088
+constfilt,Elliptic,6,tustin,25.925,38.6,0.94406088
+constfilt,Elliptic,8,tustin,34.192,29.2,0.94406088
+constfilt,Elliptic,2,zoh,10.809,92.5,0.94406088
+constfilt,Elliptic,4,zoh,18.345,54.5,0.94406088
+constfilt,Elliptic,6,zoh,25.902,38.6,0.94406088
+constfilt,Elliptic,8,zoh,34.761,28.8,0.94406088
+constfilt,Elliptic,2,matchedz,11.046,90.5,0.94406088
+constfilt,Elliptic,4,matchedz,18.812,53.2,0.94406088
+constfilt,Elliptic,6,matchedz,25.954,38.5,0.94406088
+constfilt,Elliptic,8,matchedz,34.177,29.3,0.94406088
+iir1,Butterworth,2,runtime,10.599,94.4,1.00000000
+iir1,Butterworth,4,runtime,18.020,55.5,1.00000000
+iir1,Butterworth,6,runtime,28.847,34.7,1.00000000
+iir1,Butterworth,8,runtime,39.983,25.0,1.00000000
+kfr,Butterworth,2,runtime+simd,602.752,1.7,1.00000000
+kfr,Butterworth,4,runtime+simd,386.308,2.6,1.00000000
+kfr,Butterworth,6,runtime+simd,307.883,3.2,1.00000000
+kfr,Butterworth,8,runtime+simd,307.968,3.2,1.00000000
+kfr,Elliptic,2,runtime+simd,602.510,1.7,0.94406088
+kfr,Elliptic,4,runtime+simd,386.593,2.6,0.94406088
+kfr,Elliptic,6,runtime+simd,308.161,3.2,0.94406088
+kfr,Elliptic,8,runtime+simd,308.193,3.2,0.94406088

include/constfilt/analog_filter.hpp

@@ -37,6 +37,13 @@ namespace constfilt
 //   AnalogFilter(continuous_tf, sample_rate_hz, method_tag)
 //     method_tag     - method instance; for TustinPW supply
 //                      constfilt::prewarp(warp_hz)
+//
+//   AnalogFilter(continuous_tf, factored_tf, sample_rate_hz, method_tag)
+//     factored_tf    - the poles, zeros, and gain of continuous_tf in
+//                      factored form. Required alongside continuous_tf.
+//                      For ZOH, the poles are used to build the Vandermonde
+//                      eigenvector matrix analytically. For MatchedZ, the
+//                      poles and zeros are mapped directly to the z-domain.
 template <typename T, consteig::Size N, typename Method = TustinNW,
           bool CheckStab = true>
 class AnalogFilter : public Filter<T, N + 1u, N + 1u>
@@ -79,6 +86,15 @@ class AnalogFilter : public Filter<T, N + 1u, N + 1u>
     {
     }
 
+    constexpr AnalogFilter(TransferFunction<T, N + 1u, N + 1u> continuous_tf,
+                           const FactoredTF<T, N> &factored_tf,
+                           T sample_rate_hz, BoundMethod method_tag)
+        : AnalogFilter(checked_discretize_factored(continuous_tf.b,
+                                                   continuous_tf.a, factored_tf,
+                                                   sample_rate_hz, method_tag))
+    {
+    }
+
   private:
     constexpr explicit AnalogFilter(
         TransferFunction<T, N + 1u, N + 1u> digital_tf)
@@ -120,6 +136,16 @@ class AnalogFilter : public Filter<T, N + 1u, N + 1u>
         return analog_to_digital<T, N>(
             b_c, a_c, static_cast<T>(1) / sample_rate_hz, method_tag);
     }
+
+    static constexpr TransferFunction<T, N + 1u, N + 1u>
+    checked_discretize_factored(const T (&b_c)[N + 1u], const T (&a_c)[N + 1u],
+                                const FactoredTF<T, N> &factored_tf,
+                                T sample_rate_hz, BoundMethod method_tag)
+    {
+        return discretize_with_factored<T, N>(
+            b_c, a_c, factored_tf, static_cast<T>(1) / sample_rate_hz,
+            method_tag);
+    }
 };
 
 } // namespace constfilt