Crate rhai_sci

Source
Expand description

Github CI Crates.io docs.rs

§About rhai-sci

This crate provides some basic scientific computing utilities for the Rhai scripting language, inspired by languages like MATLAB, Octave, and R. For a complete API reference, check the docs.

§Install

To use the latest released version of rhai-sci, add this to your Cargo.toml:

rhai-sci = "0.2.2"

§Usage

Using this crate is pretty simple! If you just want to evaluate a single line of Rhai, then you only need:

use rhai::INT;
use rhai_sci::eval;
let result = eval::<INT>("argmin([43, 42, -500])").unwrap();

If you need to use rhai-sci as part of a persistent Rhai scripting engine, then do this instead:

use rhai::{Engine, packages::Package, INT};
use rhai_sci::SciPackage;

// Create a new Rhai engine
let mut engine = Engine::new();

// Add the rhai-sci package to the new engine
engine.register_global_module(SciPackage::new().as_shared_module());

// Now run your code
let value = engine.eval::<INT>("argmin([43, 42, -500])").unwrap();

§Features

FeatureDefaultDescription
metadataDisabledEnables exporting function metadata and is necessary for running doc-tests on Rhai examples.
ioEnabledEnables the read_matrix function but pulls in several additional dependencies (polars, url, temp-file, csv-sniffer, minreq).
nalgebraEnabledEnables several functions (regress, inv, mtimes, horzcat, vertcat, repmat, svd, hessenberg, and qr) but brings in the nalgebra and linregress crates.
randEnabledEnables the rand function for generating random FLOAT values and random matrices, but brings in the rand crate.

§API

This package provides a large variety of functions to help with scientific computing: acosd   acoshd   acot   acotd   acoth   acothd   acsc   acscd   acsch   acschd   argmax   argmin   asec   asecd   asech   asechd   asind   asinhd   assert   assert_approx_eq   assert_eq   assert_ne   atand   atanhd   bounds   cart2pol   cart2sph   cosd   coshd   cot   cotd   coth   cothd   csc   cscd   csch   cschd   cummax   cummin   cumprod   cumsum   cumtrapz   deg2rad   diag   diff   eigs   eye   flatten   fliplr   flipud   hessenberg   horzcat   hypot   interp1   intersect   inv   iqr   is_column_vector   is_float_list   is_int_list   is_list   is_matrix   is_numeric_array   is_numeric_list   is_row_vector   linspace   logspace   mad   max   maxk   mean   median   meshgrid   min   mink   mode   movmad   movmax   movmean   movmedian   movmin   movprod   movstd   movsum   movvar   mtimes   ndims   nnz   numel   ones   pol2cart   prctile   prod   qr   rad2deg   rand   read_matrix   regress   repmat   rms   rot90   sec   secd   sech   sechd   sind   sinhd   size   sph2cart   std   sum   svd   tand   tanhd   transpose   trapz   union   unique   variance   vertcat   zeros  

§physical constants

Physical constants useful for science.

§pi: FLOAT

The ratio of a circle’s circumference to its diameter (non-dimensional).

assert_eq(pi, 3.14159265358979323846264338327950288);

§c: FLOAT

The speed of light in meters per second (m/s).

assert_eq(c, 299792458.0);

§e: FLOAT

Euler’s number (non-dimensional).

assert_eq(e, 2.71828182845904523536028747135266250);

§g: FLOAT

The acceleration due to gravity on Earth in meters per second per second (m/s^2).

assert_eq(g, 9.80665);

§h: FLOAT

The Planck constant in Joules per Hertz (J/Hz).

assert_eq(h, 6.62607015e-34);

§phi: FLOAT

The golden ratio (non-dimensional).

assert_eq(phi, 1.61803398874989484820);

§G: FLOAT

The Newtonian gravitational constant (non-dimensional).

assert_eq(G, 6.6743015e-11);

§acosd(x: f64) -> f64

Returns the inverse cosine in degrees

assert_eq(acosd(-1.0), 180.0);
assert_eq(acosd(0.0), 90.0);
assert_eq(acosd(1.0), 0.0);

§acoshd(x: f64) -> f64

Returns the inverse hyperbolic cosine in degrees

assert_eq(acoshd(1.0), 0.0)
assert_eq(acoshd(10.0), 180.0/pi*acosh(10.0))

§acot(x: f64) -> f64

Returns the inverse of the cotangent in radians

assert_eq(acot(1.0), pi/4);
assert_eq(acot(-1.0), -pi/4);
assert_eq(acot(0.0), pi/2);

§acotd(x: f64) -> f64

Returns the inverse of the cotangent in degrees

assert_eq(acotd(1.0), 45.0);
assert_eq(acotd(-1.0), -45.0);
assert_eq(acotd(0.0), 90.0);

§acoth(x: f64) -> f64

Returns the inverse hyperbolic cotangent of the argument in radians

assert_eq(acoth(1.0), atanh(1.0));
assert_eq(acoth(-1.0), atanh(-1.0));

§acothd(x: f64) -> f64

Returns the inverse hyperbolic cotangent of the argument in degrees

assert_eq(acothd(1.0), atanhd(1.0));
assert_eq(acothd(-1.0), atanhd(-1.0));

§acsc(x: f64) -> f64

Returns the inverse cosecant in radians

assert_eq(acsc(-1.0), -pi/2)
assert_eq(acsc(inf), 0.0)
assert_eq(acsc(1.0), pi/2)

§acscd(x: f64) -> f64

