owlbarn / owl Goto Github PK

Liang,

In this example http://members.logical.net/~marshall/foo05to06.pdf , the figures on the right were generated by Plot.plot. Notice that the y axis tick labels are rendered in a kind of scientific notation rather than using regular decimal numbers. That's not what I want.

I don't fully undesrstand why this happens. This has something to do with the many decimal points of precision in the y values, I think. (I'm using the default Mat module, and all of the numbers are between 0 and 1.)

You can control this behavior with plsxax, plsyax, and plszax for the three different axes. plprec is similar, but I think that applies to all three axes, and the description of its arguments is different.

What I've done in my fork marshalls_current branch is to add {x,y,z}digits fields to page, corresponding set_ functions, and then in _prepare_page I call plsxax, plsyax and plszax. This has to happen before running some of plenv, pladv, plbox3, or it won't apply to the first plot on a multi-plot page. I haven't tested this completely thoroughly, but it fixed my problem: Now I don't get scientific notation.

However, I'm not sure how to get the default treatment of tick labels now. You can make the tick labels switch to scientific notation when the data (?) exceeds the number of decimal points that you pass to pls?ax, or leave the value at the default I set, zero, in which case you'll never get scientific notation. I haven't figured out what the internal default number of decimal points is for plplot is for this operation. Maybe the behavior I've set up is OK, since a user can always override it by calling set_?digits.

If you want, feel free to request a PR any time, or to copy the code. I can probably pull in the current master and insert my code pretty quickly, although it might take a few days if I get too busy. Or maybe you'd prefer to wait until I understand this better, unless you have more insight about it.

(I know I said I wasn't going to be contributing for a couple of weeks because of my upcoming conference. For the most part I'm just adding comments about things that I have done for the sake of my conference presentation. This is one of those things.)

Here's a compare function written using fold2 (see issue #39). This can be used, for example, to allow matrices to be used with Jane Street's Interval module.

(** A compare function for matrices that returns zero if all elements of
    both matrices are equal, -1 if each element of the first is less
    than or equal to the corresponding element of the second, or 1 if
    each element of the first is greater than or equal to each element
    of the second.  Raises an exception otherwise. *)
let compare m1 m2 =
  let f acc e1 e2 =
    match acc with
         | -1 -> if e1 <= e2 then -1 else failwith "compare: incomparable"
         |  1 -> if e1 >= e2 then  1 else failwith "compare: incomparable"
         |  0 -> if e1 = e2 then 0
                 else if e1 < e2 then -1
                 else 1
         |  _ -> failwith "bug: acc is not -1, 0, or 1" (* avoid match warning *)
  in fold2 f 0 m1 m2

Treating a matrix as less than (greater than) another when every element is less than (greater than) or equal to corresponding elements is a choice. Perhaps someone might want a function that treated matrices as incomparable if some pairs were less than (greater than) while others were equal. Should there be two variants of compare?

Hi, I searched this repository for a function to the solve the eigenvalue problem, but could not find a suitable one.
Is there any way to solve it using Owl? If not, will Owl provide it in future? Or Do I miss somthing?

hi, it would be great to have the equivalent of numpy's array[3:10, 2: -4, ::2, ::-1]. afaics so far we have only array[:,3].

extending the functionality would probably require a richer variant type for the slice indices than int option. something like All | I of int | R of int * int | R' of int * int * int where R and R' are ranges with increment 1 or with given increment, respectively, and I stands for a single index. About the negative indices: i think jane street uses those as well for arrays, may be worth checking that.

it would be good to have some index that means end. but i think it's too cumbersome to replace the int's above by another variant End | Interior of int, because it requires typing Interior (or something shorter) all the time. so using -1 for the end as in numpy or potentially 0 when dealing with fortran arrays seems better. (anyway, i guess in owl we always have c-layout, no?)

hi, some comments about the infix operators:

@@ overlaps with the builtin application f @@ g x == f (g x) which i think may be a bit intrusive?

precedence:
https://caml.inria.fr/pub/docs/manual-ocaml/expr.html
i think $@ for dot product has lower precedence than summation which it probably should not have. you probably want something that starts with * or % for instance *.@. there may be other precedence issues, not sure.

why are there /$ and $/ both? actually, i suspect that all of the $... combinations are problematic wrt precedence?

finally, one could think of adding prefix operators for negative and inverse? e.g. ~-@ and ~/@.

generally, how will this work out when vectors also enter the picture? will there be a symbol do denote vector-vector operations and another one for matrix-vector? i suspect modular implicits would allow overloading those, but is it even sure that they will come?

I've been enjoying playing around with owl - thank you very much for creating. When experimenting I've noticed that matrix multiplication is quite slow for large matrices. For example if I create two random 5000x5000 matrices using Mat.uniform then multiply together it takes just under two minutes. If I do the same with Julia it takes just under two seconds (both owl and Julia on bash on Ubuntu on windows). Any obvious reason for this?

One other minor comment: it would be useful to put Mat. of_array command more prominently in the documentation as it is handy for creating small matrices.

on the slicing tutorial i saw that [3;1] will be expanded into [3;1;-1]. while i like the other slicing rules, i doubt that this particular automatic backward slicing expansion is a good idea.

reason: it's often useful in algorithms to have backwards ranges yield an empty result. for example, if the first slice bound is a running variable [i; fixed] then as soon as i exceeds fixed, we get an empty return which is often what is needed. in contrast, iterating backwards over an array is more rare in my experience, and if it is desired, it does not hurt to be forced to explicitly give the -1 as a step in these cases. (fwiw, numpy gives empty as well)

[edit]: after trying the slicing functions out, i still think it would be feasible to return a matrix with shape 0 in all dimensions where the given range specification goes backwards without an explicit negative increment. admittedly, such matrices are normally not useful on their own, but could be useful corner cases in algorithms. e.g. iterating over the 0 dimension does nothing, folding gives the initial value, etc.

in light of the recent slicing additions here is another feature from numpy that could be potentially useful: boolean indexing.

e.g. a[a>5] where a is a numpy array, and a>5 is a boolean array that gives elementwise truth values of the comparison and has the same shape as a.
then the indexing returns a flat, one-dimensional array of only the entries that are True. in numpy this is used a lot like this: a[a>5] = 5, i.e. for assignment. but in image processing, bitmasks such as a>5 itself are also useful.

i'm not sure it's necessarily a good idea to try and copy the assignment idiom in owl. for one, i don't think there is an appropriate boolean bigarray version. then, to use it for assignment, the array would need to be a view. it might be better instead to construct some kind of an iterator function that does this conveniently. currently, what i would try would be Mat.iteri (fun i j v -> if v > 5. then Mat.{i,j} <- 5.) a which seems rather long?

