05_curve_fitting

ME 3255 Homework 5

Problem 2

Least Squares function

function [a,fx,r2] = least_squares(Z,y);
  % inputs: Z matrix and measurements y
  % outputs:
  % a   - best fit constant
  % fx  - function at xi
  % r2  - coefficient of determination
  a = Z/y;
  fx = Z*a;
  e = y-fx;
  st = std(y);
  sr = std(e);
  r2 = (st-sr)/st;
end

set defaults
% Part A
xa = [1 2 3 4 5]';
ya = [2.2 2.8 3.6 4.5 5.5]';
Z = [xa.^0,xa.^1,1./xa];
[a,fx,r2] = least_squares(Z,ya)
plot(xa,fx,xa,ya,'o','Linewidth',2)

% Part B
xb = [0 2 4 6 8 10 12 14 16 18]';
yb = [21.5 20.84 23.19 22.69 30.27 40.11 43.31 54.79 70.88 89.48]';
Z = [xb.^0,-xb.^1,xb.^2];
[a,fx,r2] = least_squares(Z,yb)
plot(xb,fx,xb,yb,'o','Linewidth',2)

% Part C
xc = [0.5 1 2 3 4 5 6 7 9]';
yc = [6 4.4 3.2 2.7 2.2 1.9 1.7 1.4 1.1]';
Z = [exp(-1.5*xc),exp(-0.3*xc),exp(-0.05*xc)];
[a,fx,r2] = least_squares(Z,yc)
plot(xc,fx,xc,yc,'o','Linewidth',2)

% Part D
xd = [0.00000000e+00 1.26933037e-01 2.53866073e-01 3.80799110e-01 5.07732146e-01 6.34665183e-01 7.61598219e-01 8.88531256e-01 1.01546429e+00 1.14239733e+00 1.26933037e+00 1.39626340e+00 1.52319644e+00 1.65012947e+00 1.77706251e+00 1.90399555e+00 2.03092858e+00 2.15786162e+00 2.28479466e+00 2.41172769e+00 2.53866073e+00 2.66559377e+00 2.79252680e+00 2.91945984e+00 3.04639288e+00 3.17332591e+00 3.30025895e+00 3.42719199e+00 3.55412502e+00 3.68105806e+00 3.80799110e+00 3.93492413e+00 4.06185717e+00 4.18879020e+00 4.31572324e+00 4.44265628e+00 4.56958931e+00 4.69652235e+00 4.82345539e+00 4.95038842e+00 5.07732146e+00 5.20425450e+00 5.33118753e+00 5.45812057e+00 5.58505361e+00 5.71198664e+00 5.83891968e+00 5.96585272e+00 6.09278575e+00 6.21971879e+00 6.34665183e+00 6.47358486e+00 6.60051790e+00 6.72745093e+00 6.85438397e+00 6.98131701e+00 7.10825004e+00 7.23518308e+00 7.36211612e+00 7.48904915e+00 7.61598219e+00 7.74291523e+00 7.86984826e+00 7.99678130e+00 8.12371434e+00 8.25064737e+00 8.37758041e+00 8.50451345e+00 8.63144648e+00 8.75837952e+00 8.88531256e+00 9.01224559e+00 9.13917863e+00 9.26611167e+00 9.39304470e+00 9.51997774e+00 9.64691077e+00 9.77384381e+00 9.90077685e+00 1.00277099e+01 1.01546429e+01 1.02815760e+01 1.04085090e+01 1.05354420e+01 1.06623751e+01 1.07893081e+01 1.09162411e+01 1.10431742e+01 1.11701072e+01 1.12970402e+01 1.14239733e+01 1.15509063e+01 1.16778394e+01 1.18047724e+01 1.19317054e+01 1.20586385e+01 1.21855715e+01 1.23125045e+01 1.24394376e+01 1.25663706e+01]';
yd=[9.15756288e-02 3.39393873e-01 6.28875306e-01 7.67713096e-01 1.05094584e+00 9.70887288e-01 9.84265740e-01 1.02589034e+00 8.53218113e-01 6.90197665e-01 5.51277193e-01 5.01564914e-01 5.25455797e-01 5.87052838e-01 5.41394658e-01 7.12365594e-01 8.14839678e-01 9.80181855e-01 9.44430709e-01 1.06728057e+00 1.15166322e+00 8.99464065e-01 7.77225453e-01 5.92618124e-01 3.08822183e-01 -1.07884730e-03 -3.46563271e-01 -5.64836023e-01 -8.11931510e-01 -1.05925186e+00 -1.13323611e+00 -1.11986890e+00 -8.88336727e-01 -9.54113139e-01 -6.81378679e-01 -6.02369117e-01 -4.78684439e-01 -5.88160325e-01 -4.93580777e-01 -5.68747320e-01 -7.51641934e-01 -8.14672884e-01 -9.53191554e-01 -9.55337518e-01 -9.85995556e-01 -9.63373597e-01 -1.01511061e+00 -7.56467517e-01 -4.17379564e-01 -1.22340361e-01 2.16273929e-01 5.16909714e-01 7.77031694e-01 1.00653798e+00 9.35718089e-01 1.00660116e+00 1.11177057e+00 9.85485116e-01 8.54344900e-01 6.26444042e-01 6.28124048e-01 4.27764254e-01 5.93991751e-01 4.79248018e-01 7.17522492e-01 7.35927848e-01 9.08802925e-01 9.38646871e-01 1.13125860e+00 1.07247935e+00 1.05198782e+00 9.41647332e-01 6.98801244e-01 4.03193543e-01 1.37009682e-01 -1.43203880e-01 -4.64369445e-01 -6.94978252e-01 -1.03483196e+00 -1.10261288e+00 -1.12892727e+00 -1.03902484e+00 -8.53573083e-01 -7.01815315e-01 -6.84745997e-01 -6.14189417e-01 -4.70090797e-01 -5.95052432e-01 -5.96497000e-01 -5.66861911e-01 -7.18239679e-01 -9.52873043e-01 -9.37512847e-01 -1.15782985e+00 -1.03858206e+00 -1.03182712e+00 -8.45121554e-01 -5.61821980e-01 -2.83427014e-01 -8.27056140e-02]';
Z = [sin(xd),sin(3*xd)];
[a,fx,r2] = least_squares(Z,yd)
plot(xd,fx,xd,yd,'o','Linewidth',2)

