ME3255_FInalProject

I noticed that I was never assigned to a group so I did this project solo

Part A

Problem 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); The output w should be a vector, but the solution represents a 2D data set w(x,y).

Approach First off for this part I created a matrix using analysis of the interior nodes of a 5x5 matrix. After that I created a vector for the Y direction, after which it is converted into a vector for the displacement of the nodes. The vector is then transformed to be the surface of the new matrix.

function [w] = membrane_solution3(T,P)


    od = ones(8,1);
    od(3:3:end) = 0;
    a = -4*diag(ones(9,1))+diag(ones(9-3,1),3)+diag(ones(9-3,1),-3)+diag(od,1)+diag(od,-1);


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

grid on
    [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]);
    surf(x,y,z)
    title('Membrane Displacement')
    xlabel('X (\muM)')
    ylabel('Y (\muM)')
    zlabel('Displacement (\muM)')
    colormap jet
    shading interp
    disp(w)

end

Part B

Problem 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.

membrane_solution3(0.006,0.001)

Part C

Problem 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); The output w should be a vector, but the solution represents a 2D data set w(x,y).

Approach This part of the assignment was handled similarly to part A, as the only thing that changed from the two scripts was the ability of the user to input their own dimensions, which was done by altering the indexing in the function.

function [w] = membrane_solution3(T,P)


    od = ones(8,1);
    od(3:3:end) = 0;
    a = -4*diag(ones(9,1))+diag(ones(9-3,1),3)+diag(ones(9-3,1),-3)+diag(od,1)+diag(od,-1);


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

grid on
    [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]);
    surf(x,y,z)
    title('Membrane Displacement')
    xlabel('X (\muM)')
    ylabel('Y (\muM)')
    zlabel('Displacement (\muM)')
    colormap jet
    shading interp
    disp(w)

end

Part D

Problem 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. Include the graph in your README.

membrane_solution(0.006,0.001,10)

Part E

Problem 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 The first thing I did was pull the solution from membrane_solution.m, after which the average displacement of each node was calculated. This was done through taking the displacement at each corner, after which it was all averaged. The same was done for the left and right points (boundary conditions).

function [pw_se,w]=SE_diff(T,P,n)

E = 1e6;
v = .31;
t = .0003;
h = 10/(n+1);
w = membrane_solution(T,P,n);
z = zeros(n+2);
z(2:end-1,2:end-1) = reshape(w,[n n]);
nt = n + 1;
wn = zeros(nt);
for i = 1:nt
    for j = 1:nt
        wn(i,j) = mean([z(i,j),z(i+1,j),z(i,j+1),z(i+1,j+1)]);
    end
end
pw = sum(sum(wn.*h^2.*P));
dwdx = zeros(nt);
dwdy = zeros(nt);
for i = 1:nt
    for j = 1:nt
        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
se = E*t*h^2/(2*(1-v^2))*sum(sum(0.25.*dwdx.^4+.25.*dwdy.^4+0.5.*(dwdx.*dwdy).^2));
pw_se = pw-se;

Part F

Problem 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 Using previous knowledge of bisecting functions, as well as the outputs of the SE_diff function, the script ran every single iteration for the different amounts of nodes. The script then saves the values as a vector T, from which it is able to compare and output error from rel_error.

function [root,fx,ea,iter]=bisect(func,xl,xu,es,maxit,varargin)
if nargin<3,error('at least 3 input arguments required'),end
test = func(xl,varargin{:})*func(xu,varargin{:});
if nargin<4||isempty(es), es=0.0001;end
if nargin<5||isempty(maxit), maxit=50;end
iter = 0; xr = xl; ea = 100;
while (1)
  xr_old = xr;
  xr = (xl + xu)/2;
  iter = iter + 1;
  if xr ~= 0,ea = abs((xr - xr_old)/xr) * 100;end
  test = func(xl,varargin{:})*func(xr,varargin{:});
  if test < 0
    xu = xr;
  elseif test > 0
    xl = xr;
  else
    ea = 0;
  end
  if ea <= es || iter >= maxit,break,end
end
root = xr; fx = func(xr, varargin{:});

function [T,ea] = tension_sol(P,n)
y =@(T) SE_diff(T,P,n);
[T,fx,ea,iter]=bisect(y,.01,1);

n=[3,20:5:40];
P=0.001;
T = zeros(1,length(n));
ea = zeros(1,length(n));

 for i = 1:length(n)
[T(i), ea(i)] = tension_sol(P,n(i));
 end

number of nodes	Tension (uN/um)	rel. error
3	0.049	n/a
20	0.0599	20.5%
25	0.0601	0.27%
30	0.0602	0.16%
35	0.0602	0.09%
40	0.0603	0.05%

Part G

Problem 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. Plot the data and best-fit curve in your README.

Approach Utilizing the previous tension_sol results, the script runs all the iterations of pressure and tension. After this it also calculates the maximum deflection of every single iteration. The results are plotted as a linear regression, and plotted on the chart pictured below.

No = linspace(.001,.01,10);
n = 20;

T = zeros(1,length(No));
wmax = zeros(1,length(No));

for i = 1:length(No)
    T(i) = tension_sol(No(i),n);
    w = membrane_solution(T(i),No(i),n);
    wmax(i) = max(w);
end

clf
x = wmax';
y = No';
Z=x.^3;
a=Z\y;
x_fcn=linspace(min(x),max(x));
plot(x,y,'o',x_fcn,a*x_fcn.^3)
title('Pressure vs Maximum Deflection')
xlabel('Maximum Deflection (\muM)')
ylabel('Pressure (MPa)')

rcc02007/ME3255_FInalProject

About

Resources

Stars

Watchers

Forks

Releases

Languages