Returns the inverse cosecant in degrees

assert_eq(acscd(-1.0), -90.0)
assert_eq(acscd(inf), 0.0)
assert_eq(acscd(1.0), 90.0)

§acsch(x: f64) -> f64

Returns the inverse hyperbolic cosecant in radians

assert_eq(acsch(inf), 0.0)
assert_eq(acsch(1.0), asinh(1.0))
assert_eq(acsch(-1.0), asinh(-1.0))

§acschd(x: f64) -> f64

Returns the inverse hyperbolic cosecant in degrees

assert_eq(acschd(inf), 0.0)
assert_eq(acschd(1.0), asinhd(1.0))
assert_eq(acschd(-1.0), asinhd(-1.0))

§argmax(arr: Array) -> Dynamic

Return the index of the largest array element. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [1, 2, 3];
let m = argmax(data);
assert_eq(m, 2);

§argmin(arr: Array) -> Dynamic

Return the index of the smallest array element. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [1, 2, 3];
let m = argmin(data);
assert_eq(m, 0);

§asec(x: f64) -> f64

Returns the inverse secant in radians

assert_eq(asec(1.0), 0.0);
assert_eq(asec(-1.0), pi);
assert_eq(asec(0.5), acos(2.0));

§asecd(x: f64) -> f64

Returns the inverse secant in degrees

assert_eq(asecd(1.0), 0.0);
assert_eq(asecd(-1.0), 180.0);
assert_eq(asecd(0.5), acosd(2.0));

§asech(x: f64) -> f64

Returns the inverse hyperbolic secant in radians

assert_eq(asech(1.0), 0.0);
assert_eq(asech(0.5), acosh(2.0));
assert_eq(asech(0.1), acosh(10.0));

§asechd(x: f64) -> f64

Returns the inverse hyperbolic secant of the argument in degrees

assert_eq(asechd(1.0), 0.0);

§asind(x: f64) -> f64

Returns the inverse sine in degrees

assert_eq(asind(-1.0), -90.0);
assert_eq(asind(0.0), 0.0);
assert_eq(asind(1.0), 90.0);

§asinhd(x: f64) -> f64

Returns the inverse hyperbolic sine in degrees

assert_eq(asinhd(0.0), 0.0)
assert_eq(asinhd(10.0), 180.0/pi*asinh(10.0))

§assert(comparison: bool) -> bool

Assert that a statement is true and throw an error if it is not.

assert(2==2);

§assert_approx_eq

§assert_approx_eq(lhs: Array, rhs: Array) -> bool

Assert that two arrays are approximately equal and return an error if they are not. Use the default tolerance of 1e-10 for the comparison.

assert_approx_eq([2.0, 2.0], [2.0, 2.000000000000000001]);

§assert_approx_eq(lhs: f64, rhs: f64) -> bool

Assert that two floats are approximately equal and return an error if they are not. Use the default tolerance of 1e-10 for the comparison.

assert_approx_eq(2.0, 2.000000000000000001);

§assert_approx_eq(lhs: Array, rhs: Array, eps: f64) -> bool

Assert that two arrays are approximately equal (within eps) and return an error if they are not.

assert_approx_eq([2.0, 2.0], [2.0, 2.000000000000000001], 1e-10);

§assert_approx_eq(lhs: f64, rhs: f64, eps: f64) -> bool

Assert that two floats are approximately equal (within eps) and return an error if they are not.

assert_approx_eq(2.0, 2.000000000000000001, 1e-10);

§assert_eq(lhs: Dynamic, rhs: Dynamic) -> bool

Assert that two arguments are equal and throw an error if they are not.

assert_eq(2, 2);

§assert_ne(lhs: Dynamic, rhs: Dynamic) -> bool

Assert that two arguments are unequal and throw an error if they are not.

assert_ne(2, 1);

§atand

§atand(x: f64) -> f64

Returns the inverse tangent in degrees

assert_approx_eq(atand(-1.0), -45.0);
assert_eq(atand(0.0), 0.0);
assert_approx_eq(atand(1.0), 45.0);

§atand(x: f64, y: f64) -> f64

Returns the inverse tangent in degrees , taking two arguments as input.

assert_approx_eq(atand(-1.0, 1.0), -45.0);
assert_eq(atand(0.0, 1.0), 0.0);
assert_approx_eq(atand(1.0, 1.0), 45.0);

§atanhd(x: f64) -> f64

Returns the inverse hyperbolic tangent in degrees

assert_eq(atanhd(0.0), 0.0)
assert_eq(atanhd(10.0), 180.0/pi*atanh(10.0))

§bounds(arr: Array) -> Array

Return the highest value from an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let high_and_low = bounds([2, 3, 4, 5]);
assert_eq(high_and_low, [2, 5]);

§cart2pol

§cart2pol(x: f64, y: f64) -> Array

Convert the argument from 2D Cartesian coordinates to polar coordinates.

assert_eq(cart2pol(1.0, 1.0), [pi/4, sqrt(2.0)])

§cart2pol(x: f64, y: f64, z: f64) -> Array

Convert the argument from 3D Cartesian coordinates to polar coordinates.

assert_eq(cart2pol(1.0, 1.0, 1.0), [pi/4, sqrt(2.0), 1.0])

§cart2sph(x: f64, y: f64, z: f64) -> Array

Convert the argument from 3D Cartesian coordinates to spherical coordinates.

assert_approx_eq(cart2sph(1.0, 0.0, 1.0), [0.0, pi/4, sqrt(2.0)])

