{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Leakage characterization using GST\n",
    "This tutorial demonstrates how to perform GST on a \"leaky-qubit\" described by a 3-level (instead of the desired 2-level) system. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pygsti\n",
    "import pygsti.construction as pc\n",
    "import pygsti.modelpacks.legacy.std1Q_XYI as std1Q\n",
    "import numpy as np\n",
    "import scipy.linalg as sla\n",
    "#import pickle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def to_3level_unitary(U_2level):\n",
    "    U_3level = np.zeros((3,3),complex)\n",
    "    U_3level[0:2,0:2] = U_2level\n",
    "    U_3level[2,2] = 1.0\n",
    "    return U_3level\n",
    "\n",
    "def unitary_to_gmgate(U):\n",
    "    return pygsti.tools.change_basis( \n",
    "        pygsti.tools.unitary_to_process_mx(U), 'std','gm')\n",
    "\n",
    "def state_to_gmvec(state):\n",
    "    pygsti.tools.stdmx_to_gmvec\n",
    "\n",
    "Us = pygsti.tools.internalgates.get_standard_gatename_unitaries()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "mdl_2level_ideal = std1Q.target_model()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "rho0 = np.array( [[1,0,0],\n",
    "                  [0,0,0],\n",
    "                  [0,0,0]], complex)\n",
    "E0 = rho0\n",
    "E1 = np.array( [[0,0,0],\n",
    "                [0,1,0],\n",
    "                [0,0,1]], complex)\n",
    "\n",
    "sslbls = pygsti.obj.StateSpaceLabels(['Qubit+Leakage'],[9])\n",
    "mdl_3level_ideal = pygsti.obj.ExplicitOpModel(sslbls, 'gm')\n",
    "mdl_3level_ideal['rho0'] = pygsti.tools.stdmx_to_gmvec(rho0)\n",
    "mdl_3level_ideal['Mdefault'] = pygsti.obj.TPPOVM([('0',pygsti.tools.stdmx_to_gmvec(E0)),\n",
    "                                                    ('1',pygsti.tools.stdmx_to_gmvec(E1))])\n",
    "\n",
    "mdl_3level_ideal['Gi'] = unitary_to_gmgate( to_3level_unitary(Us['Gi']))\n",
    "mdl_3level_ideal['Gx'] = unitary_to_gmgate( to_3level_unitary(Us['Gxpi2']))\n",
    "mdl_3level_ideal['Gy'] = unitary_to_gmgate( to_3level_unitary(Us['Gypi2']))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigmaX = np.array([[0,1],[1,0]],complex)\n",
    "rot = sla.expm(1j * 0.1 * sigmaX)\n",
    "Uleakage = np.identity(3,complex)\n",
    "Uleakage[1:3,1:3] = rot\n",
    "leakageOp = unitary_to_gmgate(Uleakage)\n",
    "#print(Uleakage)\n",
    "\n",
    "#Guess of a model w/just unitary leakage\n",
    "mdl_3level_guess = mdl_3level_ideal.copy()\n",
    "mdl_3level_guess['Gi'] = np.dot(leakageOp, mdl_3level_guess['Gi'])\n",
    "#mdl_3level_guess['Gx'] = np.dot(leakageOp, mdl_3level_guess['Gx'])\n",
    "#mdl_3level_guess['Gy'] = np.dot(leakageOp, mdl_3level_guess['Gy'])\n",
    "\n",
    "#Actual model used for data generation (some depolarization too)\n",
    "mdl_3level_noisy = mdl_3level_ideal.depolarize(op_noise=0.005, spam_noise=0.01)\n",
    "mdl_3level_noisy['Gi'] = np.dot(leakageOp, mdl_3level_noisy['Gi'])\n",
    "#mdl_3level_noisy['Gx'] = np.dot(leakageOp, mdl_3level_noisy['Gx'])\n",
    "#mdl_3level_noisy['Gy'] = np.dot(leakageOp, mdl_3level_noisy['Gy'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "#print(mdl_3level_guess)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get sequences using expected model\n",
    "generate_fiducials = False\n",
    "\n",
    "if generate_fiducials:\n",
    "    prepfids, measfids = pygsti.algorithms.generate_fiducials(\n",
    "        mdl_3level_guess, omitIdentity=False, maxFidLength=4, verbosity=4)\n",
    "    pygsti.io.write_circuit_list(\"example_files/leakage_prepfids.txt\", prepfids)\n",
    "    pygsti.io.write_circuit_list(\"example_files/leakage_measfids.txt\", measfids)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "prepfids = pygsti.io.load_circuit_list(\"example_files/leakage_prepfids.txt\")\n",
    "measfids = pygsti.io.load_circuit_list(\"example_files/leakage_measfids.txt\")\n",
    "germs = std1Q.germs\n",
    "maxLengths = [1,]\n",
    "expList = pc.make_lsgst_experiment_list(mdl_3level_noisy, prepfids, measfids, germs, maxLengths)\n",
    "ds = pc.generate_fake_data(mdl_3level_noisy, expList, 1000, 'binomial', seed=1234)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--- Circuit Creation ---\n",
      "   280 sequences created\n",
      "-- Std Practice:  Iter 1 of 1  (CPTP) --: \n",
      "  --- Iterative MLGST: Iter 1 of 1  280 operation sequences ---: \n",
      "    --- Minimum Chi^2 GST ---\n",
      "    Sum of Chi^2 = 801.665 (280 data params - 31 model params = expected mean of 249; p-value = 0)\n",
      "    Completed in 8.6s\n",
      "    2*Delta(log(L)) = 809.696\n",
      "    Iteration 1 took 8.7s\n",
      "    \n",
      "    Switching to ML objective (last iteration)\n",
      "    --- MLGST ---\n",
      "      Maximum log(L) = 398.424 below upper bound of -421916\n",
      "        2*Delta(log(L)) = 796.848 (280 data params - 31 model params = expected mean of 249; p-value = 0)\n",
      "      Completed in 1.3s\n",
      "    2*Delta(log(L)) = 796.848\n",
      "    Final MLGST took 1.3s\n",
      "    \n",
      "  Iterative MLGST Total Time: 10.0s\n",
      "--- Circuit Creation ---\n",
      "   280 sequences created\n",
      "-- Std Practice:  Iter 1 of 2  (CPTP) --: \n",
      "  --- Iterative MLGST: Iter 1 of 1  280 operation sequences ---: \n",
      "    --- Minimum Chi^2 GST ---\n",
      "      bulk_evaltree: created initial tree (280 strs) in 0s\n",
      "      bulk_evaltree: split tree (1 subtrees) in 0s\n",
      "      Created evaluation tree with 1 subtrees.  Will divide 1 procs into 1 (subtree-processing)\n",
      "       groups of ~1 procs each, to distribute over 360 params (taken as 1 param groups of ~360 params).\n",
      "      --- Outer Iter 0: norm_f = 1.53652e+06, mu=1, |x|=6.32456e-08, |J|=1918.56\n",
      "      --- Outer Iter 1: norm_f = 5099.26, mu=285.685, |x|=0.119235, |J|=2823.13\n",
      "      --- Outer Iter 2: norm_f = 2609.34, mu=2283.99, |x|=0.213143, |J|=2042.32\n",
      "      --- Outer Iter 3: norm_f = 552.948, mu=1173.1, |x|=0.201959, |J|=2169.09\n",
      "      --- Outer Iter 4: norm_f = 304.854, mu=774.594, |x|=0.212428, |J|=2238.79\n",
      "      --- Outer Iter 5: norm_f = 280.095, mu=1549.12, |x|=0.216086, |J|=2252.51\n",
      "      --- Outer Iter 6: norm_f = 265.947, mu=1548.86, |x|=0.223266, |J|=2286.38\n",
      "      --- Outer Iter 7: norm_f = 263.391, mu=2015.59, |x|=0.230447, |J|=2305.34\n",
      "      --- Outer Iter 8: norm_f = 257.986, mu=2037.5, |x|=0.235554, |J|=2314.3\n",
      "      --- Outer Iter 9: norm_f = 254.684, mu=2065.78, |x|=0.240654, |J|=2322.02\n",
      "      --- Outer Iter 10: norm_f = 253.632, mu=2658.14, |x|=0.246082, |J|=2322.53\n",
      "      --- Outer Iter 11: norm_f = 251.262, mu=2659.7, |x|=0.24992, |J|=2326.78\n",
      "      --- Outer Iter 12: norm_f = 249.992, mu=2659.7, |x|=0.253678, |J|=2328.16\n",
      "      Least squares message = Both actual and predicted relative reductions in the sum of squares are at most 0.01\n",
      "    Sum of Chi^2 = 249.992 (280 data params - 161 model params = expected mean of 119; p-value = 2.49065e-11)\n",
      "    Completed in 18.6s\n",
      "    2*Delta(log(L)) = 250.944\n",
      "    Iteration 1 took 18.7s\n",
      "    \n",
      "    Switching to ML objective (last iteration)\n",
      "    --- MLGST ---\n",
      "      --- Outer Iter 0: norm_f = 125.472, mu=1, |x|=0.253678, |J|=1648.25\n",
      "      Least squares message = Both actual and predicted relative reductions in the sum of squares are at most 0.01\n",
      "      Maximum log(L) = 125.472 below upper bound of -421916\n",
      "        2*Delta(log(L)) = 250.944 (280 data params - 161 model params = expected mean of 119; p-value = 1.92685e-11)\n",
      "      Completed in 1.0s\n",
      "    2*Delta(log(L)) = 250.944\n",
      "    Final MLGST took 1.0s\n",
      "    \n",
      "  Iterative MLGST Total Time: 19.6s\n",
      "      -- Performing 'stdgaugeopt' gauge optimization on CPTP estimate --\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n",
      "WARNING: MLGST failed to improve logl: retaining chi2-objective estimate\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-- Std Practice:  Iter 2 of 2  (True) --: \n",
      "      -- Performing 'stdgaugeopt' gauge optimization on True estimate --\n"
     ]
    }
   ],
   "source": [
    "results_2level = pygsti.do_stdpractice_gst(ds, mdl_2level_ideal, prepfids, measfids,\n",
    "                                           germs, maxLengths, modes=\"CPTP\", verbosity=3)\n",
    "results_3level = pygsti.do_stdpractice_gst(ds, mdl_3level_ideal, prepfids, measfids,\n",
    "                                           germs, maxLengths, modes=\"CPTP,True\",\n",
    "                                           modelsToTest={'True': mdl_3level_noisy}, \n",
    "                                           verbosity=4, advancedOptions={'all': {'tolerance': 1e-2}})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Running idle tomography\n",
      "Computing switchable properties\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1236: UserWarning:\n",
      "\n",
      "Idle tomography failed:\n",
      "operands could not be broadcast together with shapes (9,) (16,) \n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Found standard clifford compilation from std1Q_XYI\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1303: UserWarning:\n",
      "\n",
      "PDF output is not available for reports with multiple datasets\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n"
     ]
    }
   ],
   "source": [
    "pygsti.report.construct_standard_report(\n",
    "    {'two-level': results_2level, 'three-level': results_3level},\n",
    "    \"Leakage Example Report\"\n",
    ").write_html('example_files/leakage_report')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "#try a different basis:\n",
    "gm_basis = pygsti.obj.Basis.cast('gm',9)\n",
    "   \n",
    "leakage_basis_mxs = [ np.sqrt(2)/3*(np.sqrt(3)*gm_basis[0] + 0.5*np.sqrt(6)*gm_basis[8]),\n",
    "                      gm_basis[1], gm_basis[4], gm_basis[7],\n",
    "                     gm_basis[2], gm_basis[3], gm_basis[5], gm_basis[6],\n",
    "                     1/3*(np.sqrt(3)*gm_basis[0] - np.sqrt(6)*gm_basis[8]) ]\n",
    "#for mx in leakage_basis_mxs:\n",
    "#    pygsti.tools.print_mx(mx)\n",
    "\n",
    "check = np.zeros( (9,9), complex)\n",
    "for i,m1 in enumerate(leakage_basis_mxs):\n",
    "    for j,m2 in enumerate(leakage_basis_mxs):\n",
    "        check[i,j] = np.trace(np.dot(m1,m2))\n",
    "assert(np.allclose(check, np.identity(9,complex)))\n",
    "\n",
    "leakage_basis = pygsti.obj.ExplicitBasis(leakage_basis_mxs, name=\"LeakageBasis\",  \n",
    "                                        longname=\"2+1 level leakage basis\", real=True,\n",
    "                                        labels=['I','X','Y','Z','LX0','LX1','LY0','LY1','L'])\n",
    "\n",
    "def changebasis_3level_model(mdl):\n",
    "    new_mdl = mdl.copy()\n",
    "    new_mdl.preps['rho0'] = pygsti.obj.FullSPAMVec(\n",
    "        pygsti.tools.change_basis(mdl.preps['rho0'].todense(), gm_basis, leakage_basis))\n",
    "    new_mdl.povms['Mdefault'] = pygsti.obj.UnconstrainedPOVM(\n",
    "        [('0', pygsti.tools.change_basis(mdl.povms['Mdefault']['0'].todense(), gm_basis, leakage_basis)),\n",
    "         ('1', pygsti.tools.change_basis(mdl.povms['Mdefault']['1'].todense(), gm_basis, leakage_basis))])\n",
    "    \n",
    "    for lbl,op in mdl.operations.items():\n",
    "        new_mdl.operations[lbl] = pygsti.obj.FullDenseOp(\n",
    "            pygsti.tools.change_basis(op.todense(), gm_basis, leakage_basis))\n",
    "    new_mdl.basis = leakage_basis\n",
    "    return new_mdl\n",
    "\n",
    "def changebasis_3level_results(results):\n",
    "    new_results = results.copy()\n",
    "    for estlbl,est in results.estimates.items():\n",
    "        for mlbl,mdl in est.models.items():\n",
    "            if isinstance(mdl,(list,tuple)): #assume a list/tuple of models\n",
    "                new_results.estimates[estlbl].models[mlbl] = \\\n",
    "                    [ changebasis_3level_model(m) for m in mdl ]\n",
    "            else:\n",
    "                new_results.estimates[estlbl].models[mlbl] = changebasis_3level_model(mdl)\n",
    "    return new_results\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "results_3level_leakage_basis = changebasis_3level_results( results_3level )    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Running idle tomography\n",
      "Computing switchable properties\n",
      "Found standard clifford compilation from std1Q_XYI\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1236: UserWarning:\n",
      "\n",
      "Idle tomography failed:\n",
      "operands could not be broadcast together with shapes (9,) (16,) \n",
      "\n",
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1303: UserWarning:\n",
      "\n",
      "PDF output is not available for reports with multiple datasets\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n"
     ]
    }
   ],
   "source": [
    "pygsti.report.construct_standard_report(\n",
    "    {'two-level': results_2level, 'three-level': results_3level_leakage_basis},\n",
    "    \"Leakage Example Report\"\n",
    ").write_html('example_files/leakage_report')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Open the report [here](example_files/leakage_report/main.html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# use \"kite\" density-matrix structure\n",
    "def to_2plus1_superop(superop_2level):\n",
    "    ret = np.zeros((5,5),'d')\n",
    "    ret[0:4,0:4] = superop_2level\n",
    "    ret[4,4] = 1.0 #leave leakage population where it is\n",
    "    return ret\n",
    "\n",
    "#Tack on a single extra \"0\" for the 5-th dimension corresponding\n",
    "# to the classical leakage level population.\n",
    "rho0 = np.concatenate( (mdl_2level_ideal.preps['rho0'],[[0]]), axis=0)\n",
    "E0 = np.concatenate( (mdl_2level_ideal.povms['Mdefault']['0'],[[0]]), axis=0)\n",
    "E1 = np.concatenate( (mdl_2level_ideal.povms['Mdefault']['1'],[[0]]), axis=0)\n",
    "\n",
    "\n",
    "sslbls = pygsti.obj.StateSpaceLabels([('Qubit',),('Leakage',)],[(4,),(1,)])\n",
    "mdl_2plus1_ideal = pygsti.obj.ExplicitOpModel(sslbls, 'gm')\n",
    "mdl_2plus1_ideal['rho0'] = rho0\n",
    "mdl_2plus1_ideal['Mdefault'] = pygsti.obj.UnconstrainedPOVM([('0',E0),('1',E1)])\n",
    "\n",
    "mdl_2plus1_ideal['Gi'] = to_2plus1_superop(mdl_2level_ideal['Gi'])\n",
    "mdl_2plus1_ideal['Gx'] = to_2plus1_superop(mdl_2level_ideal['Gi'])\n",
    "mdl_2plus1_ideal['Gy'] = to_2plus1_superop(mdl_2level_ideal['Gi'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--- Circuit Creation ---\n",
      "   280 sequences created\n",
      "--- Iterative MLGST: Iter 1 of 1  280 operation sequences ---: \n",
      "  --- Minimum Chi^2 GST ---\n",
      "    bulk_evaltree: created initial tree (280 strs) in 0s\n",
      "    bulk_evaltree: split tree (1 subtrees) in 0s\n",
      "    Created evaluation tree with 1 subtrees.  Will divide 1 procs into 1 (subtree-processing)\n",
      "     groups of ~1 procs each, to distribute over 90 params (taken as 1 param groups of ~90 params).\n",
      "    --- Outer Iter 0: norm_f = 1.10463e+09, mu=1, |x|=4.24264, |J|=246913\n",
      "    --- Outer Iter 1: norm_f = 106758, mu=2.30523e+06, |x|=4.23094, |J|=4045.54\n",
      "    --- Outer Iter 2: norm_f = 76699.5, mu=768410, |x|=4.21158, |J|=3590.63\n",
      "    --- Outer Iter 3: norm_f = 69764.7, mu=256137, |x|=4.2006, |J|=3504.09\n",
      "    --- Outer Iter 4: norm_f = 68704, mu=85378.9, |x|=4.19376, |J|=3476.38\n",
      "    --- Outer Iter 5: norm_f = 67412.7, mu=28459.6, |x|=4.17968, |J|=3422.34\n",
      "    --- Outer Iter 6: norm_f = 66403.8, mu=11408.7, |x|=4.166, |J|=3376.7\n",
      "    --- Outer Iter 7: norm_f = 66147.1, mu=6650.94, |x|=4.16608, |J|=3369.15\n",
      "    --- Outer Iter 8: norm_f = 66093.8, mu=2216.98, |x|=4.1746, |J|=3379.62\n",
      "    --- Outer Iter 9: norm_f = 66019, mu=738.994, |x|=4.19649, |J|=3445.39\n",
      "    --- Outer Iter 10: norm_f = 65902.5, mu=246.331, |x|=4.23102, |J|=3579.46\n",
      "    --- Outer Iter 11: norm_f = 65735.9, mu=82.1104, |x|=4.27791, |J|=3821.13\n",
      "    --- Outer Iter 12: norm_f = 64623.4, mu=27.3701, |x|=4.33488, |J|=4241.65\n",
      "    --- Outer Iter 13: norm_f = 64142, mu=2009.45, |x|=4.39261, |J|=5860.49\n",
      "    --- Outer Iter 14: norm_f = 61597.7, mu=1620.53, |x|=4.358, |J|=4917.3\n",
      "    --- Outer Iter 15: norm_f = 60952.1, mu=1631.1, |x|=4.41398, |J|=5481.58\n",
      "    --- Outer Iter 16: norm_f = 59217.3, mu=34796.7, |x|=4.41967, |J|=5771.16\n",
      "    --- Outer Iter 17: norm_f = 58825.1, mu=11598.9, |x|=4.43114, |J|=6106.27\n",
      "    --- Outer Iter 18: norm_f = 58221.7, mu=3866.3, |x|=4.45154, |J|=6390.61\n",
      "    --- Outer Iter 19: norm_f = 57843.9, mu=3807.24, |x|=4.47813, |J|=55133.6\n",
      "    --- Outer Iter 20: norm_f = 57288.3, mu=2123.74, |x|=4.48377, |J|=32733.1\n",
      "    --- Outer Iter 21: norm_f = 56947.3, mu=707.912, |x|=4.467, |J|=13187.7\n",
      "    --- Outer Iter 22: norm_f = 56610.9, mu=235.971, |x|=4.42759, |J|=10147.5\n",
      "    --- Outer Iter 23: norm_f = 56271.9, mu=78.6569, |x|=4.37557, |J|=8993.2\n",
      "    --- Outer Iter 24: norm_f = 56029.3, mu=26.219, |x|=4.33367, |J|=12273.5\n",
      "    --- Outer Iter 25: norm_f = 55968.6, mu=8949.41, |x|=4.33517, |J|=10190.1\n",
      "    --- Outer Iter 26: norm_f = 55717.8, mu=5966.27, |x|=4.33879, |J|=16582.5\n",
      "    --- Outer Iter 27: norm_f = 55091.9, mu=15910.1, |x|=4.34563, |J|=14702.6\n",
      "    --- Outer Iter 28: norm_f = 53733.6, mu=5303.35, |x|=4.37322, |J|=6980.79\n",
      "    --- Outer Iter 29: norm_f = 53575.5, mu=1.81021e+06, |x|=4.37335, |J|=12537.6\n",
      "    --- Outer Iter 30: norm_f = 53541.4, mu=1.0088e+06, |x|=4.3736, |J|=11288.9\n",
      "    --- Outer Iter 31: norm_f = 53517.7, mu=336267, |x|=4.37396, |J|=11721.2\n",
      "    --- Outer Iter 32: norm_f = 53456.6, mu=112089, |x|=4.37489, |J|=11347.8\n",
      "    --- Outer Iter 33: norm_f = 53308.1, mu=37363, |x|=4.37748, |J|=10869\n",
      "    --- Outer Iter 34: norm_f = 52988.4, mu=12454.3, |x|=4.38471, |J|=9499.89\n",
      "    --- Outer Iter 35: norm_f = 52265.9, mu=4151.44, |x|=4.40602, |J|=7563.08\n",
      "    --- Outer Iter 36: norm_f = 50059.7, mu=1383.81, |x|=4.47507, |J|=5821.82\n",
      "    --- Outer Iter 37: norm_f = 49913.7, mu=1.51149e+07, |x|=4.47508, |J|=6066.77\n",
      "    --- Outer Iter 38: norm_f = 49681.8, mu=5.03831e+06, |x|=4.47516, |J|=9698.72\n",
      "    --- Outer Iter 39: norm_f = 49544.5, mu=1.67944e+06, |x|=4.47557, |J|=13114.7\n",
      "    --- Outer Iter 40: norm_f = 49363.6, mu=559812, |x|=4.47658, |J|=10953.2\n",
      "    --- Outer Iter 41: norm_f = 49060.3, mu=186604, |x|=4.47909, |J|=11919\n",
      "    --- Outer Iter 42: norm_f = 48497.6, mu=62201.4, |x|=4.4844, |J|=10073.4\n",
      "    --- Outer Iter 43: norm_f = 47036.6, mu=20733.8, |x|=4.4955, |J|=8287.01\n",
      "    --- Outer Iter 44: norm_f = 43664.5, mu=6911.26, |x|=4.51367, |J|=18627.6\n",
      "    --- Outer Iter 45: norm_f = 42439.7, mu=18430, |x|=4.51077, |J|=18349.8\n",
      "    --- Outer Iter 46: norm_f = 40690.4, mu=6143.34, |x|=4.50065, |J|=14519.8\n",
      "    --- Outer Iter 47: norm_f = 39258.6, mu=4095.56, |x|=4.49322, |J|=15991.5\n",
      "    --- Outer Iter 48: norm_f = 39137.7, mu=4.08175e+06, |x|=4.49239, |J|=27650.7\n",
      "    --- Outer Iter 49: norm_f = 39077.6, mu=3.20948e+06, |x|=4.49197, |J|=15364.5\n",
      "    --- Outer Iter 50: norm_f = 39060.6, mu=1.06983e+06, |x|=4.49167, |J|=14759.3\n",
      "    --- Outer Iter 51: norm_f = 39026.2, mu=356608, |x|=4.49133, |J|=14753.7\n",
      "    --- Outer Iter 52: norm_f = 38948.5, mu=118869, |x|=4.49111, |J|=14667.9\n",
      "    --- Outer Iter 53: norm_f = 38746.7, mu=39623.2, |x|=4.49069, |J|=14725.2\n",
      "    --- Outer Iter 54: norm_f = 38191.5, mu=13207.7, |x|=4.48815, |J|=14805.5\n",
      "    --- Outer Iter 55: norm_f = 37334.8, mu=8805.15, |x|=4.48211, |J|=15766.5\n",
      "    --- Outer Iter 56: norm_f = 36044, mu=5870.1, |x|=4.46811, |J|=16284.8\n",
      "    --- Outer Iter 57: norm_f = 33558, mu=1956.7, |x|=4.40536, |J|=9401.41\n",
      "    --- Outer Iter 58: norm_f = 33362.2, mu=2.13724e+07, |x|=4.40534, |J|=10788.6\n",
      "    --- Outer Iter 59: norm_f = 33226.9, mu=7.12413e+06, |x|=4.40534, |J|=14217.4\n",
      "    --- Outer Iter 60: norm_f = 33107.4, mu=2.37471e+06, |x|=4.40557, |J|=14789.5\n",
      "    --- Outer Iter 61: norm_f = 32975.4, mu=791569, |x|=4.40624, |J|=15911.7\n",
      "    --- Outer Iter 62: norm_f = 32872.5, mu=263856, |x|=4.40635, |J|=16506.5\n",
      "    --- Outer Iter 63: norm_f = 32731, mu=87952.2, |x|=4.40395, |J|=16222\n",
      "    --- Outer Iter 64: norm_f = 32395.3, mu=29317.4, |x|=4.39537, |J|=15947.3\n",
      "    --- Outer Iter 65: norm_f = 31534.2, mu=9772.46, |x|=4.36765, |J|=14973.7\n",
      "    --- Outer Iter 66: norm_f = 29270.6, mu=3257.49, |x|=4.28625, |J|=26456.8\n",
      "    --- Outer Iter 67: norm_f = 29216.4, mu=4.2735e+06, |x|=4.28698, |J|=35768.3\n",
      "    --- Outer Iter 68: norm_f = 29010.6, mu=3.73779e+06, |x|=4.28722, |J|=15713.7\n",
      "    --- Outer Iter 69: norm_f = 28977.7, mu=2.95456e+06, |x|=4.28728, |J|=13925.3\n",
      "    --- Outer Iter 70: norm_f = 28955.7, mu=984852, |x|=4.28731, |J|=14223.8\n",
      "    --- Outer Iter 71: norm_f = 28899.2, mu=328284, |x|=4.28696, |J|=14106.5\n",
      "    --- Outer Iter 72: norm_f = 28753.2, mu=109428, |x|=4.28512, |J|=14081.2\n",
      "    --- Outer Iter 73: norm_f = 28356.5, mu=36476, |x|=4.27824, |J|=13764.1\n",
      "    --- Outer Iter 74: norm_f = 27225.4, mu=12158.7, |x|=4.25627, |J|=13902.6\n",
      "    --- Outer Iter 75: norm_f = 25349.9, mu=8105.78, |x|=4.22613, |J|=15683\n",
      "    --- Outer Iter 76: norm_f = 23846.6, mu=16203.5, |x|=4.19325, |J|=47491.6\n",
      "    --- Outer Iter 77: norm_f = 21781.4, mu=16393.1, |x|=4.18138, |J|=25426.8\n",
      "    --- Outer Iter 78: norm_f = 14887.8, mu=10881.9, |x|=4.18755, |J|=32961.5\n",
      "    --- Outer Iter 79: norm_f = 11140.5, mu=29018.5, |x|=4.19924, |J|=13944.4\n",
      "    --- Outer Iter 80: norm_f = 8733.16, mu=15541.7, |x|=4.20191, |J|=14970.4\n",
      "    --- Outer Iter 81: norm_f = 5128.43, mu=5996.24, |x|=4.21039, |J|=16179.1\n",
      "    --- Outer Iter 82: norm_f = 2758.87, mu=4655.2, |x|=4.21444, |J|=96415.2\n",
      "    --- Outer Iter 83: norm_f = 2743.64, mu=17885.1, |x|=4.21137, |J|=15328.5\n",
      "    --- Outer Iter 84: norm_f = 1352.67, mu=5961.7, |x|=4.20794, |J|=11149.5\n",
      "    --- Outer Iter 85: norm_f = 977.541, mu=1987.23, |x|=4.20464, |J|=12490.5\n",
      "    --- Outer Iter 86: norm_f = 846.579, mu=5299.29, |x|=4.20339, |J|=11201.8\n",
      "    --- Outer Iter 87: norm_f = 733.459, mu=1766.43, |x|=4.20559, |J|=11365.5\n",
      "    --- Outer Iter 88: norm_f = 656.879, mu=588.81, |x|=4.21246, |J|=11168.4\n",
      "    --- Outer Iter 89: norm_f = 645.142, mu=196.27, |x|=4.21567, |J|=11058.9\n",
      "    --- Outer Iter 90: norm_f = 643.842, mu=65.4233, |x|=4.21665, |J|=11053.9\n",
      "    --- Outer Iter 91: norm_f = 642.964, mu=21.8078, |x|=4.22061, |J|=11038.7\n",
      "    --- Outer Iter 92: norm_f = 641.277, mu=7.26925, |x|=4.23584, |J|=11018.3\n",
      "    --- Outer Iter 93: norm_f = 639.658, mu=7.23618, |x|=4.27663, |J|=10871.9\n",
      "    --- Outer Iter 94: norm_f = 637.107, mu=7.28484, |x|=4.31317, |J|=10649.4\n",
      "    --- Outer Iter 95: norm_f = 630.214, mu=11.6824, |x|=4.33987, |J|=10484.5\n",
      "    --- Outer Iter 96: norm_f = 626.735, mu=9.51284, |x|=4.34718, |J|=10472.6\n",
      "    --- Outer Iter 97: norm_f = 625.322, mu=8.10317, |x|=4.34505, |J|=10398.2\n",
      "    --- Outer Iter 98: norm_f = 624.68, mu=7.5562, |x|=4.33819, |J|=10333.6\n",
      "    --- Outer Iter 99: norm_f = 624.289, mu=7.11888, |x|=4.33065, |J|=10283.2\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    --- Outer Iter 100: norm_f = 624.028, mu=6.93394, |x|=4.32467, |J|=10246.6\n",
      "    --- Outer Iter 101: norm_f = 623.818, mu=6.8888, |x|=4.32114, |J|=10218.6\n",
      "    --- Outer Iter 102: norm_f = 623.629, mu=6.88769, |x|=4.32045, |J|=10195.9\n",
      "    --- Outer Iter 103: norm_f = 623.459, mu=6.9114, |x|=4.32287, |J|=10173.4\n",
      "    --- Outer Iter 104: norm_f = 623.308, mu=7.31403, |x|=4.32875, |J|=10148.8\n",
      "    --- Outer Iter 105: norm_f = 623.125, mu=7.90122, |x|=4.3379, |J|=10119.4\n",
      "    --- Outer Iter 106: norm_f = 622.938, mu=9.174, |x|=4.35033, |J|=10086.6\n",
      "    --- Outer Iter 107: norm_f = 622.607, mu=9.64888, |x|=4.36504, |J|=10058.8\n",
      "    --- Outer Iter 108: norm_f = 622.501, mu=15.1712, |x|=4.3845, |J|=10029.4\n",
      "    --- Outer Iter 109: norm_f = 621.23, mu=13.3092, |x|=4.39997, |J|=10069.7\n",
      "    --- Outer Iter 110: norm_f = 620.384, mu=9.13225, |x|=4.4113, |J|=10119.6\n",
      "    --- Outer Iter 111: norm_f = 620.268, mu=23.4506, |x|=4.43217, |J|=10096.4\n",
      "    --- Outer Iter 112: norm_f = 619.264, mu=14.5716, |x|=4.45056, |J|=10186.6\n",
      "    --- Outer Iter 113: norm_f = 618.615, mu=11.9346, |x|=4.46778, |J|=10253.6\n",
      "    --- Outer Iter 114: norm_f = 618.129, mu=20.5586, |x|=4.49041, |J|=10305.2\n",
      "    --- Outer Iter 115: norm_f = 617.563, mu=17.4314, |x|=4.51711, |J|=10394.6\n",
      "    --- Outer Iter 116: norm_f = 617.054, mu=16.2147, |x|=4.54767, |J|=10506.3\n",
      "    --- Outer Iter 117: norm_f = 616.572, mu=14.439, |x|=4.57815, |J|=10643.3\n",
      "    --- Outer Iter 118: norm_f = 616.234, mu=13.7529, |x|=4.60851, |J|=10773.3\n",
      "    --- Outer Iter 119: norm_f = 615.955, mu=12.282, |x|=4.6361, |J|=10898.2\n",
      "    --- Outer Iter 120: norm_f = 615.776, mu=11.5989, |x|=4.6626, |J|=10997.5\n",
      "    --- Outer Iter 121: norm_f = 615.602, mu=9.78855, |x|=4.6876, |J|=11077.5\n",
      "    --- Outer Iter 122: norm_f = 615.382, mu=7.97086, |x|=4.71647, |J|=11114.4\n",
      "    --- Outer Iter 123: norm_f = 615.083, mu=8.21206, |x|=4.76099, |J|=11091.5\n",
      "    --- Outer Iter 124: norm_f = 610.816, mu=34.0813, |x|=4.77627, |J|=11538.9\n",
      "    --- Outer Iter 125: norm_f = 608.377, mu=24.7929, |x|=4.78397, |J|=11214\n",
      "    --- Outer Iter 126: norm_f = 607.645, mu=26.8509, |x|=4.7967, |J|=11426.2\n",
      "    --- Outer Iter 127: norm_f = 605.254, mu=8.95029, |x|=4.80474, |J|=11630.8\n",
      "    --- Outer Iter 128: norm_f = 605.037, mu=8.93571, |x|=4.83075, |J|=11812.2\n",
      "    --- Outer Iter 129: norm_f = 604.797, mu=2.97857, |x|=4.84959, |J|=11861.7\n",
      "    --- Outer Iter 130: norm_f = 604.766, mu=3.54779, |x|=4.87109, |J|=11953.4\n",
      "    --- Outer Iter 131: norm_f = 604.751, mu=3.54567, |x|=4.89697, |J|=12041.5\n",
      "    --- Outer Iter 132: norm_f = 604.725, mu=2.42057, |x|=4.9179, |J|=12112.3\n",
      "    --- Outer Iter 133: norm_f = 604.715, mu=2.43974, |x|=4.94649, |J|=12184\n",
      "    --- Outer Iter 134: norm_f = 604.686, mu=2.07631, |x|=4.97295, |J|=12248.9\n",
      "    --- Outer Iter 135: norm_f = 604.672, mu=2.08699, |x|=5.00495, |J|=12320.3\n",
      "    --- Outer Iter 136: norm_f = 604.643, mu=2.0358, |x|=5.03653, |J|=12394.1\n",
      "    --- Outer Iter 137: norm_f = 604.622, mu=2.03575, |x|=5.07018, |J|=12472.4\n",
      "    --- Outer Iter 138: norm_f = 604.597, mu=2.03155, |x|=5.10442, |J|=12556.7\n",
      "    --- Outer Iter 139: norm_f = 604.573, mu=2.03116, |x|=5.13967, |J|=12646.7\n",
      "    --- Outer Iter 140: norm_f = 604.549, mu=2.03069, |x|=5.17567, |J|=12743.2\n",
      "    --- Outer Iter 141: norm_f = 604.525, mu=2.03052, |x|=5.21242, |J|=12846.3\n",
      "    --- Outer Iter 142: norm_f = 604.501, mu=2.03038, |x|=5.24981, |J|=12956.3\n",
      "    --- Outer Iter 143: norm_f = 604.477, mu=2.03025, |x|=5.28778, |J|=13073.3\n",
      "    --- Outer Iter 144: norm_f = 604.453, mu=2.03009, |x|=5.32625, |J|=13197.3\n",
      "    --- Outer Iter 145: norm_f = 604.429, mu=2.02983, |x|=5.36517, |J|=13328.3\n",
      "    --- Outer Iter 146: norm_f = 604.406, mu=2.02939, |x|=5.40444, |J|=13466.4\n",
      "    --- Outer Iter 147: norm_f = 604.384, mu=2.02862, |x|=5.44402, |J|=13611.5\n",
      "    --- Outer Iter 148: norm_f = 604.361, mu=2.02732, |x|=5.48385, |J|=13763.5\n",
      "    --- Outer Iter 149: norm_f = 604.339, mu=2.02517, |x|=5.52386, |J|=13922.3\n",
      "    --- Outer Iter 150: norm_f = 604.318, mu=2.02183, |x|=5.56403, |J|=14087.8\n",
      "    --- Outer Iter 151: norm_f = 604.297, mu=2.01688, |x|=5.60433, |J|=14260\n",
      "    --- Outer Iter 152: norm_f = 604.277, mu=2.00991, |x|=5.64472, |J|=14438.7\n",
      "    --- Outer Iter 153: norm_f = 604.257, mu=2.00055, |x|=5.6852, |J|=14624\n",
      "    --- Outer Iter 154: norm_f = 604.238, mu=1.98851, |x|=5.72578, |J|=14815.8\n",
      "    --- Outer Iter 155: norm_f = 604.219, mu=1.97361, |x|=5.76647, |J|=15014.3\n",
      "    --- Outer Iter 156: norm_f = 604.201, mu=1.95581, |x|=5.80728, |J|=15219.6\n",
      "    --- Outer Iter 157: norm_f = 604.184, mu=1.93518, |x|=5.84824, |J|=15431.8\n",
      "    --- Outer Iter 158: norm_f = 604.167, mu=1.91188, |x|=5.88938, |J|=15651\n",
      "    --- Outer Iter 159: norm_f = 604.15, mu=1.88615, |x|=5.93071, |J|=15877.5\n",
      "    --- Outer Iter 160: norm_f = 604.134, mu=1.85827, |x|=5.97227, |J|=16111.3\n",
      "    --- Outer Iter 161: norm_f = 604.118, mu=1.82855, |x|=6.01407, |J|=16352.8\n",
      "    --- Outer Iter 162: norm_f = 604.103, mu=1.79727, |x|=6.05613, |J|=16602\n",
      "    --- Outer Iter 163: norm_f = 604.088, mu=1.76473, |x|=6.09848, |J|=16859.2\n",
      "    --- Outer Iter 164: norm_f = 604.073, mu=1.73119, |x|=6.14113, |J|=17124.5\n",
      "    --- Outer Iter 165: norm_f = 604.058, mu=1.69687, |x|=6.18409, |J|=17397.9\n",
      "    --- Outer Iter 166: norm_f = 604.044, mu=1.662, |x|=6.22736, |J|=17679.8\n",
      "    --- Outer Iter 167: norm_f = 604.03, mu=1.62674, |x|=6.27096, |J|=17970.1\n",
      "    --- Outer Iter 168: norm_f = 604.017, mu=1.59127, |x|=6.31489, |J|=18269\n",
      "    --- Outer Iter 169: norm_f = 604.004, mu=1.55571, |x|=6.35914, |J|=18576.5\n",
      "    --- Outer Iter 170: norm_f = 603.991, mu=1.52017, |x|=6.40373, |J|=18892.8\n",
      "    --- Outer Iter 171: norm_f = 603.978, mu=1.48476, |x|=6.44866, |J|=19217.9\n",
      "    --- Outer Iter 172: norm_f = 603.966, mu=1.44955, |x|=6.49391, |J|=19551.9\n",
      "    --- Outer Iter 173: norm_f = 603.954, mu=1.41462, |x|=6.53949, |J|=19894.8\n",
      "    --- Outer Iter 174: norm_f = 603.942, mu=1.38003, |x|=6.58539, |J|=20246.7\n",
      "    --- Outer Iter 175: norm_f = 603.93, mu=1.34582, |x|=6.63161, |J|=20607.6\n",
      "    --- Outer Iter 176: norm_f = 603.919, mu=1.31203, |x|=6.67814, |J|=20977.5\n",
      "    --- Outer Iter 177: norm_f = 603.908, mu=1.2787, |x|=6.72498, |J|=21356.4\n",
      "    --- Outer Iter 178: norm_f = 603.897, mu=1.24587, |x|=6.77213, |J|=21744.4\n",
      "    --- Outer Iter 179: norm_f = 603.886, mu=1.21354, |x|=6.81956, |J|=22141.5\n",
      "    --- Outer Iter 180: norm_f = 603.876, mu=1.18175, |x|=6.86729, |J|=22547.6\n",
      "    --- Outer Iter 181: norm_f = 603.866, mu=1.1505, |x|=6.91529, |J|=22962.8\n",
      "    --- Outer Iter 182: norm_f = 603.856, mu=1.11981, |x|=6.96356, |J|=23386.9\n",
      "    --- Outer Iter 183: norm_f = 603.847, mu=1.0897, |x|=7.0121, |J|=23820.1\n",
      "    --- Outer Iter 184: norm_f = 603.837, mu=1.06016, |x|=7.06088, |J|=24262.2\n",
      "    --- Outer Iter 185: norm_f = 603.828, mu=1.0312, |x|=7.10992, |J|=24713.3\n",
      "    --- Outer Iter 186: norm_f = 603.819, mu=1.00283, |x|=7.15918, |J|=25173.3\n",
      "    --- Outer Iter 187: norm_f = 603.811, mu=0.97505, |x|=7.20868, |J|=25642.1\n",
      "    --- Outer Iter 188: norm_f = 603.802, mu=0.947853, |x|=7.25838, |J|=26119.7\n",
      "    --- Outer Iter 189: norm_f = 603.794, mu=0.921243, |x|=7.3083, |J|=26606\n",
      "    --- Outer Iter 190: norm_f = 603.786, mu=0.895216, |x|=7.35841, |J|=27101\n",
      "    --- Outer Iter 191: norm_f = 603.779, mu=0.869769, |x|=7.4087, |J|=27604.6\n",
      "    --- Outer Iter 192: norm_f = 603.771, mu=0.844897, |x|=7.45917, |J|=28116.7\n",
      "    --- Outer Iter 193: norm_f = 603.764, mu=0.820596, |x|=7.5098, |J|=28637.2\n",
      "    --- Outer Iter 194: norm_f = 603.757, mu=0.796859, |x|=7.56059, |J|=29166\n",
      "    --- Outer Iter 195: norm_f = 603.75, mu=0.773679, |x|=7.61151, |J|=29703.2\n",
      "    --- Outer Iter 196: norm_f = 603.743, mu=0.75105, |x|=7.66257, |J|=30248.4\n",
      "    --- Outer Iter 197: norm_f = 603.737, mu=0.728965, |x|=7.71374, |J|=30801.7\n",
      "    --- Outer Iter 198: norm_f = 603.731, mu=0.707414, |x|=7.76503, |J|=31362.9\n",
      "    --- Outer Iter 199: norm_f = 603.725, mu=0.686389, |x|=7.81641, |J|=31932\n",
      "    --- Outer Iter 200: norm_f = 603.719, mu=0.665883, |x|=7.86787, |J|=32508.7\n",
      "    --- Outer Iter 201: norm_f = 603.713, mu=0.645885, |x|=7.91941, |J|=33092.9\n",
      "    --- Outer Iter 202: norm_f = 603.707, mu=0.626386, |x|=7.971, |J|=33684.6\n",
      "    --- Outer Iter 203: norm_f = 603.702, mu=0.607378, |x|=8.02265, |J|=34283.6\n",
      "    --- Outer Iter 204: norm_f = 603.697, mu=0.58885, |x|=8.07433, |J|=34889.7\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    --- Outer Iter 205: norm_f = 603.692, mu=0.570793, |x|=8.12603, |J|=35502.8\n",
      "    --- Outer Iter 206: norm_f = 603.687, mu=0.553197, |x|=8.17775, |J|=36122.7\n",
      "    --- Outer Iter 207: norm_f = 603.683, mu=0.536052, |x|=8.22946, |J|=36749.3\n",
      "    --- Outer Iter 208: norm_f = 603.678, mu=0.519349, |x|=8.28117, |J|=37382.4\n",
      "    --- Outer Iter 209: norm_f = 603.674, mu=0.503077, |x|=8.33284, |J|=38021.7\n",
      "    --- Outer Iter 210: norm_f = 603.67, mu=0.487227, |x|=8.38448, |J|=38667.2\n",
      "    --- Outer Iter 211: norm_f = 603.666, mu=0.471789, |x|=8.43607, |J|=39318.6\n",
      "    --- Outer Iter 212: norm_f = 603.662, mu=0.456753, |x|=8.48759, |J|=39975.8\n",
      "    --- Outer Iter 213: norm_f = 603.658, mu=0.44211, |x|=8.53903, |J|=40638.4\n",
      "    --- Outer Iter 214: norm_f = 603.655, mu=0.427849, |x|=8.59039, |J|=41306.4\n",
      "    --- Outer Iter 215: norm_f = 603.651, mu=0.413962, |x|=8.64164, |J|=41979.5\n",
      "    --- Outer Iter 216: norm_f = 603.648, mu=0.400438, |x|=8.69277, |J|=42657.5\n",
      "    --- Outer Iter 217: norm_f = 603.645, mu=0.387269, |x|=8.74378, |J|=43340.1\n",
      "    --- Outer Iter 218: norm_f = 603.642, mu=0.374445, |x|=8.79464, |J|=44027.2\n",
      "    --- Outer Iter 219: norm_f = 603.639, mu=0.361957, |x|=8.84534, |J|=44718.3\n",
      "    --- Outer Iter 220: norm_f = 603.636, mu=0.349796, |x|=8.89588, |J|=45413.5\n",
      "    --- Outer Iter 221: norm_f = 603.634, mu=0.337954, |x|=8.94623, |J|=46112.2\n",
      "    --- Outer Iter 222: norm_f = 603.631, mu=0.326421, |x|=8.99638, |J|=46814.4\n",
      "    --- Outer Iter 223: norm_f = 603.629, mu=0.315189, |x|=9.04632, |J|=47519.7\n",
      "    --- Outer Iter 224: norm_f = 603.626, mu=0.30425, |x|=9.09604, |J|=48227.8\n",
      "    --- Outer Iter 225: norm_f = 603.624, mu=0.293595, |x|=9.14552, |J|=48938.4\n",
      "    --- Outer Iter 226: norm_f = 603.622, mu=0.283216, |x|=9.19474, |J|=49651.4\n",
      "    --- Outer Iter 227: norm_f = 603.62, mu=0.273106, |x|=9.24369, |J|=50366.2\n",
      "    --- Outer Iter 228: norm_f = 603.618, mu=0.263256, |x|=9.29237, |J|=51082.7\n",
      "    --- Outer Iter 229: norm_f = 603.616, mu=0.253659, |x|=9.34075, |J|=51800.6\n",
      "    --- Outer Iter 230: norm_f = 603.614, mu=0.244307, |x|=9.38881, |J|=52519.4\n",
      "    --- Outer Iter 231: norm_f = 603.613, mu=0.235194, |x|=9.43655, |J|=53238.9\n",
      "    --- Outer Iter 232: norm_f = 603.611, mu=0.226311, |x|=9.48395, |J|=53958.7\n",
      "    --- Outer Iter 233: norm_f = 603.61, mu=0.217652, |x|=9.53099, |J|=54678.5\n",
      "    --- Outer Iter 234: norm_f = 603.608, mu=0.209209, |x|=9.57766, |J|=55397.9\n",
      "    --- Outer Iter 235: norm_f = 603.607, mu=0.200977, |x|=9.62394, |J|=56116.5\n",
      "    --- Outer Iter 236: norm_f = 603.606, mu=0.192948, |x|=9.66982, |J|=56834\n",
      "    --- Outer Iter 237: norm_f = 603.605, mu=0.185116, |x|=9.71528, |J|=57549.9\n",
      "    --- Outer Iter 238: norm_f = 603.603, mu=0.177474, |x|=9.76031, |J|=58263.9\n",
      "    --- Outer Iter 239: norm_f = 603.602, mu=0.170017, |x|=9.80488, |J|=58975.6\n",
      "    --- Outer Iter 240: norm_f = 603.601, mu=0.162738, |x|=9.84899, |J|=59684.5\n",
      "    --- Outer Iter 241: norm_f = 603.6, mu=0.155631, |x|=9.89261, |J|=60390.3\n",
      "    --- Outer Iter 242: norm_f = 603.6, mu=0.14869, |x|=9.93573, |J|=61092.4\n",
      "    --- Outer Iter 243: norm_f = 603.599, mu=0.14191, |x|=9.97832, |J|=61790.5\n",
      "    --- Outer Iter 244: norm_f = 603.598, mu=0.135285, |x|=10.0204, |J|=62484\n",
      "    --- Outer Iter 245: norm_f = 603.597, mu=0.128809, |x|=10.0619, |J|=63172.5\n",
      "    --- Outer Iter 246: norm_f = 603.597, mu=0.122477, |x|=10.1028, |J|=63855.5\n",
      "    --- Outer Iter 247: norm_f = 603.596, mu=0.116283, |x|=10.1431, |J|=64532.6\n",
      "    --- Outer Iter 248: norm_f = 603.595, mu=0.110222, |x|=10.1829, |J|=65203.1\n",
      "    --- Outer Iter 249: norm_f = 603.595, mu=0.104289, |x|=10.2219, |J|=65866.5\n",
      "    --- Outer Iter 250: norm_f = 603.594, mu=0.0984776, |x|=10.2604, |J|=66522.4\n",
      "    --- Outer Iter 251: norm_f = 603.594, mu=0.0927839, |x|=10.2981, |J|=67170.1\n",
      "    Least squares message = Both actual and predicted relative reductions in the sum of squares are at most 1e-06\n",
      "  Sum of Chi^2 = 603.594 (280 data params - 65 model params = expected mean of 215; p-value = 0)\n",
      "  Completed in 23.5s\n",
      "  2*Delta(log(L)) = 613.182\n",
      "  Iteration 1 took 23.6s\n",
      "  \n",
      "  Switching to ML objective (last iteration)\n",
      "  --- MLGST ---\n",
      "    --- Outer Iter 0: norm_f = 306.591, mu=1, |x|=10.2981, |J|=47563.3\n",
      "    --- Outer Iter 1: norm_f = 302.819, mu=193418, |x|=10.2982, |J|=47623\n",
      "    --- Outer Iter 2: norm_f = 301.967, mu=64472.5, |x|=10.2982, |J|=47611.6\n",
      "    --- Outer Iter 3: norm_f = 301.659, mu=21490.8, |x|=10.2983, |J|=47635.2\n",
      "    --- Outer Iter 4: norm_f = 301.541, mu=7163.61, |x|=10.2983, |J|=47682.2\n",
      "    --- Outer Iter 5: norm_f = 301.465, mu=2387.87, |x|=10.2983, |J|=47747.8\n",
      "    --- Outer Iter 6: norm_f = 301.407, mu=795.957, |x|=10.2983, |J|=47826\n",
      "    --- Outer Iter 7: norm_f = 301.374, mu=265.319, |x|=10.2983, |J|=47885.2\n",
      "    --- Outer Iter 8: norm_f = 301.344, mu=88.4396, |x|=10.2983, |J|=47898.7\n",
      "    --- Outer Iter 9: norm_f = 301.293, mu=29.4799, |x|=10.2986, |J|=47914.9\n",
      "    --- Outer Iter 10: norm_f = 301.241, mu=15.7421, |x|=10.2992, |J|=48107.2\n",
      "    --- Outer Iter 11: norm_f = 301.217, mu=10.5499, |x|=10.3, |J|=48492.9\n",
      "    --- Outer Iter 12: norm_f = 301.208, mu=5.33888, |x|=10.3003, |J|=48828.2\n",
      "    --- Outer Iter 13: norm_f = 301.206, mu=4.69523, |x|=10.2998, |J|=49119.2\n",
      "    --- Outer Iter 14: norm_f = 301.204, mu=1.56508, |x|=10.2988, |J|=49274.6\n",
      "    --- Outer Iter 15: norm_f = 301.203, mu=0.521692, |x|=10.2946, |J|=49371.3\n",
      "    --- Outer Iter 16: norm_f = 301.202, mu=0.173897, |x|=10.2809, |J|=49291\n",
      "    --- Outer Iter 17: norm_f = 301.202, mu=0.129846, |x|=10.2385, |J|=48801.5\n",
      "    --- Outer Iter 18: norm_f = 301.201, mu=0.129846, |x|=10.1819, |J|=48116.9\n",
      "    --- Outer Iter 19: norm_f = 301.2, mu=0.12793, |x|=10.1256, |J|=47426.7\n",
      "    --- Outer Iter 20: norm_f = 301.199, mu=0.126963, |x|=10.0689, |J|=46731.8\n",
      "    --- Outer Iter 21: norm_f = 301.198, mu=0.126124, |x|=10.0125, |J|=46041.1\n",
      "    --- Outer Iter 22: norm_f = 301.196, mu=0.125355, |x|=9.9567, |J|=45359.9\n",
      "    --- Outer Iter 23: norm_f = 301.195, mu=0.124608, |x|=9.90153, |J|=44690.5\n",
      "    --- Outer Iter 24: norm_f = 301.194, mu=0.123843, |x|=9.84716, |J|=44034.9\n",
      "    --- Outer Iter 25: norm_f = 301.193, mu=0.123041, |x|=9.79365, |J|=43394.2\n",
      "    --- Outer Iter 26: norm_f = 301.192, mu=0.122186, |x|=9.74104, |J|=42768.8\n",
      "    --- Outer Iter 27: norm_f = 301.192, mu=0.121266, |x|=9.68936, |J|=42159.1\n",
      "    --- Outer Iter 28: norm_f = 301.191, mu=0.120272, |x|=9.63863, |J|=41565.1\n",
      "    --- Outer Iter 29: norm_f = 301.19, mu=0.119198, |x|=9.58883, |J|=40986.5\n",
      "    --- Outer Iter 30: norm_f = 301.189, mu=0.118037, |x|=9.53998, |J|=40423.3\n",
      "    --- Outer Iter 31: norm_f = 301.188, mu=0.116786, |x|=9.49206, |J|=39875\n",
      "    --- Outer Iter 32: norm_f = 301.187, mu=0.115443, |x|=9.44507, |J|=39341.4\n",
      "    --- Outer Iter 33: norm_f = 301.187, mu=0.114005, |x|=9.39899, |J|=38822.1\n",
      "    --- Outer Iter 34: norm_f = 301.186, mu=0.112471, |x|=9.35382, |J|=38316.8\n",
      "    --- Outer Iter 35: norm_f = 301.185, mu=0.110842, |x|=9.30953, |J|=37825.1\n",
      "    --- Outer Iter 36: norm_f = 301.185, mu=0.109116, |x|=9.26613, |J|=37346.7\n",
      "    --- Outer Iter 37: norm_f = 301.184, mu=0.107295, |x|=9.2236, |J|=36881.3\n",
      "    --- Outer Iter 38: norm_f = 301.184, mu=0.105378, |x|=9.18193, |J|=36428.5\n",
      "    --- Outer Iter 39: norm_f = 301.183, mu=0.103367, |x|=9.14111, |J|=35988.1\n",
      "    --- Outer Iter 40: norm_f = 301.183, mu=0.101263, |x|=9.10113, |J|=35559.7\n",
      "    --- Outer Iter 41: norm_f = 301.182, mu=0.099066, |x|=9.06198, |J|=35143.1\n",
      "    --- Outer Iter 42: norm_f = 301.182, mu=0.0967775, |x|=9.02366, |J|=34738.1\n",
      "    --- Outer Iter 43: norm_f = 301.181, mu=0.0943984, |x|=8.98616, |J|=34344.4\n",
      "    --- Outer Iter 44: norm_f = 301.181, mu=0.0919296, |x|=8.94947, |J|=33961.7\n",
      "    --- Outer Iter 45: norm_f = 301.18, mu=0.0893721, |x|=8.91359, |J|=33589.9\n",
      "    --- Outer Iter 46: norm_f = 301.18, mu=0.0867267, |x|=8.87852, |J|=33228.8\n",
      "    --- Outer Iter 47: norm_f = 301.18, mu=0.0839945, |x|=8.84424, |J|=32878\n",
      "    --- Outer Iter 48: norm_f = 301.179, mu=0.0811763, |x|=8.81076, |J|=32537.6\n",
      "    --- Outer Iter 49: norm_f = 301.179, mu=0.0782729, |x|=8.77807, |J|=32207.2\n",
      "    --- Outer Iter 50: norm_f = 301.179, mu=0.0752852, |x|=8.74618, |J|=31886.8\n",
      "    --- Outer Iter 51: norm_f = 301.179, mu=0.0722142, |x|=8.71507, |J|=31576.2\n",
      "    --- Outer Iter 52: norm_f = 301.178, mu=0.0690606, |x|=8.68476, |J|=31275.2\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    --- Outer Iter 53: norm_f = 301.178, mu=0.0658253, |x|=8.65523, |J|=30983.7\n",
      "    Least squares message = Both actual and predicted relative reductions in the sum of squares are at most 1e-06\n",
      "    Maximum log(L) = 301.178 below upper bound of -421916\n",
      "      2*Delta(log(L)) = 602.356 (280 data params - 65 model params = expected mean of 215; p-value = 0)\n",
      "    Completed in 5.1s\n",
      "  2*Delta(log(L)) = 602.356\n",
      "  Final MLGST took 5.1s\n",
      "  \n",
      "Iterative MLGST Total Time: 28.7s\n",
      "    -- Performing 'go0' gauge optimization on default estimate --\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/algorithms/gaugeopt.py:242: UserWarning:\n",
      "\n",
      "No gauge group specified, so no gauge optimization performed.\n",
      "\n"
     ]
    }
   ],
   "source": [
    "results_2plus1 = pygsti.do_long_sequence_gst(ds, mdl_2plus1_ideal, prepfids, measfids,\n",
    "                                             germs, maxLengths, verbosity=3,\n",
    "                                             advancedOptions={\"starting point\": \"target\"})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Running idle tomography\n",
      "Computing switchable properties\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1236: UserWarning:\n",
      "\n",
      "Idle tomography failed:\n",
      "operands could not be broadcast together with shapes (9,) (16,) \n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Found standard clifford compilation from std1Q_XYI\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/enielse/pyGSTi/pygsti/report/factory.py:1303: UserWarning:\n",
      "\n",
      "PDF output is not available for reports with multiple datasets\n",
      "\n",
      "/Users/enielse/pyGSTi/pygsti/tools/matrixtools.py:613: UserWarning:\n",
      "\n",
      "Cannot construct a real log: unpaired negative real eigenvalues: [(-1.0084172062338246+0j)]\n",
      "\n",
      "/Users/enielse/pyGSTi/pygsti/report/reportables.py:1090: UserWarning:\n",
      "\n",
      "Truncating imaginary logarithm!\n",
      "\n",
      "/Users/enielse/pyGSTi/pygsti/tools/optools.py:1239: UserWarning:\n",
      "\n",
      "Near-identity matrix log failed; falling back to approximate log for logGTi error generator:\n",
      "Failed to construct a real logarithm! This is probably because M is not near the identity.\n",
      "Its eigenvalues are: [-1.00841721+0.j          0.01252335+0.97927003j  0.01252335-0.97927003j\n",
      "  1.00037875+0.j          1.00405493+0.j        ]\n",
      "\n",
      "/Users/enielse/pyGSTi/pygsti/tools/optools.py:1243: UserWarning:\n",
      "\n",
      "Falling back to approximate log for logGTi error generator\n",
      "\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n",
      "Statistical hypothesis tests did NOT find inconsistency between the datasets at 5.00% significance.\n"
     ]
    }
   ],
   "source": [
    "pygsti.report.construct_standard_report(\n",
    "    {'two-level': results_2level, 'three-level': results_3level_leakage_basis,\n",
    "     'two+one level': results_2plus1},\n",
    "    \"Leakage Example Report\"\n",
    ").write_html('example_files/leakage_report', autosize='none')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Open the report [here](example_files/leakage_report/main.html)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}