Compare commits

..

1 commit
main ... reps

Author SHA1 Message Date
Kiana Sheibani 4eb0c27203
Update ipkg 2023-09-25 23:49:21 -04:00
25 changed files with 289 additions and 776 deletions

View file

@ -27,22 +27,14 @@ combine the most useful features of each library.
## Documentation ## Documentation
Most of the exported utility functions have docstrings. There is also Most of the exported utility functions have docstrings. There are also plans for
a (currently unfinished) guide on how to use NumIdr [here](docs/Contents.md). an in-depth tutorial on NumIdr's features, though that is not available yet.
## Usage ## Usage
To install using idris2 directly: Basic package install:
``` sh ``` shell
git clone https://github.com/kiana-S/numidr
cd numidr
idris2 --install numidr.ipkg idris2 --install numidr.ipkg
``` ```
Or you can install using pack:
``` sh
pack install numidr
```

View file

@ -1,19 +0,0 @@
# NumIdr Guide
Welcome to NumIdr's guide! This guide will serve as a basic introduction to the core ideas of NumIdr, as well as a reference for common functions and utilities that you might otherwise miss.
If you're familiar with the Python library [NumPy](https://numpy.org/) or anything based on it, then a lot of the early content will be familiar, as NumIdr is based on the same array type.
> [!NOTE]
> This guide assumes a basic understanding of linear algebra, and it will not explain any linear algebra concepts necessary for understanding NumIdr's structure.
### Primary Features
1. [Fundamental Data Types](DataTypes.md)
2. [Basic Operations on Arrays](Operations.md)
3. [Working with Vectors and Matrices](VectorsMatrices.md)
4. [Transforms](Transforms.md)
### Advanced Concepts
5. Array Representations (coming soon!)

View file

@ -1,133 +0,0 @@
# Fundamental Data Types
NumIdr exports a number of different datatypes. The most important type and the cornerstone of the library is the _array_, but there are other useful types as well.
## Arrays
### What is an Array?
In most programming languages, the word "array" is used to mean a one-dimensional list of values that is contiguous in memory. A typical array of integers may be written in list form like this:
```
[1, 4, 10, 2, -5, 18]
```
In this kind of array, elements are indexed by a single integer, starting at zero and increasing from left to right.
NumIdr, however, uses the word a bit more generally: a NumIdr array is a multi-dimensional structure that can be indexed by any number of integers. NumIdr arrays are written as nested lists:
```
[[4, -9, -2],
[5, -6, 1]]
```
Unlike in other languages, however, this is not a nested structure. The above is a single array, and it is always manipulated as one object.
### Properties of Arrays
The `Array` datatype has the following parameters:
```idris
Array : (s : Vect rk Nat) -> (a : Type) -> Type
```
The first parameter is the _shape_, a list of numbers (the _dimensions_) where each dimension is the length of a particular _axis_ of the array. The second parameter is the _element type_, the type of the values inside the array.
Let's return to the array example from earlier:
```
[[4, -9, -2],
[5, -6, 1]]
```
This is a rank-2 array, meaning that it has two axes. Rank-2 arrays are typically called matrices. To determine the dimensions of the array, we count the size of each nested list from the outside in, which in the case of matrices means the row axis comes before the column axis. This matrix has 2 rows and 3 columns, making its shape `[2, 3]`. Thus, a possible type for this array could be `Array [2, 3] Int`.
When determining the index of a value inside the array, the order of the indices is the same as the order of the dimensions, and each index number counts from zero. For example, the index `[1, 0]` indicates the second row and first column, which contains `5`.
> [!NOTE]
> The word "dimensions" is often ambiguously used to either refer to the rank of an array
> (as in "multi-dimensional array" in the previous section), or to the lengths of its
> axes. Conventionally, NumIdr reserves "dimension" for the second meaning, and uses
> "rank" for the first meaning.
## Types of Arrays
Arrays are loosely divided into multiple subtypes mostly based on their rank. Each array subtype has an alias for convenience.
### Scalars
A scalar is a rank-0 array, meaning that it is indexed by 0 integers. Its alias is `Scalar`:
```idris
Scalar : (a : Type) -> Type
Scalar a = Array [] a
```
A scalar has exactly one index, the empty list `[]`. This means that it is exactly the same as a single value and as such is largely pointless, but NumIdr still provides an alias for it just in case you need it.
### Vectors
A vector is a rank-1 array:
```idris
Vector : (n : Nat) -> (a : Type) -> Type
Vector n a = Array [n] a
```
A vector's type signature and stored data is effectively identical to that of the standard library type `Vect`, whose elements are confusingly also called "vectors"; we often refer to those as "vects" to differentiate.
### Matrices
As mentioned before, a matrix is a rank-2 array:
```idris
Matrix : (m, n : Nat) -> (a : Type) -> Type
Matrix m n a = Array [m, n] a
```
There is also an alias `Matrix'` for square matrices.
```idris
Matrix' : (n : Nat) -> (a : Type) -> Type
Matrix' n a = Array [n, n] a
```
As a linear algebra library, the majority of the operations in NumIdr revolve around matrices.
#### Homogeneous Matrices
NumIdr also provides aliases for homogeneous matrices:
```idris
HMatrix : (m, n : Nat) -> (a : Type) -> Type
HMatrix m n a = Array [S m, S n] a
HMatrix' : (n : Nat) -> (a : Type) -> Type
HMatrix' n a = Array [S n, S n] a
-- To use with homogeneous matrices
HVector : (n : Nat) -> (a : Type) -> Type
HVector n a = Array [S n] a
```
These are useful for clarity when working with both homogeneous and non-homogeneous matrices.
## Other Datatypes
### Transforms
A transform is a wrapper type for a square matrix with certain properties that can be used to transform points in space.
```idris
Transform : (ty : TransType) -> (n : Nat) -> (a : Type) -> Type
```
The `TransType` parameter dictates what kind of transform it is. More information on transforms and their operations can be found in the [Transforms](Transforms.md) chapter.
### Permutations
The type `Permutation n` represents a permutation of `n` elements. Permutations are mostly used internally for various algorithms, but they are also an input in some operations, such as permuting the rows or columns of a matrix.
Permutations can be composed using `(*.)`, and a permutation can be converted into a matrix using `permuteM`.
[Contents](Contents.md) | [Next](Operations.md)

View file

