%(1) decide upper bound of R
% pi*rmax^2 / 2 <= 1600
rmax = sqrt(2*1600/pi);

%(2) Define R vector
R=[0.01:0.01:rmax];

%(3) Calculate L vector  (use element-by-element operations, p85)
L=(1600 - pi*R.^2/2)./(2*R);

%(4) compute the array of cost
Cost = 2*(R+L)*30 + pi*R*40;

%(5) get the min cost (use the min function, p77)
% x is the min value, and k is the indice in the Cost vector
[x, k] = min(Cost);
mincost = x

%(6) get the corresponding R and L value
% k is the indice of that value R and L in their vectors
l = L(k)
r = R(k)