Executing the following code a couple of times in utop may cause segmentation fault. I have not figured out the reason for this.

let x = Dense.linspace (-2.5) 2.5 100 in
let y = Dense.linspace (-2.5) 2.5 50 in
let u, v = Dense.meshup x y in
let z = Dense.(Maths.sin @@ ((u *@ 2.) +@ (v *@ 2.))) in
let h = Plot.create "" in
let _ = Plot.contour ~h u v z in
Plot.output h;;

Passing a non-square matrix to inv aborts utop in addition to the assert failure.

Welcome to utop version 2.0.1 (using OCaml version 4.04.2)
...
# Owl.Mat.inv (Owl.Mat.sequential 2 3);;
Assertion failed: (rows() == cols()), function inverse, file ./Eigen/src/LU/InverseImpl.h, line 338.
Abort trap: 6
~$

I used the latest version of the master branch of Owl with the version of the OCaml Eigen lib currently available in opam for OCaml 4.04.2:

~$ opam list | grep -i eigen
eigen                           0.0.4  Owl's OCaml interface to Eigen3 C++ libra

Let me know if you need other info about my configuration. Thanks-

This might end up being a question about opam that I should ask about on discuss.ocaml.org or StackOverflow, but it could be that the difficulty I'm experiencing has to do with how the current versions of Owl and Eigen work together, so I thought I'd ask here first.

The latest version of Owl won't compile with the opam version of Eigen, producing this error:

# Error: Unbound value Eigen.Dense.S.gemm

I saw in the Eigen git history that dot was recently renamed to gemm, so I assumed that I just needed a more recent version of Eigen. I got Eigen from github and was able to compile and install it using 'opam pin add eigen <my source dir>. I don't quite understand all of the command variants in the Eigen source, but I think the appropriate version of gemm is in there (in eigen_dismat_s.ml?).

However, when try to pin my local owl git tree, I get the same error, ie.. that gemm is unbound (see below). I don't understand this.

On the other hand, if I then compile a file to a native executable (rather than running utop), Owl seems to be available. I have no problem compiling code that uses Owl.

Any thoughts about what might be going wrong?

Thanks-
Marshall

/docs/src/ml/owl(master)*$ opam pin add owl `pwd`
[NOTE] Package owl is already path-pinned to /Users/marshall/docs/src/ml/owl.
       This will erase any previous custom definition.
Proceed ? [Y/n] y

[owl] /Users/marshall/docs/src/ml/owl/ synchronized
[owl] Installing new package description from /Users/marshall/docs/src/ml/owl

owl needs to be installed.
The following actions will be performed:
  ∗  install owl 0.2.4*
Do you want to continue ? [Y/n] y

=-=- Gathering sources =-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
[owl.0.2.4] /Users/marshall/docs/src/ml/owl/ already up-to-date

=-=- Processing actions -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
[ERROR] The compilation of owl failed at "ocaml setup.ml -build".
Processing  1/1: [owl: ocaml]
#=== ERROR while installing owl.0.2.4 =========================================#
# opam-version 1.2.2
# os           darwin
# command      ocaml setup.ml -build
# path         /Users/marshall/.opam/yo/build/owl.0.2.4
# compiler     4.04.1
# exit-code    1
# env-file     /Users/marshall/.opam/yo/build/owl.0.2.4/owl-64484-d2d111.env
# stdout-file  /Users/marshall/.opam/yo/build/owl.0.2.4/owl-64484-d2d111.out
# stderr-file  /Users/marshall/.opam/yo/build/owl.0.2.4/owl-64484-d2d111.err
### stdout ###
# [...]
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamldep -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -modules lib/owl_stats.ml > lib/owl_stats.ml.depends
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamldep -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -modules lib/owl_dense_common.ml > lib/owl_dense_common.ml.depends
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamldep -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -modules lib/owl_pretty.mli > lib/owl_pretty.mli.depends
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamlc -c -g -annot -bin-annot -ccopt -I/Users/marshall/.opam/4.04.1/lib/ctypes -ccopt -O3 -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -I lib -I lib/backend -I lib/ext -I lib/nlp -o lib/owl_pretty.cmi lib/owl_pretty.mli
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamldep -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -modules lib/owl_dense_ndarray_a.ml > lib/owl_dense_ndarray_a.ml.depends
# /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamlc -c -g -annot -bin-annot -ccopt -I/Users/marshall/.opam/4.04.1/lib/ctypes -ccopt -O3 -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -I lib -I lib/backend -I lib/ext -I lib/nlp -o lib/owl_dense_common.cmo lib/owl_dense_common.ml
# + /Users/marshall/.opam/4.04.1/bin/ocamlfind ocamlc -c -g -annot -bin-annot -ccopt -I/Users/marshall/.opam/4.04.1/lib/ctypes -ccopt -O3 -package str -package plplot -package gsl -package eigen -package dolog -package ctypes.stubs -package ctypes -I lib -I lib/backend -I lib/ext -I lib/nlp -o lib/owl_dense_common.cmo lib/owl_dense_common.ml
# File "lib/owl_dense_common.ml", line 202, characters 17-35:
# Error: Unbound value Eigen.Dense.S.gemm
# Command exited with code 2.
### stderr ###
# E: Failure("Command ''/Users/marshall/.opam/4.04.1/bin/ocamlbuild' lib/libowl_stubs.a lib/dllowl_stubs.so lib/owl.cma lib/owl.cmxa lib/owl.a lib/owl.cmxs lib/neural/owl_neural.cma lib/neural/owl_neural.cmxa lib/neural/owl_neural.a lib/neural/owl_neural.cmxs perftest/perf_common.cma perftest/perf_common.cmxa perftest/perf_common.a perftest/perf_common.cmxs perftest/perf_dense_real.native perftest/perf_sparse_real.native perftest/perf_sparse_complex.native perftest/perf_dense_ndarray.native perftest/perf_dense_matrix.native perftest/perf_sparse_ndarray.native perftest/perf_sparse_matrix.native perftest/perf_specific_fun.native examples/test_neural.native -use-ocamlfind -tag debug' terminated with error code 10")



=-=- Error report -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
The following actions failed
  ∗  install owl 0.2.4
No changes have been performed
[NOTE] Pinning command successful, but your installed packages may be out of sync.

This might be planned anyway--I know that work on the plot routines are in progress--but it would be helpful to be able to specify the orientation of 3D plots. Thanks for including the plot routines. I just took a look at plplot, and I appreciate the way in which Owl greatly simplifies usage, as Tutorial 4 notes.