§cosd(degrees: f64) -> f64

Returns the cosine of an argument given in degrees

assert_eq(cosd(0.0), 1.0);
assert_approx_eq(cosd(90.0), 0.0);
assert_eq(cosd(180.0), -1.0);
assert_approx_eq(cosd(270.0), 0.0);

§coshd(degrees: f64) -> f64

Returns the hyperbolic cosine of the argument given in degrees

assert_eq(coshd(0.0), 1.0)
assert_eq(coshd(10.0), cosh(10.0*pi/180.0))

§cot(radians: f64) -> f64

Returns the cotangent of the argument given in radians

assert_approx_eq(cot(pi/4), 1.0, 1e-10);
assert_approx_eq(cot(pi/2), 0.0, 1e-10);
assert_approx_eq(cot(3*pi/4), -1.0, 1e-10);

§cotd(degrees: f64) -> f64

Returns the cotangent of the argument given in degrees

assert_approx_eq(cotd(45.0), 1.0, 1e-10);
assert_approx_eq(cotd(90.0), 0.0, 1e-10);
assert_approx_eq(cotd(135.0), -1.0, 1e-10);

§coth(radians: f64) -> f64

Returns the hyperbolic cotangent of the argument given in radians

assert_approx_eq(coth(1.0), cosh(1.0)/sinh(1.0), 1e-10);
assert_approx_eq(coth(0.5), cosh(0.5)/sinh(0.5), 1e-10);
assert_approx_eq(coth(0.1), cosh(0.1)/sinh(0.1), 1e-10);

§cothd(degrees: f64) -> f64

Returns the hyperbolic cotangent of the argument given in degrees

assert_approx_eq(cothd(1.0), coshd(1.0)/sinhd(1.0), 1e-10);
assert_approx_eq(cothd(0.5), coshd(0.5)/sinhd(0.5), 1e-10);
assert_approx_eq(cothd(0.1), coshd(0.1)/sinhd(0.1), 1e-10);

§csc(radians: f64) -> f64

Returns the cosecant of the argument given in radians

assert_eq(csc(-pi/2), -1.0)
assert_eq(csc(0.0), inf)
assert_eq(csc(pi/2), 1.0)

§cscd(degrees: f64) -> f64

Returns the cosecant of the argument given in degrees

assert_eq(cscd(-90.0), -1.0)
assert_eq(cscd(0.0), inf)
assert_eq(cscd(90.0), 1.0)

§csch(radians: f64) -> f64

Returns the hyperbolic cosecant of the argument given in radians

assert_eq(csch(0.0), inf)
assert_eq(csch(10.0), 1.0/sinh(10.0))
assert_eq(csch(pi/2), 1.0/sinh(pi/2))

§cschd(degrees: f64) -> f64

Returns the hyperbolic cosecant of the argument given in degrees

assert_eq(cschd(0.0), inf)
assert_eq(cschd(10.0), 1.0/sinhd(10.0))
assert_eq(cschd(90.0), 1.0/sinhd(90.0))

§cummax(arr: Array) -> Array

Returns an array representing the cumulative maximum of a 1-D array.

let arr = [1, 4, 5, 3, 9, 8];
let c = cummax(arr);
assert_eq(c, [1, 4, 5, 5, 9, 9]);

§cummin(arr: Array) -> Array

Returns an array representing the cumulative minimum of a 1-D array.

let arr = [8, 9, 3, 5, 4, 1];
let c = cummin(arr);
assert_eq(c, [8, 8, 3, 3, 3, 1]);

§cumprod(arr: Array) -> Array

Returns an array representing the cumulative product of a 1-D array.

let arr = [1, 2, 3, 4, 5];
let c = cumprod(arr);
assert_eq(c, [1, 2, 6, 24, 120]);

§cumsum(arr: Array) -> Array

Returns an array representing the cumulative product of a 1-D array.

let arr = [1.1, 2.5, 3.4];
let c = cumsum(arr);
assert_eq(c, [1.1, 3.6, 7.0]);

§cumtrapz

§cumtrapz(y: Array) -> Array

Returns the cumulative approximate integral of the curve defined by Y and x using the trapezoidal method. Assumes unit spacing in the x direction.

let y = [1, 2, 3];
let c = cumtrapz(y);
assert_eq(c, [0.0, 1.5, 4.0]);    

§cumtrapz(x: Array, y: Array) -> Array

Returns the cumulative approximate integral of the curve defined by Y and x using the trapezoidal method.

let y = [1, 2, 3];
let x = [1, 2, 3];
let c = cumtrapz(x, y);
assert_eq(c, [0.0, 1.5, 4.0]);    

§deg2rad(degrees: f64) -> f64

Converts the argument from degrees to radians

assert_eq(deg2rad(180.0), pi);

§diag(matrix: Array) -> Array

This function can be used in two distinct ways.

  1. If the argument is an 2-D array, diag returns an array containing the diagonal of the array.
  2. If the argument is a 1-D array, diag returns a matrix containing the argument along the diagonal and zeros elsewhere.
 let matrix = [[1, 2, 3],
               [4, 5, 6],
               [7, 8, 9]];
 let d = diag(matrix);
 assert_eq(d, [1, 5, 9]);
 let diagonal = [1.0, 2.0, 3.0];
 let matrix = diag(diagonal);
 assert_eq(matrix, [[1.0, 0.0, 0.0],
                    [0.0, 2.0, 0.0],
                    [0.0, 0.0, 3.0]]);

