MathWorks Blogs

The final MATLAB example of "Inplace modification" is not correct.

function f(out, x)
out = x.^2
end
x = rand(10)
y = zeros(length(x), 1)
f(y, x)

What happens here when you call f(y, x) is:

1. The arrays x and y are passed to f (no memory copying done yet, since MATLAB uses copy on write)

2. When we have `out = x.^2`, this will allocate a new array in memory and store in it the result of `x.^2`, and will call this `out`. The original `out` which was passed into the function can now be garbage collected (although it won't be, because it's still in the parent scope as 'y')

3. When the function exits, the new `out` goes out of scope and can be garbage collected.

So all this example does is allocate a new array, assign x.^2 to it, and then throw that result away again. There's no in-place modification.

You can't really pass a matrix by reference in MATLAB like you can in Python, however in some cases you can write your function in a specific way which will ensure in place operations are done and prevent a huge matrix copy. For example (pauses are just there so that you can watch task manager memory usage):

function x = inPlace()
x = rand(100000000, 1);
pause(3);
x = doSomething(x);
pause(3)
end

function y = doSomething(y)
y = y.^2;
end

Conditions required are:

* Input and output variables must have same name in caller

* Input and output variables must have same name in callee

* Callee must be inside a function,…

Cleve Moler was one of the big names in numerical methods, and participated in creation of canonical FORTRAN libraries for solving linear equations, and matrix algorithms more generally.

To teach this more conveniently to his students, he wrote the original version of MATrixLABoratory to allow interactive exploration of the library functions without having to compile FORTRAN code. The original version was about 2000 lines of code in FORTRAN.

Engineering students loved it so much that he decided to make a company around this product. His buddy expanded and rewrote the interpreter in C, for a PC, and the rest is history:

"In 1983 Jack Little suggested the creation of a commercial product based on MATLAB. I said I thought that was a good idea, but I didn't join him initially. The IBM PC had been introduced only two years earlier and was barely powerful enough to run something like MATLAB, but Little anticipated its evolution. He left his job, bought a Compaq PC clone at Sears, moved into the hills behind Stanford, and, with my encouragement, spent a year and a half creating a new and extended version of MATLAB written in C. A friend, Steve Bangert, joined the project and worked on the new MATLAB in his spare time."

User guide for the original version of MATLAB:
https://blogs.mathworks.com/cleve/2018/02/05/the-historic-matlab-users-guide/

The source code of the very early (1982?) FORTRAN version of MATLAB:
https://github.com/johnsonjh/matlab

The origins of the first PC…

The paper by Lyness and Moler is nice, but if it's really "about the idea in this post" then the fact is somewhat hidden. It's mostly about using the Cauchy integral formula to turn derivatives into integrals, using the Poisson summation formula to relate those to finite sums, and using the Moebius inversion formula to get those relations into a form from which you can extract the derivatives.

Maybe taking a very short partial sum of one of the series in Theorem 2 gives you the result in the blog post here, but that seems rather sledgehammer-to-crack-a-nut-ish when (as the blog post says) all you need is the Cauchy-Riemann equations.

I think it's true that the Lyness-Moler idea of "numerical differentiation by numerical complex integration" was, historically, part of the chain of ideas that led to the very simple "complex-step approximation" discussed in the blog post. But that's a far cry from saying that their paper is "about the idea in this post"!

There's some discussion of the history here: http://blogs.mathworks.com/cleve/2013/10/14/complex-step-differentiation/ where Moler (as in "Lyness and Moler") says "The complex step algorithm that I'm about to describe is much simpler than anything described in that 1967 paper. So, although we have some claim of precedence for the idea of using complex arithmetic on real functions, we certainly cannot claim to have invented today's algorithm."

And yes, Chebfun is very nice.

That manual is great! I highly recommend reading through it for anyone interested in computing history - and anyone interested in how to communicate clearly in technical writing.

Here's a fun report from someone who had the 205 as their first personal computer:

https://blogs.mathworks.com/cleve/2017/10/09/my-first-computer-the-burroughs-205-datatron/

I loved this story at the end:

> Skip on Minus

> One of the programs that I wrote during the summer at the University of Utah involved a loop that was exactly 20 instructions long. That meant the entire loop could be loaded into the high speed memory. It didn't even need the jump instruction from the end of the loop back to the beginning. It would just continue execution from the duplicate of the code in the next block on the drum. The trouble was there was no way to exit the loop.

> I needed a "Skip on Minus" instruction. There was a "Skip on Plus", but no "Skip on Minus". The EE grad student who maintained the machine was a guy named Joe. I wish I could remember his last name. Anyway, I told him I needed a "Skip on Minus" instruction. He said, "I'll put one in". He disappeared behind the cabinets with circuit diagrams and a soldering iron. He called out, "the op code will be 51".

> I wrote my program with a "51" instruction. It worked perfectly. And so the University of Utah had the only Burroughs 205 Datatron in the world with a "Skip on Minus" instruction.

I'm well aware of how and when it works. It is the recognized idiom to perform in-place operations in practice in Matlab code [1].

It would be nice to take advantage of the same idea in mex code where possible. For example an API function that indicates if it is safe to over-write an input array and return it as the corresponding output or not. I don't think Mathworks currently provide documented and future-safe means to do so however. (If one controls all the code, updating things in-place is possible in practice if careful, using the A(1)=A(1); trick I gave in a parent post. Such code is just not guaranteed to be future-proof.)

[1] http://blogs.mathworks.com/loren/2007/03/22/in-place-operations-on-data/

:) Personally I've never used any serious complex analysis in my job (I'm very grateful too, because I always struggled a bit with it). The closest thing I've seen, which I did run into recently, is the use of complex numbers to compute very accurate finite differences. It's one of those delightful tricks that is both elegant and useful: https://blogs.mathworks.com/cleve/2013/10/14/complex-step-differentiation/

I've been working in golang, which fortunately has built-in complex128 types, so it's proved very helpful in a project!