Source file vec2_D.ml

(* File: vec_SD.ml

   Copyright (C) 2001-

     Markus Mottl
     email: markus.mottl@gmail.com
     WWW: http://www.ocaml.info

     Christophe Troestler
     email: Christophe.Troestler@umons.ac.be
     WWW: http://math.umons.ac.be/anum/

   This library is free software; you can redistribute it and/or
   modify it under the terms of the GNU Lesser General Public
   License as published by the Free Software Foundation; either
   version 2.1 of the License, or (at your option) any later version.

   This library is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
   Lesser General Public License for more details.

   You should have received a copy of the GNU Lesser General Public
   License along with this library; if not, write to the Free Software
   Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
*)

open Vec4_D
open Utils
open Float64

let random ?rnd_state ?(from = -1.) ?(range = 2.) n =
  let vec = create n in
  let state =
    match rnd_state with
    | None -> Random.get_state ()
    | Some state -> state in
  for row = 1 to n do
    vec.{row} <- Random.State.float state range +. from
  done;
  if rnd_state = None then Random.set_state state;
  vec

let get_y_vec ~loc ~ofsy ~incy ~n y = get_vec loc y_str y ofsy incy n create
let get_z_vec ~loc ~ofsz ~incz ~n z = get_vec loc z_str z ofsz incz n create


(* Unary vector operations *)

let int_abs = abs

let unop direct loc =
  let loc = "Lacaml.D.Vec." ^ loc in
  fun ?n ?ofsy ?incy ?y ?ofsx ?incx x ->
    let ofsx, incx = get_vec_geom loc x_str ofsx incx in
    let ofsy, incy = get_vec_geom loc y_str ofsy incy in
    let n = get_dim_vec loc x_str ofsx incx x n_str n in
    let y = get_y_vec ~loc ~ofsy ~incy ~n y in
    direct ~n ~ofsy ~incy ~y ~ofsx ~incx ~x;
    y

external direct_abs :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dabs_stub_bc" "lacaml_Dabs_stub"

let abs = unop direct_abs "abs"

external direct_signum :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsignum_stub_bc" "lacaml_Dsignum_stub"

let signum = unop direct_signum "signum"

external direct_sqr :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsqr_stub_bc" "lacaml_Dsqr_stub"

let sqr = unop direct_sqr "sqr"

external direct_sqrt :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsqrt_stub_bc" "lacaml_Dsqrt_stub"

let sqrt = unop direct_sqrt "sqrt"

external direct_cbrt :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dcbrt_stub_bc" "lacaml_Dcbrt_stub"

let cbrt = unop direct_cbrt "cbrt"

external direct_exp :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dexp_stub_bc" "lacaml_Dexp_stub"

let exp = unop direct_exp "exp"

external direct_exp2 :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dexp2_stub_bc" "lacaml_Dexp2_stub"

let exp2 = unop direct_exp2 "exp2"

external direct_expm1 :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dexpm1_stub_bc" "lacaml_Dexpm1_stub"

let expm1 = unop direct_expm1 "expm1"

external direct_log :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dlog_stub_bc" "lacaml_Dlog_stub"

let log = unop direct_log "log"

external direct_log10 :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dlog10_stub_bc" "lacaml_Dlog10_stub"

let log10 = unop direct_log10 "log10"

external direct_log2 :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dlog2_stub_bc" "lacaml_Dlog2_stub"

let log2 = unop direct_log2 "log2"

external direct_log1p :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dlog1p_stub_bc" "lacaml_Dlog1p_stub"

let log1p = unop direct_log1p "log1p"

external direct_sin :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsin_stub_bc" "lacaml_Dsin_stub"

let sin = unop direct_sin "sin"

external direct_cos :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dcos_stub_bc" "lacaml_Dcos_stub"

let cos = unop direct_cos "cos"

external direct_tan :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dtan_stub_bc" "lacaml_Dtan_stub"

let tan = unop direct_tan "tan"

external direct_asin :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dasin_stub_bc" "lacaml_Dasin_stub"

let asin = unop direct_asin "asin"

external direct_acos :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dacos_stub_bc" "lacaml_Dacos_stub"

let acos = unop direct_acos "acos"

external direct_atan :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Datan_stub_bc" "lacaml_Datan_stub"

let atan = unop direct_atan "atan"

external direct_sinh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsinh_stub_bc" "lacaml_Dsinh_stub"

let sinh = unop direct_sinh "sinh"

external direct_cosh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dcosh_stub_bc" "lacaml_Dcosh_stub"

let cosh = unop direct_cosh "cosh"

external direct_tanh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dtanh_stub_bc" "lacaml_Dtanh_stub"

let tanh = unop direct_tanh "tanh"

external direct_asinh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dasinh_stub_bc" "lacaml_Dasinh_stub"

let asinh = unop direct_asinh "asinh"

