diff --git a/HW5/README.html b/HW5/README.html new file mode 100644 index 0000000..8b9c93d --- /dev/null +++ b/HW5/README.html @@ -0,0 +1,186 @@ + + + + + + + + + + + + +

Homework #5

+

due 4/14 by 11:59pm

+

1. Create a new github repository called ‘05_curve_fitting’.

+
    +

  1. Add rcc02007 and zhs15101 as collaborators.

    2. +

  2. Clone the repository to your computer.

    3. +
+

2. Create a least-squares function called least_squares.m that accepts a Z-matrix and dependent variable y as input and returns the vector of best-fit constants, a, the best-fit function evaluated at each point \(f(x_{i})\), and the coefficient of determination, r2.

+
[a,fx,r2]=least_squares(Z,y);
+

Test your function on the sets of data in script problem_2_data.m and show that the following functions are the best fit lines. Report the coefficient of determination in your README.

+
    +

  1. y=0.3745+0.98644x+0.84564/x

    2. +

  2. y= 22.47-1.36x+0.28x^2

    3. +

  3. y=4.0046e(-1.5x)+2.9213e(-0.3x)+1.5647e^(-0.05x)

    4. +

  4. y=0.99sin(t)+0.5sin(3t)

    5. +
+

3. Load the data from the class dart experiment. Show your work in a script with filename dart_statistics.m in your repository.

+
    +

  1. Which user (0-32) was the most accurate dart thrower e.g. mean(x)+mean(y) was closest to 0 cm?

    2. +

  2. Which user (0-32) was the most precise dart thrower e.g. std(x)+std(y) was closest to 0 cm?

    3. +
+

4. 100 steel rods are going to be used to support a 1000 kg structure. The rods will buckle when the load in any rod exceeds the critical buckling load

+

\(P_{cr}=\frac{\pi^3 Er^4}{16L^2}\) buckle

+

where E=200e9 Pa, r=0.01 m +/-0.001 m, and L is the length of the rod which can only be controlled to 1% tolerance. Create a Monte Carlo model buckle_monte_carlo.m that predicts the mean and standard deviation of the buckling load for 100 samples with normally distributed dimensions r and L.

+

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

+
    +

  1. What is the mean_buckle_load and std_buckle_load for L=5 m?

    2. +

  2. What length, L, should the beams be so that only 2.5% will reach the critical buckling load?

    3. +
+

5. The drag coefficient for spheres such as sporting balls is known to vary as a function of the Reynolds number Re, a dimensionless number that gives a measure of the ratio of inertial forces to viscous forces:

+

\(Re = \frac{\rho V D}{\mu}\) Reynolds

+

where ρ= the fluid’s density (kg/m3), V= its velocity (m/s), D= diameter (m), and μ= dynamic viscosity (N⋅s/m2). The following table provides values for a smooth spherical ball:

+ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
Re (x1e-4)C_D
20.52
5.80.52
16.80.52
27.20.5
29.90.49
33.90.44
36.30.18
400.074
460.067
600.08
1000.12
2000.16
4000.19
+
    +
  1. Create a function sphere_drag.m that outputs the drag coefficient based on the given table and an input Reynolds number using a spline interpolation of either linear (‘linear’), piecewise cubic (‘pchip’), or cubic (‘cubic’):
    2. +
+

[Cd_out]=sphere_drag(Re_in,spline_type)

+
    +
  1. Use the following physical constants to plot the drag force vs velocity for a baseball: ρ= 1.3 kg/m3, V= 4 - 40 (m/s), D= 23.5 cm, and μ=1.78e-5 Pa-s. Plot all three interpolation methods on a single plot. Show the plot in your README.
    2. +
+

6. Evaluate the integral of the following function:

+

\(\int_{2}^{3}f(x) = \int_{2}^{3} 1/6x^3 + 1/2x^2+x dx\)

+
    +

  1. analytically

    2. +

  2. with 1-point Gauss quadrature

    3. +

  3. with 2-point Gauss quadrature

    4. +

  4. with 3-point Gauss quadrature

    5. +

  5. include the results in a table in your README

    6. +
