%
% MorrisLecar.m
%
% MATLAB file to be used with MATCONT.
%
% This file was generated by the program VFGEN (Version:2.4.0)
% Generated on 10-Jul-2008 at 14:18
%
%
function out = MorrisLecar
    out{1} = @MorrisLecar_init;
    out{2} = @MorrisLecar_vf;
    out{3} = @MorrisLecar_jac;
    out{4} = @MorrisLecar_jacp;
    out{5} = @MorrisLecar_hess;
    out{6} = @MorrisLecar_hessp;
    out{7} = @MorrisLecar_der3;
    out{8} = [];
    out{9} = [];

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

function [tspan, y0, options] = MorrisLecar_init

    tspan = [0; 10];

    % Get the constants and default parameters, because the default
    % initial conditions can depend on them.
    gca = 5.5;
    gk = 8.0;
    gl = 2.0;
    vca = 115.0;
    vk = -84.0;
    vl = -55.0;
    c = 20.0;
    phi = 0.22;
    ic = 90.0;
    v1 = -1.2;
    v2 = 18.0;
    v3 = 2.0;
    v4 = 30.0;
    % Set the default initial point.
    y0 = zeros(2,1);
    y0(1) = 0;
    y0(2) = 0;

    handles = feval(@MorrisLecar);
    options = odeset('Jacobian',handles(3),'JacobianP',handles(4), ...
                     'Hessians',handles(5),'HessiansP',handles(6), ...
                     'Der3',handles(7) );

%--------------------------------------------------
%
% The vector field
%
function vf_ = MorrisLecar_vf(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    minf = 1/2+1/2*tanh(-(v1-v)*v2^(-1));
    winf = 1/2+1/2*tanh(-(v3-v)*v4^(-1));
    tauw = cosh(-1/2*(v3-v)*v4^(-1))^(-1);

    vf_ = zeros(2,1);
    vf_(1) = -c^(-1)*(w*(v-vk)*gk+gl*(v-vl)-gca*minf*(vca-v)-ic);
    vf_(2) = -tauw^(-1)*(w-winf)*phi;

%--------------------------------------------------
%
% The Jacobian of the vector field with respect to the variables
%
function jac_ = MorrisLecar_jac(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    jac_ = zeros(2,2);
    jac_(1,1) = -1/2*c^(-1)*(2*gl+2*w*gk+gca*(1+tanh(-(v1-v)*v2^(-1)))+gca*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1));
    jac_(1,2) = -c^(-1)*(v-vk)*gk;
    jac_(2,1) = -1/2*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-1)*cosh(-1/2*(v3-v)*v4^(-1))*phi-1/4*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-1)*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi;
    jac_(2,2) = -cosh(-1/2*(v3-v)*v4^(-1))*phi;

%--------------------------------------------------
%
% The Jacobian of the vector field with respect to the parameters
%
function jacp_ = MorrisLecar_jacp(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    jacp_ = zeros(2,13);
    jacp_(1,1) = 1/2*c^(-1)*(1+tanh(-(v1-v)*v2^(-1)))*(vca-v);
    jacp_(1,2) = -c^(-1)*w*(v-vk);
    jacp_(1,3) = -c^(-1)*(v-vl);
    jacp_(1,4) = 1/2*c^(-1)*gca*(1+tanh(-(v1-v)*v2^(-1)));
    jacp_(1,5) = c^(-1)*w*gk;
    jacp_(1,6) = c^(-1)*gl;
    jacp_(1,7) = 1/2*c^(-2)*(2*w*(v-vk)*gk+2*gl*(v-vl)-gca*(1+tanh(-(v1-v)*v2^(-1)))*(vca-v)-2*ic);
    jacp_(1,9) = c^(-1);
    jacp_(1,10) = 1/2*c^(-1)*gca*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1);
    jacp_(1,11) = -1/2*c^(-1)*gca*(v1-v)*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2);
    jacp_(2,8) = -1/2*cosh(-1/2*(v3-v)*v4^(-1))*(-1+2*w-tanh(-(v3-v)*v4^(-1)));
    jacp_(2,12) = 1/2*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-1)*cosh(-1/2*(v3-v)*v4^(-1))*phi+1/4*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-1)*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi;
    jacp_(2,13) = -1/2*(v3-v)*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*phi-1/4*sinh(-1/2*(v3-v)*v4^(-1))*(v3-v)*v4^(-2)*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi;

%--------------------------------------------------
%
% The Hessians function.
% This function returns a 3D matrix.
% hess_(n,:,:) is the Hessian of the n-th component of the vector field. That is,
% hess_(n,i,j) is the second partial derivative of the n-th component
% of the vector field, taken with respect to the i-th and j-th variables.
%
function hess_ = MorrisLecar_hess(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    hess_ = zeros(2,2,2);

    hess_(1,1,1) = c^(-1)*(gca*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1)+gca*tanh(-(v1-v)*v2^(-1))*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2));
    hess_(1,1,2) = -c^(-1)*gk;
    hess_(1,2,1) = hess_(1,1,2);

    hess_(2,1,1) = -1/2*sinh(-1/2*(v3-v)*v4^(-1))*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*phi+(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*phi*tanh(-(v3-v)*v4^(-1))-1/8*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi;
    hess_(2,1,2) = -1/2*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-1)*phi;
    hess_(2,2,1) = hess_(2,1,2);

