在 Matlab 中计算矩阵幂（矩阵多项式）的加权求和

Question

给定一个 nxn 矩阵 A_k 和一个 nx1 向量 x，是否有任何聪明的方法来计算

使用 Matlab？ x_i是向量x的元素，因此J是矩阵之和。到目前为止，我使用了 for 循环，但我想知道是否有更聪明的方法。

Answer 1

简答：您可以使用内置的 matlab 函数 polyvalm 进行矩阵多项式计算，如下所示：

x = x(end:-1:1); % flip the order of the elements
x(end+1) = 0; % append 0
J = polyvalm(x, A);

长答案：Matlab 在内部使用循环。因此，如果您优化自己的实现（请参阅我的 calcJ_loopOptimised 函数），您并没有获得那么多，或者表现更差：

% construct random input
n = 100;
A = rand(n);
x = rand(n, 1);

% calculate the result using different methods
Jbuiltin = calcJ_builtin(A, x);
Jloop = calcJ_loop(A, x);
JloopOptimised = calcJ_loopOptimised(A, x);

% check if the functions are mathematically equivalent (should be in the order of `eps`)
relativeError1 = max(max(abs(Jbuiltin - Jloop)))/max(max(Jbuiltin))
relativeError2 = max(max(abs(Jloop - JloopOptimised)))/max(max(Jloop))

% measure the execution time
t_loopOptimised = timeit(@() calcJ_loopOptimised(A, x))
t_builtin = timeit(@() calcJ_builtin(A, x))
t_loop = timeit(@() calcJ_loop(A, x))

% check if builtin function is faster
builtinFaster = t_builtin < t_loopOptimised

% calculate J using Matlab builtin function
function J = calcJ_builtin(A, x)
  x = x(end:-1:1);
  x(end+1) = 0;
  J = polyvalm(x, A);
end

% naive loop implementation
function J = calcJ_loop(A, x)
  n = size(A, 1);
  J = zeros(n,n);
  for i=1:n
    J = J + A^i * x(i);
  end
end

% optimised loop implementation (cache result of matrix power)
function J = calcJ_loopOptimised(A, x)
  n = size(A, 1);
  J = zeros(n,n);
  A_ = eye(n);
  for i=1:n
    A_ = A_*A;
    J = J + A_ * x(i);
  end
end

对于 n=100，我得到以下信息：

t_loopOptimised = 0.0077
t_builtin       = 0.0084
t_loop          = 0.0295

对于 n=5，我得到以下信息：

t_loopOptimised = 7.4425e-06
t_builtin       = 4.7399e-05
t_loop          = 1.0496e-04

请注意，我的时间在不同的运行之间略有波动，但优化的循环几乎总是比内置函数快（小 n 高达 6 倍）。

在 Matlab 中计算矩阵幂（矩阵多项式）的加权求和

Compute weighted summation of matrix power (matrix polynomial) in Matlab

matlab

matrix-multiplication

polynomials