@ -1,179 +0,0 @@
# Basic Operations on Arrays
> [!CAUTION]
> Arrays and their associated functions are not intended to be evaluated at compile-time. If you try to compute an array in the REPL, you will not get the output you expect!
>
> If you really need to use the REPL to test array code, use `:exec`.
## Constructing Arrays
The most important array constructor is `array`, which returns an array of the specified values:
```idris
array [[1, 2, 3], [4, 5, 6]]
```
Scalars, vectors and matrices have their own constructors, used in exactly the same way (`scalar`, `vector`, and `matrix`). These should be used instead of `array` wherever possible, as they provide more information to the type-checker.
There are also a few other, more basic constructors for convenience.
```idris
-- A 2x2x3 array filled with zeros
repeat 0 [2, 2, 3]
-- Same as previous
zeros [2, 2, 3]
-- A 2x2x3 array filled with ones
ones [2, 2, 3]
```
## Accessing Array Properties
There are a few simple functions for accessing basic properties of arrays: `shape` and `rank`, which are self-explanatory, and `size`, which returns the total number of elements in the array.
> [!TIP]
> The `shape` accessor is sufficient for most uses, but it can cause problems with the type-checker, as for an array `arr : Array s a` the type checker does not know that `shape arr` and `s` are equal. To solve this problem, you can use the `ShapeView`:
>
> ```idris
> example {s} arr with (viewShape arr)
> _ | Shape s = ...
> ```
## Indexing Arrays
NumIdr provides multiple different indexing functions for different purposes. These functions can be grouped based on these categories:
**Operation**
- **Access** - Accesses and returns elements from the array.
- **Update** - Returns a new array with the specified element or range updated by a function. These are indicated by an `Update` suffix.
- **Set** - Returns a new array with the specified element or range set to new values. These are indicated by a `Set` suffix.
**Range**
- **Default** - Operates on a single array element.
- **Range** - Operates on multiple elements at once. These are indicated by a `Range` suffix.
**Safety**
- **Safe** - Guarantees through its type that the index is within range by requiring each index to be a `Fin` value.
- **Non-Bounded** - Does not guarantee through its type that the index is within range, and returns `Nothing` if the provided index is out of bounds. These are indicated by an `NB` suffix.
- **Unsafe** - Does not perform any bounds checks at all. These are indicated by an `Unsafe` suffix. **Only use these if you really know what you are doing!**
Not all combinations of these categories are defined by the library. Here are the currently provided indexing functions:
| | Safe | Ranged Safe | Non-Bounded | Ranged Non-Bounded | Unsafe | Ranged Unsafe |
|------------|-----------------|------------------------|-------------------|--------------------------|-----------------------|------------------------------|
| **Access** | `index`, `(!!)` | `indexRange`, `(!!..)` | `indexNB`, `(!?)` | `indexRangeNB`, `(!?..)` | `indexUnsafe`, `(!#)` | `indexRangeUnsafe`, `(!#..)` |
| **Update** | `indexUpdate` | `indexUpdateRange` | `indexUpdateNB` | | | |
| **Set** | `indexSet` | `indexSetRange` | `indexSetNB` | | | |
The accessor functions have operator forms for convenience.
### Specifying Coordinates
When indexing a single element, a coordinate is specified as a list of numbers, where each number is either a `Fin` value for safe indexing or a `Nat` value for non-bounded and unsafe indexing.
```idris
index [1, 0] arr
-- Equivalent operator form
arr !! [1, 0]
```
With ranged indexing, a sub-array of the original array is accessed or modified. This sub-array is given by a list of _range specifiers_, one for each axis, which can be one of the following:
- `Bounds x y` - Every index from `x` (inclusive) to `y` (exclusive)
- `StartBound x` - Every index from `x` to the end of the axis
- `EndBound y` - Every index from the start of the axis to `y`
- `One i`, `One' i` - The single index `i`
- `All` - Every index in the axis
- `Indices xs` - A list of indices in the axis
- `Filter p` - A predicate specifying which indices to access or modify
The range specifier `One` is special, as it is the only range specifier that decreases the rank of the resulting array. For example, the following matrix access returns the first column as a vector due to `One` being used:
```idris
matrix [[2, 0], [1, 2]] !!.. [All, One 0]
== vector [2, 1]
```
`One'` does not decrease the rank of the result in this way.
## Standard Interface Methods
The array type implements a number of standard interfaces. Most of these implementations are unremarkable, but a few have caveats that are worth noting.
### Numeric Operations
Arrays implement the numeric interfaces (`Num`, `Neg`, `Fractional`), as well as `Semigroup` and `Monoid`, if their element type supports those operations. These functions are computed elementwise.
```idris
matrix [[1, 1], [2, 5]] + matrix [[2, 3], [-1, 3]]
== matrix [[3, 4], [1, 8]]
```
This elementwise behavior holds for `(+)`, `(*)`, `(-)`, `(/)`, and `(<+>)`. **`(*)` is not matrix multiplication!** For the generalized multiplication operator, which includes matrix multiplication, see [the next chapter](VectorsMatrices.md).
> [!WARNING]
> Due to unfortunate restrictions in Idris's standard `Num` interface, `(+)` and `(*)` can only be used when the array's shape is available at run-time. If this is not the case, you must use `zipwith (+)` or `zipWith (*)` instead.
### `Foldable` and `Traversable`
When folding or traversing the elements of an array, these elements are ordered in row-major or "C-style" order, which corresponds to the order in which elements are written and displayed. This behavior should not be depended on, however, as it can change based on the internal [array representation](Representations.md); use the `elements` function if you specifically want row-major order.
## Other Common Operations
### Concatenation and Stacking
Two arrays can be concatenated along an axis (`concat`), so long as all other axes have the same dimensions. Two matrices being concatenated along the row axis requires that they must have the same number of columns.
```idris
-- 0 is the first axis i.e. the row axis
concat 0 (matrix [[1, 2], [3, 4]]) (matrix [[5, 6], [7, 8]])
== matrix [[1, 2],
[3, 4],
[5, 6],
[7, 8]]
```
Stacking (`stack`) is similar to concatenation, but slightly different. Stacking combines arrays with the exact same shape into a single array that is one rank higher. For example, vectors can be stacked along the row axis to obtain a matrix whose rows are the original vectors.
```idris
stack 0 [vector [1, 2], vector [3, 4]]
== matrix [[1, 2],
[3, 4]]
```
There are also specialized functions for operating on vectors and matrices: `(++)` for concatenating vectors; `vconcat` and `hconcat` for concatenating two matrices vertically (by row) and horizontally (by column) respectively, and `vstack` and `hstack` for stacking row vectors and column vectors respectively.
### Reshaping and Resizing
An array can be "reshaped" into any other shape, so long as the total number of elements is the same. This reshaping is done by arranging the elements into a linear order before inserting them into the new array. As with folding, the default order is row-major.
```idris
reshape [3, 2] (vector [1, 2, 3, 4, 5, 6])
== matrix [[1, 2],
[3, 4],
[5, 6]]
```
Arrays can also be resized with `resize`, which changes their shape while keeping every element at the same index. A default element must be provided to fill any indices that did not exist in the original array.
```idris
resize [2, 4] 10 (matrix [[1, 2],
[3, 4],
[5, 6]])
== matrix [[1, 2, 10, 10],
[3, 4, 10, 10]]
```
Instead of the `resize` function, one can also use the `resizeLTE` function, which does not require a default element, but only works if the new array would be provably smaller than the original one.
### Transpose
The `transpose` function reverses the axis order of an array: For `arr : Array [3,2,4] Int`, we have `transpose arr : Array [4,2,3] Int`. For matrices, this corresponds to the usual definition of switching rows and columns. There is also a postfix form `(.T)`.
For more fine-grained control when rearranging arrays, there are the `swapAxes` and `permuteAxes` functions, where the first swaps only two axes and the second takes an arbitrary [permutation](DataTypes.md#Permutations). There are also `swapInAxis` and `permuteInAxis`, which permute inside an axis, e.g. swapping rows or columns in a matrix.
Like with concatenation and stacking, the swap and permute functions have forms specific to vectors and matrices: `swapCoords` and `permuteCoords` for vectors, and `swapRows`, `permuteRows`, `swapColumns`, and `permuteColumns` for matrices.
[Previous](DataTypes.md) | [Contents](Contents.md) | [Next](VectorsMatrices.md)

View file

@ -1,2 +0,0 @@
# Array Representations

View file

@ -1,108 +0,0 @@
# Transforms
In code that works with linear algebra, it is common to use matrices as _transformations_, functions that take in a vector and output a new vector. These matrices can often be divided into categories based on the operations they perform, such as a rotation matrix or an affine transformation matrix.
The `Transform` wrapper exists to encode these differing properties on the type level, as well as to provide extra utilities for working with matrices in this fashion.
## Types of Transforms
`Transform` has the following type signature:
```idris
Transform : (ty : TransType) -> (n : Nat) -> (a : Type) -> Type
```
A `Transform ty n a` is a wrapper over an `HMatrix' n a`. The `ty` parameter is the transform type, which dictates what properties the transform has. These eight options are currently available:
**Affine Types:**
- `TAffine`
- `TIsometry` (rotation + reflection + translation)
- `TRigid` (rotation + translation)
- `TTranslation`
**Linear Types:**
- `TLinear`
- `TOrthonormal` (rotation + reflection)
- `TRotation`
- `TTrivial` (always the `identity`)
The capital T at the beginning of each of these names identifies it as a `TransType` value. To make working with transforms smoother, NumIdr provides synonyms for transforms of each type: for example, `Isometry n a` is a synonym for `Transform TIsometry n a`.
### Linear and Affine
Transform types are divided into linear and affine types. Linear transforms must preserve the origin point, whereas affine transforms do not have this restriction.
Linear and affine transform types are in a one-to-one correspondence: a linear transform can be converted to and from an affine transform by adding or removing a translation component.
```
Linear <-> Affine
Orthonormal <-> Isometry
Rotation <-> Rigid
Trivial <-> Translation
```
The `setTranslation` and `linearize` functions perform these conversions.
For simplicity, both categories of transform are wrappers over homogeneous matrices, even though linear transforms could be represented by non-homogeneous matrices.
### Transform Type Casting
Some transform types can be cast into other types. For example, a `Rotation` can be cast into an `Orthonormal`, as all rotation matrices are orthonormal.
```idris
rt : Rotation 3 Double
cast rt : Orthonormal 3 Double
```
In the diagram from the previous section, lower types can be cast into types higher up. Each linear type (on the left) can also be cast into the corresponding affine type (on the right).
## Constructing Transforms
There are multiple ways to construct transforms, either by wrapping a matrix or directly through constructor functions.
For each transform type, `fromHMatrix` can be used to test if a homogeneous matrix satisfies the right properties, and converts it into a transform if it does. The `Rotation`, `Orthonormal` and `Linear` types also have `fromMatrix` for non-homogeneous matrices.
> [!WARNING]
> The `fromHMatrix` and `fromMatrix` constructors use exact equality comparisons when testing matrices, which can be a problem if your element type is `Double` or a similar inexact number type. To fix this issue, NumIdr provides a `WithEpsilon` named implementation that defines equality approximately and allows you to specify an epsilon value.
>
> ```idris
> fromHMatrix @{WithEpsilon 1.0e-6} mat
> ```
There are also direct constructors for transforms, which are often more convenient as they do not have the possibility of failing. There are too many of these constructors to exhaustively list here, so I encourage you to look through the functions in the `Data.NumIdr.Transform.*` modules to see what is available.
## Multiplication with Transforms
Like most objects in NumIdr, transforms multiply with the generalized multiplication operator `(*.)`, and `identity` and `inverse` can also be used with transforms.
> [!NOTE]
> There is no `tryInverse` function for transforms. This is because all transforms are required to be invertible, so there isn't any need for it; you can safely use `inverse` in all circumstances.
Transforms of any types can be multiplied. When two transforms of different types are multiplied, the resulting transform type is determined by taking the most specific type that both original types can be cast to. For example, an `Orthonormal` transform multiplied by a `Translation` returns an `Isometry`.
### The Point Type
Transforms behave differently from regular matrices when applied to vectors. When an affine transform is applied in this way, it is first linearized, so that vectors only have linear transforms applied to them. **This is not a bug!**
In order to properly apply affine transforms, the `Point` type must be used, which is a wrapper around the `Vector` type that supports these transforms. A point can be constructed with the `point` function, which is used exactly the same as the `vector` constructor.
```idris
point [4, 3, 6]
```
Points support most basic operations that vectors do, including indexing operations and standard library methods. However, a point cannot be added to another point. Instead, a vector must be added to a point:
```idris
(+.) : Vector n a -> Point n a -> Point n a
(.+) : Point n a -> Vector n a -> Point n a
(-.) : Point n a -> Point n a -> Vector n a
```
To remember the distinction between the two addition operators, the dot is always on the side of the point, not the vector.
This separation between points and vectors is intended to make working with affine transformations more convenient, as it mirrors the separation between points and vectors in affine algebra. These may feel like arbitrary restrictions, but you might be surprised by how convenient they are to work with!
[Previous](VectorsMatrices.md) | [Contents](Contents.md)

View file

@ -1,114 +0,0 @@
# Working with Vectors and Matrices
As linear algebra is one of the main concerns of NumIdr, most of its provided functions are dedicated to vectors (rank-1 arrays) and matrices (rank-2 arrays).
## The Generalized Multiplication Operator
A linear algebra library wouldn't be very useful without matrix multiplication! Since `(*)` is already used for element-wise multiplication, NumIdr defines a new interface `Mult` that can accept and return values of different types:
```idris
interface Mult a b c where
(*.) : a -> b -> c
```
The generalized multiplication operator `(*.)` covers matrix multiplication, scalar-vector multiplication, and any other linear algebra operation that's vaguely multiplication-like.
## Vectors
### Algebraic Operations
Vectors can be added together with element-wise addition `(+)`. Scalar-vector multiplication is done with the generalized multiplication operator `(*.)`.
```idris
2 *. (vector [1, 1] + vector [2, 3])
== vector [6, 8]
```
A few other basic linear algebra operations are available:
- `dot`, The dot product
- `cross`, The cross product
- `perp`, The perpendicular product (sometimes called the 2D cross product)
- `triple`, The scalar triple product
### Indexing
NumIdr provides special versions of `index` and `indexNB` and their infix forms `(!!)` and `(!?)` for use with vectors. These take a single numeric index instead of a list.
```idris
Vector.index 2 v == index [2] v
v !! 2 == v !! [2]
```
For convenience, when working with two- or three-dimensional vectors, there are postfix accessors `(.x)`, `(.y)`, and `(.z)`:
```idris
v = vector [5, 6, 2]
v.x == 5
v.y == 6
v.z == 2
```
### Other Operations
- `toVect` - Convert a vector into a `Vect`
- `dim` - Returns the vector's length
- `(++)` - Concatenate two vectors
## Matrices
### Arithmetic Operations
Like vectors, matrices can be added together using `(+)`. Matrix multiplication, as well as matrix-vector and matrix-scalar multiplication, are performed using `(*.)`.
For the purposes of working with matrices and matrix-like objects, the sub-interfaces `MultMonoid` and `MultGroup` are defined:
```idris
interface Mult' a => MultMonoid a where
identity : a
interface MultMonoid a => MultGroup a where
inverse : a -> a
```
The `identity` function returns an identity matrix, and `inverse` calculates a matrix's inverse. Note that `inverse` cannot tell you if an inverse of your matrix does not exist; if you want to handle that possibility, use `tryInverse` instead.
```idris
tryInverse : FieldCmp a => Matrix' n a -> Maybe (Matrix' n a)
```
You can also use the `invertible` predicate to test if a matrix has an inverse.
#### LU and LUP Decomposition
The functions `decompLU` and `decompLUP` compute LU and LUP decomposition on a matrix.
```idris
decompLU : Field a => (mat : Matrix m n a) -> Maybe (DecompLU mat)
decompLUP : FieldCmp a => (mat : Matrix m n a) -> DecompLUP mat
```
`DecompLU` and `DecompLUP` are record types holding the results of the corresponding decomposition. The accessors `lower`, `upper` and `permute` can be applied to get each component of the decomposition; `lower` and `upper` return matrices, and `permute` returns a `Permutation` value.
#### Other Algebraic Operations
- `trace` - The sum of the matrix's diagonal
- `outer` - The matrix-valued outer product (or tensor product) of two vectors
- `det` - Determinant of the matrix
- `solve` - Apply an inverse matrix to a vector, useful for solving linear equations
> [!TIP]
> The `det` and `solve` operations require computing an LUP decomposition, which can be expensive. The variants `detWithLUP` and `solveWithLUP` allow an existing LUP decomposition to be passed in, which can make your code more efficient.
### Indexing
Aside from the usual array indexing functions, there are a few functions specialized to matrix indexing:
- `getRow` and `getColumn` - Returns a specific row or column of the matrix
- `diagonal` - Returns the diagonal elements of the matrix as a vector
- `minor` - Removes a single row and column from the matrix
[Previous](Operations.md) | [Contents](Contents.md) | [Next](Transforms.md)

View file

@ -1,6 +1,5 @@
package numidr package numidr
version = 0.3.0 version = 2.0.0
brief = "Linear algebra and data science library"
authors = "Kiana Sheibani" authors = "Kiana Sheibani"
license = "MIT" license = "MIT"
@ -10,8 +9,14 @@ langversion >= 0.6.0
sourcedir = "src" sourcedir = "src"
readme = "README.md" readme = "README.md"
modules = Data.Permutation, modules = Data.NP,
Data.Permutation,
Data.NumIdr, Data.NumIdr,
Data.NumIdr.PrimArray,
Data.NumIdr.PrimArray.Bytes,
Data.NumIdr.PrimArray.Boxed,
Data.NumIdr.PrimArray.Linked,
Data.NumIdr.PrimArray.Delayed,
Data.NumIdr.Array, Data.NumIdr.Array,
Data.NumIdr.Array.Array, Data.NumIdr.Array.Array,
Data.NumIdr.Array.Coords, Data.NumIdr.Array.Coords,

39
sirdi.json Normal file
View file

@ -0,0 +1,39 @@
[
{
"name": "numidr",
"version": "1.0.0",
"deps": [],
"modules": [
"Data.NP",
"Data.Permutation",
"Data.NumIdr",
"Data.NumIdr.Array",
"Data.NumIdr.Array.Array",
"Data.NumIdr.Array.Coords",
"Data.NumIdr.Array.Order",
"Data.NumIdr.Homogeneous",
"Data.NumIdr.Matrix",
"Data.NumIdr.Interfaces",
"Data.NumIdr.PrimArray",
"Data.NumIdr.Scalar",
"Data.NumIdr.Vector",
"Data.NumIdr.LArray",
"Data.NumIdr.Transform",
"Data.NumIdr.Transform.Affine",
"Data.NumIdr.Transform.Isometry",
"Data.NumIdr.Transform.Linear",
"Data.NumIdr.Transform.Orthonormal",
"Data.NumIdr.Transform.Point",
"Data.NumIdr.Transform.Rigid",
"Data.NumIdr.Transform.Rotation",
"Data.NumIdr.Transform.Transform",
"Data.NumIdr.Transform.Translation",
"Data.NumIdr.Transform.Trivial"
],
"passthru": {
"authors": "Kiana Sheibani",
"license": "MIT",
"sourceloc": "https://github.com/kiana-S/numidr"
}
}
]

43
src/Data/NP.idr Normal file
View file

@ -0,0 +1,43 @@
module Data.NP
import Data.Vect
%default total
||| A type constructor for heterogenous lists.
public export
data NP : (f : k -> Type) -> (ts : Vect n k) -> Type where
Nil : NP {k} f []
(::) : f t -> NP {k} f ts -> NP {k} f (t :: ts)
public export
(all : NP (Eq . f) ts) => Eq (NP f ts) where
(==) {all=[]} [] [] = True
(==) {all=_::_} (x :: xs) (y :: ys) = (x == y) && (xs == ys)
public export
(all : NP (Show . f) ts) => Show (NP f ts) where
show xs = "[" ++ elems xs ++ "]"
where
elems : {0 f,ts : _} -> (all : NP (Show . f) ts) => NP f ts -> String
elems {all=[]} [] = ""
elems {all=[_]} [x] = show x
elems {all=_::_} (x :: xs) = show x ++ ", " ++ elems xs
-- use this when idris can't figure out the constraint above
public export
eqNP : forall f. (forall t. Eq (f t)) => (x, y : NP f ts) -> Bool
eqNP [] [] = True
eqNP (x :: xs) (y :: ys) = x == y && eqNP xs ys
public export
head : NP f (t :: ts) -> f t
head (x :: _) = x
public export
index : (i : Fin n) -> NP {n} f ts -> f (index i ts)
index FZ (x :: xs) = x
index (FS i) (x :: xs) = index i xs

View file

@ -1,8 +1,8 @@
module Data.NumIdr module Data.NumIdr
import public Data.NP
import public Data.Permutation import public Data.Permutation
import Data.NumIdr.PrimArray
import public Data.NumIdr.Interfaces import public Data.NumIdr.Interfaces
import public Data.NumIdr.Array import public Data.NumIdr.Array
import public Data.NumIdr.Scalar import public Data.NumIdr.Scalar

View file

@ -1,5 +1,6 @@
module Data.NumIdr.Array module Data.NumIdr.Array
import public Data.NP
import public Data.NumIdr.Array.Array import public Data.NumIdr.Array.Array
import public Data.NumIdr.Array.Coords import public Data.NumIdr.Array.Coords
import public Data.NumIdr.Array.Rep import public Data.NumIdr.Array.Rep

View file

@ -3,6 +3,7 @@ module Data.NumIdr.Array.Array
import Data.List import Data.List
import Data.Vect import Data.Vect
import Data.Zippable import Data.Zippable
import Data.NP
import Data.Permutation import Data.Permutation
import Data.NumIdr.Interfaces import Data.NumIdr.Interfaces
import Data.NumIdr.PrimArray import Data.NumIdr.PrimArray
@ -12,8 +13,6 @@ import Data.NumIdr.Array.Coords
%default total %default total
||| The type of an array.
|||
||| Arrays are the central data structure of NumIdr. They are an `n`-dimensional ||| Arrays are the central data structure of NumIdr. They are an `n`-dimensional
||| grid of values, where `n` is a value known as the *rank* of the array. Arrays ||| grid of values, where `n` is a value known as the *rank* of the array. Arrays
||| of rank 0 are single values, arrays of rank 1 are vectors, and arrays of rank ||| of rank 0 are single values, arrays of rank 1 are vectors, and arrays of rank
@ -24,18 +23,22 @@ import Data.NumIdr.Array.Coords
||| array's total size. ||| array's total size.
||| |||
||| Arrays are indexed by row first, as in the standard mathematical notation for ||| Arrays are indexed by row first, as in the standard mathematical notation for
||| matrices. ||| matrices. This is independent of the actual order in which they are stored; the
||| default order is row-major, but this is configurable.
||| |||
||| @ rk The rank of the array ||| @ rk The rank of the array
||| @ s The shape of the array ||| @ s The shape of the array
||| @ a The type of the array's elements ||| @ a The type of the array's elements
export export
data Array : (s : Vect rk Nat) -> (a : Type) -> Type where data Array : (s : Vect rk Nat) -> (a : Type) -> Type where
||| Internally, arrays are stored via one of a handful of representations. ||| Internally, arrays are stored using Idris's primitive array type. This is
||| stored along with the array's shape, and a vector of values called the
||| *strides*, which determine how indexes into the internal array should be
||| performed. This is how the order of the array is configurable.
||| |||
||| @ ord The order of the elements of the array
||| @ sts The strides of the array
||| @ s The shape of the array ||| @ s The shape of the array
||| @ rep The internal representation of the array
||| @ rc A witness that the element type satisfies the representation constraint
MkArray : (rep : Rep) -> (rc : RepConstraint rep a) => (s : Vect rk Nat) -> MkArray : (rep : Rep) -> (rc : RepConstraint rep a) => (s : Vect rk Nat) ->
PrimArray rep s a @{rc} -> Array s a PrimArray rep s a @{rc} -> Array s a
@ -53,32 +56,24 @@ unsafeMkArray = MkArray
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| The shape of the array. ||| The shape of the array
export export
shape : Array {rk} s a -> Vect rk Nat shape : Array {rk} s a -> Vect rk Nat
shape (MkArray _ s _) = s shape (MkArray _ s _) = s
||| The size of the array, i.e. the total number of elements. ||| The rank of the array
export
size : Array s a -> Nat
size = product . shape
||| The rank of the array.
export export
rank : Array s a -> Nat rank : Array s a -> Nat
rank = length . shape rank = length . shape
||| The internal representation of the array.
export export
getRep : Array s a -> Rep getRep : Array s a -> Rep
getRep (MkArray rep _ _) = rep getRep (MkArray rep _ _) = rep
||| The representation constraint of the array.
export export
getRepC : (arr : Array s a) -> RepConstraint (getRep arr) a getRepC : (arr : Array s a) -> RepConstraint (getRep arr) a
getRepC (MkArray _ @{rc} _ _) = rc getRepC (MkArray _ @{rc} _ _) = rc
||| Extract the primitive backend array.
export export
getPrim : (arr : Array s a) -> PrimArray (getRep arr) s a @{getRepC arr} getPrim : (arr : Array s a) -> PrimArray (getRep arr) s a @{getRepC arr}
getPrim (MkArray _ _ pr) = pr getPrim (MkArray _ _ pr) = pr
@ -176,12 +171,22 @@ fromFunction : {default B rep : Rep} -> RepConstraint rep a =>
(s : Vect rk Nat) -> (Coords s -> a) -> Array s a (s : Vect rk Nat) -> (Coords s -> a) -> Array s a
fromFunction s f = MkArray rep s (PrimArray.fromFunction s f) fromFunction s f = MkArray rep s (PrimArray.fromFunction s f)
||| Construct an array using a structure of nested vectors. The elements are arranged
||| to the specified order before being written.
|||
||| @ s The shape of the constructed array
||| @ ord The order of the constructed array
export
array' : {default B rep : Rep} -> RepConstraint rep a =>
(s : Vect rk Nat) -> Vects s a -> Array s a
array' s v = MkArray rep s (fromVects s v)
||| Construct an array using a structure of nested vectors. ||| Construct an array using a structure of nested vectors.
||| To explicitly specify the shape and order of the array, use `array'`. ||| To explicitly specify the shape and order of the array, use `array'`.
export export
array : {default B rep : Rep} -> RepConstraint rep a => array : {default B rep : Rep} -> RepConstraint rep a =>
{s : Vect rk Nat} -> Vects s a -> Array s a {s : Vect rk Nat} -> Vects s a -> Array s a
array v = MkArray rep s (fromVects s v) array {rep} = array' {rep} _
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
@ -335,26 +340,14 @@ arr !#.. is = indexRangeUnsafe is arr
-- Operations on arrays -- Operations on arrays
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| Map a function over an array.
|||
||| You should almost always use `map` instead; only use this function if you
||| know what you are doing!
export export
mapArray' : (a -> a) -> Array s a -> Array s a mapArray' : (a -> a) -> Array s a -> Array s a
mapArray' f (MkArray rep _ arr) = MkArray rep _ (mapPrim f arr) mapArray' f (MkArray rep _ arr) = MkArray rep _ (mapPrim f arr)
||| Map a function over an array.
|||
||| You should almost always use `map` instead; only use this function if you
||| know what you are doing!
export export
mapArray : (a -> b) -> (arr : Array s a) -> RepConstraint (getRep arr) b => Array s b mapArray : (a -> b) -> (arr : Array s a) -> RepConstraint (getRep arr) b => Array s b
mapArray f (MkArray rep _ arr) @{rc} = MkArray rep @{rc} _ (mapPrim f arr) mapArray f (MkArray rep _ arr) @{rc} = MkArray rep @{rc} _ (mapPrim f arr)
||| Combine two arrays of the same element type using a binary function.
|||
||| You should almost always use `zipWith` instead; only use this function if
||| you know what you are doing!
export export
zipWithArray' : (a -> a -> a) -> Array s a -> Array s a -> Array s a zipWithArray' : (a -> a -> a) -> Array s a -> Array s a -> Array s a
zipWithArray' {s} f a b with (viewShape a) zipWithArray' {s} f a b with (viewShape a)
@ -363,10 +356,6 @@ zipWithArray' {s} f a b with (viewShape a)
$ PrimArray.fromFunctionNB @{mergeRepConstraint (getRepC a) (getRepC b)} _ $ PrimArray.fromFunctionNB @{mergeRepConstraint (getRepC a) (getRepC b)} _
(\is => f (a !# is) (b !# is)) (\is => f (a !# is) (b !# is))
||| Combine two arrays using a binary function.
|||
||| You should almost always use `zipWith` instead; only use this function if
||| you know what you are doing!
export export
zipWithArray : (a -> b -> c) -> (arr : Array s a) -> (arr' : Array s b) -> zipWithArray : (a -> b -> c) -> (arr : Array s a) -> (arr' : Array s b) ->
RepConstraint (mergeRep (getRep arr) (getRep arr')) c => Array s c RepConstraint (mergeRep (getRep arr) (getRep arr')) c => Array s c
@ -388,12 +377,6 @@ export
convertRep : (rep : Rep) -> RepConstraint rep a => Array s a -> Array s a convertRep : (rep : Rep) -> RepConstraint rep a => Array s a -> Array s a
convertRep rep (MkArray _ s arr) = MkArray rep s (convertRepPrim arr) convertRep rep (MkArray _ s arr) = MkArray rep s (convertRepPrim arr)
||| Temporarily convert an array to a delayed representation to make modifying
||| it more efficient, then convert it back to its original representation.
export
delayed : (Array s a -> Array s' a) -> Array s a -> Array s' a
delayed f arr = convertRep (getRep arr) @{getRepC arr} $ f $ convertRep Delayed arr
||| Resize the array to a new shape, preserving the coordinates of the original ||| Resize the array to a new shape, preserving the coordinates of the original
||| elements. New coordinates are filled with a default value. ||| elements. New coordinates are filled with a default value.
||| |||
@ -411,7 +394,7 @@ resize s' def arr = fromFunction {rep=getRep arr} @{getRepC arr} s' (fromMaybe d
export export
-- HACK: Come up with a solution that doesn't use `believe_me` or trip over some -- HACK: Come up with a solution that doesn't use `believe_me` or trip over some
-- weird bug in the type-checker -- weird bug in the type-checker
resizeLTE : (s' : Vect rk Nat) -> (0 ok : All Prelude.id (zipWith LTE s' s)) => resizeLTE : (s' : Vect rk Nat) -> (0 ok : NP Prelude.id (zipWith LTE s' s)) =>
Array {rk} s a -> Array s' a Array {rk} s a -> Array s' a
resizeLTE s' arr = resize s' (believe_me ()) arr resizeLTE s' arr = resize s' (believe_me ()) arr
@ -463,7 +446,6 @@ stack axis arrs = rewrite sym (lengthCorrect arrs) in
getAxisInd FZ (i :: is) = (i, is) getAxisInd FZ (i :: is) = (i, is)
getAxisInd {s=_::_} (FS ax) (i :: is) = mapSnd (i::) (getAxisInd ax is) getAxisInd {s=_::_} (FS ax) (i :: is) = mapSnd (i::) (getAxisInd ax is)
||| Join the axes of a nested array structure to form a single array.
export export
joinAxes : {s' : _} -> Array s (Array s' a) -> Array (s ++ s') a joinAxes : {s' : _} -> Array s (Array s' a) -> Array (s ++ s') a
joinAxes {s} arr with (viewShape arr) joinAxes {s} arr with (viewShape arr)
@ -478,7 +460,6 @@ joinAxes {s} arr with (viewShape arr)
dropUpTo (x::xs) (y::ys) = dropUpTo xs ys dropUpTo (x::xs) (y::ys) = dropUpTo xs ys
||| Split an array into a nested array structure along the specified axes.
export export
splitAxes : (rk : Nat) -> {0 rk' : Nat} -> {s : _} -> splitAxes : (rk : Nat) -> {0 rk' : Nat} -> {s : _} ->
Array {rk=rk+rk'} s a -> Array (take {m=rk'} rk s) (Array (drop {m=rk'} rk s) a) Array {rk=rk+rk'} s a -> Array (take {m=rk'} rk s) (Array (drop {m=rk'} rk s) a)
@ -596,8 +577,7 @@ export
-- Foldable and Traversable operate on the primitive array directly. This means -- Foldable and Traversable operate on the primitive array directly. This means
-- that their operation is dependent on the internal representation of the -- that their operation is dependent on the order of the array.
-- array.
export export
Foldable (Array s) where Foldable (Array s) where

View file

@ -4,8 +4,7 @@ import Data.Either
import Data.List import Data.List
import Data.List1 import Data.List1
import Data.Vect import Data.Vect
import Data.NP
import public Data.Vect.Quantifiers
%default total %default total
@ -26,9 +25,6 @@ export
rangeLenZ : (x : Nat) -> length (range 0 x) = x rangeLenZ : (x : Nat) -> length (range 0 x) = x
rangeLenZ x = rangeLen 0 x `trans` minusZeroRight x rangeLenZ x = rangeLen 0 x `trans` minusZeroRight x
export %unsafe
assertFin : Nat -> Fin n
assertFin n = natToFinLt n @{believe_me Oh}
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
-- Array coordinate types -- Array coordinate types
@ -39,16 +35,7 @@ assertFin n = natToFinLt n @{believe_me Oh}
||| values of `Fin dim`, where `dim` is the dimension of each axis. ||| values of `Fin dim`, where `dim` is the dimension of each axis.
public export public export
Coords : (s : Vect rk Nat) -> Type Coords : (s : Vect rk Nat) -> Type
Coords = All Fin Coords = NP Fin
-- Occasionally necessary for reasons of Idris not being great at
-- resolving interface constraints
public export
[eqCoords] Eq (Coords s) where
[] == [] = True
(x :: xs) == (y :: ys) = x == y && xs == ys
||| Forget the shape of the array by converting each index to type `Nat`. ||| Forget the shape of the array by converting each index to type `Nat`.
export export
@ -66,7 +53,6 @@ namespace Strict
public export public export
data CRange : Nat -> Type where data CRange : Nat -> Type where
One : Fin n -> CRange n One : Fin n -> CRange n
One' : Fin n -> CRange n
All : CRange n All : CRange n
StartBound : Fin (S n) -> CRange n StartBound : Fin (S n) -> CRange n
EndBound : Fin (S n) -> CRange n EndBound : Fin (S n) -> CRange n
@ -77,14 +63,12 @@ namespace Strict
public export public export
CoordsRange : (s : Vect rk Nat) -> Type CoordsRange : (s : Vect rk Nat) -> Type
CoordsRange = All CRange CoordsRange = NP CRange
namespace NB namespace NB
public export public export
data CRangeNB : Type where data CRangeNB : Type where
One : Nat -> CRangeNB
One' : Nat -> CRangeNB
All : CRangeNB All : CRangeNB
StartBound : Nat -> CRangeNB StartBound : Nat -> CRangeNB
EndBound : Nat -> CRangeNB EndBound : Nat -> CRangeNB
@ -134,7 +118,6 @@ namespace Strict
public export public export
cRangeToList : {n : Nat} -> CRange n -> Either Nat (List Nat) cRangeToList : {n : Nat} -> CRange n -> Either Nat (List Nat)
cRangeToList (One x) = Left (cast x) cRangeToList (One x) = Left (cast x)
cRangeToList (One' x) = Right [cast x]
cRangeToList All = Right $ range 0 n cRangeToList All = Right $ range 0 n
cRangeToList (StartBound x) = Right $ range (cast x) n cRangeToList (StartBound x) = Right $ range (cast x) n
cRangeToList (EndBound x) = Right $ range 0 (cast x) cRangeToList (EndBound x) = Right $ range 0 (cast x)
@ -155,8 +138,8 @@ namespace Strict
newShape : {s : _} -> (rs : CoordsRange s) -> Vect (newRank rs) Nat newShape : {s : _} -> (rs : CoordsRange s) -> Vect (newRank rs) Nat
newShape [] = [] newShape [] = []
newShape (r :: rs) with (cRangeToList r) newShape (r :: rs) with (cRangeToList r)
_ | Left _ = newShape rs newShape (r :: rs) | Left _ = newShape rs
_ | Right xs = length xs :: newShape rs newShape (r :: rs) | Right xs = length xs :: newShape rs
getNewPos : {s : _} -> (rs : CoordsRange {rk} s) -> Vect rk Nat -> Vect (newRank rs) Nat getNewPos : {s : _} -> (rs : CoordsRange {rk} s) -> Vect rk Nat -> Vect (newRank rs) Nat
@ -182,14 +165,6 @@ namespace NB
validateCRange (d :: s) (r :: rs) = [| validate' d r :: validateCRange s rs |] validateCRange (d :: s) (r :: rs) = [| validate' d r :: validateCRange s rs |]
where where
validate' : (n : Nat) -> CRangeNB -> Maybe (CRange n) validate' : (n : Nat) -> CRangeNB -> Maybe (CRange n)
validate' n (One i) =
case isLT i n of
Yes _ => Just (One (natToFinLT i))
_ => Nothing
validate' n (One' i) =
case isLT i n of
Yes _ => Just (One' (natToFinLT i))
_ => Nothing
validate' n All = Just All validate' n All = Just All
validate' n (StartBound x) = validate' n (StartBound x) =
case isLTE x n of case isLTE x n of
@ -212,44 +187,29 @@ namespace NB
export %unsafe export %unsafe
assertCRange : (s : Vect rk Nat) -> Vect rk CRangeNB -> CoordsRange s assertCRange : (s : Vect rk Nat) -> Vect rk CRangeNB -> CoordsRange s
assertCRange [] [] = [] assertCRange [] [] = []
assertCRange (d :: s) (r :: rs) = assert' r :: assertCRange s rs assertCRange (d :: s) (r :: rs) = assert' d r :: assertCRange s rs
where where
assert' : forall n. CRangeNB -> CRange n assert' : (n : Nat) -> CRangeNB -> CRange n
assert' (One i) = One (assertFin i) assert' n All = All
assert' (One' i) = One' (assertFin i) assert' n (StartBound x) = StartBound (believe_me x)
assert' All = All assert' n (EndBound x) = EndBound (believe_me x)
assert' (StartBound x) = StartBound (assertFin x) assert' n (Bounds x y) = Bounds (believe_me x) (believe_me y)
assert' (EndBound x) = EndBound (assertFin x) assert' n (Indices xs) = Indices (believe_me <$> xs)
assert' (Bounds x y) = Bounds (assertFin x) (assertFin y) assert' n (Filter f) = Filter (f . finToNat)
assert' (Indices xs) = Indices (assertFin <$> xs)
assert' (Filter f) = Filter (f . finToNat)
public export public export
cRangeNBToList : Nat -> CRangeNB -> Either Nat (List Nat) cRangeNBToList : Nat -> CRangeNB -> List Nat
cRangeNBToList s (One i) = Left i cRangeNBToList s All = range 0 s
cRangeNBToList s (One' i) = Right [i] cRangeNBToList s (StartBound x) = range x s
cRangeNBToList s All = Right $ range 0 s cRangeNBToList s (EndBound x) = range 0 x
cRangeNBToList s (StartBound x) = Right $ range x s cRangeNBToList s (Bounds x y) = range x y
cRangeNBToList s (EndBound x) = Right $ range 0 x cRangeNBToList s (Indices xs) = nub xs
cRangeNBToList s (Bounds x y) = Right $ range x y cRangeNBToList s (Filter p) = filter p $ range 0 s
cRangeNBToList s (Indices xs) = Right $ nub xs
cRangeNBToList s (Filter p) = Right $ filter p $ range 0 s
public export
newRank : Vect rk Nat -> Vect rk CRangeNB -> Nat
newRank _ [] = 0
newRank (d :: s) (r :: rs) =
case cRangeNBToList d r of
Left _ => newRank s rs
Right _ => S (newRank s rs)
||| Calculate the new shape given by a coordinate range. ||| Calculate the new shape given by a coordinate range.
public export public export
newShape : (s : Vect rk Nat) -> (is : Vect rk CRangeNB) -> Vect (newRank s is) Nat newShape : Vect rk Nat -> Vect rk CRangeNB -> Vect rk Nat
newShape [] [] = [] newShape = zipWith (length .: cRangeNBToList)
newShape (d :: s) (r :: rs) with (cRangeNBToList d r)
_ | Left _ = newShape s rs
_ | Right xs = length xs :: newShape s rs
export export
getAllCoords' : Vect rk Nat -> List (Vect rk Nat) getAllCoords' : Vect rk Nat -> List (Vect rk Nat)

View file

@ -12,10 +12,10 @@ import Data.Buffer
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| An order is an abstract representation of the way in which contiguous array ||| An order is an abstract representation of the way in which array
||| elements are stored in memory. Orders are used to calculate strides, which ||| elements are stored in memory. Orders are used to calculate strides,
||| provide a method of converting an array coordinate into a linear memory ||| which provide a method of converting an array coordinate into a linear
||| location. ||| memory location.
public export public export
data Order : Type where data Order : Type where
||| C-like order, or contiguous order. This order stores elements in a ||| C-like order, or contiguous order. This order stores elements in a
@ -49,57 +49,24 @@ calcStrides FOrder v@(_::_) = scanl (*) 1 $ init v
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| An internal representation of an array.
|||
||| Each array internally stores its values based on one of these
||| representations.
public export public export
data Rep : Type where data Rep : Type where
||| Byte-buffer array representation. Bytes : Order -> Rep
||| Boxed : Order -> Rep
||| This representations stores elements by converting them into byte values
||| storing them in a `Buffer`. Use of this representation is only valid if
||| the element type implements `ByteRep`.
|||
||| @ o The order to store array elements in
Bytes : (o : Order) -> Rep
||| Boxed array representation.
|||
||| This representation uses a boxed array type to store its elements.
|||
||| @ o The order to store array elements in
Boxed : (o : Order) -> Rep
||| Linked list array representation.
|||
||| This representation uses Idris's standard linked list types to store its
||| elements.
Linked : Rep Linked : Rep
||| Delayed/Lazy array representation.
|||
||| This representation delays the computation of the array's elements, which
||| may be useful in writing efficient algorithms.
Delayed : Rep Delayed : Rep
%name Rep rep %name Rep rep
||| An alias for the representation `Boxed COrder`, the C-like boxed array
||| representation.
|||
||| This representation is the default for all newly created arrays.
public export public export
B : Rep B : Rep
B = Boxed COrder B = Boxed COrder
||| An alias for the representation `Linked COrder`.
public export public export
L : Rep L : Rep
L = Linked L = Linked
||| An alias for the representation `Delayed`.
public export public export
D : Rep D : Rep
D = Delayed D = Delayed
@ -113,8 +80,6 @@ Eq Rep where
Delayed == Delayed = True Delayed == Delayed = True
_ == _ = False _ == _ = False
||| A predicate on representations for those that store their elements
||| linearly.
public export public export
data LinearRep : Rep -> Type where data LinearRep : Rep -> Type where
BytesIsL : LinearRep (Bytes o) BytesIsL : LinearRep (Bytes o)
@ -175,19 +140,11 @@ forceRepNCCorrect Delayed = DelayedNC
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| An interface for elements that can be converted into raw bytes.
|||
||| An implementation of this interface is required to use the `Bytes` array
||| representation.
public export public export
interface ByteRep a where interface ByteRep a where
||| The number of bytes used to store each value.
bytes : Nat bytes : Nat
||| Convert a value into a list of bytes.
toBytes : a -> Vect bytes Bits8 toBytes : a -> Vect bytes Bits8
||| Convert a list of bytes into a value.
fromBytes : Vect bytes Bits8 -> a fromBytes : Vect bytes Bits8 -> a
export export
@ -337,11 +294,6 @@ export
in fromBytes bs1 :: fromBytes bs2 in fromBytes bs1 :: fromBytes bs2
||| The constraint that each array representation requires.
|||
||| Currently, only the `Bytes` representation has a constraint, requiring an
||| implementation of `ByteRep`. All other representations can be used without
||| constraint.
public export public export
RepConstraint : Rep -> Type -> Type RepConstraint : Rep -> Type -> Type
RepConstraint (Bytes _) a = ByteRep a RepConstraint (Bytes _) a = ByteRep a

View file

@ -31,16 +31,16 @@ export
FieldCmp Double where FieldCmp Double where
abslt = (<) `on` abs abslt = (<) `on` abs
-- Alternative implementations of `Eq` and `FieldCmp` that compare approximately, -- Alternative implementations of `Eq` and `FieldCmp` that compare floating
-- useful when working with flating point numbers -- point numbers approximately, useful when working with transforms
namespace Eq namespace Eq
export export
WithEpsilon : (Neg a, Abs a, Ord a) => a -> Eq a WithEpsilon : Double -> Eq Double
WithEpsilon ep = MkEq (\x,y => abs (x - y) < ep) (\x,y => abs (x - y) >= ep) WithEpsilon ep = MkEq (\x,y => abs (x - y) < ep) (\x,y => abs (x - y) >= ep)
namespace FieldCmp namespace FieldCmp
export export
WithEpsilon : (Neg a, Fractional a, Abs a, Ord a) => a -> FieldCmp a WithEpsilon : Double -> FieldCmp Double
WithEpsilon ep = MkFieldCmp @{WithEpsilon ep} ((<) `on` abs) WithEpsilon ep = MkFieldCmp @{WithEpsilon ep} ((<) `on` abs)
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------

View file

@ -33,7 +33,7 @@ Matrix' n = Array [n,n]
export export
matrix : {default B rep : Rep} -> RepConstraint rep a => {m, n : _} -> matrix : {default B rep : Rep} -> RepConstraint rep a => {m, n : _} ->
Vect m (Vect n a) -> Matrix m n a Vect m (Vect n a) -> Matrix m n a
matrix x = array {rep, s=[m,n]} x matrix x = array' {rep} [m,n] x
||| Construct a matrix with a specific value along the diagonal. ||| Construct a matrix with a specific value along the diagonal.
@ -168,8 +168,7 @@ diagonal {n} mat with (viewShape mat)
export export
-- TODO: throw an actual proof in here to avoid the unsafety -- TODO: throw an actual proof in here to avoid the unsafety
minor : Fin (S m) -> Fin (S n) -> Matrix (S m) (S n) a -> Matrix m n a minor : Fin (S m) -> Fin (S n) -> Matrix (S m) (S n) a -> Matrix m n a
minor i j mat = replace {p = flip Array a} (believe_me $ Refl {x = ()}) minor i j mat = believe_me $ mat!!..[Filter (/=i), Filter (/=j)]
$ mat!!..[Filter (/=i), Filter (/=j)]
filterInd : Num a => (Nat -> Nat -> Bool) -> Matrix m n a -> Matrix m n a filterInd : Num a => (Nat -> Nat -> Bool) -> Matrix m n a -> Matrix m n a
@ -404,12 +403,12 @@ gaussStep i lu =
export export
decompLU : Field a => (mat : Matrix m n a) -> Maybe (DecompLU mat) decompLU : Field a => (mat : Matrix m n a) -> Maybe (DecompLU mat)
decompLU {m,n} mat with (viewShape mat) decompLU {m,n} mat with (viewShape mat)
_ | Shape [m,n] = map MkLU _ | Shape [m,n] = map (MkLU . convertRep _ @{getRepC mat})
$ iterateN (minimum m n) (\i => (>>= gaussStepMaybe i)) (Just mat) $ iterateN (minimum m n) (\i => (>>= gaussStepMaybe i)) (Just (convertRep Delayed mat))
where where
gaussStepMaybe : Fin (minimum m n) -> Matrix m n a -> Maybe (Matrix m n a) gaussStepMaybe : Fin (minimum m n) -> Matrix m n a -> Maybe (Matrix m n a)
gaussStepMaybe i mat = if mat!#[cast i,cast i] == 0 then Nothing gaussStepMaybe i mat = if mat!#[cast i,cast i] == 0 then Nothing
else Just $ gaussStep i mat else Just (gaussStep i mat)
||| The LUP decomposition of a matrix. ||| The LUP decomposition of a matrix.
@ -520,7 +519,7 @@ decompLUP {m,n} mat with (viewShape mat)
maxIndex x [] = x maxIndex x [] = x
maxIndex _ [x] = x maxIndex _ [x] = x
maxIndex x ((a,b)::(c,d)::xs) = maxIndex x ((a,b)::(c,d)::xs) =
if abslt b d then maxIndex x ((c,d)::xs) if abslt b d then assert_total $ maxIndex x ((c,d)::xs)
else assert_total $ maxIndex x ((a,b)::xs) else assert_total $ maxIndex x ((a,b)::xs)
gaussStepSwap : Fin (S $ minimum m n) -> DecompLUP mat -> DecompLUP mat gaussStepSwap : Fin (S $ minimum m n) -> DecompLUP mat -> DecompLUP mat
@ -568,8 +567,9 @@ solveLowerTri' {n} mat b with (viewShape b)
construct [] = [] construct [] = []
construct {i=S i} (b :: bs) = construct {i=S i} (b :: bs) =
let xs = construct bs let xs = construct bs
in (b - sum (zipWith (*) xs (reverse $ toVect $ replace {p = flip Array a} (believe_me $ Refl {x=()}) $ i' = assert_total $ case natToFin i n of Just i' => i'
mat !#.. [One i, EndBound i]))) / mat!#[i,i] :: xs in (b - sum (zipWith (*) xs (reverse $ toVect $ believe_me $
mat !!.. [One i', EndBound (weaken i')]))) / mat!#[i,i] :: xs
solveUpperTri' : Field a => Matrix' n a -> Vector n a -> Vector n a solveUpperTri' : Field a => Matrix' n a -> Vector n a -> Vector n a
@ -580,8 +580,9 @@ solveUpperTri' {n} mat b with (viewShape b)
construct _ [] = [] construct _ [] = []
construct i (b :: bs) = construct i (b :: bs) =
let xs = construct (S i) bs let xs = construct (S i) bs
in (b - sum (zipWith (*) xs (toVect $ replace {p = flip Array a} (believe_me $ Refl {x=()}) $ i' = assert_total $ case natToFin i n of Just i' => i'
mat !#.. [One i, StartBound (S i)]))) / mat!#[i,i] :: xs in (b - sum (zipWith (*) xs (toVect $ believe_me $
mat !!.. [One i', StartBound (FS i')]))) / mat!#[i,i] :: xs
||| Solve a linear equation, assuming the matrix is lower triangular. ||| Solve a linear equation, assuming the matrix is lower triangular.
@ -626,7 +627,6 @@ solve mat = solveWithLUP mat (decompLUP mat)
-------------------------------------------------------------------------------- --------------------------------------------------------------------------------
||| Determine whether a matrix has an inverse.
export export
invertible : FieldCmp a => Matrix' n a -> Bool invertible : FieldCmp a => Matrix' n a -> Bool
invertible {n} mat with (viewShape mat) invertible {n} mat with (viewShape mat)

View file

@ -3,6 +3,7 @@ module Data.NumIdr.PrimArray
import Data.Buffer import Data.Buffer
import Data.Vect import Data.Vect
import Data.Fin import Data.Fin
import Data.NP
import Data.NumIdr.Array.Rep import Data.NumIdr.Array.Rep
import Data.NumIdr.Array.Coords import Data.NumIdr.Array.Coords
import public Data.NumIdr.PrimArray.Bytes import public Data.NumIdr.PrimArray.Bytes
@ -97,7 +98,7 @@ indexUpdate {rep = Linked} is f arr = update is f arr
indexUpdate {rep = Delayed} is f arr = indexUpdate {rep = Delayed} is f arr =
\is' => \is' =>
let x = arr is' let x = arr is'
in if (==) @{eqCoords} is is' then f x else x in if eqNP is is' then f x else x
export export
indexRange : {rep,s : _} -> RepConstraint rep a => (rs : CoordsRange s) -> PrimArray rep s a indexRange : {rep,s : _} -> RepConstraint rep a => (rs : CoordsRange s) -> PrimArray rep s a

View file

@ -62,7 +62,7 @@ unsafeFromIns s ins = arrayAction s $ \arr =>
export export
create : {o : _} -> (s : Vect rk Nat) -> (Nat -> a) -> PrimArrayBoxed o s a create : {o : _} -> (s : Vect rk Nat) -> (Nat -> a) -> PrimArrayBoxed o s a
create s f = arrayAction s $ \arr => create s f = arrayAction s $ \arr =>
for_ (range 0 (product s)) $ \i => arrayDataSet i (f i) arr for_ [0..pred (product s)] $ \i => arrayDataSet i (f i) arr
||| Index into a primitive array. This function is unsafe, as it performs no ||| Index into a primitive array. This function is unsafe, as it performs no
||| boundary check on the index given. ||| boundary check on the index given.
@ -88,7 +88,7 @@ export
bytesToBoxed : {s : _} -> ByteRep a => PrimArrayBytes o s -> PrimArrayBoxed o s a bytesToBoxed : {s : _} -> ByteRep a => PrimArrayBytes o s -> PrimArrayBoxed o s a
bytesToBoxed (MkPABytes sts buf) = MkPABoxed sts $ unsafePerformIO $ do bytesToBoxed (MkPABytes sts buf) = MkPABoxed sts $ unsafePerformIO $ do
arr <- newArrayData (product s) (believe_me ()) arr <- newArrayData (product s) (believe_me ())
for_ (range 0 (product s)) $ \i => arrayDataSet i !(readBuffer i buf) arr for_ [0..pred (product s)] $ \i => arrayDataSet i !(readBuffer i buf) arr
pure arr pure arr
export export
@ -96,7 +96,7 @@ boxedToBytes : {s : _} -> ByteRep a => PrimArrayBoxed o s a -> PrimArrayBytes o
boxedToBytes @{br} (MkPABoxed sts arr) = MkPABytes sts $ unsafePerformIO $ do boxedToBytes @{br} (MkPABoxed sts arr) = MkPABytes sts $ unsafePerformIO $ do
Just buf <- newBuffer $ cast (product s * bytes @{br}) Just buf <- newBuffer $ cast (product s * bytes @{br})
| Nothing => die "Cannot create array" | Nothing => die "Cannot create array"
for_ (range 0 (product s)) $ \i => writeBuffer i !(arrayDataGet i arr) buf for_ [0..pred (product s)] $ \i => writeBuffer i !(arrayDataGet i arr) buf
pure buf pure buf
@ -115,7 +115,7 @@ foldl f z (MkPABoxed sts arr) =
if product s == 0 then z if product s == 0 then z
else unsafePerformIO $ do else unsafePerformIO $ do
ref <- newIORef z ref <- newIORef z
for_ (range 0 (product s)) $ \n => do for_ [0..pred $ product s] $ \n => do
x <- readIORef ref x <- readIORef ref
y <- arrayDataGet n arr y <- arrayDataGet n arr
writeIORef ref (f x y) writeIORef ref (f x y)
@ -127,7 +127,7 @@ foldr f z (MkPABoxed sts arr) =
if product s == 0 then z if product s == 0 then z
else unsafePerformIO $ do else unsafePerformIO $ do
ref <- newIORef z ref <- newIORef z
for_ (range 0 (product s)) $ \n => do for_ [pred $ product s..0] $ \n => do
x <- arrayDataGet n arr x <- arrayDataGet n arr
y <- readIORef ref y <- readIORef ref
writeIORef ref (f x y) writeIORef ref (f x y)
@ -137,3 +137,99 @@ export
traverse : {o,s : _} -> Applicative f => (a -> f b) -> PrimArrayBoxed o s a -> f (PrimArrayBoxed o s b) traverse : {o,s : _} -> Applicative f => (a -> f b) -> PrimArrayBoxed o s a -> f (PrimArrayBoxed o s b)
traverse f = map (Boxed.fromList _) . traverse f . foldr (::) [] traverse f = map (Boxed.fromList _) . traverse f . foldr (::) []
{-
export
updateAt : Nat -> (a -> a) -> PrimArrayBoxed o s a -> PrimArrayBoxed o s a
updateAt n f arr = if n >= length arr then arr else
unsafePerformIO $ do
let cpy = copy arr
x <- arrayDataGet n cpy.content
arrayDataSet n (f x) cpy.content
pure cpy
export
unsafeUpdateInPlace : Nat -> (a -> a) -> PrimArrayBoxed a -> PrimArrayBoxed a
unsafeUpdateInPlace n f arr = unsafePerformIO $ do
x <- arrayDataGet n arr.content
arrayDataSet n (f x) arr.content
pure arr
||| Convert a primitive array to a list.
export
toList : PrimArrayBoxed a -> List a
toList arr = iter (length arr) []
where
iter : Nat -> List a -> List a
iter Z acc = acc
iter (S n) acc = let el = index n arr
in iter n (el :: acc)
||| Construct a primitive array from a list.
export
fromList : List a -> PrimArrayBoxed a
fromList xs = create (length xs)
(\n => assert_total $ fromJust $ getAt n xs)
where
partial
fromJust : Maybe a -> a
fromJust (Just x) = x
export
unsafeZipWith : (a -> b -> c) -> PrimArrayBoxed a -> PrimArrayBoxed b -> PrimArrayBoxed c
unsafeZipWith f a b = create (length a) (\n => f (index n a) (index n b))
export
unsafeZipWith3 : (a -> b -> c -> d) ->
PrimArrayBoxed a -> PrimArrayBoxed b -> PrimArrayBoxed c -> PrimArrayBoxed d
unsafeZipWith3 f a b c = create (length a) (\n => f (index n a) (index n b) (index n c))
export
unzipWith : (a -> (b, c)) -> PrimArrayBoxed a -> (PrimArrayBoxed b, PrimArrayBoxed c)
unzipWith f arr = (map (fst . f) arr, map (snd . f) arr)
export
unzipWith3 : (a -> (b, c, d)) -> PrimArrayBoxed a -> (PrimArrayBoxed b, PrimArrayBoxed c, PrimArrayBoxed d)
unzipWith3 f arr = (map ((\(x,_,_) => x) . f) arr,
map ((\(_,y,_) => y) . f) arr,
map ((\(_,_,z) => z) . f) arr)
export
foldl : (b -> a -> b) -> b -> PrimArrayBoxed a -> b
foldl f z (MkPABoxed size arr) =
if size == 0 then z
else unsafePerformIO $ do
ref <- newIORef z
for_ [0..pred size] $ \n => do
x <- readIORef ref
y <- arrayDataGet n arr
writeIORef ref (f x y)
readIORef ref
export
foldr : (a -> b -> b) -> b -> PrimArrayBoxed a -> b
foldr f z (MkPABoxed size arr) =
if size == 0 then z
else unsafePerformIO $ do
ref <- newIORef z
for_ [pred size..0] $ \n => do
x <- arrayDataGet n arr
y <- readIORef ref
writeIORef ref (f x y)
readIORef ref
export
traverse : Applicative f => (a -> f b) -> PrimArrayBoxed a -> f (PrimArrayBoxed b)
traverse f = map fromList . traverse f . toList
||| Compares two primitive arrays for equal elements. This function assumes the
||| arrays have the same length; it must not be used in any other case.
export
unsafeEq : Eq a => PrimArrayBoxed a -> PrimArrayBoxed a -> Bool
unsafeEq a b = unsafePerformIO $
map (concat @{All}) $ for [0..pred (arraySize a)] $
\n => (==) <$> arrayDataGet n (content a) <*> arrayDataGet n (content b)
-}

View file

@ -48,7 +48,8 @@ bufferAction @{br} s action = fromBuffer s $ unsafePerformIO $ do
export export
constant : {o : _} -> ByteRep a => (s : Vect rk Nat) -> a -> PrimArrayBytes o s constant : {o : _} -> ByteRep a => (s : Vect rk Nat) -> a -> PrimArrayBytes o s
constant @{br} s x = bufferAction @{br} s $ \buf => do constant @{br} s x = bufferAction @{br} s $ \buf => do
for_ (range 0 (product s)) $ \i => writeBuffer i x buf let size = product s
for_ [0..pred size] $ \i => writeBuffer i x buf
export export
unsafeDoIns : ByteRep a => List (Nat, a) -> PrimArrayBytes o s -> IO () unsafeDoIns : ByteRep a => List (Nat, a) -> PrimArrayBytes o s -> IO ()
@ -62,7 +63,8 @@ unsafeFromIns @{br} s ins = bufferAction @{br} s $ \buf =>
export export
create : {o : _} -> ByteRep a => (s : Vect rk Nat) -> (Nat -> a) -> PrimArrayBytes o s create : {o : _} -> ByteRep a => (s : Vect rk Nat) -> (Nat -> a) -> PrimArrayBytes o s
create @{br} s f = bufferAction @{br} s $ \buf => do create @{br} s f = bufferAction @{br} s $ \buf => do
for_ (range 0 (product s)) $ \i => writeBuffer i (f i) buf let size = product s
for_ [0..pred size] $ \i => writeBuffer i (f i) buf
export export
index : ByteRep a => Nat -> PrimArrayBytes o s -> a index : ByteRep a => Nat -> PrimArrayBytes o s -> a
@ -100,7 +102,7 @@ foldl f z (MkPABytes sts buf) =
if product s == 0 then z if product s == 0 then z
else unsafePerformIO $ do else unsafePerformIO $ do
ref <- newIORef z ref <- newIORef z
for_ (range 0 (product s)) $ \n => do for_ [0..pred $ product s] $ \n => do
x <- readIORef ref x <- readIORef ref
y <- readBuffer n buf y <- readBuffer n buf
writeIORef ref (f x y) writeIORef ref (f x y)

View file

@ -2,9 +2,9 @@ module Data.NumIdr.PrimArray.Delayed
import Data.List import Data.List
import Data.Vect import Data.Vect
import Data.NP
import Data.NumIdr.Array.Rep import Data.NumIdr.Array.Rep
import Data.NumIdr.Array.Coords import Data.NumIdr.Array.Coords
import Data.NumIdr.PrimArray.Linked
%default total %default total
@ -22,8 +22,8 @@ export
indexRange : {s : _} -> (rs : CoordsRange s) -> PrimArrayDelayed s a -> PrimArrayDelayed (newShape rs) a indexRange : {s : _} -> (rs : CoordsRange s) -> PrimArrayDelayed s a -> PrimArrayDelayed (newShape rs) a
indexRange [] v = v indexRange [] v = v
indexRange (r :: rs) v with (cRangeToList r) indexRange (r :: rs) v with (cRangeToList r)
_ | Left i = indexRange rs (\is' => v (assertFin i :: is')) _ | Left i = indexRange rs (\is' => v (believe_me i :: is'))
_ | Right is = \(i::is') => indexRange rs (\is'' => v (assertFin (Vect.index i (fromList is)) :: is'')) is' _ | Right is = \(i::is') => indexRange rs (\is'' => v (believe_me (Vect.index i (fromList is)) :: is'')) is'
export export
indexSetRange : {s : _} -> (rs : CoordsRange s) -> PrimArrayDelayed (newShape rs) a indexSetRange : {s : _} -> (rs : CoordsRange s) -> PrimArrayDelayed (newShape rs) a

View file

@ -2,6 +2,7 @@ module Data.NumIdr.PrimArray.Linked
import Data.Vect import Data.Vect
import Data.DPair import Data.DPair
import Data.NP
import Data.NumIdr.Array.Rep import Data.NumIdr.Array.Rep
import Data.NumIdr.Array.Coords import Data.NumIdr.Array.Coords
@ -27,17 +28,15 @@ export
indexRange : {s : _} -> (rs : CoordsRange s) -> Vects s a -> Vects (newShape rs) a indexRange : {s : _} -> (rs : CoordsRange s) -> Vects s a -> Vects (newShape rs) a
indexRange [] v = v indexRange [] v = v
indexRange (r :: rs) v with (cRangeToList r) indexRange (r :: rs) v with (cRangeToList r)
_ | Left i = indexRange rs (Vect.index (assertFin i) v) _ | Left i = indexRange rs (Vect.index (believe_me i) v)
_ | Right is = assert_total $ _ | Right is = believe_me $ map (\i => indexRange rs (Vect.index (believe_me i) v)) is
case toVect _ (map (\i => indexRange rs (Vect.index (assertFin i) v)) is) of
Just v => v
export export
indexSetRange : {s : _} -> (rs : CoordsRange s) -> Vects (newShape rs) a -> Vects s a -> Vects s a indexSetRange : {s : _} -> (rs : CoordsRange s) -> Vects (newShape rs) a -> Vects s a -> Vects s a
indexSetRange {s=[]} [] rv _ = rv indexSetRange {s=[]} [] rv _ = rv
indexSetRange {s=_::_} (r :: rs) rv v with (cRangeToList r) indexSetRange {s=_::_} (r :: rs) rv v with (cRangeToList r)
_ | Left i = updateAt (assertFin i) (indexSetRange rs rv) v _ | Left i = updateAt (believe_me i) (indexSetRange rs rv) v
_ | Right is = foldl (\v,i => updateAt (assertFin i) (indexSetRange rs (Vect.index (assertFin i) rv)) v) v is _ | Right is = foldl (\v,i => updateAt (believe_me i) (indexSetRange rs (Vect.index (believe_me i) rv)) v) v is
export export

View file

@ -176,7 +176,7 @@ Traversable (Point n) where
export export
Show a => Show (Point n a) where Show a => Show (Point n a) where
showPrec d (MkPoint v) = showCon d "point" $ showArg $ toVect v showPrec d (MkPoint v) = showCon d "point" $ showArg $ elements v
export export
Cast a b => Cast (Point n a) (Point n b) where Cast a b => Cast (Point n a) (Point n b) where

View file

@ -25,13 +25,12 @@ isTranslation : Eq a => Num a => HMatrix' n a -> Bool
isTranslation {n} mat with (viewShape mat) isTranslation {n} mat with (viewShape mat)
_ | Shape [S n,S n] = isHMatrix mat && getMatrix mat == identity _ | Shape [S n,S n] = isHMatrix mat && getMatrix mat == identity
||| Try to construct a translation from a homogeneous matrix.
export
fromHMatrix : Eq a => Num a => HMatrix' n a -> Maybe (Translation n a)
fromHMatrix mat = if isTranslation mat then Just (unsafeMkTrans mat) else Nothing
||| Construct a translation given a vector. ||| Construct a translation given a vector.
export export
translate : Num a => Vector n a -> Translation n a translate : Num a => Vector n a -> Translation n a
translate v = unsafeMkTrans (translationH v) translate v = unsafeMkTrans (translationH v)
||| Try to construct a translation from a homogeneous matrix.
export
fromHMatrix : Eq a => Num a => HMatrix' n a -> Maybe (Translation n a)
fromHMatrix mat = if isTranslation mat then Just (unsafeMkTrans mat) else Nothing

View file

@ -40,8 +40,7 @@ vector v = rewrite sym (lengthCorrect v)
||| Convert a vector into a `Vect`. ||| Convert a vector into a `Vect`.
export export
toVect : Vector n a -> Vect n a toVect : Vector n a -> Vect n a
toVect v = replace {p = flip Vect a} (believe_me $ Refl {x=()}) $ toVect v = believe_me $ Vect.fromList $ Prelude.toList v
Vect.fromList $ Prelude.toList v
||| Return the `i`-th basis vector. ||| Return the `i`-th basis vector.
export export
@ -186,8 +185,8 @@ perp a b = a.x * b.y - a.y * b.x
||| Calculate the cross product of the two vectors. ||| Calculate the cross product of the two vectors.
export export
cross : Neg a => Vector 3 a -> Vector 3 a -> Vector 3 a cross : Neg a => Vector 3 a -> Vector 3 a -> Vector 3 a
cross v1 v2 = let [a, b, c] = toVect v1 cross v1 v2 = let [a, b, c] = elements v1
[x, y, z] = toVect v2 [x, y, z] = elements v2
in vector [b*z - c*y, c*x - a*z, a*y - b*x] in vector [b*z - c*y, c*x - a*z, a*y - b*x]
||| Calculate the triple product of the three vectors. ||| Calculate the triple product of the three vectors.