+ + + + + + + + + + + + + + + + + + + + + + + + + +
methodvalueerror
analytical…0%
1 Gauss point…..%
………
+ + diff --git a/HW5/README.md b/HW5/README.md new file mode 100644 index 0000000..bb12d64 --- /dev/null +++ b/HW5/README.md @@ -0,0 +1,119 @@ +# Homework #5 +## due 11/17 by 11:59pm + + +### Problem 1 due 11/9 + +**1\.** Create a new github repository called '05_curve_fitting'. + +a. Add rcc02007 and zhs15101 as collaborators. + +b. Clone the repository to your computer. + + +**2\.** Create a least-squares function called `least_squares.m` that accepts a Z-matrix and +dependent variable y as input and returns the vector of best-fit constants, a, the +best-fit function evaluated at each point $f(x_{i})$, and the coefficient of +determination, r2. + +```matlab +[a,fx,r2]=least_squares(Z,y); +``` + +Test your function on the sets of data in script `problem_2_data.m` and show that the +following functions are the best fit lines. Report the coefficient of determination in your README. + +a. y=0.3745+0.98644x+0.84564/x + +b. y= 22.47-1.36x+0.28x^2 + +c. y=4.0046e^(-1.5x)+2.9213e^(-0.3x)+1.5647e^(-0.05x) + +d. y=0.99sin(t)+0.5sin(3t) + +**3\.** Load the data from the [class dart +experiment](https://github.com/cooperrc/stats_and_linear_regression/blob/master/compiled_data.csv). +Show your work in a script with filename `dart_statistics.m` in your repository. + +a. Which user (0-32) was the most accurate dart thrower e.g. mean(x)+mean(y) was closest to 0 cm? + +b. Which user (0-32) was the most precise dart thrower e.g. std(x)+std(y) was closest to 0 cm? + +**4\.** 100 steel rods are going to be used to support a 1000 kg structure. The +rods will buckle when the load in any rod exceeds the [critical buckling +load](https://en.wikipedia.org/wiki/Euler%27s_critical_load) + +$P_{cr}=\frac{\pi^3 Er^4}{16L^2}$ +![buckle](./equations/buckle.png) + +where E=200e9 Pa, r=0.01 m +/-0.001 m, and L is the length of the rod which can only be +controlled to 1% tolerance. Create a Monte +Carlo model `buckle_monte_carlo.m` that predicts the mean and standard deviation of the buckling load for 100 +samples with normally distributed dimensions r and L. + +`[mean_buckle_load,std_buckle_load]=buckle_monte_carlo(E,r_mean,r_std,L_mean,L_std)` + +a. What is the mean_buckle_load and std_buckle_load for L=5 m? + +b. What length, L, should the beams be so that only 2.5% will reach the critical buckling +load? + +**5\.** The drag coefficient for spheres such as sporting balls is known to vary as a +function of the Reynolds number Re, a dimensionless number that gives a measure of the +ratio of inertial forces to viscous forces: + +$Re = \frac{\rho V D}{\mu}$ +![Reynolds](./equations/reynolds.png) + +where ρ= the fluid’s density (kg/m3), V= its velocity (m/s), D= diameter (m), and μ= +dynamic viscosity (N⋅s/m2). The +following table provides values for a smooth spherical ball: + +|Re (x1e-4) | C_D | +|---|---| +|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| + +a. Create a function `sphere_drag.m` that outputs the drag coefficient based on the given +table and an input Reynolds number using a spline interpolation of either linear +('linear'), +piecewise cubic ('pchip'), or cubic ('cubic'): + +`[Cd_out]=sphere_drag(Re_in,spline_type)` + +b. Use the following physical constants to plot the drag force vs velocity for a baseball: +ρ= 1.3 kg/m3, V= 4 - 40 (m/s), D= 23.5 cm, and μ=1.78e-5 Pa-s. Plot all three +interpolation methods on a single plot. Show the plot in your README. + +**6\.** Evaluate the integral of the following function: + +$\int_{2}^{3}f(x) = \int_{2}^{3} 1/6x^3 + 1/2x^2+x dx$ + +a. analytically + +b. with 1-point Gauss quadrature + +c. with 2-point Gauss quadrature + +d. with 3-point Gauss quadrature + +e. include the results in a table in your README + +|method|value|error| +|---|---|---| +|analytical| ... | 0%| +|1 Gauss point | ... | ..% | +|... | ... | ... | + + diff --git a/HW5/equations/buckle.png b/HW5/equations/buckle.png new file mode 100644 index 0000000..0031a09 Binary files /dev/null and b/HW5/equations/buckle.png differ diff --git a/HW5/equations/buckle.tex b/HW5/equations/buckle.tex new file mode 100644 index 0000000..37e36bd --- /dev/null +++ b/HW5/equations/buckle.tex @@ -0,0 +1,2 @@ +P_{cr}=\frac{\pi^2 EI}{4L^2} + diff --git a/HW5/equations/fx.png b/HW5/equations/fx.png new file mode 100644 index 0000000..d4d0b17 Binary files /dev/null and b/HW5/equations/fx.png differ diff --git a/HW5/equations/fx.tex b/HW5/equations/fx.tex new file mode 100644 index 0000000..3f14aac --- /dev/null +++ b/HW5/equations/fx.tex @@ -0,0 +1,2 @@ +\int_{2}^{3}f(x) = \int_{2}^{3} 1/6x^3 + 1/2x^2+x dx + diff --git a/HW5/equations/reynolds.png b/HW5/equations/reynolds.png new file mode 100644 index 0000000..2b53d3c Binary files /dev/null and b/HW5/equations/reynolds.png differ diff --git a/HW5/equations/reynolds.tex b/HW5/equations/reynolds.tex new file mode 100644 index 0000000..7414d0e --- /dev/null +++ b/HW5/equations/reynolds.tex @@ -0,0 +1,2 @@ +Re = \frac{\rho V D}{\mu} + diff --git a/HW5/problem_2_data.m b/HW5/problem_2_data.m new file mode 100644 index 0000000..558debc --- /dev/null +++ b/HW5/problem_2_data.m @@ -0,0 +1,19 @@ + +% part a +xa=[1 2 3 4 5]'; +ya=[2.2 2.8 3.6 4.5 5.5]'; + +% 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]'; + +% 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]'; + +% 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]'; diff --git a/extra_credit/README.md b/extra_credit/README.md index 0ac65ae..1234d6d 100644 --- a/extra_credit/README.md +++ b/extra_credit/README.md @@ -21,3 +21,45 @@ columns with your netid on each row as such, |rcc02007 | 1 | 30 | |rcc02007 | ...| ... | |rcc02007 | ...| ... | + +# Extra Credit Assignment \#3 +## Due 11/17 by 11:59pm + +**Nonlinear Regression - Logistic Regression** + +![logistic regression of Challenger O-ring failure](http://www.stat.ufl.edu/~winner/cases/challenger.ppt) + +Use the Temperature and failure data from the Challenger O-rings +[challenger_oring.csv](./challenger_oring.csv). Your independent variable is temperature and your dependent +variable is failure (1=fail, 0=pass). Create a function called `cost_logistic.m` that +takes the vector `a`, and independent variable `x` and dependent variable `y`. Use the +function, $\sigma(t)=\frac{1}{1+e^{-t}}$ where $t=a_{0}+a_{1}x$. Use the cost function, + +$J(a_{0},a_{1})=1/m\sum_{i=1}^{n}\left[-y_{i}\log(\sigma(t_{i}))-(1-y_{i})\log((1-\sigma(t_{i})))\right]$ + +and gradient + +$\frac{\partial J}{\partial a_{i}}= +1/m\sum_{k=1}^{N}\left(\sigma(t_{k})-y_{k}\right)x_{k}^{i}$ + +where $x_{k}^{i} is the k-th value of temperature raised to the i-th power (0, and 1) + +a. edit `cost_logistic.m` so that the output is `[J,grad]` or [cost, gradient] + +b. use the following code to solve for a0 and a1 + +```matlab +% Set options for fminunc +options = optimset('GradObj', 'on', 'MaxIter', 400); +% Run fminunc to obtain the optimal theta +% This function will return theta and the cost +[theta, cost] = ... +fminunc(@(a)(costFunction(a, x, y)), initial_a, options); +``` + +c. plot the data and the best-fit logistic regression model + +```matlab +plot(x,y, x, sigma(a(1)+a(2)*x)) +``` +