From c6d31b809f5b55bdd1fe459105f3f01963519364 Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Mon, 25 Nov 2019 18:44:46 +0100 Subject: [PATCH] Added SpVec notebook --- Notebooks/SpVec.ipynb | 1308 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 1308 insertions(+) create mode 100644 Notebooks/SpVec.ipynb diff --git a/Notebooks/SpVec.ipynb b/Notebooks/SpVec.ipynb new file mode 100644 index 0000000..3e0623c --- /dev/null +++ b/Notebooks/SpVec.ipynb @@ -0,0 +1,1308 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "collapsed": true + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "- : unit = ()\n", + "Findlib has been successfully loaded. Additional directives:\n", + " #require \"package\";; to load a package\n", + " #list;; to list the available packages\n", + " #camlp4o;; to load camlp4 (standard syntax)\n", + " #camlp4r;; to load camlp4 (revised syntax)\n", + " #predicates \"p,q,...\";; to set these predicates\n", + " Topfind.reset();; to force that packages will be reloaded\n", + " #thread;; to enable threads\n", + "\n", + "- : unit = ()\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + "/home/scemama/qp2/external/opam/default/lib/bytes: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/base64: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/base64/base64.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/compiler-libs: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/compiler-libs/ocamlcommon.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/result: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/result/result.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ppx_deriving/runtime: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/ppx_deriving/runtime/ppx_deriving_runtime.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ppx_deriving_yojson/runtime: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/ppx_deriving_yojson/runtime/ppx_deriving_yojson_runtime.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/unix.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/uuidm: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/uuidm/uuidm.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/easy-format: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/easy-format/easy_format.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/biniou: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/biniou/biniou.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/yojson: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/yojson/yojson.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/jupyter: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/jupyter/jupyter.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/jupyter/notebook: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/jupyter/notebook/jupyter_notebook.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/compiler-libs/ocamlbytecomp.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/compiler-libs/ocamltoplevel.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/ocaml/bigarray.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/lacaml: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/lacaml/lacaml.cma: loaded\n", + "/home/scemama/qp2/external/opam/default/lib/lacaml/top: added to search path\n", + "/home/scemama/qp2/external/opam/default/lib/lacaml/top/lacaml_top.cma: loaded\n" + ] + } + ], + "source": [ + "#use \"topfind\";;\n", + "#require \"jupyter.notebook\";;\n", + "#require \"lacaml.top\" ;;" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Sparse Vector module\n", + "----------------" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "A sparse vector is a structure made of:\n", + "* The dimension of the vector space\n", + "* The number of non-zeros\n", + "* An array of indices\n", + "* An array of values\n", + "\n", + "The indices are stored in an ``int Bigarray`` and the values are stored in a ``Lacaml.Vec.t``\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Types" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "module L = Lacaml.D\n" + ] + }, + "execution_count": 2, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "type t = {\n", + " dim : int;\n", + " nnz : int;\n", + " indices :\n", + " (int, Bigarray.int_elt, Bigarray.fortran_layout) Bigarray.Array1.t;\n", + " values : L.Vec.t;\n", + "}\n" + ] + }, + "execution_count": 2, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "module L = Lacaml.D\n", + "\n", + "type t =\n", + " {\n", + " dim: int ;\n", + " nnz: int ;\n", + " indices: (int, Bigarray.int_elt, Bigarray.fortran_layout) Bigarray.Array1.t ; (* Indices *)\n", + " values: L.Vec.t\n", + " }\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Printers" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val pp : Format.formatter -> t -> unit = \n" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let pp ppf t =\n", + " let pp_data ppf t =\n", + " for i=1 to t.nnz do\n", + " Format.fprintf ppf \"@[(%d,@ %f)@]@;\" t.indices.{i} t.values.{i}\n", + " done\n", + " in\n", + " Format.fprintf ppf \"@[{@[dim:@ %d@]@;@[%a@]}@]\" t.dim pp_data t" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Creators" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val of_vec : ?threshold:float -> Lacaml.D.vec -> t = \n" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let of_vec ?(threshold=0.) v =\n", + " let dim = L.Vec.dim v in\n", + " let buffer_idx = Bigarray.(Array1.create int fortran_layout) dim in\n", + " let buffer_val = Bigarray.(Array1.create float64 fortran_layout) dim in\n", + " let check = \n", + " if threshold = 0. then\n", + " fun x -> x <> 0.\n", + " else\n", + " fun x -> (abs_float x) > 0.\n", + " in\n", + " let rec aux k i =\n", + " if i > dim then\n", + " k-1\n", + " else if check v.{i} then\n", + " ( buffer_idx.{k} <- i ;\n", + " buffer_val.{k} <- v.{i} ;\n", + " aux (k+1) (i+1)\n", + " )\n", + " else\n", + " aux k (i+1)\n", + " in\n", + " let nnz = aux 1 1 in\n", + " let indices = Bigarray.(Array1.create int fortran_layout) nnz in\n", + " let values = L.Vec.create nnz in\n", + " for i=1 to nnz do\n", + " indices.{i} <- buffer_idx.{i};\n", + " values.{i} <- buffer_val.{i};\n", + " done ;\n", + " { dim ; nnz ; indices ; values }\n", + "\n" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val make0 : int -> t = \n" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let make0 dim =\n", + " { dim ;\n", + " nnz = 0;\n", + " indices = Bigarray.(Array1.create int fortran_layout) 32 ;\n", + " values = L.Vec.create 32;\n", + " }" + ] + }, + { + "cell_type": "code", + "execution_count": 77, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val of_vec : ?threshold:float -> Lacaml.D.vec -> t = \n" + ] + }, + "execution_count": 77, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val of_array : ?threshold:float -> float array -> t = \n" + ] + }, + "execution_count": 77, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let of_vec ?(threshold=0.) v =\n", + " let dim = L.Vec.dim v in\n", + " let buffer_idx = Bigarray.(Array1.create int fortran_layout) dim in\n", + " let buffer_val = Bigarray.(Array1.create float64 fortran_layout) dim in\n", + " let check = \n", + " if threshold = 0. then\n", + " fun x -> x <> 0.\n", + " else\n", + " fun x -> (abs_float x) > 0.\n", + " in\n", + " let rec aux k i =\n", + " if i > dim then\n", + " k-1\n", + " else if check v.{i} then\n", + " ( buffer_idx.{k} <- i ;\n", + " buffer_val.{k} <- v.{i} ;\n", + " aux (k+1) (i+1)\n", + " )\n", + " else\n", + " aux k (i+1)\n", + " in\n", + " let nnz = aux 1 1 in\n", + " let indices = Bigarray.(Array1.create int fortran_layout) nnz in\n", + " let values = L.Vec.create nnz in\n", + " for i=1 to nnz do\n", + " indices.{i} <- buffer_idx.{i};\n", + " values.{i} <- buffer_val.{i};\n", + " done ;\n", + " { dim ; nnz ; indices ; values }\n", + "\n", + "\n", + "let of_array ?(threshold=0.) a =\n", + " L.Vec.of_array a\n", + " |> of_vec ~threshold " + ] + }, + { + "cell_type": "code", + "execution_count": 78, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val copy : t -> t = \n" + ] + }, + "execution_count": 78, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let copy t =\n", + " let indices = \n", + " Bigarray.(Array1.create int fortran_layout) t.nnz\n", + " in\n", + " Bigarray.Array1.blit t.indices indices ;\n", + " let values = L.copy t.values in\n", + " { dim = t.dim ;\n", + " nnz = t.nnz ;\n", + " indices ; values }" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Test" + ] + }, + { + "cell_type": "code", + "execution_count": 79, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val x : t =\n", + " {dim = 10; nnz = 0; indices = ;\n", + " values =\n", + " R1 R2 R3 R30 R31 R32\n", + " 6.9331E-310 6.9331E-310 0 ... 6.9331E-310 6.9331E-310 6.9331E-310}\n" + ] + }, + "execution_count": 79, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val dense_a : t =\n", + " {dim = 9; nnz = 5; indices = ;\n", + " values = R1 R2 R3 R4 R5\n", + " 1 -2 0.5 1E-08 3}\n" + ] + }, + "execution_count": 79, + "metadata": {}, + "output_type": "execute_result" + }, + { + "ename": "error", + "evalue": "compile_error", + "output_type": "error", + "traceback": [ + "\u001b[32mFile \"[79]\", line 6, characters 22-29:\n\u001b[31mError: This expression has type t but an expression was expected of type\n Lacaml.D.vec =\n (float, Bigarray.float64_elt, Bigarray.fortran_layout)\n Bigarray.Array1.t\n\u001b[36m 5: \u001b[30m \n\u001b[36m 6: \u001b[30mlet sparse_a = of_vec \u001b[4mdense_a\u001b[0m\u001b[30m \n\u001b[36m 7: \u001b[30m\u001b[0m\n" + ] + } + ], + "source": [ + "let x = make0 10 \n", + "\n", + "let dense_a = \n", + " of_array [| 1. ; -2. ; 0. ; 0. ; 0.5 ; 1.e-8 ; 0. ; 3. ; 0. |]\n", + " \n", + "let sparse_a = of_vec dense_a \n", + "\n", + "let _ =\n", + " copy sparse_a = sparse_a\n", + " \n", + "let _ =\n", + " copy sparse_a == sparse_a\n", + "\n", + "let () = \n", + " Format.printf \"@.@[%a@]@.\" pp sparse_a" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# aX + Y\n", + "\n", + "Run along all the entries of `X` and `Y` simultaneously with indices `k` and `l`. `m` is the index of the new array.\n", + "\n", + "if `kl`, update using `y[l]`.\n", + "\n", + "if `k=l`, update using `a*x[k] + y[l]`." + ] + }, + { + "cell_type": "code", + "execution_count": 105, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val axpy : ?threshold:float -> ?alpha:float -> t -> t -> t = \n" + ] + }, + "execution_count": 105, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let axpy ?(threshold=0.) ?(alpha=1.) x y =\n", + "\n", + " if dim x <> dim y then \n", + " invalid_arg \"Inconsistent dimensions\";\n", + "\n", + " let check = (* Test if value should be added wrt threshold *)\n", + " if threshold = 0. then\n", + " fun x -> x <> 0.\n", + " else\n", + " fun x -> abs_float x > 0.\n", + " in\n", + " \n", + " let f = (* if a=1 in ax+y, then do x+y. If a=0 then do y *)\n", + " if alpha = 1. then\n", + " fun x y -> x +. y\n", + " else if alpha = 0. then\n", + " fun _ y -> y\n", + " else\n", + " fun x y -> alpha *. x +. y\n", + " in\n", + " \n", + " let dim = dim x in\n", + " let nnz = x.nnz + y.nnz in\n", + " let new_indices = Bigarray.(Array1.create int fortran_layout) nnz in\n", + " let new_values = L.Vec.create nnz in\n", + " \n", + " let rec aux k l m =\n", + " match k <= x.nnz, l <= y.nnz with\n", + " | true , true -> (* Both arrays are running *)\n", + " begin\n", + " if x.indices.{k} < y.indices.{l} then (\n", + " let w = f x.values.{k} 0. in\n", + " if check w then (\n", + " new_indices.{m} <- x.indices.{k};\n", + " new_values.{m} <- w\n", + " );\n", + " (aux [@tailcall]) (k+1) l (m+1)\n", + " )\n", + " else if x.indices.{k} > y.indices.{l} then (\n", + " let w = y.values.{l} in\n", + " if check w then (\n", + " new_indices.{m} <- y.indices.{l};\n", + " new_values.{m} <- w\n", + " );\n", + " (aux [@tailcall]) k (l+1) (m+1)\n", + " )\n", + " else (\n", + " let w = f x.values.{k} y.values.{l} in\n", + " if check w then (\n", + " new_indices.{m} <- x.indices.{k};\n", + " new_values.{m} <- w\n", + " );\n", + " (aux [@tailcall]) (k+1) (l+1) (m+1)\n", + " )\n", + " end\n", + " | false, true -> (* Array x is done running *)\n", + " begin\n", + " let m = ref m in\n", + " for i=l to y.nnz do\n", + " let w = y.values.{i} in\n", + " if check w then (\n", + " new_indices.{!m} <- y.indices.{i};\n", + " new_values.{!m} <- w;\n", + " incr m;\n", + " )\n", + " done; !m\n", + " end\n", + " | true, false -> (* Array y is done running *)\n", + " begin\n", + " let m = ref m in\n", + " for i=k to x.nnz do\n", + " let w = alpha *. x.values.{i} in\n", + " if check w then (\n", + " new_indices.{!m} <- x.indices.{i};\n", + " new_values.{!m} <- w;\n", + " incr m;\n", + " )\n", + " done; !m\n", + " end\n", + " | false, false -> (* Both arrays are done *)\n", + " m\n", + " in\n", + " let nnz = (aux 1 1 1) - 1 in\n", + " { dim ; nnz ;\n", + " indices = new_indices ; values = new_values }\n", + " " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Test" + ] + }, + { + "cell_type": "code", + "execution_count": 108, + "metadata": { + "collapsed": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "val m_A : Lacaml.D.vec = R1 R2 R3 R8 R9 R10\n", + " 1 2 0 ... -0.001 0 0\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_B : Lacaml.D.vec = R1 R2 R3 R8 R9 R10\n", + " 0 1 2 ... 0 -0.001 2\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_As : t =\n", + " {dim = 10; nnz = 5; indices = ;\n", + " values = R1 R2 R3 R4 R5\n", + " 1 2 0.01 -2 -0.001}\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_Bs : t =\n", + " {dim = 10; nnz = 6; indices = ;\n", + " values = R1 R2 R3 R4 R5 R6\n", + " 1 2 0.01 -2 -0.001 2}\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_C : L.vec = R1 R2 R3 R8 R9 R10\n", + " 0 1 2 ... 0 -0.001 2\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : L.vec = R1 R2 R3 R8 R9 R10\n", + " 2 5 2 ... -0.002 -0.001 2\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_D : L.vec = R1 R2 R3 R8 R9 R10\n", + " 1 2 0 ... -0.001 0 0\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : L.vec = R1 R2 R3 R8 R9 R10\n", + " 1 4 4 ... -0.001 -0.002 4\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_Cs : t =\n", + " {dim = 10; nnz = 9; indices = ;\n", + " values =\n", + " R1 R2 R3 R9 R10 R11\n", + " 2 5 2 ... 2 1.04136E-71 8.95166E+271}\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : bool = false\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val m_Ds : t =\n", + " {dim = 10; nnz = 9; indices = ;\n", + " values =\n", + " R1 R2 R3 R9 R10 R11\n", + " 1 4 4 ... 4 6.82475E-38 8.26993E-72}\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : bool = false\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "1 2.000000 2.000000\n", + "2 5.000000 5.000000\n", + "3 2.000000 2.000000\n", + "4 0.000000 0.000000\n", + "5 0.020000 0.020000\n", + "6 -3.990000 -3.990000\n", + "7 -2.000000 -2.000000\n", + "8 -0.002000 -0.002000\n", + "9 -0.001000 -0.001000\n", + "10 2.000000 2.000000\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "1 1.000000 1.000000\n", + "2 4.000000 4.000000\n", + "3 4.000000 4.000000\n", + "4 0.000000 0.000000\n", + "5 0.010000 0.010000\n", + "6 -1.980000 -1.980000\n", + "7 -4.000000 -4.000000\n", + "8 -0.001000 -0.001000\n", + "9 -0.002000 -0.002000\n", + "10 4.000000 4.000000\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 108, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let m_A = L.Vec.of_array [| 1. ; 2. ; 0. ; 0. ; 0.01 ; -2. ; 0. ; -1.e-3 ; 0. ; 0.|]\n", + "let m_B = L.Vec.of_array [| 0. ; 1. ; 2. ; 0. ; 0. ; 0.01 ; -2. ; 0. ; -1.e-3 ; 2. |]\n", + "\n", + "let m_As = of_vec m_A\n", + "let m_Bs = of_vec m_B\n", + "\n", + "let m_C = L.copy m_B ;;\n", + "L.axpy ~alpha:2. m_A m_C;;\n", + "m_C;;\n", + "\n", + "let m_D = L.copy m_A;;\n", + "L.axpy ~alpha:2. m_B m_D;;\n", + "m_D;;\n", + "\n", + "let m_Cs = axpy ~alpha:2. m_As m_Bs\n", + "\n", + "\n", + "let _ = of_vec m_C = m_Cs;;\n", + "\n", + "let m_Ds = axpy ~alpha:2. m_Bs m_As\n", + "\n", + "let _ = of_vec m_D = m_Ds;;\n", + "\n", + "L.Vec.iteri (fun i x -> Format.printf \"%d %f %f\\n%!\" i x (get m_Cs i)) m_C;;\n", + "L.Vec.iteri (fun i x -> Format.printf \"%d %f %f\\n%!\" i x (get m_Ds i)) m_D;;\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Accessors" + ] + }, + { + "cell_type": "code", + "execution_count": 24, + "metadata": { + "collapsed": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "val dim : t -> int = \n" + ] + }, + "execution_count": 24, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val nnz : t -> int = \n" + ] + }, + "execution_count": 24, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val indices :\n", + " t -> (int, Bigarray.int_elt, Bigarray.fortran_layout) Bigarray.Array1.t =\n", + " \n" + ] + }, + "execution_count": 24, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val values : t -> L.Vec.t = \n" + ] + }, + "execution_count": 24, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let dim t = t.dim\n", + "let nnz t = t.nnz\n", + "let indices t = t.indices\n", + "let values t = t.values" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val get : t -> int -> float = \n" + ] + }, + "execution_count": 25, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let get t i =\n", + " if i < 1 || i > dim t then invalid_arg \"index out of bounds\";\n", + " \n", + " let rec binary_search index value low high =\n", + " if high = low then\n", + " if index.{low} = value then\n", + " low\n", + " else\n", + " raise Not_found\n", + " else let mid = (low + high) / 2 in\n", + " if index.{mid} > value then\n", + " binary_search index value low (mid - 1)\n", + " else if index.{mid} < value then\n", + " binary_search index value (mid + 1) high\n", + " else\n", + " mid\n", + " in\n", + " try\n", + " let k = \n", + " let id = indices t in\n", + " binary_search id i id.{1} (nnz t)\n", + " in\n", + " t.values.{k}\n", + " with Not_found -> 0." + ] + }, + { + "cell_type": "code", + "execution_count": 26, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val iter : (int -> float -> 'a) -> t -> unit = \n" + ] + }, + "execution_count": 26, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let iter f t =\n", + " for k=1 to nnz t do\n", + " f t.indices.{k} t.values.{k}\n", + " done" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Test " + ] + }, + { + "cell_type": "code", + "execution_count": 28, + "metadata": { + "collapsed": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "- : bool = true\n" + ] + }, + "execution_count": 28, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "1 1.000000\n", + "2 -2.000000\n", + "5 0.500000\n", + "6 0.000000\n", + "8 3.000000\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 28, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "dense_a = L.Vec.init (dim sparse_a) (get sparse_a);;\n", + "iter (fun i v -> Printf.printf \"%d %f\\n%!\" i v) sparse_a;;" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Converters" + ] + }, + { + "cell_type": "code", + "execution_count": 30, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val to_assoc_list : t -> (int * float) list = \n" + ] + }, + "execution_count": 30, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val to_vec : t -> Lacaml.D.vec = \n" + ] + }, + "execution_count": 30, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let to_assoc_list t = \n", + " let rec aux k accu =\n", + " if k = 0 then\n", + " accu\n", + " else\n", + " aux (k-1) ( (t.indices.{k}, t.values.{k})::accu )\n", + " in\n", + " aux (nnz t) []\n", + " \n", + " \n", + "let to_vec t =\n", + " let result = L.Vec.make0 (dim t) in\n", + " iter (fun k v -> result.{k} <- v) t;\n", + " result" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Test" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "metadata": { + "collapsed": true + }, + "outputs": [ + { + "data": { + "text/plain": [ + "- : (int * float) list = [(1, 1.); (2, -2.); (5, 0.5); (6, 1e-08); (8, 3.)]\n" + ] + }, + "execution_count": 31, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : bool = true\n" + ] + }, + "execution_count": 31, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "to_assoc_list sparse_a;;\n", + "to_vec sparse_a = dense_a;;" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Operations" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## One-vector operations" + ] + }, + { + "cell_type": "code", + "execution_count": 64, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val immutable : (t -> 'a) -> t -> t = \n" + ] + }, + "execution_count": 64, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val scale_mut : float -> t -> unit = \n" + ] + }, + "execution_count": 64, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val scale : float -> t -> t = \n" + ] + }, + "execution_count": 64, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val neg : t -> t = \n" + ] + }, + "execution_count": 64, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let immutable f t =\n", + " let result = copy t in\n", + " f result;\n", + " result\n", + " \n", + "let scale_mut x t = \n", + " L.scal x t.values \n", + "\n", + "let scale x = immutable @@ scale_mut x\n", + " \n", + "let neg t =\n", + " { t with values = L.Vec.neg t.values }\n", + " " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Test" + ] + }, + { + "cell_type": "code", + "execution_count": 67, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "val sparse_b : t =\n", + " {dim = 9; nnz = 5; indices = ;\n", + " values = R1 R2 R3 R4 R5\n", + " 1 -2 0.5 1E-08 3}\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "val sparse_c : t =\n", + " {dim = 9; nnz = 5; indices = ;\n", + " values = R1 R2 R3 R4 R5\n", + " 0.5 -1 0.25 5E-09 1.5}\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{dim: 9\n", + "(1, 1.000000) (2, -2.000000) (5, 0.500000) (6, 0.000000) (8, 3.000000) }\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{dim: 9\n", + "(1, 0.500000) (2, -1.000000) (5, 0.250000) (6, 0.000000) (8, 1.500000) }\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{dim: 9\n", + "(1, 0.500000) (2, -1.000000) (5, 0.250000) (6, 0.000000) (8, 1.500000) }\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "{dim: 9\n", + "(1, -1.000000) (2, 2.000000) (5, -0.500000) (6, -0.000000) (8, -3.000000) }\n" + ] + }, + { + "data": { + "text/plain": [ + "- : unit = ()\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + }, + { + "data": { + "text/plain": [ + "- : bool = true\n" + ] + }, + "execution_count": 67, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "let sparse_b = copy sparse_a;;\n", + "scale_mut 0.5 sparse_b;;\n", + "let sparse_c = scale 0.5 sparse_a;;\n", + "Format.printf \"%a@.\" pp sparse_a;;\n", + "Format.printf \"%a@.\" pp sparse_b;;\n", + "Format.printf \"%a@.\" pp sparse_c;;\n", + "Format.printf \"%a@.\" pp (neg sparse_a);;\n", + " \n", + "let _ =\n", + " let n1 = neg sparse_a in\n", + " neg n1 = sparse_a" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "celltoolbar": "Raw Cell Format", + "kernelspec": { + "display_name": "OCaml default", + "language": "OCaml", + "name": "ocaml-jupyter" + }, + "language_info": { + "codemirror_mode": "text/x-ocaml", + "file_extension": ".ml", + "mimetype": "text/x-ocaml", + "name": "OCaml", + "nbconverter_exporter": null, + "pygments_lexer": "OCaml", + "version": "4.07.1" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +}