lax/
opnorm.rs

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
//! Operator norm

use super::{AsPtr, NormType};
use crate::{layout::MatrixLayout, *};
use cauchy::*;

pub struct OperatorNormWork<T: Scalar> {
    pub ty: NormType,
    pub layout: MatrixLayout,
    pub work: Vec<MaybeUninit<T::Real>>,
}

pub trait OperatorNormWorkImpl {
    type Elem: Scalar;
    fn new(t: NormType, l: MatrixLayout) -> Self;
    fn calc(&mut self, a: &[Self::Elem]) -> <Self::Elem as Scalar>::Real;
}

macro_rules! impl_operator_norm {
    ($s:ty, $lange:path) => {
        impl OperatorNormWorkImpl for OperatorNormWork<$s> {
            type Elem = $s;

            fn new(ty: NormType, layout: MatrixLayout) -> Self {
                let m = layout.lda();
                let work = match (ty, layout) {
                    (NormType::Infinity, MatrixLayout::F { .. })
                    | (NormType::One, MatrixLayout::C { .. }) => vec_uninit(m as usize),
                    _ => Vec::new(),
                };
                OperatorNormWork { ty, layout, work }
            }

            fn calc(&mut self, a: &[Self::Elem]) -> <Self::Elem as Scalar>::Real {
                let m = self.layout.lda();
                let n = self.layout.len();
                let t = match self.layout {
                    MatrixLayout::F { .. } => self.ty,
                    MatrixLayout::C { .. } => self.ty.transpose(),
                };
                unsafe {
                    $lange(
                        t.as_ptr(),
                        &m,
                        &n,
                        AsPtr::as_ptr(a),
                        &m,
                        AsPtr::as_mut_ptr(&mut self.work),
                    )
                }
            }
        }
    };
}
impl_operator_norm!(c64, lapack_sys::zlange_);
impl_operator_norm!(c32, lapack_sys::clange_);
impl_operator_norm!(f64, lapack_sys::dlange_);
impl_operator_norm!(f32, lapack_sys::slange_);