Double Integral Using Riemann Sums

%% Double Integral using Riemann Sums
% READ THROUGH EVERYTHING

% Using Reimann sum S#,# to estimate the double integral of *function* over R =
% [matrix1] x [matrix2] using regular partition and upper-right vertices
% of the sub-rectangles as sample points.

% USED FOR f(x,y) = #xy

S = input('Riemann sum divisions');
R1 = input('matrix 1 in brackets');
R2 = input('matrix 2 in brackets');
delta_x = (R1(2)-R1(1))/S(1);
delta_y = (R2(2)-R2(1))/S(2);
delta_A = delta_x * delta_y;
a = R1(1):delta_x:R1(2);
b = R2(1):delta_y:R2(2);

sum = 0;
 for j = 2:length(b);
    for i = 2:length(a);
        c = a(i)*b(j);
        sum  = sum + c;
    end
end

sum*delta_A

%  ANSWER IS sum*delta_A* # from the f(x,y) function