eigw.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.002408201 @rate -1 vs 0.11626212 @rate -1\n",
      "4.61e-5 @rate -5.707045250487349 vs 5.2525e-5 @rate -11.11208924362413\n",
      "8.0361e-5 @rate 0.8017287662048508 vs 6.8902e-5 @rate 0.39154160336367844\n",
      "8.6112e-5 @rate 0.09971877923150411 vs 8.8803e-5 @rate 0.3660625551874217\n",
      "0.0004025 @rate 2.224702582023649 vs 0.000226741 @rate 1.3523649666280972\n",
      "0.001124501 @rate 1.4822242557410659 vs 0.000605482 @rate 1.4170387843781134\n",
      "0.003957307 @rate 1.8152340469807735 vs 0.002082158 @rate 1.7819235718640394\n",
      "0.015944531 @rate 2.0104707659234657 vs 0.007578031 @rate 1.8637434930507895\n",
      "0.066811255 @rate 2.0670294969626295 vs 0.026756799 @rate 1.8200105854499318\n",
      "0.475548726 @rate 2.8314301068960006 vs 0.12148153 @rate 2.1827595464945193\n",
      "3.712738166 @rate 2.9648184991959483 vs 0.538771185 @rate 2.1489357097050266\n",
      "29.617420932 @rate 2.9958905365506157 vs 2.387525886 @rate 2.147771776618192\n"
     ]
    }
   ],
   "source": [
    "# Performance difference between unoptimized and optimized\n",
    "tfull0, tsym0 = -1, -1\n",
    "for i in 1:12\n",
    "    ncells = 2^i\n",
    "    h = 1/ncells\n",
    "\n",
    "    dA, uA = fill(2/h, ncells+1), fill(-1/h, ncells)\n",
    "    A = SymTridiagonal(dA, uA)\n",
    "\n",
    "    tic(); eig(full(A)); tfull = toq()\n",
    "    tic(); eig(A); tsym = toq()\n",
    "\n",
    "    if tfull0 == -1\n",
    "        ratefull, ratesym = -1, -1\n",
    "    else\n",
    "        ratefull = log(tfull/tfull0)/log(2)\n",
    "        ratesym = log(tsym/tsym0)/log(2)\n",
    "    end\n",
    "    println(\"$(tfull) @rate $(ratefull) vs $(tsym) @rate $(ratesym)\")\n",
    "\n",
    "    tfull0, tsym0 = tfull, tsym\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.4444106721420094 vs 1.5278640450004204, 17.338198023510163 vs 10.472135954999587, 1.6556505853261574\n",
      "2.5796859100729246 vs 2.6350279389150306, 130.1756467559433 vs 55.31385720702505, 2.3534002748846548\n",
      "4.472658485694915 vs 4.49823217766252, 687.7477709395416 vs 246.42260458071084, 2.7909280973218653\n",
      "6.775048844637732 vs 6.7723958993493, 2985.4686949113047 vs 1014.2022651946505, 2.943662026172185\n",
      "8.40452647785413 vs 8.3975538482612, 12199.757130300295 vs 4086.1481746950753, 2.985637477821198\n",
      "9.189856976904816 vs 9.186977044310241, 49063.32000768427 vs 16374.134814440724, 2.9963915995374735\n",
      "9.545182348031528 vs 9.544321857695731, 196519.210198246 vs 65526.13149692229, 2.999096783356978\n",
      "9.711378970747527 vs 9.711147156553004, 786343.1827203796 vs 262134.13067049734, 2.9997741259752746\n",
      "9.7914943560518 vs 9.791434387902354, 3.1456391758502815e6 vs 1.0485661304642688e6, 2.9999435271264208\n",
      "9.830800318691868 vs 9.830785078921423, 1.258282317413284e7 vs 4.194294130412761e6, 2.9999858815086395\n"
     ]
    }
   ],
   "source": [
    "# How  A, M compared with A, lumped M\n",
    "for i in 1:10\n",
    "    ncells = 2^i\n",
    "    h = 1/ncells\n",
    "\n",
    "    dA, uA = fill(2/h, ncells+1), fill(-1/h, ncells)\n",
    "    A = SymTridiagonal(dA, uA)\n",
    "    dM, uM = [1.; fill(4h/6, ncells-1); 1.], fill(h/6, ncells)\n",
    "    M = SymTridiagonal(dM, uM)\n",
    "    eigw, _ = eig(full(A), full(M))\n",
    "\n",
    "    \n",
    "    idMl = 1./[1; fill(4h/6, ncells-1) + 2*fill(h/6, ncells-1); 1]\n",
    "    invM = Diagonal(idMl)\n",
    "    B = invM*A\n",
    "    eigwl, _ = eig(B)\n",
    "\n",
    "    lmin, lmax = minimum(eigw), maximum(eigw)\n",
    "    lminl, lmaxl = minimum(eigwl), maximum(eigwl)\n",
    "\n",
    "    println(\"$(lmin) vs $(lminl), $(lmax) vs $(lmaxl), $(lmax/lmaxl)\")\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.4.6-pre",
   "language": "julia",
   "name": "julia-0.4"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.4.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}