I couldn't find sorting functions for matrices or vectors. it would be good to have

• sorting for vectors with customizable compare function, with a standard default compare

• sorting for matrices by row or by column (i guess here a compare function acting on the row/column has to be given always as there is no standard way),

• potentially also for Nd arrays where the sort dimension can be given and the compare function would take a slice over all non-sorted dimensions. this seems a little less standard

it might be that in-place versions are also needed? ideally both...

some questions that did not get mentioned in the tutorial:

[edited after actually trying it out... sorry for unnecessary questions before]

a) it might be a good feature to have views as well, which always share data?

b) what about increasing/decreasing the rank of an array? something like a[:,2,:] or a[:, np.newaxis] in numpy?
a non-automatic version of rank reduction would be a separate squeeze function which takes all slice dimensions of length 1 and eliminates that dimension from the array, keeping all the same data. i.e. if the shape is [|3;1;4|] it will be [|3;4|] afterwards. for typing reasons this function could only be Ndarray -> Ndarray. analogous functions could exists from matrices to vectors.
alternatively, the Ndarray version of slice could get an optional boolean ?squeeze argument.
for rank increase one could have something like increase_rank 1 a which makes a shape [|3;1;4|] out of [|3;4|]. again, only for Ndarray. for vectors to matrices the functionality is probably already there.

i believe increasing/decreasing the rank is important especially in the absence of general broadcasting, so that arrays can be brought into the same shape easily before combining them elementwise.

c) how about indexing a general list of elements, e.g. indices 1,2,4,4,3 along a dimension of an array?

this seems to require a different syntax; the list syntax is already taken for slicing.
would it be feasible/desirable to accommodate this list indexing in some way? when trying to build this into the slice function, one could use an array [|1;2;4;4;3|]. but i guess it's difficult to allow for either a slice specification as an int list or a direct index set as an int array as an argument to the slice function, without having to put them into a variant type like Slice_spec of int list | Index_set of int array. so that solution would require very verbose syntax i fear. (for comparison, numpy has the colon syntax for slice specs and also allows a list syntax for explicit indices).

another option would be to have a separate function for index lists, in addition to slice for index ranges. this would probably mostly work, but it would be impossible to mix slicing and set indexing along different dimensions of the same array. (what a[:,[1,2,1]] in numpy does). maybe it's acceptable to have to do it separately? although it would probably allocate intermediate arrays.

potentially the following would be a concise enough syntax: only one slice function, but a constructor S with two arguments for regular slices, one Ss with three arguments for strided slices and a constructor Il with a list argument for index lists? something like general_slice [S (1, 3); Ss (0, 4, 2); Il [1;2;4;4;3]] x is relatively cumbersome to type compared to numpy; but maybe it would be tolerable?
one could reserve short function names which just apply the constructors, to be able to write general_slice [s 1 3; s 0 4 ~by:2; il [1;2;4;4;3]]; this may be acceptably concise? but it is a bit intrusive, claiming the symbol s for example seems excessive.

just a collection of thoughts... i hope they're helpful

Cool project! It's always nice to see more people doing good work in numericcal computing.

I noticed that the Julia benchmark code here uses the non-constant global variable shape inside of functions, which is a performance no-no. The simplest fix is to just change that line to const shape = (10,1000,10000) since it never changes.

It is also a bit unclear from the benchmark code whether the standard recommendation not timing the first invocation of functions to avoid timing JIT is followed – only the function definitions are included in the gist, not their invocation. For real-world code, this doesn't matter since for long-running or repeated computations (i.e. the ones whose speed one cares about), JIT is a negligible part of the time, but benchmarks tend to be short, so JIT ends up looking substantial.

Is there a reason Owl doesn't use DevFile plugin in _oasis? This makes the usual clone&make workflow a bit painful, because all target does not generate setup.data but rather expects it to be already generated. Therefore one has to a) read the Makefile and b) do make oasis all instead.

This is just an idea, but I think it would be nice to provide a vector/matrix interface that abstracts away the concrete encoding, spare or dense. The abstraction layer would automatically convert between encodings, transparently to the user. For example,

• multiplying two sparse matrices would yield a sparse matrix
• multiplying two dense matrices would yield a dense matrix
• multiplying a dense and a sparse matrix would either
  • if the dense matrix is spare, convert it to a sparse matrix and return a sparse matrix
  • otherwise, convert the sparse matrix to dense one and return a dense matrix
    etc.

Since sparsity checks may be expensive, an important consideration would be when and where to perform them. There could be both an automatic mode, where sparsity checks & conversions are performed fully automatically, and a manual mode, where the user could say explicitly through the API when to perform sparsity checks.

I think such an abstraction layer could be super useful for prototyping.