§diff(arr: Array) -> Array

Returns the difference between successive elements of a 1-D array.

let arr = [2, 5, 1, 7, 8];
let d = diff(arr);
assert_eq(d, [3, -4, 6, 1]);

§eigs(matrix: Array) -> Map

Calculate the eigenvalues and eigenvectors for a matrix. Specifically, the output is an object map with entries for real_eigenvalues, imaginary_eigenvalues, eigenvectors, and residuals.

let matrix = eye(5);
let eig = eigs(matrix);
assert(sum(eig.residuals) < 0.000001);
let matrix = [[ 0.0,  1.0],
              [-2.0, -3.0]];
let eig = eigs(matrix);
assert(sum(eig.residuals) < 0.000001);

§eye

§eye(n: Dynamic) -> Array

Returns an identity matrix. If argument is a single number, then the output is a square matrix. The argument can also be an array specifying the dimensions separately.

let matrix = eye(3);
assert_eq(matrix, [[1.0, 0.0, 0.0],
                   [0.0, 1.0, 0.0],
                   [0.0, 0.0, 1.0]]);
let matrix = eye([3, 4]);
assert_eq(matrix, [[1.0, 0.0, 0.0, 0.0],
                   [0.0, 1.0, 0.0, 0.0],
                   [0.0, 0.0, 1.0, 0.0]]);

§eye(nx: i64, ny: i64) -> Array

Returns the identity matrix, specifying the number of rows and columns separately.

let matrix = eye(3, 4);
assert_eq(matrix, [[1.0, 0.0, 0.0, 0.0],
                   [0.0, 1.0, 0.0, 0.0],
                   [0.0, 0.0, 1.0, 0.0]]);

§flatten(matrix: Array) -> Array

Returns the contents of a multidimensional array as a 1-D array.

let matrix = ones(3, 5);
let flat = flatten(matrix);
assert_eq(len(flat), 15);
let matrix = [[1.0, 2.0, 3.0], [1.0]];
let flat = flatten(matrix);
assert_eq(len(flat), 4);

§fliplr(matrix: Array) -> Array

Flip a matrix left-to-right

let matrix = fliplr([[1.0, 0.0],
                     [0.0, 2.0]]);
assert_eq(matrix, [[0.0, 1.0],
                   [2.0, 0.0]]);

§flipud(matrix: Array) -> Array

Flip a matrix up-down

let matrix = flipud([[1.0, 0.0],
                     [0.0, 2.0]]);
assert_eq(matrix, [[0.0, 2.0],
                   [1.0, 0.0]]);

§hessenberg(matrix: Array) -> Map

Calculates the QR decomposition of a matrix

let matrix = eye(5);
let h_results = hessenberg(matrix);
assert_eq(h_results, #{"h": eye(5), "q": eye(5)});

§horzcat(matrix1: Array, matrix2: Array) -> Array

Concatenate two arrays horizontally.

let arr1 = eye(3);
let arr2 = eye(3);
let combined = horzcat(arr1, arr2);
assert_eq(size(combined), [3, 6]);

§hypot(x: f64, y: f64, z: f64) -> f64

Extends the built-in hypot function to compute distance in 3D cartesian space

assert_eq(hypot(2.0, 3.0, 6.0), 7.0);

§interp1(x: Array, y: Array, xq: Dynamic) -> f64

Given reference data, perform linear interpolation.

Both arrays must be sorted and have the same length.

Out-of-bound xq values are clamped to the minimum and maximum values of y respectively.

let x = [0, 1];
let y = [1, 2];
let xq = 0.5;
let yq = interp1(x, y, xq);
assert_eq(yq, 1.5);

§intersect(arr1: Array, arr2: Array) -> Array

Performs set intersection of two arrays

 let set1 = [7, 1, 7, 7, 4];
 let set2 = [7, 0, 4, 4, 0];
let x = intersect(set1, set2);
assert_eq(x, [4, 7]);

§inv(matrix: Array) -> Array

Calculates the inverse of a matrix. Fails if the matrix if not invertible, or if the elements of the matrix aren’t FLOAT or INT.

let x = [[ 1.0,  0.0,  2.0],
         [-1.0,  5.0,  0.0],
         [ 0.0,  3.0, -9.0]];
let x_inverted = inv(x);
assert_eq(x_inverted, [[0.8823529411764706,  -0.11764705882352941,   0.19607843137254902],
                       [0.17647058823529413,  0.17647058823529413,   0.0392156862745098 ],
                       [0.058823529411764705, 0.058823529411764705, -0.09803921568627451]]
);
let x = [[1, 2],
         [3, 4]];
let x_inverted = inv(x);
assert_eq(x_inverted, [[-2.0, 1.0],
                       [1.5, -0.5]]
);

§iqr(arr: Array) -> f64

Returns the inter-quartile range for a 1-D array.

let data = [1, 1, 1, 1, 1, 1, 1, 5, 6, 9, 9, 9, 9, 9, 9, 9, 9];
let inter_quartile_range = iqr(data);
assert_eq(inter_quartile_range, 8.0);

§is_column_vector(arr: Array) -> bool

Tests whether the input is a column vector

let x = ones([5, 1]);
assert_eq(is_column_vector(x), true)
let x = ones([5, 5]);
assert_eq(is_column_vector(x), false)

§is_float_list(arr: Array) -> bool

Tests whether the input in a simple list array composed of floating point values.

let x = [1.0, 2.0, 3.0, 4.0];
assert_eq(is_float_list(x), true)
let x = [1, 2, 3, 4];
assert_eq(is_float_list(x), false)

§is_int_list(arr: Array) -> bool

Tests whether the input in a simple list array composed of integer values.

let x = [1.0, 2.0, 3.0, 4.0];
assert_eq(is_int_list(x), false)
let x = [1, 2, 3, 4];
assert_eq(is_int_list(x), true)

§is_list(arr: Array) -> bool

Tests whether the input in a simple list array

let x = [1, 2, 3, 4];
assert_eq(is_list(x), true);
let x = [[[1, 2], [3, 4]]];
assert_eq(is_list(x), false);

§is_matrix(arr: Array) -> bool

Tests whether the input is a matrix

let x = ones([3, 5]);
assert_eq(is_matrix(x), true)
let x = ones([5, 5, 5]);
assert_eq(is_matrix(x), false)

§is_numeric_array(arr: Array) -> bool

Determines if the entire array is numeric (ints or floats).

let x = [1, 2, 3.0, 5.0];
assert_eq(is_numeric_array(x), true);
let x = [1, 2, 3.0, 5.0, "a"];
assert_eq(is_numeric_array(x), false);

§is_numeric_list(arr: Array) -> bool

Tests whether the input in a simple list array composed of either floating point or integer values.

let x = [1.0, 2.0, 3.0, 4.0];
assert_eq(is_numeric_list(x), true)
let x = [1, 2, 3, 4];
assert_eq(is_numeric_list(x), true)
let x = ["a", "b", "c", "d"];
assert_eq(is_numeric_list(x), false)

§is_row_vector(arr: Array) -> bool

Tests whether the input is a row vector

let x = ones([1, 5]);
assert_eq(is_row_vector(x), true)
let x = ones([5, 5]);
assert_eq(is_row_vector(x), false)

§linspace(x1: Dynamic, x2: Dynamic, n: i64) -> Array

Returns an array containing a number of elements linearly spaced between two bounds.

let x = linspace(1, 2, 5);
assert_eq(x, [1.0, 1.25, 1.5, 1.75, 2.0]);

§logspace(a: Dynamic, b: Dynamic, n: i64) -> Array

Returns an array containing a number of elements logarithmically spaced between two bounds.

let x = logspace(1, 3, 3);
assert_eq(x, [10.0, 100.0, 1000.0]);

§mad(arr: Array) -> Dynamic

Returns the median absolute deviation of a 1-D array.

let data = [1.0, 2.0, 3.0, 3.0, 4.0, 4.0, 4.0, 5.0, 5.5, 6.0, 6.0, 6.5, 7.0, 7.0, 7.5, 8.0, 9.0, 12.0, 52.0, 90.0];
let m = mad(data);
assert_eq(m, 2.0);

§max

§max(arr: Array) -> Dynamic

Return the highest value from an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let the_highest_number = max([2, 3, 4, 5]);
assert_eq(the_highest_number, 5);
let the_highest_number = max([2, 3.0, 4.12, 5]);
assert_eq(the_highest_number, 5.0);

§max(a: Dynamic, b: Dynamic) -> Dynamic

Return the highest value from a pair of numbers. Fails if the numbers are anything other than INT or FLOAT.

let the_higher_number = max(2, 3);
assert_eq(the_higher_number, 3);
let the_higher_number = max(2.0, 3.0);
assert_eq(the_higher_number, 3.0);

§maxk(arr: Array, k: i64) -> Array

Returns the k highest values from an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [32, 15, -7, 10, 1000, 41, 42];
let mk = maxk(data, 3);
assert_eq(mk, [41, 42, 1000]);
let data = [32, 15, -7.0, 10, 1000, 41.0, 42];
let mk = maxk(data, 3);
assert_eq(mk, [41.0, 42.0, 1000.0]);

§mean(arr: Array) -> Dynamic

Return the average of an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [1, 2, 3];
let m = mean(data);
assert_eq(m, 2.0);

§median(arr: Array) -> Dynamic

Returns the variance of a 1-D array.

let data = [1, 1, 1, 1, 2, 5, 6, 7, 8];
let m = median(data);
assert_eq(m, 2.0);

§meshgrid(x: Array, y: Array) -> Map

Returns an object map containing 2-D grid coordinates based on the uni-axial coordinates contained in arguments x and y.

let x = [1, 2];
let y = [3, 4];
let g = meshgrid(x, y);
assert_eq(g, #{"x": [[1, 2],
                     [1, 2]],
               "y": [[3, 3],
                     [4, 4]]});

§min

§min(arr: Array) -> Dynamic

Return the lowest value from an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let the_lowest_number = min([2, 3, 4, 5]);
assert_eq(the_lowest_number, 2);
let the_lowest_number = min([2, 3.0, 4.12, 5]);
assert_eq(the_lowest_number, 2.0);

§min(a: Dynamic, b: Dynamic) -> Dynamic

Return the lowest value from a pair of numbers. Fails if the numbers are anything other than INT or FLOAT.

let the_lower_number = min(2, 3);
assert_eq(the_lower_number, 2);
let the_lower_number = min(2.0, 3.0);
assert_eq(the_lower_number, 2.0);

§mink(arr: Array, k: i64) -> Array

Return the k lowest values in an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [32, 15, -7, 10, 1000, 41, 42];
let mk = mink(data, 3);
assert_eq(mk, [-7, 10, 15]);
let data = [32, 15.1223232, -7, 10, 1000.00000, 41, 42];
let mk = mink(data, 3);
assert_eq(mk, [-7.0, 10.0, 15.1223232]);

