{
"cells": [
 {
		   "cell_type": "markdown",
	   "metadata": {},
	   "source": [
       "# Chapter 8: Bifurcations Revisited"
	   ]
	},
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.1_1: Saddle_Node_Transcritical_Pitchfork_Bifurcations.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"clear;\n",
"clear;\n",
"clc;\n",
"close;\n",
"set(gca(),'auto_clear','off')        //hold on\n",
"\n",
"a=0.4;\n",
"b=0.8;\n",
"for x=0:0.1:3\n",
"    y1=a*x;\n",
"    y2=(x^2)/(b*(1+x^2));\n",
"    plot2d(x,y1,style=-2)\n",
"    plot2d(x,y2,style=-3)\n",
"end\n",
"\n",
"// Classification of fixed points :\n",
"\n",
"A1=[-a 1;0 -b]                //Jacobian at (0,0)\n",
"T=trace(A1)                    //Trace of A\n",
"D=det(A1)                      //Determinant of A\n",
"\n",
"disp('Since, D>0, T<0 , orign is always a fixed point.')\n",
"\n",
"//Now using the arguments given in book and the figure obtained through this example, we conclude :\n",
"\n",
"disp('Middle Fixed Point lies between 0<x*<1, Thus is a Saddle Point.')\n",
"disp('The Thied fixed point is with x*>1, Thus always a stable node.')\n",
"\n",
"xtitle('Nullclines--Showing Intersection of x(dot) and y(dot)','x-Axis ( x )','y-Axis ( y )')\n",
"figure\n",
"a=0.4;\n",
"b=0.8;\n",
"  function xd=linear811(t,x)\n",
"    xd(1)=-a*x(1)+x(2);\n",
"    xd(2)=((x(1)^2)/(1+x(1)^2))-b*x(2);\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[0,0,4,4];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=20;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:20\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear811);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.1_2: Saddle_Node_Transcritical_Pitchfork_Bifurcations.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"clear;\n",
"clear;\n",
"clc;\n",
"close;\n",
"set(gca(),'auto_clear','off')        //hold on\n",
"\n",
"//f = x(dot)\n",
"//g = y(dot)\n",
"//Obtain Jacobian as :\n",
"//J=[df/dx  df/dy; dg/dx  dg/dy]\n",
"//Thus, J=[mew-3*x^2  0; 0  -1]\n",
"mew=2;\n",
"subplot(221)\n",
"for x=-2:0.1:2\n",
"    y1=(mew*x)-x^3;\n",
"    y2=0;\n",
"    plot2d(x,y1,style=-2)\n",
"    plot2d(x,y2,style=-3)\n",
"end\n",
"xtitle('Nullclines--Showing Intersection of x(dot) and y(dot) for mew > 0','x-Axis ( x )','y-Axis ( y )')\n",
"//Stabilities for mew >0\n",
"A1=[mew 0;0 -1]            //Jacobian at (0,0)\n",
"T1=trace(A1)\n",
"D1=det(A1)\n",
"A2=[-2*mew 0;0 -1]        //Jacobian at (+sqrt(mew),0) and (-sqrt(mew),0)\n",
"T2=trace(A2)\n",
"D2=det(A2)\n",
"Det = (T2^2) - 4*D2\n",
"\n",
"disp('Since D1<0, Thus (0,0) is a Saddle Point.')\n",
"disp('Since T2<0, D2>0 and Det>0, thus (+sqrt(mew),0) and (-sqrt(mew),0) are stable points.')\n",
"subplot(224)\n",
"mew=-2;\n",
"for x=-2:0.1:2\n",
"    y1=(mew*x)-x^3;\n",
"    y2=0;\n",
"    plot2d(x,y1,style=-2)\n",
"    plot2d(x,y2,style=-3)\n",
"end\n",
"xtitle('Nullclines--Showing Intersection of x(dot) and y(dot) for mew < 0','x-Axis ( x )','y-Axis ( y )')\n",
"// Classification of fixed points :\n",
"\n",
"A3=[mew 0;0 -1]                //Jacobian at (0,0)\n",
"T3=trace(A3)                    //Trace of A\n",
"D3=det(A3)                      //Determinant of A\n",
"Det3=T3^2 - 4*D3\n",
"disp('Since, D>0, T<0 , and Det3>0 orign is a stable fixed point.')\n",
"\n",
"\n",
"//Note Stabilities can be deduced from the figures also as done in previous chapters.\n",
"figure\n",
"mew=-4;\n",
"  function xd=linear812(t,x)\n",
"    xd(1)=mew*x(1)-(x(1)^3);\n",
"    xd(2)=-x(2);\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[-4,-4,4,14];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=10;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:10\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear812);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.1_3: Saddle_Node_Transcritical_Pitchfork_Bifurcations.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"mew=-2.1;\n",
"  function xd=linear813(t,x)\n",
"    xd(1)=mew*x(1)+x(2)+sin(x(1));\n",
"    xd(2)=x(1)-x(2);\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[-2,-2,2,2];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=20;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:20\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear813);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.2_1: Hopf_Bifurcations.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"mew=-0.2;\n",
"  function xd=linear821(t,x)\n",
"    xd(1)=mew*x(1)-x(2)+x(1)*x(2)^2;\n",
"    xd(2)=x(1)+mew*x(2)+x(2)^3;\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[-0.8,-0.8,0.8,0.8];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=35;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:35\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear821);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.3_1: Oscillating_Chemical_Reactions.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"a=10;\n",
"b=3;\n",
"  function xd=linear831(t,x)\n",
"    xd(1)=a-x(1)-((4*x(1)*x(2))/(1+x(1)^2));\n",
"    xd(2)=(b*x(1))*(1-(x(2)/(1+x(1)^2)));\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[0,3,3,7];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=25;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:25\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear831);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
,