What algorithm pagemldivide uses?

10 visualizzazioni (ultimi 30 giorni)
Bruno Luong
Bruno Luong il 18 Lug 2022
Commentato: Paul il 27 Set 2022
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.
  23 Commenti
Heiko Weichelt
Heiko Weichelt il 27 Set 2022
@Paul, now that R2022b is released, please compare the updated doc pages for MLDIVIDE and PAGEMLDIVIDE, especially the Algorithms section (which is new for PAGEMLDIVIDE):
In R2022b, we tried to make both flow-graphs, and hence the used methods for the different inputs, as equal as possible. This yields that the above example now shows the same results.
However, you will notice that MLDIVIDE has a few more special solvers that PAGEMLDIVIDE does not support, i.e., 'permuted triangular' and 'upper Hessenberg'.
Both of these cases are rarely to find in created data but are rather the result of previos calls. E.g. 2-output LU returns a permuted triangular L factor in general. Another example is HESS or SCHUR, that creates upper Hessenberg form matrices as outputs.
Since we do not have page-wise solvers for those methods, we decided to not add these cases to PAGEMLDIVIDE.
Overall, as I pointed out earlier, even if your inputs will pick up the same solver, you can, in general, not assume that we will reach exactly the same results as the threading strategy might change between page-wise threading and threading over all pages. These settings will create differences similar to comparing results using different number of threads.
Since you uses 'maxNumCompThreads(1)' in the example above, you most likely will see the same results unless your inputs hit any of these special matrix forms. However, this is not guranteed as underlying algorithms might change.
As general tip for all page functions I'd say that results obtained using pagemldivide (or any other page function) are numerically equivalent to a for-loop calculating over the same matrices. However, the two results might differ due to floating-point round-off error.
I hope this helps.
Also, check out the performance release notes for both methods:
Paul
Paul il 27 Set 2022
Thank you very much for your responses in this thread.

Accedi per commentare.

Risposte (0)

Categorie

Scopri di più su Arithmetic Operations in Help Center e File Exchange

Prodotti


Release

R2022a

Community Treasure Hunt

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

Start Hunting!

Translated by