0

enter image description here

Matrix exponential is defined as: enter image description here

Why is that the computational error decreases until around k = 75, then it stays the same? Is it because the factorial is too big, and because Matlab only uses 16 digits of precision, the round-off error is neglected?

Here's my Matlab algorithm:

format long
%     load('CA3matrix.mat');
    A = randn(500);
    terms_arr = 5:5:150; % number of terms in the exp(A)
    prev_term = 1; % store the previous terms, avoid duplicate calculation
    expAk = A^0;
    errors = [];

    for i = terms_arr
       for k = 1:i
            term = (1/k) * A * prev_term;
            prev_term = term;
            expAk = expAk + term;
            %expAk = expAk + (1/factorial(k)) * A^(k);
       end
       expA = expm(A);
       err = norm(expA - expAk)/norm(expA);
       errors = [errors err];
       expAk = A^0;
       prev_term = 1;
    end

    semilogy(terms_arr, errors, '*')
    title('Computational error versus k')
    xlabel('k','fontsize',16)
    ylabel('Computational error','fontsize',16)
    set(gca,'FontSize',14)
  • What does this instruction do: $expAk = A^0$;? – NoChance Jun 07 '19 at 23:27
  • @NoChance to generate the Identity matrix. Because the algorithm has that – Joshua Leung Jun 07 '19 at 23:31
  • @Moo, we are expected to run the algorithm in Matlab. Round-off error is expected, but we need to explain the phenomenon. I can't figure out why. – Joshua Leung Jun 07 '19 at 23:40
  • Also expm(A) will differ from the exact result by some multiple of the machine precision. But since it uses a different algorithm, it will be a different error. So the difference will always contain some random error on the machine precision scale plus the truncation error. With an increasing number of terms, the accumulated floating point noise will also increase linearly, but that is not visible in your graph. – Lutz Lehmann Jun 08 '19 at 09:10

0 Answers0