Please write all the MATLAB code, and if you can try to comment how it works or explain it below.

Problem statement

The role of the cochlea (inner ear) is to convert sounds into electrical signals sent to the brain. The cochlea consists of two fluid-filled ducts separated by a structure called basilar membrane. Because the basilar membrane is much stiffer in the zone close to the middle ear (called basal end) than at the other end (called apical end), the basilar membrane responds primarily to high frequency (high pitch) sounds at the basal end and to low frequency sounds (low pitch) at the apical end. This allows us to distinguish between different sound frequencies. In this computer project you will code a simple model of the cochlea based by a one-dimensional wave equation:

d2P 2i P(x)=0 (1)

where P is the fluid pressure, is the fluid density, i is the unit imaginary number (1i in Matlab), H is the height of the cochlear ducts, Zbm is the impedance and is the radian frequency ( = 2f where f is the sound frequency in Hz). Zbm is given by the following function of x:

Zbm(x) = iM0 + K0M0exp((1/2)x) + K0/i exp(x) (2)

where is the radian frequency, M0 is the mass, K0 is the stiffness at x = 0 and is the damping constant and is a constant parameter. The boundary conditions are:

dP/dX (0) = 1 (3)

P (L) = 0

Once the pressure is calculated, the basilar membrane velocity can be computed using:

vbm(x) = 2P(x)/Zbm(x) (5)

The parameters for the human cochlea are listed in the table below.

Parameter Value L 0.035m H 0.001m K0 1 1010N/m3 M0 0.5 kg/m2 0.03 1000kg/m3

3 102m1

For this project, all Matlab codes should be in the file named CP5.m where. Before writing any code, look at the example file named exampleFiniteDifferenceRod.m.

% -------------------------------------------------------------------------------------------------------

function exampleFiniteDifferenceRod close all clear all

% model parameters L=10; %length of the rod T1=40; % temperature at the left end of the rod T2=400; % temperatuer at the right end of the rod Ta=20; %ambient temperature h=0.01; %heat tranfer coefficient

n=10000; %number of interior nodes deltaX=L/(n+1);

% full matrices ------------------------------------------------------- % A=zeros(n,n); %create matrix full of zeros % F=zeros(n,1); % for i=1:n % C=h*deltaX^2+2; % D=h*deltaX^2*Ta; % A(i,i)=-C; % % if i1 % A(i,i-1)=1; %term before the diagonal % end % % %form RHS vector % F(i)=-D; % if i==1 % F(1)=F(1)-T1; % elseif i==n % F(n)=F(n)-T2; % end % end % T=A^-1*F; % T=A\F; %solve to find vector T (much faster than A^-1*F for large n) A=sparse(n,n); %create sparse matrix (makes computation much faster, % and requires much less memory) F=sparse(n,1); for i=1:n C=h*deltaX^2+2; D=h*deltaX^2*Ta; A(i,i)=-C; if i1 A(i,i-1)=1; %term before the diagonal end %form RHS vector F(i)=-D; if i==1 F(1)=F(1)-T1; elseif i==n F(n)=F(n)-T2; end end T=A\F; %solve to find vector T % Textended=[T1;T;T2]; %add the 2 BCs to the vector T x=linspace(0,L,n+2);<><>

plot(x,Textended,'b-') hold on xlabel('x');ylabel('T') %using bvp4c function ----------------------------------------------- solinit=bvpinit(x,zeros(2,1)); %need to initialize the solution sol=bvp4c(@dXdx,@BC,solinit); X=deval(sol,x); %evaluate the solution on the mesh T=X(1,:); % T is the 1st row of the matrix X plot(x,T,'g--','LineWidth',2) legend('Finite difference code','bvp4c') % end

% function for the ODE for heat transfer in rod usng bvp4c % inputs: % - x: position (it is not used for this example; but it needs to be % included, or MATLAB bvp4c function will return an error) % - vector X=[T,Z]^T where Z=dT/dx % outputs: vector dT=[dx1/dx,dx2/dx]^T function [dX] = dXdx(x,X) Ta=20; %ambient temperature h=0.01; %heat tranfer coefficient T=X(1,:); Z=X(2,:); dT=Z; dZ=h*(T-Ta); dX=[dT;dZ]; end

% function for the BCs for heat transfer example usng bvp4c %inputs: % - X_0=[T(0),Z(0)] % - X_L=[T(L),Z(L)] % outputs % - res=residual vector (need to be equal to 0 to satisfy the BCs) function [res] = BC(X_0,X_L) T1=40; % temperature at the left end of the rod T2=400; % temperature at the right end of the rod T_0=X_0(1); T_L=X_L(1); res(1)=T_0-T1; % to enforce T(0)=T1 res(2)=T_L-T2; % to enforce T(L)=T2 end