§mode(arr: Array) -> Dynamic

Returns the mode of a 1-D array.

let data = [1, 2, 2, 2, 2, 3];
let m = mode(data);
assert_eq(m, 2);
let data = [1.0, 2.0, 2.0, 2.0, 2.0, 3.0];
let m = mode(data);
assert_eq(m, 2.0);

§movmad(arr: Array, k: i64) -> Array

Returns an array of the moving maximum absolute deviation (with a given width) across the input array.

let data = [1, 2, 4, -1, -2, -3, -1, 3, 2, 1];
let m = movmad(data, 3);
assert_eq(m, [0.5, 1.0, 2.0, 1.0, 1.0, 1.0, 2.0, 1.0, 1.0, 0.5]);

§movmax(arr: Array, k: i64) -> Array

Returns an array of the moving maximum (with a given width) across the input array.

let data = [1, 2, 4, -1, -2, -3, -1, 3, 2, 1];
let m = movmax(data, 3);
assert_eq(m, [2, 4, 4, 4, -1, -1, 3, 3, 3, 2]);

§movmean(arr: Array, k: i64) -> Array

Returns an array of the moving average (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movmean(data, 3);
assert_eq(m, [1.5, 2.0, 3.0, 4.0, 5.0, 5.5]);

§movmedian(arr: Array, k: i64) -> Array

Returns an array of the moving median (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movmedian(data, 3);
assert_eq(m, [1.5, 2.0, 3.0, 4.0, 5.0, 5.5]);

§movmin(arr: Array, k: i64) -> Array

Returns an array of the moving minimum (with a given width) across the input array.

let data = [1, 2, 4, -1, -2, -3, -1, 3, 2, 1];
let m = movmin(data, 3);
assert_eq(m, [1, 1, -1, -2, -3, -3, -3, -1, 1, 1]);

§movprod(arr: Array, k: i64) -> Array

Returns an array of the moving product (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movprod(data, 3);
assert_eq(m, [2, 6, 24, 60, 120, 30]);

§movstd(arr: Array, k: i64) -> Array

Returns an array of the moving standard deviation (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movstd(data, 3);
assert_eq(m, [0.7071067811865476, 1.0, 1.0, 1.0, 1.0, 0.7071067811865476]);

§movsum(arr: Array, k: i64) -> Array

Returns an array of the moving sum (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movsum(data, 3);
assert_eq(m, [3, 6, 9, 12, 15, 11]);

§movvar(arr: Array, k: i64) -> Array

Returns an array of the moving variance (with a given width) across the input array.

let data = [1, 2, 3, 4, 5, 6];
let m = movvar(data, 3);
assert_eq(m, [0.5, 1.0, 1.0, 1.0, 1.0, 0.5]);

§mtimes(matrix1: Array, matrix2: Array) -> Array

Perform matrix multiplication.

let a = eye(3);
let b = ones(3);
let c = mtimes(a, b);
assert_eq(b, c);

§ndims(matrix: Array) -> i64

Return the number of dimensions in matrix, passed by reference.

let matrix = ones(4, 6);
let n = ndims(matrix);
assert_eq(n, 2);

§nnz(matrix: Array) -> i64

Returns the number of non-zero elements in a matrix, passed by reference.

let matrix = ones(4, 6);
let n = nnz(matrix);
assert_eq(n, 24);
let matrix = eye(4);
let n = nnz(matrix);
assert_eq(n, 4);

§numel(matrix: Array) -> i64

Returns the number of elements in a matrix, passed by reference.

let matrix = ones(4, 6);
let n = numel(matrix);
assert_eq(n, 24);
let matrix = [1, [1, 2, 3]];
let n = numel(matrix);
assert_eq(n, 4);

§ones

§ones(n: Dynamic) -> Array

Return a matrix of ones. Can be called with a single integer argument (indicating the square matrix of that size) or with an array argument (indicating the size for each dimension).

let matrix = ones(3);
assert_eq(matrix, [[1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0]]);
let matrix = ones([3, 3]);
assert_eq(matrix, [[1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0]]);
let matrix = ones([3]);
assert_eq(matrix, [1.0, 1.0, 1.0]);
let matrix = ones([3, 3, 3]);
assert_eq(matrix, [[[1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0]],
                   [[1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0]],
                   [[1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0],
                    [1.0, 1.0, 1.0]]]);

§ones(nx: i64, ny: i64) -> Array

Return a matrix of ones. Arguments indicate the number of rows and columns in the matrix.

let matrix = ones(3, 3);
assert_eq(matrix, [[1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0],
                   [1.0, 1.0, 1.0]]);

§pol2cart

§pol2cart(theta: f64, r: f64) -> Array

Convert the argument from 2D polar coordinates to Cartesian coordinates.

assert_approx_eq(pol2cart(pi/4, sqrt(2.0)), [1.0, 1.0])

§pol2cart(theta: f64, r: f64, z: f64) -> Array

Convert the argument from 3D polar coordinates to Cartesian coordinates.

assert_approx_eq(pol2cart(pi/4, sqrt(2.0), 1.0), [1.0, 1.0, 1.0])

§prctile(arr: Array, p: Dynamic) -> f64

Returns a given percentile value for a 1-D array of data.

The array must not be empty.

If the percentile value is <= 0 or >= 100, returns the minimum and maximum values of the array respectively.

let data = [1, 2, 0, 3, 4];
let p = prctile(data, 0);
assert_eq(p, 0.0);
let data = [1, 2, 0, 3, 4];
let p = prctile(data, 50);
assert_eq(p, 2.0);
let data = [1, 2, 0, 3, 4];
let p = prctile(data, 100);
assert_eq(p, 4.0);

§prod(arr: Array) -> Dynamic

Compute the product of an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [1, 2, 3];
let m = prod(data);
assert_eq(m, 6);
let data = [3, 6, 10];
let m = prod(data);
assert_eq(m, 180);

§qr(matrix: Array) -> Map

Calculates the QR decomposition of a matrix

let matrix = eye(5);
let qr_results = qr(matrix);
assert_eq(qr_results, #{"q": eye(5), "r": eye(5)});

§rad2deg(radians: f64) -> f64

Converts the argument from radians to degrees

assert_eq(rad2deg(pi), 180.0);

§rand

§rand() -> f64

Returns a random number between zero and one.

let r = rand();
assert(r >= 0.0 && r <= 1.0);

§rand(n: Dynamic) -> Array

Returns a matrix of random values, each between zero and one. Can be called with a single integer argument (indicating the square matrix of that size) or with an array argument (indicating the size for each dimension).

let matrix = rand(3);
assert_eq(size(matrix), [3, 3]);
let matrix = rand([3, 3]);
assert_eq(size(matrix), [3, 3]);

§rand(nx: i64, ny: i64) -> Array

Return a matrix of random values, each between zero and one. Arguments indicate the number of rows and columns in the matrix.

let matrix = rand(3, 3);
assert_eq(size(matrix), [3, 3]);

§read_matrix(file_path: String) -> Array

Reads a numeric csv file from a url

let url = "https://raw.githubusercontent.com/plotly/datasets/master/diabetes.csv";
let x = read_matrix(url);
assert_eq(size(x), [768, 9]);

§regress(x: Array, y: Array) -> Map

Performs ordinary least squares regression and provides a statistical assessment.

let x = [[1.0, 0.0],
         [1.0, 1.0],
         [1.0, 2.0]];
let y = [[0.1],
         [0.8],
         [2.1]];
let b = regress(x, y);
assert_eq(b,  #{"parameters": [-2.220446049250313e-16, 1.0000000000000002],
                "pvalues": [1.0, 0.10918255350924745],
                "standard_errors": [0.11180339887498947, 0.17320508075688767]});

§repmat(matrix: Array, nx: i64, ny: i64) -> Array

Repeats copies of a matrix

let matrix = eye(3);
let combined = repmat(matrix, 2, 2);
assert_eq(size(combined), [6, 6]);

§rms(arr: Array) -> Dynamic

Returns the variance of a 1-D array.

let data = [1, 2, 3, 4, 5];
let r = rms(data);
assert_eq(r, 3.3166247903554);

§rot90

§rot90(matrix: Array) -> Array

Rotate a matrix counterclockwise once

let matrix = rot90([[1.0, 0.0],
                   [0.0, 2.0]]);
assert_eq(matrix, [[0.0, 2.0],
                  [1.0, 0.0]]);

§rot90(matrix: Array, k: i64) -> Array

Rotate a matrix counterclockwise k times

let matrix = rot90([[1.0, 0.0],
                    [0.0, 2.0]], 2);
assert_eq(matrix, [[2.0, 0.0],
                   [0.0, 1.0]]);

§sec(radians: f64) -> f64

Returns the secant of the argument given in radians

assert_eq(sec(0.0), 1.0);
assert_eq(sec(pi/2), 1/cos(pi/2));
assert_eq(sec(pi), -1.0);
## `secd(degrees: f64) -> f64`
Returns the secant of the argument given in degrees
```typescript
assert_eq(secd(0.0), 1.0);
assert_eq(secd(90.0), 1/cosd(90.0));
assert_eq(secd(180.0), -1.0);

§sech(radians: f64) -> f64

Returns the hyperbolic secant of the argument given in radians

assert_eq(sech(0.0), 1.0);
assert_eq(sech(10.0), 1.0/cosh(10.0));
assert_eq(sech(pi/2), 1.0/cosh(pi/2));

§sechd(degrees: f64) -> f64

Returns the hyperbolic secant of the argument given in degrees

assert_eq(sechd(0.0), 1.0);
assert_eq(sechd(10.0), 1.0/coshd(10.0));
assert_eq(sechd(90.0), 1.0/coshd(90.0));

§sind(degrees: f64) -> f64

Returns the sine of an argument given in degrees

assert_eq(sind(0.0),  0.0);
assert_eq(sind(90.0), 1.0);
assert_approx_eq(sind(180.0), 0.0);
assert_eq(sind(270.0), -1.0);

§sinhd(degrees: f64) -> f64

Returns the hyperbolic sine of the argument given in degrees

assert_eq(sinhd(0.0), 0.0)
assert_eq(sinhd(10.0), sinh(10.0*pi/180.0))

§size(matrix: Array) -> Array

Returns an array indicating the size of the matrix along each dimension, passed by reference.

let matrix = ones(3, 5);
assert_eq(size(matrix), [3, 5]);
let matrix = [[[1, 2]]];
assert_eq(size(matrix), [1, 1, 2]);

§sph2cart(azimuth: f64, elevation: f64, r: f64) -> Array

Convert the argument from spherical coordinates to 3D Cartesian coordinates.

assert_approx_eq(sph2cart(0.0, pi/4, sqrt(2.0)), [1.0, 0.0, 1.0])

§std(arr: Array) -> Dynamic

Returns the standard deviation of a 1-D array.

let data = [1, 2, 3];
let v = std(data);
assert_eq(v, 1.0);

§sum(arr: Array) -> Dynamic

Sum an array. Fails if the input is not an array, or if it is an array with elements other than INT or FLOAT.

let data = [1, 2, 3];
let m = sum(data);
assert_eq(m, 6);
let data = [1, 2.0, 3];
let m = sum(data);
assert_eq(m, 6.0);

§svd(matrix: Array) -> Map

Calculates the singular value decomposition of a matrix

let matrix = eye(5);
let svd_results = svd(matrix);
assert_eq(svd_results, #{"s": ones([5]), "u": eye(5), "v": eye(5)});

§tand(degrees: f64) -> f64

Returns the tangent of an argument given in degrees

assert_approx_eq(tand(-45.0), -1.0);
assert_eq(tand(0.0), 0.0);
assert_approx_eq(tand(45.0), 1.0);

§tanhd(degrees: f64) -> f64

Returns the hyperbolic tangent of the argument given in degrees

assert_eq(tanhd(0.0), 0.0)
assert_eq(tanhd(10.0), tanh(10.0*pi/180.0))

§transpose(matrix: Array) -> Array

Transposes a matrix.

let row = [[1, 2, 3, 4]];
let column = transpose(row);
assert_eq(column, [[1],
                   [2],
                   [3],
                   [4]]);
let matrix = transpose(eye(3));
assert_eq(matrix, eye(3));

§trapz

§trapz(arr: Array) -> Dynamic

Returns the approximate integral of the curve defined by y using the trapezoidal method. Assumes that x-values have unit spacing.

let y = [1.0, 1.5, 2.0];
let A = trapz(y);
assert_eq(A, 3.0);
let y = [1, 2, 3];
let A = trapz(y);
assert_eq(A, 4.0);

§trapz(x: Array, y: Array) -> Dynamic

Returns the approximate integral of the curve defined by y and x using the trapezoidal method.

let y = [1.0, 1.5, 2.0];
let x = [1.0, 2.0, 3.0];
let A = trapz(x, y);
assert_eq(A, 3.0);
let y = [1, 2, 3];
let x = [1, 2, 3];
let A = trapz(x, y);
assert_eq(A, 4.0);

§union(arr1: Array, arr2: Array) -> Array

Returns the set union of two arrays.

let set1 = [7, 1, 7, 7, 4];
let set2 = [7, 0, 4, 4, 0];
let u = union(set1, set2);
assert_eq(u, [0, 1, 4, 7]);

§unique(arr: Array) -> Array

Returns an array of the unique elements in an array.

let data = [1, 2, 2, 2, 5, 4, 4, 2, 5, 8];
let u = unique(data);
assert_eq(u, [1, 2, 4, 5, 8]);

§variance(arr: Array) -> Dynamic

Returns the variance of a 1-D array.

let data = [1, 2, 3];
let v = variance(data);
assert_eq(v, 1.0);

§vertcat(matrix1: Array, matrix2: Array) -> Array

Concatenates two array vertically.

let arr1 = eye(3);
let arr2 = eye(3);
let combined = vertcat(arr1, arr2);
assert_eq(size(combined), [6, 3]);

§zeros

§zeros(n: Dynamic) -> Array

Return a matrix of zeros. Can be called with a single integer argument (indicating the square matrix of that size) or with an array argument (indicating the size for each dimension).

let matrix = zeros(3);
assert_eq(matrix, [[0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0]]);
let matrix = zeros([3, 3]);
assert_eq(matrix, [[0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0]]);
let matrix = zeros([3]);
assert_eq(matrix, [0.0, 0.0, 0.0]);
let matrix = zeros([3, 3, 3]);
assert_eq(matrix, [[[0.0, 0.0, 0.0],
                    [0.0, 0.0, 0.0],
                    [0.0, 0.0, 0.0]],
                   [[0.0, 0.0, 0.0],
                    [0.0, 0.0, 0.0],
                    [0.0, 0.0, 0.0]],
                   [[0.0, 0.0, 0.0],
                    [0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0]]]);

§zeros(nx: i64, ny: i64) -> Array

Return a matrix of zeros. Arguments indicate the number of rows and columns in the matrix.

let matrix = zeros(3, 3);
assert_eq(matrix, [[0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0],
                   [0.0, 0.0, 0.0]]);

Modules§

assert_functions
constant_definitions
cum_functions
int_and_diff
matrix_functions
misc_functions
moving_functions
set_functions
stats
trig_functions
validation_functions

Structs§

SciPackage
Package for scientific computing

Enums§

FOIL
Matrix compatibility conditions

Functions§

array_to_vec_float
array_to_vec_int
eval
This provides the ability to easily evaluate a line (or lines) of code without explicitly setting up a script engine
if_int_convert_to_float_and_do
If the input is an int, convert to a float and do the function. if the input is a float already, the function is still performed.
if_int_do_else_if_array_do
if_list_convert_to_vec_float_and_do
if_list_do
Does a function if the input is a list, otherwise throws an error.
if_list_do_int_or_do_float
if_matrices_and_compatible_convert_to_vec_array_and_do
if_matrix_convert_to_vec_array_and_do
If the input is a
if_matrix_do
int_and_float_totals
omatrix_to_vec_dynamic
ovector_to_vec_dynamic