function [delU]=Quasi_1D_exit_SI_bc(delt,delx,Cv,a,To,po,gamma,pI,pI_1,UI,UI_1,sI,sI_1)

gamma_fac=(gamma-1)/(gamma+1);   %Used coefficient
rhoI=UI(1);              rhoI_1=UI_1(1);
uI=UI(2)/UI(1);          uI_1=UI_1(2)/UI(1);
eI=UI(3);                eI_1=UI_1(3);

%Get dependent variables
%Temperature
TI=To*(1-gamma_fac*uI^2/a^2);
TI_1=To*(1-gamma_fac*uI_1^2/a^2);
%Pressure
%pI=po*(1+gamma_fac*uI^2/a^2)^(gamma/(gamma-1));
%pI_1=po*(1+gamma_fac*uI_1^2/a^2)^(gamma/(gamma-1))
pI=(gamma-1)*rhoI*(eI/rhoI-.5*uI^2);
pI_1=(gamma-1)*rhoI_1*(eI_1/rhoI_1-.5*uI_1^2)

%Speed of Sound
cI=sqrt(gamma*pI/rhoI);
cI_1=sqrt(gamma*pI_1/rhoI_1);

%Average the variables
u=(uI+uI_1)/2;          p=(pI+pI_1)/2;
c=(cI+cI_1)/2;          rho=(rhoI+rhoI_1)/2;

%Eigenvalues
lam_1=u*delt/delx;
lam_2=(u+c)*delt/delx;
lam_3=(u-c)*delt/delx;
%Residuals
R1=((-lam_1)/(1+lam_1))*(rhoI-rhoI_1-(1/c^2)*(pI-pI_1));
R2=((-lam_2)/(1+lam_2))*(pI-pI_1+rho*c*(uI-uI_1))-...
   gamma*delt*p*u/(1+lam_2)*((sI-sI_1)/(sI*delx));
R3=((-lam_3)/(1+lam_3))*(pI-pI_1-rho*c*(uI-uI_1))-...
   gamma*delt*p*u/(1+lam_3)*((sI-sI_1)/(sI*delx));
%Use Residuals to obtain differences
if u/c > 1
    delp=(R2+R3)/2;
else
    delp=0;
end

delrho=R1+delp/c^2;
delu=(R2-delp)/(rho*c);

pI=pI+delp;                 
uI=uI+delu;
rhoI=rhoI+delrho;

TI=pI/((gamma-1)*Cv*rhoI);
eI=rhoI*(Cv*TI+.5*uI^2);

delU=[rhoI;rhoI*uI;eI];

