{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Numerical integration\n",
    "\n",
    "`NIntegrate` performs numerical integration. The integrand can be a Symata expression or a compiled function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "application/javascript": [
       "$(\"div.input_prompt\").css({\"color\": \"blue\"})"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "application/javascript": [
       "$(\"div.input_prompt\").css({\"color\": \"blue\"})"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "using Symata; isymata()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the gamma function for `a=3`.  The result and estimated error are returned."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$  \\left[ 2.0000000000000204,4.460650935347041e-9 \\right]  $$"
      ],
      "text/plain": [
       "L\"$$  \\left[ 2.0000000000000204,4.460650935347041e-9 \\right]  $$\""
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "( a = 3.0,\n",
    " NIntegrate(Exp(-t)*t^(a-1), [t,0,Infinity]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Define a function `gamma`, where `gamma(a)` is the gamma function at `a` and `gamma(a,x)` is the lower incomplete gamma function. The default value of the second argument is `Infinity`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "gamma(a_, x_:Infinity) := NIntegrate(Exp(-t)*t^(a-1), [t,0,x])[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Test the function on a `List` of numbers."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$  \\left[ 2.0000000000000204,6.0,23.999999999999996 \\right]  $$"
      ],
      "text/plain": [
       "L\"$$  \\left[ 2.0000000000000204,6.0,23.999999999999996 \\right]  $$\""
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Map(gamma, [3.0, 4.0, 5.0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The integrand in the function `gamma` is a Symata expression that is re-evaluated every time the integration routine asks for a value. Using a compiled function in the integrand is more efficient.\n",
    "\n",
    "Create the compiled function like this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    " gammaint = Compile([t] , Exp(-t)*t^(a-1));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`gammaint` refers to a Julia function of one variable, `x`. The variable `a` is free. To use `gammaint`, we have to set the value of the Julia variable `a` using `SetJ`.\n",
    "\n",
    "The alternative gamma function is"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "gamma1(a1_, x_:Infinity) := (SetJ(a,a1), NIntegrate(gammaint, [0,x]))[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$ 2.0000000000000213 $$"
      ],
      "text/plain": [
       "L\"$$ 2.0000000000000213 $$\""
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gamma1(3.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Timing shows that `gamma1` is about 50 times faster than `gamma`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$  \\left[ 0.010876366,2.0000000000000204 \\right]  $$"
      ],
      "text/plain": [
       "L\"$$  \\left[ 0.010876366,2.0000000000000204 \\right]  $$\""
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Timing(gamma(3.0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$  \\left[ 0.000649037,2.0000000000000213 \\right]  $$"
      ],
      "text/plain": [
       "L\"$$  \\left[ 0.000649037,2.0000000000000213 \\right]  $$\""
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Timing(gamma1(3.0))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### Version and date"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "symata version 0.3.0-dev.7\n",
      "julia version  0.6.0-dev.435\n",
      "python version 2.7.12\n",
      "sympy version  1.0\n"
     ]
    }
   ],
   "source": [
    "VersionInfo()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$ 2016-11-28T22:48:19.257 $$"
      ],
      "text/plain": [
       "L\"$$ 2016-11-28T22:48:19.257 $$\""
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Now()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.6.0-dev",
   "language": "julia",
   "name": "julia-0.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}