Approximating double integral with only sum of sum and no for loops

I've stumbled across this problem which I can't seem to solve. I need to calculate an approximation of an integral.

The constriction is that I can't use for loops or create matrixes (such as meshgrid). I have to used a nested sum of sum approach.

How am I to do this? This is a working example of a single variable integral using sum:

if true
n=100;
a=0; b=1;
f=@(x)x.^2;
x=linspace(a,b,n+1);
h=(b-a)/n;
q=sum(h*f(x(1:n)));  
end

And this is my non-working attempt below. If I somehow could specify in what order to do x or y (with a for loop which I can't use) this would be simple.

Does anyone have a smart solution?

if true
n=100;
m=100;
a=0; b=2; c=0; d=2; 
f=@(x,y) x.^2 + y.^2;
x=linspace(a,b,n+1);
y=linspace(c,d,m+1);
h=(b-a)/n;
l=(d-c)/m;
q=sum ( ( sum(h*l*f(x(1:n),y(1:m)) ) ));
end

Answer 1

[David Young]

How about integral2? There's also a recent blog on double integration which might help you get the most accurate results: http://blogs.mathworks.com/loren/2014/02/12/double-integration-in-matlab-methods-and-handling-discontinuities-singularities-and-more/

Why the restriction on meshgrid if it helps to solve the problem?

[Joel]

As it turned out the only possible way was to use meshgrid and then sum of sum on it. I guess I misunderstood the assignment.