What algorithm pagemldivide uses?

29 views (last 30 days)
Bruno Luong
Bruno Luong on 18 Jul 2022
Commented: Paul about 6 hours ago
The question starts here from Paul remark that pagemldivide does not return the same numerical solution than backslash "\" applied on each page.
Further test seem to show that happens only for size(A) < 10 in case of square matrices.
nmax = 100;
errx = zeros(1,nmax);
for n=1:nmax
A = rand(n,n,10);
B = rand(n,n,10);
X = pagemldivide(A,B);
errx(n) = norm(X(:,:,end)-A(:,:,end)\B(:,:,end));
end
% Check if solutions match backskash for n >= 10
all(errx(n>=10)==0)
ans = logical
1
plot(1:nmax, errx);
xlabel('Matrix size');
ylabel('Error between pagemldivide and \\');
The question is: what algorithm pagemldivide uses (for n < 10) ?
Note that it seems not to be LU-based either (test not showed).
The document implies the output of pagemldivide should match backslash.
  21 Comments
Paul
Paul about 5 hours ago
Rerunning example code in 2022b
maxNumCompThreads(1);
rng(100)
nmax = 100;
[errmld,errmdrd,errmt,errinv,errsvd] = deal(zeros(1,nmax));
for n=1:nmax
A = rand(n,n,10);
B = rand(n,n,10);
X = pagemldivide(A,B);
errmld(n) = norm(X(:,:,end)-A(:,:,end)\B(:,:,end));
X = pagemrdivide(A,B);
errmrd(n) = norm(X(:,:,end)-A(:,:,end)/B(:,:,end));
X = pagemtimes(A,B);
errmt(n) = norm(X(:,:,end)-A(:,:,end)*B(:,:,end));
X = pageinv(A);
errinv(n) = norm(X(:,:,end)-inv(A(:,:,end)));
X = pagesvd(A);
errsvd(n) = norm(X(:,:,end)-svd(A(:,:,end)));
end
% Check if solutions match backskash for n >= 10
%all(errx(n>=10)==0)
figure;
subplot(5,1,1);plot(1:nmax,errmld),ylabel('mldivide')
subplot(5,1,2);plot(1:nmax,errmrd),ylabel('mrdivide')
subplot(5,1,3);plot(1:nmax,errmt),ylabel('mtimes')
subplot(5,1,4);plot(1:nmax,errinv),ylabel('inv')
subplot(5,1,5);plot(1:nmax,errsvd),ylabel('svd')
xlabel('Matrix size');
Will these new results for (page)ml(r)divide hold for all inputs or only certain types of inputs? Any other insights you can offer now that 2022b is released may be of interest.

Sign in to comment.

Answers (0)

Products


Release

R2022a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!