%% Multiple Cumulative Integrations
% Perform nested integrations in the x and y directions. Plot the results to visualize the cumulative integral value in both directions.
%
% Create a grid of values for the domain.

x = -2:0.1:2;
y = -2:0.2:2;
[X,Y] = meshgrid(x,y);

%
% Calculate the function f\left(x,y\right)=10x^2 +20y^2 on the grid.

F = 10*X.^2 + 20*Y.^2;

% cumtrapz integrates numeric data rather than functional expressions, so in general the underlying function does not need to be known to use cumtrapz on a matrix of data. In cases where the functional expression is known, you can instead use integral, integral2, or integral3.
%
% Use cumtrapz to approximate the double integral
% I\left(a,b\right) = \int_{-2}^{b} \int_{-2}^{a} \left( 10x^2+20y^2 \right) \, dx \,dy.
% To perform this double integration, use nested function calls to cumtrapz. The inner call first integrates the rows of data, then the outer call integrates the columns.

I = cumtrapz(y,cumtrapz(x,F,2));

%
% Plot the surface representing the original function as well as the surface representing the cumulative integration. Each point on the surface of the cumulative integration gives an intermediate value of the double integral. The last value in I gives the overall approximation of the double integral, I(end) = 642.4. Mark this point in the plot with a red star.

surf(X,Y,F,'EdgeColor','none')
xlabel('X')
ylabel('Y')
hold on
surf(X,Y,I,'FaceAlpha',0.5,'EdgeColor','none')
plot3(X(end),Y(end),I(end),'r*')
hold off

% Copyright 2012 The MathWorks, Inc.