Permalink
Cannot retrieve contributors at this time
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?
05_curve_fitting/README.md
Go to fileThis commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
347 lines (308 sloc)
9.59 KB
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# 05_Curve_Fitting | |
ME 3255 Homework 5 | |
# Problem 2 | |
##Least Squares Function | |
``` | |
function [a,fx,r2] = least_squares_2(Z,y) | |
%least squares solution for input matrix Z | |
%and measurements y | |
%outputs are: | |
%a: constants for best fit function | |
%fx: function evaluations at xi | |
%r2: coefficient of determination | |
a = Z\y; | |
fx = Z*a; | |
e=y-fx; | |
St=std(y); | |
Sr=std(e); | |
r2 = 1-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=[ones(size(xa)) xa xa.^-1]; | |
[a_a,fx_a,r2_a] = least_squares_2 (Z,ya) | |
xa_fcn = linspace(min(xa),max(xa)); | |
a = figure(1); | |
plot(xa,ya,'o',xa_fcn,a_a(1)+a_a(2)*xa_fcn+a_a(3)*xa_fcn.^-1) | |
saveas(a, 'figure_01.png') | |
``` | |
## 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=[ones(size(xb)) xb xb.^2 xb.^3]; | |
[a_b,fx_b,r2_b] = least_squares_2 (Z,yb) | |
xb_fcn = linspace(min(xb),max(xb)); | |
b = figure(2); | |
plot(xb,yb,'o',xb_fcn,a_b(1)+a_b(2)*xb_fcn.^1+a_b(3)*xb_fcn.^2+a_b(4)*xb_fcn.^3) | |
saveas(b, 'figure_02.png') | |
``` | |
##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_c,fx_c,r2_c] = least_squares_2 (Z,yc) | |
xc_fcn = linspace(min(xc),max(xc)); | |
c = figure(3); | |
plot(xc,yc,'o',xc_fcn,a_c(1)*exp(-1.5.*xc_fcn) + a_c(2)*exp(-0.3.*xc_fcn) + a_c(3)*exp(-0.05.*xc_fcn)); | |
saveas(c, 'figure_03.png') | |
``` | |
## 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_d,fx_d,r2_d] = least_squares_2 (Z,yd) | |
xd_fcn = linspace(min(xd),max(xd)); | |
d = figure(4); | |
plot(xd,yd,'o',xd_fcn,a_d(1)*sin(xd_fcn)+a_d(2)*sin(3*xd_fcn)); | |
saveas(d, 'figure_04.png') | |
``` | |
##Output: | |
``` | |
Part A: r2 = 0.9801 | |
Part B: r2 = 0.9232 | |
Part C: r2 = 0.9462 | |
Part D: r2 = 0.9219 | |
``` | |
Part A | |
![Part A](./Figures/figure_10.png) | |
Part B | |
![Part B](./Figures/figure_11.png) | |
Part C | |
![Part C](./Figures/figure_12.png) | |
Part D | |
![Part D](./Figures/figure_13.png) | |
# Problem 3 | |
``` | |
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; | |
abs_accuracy = abs(accuracy-val); | |
[~, index] = min(abs_accuracy) | |
closest_value = accuracy(index) | |
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; | |
abs_accuracy = abs(precision-val); | |
[~, index] = min(abs_accuracy) | |
closest_value = precision(index) | |
``` | |
Part A | |
``` | |
User number 32 was closest with a value of 0.0044. | |
``` | |
Part B | |
```matlab | |
User number 33 was the most precise with a standard deviation of 3.4073 | |
``` | |
# Problem 4 | |
``` | |
Part A | |
```matlab | |
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 | |
``` | |
``` | |
```matlab | |
Mean Buckle Load = 164.4154 N | |
Standard Deviation = 63.6238 N | |
``` | |
Part B | |
``` | |
```matlab | |
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 | |
``` | |
```matlab | |
L = 5.0108 m | |
``` | |
#Problem 5 | |
Part A | |
``` | |
```matlab | |
function [Cd_out] = sphere_drag(Re_in,spline_type) | |
data = [2 0.52 | |
5.8 0.52 | |
16.8 0.52 | |
27.2 0.5 | |
29.9 0.49 | |
33.9 0.44 | |
36.3 0.18 | |
40 0.074 | |
46 0.067 | |
60 0.08 | |
100 0.12 | |
200 0.16 | |
400 0.19]; | |
Cd_out=interp1(data(:,1),data(:,2),Re_in,spline_type); | |
rho=1.3;V=linspace(4,40);D=23.5e-2;mu=1.78e-5; | |
Re=rho*V*D/mu*1e-4; | |
end | |
``` | |
``` | |
```matlab | |
Cd_out1 = sphere_drag(150,'linear') | |
Cd_out2 = sphere_drag(150,'spline') | |
Cd_out3 = sphere_drag(150,'pchip') | |
``` | |
Output | |
``` | |
```matlab | |
Linear = 0.1400 | |
Spline = 0.1428 | |
Pchip = 0.1447 | |
``` | |
Part B Code | |
``` | |
```matlab | |
function [Cd_out] = sphere_drag(Re_in,spline_type) | |
data = [2 0.52 | |
5.8 0.52 | |
16.8 0.52 | |
27.2 0.5 | |
29.9 0.49 | |
33.9 0.44 | |
36.3 0.18 | |
40 0.074 | |
46 0.067 | |
60 0.08 | |
100 0.12 | |
200 0.16 | |
400 0.19]; | |
Cd_out=interp1(data(:,1),data(:,2),Re_in,spline_type) | |
Re=linspace(data(1,1),data(end,1),1000); | |
rho=1.3;V=linspace(4,40);D=23.5e-2;mu=1.78e-5; | |
Re=rho*V*D/mu*1e-4; | |
end | |
``` | |
Part B Script | |
``` | |
```matlab | |
Re=linspace(data(1,1),data(end,1),1000); | |
rho=1.3;V=linspace(4,40);D=23.5e-2;mu=1.78e-5; | |
Re=rho*V*D/mu*1e-4 | |
a=figure(1); | |
cd_interp=sphere_drag(linspace(data(1,1),data(end,1),100),'linear'); | |
plot(data(:,1),data(:,2),'o',Re,cd_interp,'b-','Linewidth',4,'Markersize',10) | |
hold on; | |
cd_interp=sphere_drag(linspace(data(1,1),data(end,1),100),'spline'); | |
plot(Re,cd_interp,'r-', 'Linewidth',4,'Markersize',10) | |
cd_interp=sphere_drag(linspace(data(1,1),data(end,1),100),'pchip'); | |
plot(Re,cd_interp,'g-', 'Linewidth',4,'Markersize',10) | |
legend('data','linear','spline','pchip') | |
saveas(a,'figure_14.png') | |
``` | |
Script | |
``` | |
```matlab | |
Cd_out1 = sphere_drag(10,'linear') | |
Cd_out2 = sphere_drag(10,'pchip') | |
Cd_out3 = sphere_drag(10,'spline') | |
plot(data(:,1),data(:,2),'o',Re,Cd_out,'b-','Linewidth',4,'Markersize',10) | |
plot(data(:,1),data(:,2),Re,Cd_out2,'r-', 'Linewidth',4,'Markersize',10) | |
plot(data(:,1),data(:,2),Re,Cd_out3,'g-', 'Linewidth',4,'Markersize',10) | |
legend('data','linear','pchip','spline') | |
``` | |
![Part B](./Figures/figure_14.png) | |
# Problem 6 | |
##Analytical Solution | |
``` | |
```matlab | |
fun = @(x) (1/6).*(x).^3+(1/2).*(x).^2+(x); | |
q = integral(fun,2,3) | |
``` | |
##Gauss Quadtrature | |
``` | |
```matlab | |
%One 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) | |
``` | |
##Two Point | |
``` | |
%Two 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) | |
``` | |
%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) | |
``` | |
``` | |
matlab | |
| 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% | | |
``` |