Before the lab

• Complete the Making 3D plots part of last lab.

Background

• GOBG Chapter 2 (Cell arrays)
• Index expressions
• Variable-length argument lists
• Optional Arguments

Using arrayfun

• Start with your solution to Making 3D plots part of last lab or use mine (delta)

• Use an anonymous function to define f so that the following octave code produces the same plot as before. The nicer (although perhaps slightly slower) solution doesn't involve copy pasting the definition of delta.

a = [4;4];
beta = 7.5;
range = [-4:0.1:8];
[X Y] = meshgrid(range,range);

f=
Z=arrayfun(f,X,Y);

surf(X,Y,Z);

• Use the techniques from L18 to measure the speedup from the for loop version. You probably want to comment out the actual plotting step and just measure the matrix computation.

• How sensitive are both versions to the number of elements in range? If you double the size (half the step size), does the time double?

Fully vectorizing

If we expand out a single call to delta in our for loop, we get

Z(i,j) = beta * (a(1)*X(i,j) + a(2)*Y(i,j)) - (beta-1) * (X(i,j)^2 + Y(i,j)^2)

• Use this equation to complete the vectorized function arraydelta that satisfies the following test.

function ret = arraydelta(beta,a,X,Y)

  ret =
endfunction

%!test
%! a = rand(2,1)
%! X = [2,2;2,2]
%! Y = [-1,-1;-1,-1]
%! assert (arraydelta(0, a, X,Y),[5,5; 5,5],eps)

• Use the anonymous function trick from above to complete the following test

%!test
%! a = rand(2,1);
%! beta = rand*6;
%! X = rand(10,10);
%! Y = rand(10,10);
%! assert (arraydelta(beta, a, X,Y), arrayfun(                            , X, Y), eps)

• Finally measure the speedup for new version. How sensitive is it to the size of range?

a = [4;4];
beta = 7.5;

range = [-4:0.1:8];
[X Y] = meshgrid(range,range);

Z=arraydelta(beta,a,X,Y);

surf(X,Y,Z);

Cell Arrays

So far we have mainly concentrated on arrays, which in Octave consist entirely of numbers. These have efficiency advantages for certain computations, but they're not very flexible. Roughly speaking, cell arrays are like multidimensional versions of Python's lists. They support elements of different types, along with random access.

One slightly surprising, but useful application of cell arrays is to use them as lists of vectors.

data = [0,0,0; 0,0,1; 0,1,0; 1,1,1]
cells = num2cell(data,2)
cellfun(@sum,cells)

• How can we change this example to find the column sums?
• Generalize the example to take the function and matrix as a parameter

function out=rowfun(fun,mat)


endfunction

%!test
%!shared data
%! data = [0,0,0; 0,0,1; 0,1,0; 1,1,1];
%! assert(rowfun(@norm,data), [0;1;1;sqrt(3)], eps)

%!assert(rowfun(@(v)(dot(v,v)),data), [0;1;1;3], eps)

%!assert(rowfun(@prod,data), [0;0;0;1], eps)

Variable-length argument lists

It turns out that the function timeit we used in Lab 20 is not quite general enough, since it works only for single argument functions. Octave provides a special input argument varargin, which collects all remaining arguments into a cell array, providing similar functionality to the . rest argument in racket, or the ...rest argument in JavaScript ES2015.

Use varargin to complete the following code for an updated timeit function

# Based on an example from the Julia microbenchmark suite.

function timeit(reps, func,         )
    times = zeros(reps, 1);

    for i=1:reps
      tic(); func(             ); times(i) = toc();
    end

    times = sort(times);
    fprintf ('%s\tmedian=%.3fms median=%.3fms total=%.3fms\n',func2str(func), median(times)*1000,
             mean(times)*1000, sum(times)*1000);
endfunction

%!test "nullary function"
%! timeit(10000,@rand)

%!test "unary function"
%! timeit(10000,@rand,1)

%!test "Binary function"
%! timeit(10000,@plus,1,2)

%!test "Ternery function"
%! timeit(10000,@plus,1,2,3)