% Class 25 
% Coded by NFR on 4.16.2019
%
% Example of Romberg method using Chapra algroithm
% Book problem 20.23
g = 9.81;
m = 80;
Cd = 0.2;
velfunc = @(t) sqrt(g*m/Cd)*tanh(sqrt(g*Cd/m)*t)
Distance = romberg(velfunc,0,8,1)
%
% Book problem 20.3
% Evaluate with 2 point Guass Quadrature
F1 = @(x) (1.5+1.5*x)*exp(2*(1.5+1.5*x))*1.5
% Look up on table 20.1 the coefficients and
% x values to use
co = 1;
c1 = 1;
xo = -1/sqrt(3);
x1 = 1/sqrt(3);
I_GQ = co*F1(xo)+c1*F1(x1)
% Verify with Romberg method
Freg = @(x) x.*exp(2.*x);
[I_ROM,ea,Iter] = romberg(Freg,0,3,0.5)
% Matlab has built in functions
% Most common are the integral, integral2, integral3 series
I_mat = integral(Freg,0,3)
% Matlab does have adaptive quadrature
I_aq = quadgk(Freg,0,3)
% 2 dimensional integrals
% This shows how to plot and perform integration with polar coordinates
%
Rmax = @(theta) 2+cos(theta)-theta./4.*sin(theta);
Rmin = 1;
Theta_vec = linspace(0,2*pi,100);
Rmax_vec = Rmax(Theta_vec);
Rmin_vec = ones(100,1);
polarplot(Theta_vec,Rmin_vec,Theta_vec,Rmax_vec)
I = integral2(@(theta,r) r,0,2*pi,Rmin,Rmax)
F1 = @(theta,r) r;
I2 = integral2(F1,0,2*pi,Rmin,Rmax)
% Using the integral function on a 'local function'
% Problem 20.19 book
% Solving for the distance traveled, we
% simply use the subfunction
distance_rocket = romberg(@subfunc,0,30,1)

function v = subfunc(t)
% Enter the piecewise function
if t >= 0 && t <= 10
    v = 11*t^2-5*t;
elseif t > 10 && t <= 20
    v = 1100-5*t;
elseif t > 20 && t <= 30
    v = 50*t+2*(t-20)^2;
end
end