Non-scalar in Uniform output error in arrayfun. Ho

Question

I'm not sure what it means in this context. I tried adding " 'uniformoutput',false " to the end of arrayfun, but then it got upset with the "+" operator saying "Undefined operator '+' for input arguments of type 'cell'." I changed it to ".+" but got a parse error :( What am doing wrong?

Here is an image of the part that is broken and the error. The entire code is below in case someone would like to try running it or copy the broken part.

the whole code:

    function gbp2(zi,zf)

global beam xlist ylist wi qi Ri wf qf Rf Psii Psif x n
beam = struct('E',[],'lambda',[],'w',[],'R',[],'q',[],'a',[]);

E = 1;                              % electric field
lambda = 1064*10^-9;                % wavelength
k = 2*pi/lambda;                    % wave number
wi = 10^-3;                         % initial waist width (minimum spot size)
zr = (pi*wi^2)/lambda;              % Rayleigh range
Ri = zi + zr^2/zi;
qi = 1/(1/Ri-1i*lambda/(pi*wi^2));  % initial complex beam parameter
Psii = atan(real(qi)/imag(qi));     % Gouy phase
mat = [1 zf; 0 1];                  % transformation matrix
A = mat(1,1); B = mat(1,2); C = mat(2,1); D = mat(2,2);
qf = (A*qi + B)/(C*qi + D);
wf = sqrt(-lambda/pi*(1/imag(1/qf)));
Rf = 1/real(1/qf);

u = @(z, coor, mode, w, R, Psi) (2/pi)^(1/4)*sqrt(exp(1i*(2*mode+1)*Psi)/(2^mode*factorial(mode)*w))*...
            hermiteH(mode,sqrt(2)*coor/w)*exp(-coor^2*(1/w^2+1i*k/(2*R))-1i*k*z);

% -------------------- ERROR IN THIS PIECE (below) ----------------------------
xlist = containers.Map('KeyType','double','ValueType','any');
ylist = containers.Map('KeyType','double','ValueType','any');

    function pts(z, w, R, Psi)
        xlist(z) = -2*w:10^-4:2*w;
        ylist(z) = zeros(1,size(xlist(z),2));
        for mode = 0:2:10
            ylist(z) = ylist(z) + arrayfun(@(coor) u(z, coor, mode, w, R, Psi),xlist(z),'uniformoutput',false);
        end
    end

pts(zi,wi,Ri,Psii)
pts(zf,wf,Rf,Psif)
plot(xlist(zi),ylist(zi),xlist(zf),ylist(zf)) end

I tried writing a similar but simpler function and it seems to work just fine:

function test(zi, zf)

u = @(z,coor,mode) z*mode + integral(@(coor)coor,0,1);
xlist = containers.Map('KeyType','double','ValueType','any');
ylist = containers.Map('KeyType','double','ValueType','any');

    function pts(z)
        xlist(z) = -5:5;
        ylist(z) = zeros(1,size(xlist(z),2));
        for mode = 0:2:10
            ylist(z) = ylist(z) + arrayfun(@(coor) u(z,coor,mode),xlist(z));
        end
    end

pts(zi)
pts(zf)
plot(xlist(zi),ylist(zi),xlist(zf),ylist(zf))

end

so I'm not sure what the problem is with my code.

Answer 1

The error message give you a big hint where to look:

Undefined operator '+' for input arguments of type 'cell'.

Error in gbp2/pts (line 36)
                   ylist(z) = ylist(z) + arrayfun(@(coor) u(z, coor, mode, w, R,
                   Psi),xlist(z).','uniformoutput',false);

From the documentation for arrayfun and the 'UniformOutput' option:

Requests that the arrayfun function combine the outputs into cell arrays B1,...,Bm. The outputs of function func can be of any size or type.

Indeed, if you check, the output of this line is a cell array:

arrayfun(@(coor) u(z, coor, mode, w, R, Psi),xlist(z),'uniformoutput',false)

You can't sum the values from the cell array directly. Here's one of several ways you can do this:

v = arrayfun(@(coor) u(z, coor, mode, w, R, Psi),xlist(z),'uniformoutput',false);
ylist(z) = ylist(z) + [v{:}];

However, I don't see why you need to used the 'UniformOutput' option or even the slow arrayfun at all. Just vectorize your function with respect to coor:

u = @(z, coor, mode, w, R, Psi)(2/pi)^(1/4)*sqrt(exp(1i*(2*mode+1)*Psi)/(2^mode*factorial(mode)*w))*...
        hermiteH(mode,sqrt(2)*coor/w).*exp(-coor.^2*(1/w^2+1i*k/(2*R))-1i*k*z);

Now

ylist(z) = ylist(z) + u(z, xlist(z), mode, w, R, Psi);

Some additional suggestions: Don't use global variables – they're inefficient and there are almost always better solutions. If pts is meant to be a nested function, you're missing a closing end for the main gbp2 function. It might be a good idea to rename your variable mode so that it doesn't overload the built-in function of the same name. Psif isn't defined. And zeros(1,size(xlist(z),2)) can be written simply as zeros(size(xlist(z))).