I found that Math.combination gives odd results when the first argument is large, but realized that this is probably expected from the conversion to OCaml ints. My problem was solved by seeing that combination wraps Gsl.Sf.choose and using the latter instead. The result will get multiplied by a float, anyway, so using the GSL function is better. (It's good that I experienced the problem, or I would have continued to convert the int output of combination back to a float!)

This made me wonder whether it would be worth including float versions of combination in Owl.Maths. I don't know the overall goals of the module, and I now know how to call GSL directly, so I don't have to have a float version of combination in Maths, but I thought I'd raise the question.

(My use case: In population genetics, a Wright-Fisher model of the evolution of a population uses binomially distributed transition probabilities from one state of the population, in which there are i copies of a gene, to a state with j copies. This means that a transition matrix for a population of N organisms with up to 2N copies of a gene needs to compute combinations in which 2N+1 is the upper term. That produces some large numbers for even a small population of reasonable size.)

(I hope I'm not being annoying with too many minor issues and questions.)

@jzstark In the example below, the axes become dotted lines. I guess you forgot to restore the pen and line style after drawing the line extension.

let y = Mat.(gaussian 100 1 *$ 10.) in
let x = Mat.gaussian 200 1 in
let h = Plot.create ~m:2 ~n:2 "plot_025.png" in

Plot.subplot h 0 0;
Plot.set_title h "Gaussian vs. Gaussian Sample";
Plot.set_ylabel h "Quantiles of Input Sample";
Plot.set_xlabel h "Normal Distribution Quantiles";
Plot.qqplot ~h y ~x:x;

Plot.subplot h 0 1;
Plot.set_title h "Gaussian vs. Default Dist";
Plot.set_ylabel h "Quantiles of Input Sample";
Plot.set_xlabel h "Normal Distribution Quantiles";
Plot.qqplot ~h y ~spec:[RGB (0,128,255)];

Plot.subplot h 1 0;
Plot.set_title h "Gaussian vs. Rayleigh Dist";
Plot.set_ylabel h "Quantiles of Input Sample";
Plot.set_xlabel h "Rayleigh Distribution (sigma=0.5) Quantiles";
Plot.qqplot ~h y ~pd:(fun p -> Stats.Cdf.rayleigh_Pinv p 0.5);

Plot.subplot h 1 1;
Plot.set_title h "Gaussian vs. Chi-Square Dist";
Plot.set_ylabel h "Quantiles of Input Sample";
Plot.set_xlabel h "Chi-Square Distribution (k=10) Quantiles";
Plot.qqplot ~h y ~pd:(fun p ->  Stats.Cdf.chisq_Pinv p 10.);

Plot.output h;;

This has the effect that I cannot use it from my utop or ocaml top level:

utop # #require "Owl";;
The files /home/flux/ocaml/.opam/4.03.0+fp+flambda/lib/Owl/Owl.cma and
/home/flux/ocaml/.opam/4.03.0+fp+flambda/lib/ocaml/compiler-libs/toploop.cmi
disagree over interface Types

This is of course unfortunate that the collision can occur, but as long as OCaml doesn't have namespaces ;) it would be nice if Owl had all its modules prefixed with Owl_, possibly with an additional module Owl that has module aliases for the current names (module Owl = struct module Types = Owl_types etc). This would avoid introducing the module names that might conflict with other projects.

I can foresee that "utils" is also going to collide some projects..

I'm suggesting this instead of putting all inside the module Owl to avoid the binary size bloat that can occur with packing.

I have problems using owl from the toplevel (also with utop). Both when compiling owl from source, and when installing via opam, I get the following error message:

# #require "owl";;
/home/pkl/.opam/4.04.0/lib/ocaml/str.cma: loaded
/home/pkl/.opam/4.04.0/lib/ocaml/unix.cma: loaded
/home/pkl/.opam/4.04.0/lib/ocaml/bigarray.cma: loaded
/home/pkl/.opam/4.04.0/lib/bytes: added to search path
/home/pkl/.opam/4.04.0/lib/integers: added to search path
/home/pkl/.opam/4.04.0/lib/integers/integers.cma: loaded
/home/pkl/.opam/4.04.0/lib/integers/integer_printers.cma: loaded
/home/pkl/.opam/4.04.0/lib/ctypes: added to search path
/home/pkl/.opam/4.04.0/lib/ctypes/ctypes.cma: loaded
/home/pkl/.opam/4.04.0/lib/ctypes/ctypes-top.cma: loaded
/home/pkl/.opam/4.04.0/lib/ctypes/cstubs.cma: loaded
/home/pkl/.opam/4.04.0/lib/plplot: added to search path
/home/pkl/.opam/4.04.0/lib/plplot/plplot.cma: loaded
/home/pkl/.opam/4.04.0/lib/dolog: added to search path
/home/pkl/.opam/4.04.0/lib/dolog/dolog.cma: loaded
/home/pkl/.opam/4.04.0/lib/eigen: added to search path
/home/pkl/.opam/4.04.0/lib/eigen/eigen.cma: loaded
/home/pkl/.opam/4.04.0/lib/gsl: added to search path
/home/pkl/.opam/4.04.0/lib/gsl/gsl.cma: loaded
/home/pkl/.opam/4.04.0/lib/owl: added to search path
/home/pkl/.opam/4.04.0/lib/owl/owl.cma: loaded
Cannot load required shared library dllowl_stubs.
Reason: /home/pkl/.opam/4.04.0/lib/stublibs/dllowl_stubs.so: /home/pkl/.opam/4.04.0/lib/stublibs/dllowl_stubs.so: undefined symbol: _ZGVbN4v___expf_finite.

It seems owl only has examples at the moment. These are no doubt nice to have, but running and inspecting them manually after each change could be less effective than automated testing.

Have you considered adding a layer of unit/property tests?

I'm getting an error trying to create a 3D plot with mesh or surf:

plot3dcl: Y array must be strictly increasing, aborting operation

I'm able to run the first 3D example from the plotting tutorial, and when I compare the direction in which values increase in the y matrix there, it's the same direction in which values increase in my y (ys) matrix.

In the example below, the matrices are small for the sake of display, but I've tried larger matrices, I've tried dividing the values to create smaller increments, I've tried starting the x and y values at a negative number, and I've even generated the same x and y matrices using meshgrid (I'm not sure why that would make a difference, but I did it).

I can't figure out what I'm doing wrong. Maybe I'm just being blind, or maybe there is something odd going on in the plplot function, but I wondered whether the transformations of the matrices that occur within mesh and surf might be creating a problem for me. I also wonder whether the plplot error message is misleading. Maybe the problem isn't with the y matrix.

Oh--this just occurred to me. Is it possible that the z values aren't close enough together, and that's the problem? contour and heatmap work with the data below, btw.

Thanks again.

# xs, ys, zs;;

   C0 C1 C2 C3 C4 C5 C6 C7 C8
R0  0  1  2  3  4  5  6  7  8
R1  0  1  2  3  4  5  6  7  8
R2  0  1  2  3  4  5  6  7  8
R3  0  1  2  3  4  5  6  7  8
R4  0  1  2  3  4  5  6  7  8
R5  0  1  2  3  4  5  6  7  8
R6  0  1  2  3  4  5  6  7  8
R7  0  1  2  3  4  5  6  7  8


   C0 C1 C2 C3 C4 C5 C6 C7 C8
R0  0  0  0  0  0  0  0  0  0
R1  1  1  1  1  1  1  1  1  1
R2  2  2  2  2  2  2  2  2  2
R3  3  3  3  3  3  3  3  3  3
R4  4  4  4  4  4  4  4  4  4
R5  5  5  5  5  5  5  5  5  5
R6  6  6  6  6  6  6  6  6  6
R7  7  7  7  7  7  7  7  7  7


           C0        C1        C2        C3        C4        C5        C6        C7         C8
R0          0         0         0         0         1         0         0         0          0
R1 0.00390625   0.03125  0.109375   0.21875  0.273438   0.21875  0.109375   0.03125 0.00390625
R2  0.0318424 0.0756443  0.128525  0.170623   0.18673  0.170623  0.128525 0.0756443  0.0318424
R3  0.0754722 0.0898582  0.121041  0.140182  0.146893  0.140182  0.121041 0.0898582  0.0754722
R4    0.12236 0.0886335  0.109046  0.118895   0.12213  0.118895  0.109046 0.0886335    0.12236
R5   0.167026 0.0816821 0.0966715  0.102487  0.104266  0.102487 0.0966715 0.0816821   0.167026
R6   0.207606 0.0731666 0.0851005 0.0890421 0.0901699 0.0890421 0.0851005 0.0731666   0.207606
R7   0.243728 0.0647161 0.0746732 0.0776517 0.0784627 0.0776517 0.0746732 0.0647161   0.243728

val xs : Mat.mat =
val ys : Mat.mat =
val zs : Mat.mat =

# let h = Plot.create "yo.pdf";;
val h : Owl_plot.handle = <abstr>

# Plot.mesh ~h xs ys zs;;
- : unit = ()

# Plot.output h;;

*** PLPLOT ERROR, ABORTING OPERATION ***
plot3dcl: Y array must be strictly increasing, aborting operation
- : unit = ()

Have you considered swapping GSL for Cephes? It lacks some functionality compared to GSL, but is available under a much less restrictive license, see SciPy repo.

set_xticklabels, set_yticklabels, and set_zticklabels don't seem to work in 3D plots. The first two do seem to work in 2D plots:

let f x = Maths.sin x /. x in
let h = Plot.create "foo2d.pdf" in
Plot.set_xticklabels h [(8., "eight")];
Plot.plot_fun ~h f 1. 15.;
Plot.output h

but not in 3D plots:

let x, y = Mat.meshgrid (-2.5) 2.5 (-2.5) 2.5 100 100 in
let z0 = Mat.(sin ((pow1 x 2.) + (pow1 y 2.))) in
let h = Plot.create "foo3d.pdf" in
Plot.set_yticklabels h [(0., "zero")];
Plot.surf ~h x y z0;
Plot.output h

The y border is front right, and I would expect the "0" in the middle of that edge to be replaced by "zero", but the "0" remains (and no other axis is affected).

Thanks again-

Marshall

I wrote a compare function in terms of fold2, which seems generally useful. It would probably belong in owl_dense_matrix_generic.ml.

let fold2 f init m1 m2 =
  let rows, cols as dims = shape m1 in
  if dims <> shape m2 then failwith "fold2: matrices have different shapes"
  ;
  let last_col = cols - 1 in
  let apply_f acc i j = 
    f acc (get m1 i j) (get m2 i j)
  in
  let rec loop acc i j =
    if i < rows
    then loop (apply_f acc i j) (i + 1) j
    else if j < last_col         (* don't start on next col if at final col *)
         then loop acc 0 (j + 1) (* start over on next col *)
         else acc
  in loop init 0 0

It's not pretty, and the style probably shows my inexperience, but it's tail recursive. I didn't succeed in coming up with a more elegant purely functional definition. However, I've noticed that existing fold functions in owl use imperative methods internally. I assume that's more efficient, so perhaps you wouldn't want to a version of fold2 similar to this one. I'm not sure what the tradeoffs are between different iter function variants (iteri is defined in owl_dense_matrix_generic.ml, but iter there is defined by the iter in owl_dense_ndarray_generic.ml, which I haven't sorted out fully).

Could you list the tests which are yet to be implemented?

I have some tests implemented in pareto, maybe some of that could be useful for owl as well.

Hello,
I am using OCaml 4.02.1 with Opam, and I get an error when running make. Here is the output:

:-) -- scemama@lpqdh82 -- /home/scemama/owl --  master  --  -- 14:51:45 --
[0] $ make oasis
oasis setup
ocaml setup.ml -configure

Configuration: 
ocamlfind: ........................................... /home/scemama/.opam/4.02.1/bin/ocamlfind
ocamlc: .............................................. /home/scemama/.opam/4.02.1/bin/ocamlc.opt
ocamlopt: ............................................ /home/scemama/.opam/4.02.1/bin/ocamlopt.opt
ocamlbuild: .......................................... /home/scemama/.opam/4.02.1/bin/ocamlbuild
Package name: ........................................ Owl
Package version: ..................................... 0.1
os_type: ............................................. Unix
system: .............................................. linux
architecture: ........................................ amd64
ccomp_type: .......................................... cc
ocaml_version: ....................................... 4.02.1
standard_library_default: ............................ /home/scemama/.opam/4.02.1/lib/ocaml
standard_library: .................................... /home/scemama/.opam/4.02.1/lib/ocaml
standard_runtime: .................................... /home/scemama/.opam/4.02.1/bin/ocamlrun
bytecomp_c_compiler: ................................. gcc -fno-defer-pop -Wall -D_FILE_OFFSET_BITS=64 -D_REENTRANT -fPIC
native_c_compiler: ................................... gcc -Wall -D_FILE_OFFSET_BITS=64 -D_REENTRANT
model: ............................................... default
ext_obj: ............................................. .o
ext_asm: ............................................. .s
ext_lib: ............................................. .a
ext_dll: ............................................. .so
default_executable_name: ............................. a.out
systhread_supported: ................................. true
Install architecture-independent files dir: .......... /usr/local
Install architecture-dependent files in dir: ......... $prefix
User executables: .................................... $exec_prefix/bin
System admin executables: ............................ $exec_prefix/sbin
Program executables: ................................. $exec_prefix/libexec
Read-only single-machine data: ....................... $prefix/etc
Modifiable architecture-independent data: ............ $prefix/com
Modifiable single-machine data: ...................... $prefix/var
Object code libraries: ............................... $exec_prefix/lib
Read-only arch-independent data root: ................ $prefix/share
Read-only architecture-independent data: ............. $datarootdir
Info documentation: .................................. $datarootdir/info
Locale-dependent data: ............................... $datarootdir/locale
Man documentation: ................................... $datarootdir/man
Documentation root: .................................. $datarootdir/doc/$pkg_name
HTML documentation: .................................. $docdir
DVI documentation: ................................... $docdir
PDF documentation: ................................... $docdir
PS documentation: .................................... $docdir
findlib_version: ..................................... 1.6.2
is_native: ........................................... true
suffix_program: ...................................... 
Remove a file.: ...................................... rm -f
Remove a directory.: ................................. rm -rf
Turn ocaml debug flag on: ............................ true
Turn ocaml profile flag on: .......................... false
Compiler support generation of .cmxs.: ............... true
OCamlbuild additional flags: ......................... 
Create documentations: ............................... true
Compile tests executable and library and run them: ... false
ocamldoc: ............................................ /home/scemama/.opam/4.02.1/bin/ocamldoc
pkg_str: ............................................. /home/scemama/.opam/4.02.1/lib/ocaml
pkg_gsl: ............................................. /home/scemama/.opam/4.02.1/lib/gsl
pkg_ctypes: .......................................... /home/scemama/.opam/4.02.1/lib/ctypes
pkg_ctypes_foreign: .................................. /home/scemama/.opam/4.02.1/lib/ctypes
pkg_plplot: .......................................... /home/scemama/.opam/4.02.1/lib/plplot
pkg_dolog: ........................................... /home/scemama/.opam/4.02.1/lib/dolog


:-) -- scemama@lpqdh82 -- /home/scemama/owl --  master  --  -- 14:51:49 --
[0] $ make
ocaml setup.ml -build
Finished, 1 target (0 cached) in 00:00:00.
+ /home/scemama/.opam/4.02.1/bin/ocamlfind ocamldep -package ctypes -package ctypes.foreign -package dolog -package gsl -package plplot -package str -modules lib/types.ml > lib/types.ml.depends
ocamlfind: Package `ctypes.foreign.unthreaded' not found - required by `ctypes.foreign'
Command exited with code 2.
Compilation unsuccessful after building 1 target (0 cached) in 00:00:00.
E: Failure("Command ''/home/scemama/.opam/4.02.1/bin/ocamlbuild' ./lib/Owl.cma ./lib/Owl.cmxa ./lib/Owl.a ./lib/Owl.cmxs ./examples/test_sgd.byte ./examples/test_kmeans.byte ./examples/test_log.byte ./examples/test_svm.byte ./examples/test_operation.byte ./examples/test_dense.byte ./examples/test_sparse.byte -tag debug' terminated with error code 10")
Makefile:2: recipe for target 'all' failed
make: *** [all] Error 1

:-( -- scemama@lpqdh82 -- /home/scemama/owl --  master  --  -- 14:52:02 --
[0] $ opam install ctypes.foreign
[ERROR] No package matches ctypes.foreign.

what should I do?

It looks like building the current version of Owl requires that openblas be installed, specifically under /usr/local/opt/openblas:

$ make
ocaml setup.ml -build
Start compiling Owl ...
Finished, 1 target (0 cached) in 00:00:01.
+ /Users/marshall/.opam/4.04.2/bin/ocamlfind ocamlmklib -o lib/owl_stubs -L/usr/local/opt/openblas/lib -L/usr/local/lib/gcc/7 -framework Accelerate -lopenblas -lgfortran lib/owl_dense_common_c.o lib/owl_dense_common_vec_cmp.o lib/owl_dense_common_vec_map.o lib/owl_dense_common_vec_fold.o lib/owl_dense_common_vec_combine.o lib/owl_cblas_generated_stub.o lib/owl_lapacke_generated_stub.o
ld: warning: directory not found for option '-L/usr/local/opt/openblas/lib'
ld: library not found for -lopenblas
clang: error: linker command failed with exit code 1 (use -v to see invocation)
Command exited with code 2.
Compilation unsuccessful after building 15 targets (0 cached) in 00:00:25.
E: Failure("Command ''/Users/marshall/.opam/4.04.2/bin/ocamlbuild' lib/libowl_stubs.a lib/dllowl_stubs.so lib/owl.cma lib/owl.cmxa lib/owl.a lib/owl.cmxs perftest/perf_common.cma perftest/perf_common.cmxa perftest/perf_common.a perftest/perf_common.cmxs perftest/perf_dense_real.native perftest/perf_sparse_real.native perftest/perf_sparse_complex.native perftest/perf_dense_ndarray.native perftest/perf_dense_matrix.native perftest/perf_sparse_ndarray.native perftest/perf_sparse_matrix.native perftest/perf_specific_fun.native examples/test_neural.native -use-ocamlfind -tag debug' terminated with error code 10")
make: *** [all] Error 1

Building and installing openblas got me past this error (there doesn't seem to be a pre-built version for MacOS), but it looks like you have to install openblas using make PREFIX=/usr/local/opt/openblas install in order for the Owl build scripts to find it. (I thought PREFIX=/usr/local/opt would put openblas in the right place, but that's incorrect.)

If this is all as expected, then it would be useful to have this information in the README.

Thanks again-

I'm an OCaml newbie; Owl makes it a lot more useful for me. Thanks very much.

I like that matrices are nicely formatted by default in repls. That's very nice. (By contrast, in Clojure's core.matrix, the way that matrices are displayed depends on which matrix implementation you're using, and some are not displayed in an optimal manner. You have to pass matrices to a pretty-printing function to assure a nice display.)

However, when I put matrices in a list, what's displayed is a bit strange. Here's a list of five 2x2 matrices In a standard ocaml 4.04.0 toplevel in a terminal window on MacOS:

- : Owl.Mat.mat list =
[   C0 C1
 R0  2  4
 R1  1  0
;    C0 C1
            R0  0  6
            R1  0  3
;
    C0 C1
 R0  0  2
 R1  3  4
;    C0 C1
            R0  6  4
            R1  4  2
;
    C0 C1
 R0  0  3
 R1  3  6
]

Here is a similar list of five matrices in utop:

   C0 C1
R0  2  4
R1  1  4
   C0 C1
         R0  5  0
         R1  0  0
   C0 C1
                  R0  0  1
                  R1  1  4
   C0 C1
                           R0  0  1
                           R1  3  3
   C0 C1
                                    R0  0  3
                                    R1  3  5
- : mat list = [; ; ; ; ]

I guess there's some kind of interaction between the way that Owl is formatting matrices and the way the repls are formatting lists.

This is obviously not an urgent issue.

Thanks for the new "zoo" module. I was playing with it and following the tutorial. When I run the command owl -upload myscript, I got the error msg: sh: 1: owl_upload_gist.sh: not found. In the source code, the only mention of owl_upload_gist.sh is in upload_gist from owl_zoo.ml, but it doesn't specify what is this shellscript or where to find it.

I also run owl -list and it turned out that owl_list_gist.sh is also not found. So I suspect that it might be because of some configuration problem of system path on my machine. I can confirm that these owl_*_gist.sh scripts are not in /usr/local/bin where the owl command locates.

Thanks for any help!

There are variations on 3D plots that Owl doesn't currently provide but that would be useful for me. For example, this plot was created by using [PL_DRAW_LINEY] as the opt argument to plmesh:
foo.pdf
(The color variation in the plot is something I'll discuss in a different issue.)

I don't want to make Owl's plotting API overly complicated. Its simplicity is what makes it so nice to use. So I've thought about different ways to make mesh and surf more flexible. This is my current thinking:

Ideally it's good to give users as much flexibility as possible with mesh and surf.

(a) However, adding an optional argument for every option that can be passed (e.g. PL_DRAW_LINEY, etc. is confusing and would make the interface to mesh and surf a mess.

(b) The number of optional arguments could be reduced by adding only one optional argument that allows users to pass something like Plplot.([ PL_DRAW_LINEXY; PL_MAG_COLOR; PL_MESH ]). This is confusing, though, and interferes with the simplicity of the Owl interface.

(c) Another option is to define new 3D plot options that use different plplot option lists. However, it's also important to avoid writing many different plot functions with a lot of duplicate code. (On the other hand, I do think that a little bit of duplicated code can be easier to understand rather than abstracting out too much of the code into common functions.)

So there's no ideal solution, other just than (d) not allowing more flexibility in 3D plotting.

My inclination at this point is use option (b) to add an optional first argument to mesh and surf. This argument allows passing a list of plplot options for plmesh, plmeshc, or plsurf3d. It's not a nice interface, but then one can use this optional argument to define new 3D plotting functions as in (c).

For example, this is the definition of the function plots2d_in_3d that I used to generate the PDF linked above:

let plots2d_in_3d = mesh ~plplot_opt:Plplot.([PL_DRAW_LINEY])

(I couldn't think of a better name for this function.) The full revised code for mesh that I used (and surf, too) is here: https://github.com/mars0i/owl/blob/master/lib/owl_plot.ml#L877

What do you think? I didn't want to submit a PR before seeing what your thoughts about this were.

I noticed minor typos on this page, https://github.com/ryanrhymes/owl/wiki/Tutorial:-Matrix-Manipulation. It includes the following list:

Mat.(x +@ y);;    (* add two matrices *)
Mat.(x -@ y);;    (* subtract y from x *)
Mat.(x *@ y);;    (* element-wise multiplication *)
Mat.(x /@ y);;    (* element-wise division *)
Mat.(x $@ y);;    (* dot product of x and y *)

Element-wise multiplication should be * and the dot product should be *@. This might just be left over from an earlier version. (I'm not supposed to be able to propose edits to the wiki pages, right? I didn't see an edit icon.)

i saw there is a vector module now - i think this makes perfect sense. in particular, linspace should give a vector IMO, and functions like numpy.arange and numpy.logspace would be good additions as well.

Einstein summation generalizes most tensor operations (product, matrix multiplication, etc).

Is it possible to set the ticks of x axis ? I'd like to print Jan, Feb, Mar, ... and so on instead of the automatic 0 2 4 ....

Right now, plot titles specified by set_title don't display on 3D plots. This is because this line in owl_plot.ml:

 let _ = pllab p.xlabel p.ylabel p.title in ()

only appears under if not p.is_3d then, and not also within the then clause.

The problem is easy to fix by adding a similar line in the then clause.

It can't be moved out of the if/then block, though, because pllab can't run until the page is prepared by plenv or pladv, which are called in the if/then block.

(One thing that might seem confusing is that pllab sets the x and y labels, but in the 3D code section, there's a call to plbox3; this is needed to set the z label, and it also sets the x and y labels. However, in my experiments, there's no interference between pllab and plbox3. In both cases, the code sets the x and y labels from the xlabel and ylabel fields in the page structure. The x and y labels gets re-set, apparently, by whichever runs last, pllab or plbox3, but since the same string is used either way, the end result is the same afaics.)

I'll be happy to submit a PR when it's appropriate, but I think it's better to wait until the current refactoring is done. There's no rush.

For me it would be useful to have an optional ~color argument for mesh and surf, like the optional argument for the 2D plotting functions, so that axes and grids can be a different color from plot lines. Maybe this makes more sense if new features are implemented for issue #57, as illustrated in this example: foo.pdf. I can submit a PR if this seems OK.

Wow, that looks awesome.
This library makes me want to try ocaml instead of my current julia workflow because I miss ML/OCaml/haskell.
I have a feature request :) It looks like ploting part cannot generate latex style labels and titles. It would be nice to have.

and I have a question. would owl arrays work with Sundials/ML http://inria-parkas.github.io/sundialsml/ ? if not, maybe there it would be nice to have a wraper of ODE solver in owl from gsl or sundials :)

Could you sketch some of the future directions for owl?

I think OCaml community might benefit from a single numerical library and would like to contribute.

As I reported on another issue, I was able to build Owl from source successfully after I installed openBLAS on one of my machines ("2011"). Now I’m having trouble on the other machine ("2015"). Both run MacOS 10.11.6. I can’t build Owl from source on the 2015 machine—that was the first problem after I installed gcc 7 and openBLAS. (There's been fuss lately about a bug with Sandy Bridge vs. Ivy Bridge CPUs. I believe I have Ivy Bridge.)

Then I found out that I can’t install Owl via opam either, so I don’t think the problem is with compiling from raw source or with openBLAS. At this point, I don’t have owl installed on the 2015 machine because I’d gotten it uninstalled when trying to build Owl from source.

All or nearly all of the error messages mention eigen-related files. That's true when I try to build Owl from source, too, but there are more error messages.

I’m only including the error report from opam, since the problem doesn’t seem to be specific to building from raw source. Let me know if you want more, though.

No need to reply quickly. I can keep working on the 2011 machine. (It's slower, but in the end it's the one that really matters at present.)

If you get a moment, I'm wonding if any diagnosis or suggestions jump out at you from the error messages below. I’ve tried changing to a new switch, recompiling openBLAS, opam uninstalling eigen, ctypes, etc. and reinstalling. Hopefully I didn't do something silly and obvious that I'm overlooking.

Thanks!

~$ opam install owl
The following actions will be performed:
  ∗  install owl 0.2.6

=-=- Gathering sources =-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
[owl] Archive in cache

=-=- Processing actions -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
[ERROR] The compilation of owl failed at "ocaml setup.ml -build".
Processing  1/1: [owl: ocaml]
#=== ERROR while installing owl.0.2.6 =========================================#
# opam-version 1.2.2
# os           darwin
# command      ocaml setup.ml -build
# path         /Users/mabrams/.opam/4.04.2/build/owl.0.2.6
# compiler     4.04.2
# exit-code    1
# env-file     /Users/mabrams/.opam/4.04.2/build/owl.0.2.6/owl-98640-6d07ae.env
# stdout-file  /Users/mabrams/.opam/4.04.2/build/owl.0.2.6/owl-98640-6d07ae.out
# stderr-file  /Users/mabrams/.opam/4.04.2/build/owl.0.2.6/owl-98640-6d07ae.err
### stdout ###
# Undefined symbols for architecture x86_64:
# [...]
#       construction vtable for std::__1::basic_iostream<char, std::__1::char_traits<char> >-in-std::__1::basic_stringstream<char, std::__1::char_traits<char>, std::__1::allocator<char> > in libeigen.a(eigen_tensor.o)
#   "virtual thunk to std::__1::basic_iostream<char, std::__1::char_traits<char> >::~basic_iostream()", referenced from:
#       construction vtable for std::__1::basic_iostream<char, std::__1::char_traits<char> >-in-std::__1::basic_stringstream<char, std::__1::char_traits<char>, std::__1::allocator<char> > in libeigen.a(eigen_dsmat.o)
#       construction vtable for std::__1::basic_iostream<char, std::__1::char_traits<char> >-in-std::__1::basic_stringstream<char, std::__1::char_traits<char>, std::__1::allocator<char> > in libeigen.a(eigen_tensor.o)
# ld: symbol(s) not found for architecture x86_64
# collect2: error: ld returned 1 exit status
# File "caml_startup", line 1:
# Error: Error during linking
# Command exited with code 2.
### stderr ###
# E: Failure("Command ''/Users/mabrams/.opam/4.04.2/bin/ocamlbuild' lib/libowl_stubs.a lib/dllowl_stubs.so lib/owl.cma lib/owl.cmxa lib/owl.a lib/owl.cmxs lib/neural/owl_neural.cma lib/neural/owl_neural.cmxa lib/neural/owl_neural.a lib/neural/owl_neural.cmxs perftest/perf_common.cma perftest/perf_common.cmxa perftest/perf_common.a perftest/perf_common.cmxs perftest/perf_dense_real.native perftest/perf_sparse_real.native perftest/perf_sparse_complex.native perftest/perf_dense_ndarray.native perftest/perf_dense_matrix.native perftest/perf_sparse_ndarray.native perftest/perf_sparse_matrix.native perftest/perf_specific_fun.native examples/test_neural.native -use-ocamlfind -tag debug' terminated with error code 10")



=-=- Error report -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=  🐫
The following actions failed
  ∗  install owl 0.2.6
No changes have been performed

The new spec code in surf and mesh disables set_azimuth and set_altitude because it causes the azimuth and altitude fields in page to be ignored.

These fields and functions are recent additions. The functions set azimuth and altitude fields in page that had been given default values by _create_page. Now the fields in the page are set from within surf and mesh by

  p.altitude <- _get_altitude spec 33.;
  p.azimuth <- _get_azimuth spec 115.;

so if these values were changed by the set_ functions, the values are ignored. Maybe something like this would be better?

  p.altitude <- _get_altitude spec p.altitude;
  p.azimuth <- _get_azimuth spec p.azimuth;

_create_page is setting the defaults to 33. and 115. respectively, so the behavior shouldn't change in other respects.

I thought about removing the set_ functions, but I think it's useful to be able to set the altitude and azimuth and have the new values apply by default to a series of surf or mesh calls.

CORRECTION: That doesn't really work because each subplot or plot page has a new page structure, so the azimuth and altitude don't persist between calls to surf and mesh, unless one is writing multiple 3D plot calls on one subplot, which would be rare enough that it wouldn't matter.

So maybe set_altitude and set_azimuth should be removed, and the default values left in hardcoded surf and mesh?

a nitpick: the map operation is not idiomatic numpy; normally you would just do np.sin(array) + 1, as you're surely aware. maybe it would be instructive adding that result for comparison?

another interesting performance test may be something that involves taking slices and iterating over them along various dimensions, forcing iteration over non-contiguous memory?

@jzstark The Plot.image seems flawed. I used the MNIST dataset to print out some numbers for testing, the image seems flipped. A simple example you can test as follows:

let z0 = Mat.create 100 100 200. in
let z1 = Mat.zeros 100 100 in
let z = Mat.(z0 @|| z1) in
Plot.image z;;

The code should generate a figure (100 x 200) with left part is white and right part is black. However, the image function generates a flipped figure.

Hello!
Do you think if there is a possibility to make plplot dependency optional?
Currently, installation of plplot OCaml module (and, thus, owl) terminates with an error if the user does not have plplot software installed.
plplot pre-build packages aren't provided for many of the popular distros (for example, in Archlinux the package is provided in AUR only), so it might be troublesome to ask the users to perform an extra step.

I don't know how useful it would be in general to have this function available. At present my only use for it is that I'd like to experiment with a Moore-Penrose pseudoinverse to help me understand its role in a proof I'm reading. Presumably there are more serious use cases. I thought I would record the absence of the function as an issue, but can't claim any urgency for implementing it. (As for my need to experiment, I think Matlab has the function, and I can learn enough Matlab to experiment with it if necessary.)

I just wanted to say thank you for this library, @ryanrhymes! 👍 🥇 👍

It's amazing to see someone work on a numeric computation library for OCaml, and the initial performance benchmarks look super promising!

It'd be fantastic if this project could be published to OPAM.
Happy to get that setup if it's something you'd like done @ryanrhymes

i played around with owl today - great job on making pulling a docker image easy.

overall, i find it impressively comprehensive already! i have some comments, coming from a scientist with a background in numpy/scipy etc. i list them below. my main concern is that i'm unsure if the design is sufficiently extensible to be a general purpose numerical library.

• the filter function returns tuples, but setter functions do not accept tuples. would it make sense to abandon curried get i j in favor of get (i,j) ? one would lose partial application but it seems otherwise maybe more consistent?

• reshape in Ndarray copies - most of the time i would use it to do a flat iteration. is this handled by the N.iter function instead?

• i did not find an easy way to make 2d Ndarray slices into matrices.

• it seems a bit Matlab-style backwards to have vectors as 1-row matrices.

• ideally there would be a seamless interoperability between 1d, 2d, nd arrays, so that 1d slices are in fact of the same type as 1d arrays etc.

• can one get a Bigarray back from a Ndarray.t?

• slices with a non-1 stride would be great.

• also, slices to select subarrays would be great.

(actually i'm not sure if doing all of these numpy things efficiently would required some blocked memory layout for arrays?)

• in numpy, indexing arrays with boolean arrays is useful. one creates a boolean array by a filtering operation and then picks these elements. how would this be done here? with filter?

• int and bool arrays can serve useful purposes. is there a way to get these? it seems Bigarray does not provide for such a thing natively...

• (it seems odd that uniform_int should return actually floats - i would throw this function out)

• the plural of axis is axes.

I'm not sure if this should be a new issue. It's closely related to #21.

I was experimenting with the 3D plotting examples in Tutorial 4. That page still uses the old operator system, but that’s just a temporary problem.

The first 3D example includes this expression:
let z0 = Mat.(sin ((x **@ 2.) +@ (y **@ 2.))) in
+@ is now + in Mat. I think that **@ refers to the function now called pow1, which raises each element of a matrix to a scalar power, but I don’t think that pow1 or pow0 (similar, but with arguments in the other order) are assigned to operators in Mat now. This might be intended, but I thought I'd mention it in case it's an oversight.