diff --git a/Cargo.toml b/Cargo.toml index 19a1297..e9ca8d5 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "control_systems_torbox" -version = "0.2.0" +version = "0.2.1" edition = "2024" authors = ["Tor Børve Rasmussen "] description = "Control systems toolbox" diff --git a/build.rs b/build.rs index 7df9610..0490754 100644 --- a/build.rs +++ b/build.rs @@ -13,7 +13,7 @@ fn main() { .to_string(); let make_args: Vec = vec![ "-j8".to_string(), - format!("BUILD_DIR={}", slicot_build_dir).to_string(), + format!("BUILD_DIR={slicot_build_dir}").to_string(), ]; let status = make() @@ -32,7 +32,7 @@ fn main() { panic!("Failed to build SLICOT") } - println!("cargo:rustc-link-search={}", slicot_build_dir); + println!("cargo:rustc-link-search={slicot_build_dir}"); println!("cargo:rustc-link-lib=static=slicot"); println!("cargo:rustc-link-lib=static=lpkaux"); println!("cargo:rustc-link-lib=gfortran"); diff --git a/src/analysis/frequency_response.rs b/src/analysis/frequency_response.rs index 48ee084..d765123 100644 --- a/src/analysis/frequency_response.rs +++ b/src/analysis/frequency_response.rs @@ -254,8 +254,7 @@ fn freq_response_ss_mat( if info != 0 { return Err(format!( - "Failed to compute frequency response of state space system. Slicot TB05AD failed with error code: {}", - info + "Failed to compute frequency response of state space system. Slicot TB05AD failed with error code: {info}" )); } diff --git a/src/analysis/system_properties.rs b/src/analysis/system_properties.rs index 3f2af17..185ae5a 100644 --- a/src/analysis/system_properties.rs +++ b/src/analysis/system_properties.rs @@ -50,6 +50,17 @@ fn h2_norm( let m = b.ncols(); let p = c.nrows(); + assert!(m > 0); + assert!(p > 0); + + if n == 0 { + let d_all_zero = d.iter().all(|&x| x == 0.0); + if d_all_zero { + return Ok(0.); + } + return Err("Order of system is zero and D matrix is not zero".into()); + } + let lda = a.nrows(); let ldb = b.nrows(); let ldc = c.nrows(); @@ -66,6 +77,7 @@ fn h2_norm( let mut iwarn = -1 as c_int; let mut info = -1 as c_int; + // println!("lda: {}", lda); let h2_norm = unsafe { ab13bd_( time_domain.as_ptr(), @@ -92,13 +104,11 @@ fn h2_norm( if info != 0 { return Err(Box::new(std::io::Error::other(format!( - "SLICOT ab13bd failed with the follwing error indicator: {}", - info + "SLICOT ab13bd failed with the follwing error indicator: {info}" )))); } else if iwarn != 0 { return Err(Box::new(std::io::Error::other(format!( - "SLICOT ab13bd returned with warning about numerical stability. Iwarn = {}", - iwarn + "SLICOT ab13bd returned with warning about numerical stability. Iwarn = {iwarn}" )))); } @@ -122,6 +132,21 @@ fn hinf_norm( ) -> Result { Ss::::verify_dimensions(&a, &b, &c, &d).map_err(|e| e.to_string())?; + let n = a.nrows(); + let m = b.ncols(); + let p = c.nrows(); + + assert!(p > 0); + assert!(m > 0); + + if n == 0 { + // The gain is constant and equal to D + // H-inf norm is the largest singular value + let singular_values = d.singular_values(); + let max_singular_value = singular_values.max(); + return Ok(max_singular_value); + } + let time_domain = if std::any::TypeId::of::() == std::any::TypeId::of::() { @@ -137,10 +162,6 @@ fn hinf_norm( let equilibration = CString::new("S").unwrap(); // Perform scaling let d_is_nonzero = CString::new("D").unwrap(); - let n = a.nrows(); - let m = b.ncols(); - let p = c.nrows(); - let mut e = DMatrix::identity(n, n); let mut f_peak = [0.0, 1.0]; // initial guess not active @@ -199,8 +220,7 @@ fn hinf_norm( if info != 0 { return Err(format!( - "SLICOT ab13dd returned with error code. info = {}", - info + "SLICOT ab13dd returned with error code. info = {info}" )); } @@ -293,8 +313,7 @@ pub fn zeros( if info != 0 { return Err(format!( - "SLICOT failed to find zeros of state-space system. Info = {}", - info + "SLICOT failed to find zeros of state-space system. Info = {info}" )); } let ldwork = dwork[0] as usize; @@ -336,8 +355,7 @@ pub fn zeros( if info != 0 { return Err(format!( - "SLICOT failed to find zeros of state-space system. Info = {}", - info + "SLICOT failed to find zeros of state-space system. Info = {info}" )); } if num_inv_zeros == 0 { @@ -384,7 +402,7 @@ pub fn zeros( ); } if info != 0 { - return Err(format!("LAPACK DGGEV returned info = {}", info)); + return Err(format!("LAPACK DGGEV returned info = {info}")); } lwork = dwork[0] as c_int; assert!(lwork > 0); @@ -413,7 +431,7 @@ pub fn zeros( ); } if info != 0 { - return Err(format!("LAPACK DGGEV returned info = {}", info)); + return Err(format!("LAPACK DGGEV returned info = {info}")); } let mut zeros = Vec::with_capacity(n_gen_eigen); @@ -446,6 +464,7 @@ impl Ss { /// - `Ok(f64)`: The H2 norm of the system. /// - `Err(String)`: An error message if the computation fails. pub fn norm_h2(&self) -> Result { + // println!("before h2: {}", self); h2_norm::( self.a().clone(), self.b().clone(), @@ -510,6 +529,7 @@ impl Ss { #[cfg(test)] mod tests { use crate::{ + Continuous, Ss, systems::Tf, transformations::SsRealization::{ControllableCF, ObservableCF}, }; @@ -546,6 +566,13 @@ mod tests { epsilon = 1e-2 ); } + + let sys = Ss::::new_from_scalar(-2.0); + assert_abs_diff_eq!(sys.norm_hinf().unwrap(), 2.0); + + assert!(sys.norm_h2().is_err()); + let sys = Ss::::new_from_scalar(0.0); + assert_eq!(sys.norm_h2().unwrap(), 0.0); } #[test] @@ -592,9 +619,9 @@ mod tests { #[test] fn ss_zeros() { let tf = (Tf::s() - 1.0) * (Tf::s() + 4.0) / (Tf::s() + 2.0).powi(2); - println!("tf: \n{}", tf); + println!("tf: \n{tf}"); let zeros = tf.to_ss_method(ControllableCF).unwrap().zeros().unwrap(); - println!("zeros: {:?}", zeros); + println!("zeros: {zeros:?}"); assert_eq!(zeros.len(), 2); let tf = 1.0 / Tf::s(); @@ -643,18 +670,18 @@ mod tests { #[test] fn ss_is_stable() { let tf = 1.0 / Tf::s(); - assert_eq!(tf.to_ss_method(ObservableCF).unwrap().is_stable(), false); + assert!(!tf.to_ss_method(ObservableCF).unwrap().is_stable()); let tf = (Tf::s() + 1.0) / (Tf::s() + 1.0).powi(4); let ss = tf.to_ss_method(ObservableCF).unwrap(); - assert_eq!(ss.is_stable(), true); + assert!(ss.is_stable()); let tf = 1.0 / ((Tf::s() + 1.0) * (Tf::s() - 2.0)); let ss = tf.to_ss_method(ObservableCF).unwrap(); - assert_eq!(ss.is_stable(), false); + assert!(!ss.is_stable()); let tf = 1.0 / (Tf::s().powi(2) + 2.0 * 0.01 * Tf::s() + 1.0); let ss = tf.to_ss_method(ObservableCF).unwrap(); - assert_eq!(ss.is_stable(), true); + assert!(ss.is_stable()); } } diff --git a/src/systems/polynom.rs b/src/systems/polynom.rs index 838136c..7e0e01c 100644 --- a/src/systems/polynom.rs +++ b/src/systems/polynom.rs @@ -60,8 +60,8 @@ where let term_i = match exp { 0 => val.to_string(), - 1 => format!("{}{}", val, var_name), - _ => format!("{}{}^{}", val, var_name, exp), + 1 => format!("{val}{var_name}"), + _ => format!("{val}{var_name}^{exp}"), }; terms.push(term_i); } @@ -76,7 +76,7 @@ where { fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { let str_poly = self.to_string_variable("x"); - write!(f, "{}", str_poly) + write!(f, "{str_poly}") } } @@ -353,11 +353,7 @@ where assert!(num_str.len() <= max_len); assert!(den_str.len() <= max_len); - write!( - formatter, - "\n {}\n {}\n {}\n\n", - num_str, dash_line, den_str - ) + write!(formatter, "\n {num_str}\n {dash_line}\n {den_str}\n\n") } } @@ -649,7 +645,7 @@ mod tests { lhs_vals .iter() .zip(rhs_vals.iter()) - .all(|(x, y)| f64::abs_diff_eq(&x, &y, epsilon)) + .all(|(x, y)| f64::abs_diff_eq(x, y, epsilon)) } } @@ -712,9 +708,9 @@ mod tests { } { - println!("poly: {}", p); + println!("poly: {p}"); let p_str = p.to_string(); - println!("poly to string: {}", p_str); + println!("poly to string: {p_str}"); } } @@ -760,7 +756,7 @@ mod tests { assert_abs_diff_eq!(rf.num.coeffs[0], 1.0); assert_abs_diff_eq!(rf.num.coeffs.last().unwrap().clone(), 3.0); - println!("Rational function normalize: \n{}", rf); + println!("Rational function normalize: \n{rf}"); } #[test] diff --git a/src/systems/state_space.rs b/src/systems/state_space.rs index 436372d..b4407c0 100644 --- a/src/systems/state_space.rs +++ b/src/systems/state_space.rs @@ -612,7 +612,7 @@ impl fmt::Display for Ss { } else { "Unknown" }; - write!(f, "{}-time state-space model\n\n", time_domain) + write!(f, "{time_domain}-time state-space model\n\n") } } @@ -635,7 +635,7 @@ mod tests { rng: &mut U, max_order: usize, ) -> Tf { - let den_order = rng.random_range(1..=max_order) as usize; + let den_order = rng.random_range(1..=max_order); let num_order = rng.random_range(0..=den_order); let num: Vec = (0..=num_order) @@ -884,14 +884,14 @@ mod tests { #[test] fn ss_display() { let sys = (1.0 / Tf::s()).to_ss().unwrap(); - println!("1: {}", sys); + println!("1: {sys}"); let sys = Ss::::new_from_scalar(2.0); - println!("2: {}", sys); + println!("2: {sys}"); let sys = ((1.0 + Tf::s()) / ((Tf::s() - 1.0) * Tf::s()).powi(2)) .to_ss() .unwrap(); - println!("3: {}after", sys); + println!("3: {sys}after"); } } diff --git a/src/systems/transfer_function.rs b/src/systems/transfer_function.rs index 081ad2b..3719338 100644 --- a/src/systems/transfer_function.rs +++ b/src/systems/transfer_function.rs @@ -497,10 +497,10 @@ mod tests { assert_abs_diff_eq!(y.re, y_expect.re); assert_abs_diff_eq!(y.im, y_expect.im); } - println!("Tf cont: \n{}", tf); + println!("Tf cont: \n{tf}"); let tf_z = Tf::::new_from_rf(tf.rf); - println!("Tf discrete: \n {}", tf_z); + println!("Tf discrete: \n {tf_z}"); } #[test] @@ -553,16 +553,16 @@ mod tests { Tf::::new(&[0.0, 0., 0., 1., 2., 0.0], &[0.0, 1.]); assert_eq!(sys.degree_num_den(), (4, 1)); assert_eq!(sys.relative_degree(), 4 - 1); - assert_eq!(sys.is_proper(), false); - assert_eq!(sys.is_strictly_proper(), false); + assert!(!sys.is_proper()); + assert!(!sys.is_strictly_proper()); let sys = Tf::s() / Tf::s(); - assert_eq!(sys.is_proper(), true); - assert_eq!(sys.is_strictly_proper(), false); + assert!(sys.is_proper()); + assert!(!sys.is_strictly_proper()); let sys = 1. / Tf::s(); - assert_eq!(sys.is_proper(), true); - assert_eq!(sys.is_strictly_proper(), true); + assert!(sys.is_proper()); + assert!(sys.is_strictly_proper()); } #[test] @@ -582,11 +582,11 @@ mod tests { #[test] fn tf_display() { let tf = 1.0 / Tf::s(); - println!("1: {}", tf); + println!("1: {tf}"); let tf = 1.0 / Tf::s().powi(5) * (Tf::s() + 1.0); - println!("2: {}", tf); + println!("2: {tf}"); let tf = 2.1 / Tf::s() * (Tf::s() - 1.2).powi(5) / (Tf::s() + 1.0).powi(3); - println!("3: {}", tf); + println!("3: {tf}"); } } diff --git a/src/transformations.rs b/src/transformations.rs index 66924ed..13890c6 100644 --- a/src/transformations.rs +++ b/src/transformations.rs @@ -119,8 +119,7 @@ fn ss2tf_mat( if info != 0 { return Err(Box::new(std::io::Error::other(format!( - "SLICOT tb04ad_ failed with info code {}", - info + "SLICOT tb04ad_ failed with info code {info}" )))); } @@ -155,6 +154,17 @@ impl Tf { let (_, den_deg) = tf.degree_num_den(); let nx = den_deg; + let mut d = DMatrix::zeros(1, 1); + if !tf.is_strictly_proper() { + assert!(tf.numerator().coeffs().len() > nx); + d[(0, 0)] = tf.numerator().coeffs()[nx]; + } + + if nx == 0 { + assert_eq!(tf.relative_degree(), 0); + return Ok(Ss::new_from_scalar(d[(0, 0)])); + } + let mut a = DMatrix::zeros(nx, nx); a.view_mut((1, 0), (nx - 1, nx - 1)) .copy_from(&DMatrix::identity(nx - 1, nx - 1)); @@ -162,12 +172,6 @@ impl Tf { a[(row, nx - 1)] = -tf.denominator().coeffs()[row]; } - let mut d = DMatrix::zeros(1, 1); - if !tf.is_strictly_proper() { - assert!(tf.numerator().coeffs().len() > nx); - d[(0, 0)] = tf.numerator().coeffs()[nx]; - } - let mut num_extended = tf.numerator().coeffs().to_vec(); num_extended.resize(nx, 0.); @@ -397,8 +401,7 @@ pub fn minreal_mat_mut( if info != 0 { return Err(format!( - "Minreal failed. SLICOT TB01PD returned info = {}", - info + "Minreal failed. SLICOT TB01PD returned info = {info}" )); } @@ -586,7 +589,7 @@ mod tests { rng: &mut U, max_order: usize, ) -> Tf { - let den_order = rng.random_range(1..=max_order) as usize; + let den_order = rng.random_range(1..=max_order); let num_order = rng.random_range(0..=den_order); let num: Vec = (0..=num_order) @@ -639,7 +642,7 @@ mod tests { let tf = ss.to_tf().unwrap(); let tf_ans = 1. / Tf::s().powi(2); - println!("ss2Tf: \n{}", tf); + println!("ss2Tf: \n{tf}"); assert_abs_diff_eq!(tf, tf_ans); } @@ -659,6 +662,14 @@ mod tests { assert_abs_diff_eq!(ss.d(), &d_ans); assert_abs_diff_eq!(tf_ret.normalize(), tf.normalize()); + + let tf = Tf::<_, Continuous>::new_from_scalar(1.0); + let ss = tf.to_ss().unwrap(); + assert_eq!(ss.order(), 0); + assert_eq!(ss.d()[(0, 0)], 1.0); + + assert!(!tf.is_strictly_proper()); + assert!(tf.is_proper()); } #[test] @@ -695,7 +706,7 @@ mod tests { let ss = tf.to_ss().unwrap(); println!("a: {}, b: {}, c: {}, d: {}", ss.a(), ss.b(), ss.c(), ss.d()); let tf_minreal = ss.minreal(None).unwrap().to_tf().unwrap(); - println!("tf minreal:\n{}", tf_minreal); + println!("tf minreal:\n{tf_minreal}"); assert_abs_diff_eq!(1.0 / Tf::s(), tf_minreal, epsilon = 1e-2); let tf = Tf::s().powi(5) / Tf::s().powi(6); diff --git a/src/utils/traits.rs b/src/utils/traits.rs index 411427b..a7c20f1 100644 --- a/src/utils/traits.rs +++ b/src/utils/traits.rs @@ -196,11 +196,11 @@ mod tests { #[test] fn mag_db_covert() { - assert_abs_diff_eq!((1. as f64).mag2db(), 0.); - assert_abs_diff_eq!((1. as f32).mag2db(), 0.); + assert_abs_diff_eq!(1_f64.mag2db(), 0.); + assert_abs_diff_eq!(1_f32.mag2db(), 0.); - assert_abs_diff_eq!((20 as f64).db2mag(), 10.); - assert_abs_diff_eq!((20 as f32).db2mag(), 10.); + assert_abs_diff_eq!(20_f64.db2mag(), 10.); + assert_abs_diff_eq!(20_f32.db2mag(), 10.); let mut rng = rand::rng();