# ayr13001 / ME_3255_Final_Project

ME3255 Final Project (Group 27)
MATLAB

## Latest commit Fetching latest commit…
Cannot retrieve the latest commit at this time.

## Files

Type Name Latest commit message Commit time
Failed to load latest commit information. Figures .DS_Store Grade PartFscript.m README.md SE_diff.m bisect.m membrane_solution.m membrane_solution3.m partf.m partg.m

# ME 3255 Final Project

## Part A

#### Problem Statement

Create a central finite difference approximation of the gradient with 3-by-3 interior nodes of w for the given membrane solution in terms of P and T. `[w]=membrane_solution3(T,P);`

#### Approach

```function [w] = membrane_solution3(T,P)
% Central finite difference approximation of gradient
% with 3x3 (um^2) interior nodes in terms of P & T.
% Input:
% T = tension per unit length (uN/um)
% P = pressure (MPa)
% Output:
% w = displacement vector for interior nodes

od = ones(8,1);
od(3:3:end) = 0;
k = -4 * diag(ones((3^2),1)) + diag(ones((3^2)-3,1),3) + diag(ones((3^2)-3,1),-3) + diag(od,1) + diag(od,-1);

y = -(10/4)^2*(P/T)*ones(9,1);
w = k\y;

% Find displacement vector w
% which represents 2D data set w(x,y)
[x,y] = meshgrid(0:10/4:10,0:10/4:10);
z = zeros(size(x));
z(2:end-1,2:end-1) = reshape(w,[3 3]);

% Plot gradient of membrane displacement
surf(x,y,z)
zlabel('Displacement [um] (micrometers)')
end```

## Part B

#### Problem Statement

Solve for w given a pressure, P=0.001 MPa and tension, T=0.006 uN/um. Plot the result with `surf(X,Y,W)` where X, Y, and W are the x-, y-, and z-coordinates of each point on the membrane from 0-10um.

#### Approach

T = 0.006 uN/um P = 0.001 MPa

`[w] = membrane_solution3(0.006,0.001);` ## Part C

#### Problem Statement

Create a general central finite difference approximation of the gradient with n-by-n interior nodes of w for the given membrane solution in terms of P and T. `[w]=membrane_solution(T,P,n);`

#### Approach

```function [w] = membrane_solution(T,P,n)
% General Central finite difference approximation of
% gradient with n-by-n interior nodes in terms of P & T.
% Input:
% T = tension per unit length (uN/um)
% P = pressure (MPa)
% n = # of interior node rows/columns
% Output:
% w = displacement vector for interior nodes

od = ones(n^2-1,1);
od(n:n:end) = 0;
k = -4 * diag(ones(n^2,1)) + diag(ones((n^2)-n,1),n) + diag(ones((n^2)-n,1),-n) + diag(od,1) + diag(od,-1);

y = -(10/(n+1))^2*(P/T)*ones(n^2,1);
w = k\y;

% Find displacement vector w
% which represents 2D data set w(x,y)
[x,y] = meshgrid(0:10/(n+1):10,0:10/(n+1):10);
z = zeros(size(x));
z(2:end-1,2:end-1) = reshape(w,[n n]);

% Plot gradient of membrane displacement
surf(x,y,z)
zlabel('Displacement [um] (micrometers)')
end```

## Part D

#### Problem Statement

Solve for w given a pressure, P=0.001 MPa and tension, T=0.006 uN/um with 10 interior nodes. Plot the result with `surf(X,Y,W)` where X, Y, and W are the x-, y-, and z-coordinates of each point on the membrane from 0-10um.

#### Approach

• T = 0.006 uN/um
• P = 0.001 MPa
• n = 10 nodes
`[w] = membrane_solution(0.006,0.001,10)` ## Part E

#### Problem Statement

Create a function `SE_diff` that calculates the difference in strain energy (right hand side Eq.4) and work done by pressure (left hand side Eq. 4) for n-by-n elements. `[pw_se,w]=SE_diff(T,P,n)` Use the solution from part c to calculate w, then do a numerical integral over the elements to calculate work done and strain energy.

#### Approach

