(Matlab) Странная потеря точности при назначении Complex-матрицы локальной переменной - PullRequest
0 голосов
/ 19 февраля 2019

Я получаю эту неожиданную (и, по-видимому, неприличную) потерю точности с плавающей запятой, в частности значений Complex, при сохранении выражения в локальной переменной тот же оператор выполняется правильно в области действия «сценария»

(в настоящее время тестирование на Linux Matlab-2018a, похоже на ошибку, еще не тестировалось на других выпусках)

В вычислении хранится выражение transpose(a)*conj(a), где a является комплекснымматрица.Храня эквивалентный conj(a'*a), работающий без проблем, я хочу понять проблему, хотя, я извлек это из большей кодовой базы коллеги, которая получала неожиданные результаты.

В этом примере яиспользуя ishermitian(), чтобы увидеть, произошло ли неожиданное поведение, поскольку эти выражения по определению дают эрмитовы матрицы, а совместимая семантика округления с плавающей точкой будет сохранять эрмитово даже при потере точности.То же самое не происходит, если матрица "настоящая", кстати.

У меня есть следующий Minimal-Reproducible-Example: 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 - ненастоящие (как и ожидалось)


Последствия:

  1. Важное примечание от @ Cris Luengo 'ответ' - код 'function' в основном JIT 'd в современных движках Matlab, в то время как код' script '- вообще нет, отсюда и разница.

  2. Таким образом, этоэто не какая-то область видимости или свойство переменной (я начал пытаться выделить переменную в области видимости скрипта и передать ее функции и т. д.)

  3. В значительной степени мне хочетсясделать числовой код подобным образом на языке, где результаты вычислений проявляют такие допущения в системе типов (например, указание - это эрмитова матрица и т. д. в зависимости от алгебраической структуры скаляров), и использовать эту информацию для выбораалгоритмы и т. д.

1 Ответ

0 голосов
/ 19 февраля 2019

Вы можете упростить свой тест до:

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.

Я бы посоветовал вам отправить это как ошибку в MathWorks .(Ссылка «сообщить об ошибке» - самая правая, прямо под синей полосой, вам нужна действующая лицензия, чтобы сообщить об ошибке таким образом.)

...