Output:

Part A: 0.9801

Part B: 0.9112

Part C: 0.9462

Part D: 0.9219

Problem 3

Part A

darts = dlmread('compiled_data.csv',',',1,0);
x_darts =darts(:,2).*cosd(darts(:,3));
y_darts =darts(:,2).*sind(darts(:,3));
ur = darts(:,1);
i = 1;
for r = 0:32
    i_interest = find(ur==r);
    mx_ur(i)= mean(x_darts(i_interest));
    my_ur(i)= mean(y_darts(i_interest));
    i = i +1;
end
    accuracy = mx_ur + my_ur;
    val = 0; %value trying to be closest to
    abs_accuracy = abs(accuracy-val); %takes absolute value of accuracy vector
    [~, index] = min(abs_accuracy)
    closest_value = accuracy(index)

Most accurate dart thrower: User 31

• combined mean of x & y of 0.0044cm from zero

Part B

darts = dlmread('compiled_data.csv',',',1,0);
x_darts =darts(:,2).*cosd(darts(:,3));
y_darts =darts(:,2).*sind(darts(:,3));
ur = darts(:,1);
i = 1;
for r = 0:32
    i_interest = find(ur==r);
    stdx_ur(i)= std(x_darts(i_interest));
    stdy_ur(i)= std(y_darts(i_interest));
    i = i +1;
end
    precision = stdx_ur + stdy_ur;
    val = 0; %value trying to be closest to
    abs_accuracy = abs(precision-val); %takes absolute value of accuracy vector
    [~, index] = min(abs_accuracy)
    closest_value = precision(index)

Most precise dart thrower: User 32

• combined x & y STD of 3.4073

Problem 4

Part A

function [mean_buckle_load,std_buckle_load]=buckle_monte_carlo(E,r_mean,r_std,L_mean,L_std);

N=100;
r_mean=0.01;
r_std=.001;
L_mean=5;
L_std=L_mean*0.01;

rrand=normrnd(r_mean,r_std,[N,1]);
Lrand=normrnd(L_mean,L_std,[N,1]);
E = 200*10^9;

P_cr=(pi^3*E.*rrand.^4)./(16.*Lrand.^2);
x1 = mean(P_cr)
x2 = std(P_cr)
P_test=1000*9.81/100;
sum(P_cr<P_test)
end

Output:

• Mean_buckle_load = 164.4154 N
• STD = 63.6238 N

Part B

function [mean_buckle_load,std_buckle_load]=buckle_monte_carlo_2(E,r_mean,r_std,P_cr);

N=100;
r_mean=0.01;
r_std=.001;
P_test=1000*9.81/100;
sum(2.5<P_test)
P_cr = 164.4154;
rrand=normrnd(r_mean,r_std,[N,1]);
E = 200*10^9;
L = ((pi^3*E.*rrand.^4)./(16*P_cr)).^0.5;
mean(L)
end

Output:

• L = 5.0108 m

Problem 5

Part A

cd_out_linear = sphere_drag(300,'linear') = 0.1750
cd_out_spline = sphere_drag(300,'spline') = 0.1809
cd_out_pchip = sphere_drag(300,'pchip') = 0.1807

Part B

Problem 6

Analytical

fun = @(x) (1/6).*(x).^3+(1/2).*(x).^2+(x);
q = integral(fun,2,3)

Gauss Quadrature

1-point

x = [0];
c = [2];
a = 2;
b = 3;
tildec = (b-a)/2*c;
tildex = (b-a)/2*x + (b+a)/2;
f = (1/6).*(tildex).^3+(1/2).*(tildex).^2+(tildex);
value1 = sum(tildec.*f)

2-point

x = [-0.57735,0.57735];
c = [1,1];
a = 2;
b = 3;
tildec = (b-a)/2*c;
tildex = (b-a)/2*x + (b+a)/2;
f = (1/6).*(tildex).^3+(1/2).*(tildex).^2+(tildex);
value2 = sum(tildec.*f)

3-point

%Three Point
x = [-0.77459666,0,0.77459666];
c = [0.5555555,0.88888888,0.55555555];
a = 2;
b = 3;
tildec = (b-a)/2*c;
tildex = (b-a)/2*x + (b+a)/2;
f = (1/6).*(tildex).^3+(1/2).*(tildex).^2+(tildex);
value3 = sum(tildec.*f)

| Method | Value | Error |

| Analytical | 8.3750 | 0% |

| 1 Gauss Point | 8.2292 | 1.74% |

| 2 Gauss Point | 8.3750 | 0% |

| 3 Gauss Point | 8.3750 | 0% |