%% Exercise 4_5; Book Problem
%% Carol Lucas
%% Calls LcfullE5fun for ode simulation
mess='Program pauses after every figure.  Hit return to continue.'
VRAY=[-91:5:60];
caout=3;
gca=4;
init=[1;0];
figure(1);clf;
for k=1:length(VRAY)
V=VRAY(k)
[t,x]=ode45('LcfullE5fun',[0:400],init,[],V);
Ica=gca*caout*x(:,2)*V./(1-exp(V/13.4));
Imax(k)=min(Ica);
k1p=0.78/(1+exp(-(3+V)/10));
k1m=0.78-k1p;
k2pcad=-7.56e-4*caout*V/(1-exp(V/13.4));
k2m=0.002;

%% Calculate SS from Matrix - Compare to simulation
a=[-k1p k1m;k1p-k2m -k1m-k2pcad-k2m];
b=[0;k2m];
xss(:,k)=-inv(a)*b;
Iss(k)=gca*caout*xss(2,k)*V/(1-exp(V/13.4));
Isssim(k)=Ica(300);
plot(t,x,t,Ica,[300],Iss(k),'ro');
xlabel('time');ylabel('Open, Closed, Current, Calculated SS');
pause;
end;

figure(2),clf;plot(VRAY,xss);
xlabel('V Value');
ylabel('calculated ss values');
title('Exercise 4_5, V vs ss values of model');
pause;
figure(3),clf;plot(VRAY,Iss,VRAY,Imax,VRAY,Isssim);
xlabel('V Value');
ylabel('Ica minimum and equillibrium (calculated and simulated)');
title('Exercise 4-5, Voltage Dependent Model');