(Matlab)奇怪的精度损失,同时将复杂矩阵分配给局部变量
(Matlab) Strange precision loss, while assigning Complex-matrix to local variable
当我将表达式存储到局部变量时,相同的语句在 'script'范围
(目前正在 Matlab-2018a Linux 上测试,似乎是一个错误,尚未在其他版本上测试过)
计算正在存储表达式 transpose(a)*conj(a)
,其中 a
是一个复数矩阵。存储等效的 conj(a'*a)
可以正常工作,但我想了解这个问题,我已经从一位同事的更大的代码库中提取了它,该代码库得到了意想不到的结果。
在示例中,我使用 ishermitian()
来查看是否发生了意外行为,因为根据定义这些表达式给出了 Hermitian 矩阵,并且兼容的浮点舍入语义将使它保持 hermitian even 失去精度时。如果矩阵是 "real" 顺便说一句,同样的事情不会发生。
我有以下最小可复制示例:tst_bad()
说明了行为不端的版本
len = 16; % generate complex input:
t = 2*pi/len*(0:len-1)';
tt = t + (0:0.1:0.3);
a = hilbert(sin(tt));
% a = rand(len, 4)+1i*rand(len, 4); % alternative input, almost as good
m1 = transpose(a)*conj(a);
res = -ones(1, 4);
res(1) = ishermitian(m1);
[res(2) m2] = tst_bad(a);
[res(3) m3] = tst_good(a);
m4 = single(m1);
res(4) = ishermitian(m4)
display(res)
% 1 0 1 1
function [f, m] = tst_bad(a)
m = transpose(a)*conj(a);
f = ishermitian(m);
m_iseq = isequal(m, transpose(a)*conj(a)) % 'true' but :
% running "isequal(m, transpose(a)*conj(a))" in the REPL will return
% 'false' if running in matlab-debugger
end
function [f, m] = tst_good(a)
m = conj(a'*a);
f = ishermitian(m);
end
有没有人见过类似的行为?请注意,即使是病态 m2
矩阵的对角线成员也不是真实的(如预期的那样)
后果:
来自@Cris Luengo
回复的重要说明 - 'function' 代码在当前的 Matlab 引擎中大部分是 JIT
,而 'script' 代码- 通常不是,因此有所不同。
因此,这不是任何变量范围或 属性(我已经开始尝试在 'script-scope' 中分配一个变量并将其传递给函数等)
几乎让我想用一种语言来做这样的数字代码,在这种语言中,计算结果将这种假设体现在类型系统中(比如声明——这是一个厄米矩阵,等等,这取决于关于标量的代数结构),并使用此信息进一步选择算法等
您可以将测试简化为:
len = 16;
t = 2*pi/len*(0:len-1)';
tt = t + (0:0.1:0.3);
a = hilbert(sin(tt));
m1 = transpose(a)*conj(a);
disp(ishermitian(m1)) % true
m2 = tst_bad(a);
disp(ishermitian(m2)) % false
disp(isequal(m1,m2)) % false
function m = tst_bad(a)
m = transpose(a)*conj(a);
end
我刚刚在 MATLAB R2018b(在线版)上 运行 确认了这个问题。在函数内或函数外计算 transpose(a)*conj(a)
会导致不同的结果。这看起来像是 JIT 的问题。
我会建议你 submit this as a bug to the MathWorks。 ("report a bug" link 是蓝色条下方最右边的一个,您需要有效的许可证才能以这种方式报告错误。)
当我将表达式存储到局部变量时,相同的语句在 'script'范围
(目前正在 Matlab-2018a Linux 上测试,似乎是一个错误,尚未在其他版本上测试过)
计算正在存储表达式 transpose(a)*conj(a)
,其中 a
是一个复数矩阵。存储等效的 conj(a'*a)
可以正常工作,但我想了解这个问题,我已经从一位同事的更大的代码库中提取了它,该代码库得到了意想不到的结果。
在示例中,我使用 ishermitian()
来查看是否发生了意外行为,因为根据定义这些表达式给出了 Hermitian 矩阵,并且兼容的浮点舍入语义将使它保持 hermitian even 失去精度时。如果矩阵是 "real" 顺便说一句,同样的事情不会发生。
我有以下最小可复制示例:tst_bad()
说明了行为不端的版本
len = 16; % generate complex input:
t = 2*pi/len*(0:len-1)';
tt = t + (0:0.1:0.3);
a = hilbert(sin(tt));
% a = rand(len, 4)+1i*rand(len, 4); % alternative input, almost as good
m1 = transpose(a)*conj(a);
res = -ones(1, 4);
res(1) = ishermitian(m1);
[res(2) m2] = tst_bad(a);
[res(3) m3] = tst_good(a);
m4 = single(m1);
res(4) = ishermitian(m4)
display(res)
% 1 0 1 1
function [f, m] = tst_bad(a)
m = transpose(a)*conj(a);
f = ishermitian(m);
m_iseq = isequal(m, transpose(a)*conj(a)) % 'true' but :
% running "isequal(m, transpose(a)*conj(a))" in the REPL will return
% 'false' if running in matlab-debugger
end
function [f, m] = tst_good(a)
m = conj(a'*a);
f = ishermitian(m);
end
有没有人见过类似的行为?请注意,即使是病态 m2
矩阵的对角线成员也不是真实的(如预期的那样)
后果:
来自@
Cris Luengo
回复的重要说明 - 'function' 代码在当前的 Matlab 引擎中大部分是JIT
,而 'script' 代码- 通常不是,因此有所不同。因此,这不是任何变量范围或 属性(我已经开始尝试在 'script-scope' 中分配一个变量并将其传递给函数等)
几乎让我想用一种语言来做这样的数字代码,在这种语言中,计算结果将这种假设体现在类型系统中(比如声明——这是一个厄米矩阵,等等,这取决于关于标量的代数结构),并使用此信息进一步选择算法等
您可以将测试简化为:
len = 16;
t = 2*pi/len*(0:len-1)';
tt = t + (0:0.1:0.3);
a = hilbert(sin(tt));
m1 = transpose(a)*conj(a);
disp(ishermitian(m1)) % true
m2 = tst_bad(a);
disp(ishermitian(m2)) % false
disp(isequal(m1,m2)) % false
function m = tst_bad(a)
m = transpose(a)*conj(a);
end
我刚刚在 MATLAB R2018b(在线版)上 运行 确认了这个问题。在函数内或函数外计算 transpose(a)*conj(a)
会导致不同的结果。这看起来像是 JIT 的问题。
我会建议你 submit this as a bug to the MathWorks。 ("report a bug" link 是蓝色条下方最右边的一个,您需要有效的许可证才能以这种方式报告错误。)