diff --git a/Jupyter_notebook/KB_Exercises_with_solutions.ipynb b/Jupyter_notebook/KB_Exercises_with_solutions.ipynb
new file mode 100644
index 0000000..3eda468
--- /dev/null
+++ b/Jupyter_notebook/KB_Exercises_with_solutions.ipynb
@@ -0,0 +1,470 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "21e4894b",
+   "metadata": {},
+   "source": [
+    "Julia language\n",
+    "\n",
+    "## Exercise 1 ( 10 min )\n",
+    "\n",
+    "(a) write a function that constructs a quantum `state` in z-basis for N qbit system in vector format.\n",
+    "\n",
+    "Input : string of bits \"..000101101..\"\n",
+    "\n",
+    "(b) write a function that constructs a quantum `Hamiltonian` for N qbit system in matrix format\n",
+    "\n",
+    "$H = \\sum n_i n_{i+1} - \\mu\\sum n_i$\n",
+    "\n",
+    "where $n = \\begin{bmatrix}\n",
+    "           1 & 0 \\\\\n",
+    "           0 & 0 \\\\\n",
+    "         \\end{bmatrix}$\n",
+    "         \n",
+    "Input: number of qbits N and chemical potential $\\mu$ \n",
+    "\n",
+    "\n",
+    "hint: one can apply kronecker product iteratively with `foldl(⊗, container)`\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "fe85586b",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "Ham (generic function with 1 method)"
+      ]
+     },
+     "execution_count": 1,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "⊗(x, y) = kron(x, y)\n",
+    "\n",
+    "\n",
+    "function psi(bitstring)\n",
+    "    up = [1, 0]\n",
+    "    down = [0, 1]\n",
+    "    state = fill(up, length(bitstring))\n",
+    "    for (i, char) in enumerate(bitstring)\n",
+    "        if char == '0'\n",
+    "            state[i] = down\n",
+    "        end\n",
+    "    end \n",
+    "    return foldl(⊗, state) #reduce \n",
+    "end\n",
+    "\n",
+    "\n",
+    "\n",
+    "function Ham(N, μ)  # μ is \\mu<TAB>\n",
+    "    id = [\n",
+    "        1 0;0 1]\n",
+    "    ProjUp = [1 0; 0 0]\n",
+    "    H = zeros(2^N, 2^N)\n",
+    "    \n",
+    "    for i in 1: N\n",
+    "        one_body = fill(id, N)\n",
+    "        one_body[i] = ProjUp \n",
+    "        H += -μ * foldl(⊗, one_body)\n",
+    "    end\n",
+    "    \n",
+    "    for i in 1: N-1\n",
+    "        two_body = fill(id, N)\n",
+    "        two_body[i] = ProjUp\n",
+    "        two_body[i+1] = ProjUp\n",
+    "        H += foldl(⊗, two_body)\n",
+    "    end\n",
+    "    return H\n",
+    "end\n",
+    "\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "f8f901ad",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "4×4 Matrix{Float64}:\n",
+       " 0.8   0.0   0.0  0.0\n",
+       " 0.0  -0.1   0.0  0.0\n",
+       " 0.0   0.0  -0.1  0.0\n",
+       " 0.0   0.0   0.0  0.0"
+      ]
+     },
+     "execution_count": 2,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "psi(\"11\")\n",
+    "Ham(2, 0.1)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e2b00e64",
+   "metadata": {},
+   "source": [
+    "## Exercise 2 ( 10 min )\n",
+    "\n",
+    "write a function that returns a product state(s) with minimal energy, i.e. \n",
+    "\n",
+    "$|\\psi\\rangle : \\langle \\psi |H(N, \\mu)| \\psi \\rangle = E_{min}$\n",
+    "\n",
+    "among z-basis states (\"..000101101..\").\n",
+    "\n",
+    "$H = \\sum n_i n_{i+1} - \\mu\\sum n_i$\n",
+    "\n",
+    "\n",
+    "answer should have a following format `Pair(\"..0100100..\", psi)`\n",
+    "\n",
+    "hint: function that translates integer to bitsting is `string(i; base = 2, pad = string_length)`\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "id": "1f493ce6",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "psi_min (generic function with 1 method)"
+      ]
+     },
+     "execution_count": 3,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "function psi_min(H)\n",
+    "    s=size(H, 1) #size along dim 1\n",
+    "    N = Int(log2(s))\n",
+    "    \n",
+    "    #Number = {Int, Int32, Int64, Float, Float32, Float64... }\n",
+    "    \n",
+    "    energy = Dict{String, Number}()\n",
+    "    for i in 0 : 2^N-1 \n",
+    "        bstr = string(i; base = 2, pad = N)\n",
+    "        ψ = psi(bstr)\n",
+    "        energy[bstr] = (ψ' * H * ψ)\n",
+    "    end\n",
+    "    min_en = findmin(energy)[1] #finding a min\n",
+    "    return filter(elem->(elem.second == min_en), energy) #filter out enrgies ≠ min\n",
+    "end"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "id": "03c40ad7",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "  0.298671 seconds (96.02 k allocations: 479.650 MiB, 4.96% gc time, 7.43% compilation time)\n"
+     ]
+    },
+    {
+     "data": {
+      "text/plain": [
+       "Dict{String, Number} with 6 entries:\n",
+       "  \"1010100101\" => -2.0\n",
+       "  \"1001010101\" => -2.0\n",
+       "  \"0101010101\" => -2.0\n",
+       "  \"1010101010\" => -2.0\n",
+       "  \"1010101001\" => -2.0\n",
+       "  \"1010010101\" => -2.0"
+      ]
+     },
+     "execution_count": 4,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "@time psi_min(Ham(10, 0.4))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "00618a92",
+   "metadata": {},
+   "source": [
+    "## Exercise 3 ( 5 min )\n",
+    "\n",
+    "Write `time_evolution` function, that evolves a state in time from 0 to `T` with time step `dt`. \n",
+    "Evolution operator is $e^{-iHt}$\n",
+    "\n",
+    "\n",
+    "hint: its much cheaper to compute $e^{At} * v$ then $e^{At}$ and apply it to $v$.\n",
+    "\n",
+    "hint: one can test `time_evolution` implementation with complex random normalized vector of size `2^N`\n",
+    "\n",
+    "Use `exponentiate` by `using KrylovKit` package, which returns updated vector and convergence info `v, info`"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "id": "9626195b",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "time_evolution (generic function with 1 method)"
+      ]
+     },
+     "execution_count": 5,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "using KrylovKit: exponentiate\n",
+    "\n",
+    "function time_evolution(ψ, H, T, dt)\n",
+    "    n = Int(T/dt)\n",
+    "    for _ in 1:n #pay attention!\n",
+    "    #in Julia ranges are closed, i.e 1:5 means 1,2,3,4,5\n",
+    "        ψ, info = exponentiate(H, -dt * 1.0im, ψ)\n",
+    "    end\n",
+    "    return ψ\n",
+    "end\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 10,
+   "id": "9fe70a9a",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "32-element Vector{ComplexF64}:\n",
+       "  -0.11369655122605653 + 0.0644759408528872im\n",
+       "   0.18349136743927763 - 0.018804816011487188im\n",
+       "    0.1467577801718119 - 0.028354675851863194im\n",
+       "  0.026688228912223316 - 0.1124675819945794im\n",
+       "   0.15798822035534396 - 0.15212888374558808im\n",
+       "  -0.21275316455822024 - 0.023731124236512875im\n",
+       "  -0.12136589943180677 - 0.06088868446034836im\n",
+       "  -0.13852381519551052 - 0.0691409971171647im\n",
+       "   0.01834366369626367 - 0.11974752725919376im\n",
+       "  -0.07862415197345206 - 0.0635506699095704im\n",
+       "  -0.14791132144988767 + 0.19550749202417048im\n",
+       "  -0.08707132251266379 + 0.2099442167096762im\n",
+       "  -0.20052000333458775 - 0.06848655240761178im\n",
+       "                       ⋮\n",
+       "  -0.18249832705793287 - 0.10218628862041482im\n",
+       "   0.06809470175546742 + 0.1264850558212589im\n",
+       "   -0.0620675827036714 + 0.1621622664192347im\n",
+       " -0.020759338229685908 + 0.07421231764004241im\n",
+       "   0.10253211481695645 - 0.10792728617841092im\n",
+       "  -0.19711912650548005 - 0.07048267356892052im\n",
+       "  0.030117474287823014 + 0.12083556505161931im\n",
+       "   0.05088814494581848 + 0.039914978154047355im\n",
+       "   -0.1474374927031578 + 0.1302945516514859im\n",
+       "   0.04385013425934344 + 0.17969074529162735im\n",
+       "   0.09643366479407968 + 0.2445407924259601im\n",
+       "  0.021474246565776667 + 0.17810974803718116im"
+      ]
+     },
+     "execution_count": 10,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "# Usage example\n",
+    "\n",
+    "using LinearAlgebra\n",
+    "\n",
+    "L = 5\n",
+    "ψ = rand(Complex{Float64}, 2^L)\n",
+    "ψ /= norm(ψ,2)\n",
+    "\n",
+    "time_evolution(ψ, Ham(L,0.1), 5.0, 0.1)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "801e179a",
+   "metadata": {},
+   "source": [
+    "## Exercise 4 ( 10 min )\n",
+    "\n",
+    "compare execution time of `time_evolution` function (Ex 3) for different types of Hamiltonians `H(N = 14, μ = 0.5)`\n",
+    "\n",
+    "1. dense\n",
+    "2. sparse `using SparseArrays`\n",
+    "3. OPTIONAL : cuda dense (for GPU)\n",
+    "\n",
+    "hint: `CUDA` package has a `cu()` function that translates arrays to GPU.\n",
+    "    \n",
+    "\n",
+    "execution time can be measured with macro @time (or better @btime from `BenchmarkTools`)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 18,
+   "id": "640824a9",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      " 11.944268 seconds (4.55 k allocations: 118.975 MiB)\n",
+      "  1.511705 seconds (4.65 k allocations: 120.543 MiB)\n"
+     ]
+    },
+    {
+     "data": {
+      "text/plain": [
+       "4096-element Vector{ComplexF64}:\n",
+       "  -0.011853896873262104 + 0.01259998223099739im\n",
+       "   0.016679847637844233 + 0.004370082560263621im\n",
+       "   0.002114974256463454 - 0.0071131515450265324im\n",
+       " -0.0032054591819539137 - 0.00926459570257121im\n",
+       "  0.0011641368741503472 - 0.005278561585997654im\n",
+       "  -0.009475273607725068 - 0.005131796471141815im\n",
+       "  -0.022962912845868162 - 0.00540743868054884im\n",
+       "  -0.006543091664889295 + 0.0008557358944384494im\n",
+       "   0.004189808193965449 - 0.005722135568732727im\n",
+       "  -0.014262943252438908 - 0.013317884178562809im\n",
+       "  -0.006759735280905994 + 0.006290766876914184im\n",
+       " -0.0023149785963365744 + 0.024138828713581504im\n",
+       "   -0.01298834674972344 - 0.0033500004175022738im\n",
+       "                        ⋮\n",
+       "   -0.01402220755017167 - 0.012596739847918991im\n",
+       "   -0.00611078532243405 + 0.021092781754763045im\n",
+       "   0.005161613712920969 + 0.010110514799360763im\n",
+       "  -0.004415839482049863 + 0.013229135261676957im\n",
+       "  0.0015452517633196218 - 0.002950560010065211im\n",
+       "  -0.003299943671265516 + 0.0011637286724671902im\n",
+       " 0.00010975421394683415 + 0.002950189618318821im\n",
+       "   0.005115835436625776 + 0.023144095840165742im\n",
+       "  -0.007788043334072698 + 0.0018559592284559017im\n",
+       "  0.0053393551846635775 + 0.01349703943774104im\n",
+       "  0.0004950099364579647 + 0.014086862352750253im\n",
+       "   0.017059820216028993 + 0.006006884844123996im"
+      ]
+     },
+     "execution_count": 18,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "using SparseArrays\n",
+    "using BenchmarkTools\n",
+    "using LinearAlgebra\n",
+    "\n",
+    "\n",
+    "L = 12\n",
+    "dt = 0.1\n",
+    "T = 5.0\n",
+    "H = Ham(L, 0.1)\n",
+    "Hsp = sparse(H)\n",
+    "\n",
+    "ψ = rand(Complex{Float64}, 2^L)\n",
+    "ψ /= norm(ψ,2)\n",
+    "\n",
+    "@time time_evolution(ψ, H, T, dt)\n",
+    "@time time_evolution(ψ, Hsp, T, dt)\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 14,
+   "id": "4cf74f03",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/plain": [
+       "1024×1024 SparseMatrixCSC{Float64, Int64} with 1023 stored entries:\n",
+       "⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄⠀⠀\n",
+       "⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠀⠑⢄"
+      ]
+     },
+     "execution_count": 14,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "using SparseArrays\n",
+    "using BenchmarkTools\n",
+    "using LinearAlgebra\n",
+    "\n",
+    "L = 10\n",
+    "dt = 0.1\n",
+    "T = 5.0\n",
+    "H = Ham(L, 0.1)\n",
+    "Hsp = sparse(H)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "af7269bd",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Julia 1.8.2",
+   "language": "julia",
+   "name": "julia-1.8"
+  },
+  "language_info": {
+   "file_extension": ".jl",
+   "mimetype": "application/julia",
+   "name": "julia",
+   "version": "1.8.2"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}