%% Question 1

% initial estimate for L1 and L2
M1 = 1.991e30;
M2 = 5.98e24;
R = 149597870700;
r = R * nthroot(M2 / (3 * M1), 3);

% bracket for L1 and L2
a = 0.9 * r;
b = 1.1 * r;

% estimate for r1 and r2
r1 = bisect(@L1, a, b, 1e-3)
r2 = bisect(@L2, a, b, 1e-3)

% compare to fzero
z1 = fzero(@L1, a)
z2 = fzero(@L2, a)


%% Question 2

% area of half circle
A = pi / 2;

% estimates of xbar and ybar for n = 2, 3, 4, ..., 20
% note: we won't use xbar(1) or ybar(1) because we aren't
% asked for n = 1
xbar = zeros(1, 20);
ybar = zeros(1, 20);

for n = 2:20
    xbar(n) = 1 / A * tintegrate(@fxbar, 0, 2, n);
    ybar(n) = 1 / A * tintegrate(@fybar, 0, 2, n);
end

% true values of xbar and ybar
xbarTrue = 1;
ybarTrue = 2 / (3 * A);

% error
xbarErr = xbar(2:20) - xbarTrue;
ybarErr = ybar(2:20) - ybarTrue;

plot(2:20, xbarErr);
xlabel('number of panels');
ylabel('error (xbar - true value)');
title('error in xbar');

figure
plot(2:20, ybarErr);
xlabel('number of panels');
ylabel('error (ybar - true value)');
title('error in ybar');