```function [pw_se,w] = SE_diff(T,P,n)
% function that calculates the difference between
% strain energy and work done by pressure on the membrane.
% Input:
% T = tension per unit length (uN/um)
% P = pressure (MPa)
% n = # of interior node rows/columns
% Output:
% pw_se = Absolute value of difference between strain energy and work done by pressure
% w = displacement vector for interior nodes

E = 1e6;        % 1 TPa ~= 10^6 MPa
t = 3*10^-4;    % thickness [um]
h = 10/(n+1);   % height [um]
v = 0.31;       % Poisson's Ratio

% Displacement vector w found using Part C
w = membrane_solution(T,P,n);
z = zeros(n + 2);
z(2:end-1,2:end-1) = reshape(w,[n n]);

% Calculate average displacement, wavg, for each element by taking the displacement at each
% corner and then average the found values.
num = n + 1;
wavg = zeros(num);
for i = 1:num
for j = 1:num
wavg(i,j) = mean([z(i,j),z(i+1,j),z(i,j+1),z(i+1,j+1)]);
end
end

% final work done by pressure
pw = sum(sum(wavg.*h^2.*P))

% to find= change in displacement, find the change in displacement on
% the x-axis, dwdx, and the change in displacement on the y-axis, dwdy, and
% average the found values.
dwdx = zeros(num);
dwdy = zeros(num);
for i = 1:num
for j = 1:num
dwdx(i,j) = mean([z(i+1,j)-z(i,j),z(i+1,j+1)-z(i,j+1)]);
dwdy(i,j) = mean([z(i,j+1)-z(i,j),z(i+1,j+1)-z(i+1,j)]);
end
end

% Using dwdx and dwdy, calculate the strain energy, se.
se = (E*t*h^2)/(2*(1-v^2)) * sum(sum((1/4).*dwdx.^4+(1/4).*dwdy.^4+(1/4).*(dwdx.*dwdy).^2));

% Final value of difference between strain energy and work done by pressure, pw_se.
pw_se = pw - se;```

## Part F

#### Problem Statement

Use a root-finding method to calculate the tension in the membrane given a pressure, P=0.001 MPa, and n=[20:5:40] interior nodes. Show that the error in tension is decreasing with a table

#### Approach

Root-finding Method used to find tension in the Membrane, given:

• P = 0.001 MPa
• n = [20:5:40]
```% Part F script
% Take the results of 'SE_diff.m' and use the 'bisect.m'
% root-finding method to
% solve for the
% tension T,
% given a pressure P and
% a range of interior nodes [n]

[pw_se,w] = SE_diff(T,P,n)
[root,fx,ea,iter] = bisect(@(T) SE_diff(T, 0.001,40),0.001,1,0.1)ff(T,0.001,20),0.001,1,0.1)
% Run a total of 5 times, varying the value of n by 5
% from 20 to 40```
##### Relative Error
number of nodes Tension(uN/um) rel. error
3 0.0489 n/a
20 0.0599 22.3%
25 0.0601 0.27%
30 0.0602 0.17%
35 0.0603 0.09%
40 0.0603 0.06%

## Part G

#### Problem Statement

Plot the Pressure vs maximum deflection (P (y-axis) vs max(w) (x-axis)) for P = linspace(0.001,0.01,10). Use a root-finding method to determine tension, T, at each pressure. Use a cubic best-fit to find A, where, P(x)=A*dw^3. State how many interior nodes were used for the graph.

#### Approach

```% Create plot for wmax v pressure
clear
% range of pressures to be used to find tension T
P = linspace(0.001,0.01,10);
% find the tension using the 'bisect.m' method
for i = 1:length(P)
func = @(T) SE_diff(T,P(i),10);
[root(i),fx,ea,iter] = bisect(func,0.001,1,.1);
end

% Calculate each w using every root and pressure value
for i = 1:length(root)
w = membrane_solution(root(i),P(i),10);
w1(:,i) = w; % displays w values as column
end

% Find final wmax by taking the maximum w value from each
% column in the w vector
wmax = max(w1);
coefficients = polyfit(wmax,P,3);
Y = polyval(coefficients,wmax);

% plot the w_max vs. the Pressure and include a cubic best fit curve.
plot(wmax,P,wmax,Y,'or')
xlabel('Max Deflection (um)')
ylabel('Pressure (MPa)')
title('Pressure vs. Maximum Deflection')```
You can’t perform that action at this time.