Skip to content
Permalink
482822af46
Switch branches/tags

Name already in use

A tag already exists with the provided branch name. Many Git commands accept both tag and branch names, so creating this branch may cause unexpected behavior. Are you sure you want to create this branch?
Go to file
 
 
Cannot retrieve contributors at this time
219 lines (193 sloc) 7.41 KB

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 A

###Part B: 0.9112

Part B

###Part C: 0.9462

Part C

###Part D: 0.9219

Part D

##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

Drag Force vs. Velocity  (3 Interpolation Methods)

#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% |