%--------------------------------------------------
%
% The Hessians with respect to the parameters.
% This function returns a 3D matrix.
% hessp_(n,i,j) is the second partial derivative of the n-th component
% of the vector field, taken with respect to the i-th variable and
% the j-th parameter.
%
function hessp_ = MorrisLecar_hessp(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    hessp_ = zeros(2,2,13);

    hessp_(1,1,1) = -1/2*c^(-1)*(1+tanh(-(v1-v)*v2^(-1))+(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1));
    hessp_(1,1,2) = -c^(-1)*w;
    hessp_(1,1,3) = -c^(-1);
    hessp_(1,1,4) = -1/2*c^(-1)*gca*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1);
    hessp_(1,1,7) = 1/2*c^(-2)*(2*gl+2*w*gk+gca*(1+tanh(-(v1-v)*v2^(-1)))+gca*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1));
    hessp_(1,1,10) = -1/2*c^(-1)*(gca*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-1)+2*gca*tanh(-(v1-v)*v2^(-1))*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2));
    hessp_(1,1,11) = 1/2*c^(-1)*(gca*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2)+gca*(v1-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2)+2*gca*(v1-v)*tanh(-(v1-v)*v2^(-1))*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-3));
    hessp_(1,2,2) = -c^(-1)*(v-vk);
    hessp_(1,2,5) = c^(-1)*gk;
    hessp_(1,2,7) = c^(-2)*(v-vk)*gk;

    hessp_(2,1,8) = -1/2*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-1)*cosh(-1/2*(v3-v)*v4^(-1))-1/4*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-1)*(-1+2*w-tanh(-(v3-v)*v4^(-1)));
    hessp_(2,1,12) = 1/2*sinh(-1/2*(v3-v)*v4^(-1))*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*phi-(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*phi*tanh(-(v3-v)*v4^(-1))+1/8*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi;
    hessp_(2,1,13) = -1/2*sinh(-1/2*(v3-v)*v4^(-1))*(v3-v)*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-3)*phi+(v3-v)*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-3)*cosh(-1/2*(v3-v)*v4^(-1))*phi*tanh(-(v3-v)*v4^(-1))-1/8*(v3-v)*v4^(-3)*cosh(-1/2*(v3-v)*v4^(-1))*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi+1/4*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-2)*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi+1/2*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*phi;
    hessp_(2,2,8) = -cosh(-1/2*(v3-v)*v4^(-1));
    hessp_(2,2,12) = 1/2*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-1)*phi;
    hessp_(2,2,13) = -1/2*sinh(-1/2*(v3-v)*v4^(-1))*(v3-v)*v4^(-2)*phi;

%--------------------------------------------------
%
% Third derivatives of the vector field.
% This function returns a 4D matrix.
% der3_(n,i,j,k) is the third partial derivative of the n-th component
% of the vector field, taken with respect to the i-th, j-th and k-th variables.
%
function der3_ = MorrisLecar_der3(t,x_,gca,gk,gl,vca,vk,vl,c,phi,ic,v1,v2,v3,v4)
    v          = x_(1);
    w          = x_(2);

    der3_ = zeros(2,2,2,2);

    der3_(1,1,1,1) = -c^(-1)*(gca*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)^2*v2^(-3)+3*gca*tanh(-(v1-v)*v2^(-1))*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-2)+2*gca*tanh(-(v1-v)*v2^(-1))^2*(vca-v)*(-1+tanh(-(v1-v)*v2^(-1))^2)*v2^(-3));

    der3_(2,1,1,1) = -2*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-3)*cosh(-1/2*(v3-v)*v4^(-1))*phi*tanh(-(v3-v)*v4^(-1))^2-3/8*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-3)*cosh(-1/2*(v3-v)*v4^(-1))*phi-1/16*sinh(-1/2*(v3-v)*v4^(-1))*v4^(-3)*(-1+2*w-tanh(-(v3-v)*v4^(-1)))*phi+3/2*sinh(-1/2*(v3-v)*v4^(-1))*(-1+tanh(-(v3-v)*v4^(-1))^2)*v4^(-3)*phi*tanh(-(v3-v)*v4^(-1))-(-1+tanh(-(v3-v)*v4^(-1))^2)^2*v4^(-3)*cosh(-1/2*(v3-v)*v4^(-1))*phi;
    der3_(2,1,1,2) = -1/4*v4^(-2)*cosh(-1/2*(v3-v)*v4^(-1))*phi;
    der3_(2,1,2,1) = der3_(2,1,1,2);
    der3_(2,2,1,1) = der3_(2,1,1,2);