rust - 有没有办法写入代数矩阵的整行/整列？

Question

我正在使用 DMatrix 结构来分配动态大小的矩阵，我在其中使用 L2 范数通过归一化列向量重复覆盖每一列。

// a is some DMatrix of arbitrary size
let col_0 = a.column(0);
let norm_of_col_0 = col_0.normalize();

而不是遍历当前列中的每个单元格:

let row = a.shape().0;
let col = a.shape().1;
for col in 0..ncols {
    let norm_of_col = a.column(col).normalize();
    for row in 0..nrows {
        *a.index_mut((row, col)) = norm_of_col()[row];
    }
}

我想直接用标准化版本覆盖该列。代码在语义上应该是这样的:

*a.index_mut((_, col)) = norm_of_col();

(_, col) 表示我选择了 col 列，_ 表示整行。

更一般地说，有没有办法用相同大小和数据类型的新行或列覆盖行或列？ insert_columns 等方法仅向现有矩阵添加列。

如果是这样，这样做在计算上会更快，还是我应该只编写一个辅助方法来遍历每个单元格以更新矩阵？

Answer 1

你可以用 nalgebra 0.18.0 这样做:

use nalgebra::DMatrix;

fn main() {
    let mut m = DMatrix::from_vec(2, 3, (0 .. 6).map(|n| n as f64).collect());
    dbg!(&m);
    for mut col in m.column_iter_mut() {
        let normalized = col.normalize();
        col.copy_from(&normalized);
    }
    dbg!(&m);
}

与您的代码相比，我没有衡量此代码的性能。

请注意，copy_from 遍历项目而不在每一步检查边界，而是在循环之前只检查一次。我还没有检查优化器是否可以在您的代码中执行等效转换。这个简单的基准测试为我机器上这个答案的解决方案提供了优势(不确定它的代表性如何；通常的基准免责声明适用):

use criterion::{black_box, criterion_group, criterion_main, Benchmark, Criterion};
use nalgebra::DMatrix;

fn normalize_lib(m: &mut DMatrix<f64>) {
    for mut col in m.column_iter_mut() {
        let normalized = col.normalize();
        col.copy_from(&normalized);
    }
}

fn normalize_hand_rolled(a: &mut DMatrix<f64>) {
    let nrows = a.shape().0;
    let ncols = a.shape().1;
    for col in 0..ncols {
        let norm_of_col = a.column(col).normalize();
        for row in 0..nrows {
            *a.index_mut((row, col)) = norm_of_col[row];
        }
    }
}

fn benchmark(c: &mut Criterion) {
    let mut m0 = DMatrix::new_random(100, 100);
    let mut m1 = m0.clone();
    let bench = Benchmark::new("lib", move |b| b.iter(|| normalize_lib(black_box(&mut m0))))
        .with_function("hand_rolled", move |b| {
            b.iter(|| normalize_hand_rolled(black_box(&mut m1)))
        });
    c.bench("normalize", bench);
}

criterion_group!(benches, benchmark);
criterion_main!(benches);

normalize/lib           time:   [26.102 us 26.245 us 26.443 us]
normalize/hand_rolled   time:   [37.013 us 37.057 us 37.106 us]