%% 6.581 PS5
%% Qsn4
%% lambda2.m

function [phi, gradPhi] = lambda2(P)

%% DATA
%%
%% Observed Fixed points for two different trials
%%
obs1 = [33.71; 0.11; 1.57; 0.01];
obs2 = [0.05 ; 3.74; 0.02; 0.479];  


%% Find the fixed points (x1, x2)
%% Hint: this is what we did on last weeks homework
x1 = ; %CHANGE ME
x2 = ; %CHANGE ME

%%
%% Compute the Squared Error to the Data
%%
c   = ; %CHANGE ME
phi = ; %CHANGE ME

%%
%% Compute the Gradient of Phi with P.  This part will be several steps
%% you may want to compute the matrix quantites we will need here in
%% a seperate function.
%%
gradPhi = ; %CHANGE ME





%%=========================================================================
%% The ODE model of the Lambda Switch
%% This is implemented for you but you may want to copy and paste this code
%% as a starting point for other helper methods.
%%
%% These are the constants that we used in PS4
%% 
%% p0 = [ 1.0; %DNA_TOTAL
%%        0.5; %k_Cro2B_mCro
%%        1.0; %k_mCro_0
%%        1.0; %k_mCro_Cro
%%        1.0; %k_Cro_0
%%        10;  %k_Cro_Cro2
%%        1.0; %k_Cro2_Cro
%%        10;  %k_Cro2_Cro2B
%%        1.0; %k_Cro2B_Cro2
%%        0.5; %k_CI2B_mCI
%%        1.0; %k_mCI_0
%%        1.0; %k_mCI_CI
%%        1.0; %k_CI_0
%%        10;  %k_CI_CI2
%%        1.0; %k_CI2_CI
%%        10;  %k_CI2_CI2B
%%        1    %k_CI2B_CI2
%%        ];
%%=========================================================================
function dxdt=F(x, P)
%% Indexes into x
mCro  = 1;
Cro   = 2;
Cro2  = 3;
Cro2B = 4;
mCI   = 5;
CI    = 6;
CI2   = 7;
CI2B  = 8;

%% Indexes into P
DNA_tot      = 1;
k_Cro2B_mCro = 2;
k_mCro_0     = 3;
k_mCro_Cro   = 4;
k_Cro_0      = 5;
k_Cro_Cro2   = 6;
k_Cro2_Cro   = 7;
k_Cro2_Cro2B = 8;
k_Cro2B_Cro2 = 9;
k_CI2B_mCI   = 10;
k_mCI_0      = 11;
k_mCI_CI     = 12;
k_CI_0       = 13;
k_CI_CI2     = 14;
k_CI2_CI     = 15;
k_CI2_CI2B   = 16;
k_CI2B_CI2   = 17;

%% Don't forget to differentiate this when you are computing the Jacobian
%% and the gradients
DNA = P(DNA_tot) - x(Cro2B)- x(CI2B);

%% The Model
dxdt = [ ...
P(k_Cro2B_mCro)  * x(Cro2B)    -   P(k_mCro_0)     * x(mCro);
P(k_mCro_Cro)    * x(mCro)     + 2*P(k_Cro2_Cro)   * x(Cro2)    - 2*P(k_Cro_Cro2) * x(Cro)^2 - P(k_Cro_0)      * x(Cro);
P(k_Cro_Cro2)    * x(Cro)^2    +   P(k_Cro2B_Cro2) * x(Cro2B)   -   P(k_Cro2_Cro) * x(Cro2)  - P(k_Cro2_Cro2B) * x(Cro2)*DNA;
P(k_Cro2_Cro2B)  * x(Cro2)*DNA -   P(k_Cro2B_Cro2) * x(Cro2B);
P(k_CI2B_mCI)    * x(CI2B)     -   P(k_mCI_0)      * x(mCI);
P(k_mCI_CI)      * x(mCI)      + 2*P(k_CI2_CI)     * x(CI2)     - 2*P(k_CI_CI2)   * x(CI)^2  - P(k_CI_0)       * x(CI);
P(k_CI_CI2)      * x(CI)^2     +   P(k_CI2B_CI2)   * x(CI2B)    -   P(k_CI2_CI)   * x(CI2)   - P(k_CI2_CI2B)   * x(CI2)*DNA;
P(k_CI2_CI2B)    * x(CI2)*DNA  -   P(k_CI2B_CI2)   * x(CI2B)...
];