What kinds of operation can be speeded up by reducing significant digits using VPA/Digits

I am a new beginer of Matlab and learning about variable-precission arithmetic by using VPA & Digits.

I've found an example from "doc digits" about "Increase Speed by Reducing Precision", which is about a massive computation of Zeta function (that I don't know so far).

I think I can try someting much easy to understand, so I wrote:

digits 4;
a=vpa(0.123456789012345);
b=vpa(0.543210987654321);
tic
for i=1:10000000000
c=a+b;
end
disp(c)
toc

The result is far from my imagination, if there are no VPA in use (that will be a=0.123456......), it took just 10+ seconds to show up c & toc, but if it was running just as above, it maks me press CTRL-C to interrupt after a few minutes of waiting.

So my question is as the title of this post, any ideas? Thanks!!

[勇忠 黄]

Thank you John for your kind explanation!!

Answer 1

[John D'Errico]

I think you misunderstand. vpa is not just something that you can use to dial in your own choice of precision as easily as you want. MATLAB normally works using double precision arithmetic. All normal numbers are stored as a double. Single is an alternative, but it gains you mostly a decrease in memory, at a significant cost in precision. There is some gain in speeed using single as I recall, at least on some computations. Usully though, you just stick with double.

All MATLAB routines are written to run as fast as you can hope for, when using the standard arithmetic, when applied to doubles or singles.

However, at a significant, serious cost in speed and memory, you can use the symbolic toolbox to represent numbers. It can do higher precision arithmetic. It is not fast. And in fact, you cannot even use it in all possible computations, as some functionality potentially may not work. In fact, compared to the same computations using doubles, computations using these symbolic numbers will be comparatively as slow as mollasses on a cold day.

In some relatively rare circumstances, you can gain by the use of tools like VPA. (Or, for that matter, my own high precision tool HPF.) If you have an isolated computation where there is a problem with the dynamic range of a double preciosn number, you might decide to do that computation using the higher precision capability, and then immediately convert back into doubles as soon as possible.

What you do NOT want to do is try to work fully in your own chosen precision. That way lies death, as you have found. A slow, painful death.

Instead, this is why classes are taught. You will find courses in numerical methods, in numerical analysis, in numerical linear algebra. There you learn about how to do computations where you do not need to resort to high precision computations. You can learn the limitations of the various numerical methods, as well as how to solve most problems while never needing to resort to high precision.

In fact, even as the author of several tools in MATLAB that allow you to use high/variable precision, I would claim these tools can become more of a crutch than they are an aid. By allowing you to not need to worry about how to do your work in the precision available, they convince you this is a good idea. You are far better off understanding how to do computations using a double, and only as a last resort should you even onsider verging into higher precision. As I said, this can come from taking the proper courses, reading a few books, etc. It can come from learning to code from someone who is willing and able to mentor you. But high precision should be the last resort, not a general tool.

Having said all of that, are there cases where I would use high precision in my own work? Occasionally. For example, I use large integer computations in my search for huge Hampshire numbers that are prime. These searches usually wander into numbers with a hundred thousand decimal digits or more. In other circumstances, I've used these tools to compute numbers like e and pi to as much as a million digits. (Really, that is just because I can, as an exercise in numerical methods.) For almost all serious work, doubles work just fine for me.

[勇忠 黄]

Dear MVP John, Thank you very much for your answer! I think I understand what you've said.

Now here below is the example I mentioned in Matlab's document:

[X,Y] = meshgrid((0:0.0025:.75),(5:-0.05:0));
C = X + Y*i;
digits(10)
vpaC = vpa(C);
tic
zeta(vpaC);%or zeta(C), for comparison
toc

if we call the function Zeta by zeta(C) where C is double class, it will be about 4x slower then by zeta(vpaC), this example shows that vpa (or symbolic) data really can improve performance of certain calculation/operation. My real problem is "what kinds of calculation are there, for example, in zeta function, that can be speeded up by using "low-significant-digits" symbolic variables.".

Thanks again!

[John D'Errico]

And my answer is, there are almost no computations where you will truly gain in speed. The gain you might see here is totally illusory.

If you look at what happens when you use zeta on a double precision input, it still uses the same codes, but one of the things it does is convert the number to a symbolic one, using VPA.

So you have extra overhead converting numbers from doubles into a symbolic form, using 16 digits. And then you have extra overhead converting back. You can look at the code, if you don't believe me. I just did. For example, you will see code like this, in the call zeta(A), when A is a double.

 oldDigits = digits(16);
 args = cellfun(@(x)vpa(x), varargin, 'UniformOutput', false);
 cleanupObj = onCleanup(@() digits(oldDigits));
 res = cast(fn(args{:}),superiorfloat(varargin{:}));

In there, you will see the call to vpa itself, applied to x. The last line of that code, if you look carefully, has MATLAB applying the function "fn" to a cell array input. args is a cell array, containing a list of real numbers, but converted to vpa form. The function inside fn is @zeta. Then at the very end, it casts that back into a double precision output.

The reason zeta of a double is slower is NOT because vpa is faster, but because MATLAB used a lot of that time just doing class conversions, etc. It is not because VPA was inherantly faster for that, because it did all the work on vpa'ed numbers anyway. The code for zeta in MATLAB is purely one that works on symbolic inputs.

I'm sorry, but there really is no good reason I can think of why you would logically use vpa for a gain in speed when a double would not have had a problem because of the dynamic range of the double. You typically use vpa because the double will overflow to inf/-inf or underflow to 0, or because there are insufficient digits for the precision you need.

[勇忠 黄]

Got that! You DID help me a lot!! I was misled by the information provided by Matlab (screenshot below):

It makes me think of that VPA can be used to speed up the performance by lowering precssion of digits.

Now I see. :-)

Thank you very much!!

[John D'Errico]

If you change the precision using vpa, thus vpa with 10 versus 100 digits and you will see a difference.

But vpa versus double will be a rough comparison, since symbolic computations in ANY form are vastly slower than are similar computations using doubles.

If the computation can be done using doubles, then use doubles. As long as you are not in a realm where the double has numerical problems, you will gain speed there.

If the code uses symbolic tools, then odds are, the code will convert it to a symbolic form (i.e., vpa) anyway, and therefore you are not gaining speed.

A = vpa(randi(100,100,100));
D = double(A);

timeit(@() A*2)
ans =
            0.070772938132
            
timeit(@() D*2)
ans =
      5.50372533333333e-06

See that the time required for most arithmetic operations on a double array is VASTLY better than for a similar operation on a symbolic array. In the above case, we see what, a factor of more than 10000?

It gets worse than that. On many truly large computations, your computer can use more than one core. Mine has 8 real cores. So on some computations, MATLAB can intelligently choose to access all 8 cores. So I get roughly an 8x speedup for some computations. These mostly reduce to large scale linear algebraic computations. For example, if I try this:

A = randn(15000);
[L,U] = lu(A);

for a few seconds, MATLAB will have all 8 cores on my computer running flat out. However, even for large symbolic computations, MATLAB only uses one core. It does not farm out symbolic computations to multiple cores. (That may change one day. But this is true for now.) And calls that use vpa are symbolic toolbox calls. They do not gain from multiple cores at this point in time. Yes, I suppose you could use parfor, and gain some speed that way, IF you have the parallel computation toolbox.

So any computation that works using vpa, perhaps even only in 10 significant digits, will still never be faster than a double. There really is not much value in hoping to use vpa as a speedup, not when compared to a double.

digits 30
A = vpa(randn(100));
timeit(@() inv(A))
ans =
            4.266870108132
            
digits 20
A = vpa(randn(100));
timeit(@() inv(A))
ans =
            4.120985673132       
            
digits 15
A = vpa(randn(100));
timeit(@() inv(A))
ans =
            0.121588904132

D = double(A);
timeit(@() inv(D))
ans =
      0.000186797873666667

As you can see, decreasing the number of digits gains speed. The big gain when I went under 16 digits probably means the software is smart enough to convert to a double FIRST, then calls inv on the result as a double matrix, then goes back into symbolic form. So while the 15 digit inverse was relatively fast, it was still much slower than a double.

In any case, while using vpa for low precision is a speedup, it simply does not compare to anything you can do using a double. The only reason zeta was slow for a double really seems to be things like function call overhead, and several class conversions.

[勇忠 黄]

Much more clearer with the fact!

The reason of my learning about VPA is that I was asked to do a numerical test about the Wilkinson Iterative Refinement algorithm to solve linear systems Ax=b. One of the 3 iterative steps is said that the refinement element vector can be solved under low precission and will not affect the convergence of a higher precission result.

That's why I have to find ways of simulation of "low-precission-calculation" under my 64-bit Matlab running in a 4-core Win10 environment. Obviously, I don't concern about speed because I surely know any "simulation" things cost additional loads of CPU.

But finally I found VPA is not suitable for my work because of its "Guarding Digits" mechanism, while I need for example, vpa(0.33333,4) to be an exact value of 0.3333, and it returns me with this: vpa(0.33333,4)+vpa(0.33333,4)=0.6667, not 0.6666 that is what I really need.

So at that time I simply use single() & double() to do the simulation.