
{
		   "cell_type": "markdown",
		   "metadata": {},
		   "source": [
			"## Example 8.3_2: Oscillating_Chemical_Reactions.sce"
		   ]
		  },
  {
"cell_type": "code",
	   "execution_count": null,
	   "metadata": {
	    "collapsed": true
	   },
	   "outputs": [],
"source": [
"a=10;\n",
"b=4;\n",
"  function xd=linear832(t,x)\n",
"    xd(1)=a-x(1)-((4*x(1)*x(2))/(1+x(1)^2));\n",
"    xd(2)=(b*x(1))*(1-(x(2)/(1+x(1)^2)));\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[0,0,4,10];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=10;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:10\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear832);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')\n",
"  figure\n",
"  a=10;\n",
"b=2;\n",
"  function xd=linear8322(t,x)\n",
"    xd(1)=a-x(1)-((4*x(1)*x(2))/(1+x(1)^2));\n",
"    xd(2)=(b*x(1))*(1-(x(2)/(1+x(1)^2)));\n",
"             //x(dot);   x(2) means y.\n",
"             //y(dot);   x(1) means x.; \n",
"  endfunction\n",
"  bound=[0,0,5,8];        //Bounds of x-axis and y-axis as [xmin ymin xmax ymax], change them according to your needs.\n",
"  nrect=10;     //increase it to get more number of curves, i.e. more information will be available.\n",
"  set(gca(),'auto_clear','off')        //hold on\n",
"  x=linspace(bound(1),bound(3),nrect);\n",
"  y=linspace(bound(2),bound(4),nrect);\n",
"  x0=[];\n",
" \n",
" for i=1:10\n",
"     x0=[x(i);y(i)];\n",
"     t0=0;\n",
"     t=0:0.01:3000;\n",
"     xout=ode(x0,t0,t,linear8322);\n",
"     plot2d(xout(1,:),xout(2,:));\n",
" end\n",
"  xtitle('Phase Portrait','x-axis ( x )','y-axis ( y )')"
   ]
   }
],
"metadata": {
		  "kernelspec": {
		   "display_name": "Scilab",
		   "language": "scilab",
		   "name": "scilab"
		  },
		  "language_info": {
		   "file_extension": ".sce",
		   "help_links": [
			{
			 "text": "MetaKernel Magics",
			 "url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md"
			}
		   ],
		   "mimetype": "text/x-octave",
		   "name": "scilab",
		   "version": "0.7.1"
		  }
		 },
		 "nbformat": 4,
		 "nbformat_minor": 0
}