% -------------------------------------------------------------------------------------------------------

Tasks 1. Using the centered difference formula for d2P (x ) and the centered difference formula for

dx2 i dP (x0), transform this boundary value problem for P(x) into a linear system of equations of

dx the form [A]{P } = {F } where [A] is a matrix, {P }=[P (x0 ),P (x1 ),,...,P (xn )]T (where x0 = 0 and xn+1 = L) and {F } is a column vector.

In the function CP5, find the vector of the fluid pressure, {P} using the finite difference method for f =6000 Hz and f =500 Hz. Use n = 1001. Plot the absolute value of the pressure as a function of x in Figure 1.

Compute the basilar membrane velocity using Eq. 5 for f=6000Hz and f=500 Hz. Plot on a separate graph (Figure 2) the absolute value of the basilar membrane velocity a function of x.

Transform the 2nd order ordinary differential equation for P into a system of 1st order ordinary differential equations. Write the system in the form: dX = F (x, X ), where X =

dx

[dP ,P]T and F(x,X) is a column vector of length 2 that is a function of x and X. dx

Find the fluid pressure using Matlabs bpv4c function for f=6000Hz and f=500 Hz. You will

need to write two additional functions: one function named dXdx that returns F(x,X), and

a function named BC used to enforce the boundary conditions dP (0) = 1 and P(L) = 0. Plot dx

in Figure 1 the absolute value of the pressure that you obtain using bvp4c.

Plot in Figure 2 the absolute value of the basilar membrane velocity (for f=6000Hz and

f=500 Hz) that you obtain using Eq. 5 applied to the fluid pressure obtained using bvp4c.

Make sure to include axis labels, titles, and legends on your graphs. And please comment your codes.

Problem statement The role of the cochlea (inner ear) is to convert sounds into electrical signals sent to the brain. The cochlea consists of two fluid-filled ducts separated by a structure called basilar membrane. Because the basilar membrane is much stiffer in the zone close to the middle ear (called basal end) than at the other end (called apical end), the basilar membrane responds primarily to high frequency (high pitch) sounds at the basal end and to low frequency sounds (low pitch) at the apical end. This allows us to distinguish between different sound frequencies. In this computer project you will code a simple model of the cochlea based by a one-dimensional wave equation: d . (1) 2ipw dx2 HZP() = 0 where P is the fluid pressure, p is the fluid density, i is the unit imaginary number (li in Matlab), H is the height of the cochlear ducts, Zbm is the impedance and w is the radian frequency (w = 2rf where f is the sound frequency in Hz). Zbm is given by the following function of r: Korp(-OLI Zbm (x) = iw Mo + dV K, M,exp( -ax) + "exp(-ax) (2) 2 iw where w is the radian frequency, Mo is the mass, Ko is the stiffness at r = 0 and d is the damping constant and a is a constant parameter. The boundary conditions are: dP dr (0) = 1 P(L) = 0 Once the pressure is calculated, the basilar membrane velocity can be computed using: -2P() * Vbm (T) = Zbm () The parameters for the human cochlea are listed in the table below. H Parameter Value 0.035m H 0.001m | K. 1x 1010 N/m3 0.5 kg/m2 0.03 1000kg/m3 3 x 102m-1 M, 2 Problem statement The role of the cochlea (inner ear) is to convert sounds into electrical signals sent to the brain. The cochlea consists of two fluid-filled ducts separated by a structure called basilar membrane. Because the basilar membrane is much stiffer in the zone close to the middle ear (called basal end) than at the other end (called apical end), the basilar membrane responds primarily to high frequency (high pitch) sounds at the basal end and to low frequency sounds (low pitch) at the apical end. This allows us to distinguish between different sound frequencies. In this computer project you will code a simple model of the cochlea based by a one-dimensional wave equation: d . (1) 2ipw dx2 HZP() = 0 where P is the fluid pressure, p is the fluid density, i is the unit imaginary number (li in Matlab), H is the height of the cochlear ducts, Zbm is the impedance and w is the radian frequency (w = 2rf where f is the sound frequency in Hz). Zbm is given by the following function of r: Korp(-OLI Zbm (x) = iw Mo + dV K, M,exp( -ax) + "exp(-ax) (2) 2 iw where w is the radian frequency, Mo is the mass, Ko is the stiffness at r = 0 and d is the damping constant and a is a constant parameter. The boundary conditions are: dP dr (0) = 1 P(L) = 0 Once the pressure is calculated, the basilar membrane velocity can be computed using: -2P() * Vbm (T) = Zbm () The parameters for the human cochlea are listed in the table below. H Parameter Value 0.035m H 0.001m | K. 1x 1010 N/m3 0.5 kg/m2 0.03 1000kg/m3 3 x 102m-1 M, 2