external direct_acosh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dacosh_stub_bc" "lacaml_Dacosh_stub"

let acosh = unop direct_acosh "acosh"

external direct_atanh :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Datanh_stub_bc" "lacaml_Datanh_stub"

let atanh = unop direct_atanh "atanh"

external direct_floor :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dfloor_stub_bc" "lacaml_Dfloor_stub"

let floor = unop direct_floor "floor"

external direct_ceil :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dceil_stub_bc" "lacaml_Dceil_stub"

let ceil = unop direct_ceil "ceil"

external direct_round :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dround_stub_bc" "lacaml_Dround_stub"

let round = unop direct_round "round"

external direct_trunc :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dtrunc_stub_bc" "lacaml_Dtrunc_stub"

let trunc = unop direct_trunc "trunc"

external direct_erf :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Derf_stub_bc" "lacaml_Derf_stub"

let erf = unop direct_erf "erf"

external direct_erfc :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Derfc_stub_bc" "lacaml_Derfc_stub"

let erfc = unop direct_erfc "erfc"

external direct_logistic :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dlogistic_stub_bc" "lacaml_Dlogistic_stub"

let logistic = unop direct_logistic "logistic"

external direct_relu :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Drelu_stub_bc" "lacaml_Drelu_stub"

let relu = unop direct_relu "relu"

external direct_softplus :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsoftplus_stub_bc" "lacaml_Dsoftplus_stub"

let softplus = unop direct_softplus "softplus"

external direct_softsign :
  n : (int [@untagged]) ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  unit = "lacaml_Dsoftsign_stub_bc" "lacaml_Dsoftsign_stub"

let softsign = unop direct_softsign "softsign"


(* Binary vector operations *)

let binop direct loc =
  let loc = "Lacaml.D.Vec." ^ loc in
  fun ?n ?ofsz ?incz ?z ?ofsx ?incx x ?ofsy ?incy y ->
    let ofsz, incz = get_vec_geom loc z_str ofsz incz in
    let ofsx, incx = get_vec_geom loc x_str ofsx incx in
    let ofsy, incy = get_vec_geom loc y_str ofsy incy in
    let n = get_dim_vec loc x_str ofsx incx x n_str n in
    check_vec loc y_str y (ofsy + (n - 1) * int_abs incy);
    let z = get_z_vec ~loc ~ofsz ~incz ~n z in
    direct ~n ~ofsz ~incz ~z ~ofsx ~incx ~x ~ofsy ~incy ~y;
    z

external direct_pow :
  n : (int [@untagged]) ->
  ofsz : (int [@untagged]) ->
  incz : (int [@untagged]) ->
  z : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  unit = "lacaml_Dpow_stub_bc" "lacaml_Dpow_stub"

let pow = binop direct_pow "pow"

external direct_atan2 :
  n : (int [@untagged]) ->
  ofsz : (int [@untagged]) ->
  incz : (int [@untagged]) ->
  z : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  unit = "lacaml_Datan2_stub_bc" "lacaml_Datan2_stub"

let atan2 = binop direct_atan2 "atan2"

external direct_hypot :
  n : (int [@untagged]) ->
  ofsz : (int [@untagged]) ->
  incz : (int [@untagged]) ->
  z : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  unit = "lacaml_Dhypot_stub_bc" "lacaml_Dhypot_stub"

let hypot = binop direct_hypot "hypot"

external direct_min2 :
  n : (int [@untagged]) ->
  ofsz : (int [@untagged]) ->
  incz : (int [@untagged]) ->
  z : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  unit = "lacaml_Dmin2_stub_bc" "lacaml_Dmin2_stub"

let min2 = binop direct_min2 "min2"

external direct_max2 :
  n : (int [@untagged]) ->
  ofsz : (int [@untagged]) ->
  incz : (int [@untagged]) ->
  z : vec ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  ofsy : (int [@untagged]) ->
  incy : (int [@untagged]) ->
  y : vec ->
  unit = "lacaml_Dmax2_stub_bc" "lacaml_Dmax2_stub"

let max2 = binop direct_max2 "max2"


(* Misc functions *)

external direct_log_sum_exp :
  n : (int [@untagged]) ->
  ofsx : (int [@untagged]) ->
  incx : (int [@untagged]) ->
  x : vec ->
  (float [@unboxed])
  = "lacaml_Dlog_sum_exp_vec_stub_bc" "lacaml_Dlog_sum_exp_vec_stub"

let log_sum_exp =
  let loc = "Lacaml.D.Vec.log_sum_exp" in
  fun ?n ?ofsx ?incx x ->
    let ofsx, incx = get_vec_geom loc x_str ofsx incx in
    let n = get_dim_vec loc x_str ofsx incx x n_str n in
    direct_log_sum_exp ~n ~ofsx ~incx ~x