Homework, 15.1: 
#1, 5, 8, 12, 15
MAT329
Andy Long

#1: we approximate the integral with upper right points, and then midpoints:

f[x_,y_]:=x y
upperRightCorners=Table[f[x,y],{x,2,6,2},{y,2,4,2}]
midpoints=Table[f[x,y],{x,1,5,2},{y,1,3,2}]
deltaA=2*2
approxCorner=Total[Total[upperRightCorners]]*deltaA
approxMid=Total[Total[midpoints]]*deltaA

#5: This problem allows us to use the masks which I discussed in class. I’ll use mine, in addition....

data={
{-3,-5,-6,-4,-1},
{-1,-2,-3,-1,1},
{1,0,-1,1,4},
{2,2,1,3,7},
{3,4,2,5,9}
}
partA=(-2+2+(-1)+3)*1*2
weights00={
{1,1,1,1,0},
{1,1,1,1,0},
{1,1,1,1,0},
{1,1,1,1,0},
{0,0,0,0,0}};
partB=Total[Total[weights00*data]]*.5*1
BetterWeights={
{1,2,2,2,1},
{2,4,4,4,2},
{2,4,4,4,2},
{2,4,4,4,2},
{1,2,2,2,1}
}/4;
partB=Total[Total[BetterWeights*data]]*.5*1

#8: The midpoint rule with m=n=2 suggests that we simply need to estimate the function at the center of each square. Then, since each square has unit area, we multiply each value by 1 and add them up! I.e.,

result=(3.3+5.1+1.3+3.4) (* or so! :*)

#12: The height over this rectangular region [0,5] x [0,3] is given by 5-x -- that’s a plane, independent of y, so it’s a wedge. The volume will be average of the max and min heights times the area of the region. That is, the average of 0 and 5, 5/2, times 15.

volume=5/2*15
f[x_,y_]:=5-x
Plot3D[f[x,y],{x,0,5},{y,0,3}]
Integrate[f[x,y],{x,0,5},{y,0,3}]

#15: we need to find some nice way to evaluate the function at the midpoints, for a variable number of intervals:

f[x_,y_]:=(Sqrt[1+x*(E^(-y))])
(* here's what we should get: *)
NIntegrate[f[x,y],{x,0,1},{y,0,1}]
(* now we'll create an approximating function *)
approx[n_]:=
Module[
{m=Sqrt[n]},
x[i_]:=1/(2m)+(i-1)/m;
y[i_]:=1/(2m)+(i-1)/m;
1/n*Total[Total[Table[f[x[i],y[j]],{i,1,m},{j,1,m}]]]
]
(* And here are our approximations: *)
N[approx[1]]
N[approx[4]]
N[approx[16]]
N[approx[64]]
N[approx[256]]
N[approx[1024]]

Now here’s a more general version of that. I’m going to allow a different number of subintervals in each direction: m in the x direction, and n in the y direction.

f[x_,y_]:=(Sqrt[1+x*(E^(-y))])
(* here's what we should get: *)
N[Integrate[f[x,y],{x,0,1},{y,0,1}],9]
(* now we'll create an approximating function *)
approx[m_, n_]:=
Module[
{}, (* this is for local variables *)
x[i_]:=1/(2m)+(i-1)/m;
y[j_]:=1/(2n)+(j-1)/n;
Print[Table[x[i],{i,1,m}]];
1/(m*n)*Total[Total[Table[f[x[i],y[j]],{i,1,m},{j,1,n}]]]
]
(* And here are our approximations: *)
N[approx[1,1],9]
N[approx[2,2],9]
N[approx[4,4],9]
N[approx[8,8],9]
N[approx[16,16],9]
N[approx[32,32],9]
N[approx[64,64],9]

An even more general method would allow for different limits of integration, so here we go with that:

f[x_,y_]:=(Sqrt[1+x*(E^(-y))])
a=0
b=1
c=0
d=1
(* here's what we should get: *)
NIntegrate[f[x,y],{x,a,b},{y,c,d}]
(* now we'll create an approximating function that includes the endpoints: *)
approx[a_,b_,c_,d_,m_, n_]:=
Module[
{
xrange=b-a,
yrange=d-c
}, (* this is for local variables *)
x[i_]:=a+xrange(1/(2m)+(i-1)/m);
y[j_]:=c+yrange(1/(2n)+(j-1)/n);
xrange*yrange/(m*n)*Total[Total[Table[f[x[i],y[j]],{i,1,m},{j,1,n}]]]
]
(* And here are our approximations: *)
N[approx[a,b,c,d,1,1],9]
N[approx[a,b,c,d,2,2],9]
N[approx[a,b,c,d,4,4],9]
N[approx[a,b,c,d,8,8],9]
N[approx[a,b,c,d,16,16],9]
N[approx[a,b,c,d,32,32],9]
N[approx[a,b,c,d,64,64],9]

Finally we include the function as an argument to the procedure, and it becomes a self-contained mid-point method, for any rectangular interval, and any function:

Here is a function I call “MidpointApprox”:

f[x_,y_]:=(Sqrt[1+x*(E^(-y))])
a=0
b=1
c=0
d=1
(* here's what we should get: *)
NIntegrate[f[x,y],{x,a,b},{y,c,d}]
(* 
We create a general midpoint approximation function, that takes as arguments the integrand function, the endpoints in the x and y directions, and the number of sub-intervals.

Notice that I've done a little simplifying of the evaluations of x[i] and y[j]. This would speed up the code. Might as well take care of all algebra up front!
 *)
MidpointApprox[f_,a_,b_,c_,d_,m_, n_]:=
Module[
{
Deltax=(b-a)/m,
Deltay=(d-c)/n
}, (* this is for local variables *)
x[i_]:=a+Deltax(i-1/2);
y[j_]:=c+Deltay(j-1/2);
(* This was my error: I had 1/(m*n) as the subregion area, rather than the correct subregion area, which is Deltax*Deltay: *)
Deltax*Deltay*Total[Total[Table[f[x[i],y[j]],{i,1,m},{j,1,n}]]]
]
(* And here are our approximations: *)
N[MidpointApprox[f,a,b,c,d,64,64],9]
(* but if we define another function, call it g, we can use it directly in MidpointApprox: *)
g[x_,y_]:=x y
N[MidpointApprox[g,a,b,c,d,64,64],9]
N[MidpointApprox[g,-1,1,0,1,5,5],9]

Check these values! Let’s take a look at 2-4, p. 1005:

(* #2: *)
g[x_,y_]:=1-x y^2
a=0;
b=4;
c=-1;
d=2;
m=2;
n=3;
N[MidpointApprox[g,a,b,c,d,m,n]]
NIntegrate[g[x,y],{x,a,b},{y,c,d}]

(* #3: *)
g[x_,y_]:=x E^(-x y)
a=0;
b=2;
c=0;
d=1;
m=2;
n=2;
N[MidpointApprox[g,a,b,c,d,m,n]]
NIntegrate[g[x,y],{x,a,b},{y,c,d}]

(* #4: *)
g[x_,y_]:=1+x^2+3 y
a=1;
b=2;
c=0;
d=3;
m=2;
n=2;
N[MidpointApprox[g,a,b,c,d,m,n]]
NIntegrate[g[x,y],{x,a,b},{y,c,d}]

Now here’s how I made the change to UpperLeftApprox -- can you see what’s changed from MidpointApprox? How would you create LowerRightApprox?

UpperLeftApprox[f_,a_,b_,c_,d_,m_, n_]:=
Module[
{
Deltax=(b-a)/m,
Deltay=(d-c)/n
}, (* this is for local variables *)
x[i_]:=a+Deltax(i-1);
y[j_]:=c+Deltay*j;
(* This was my error: I had 1/(m*n) as the subregion area, rather than the correct subregion area, which is Deltax*Deltay: *)
Deltax*Deltay*Total[Total[Table[f[x[i],y[j]],{i,1,m},{j,1,n}]]]
]
(* #2 checks out: *)
g[x_,y_]:=1-x y^2
a=0;
b=4;
c=-1;
d=2;
m=2;
n=3;
N[UpperLeftApprox[g,a,b,c,d,m,n]]
NIntegrate[g[x,y],{x,a,b},{y,c,d}]