6. The Standard Functions | [ << ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
6. The Standard Functions
Functions are called by giving their name followed by a parenthesized
list of their arguments. For example, sqrt(a) computes the
square root of its argument `a'. The parentheses are required,
whether they contain any arguments or not. Functions may take more than
one number and kind of argument.
6.1 Basic Math abs, sin, sqrt, ... 6.2 Arrays fill, sparse, zero, ... 6.3 Sets union, intersection, ... 6.4 Linear Algebra solve, factor, fft, ... 6.5 Numerical Analysis ode4, brent, ... 6.6 Basic I/O read, printf, ... 6.7 Entity I/O get, put, ... 6.8 Execution source, system, ... 6.9 Special Tools plot, npsol, ... 6.10 Miscellaneous who, what, time, ...
abs function returns the absolute value of its numeric
argument. "Absolute value" is synonymous with "magnitude" for
complex types. If x is non-scalar, every element is replaced by
its absolute value.
x must be a numeric scalar, vector, or matrix.
acos function returns the arc cosine of its numeric
argument. If it's argument is complex or has a magnitude greater than
one, then the result is complex. Otherwise, the result is real. If
x is non-scalar, every element is replaced by its arc
cosine.The argument x must be a numeric scalar, vector, or matrix.
acosh function returns the hyperbolic arc cosine of its
numeric argument. If it's argument is complex or is less than 1.0, then
the result is complex. Otherwise, the result is real. If x is
non-scalar, every element is replaced by its hyperbolic arc cosine.
The argument x must be a numeric scalar, vector, or matrix. Note
that acosh(0) is not zero, so if acosh is applied to a
sparse vector or matrix, then a dense array is the result.
arg function returns the argument (phase angle) of its
numeric argument.
See also abs.
asin function returns the arc sine of its numeric argument.
If it's argument is complex or has a magnitude greater than one, then
the result is complex. Otherwise, the result is real. If x is
non-scalar, every element is replaced by its arc sine.The argument x must be a numeric scalar, vector, or matrix.
asinh function returns the hyperbolic arc sine of its numeric
argument. If x is non-scalar, every element is replaced by its
hyperbolic arc sine.The argument x must be a numeric scalar, vector, or matrix.
atan function returns the arc tangent of its numeric
argument. If it's argument is complex, then the result is complex.
Otherwise, the result is real. If x is non-scalar, every element
is replaced by its arc tangent.The argument x must be a numeric scalar, vector, or matrix.
atanh function returns the hyperbolic arc tangent of its
numeric argument. If it's argument is complex or has magnitude greater
than 1.0, then the result is complex. Otherwise, the result is real.
If x is non-scalar, every element is replaced by its hyperbolic
arc tangent.The argument x must be a numeric scalar, vector, or matrix.
atan2 function computes an angle corresponding to y and
x, the lengths of the opposite and adjacent sides. This is
similar to atan(y/x), but with these differences:
atan(y/x), because it avoids division by
zero.
The arguments must be scalars, vectors, or matrices, with either integer or real type.
ceil function returns the ceiling of its numeric argument.
The ceiling of x is the smallest integer that is not less than
x. If x is complex, the ceiling is applied to both real and
imaginary parts. If x is non-scalar, every element is replaced by
its ceiling. The type of x is not changed by ceil.x must be a numeric scalar, vector, or matrix.
conj function returns the complex conjugate of its numeric
argument. If x is non-scalar, every element is replaced by its
complex conjugate.x must be a numeric scalar, vector, or matrix.
cos function returns the cosine of its numeric argument. If
x is non-scalar, every element is replaced by its cosine.x must be a numeric scalar, vector, or matrix.
cosh function returns the hyperbolic cosine of its numeric
argument. If x is non-scalar, every element is replaced by its
hyperbolic cosine.x must be a numeric scalar, vector, or matrix.
erf function returns the error function of its numeric
argument. If x is non-scalar, every element is replaced by its
error function result.
The error function is defined as 2/sqrt(pi) times the integral
from 0 to x of exp(-t^2) dt. Among its useful
properties, the error function is twice the integral of the Gaussian
distribution with 0 mean and variance of 1/2.
Complex arguments are not yet implemented.
erfc function returns the complementary error function of
its numeric argument. If x is non-scalar, every element is
replaced by its complementary error function result.
The complementary error function erfc(x) is simply 1-erf(c),
but it is computed in a way that avoids round-off errors when x
is large.
Complex arguments are not yet implemented.
exp function returns the exponential of its numeric argument.
If x is non-scalar, every element is replaced by its
exponential.x must be a numeric scalar, vector, or matrix.
floor function returns the floor of its numeric argument.
The floor of x is the largest integer that is not greater than
x. If x is complex, the floor is applied to both real and
imaginary parts. If x is non-scalar, every element is replaced by
its floor. The type of x is not changed by floor.x must be a numeric scalar, vector, or matrix.
gcd(x).factors*x is equal
to gcd(x).
See also lcm and primef.
imag function returns the imaginary part of its numeric
argument. If x is non-scalar, every element is replaced by its
imaginary part. The value returned by imag has real type.
See also real.
integer function converts the real part of its numeric
argument to integer type, rounding it if necessary.
See also round.
See also gcd and primef.
log function returns the natural logarithm of its numeric
argument. If x is non-scalar, every element is replaced by its
logarithm.x must be a numeric scalar, vector, or matrix.
log10 function returns the base-10 logarithm of its numeric
argument. If x is non-scalar, every element is replaced by its
base-10 logarithm.x must be a numeric scalar, vector, or matrix.
See also gcd, lcm, and primes.
See also gcd, lcm, and primes.
real has real type.
See also imag.
round function rounds its numeric argument to the nearest
whole number. If x is complex, both real and imaginary parts are
rounded. If x is non-scalar, every element is rounded. The type
of x is not changed by round.If Algae has been compiled to use the system's rint(3m) function, then arguments with a fractional part of exactly 1/2 are rounded to the nearest even whole number. Otherwise, such arguments are rounded towards +infinity.
See also integer.
sin function returns the sine of its numeric argument. If
x is non-scalar, every element is replaced by its sine.sinh function returns the hyperbolic sine of its numeric
argument. If x is non-scalar, every element is replaced by its
hyperbolic sine.sqrt function returns the square root of its numeric
argument. If x is complex or negative, then the result is
complex. Otherwise, the result is real. If x is non-scalar,
every element is replaced by its square root.tan function returns the tangent of its numeric argument. If
x is non-scalar, every element is replaced by its tangent.tanh function returns the hyperbolic tangent of its numeric
argument. If x is non-scalar, every element is replaced by its
hyperbolic tangent.band function computes the upper and lower bandwidths and
profiles of the square matrix m. The row bandwidth of a
particular row is the number of elements to the left of the diagonal in
that row, but not including any consecutive zero elements at the
beginning of that row. Likewise, the column bandwidth of a particular
column is the number of elements above the diagonal in that column, but
not including any consecutive zero elements at the beginning of that
column. Notice that, for a diagonal matrix, all row and column
bandwidths are zero. The lower bandwidth is the largest of all the row
bandwidths, and the upper bandwidth is the largest of all the column
bandwidths. The lower profile is the sum of all the row bandwidths, and
the upper profile is the sum of all the column bandwidths.This function returns a real vector with four elements. In order, these elements are as follows: lower bandwidth, lower profile, upper bandwidth, and upper profile.
Profile reduction can be very effective in reducing the time and memory
required to factor a matrix. The function gpskca provides one
approach for doing this.
See also gpskca.
bdiag function mimics the diag function, but in a
"block" sense. The matrix x is taken to be comprised of
m rows and n columns of blocks that are themselves matrices.
If either m or n is 0, then the matrix returned has the
blocks of x on its diagonal and is elsewhere zero. Otherwise,
the blocks of the diagonal of x are appended together.An example of the former case is:
> M = magic (4)
[ 9 7 6 12 ]
[ 4 14 15 1 ]
[ 16 2 3 13 ]
[ 5 11 10 8 ]
> bdiag (M; 0; 2)
[ 9 7 . . ]
[ 4 14 . . ]
[ 16 2 . . ]
[ 5 11 . . ]
[ . . 6 12 ]
[ . . 15 1 ]
[ . . 3 13 ]
[ . . 10 8 ]
|
An example of the latter case is:
> bdiag ( M; 2; 2 )
[ 9 7 3 13 ]
[ 4 14 10 8 ]
|
x must be scalar, vector, or matrix. After being converted to a matrix, if necessary, the number of rows of x must be evenly divisible by m (unless m is 0) and the number of columns of x must be evenly divisible by n (unless n is 0).
See also btrans and diag.
btrans function mimics the transpose operator, but in a
"block" sense. The matrix x is taken to be comprised of m
rows and n columns of blocks that are themselves matrices. The
blocks themselves are not transposed, but are simply moved across the
diagonal.For example:
> M = magic (4)
[ 9 7 6 12 ]
[ 4 14 15 1 ]
[ 16 2 3 13 ]
[ 5 11 10 8 ]
> btrans (M; 2; 2)
[ 9 7 16 2 ]
[ 4 14 5 11 ]
[ 6 12 3 13 ]
[ 15 1 10 8 ]
|
x must be scalar, vector, or matrix. After being converted to a matrix, if necessary, the number of rows of x must be evenly divisible by m and the number of columns of x must be evenly divisible by n.
See also bdiag.
The circshift function shifts the elements of the array x.
The arg s is a vector--each of its elements specifies the shift
distance for the corresponding dimension of x. A positive shift is
down or right; a negative shift is up or left. If unspecified, the
shift distance is zero (meaning no shift). If x is a table, the
function is applied to each of its members.
For example:
> x = fill (3,3; 1:9)
[ 1 2 3 ]
[ 4 5 6 ]
[ 7 8 9 ]
> circshift (x; 1)
[ 7 8 9 ]
[ 1 2 3 ]
[ 4 5 6 ]
> circshift (x; 0,1)
[ 3 1 2 ]
[ 6 4 5 ]
[ 9 7 8 ]
> circshift (x; 1,1)
[ 9 7 8 ]
[ 3 1 2 ]
[ 6 4 5 ]
|
combine function appends the vectors u and v,
removing all but the first of any redundant elements. This is exactly
like the union function, except that the result is not sorted.shape[1] is returned. If shape has two elements,
then a matrix with shape[1] rows and shape[2] columns is
returned.
The elements of x are used, in order, to fill the new entity. For
matrices, this proceeds by rows. If the entity returned has more
elements than x, it is padded with zeros or null strings. For
example, cram( 5; 1,2 ) returns the vector
(1,2,0,0,0).
If the return value is an array, it may be either dense or sparse. The
code will try to make a reasonable choice, but avoids spending much time
in deciding. For example, cram( 1000; 1 ) returns a sparse
vector since it is immediately apparent that most of its elements are
zero. On the other hand, cram( 1000; (1:1000)<2 ) is dense even
though it would be more efficient to store it sparse; since its second
argument is a dense vector, cram would need to check most of its
elements to realize that.
The form function is identical to this function except that the
arrays it returns are always dense.
See also fill, form, and zero.
dense function converts its argument array to dense storage.
See also sparse.diag function performs two different tasks, depending on the
class of x. If x is a matrix, then its diagonal is returned
as a vector. If x is a vector, then a matrix is returned which
has x as its diagonal and is zero elsewhere. (If x is a
scalar, it is returned intact.)
See also bdiag.
dice function takes a character scalar s and returns a
character vector, each element of which is a single character from
s. For example, the expression dice("Go dawgs!")
yields
( "G", "o", " ", "d", "a", "w", "g", "s", "!" ) |
See also split.
diff function takes a numeric vector v and returns a
vector of differences between its elements. If v has n
elements, then the return vector is
v[2]-v[1], v[3]-v[2], ..., v[n]-v[n-1] |
The return vector has one less element than v. If v has labels, then the return vector contains all but the last one.
One simple use for diff is to compute forward-difference
approximations for the slope of a curve. If c is a vector
containing ordinates of the curve in its elements and abscissae in its
labels, then the expression
diff (c) / diff (unlabel (c.eid)) |
returns an approximation to its slope. (The call to unlabel is
probably unnecessary in most cases, but it avoids trouble when the
labels of c themselves have labels.)
mksparse function. A table
is returned, containing the members shape, rows, cols,
and values. The argument matrix x must be numeric. Despite
the name of the function, it need not be sparse.
See also mksparse.
shape[1] is returned. If shape has two elements,
then a matrix with shape[1] rows and shape[2] columns is
returned.
The elements of x are used, in order, to fill the new entity. For
matrices, this proceeds by rows. If the entity returned has more
elements than x, then the elements of x are reused. For
example, fill( 5; "a","b" ) returns the vector
("a","b","a","b","a").
When the return value from fill is a vector or matrix, it will
always be dense. See the cram function for sparse return arrays.
The form function differs from this function only in that it
pads with zeros rather than reusing elements of x.
See also form, cram, and zero.
find function locates the elements of b that have the
values given by a and returns a vector containing their element
numbers. For example, find( 2,3; 0,1,2,3,4 ) returns the vector
(3,4). One common use is in an expression like
A[find(77;A.rid);], which returns the row (or rows) of a
having the row label 77.
If a is a scalar, then find(a;b) returns the element
numbers of b, sorted from smallest to largest, for which the
corresponding element of b is equal to a. If a is a
vector, then find(a;b) returns the same as
find(a[1];b),find(a[2];b),....
See also grep and lose.
first function returns the index of the first "true"
element in the vector v. If none are found, 0 is returned.
See also find and last.
shape[1] is returned. If shape has two elements,
then a matrix with shape[1] rows and shape[2] columns is
returned.
The elements of x are used, in order, to fill the new entity. For
matrices, this proceeds by rows. If the entity returned has more
elements than x, it is padded with zeros or null strings. For
example, form( 5; "a","b" ) returns the vector
("a","b","","","").
When the return value from form is a vector or matrix, it will
always be dense. See the cram function for sparse return arrays.
The fill function differs from this function only in that it
reuses the elements of x rather than padding with zeros.
See also fill, cram, and zero.
full function converts a matrix stored in "sparse_upper"
form to the "sparse" form.gpskca function attempts to find a symmetric permutation of
the rows and columns of matrix m to reduce either its bandwidth or
its profile. If the flag argument is "true" (in the sense of
the test function, then the Gibbs-Poole-Stockmeyer algorithm is
used for bandwidth reduction; otherwise, the Gibbs-King algorithm is
used for profile reduction. See the description of the band
function for a definition of these terms.The return value is the permutation, given as an integer vector with length equal to the order of m. Consider the following example interaction:
> a = symmetric (magic(6) > 25);
> band (a)
( 4, 13, 4, 13 )
> v = gpskca (a);
> b = a[v;v];
> band (b)
( 2, 9, 2, 9 )
|
First, the sparse, symmetric matrix a is created. The
band function shows that it has a bandwidth of 4 and a profile of
13. Next, gpskca provides a permutation to reduce the profile.
The matrix b is set to equal a but with this permutation
applied. Finally, band shows that this permutation reduces both
the bandwidth (to 2) and the profile (to 9).
The matrix m must be either symmetric or hermitian. This function uses the GPSKCA subroutine written by John Lewis.
See also band.
grep function finds the elements of the vector v which
match the pattern expr. The pattern is an extended regular
expression which is given to the UNIX function egrep to do the
searching.The character strings in v must not contain a newline, or the results will generally be wrong.
Here are some examples (user input is preceded by the `>' prompt):
> grep ("a"; "ab", "ac", "bc")
( 1, 2 )
> grep ("9$"; 1:40)
( 9, 19, 29, 39 )
> m = magic (3)
[ 8 1 6 ]
[ 3 5 7 ]
[ 4 9 2 ]
> m.rid = "top", "middle", "bottom";
> m[ grep ("top|bottom"; m.rid); ]
[ 8 1 6 ]
[ 4 9 2 ]
|
This is a terribly inefficient implementation (maybe someday Algae will have builtin regular expressions), but maybe you'll find it useful.
See also find and lose.
If x has both row and column labels, they must match.
See also symmetric.
ident function returns an identity matrix with n rows
and columns.If v is a matrix, then a vector is returned, each element giving the row number of the greatest value in the corresponding column. Again, if more than one element qualifies, the first one is given.
If v is a table, then the imax function is applied to each
of its members.
The argument v may not be a scalar, and must not have complex
type. Note that imax(v) is not necessarily equal to
imin(-v).
See also imin, max, min, isort, and sort.
If v is a matrix, then a vector is returned, each element giving the row number of the least value in the corresponding column. Again, if more than one element qualifies, the last one is given.
If v is a table, then the imin function is applied to each
of its members.
The argument v may not be a scalar, and must not have complex
type. Note that imin(v) is not necessarily equal to
imax(-v).
See also imax, max, min, isort, and sort.
sort(20,40,10,30) returns the vector (3,1,4,2). The
expression v[isort(v)] actually returns the sorted vector
v, although the builtin function sort does this more
efficiently. Complex numbers are sorted primarily by real value and
secondarily by imaginary value.
See also sort, max, min, and set.
label function assigns labels to vectors and matrices. If
x is a vector, a is assigned as its element labels. If
x is a matrix, a and b are assigned as its row and
column labels, respectively. If x has any other class, an
exception is raised.
See also unlabel.
last function returns the index of the last "true"
element in the vector v. If none are found, 0 is returned.
See also find and first.
linspace function generates a vector of n elements,
spaced uniformly between a and b. All three arguments must
be (or be convertible to) scalars. The argument n must be (or
round to) an integer greater than one.
The vector generation operator `:' also generates uniformly spaced
vectors. With the operator, the number of elements is determined from
the specified spacing; with the linspace function, the spacing is
determined from the number of elements specified. Also notice that the
last element of the vector returned by linspace is equal to
b, which is not necessarily the case with the `:' operator.
See also logspace.
logspace function generates a vector of n elements,
spaced logarithmically between a and b. All three arguments must
be (or be convertible to) scalars. Both a and b must be
nonzero, and n must be (or round to) an integer greater than one.Logarithmic spacing means that the logarithm of the result vector is a uniformly spaced vector. This also means that the ratio between any two adjacent elements is constant. Note that even if both a and b are real, the result will be complex if they have opposite signs, and that a complex result will not necessarily form a straight line in the complex plane.
See also linspace.
lose function locates elements of b that do not have
the values given by a and returns a vector containing their
element numbers. For example, lose( 2,3; 0,1,2,3,4 ) returns the
vector (1,2,5). One common use is in an expression such as
A[lose(77;A.rid);], which returns all of the rows of a that
don't have the row label 77.
If a is a scalar, then lose(a;b) returns the element
numbers of b, sorted from smallest to largest, for which the
corresponding element of b is not equal to a. If a is
a vector, then lose(a;b) returns the intersection of the results of
lose(a[1];b), lose(a[2];b), etc.
See also find.
magic function returns a magic square of order n.
(What system would be complete without it?) The elements of a magic
square consist of all the integers 1 through n^2, arranged so
that the row, column, and the two diagonal sums are equal.n must be greater than 0. No magic square exists for order 2.
matrix
returns a real matrix with zero rows and zero columns.
See also scalar and vector.
If v is a table, then the max function is applied to each
of its members.
The argument v may not be a scalar, and must not have complex type.
See also min, imax, imin, isort, and sort.
This is demonstrated by the following interactive session:
> A = fill (3,3; 1:9)
[ 1 2 3 ]
[ 4 5 6 ]
[ 7 8 9 ]
> A.rid = A.cid = 1:3;
> B = fill (3,3; 10:90:10)
[ 10 20 30 ]
[ 40 50 60 ]
[ 70 80 90 ]
> B.rid = B.cid = 2:4;
> merge (A; B)
[ 1 2 3 0 ]
[ 4 15 26 30 ]
[ 7 48 59 60 ]
[ 0 70 80 90 ]
|
If x and y are not vectors or matrices, they are simply summed.
If v is a table, then the min function is applied to each
of its members.
The argument v may not be a scalar, and must not have complex type.
See also max, imax, imin, isort, and sort.
The following interactive session provides an example:
> mksparse ({shape=3,4; rows=1,2,3; cols=1,3,2; values=10,11,12})
[ 10 . . . ]
[ . . 11 . ]
[ . 12 . . ]
|
The values vector must have numeric type. If t is NULL, a real matrix with 0 rows and 0 columns is returned.
See also sparse and exsparse.
norm function computes the p-norm of x, where
p is 1, 2, "frobenius", or "infinity". (The latter two may be
abbreviated as "frob" and "inf".) If p is not specified, the
2-norm is used.For complex x, the 1 and "infinity" norms deal not with the magnitude of each element, but with the sum of the absolute values of the real and imaginary parts.
V [ pick (V < 0) ] |
returns all of the negative elements of V.
See also find.
product function computes the product of the elements of an
array. If x is a scalar, it is returned unchanged. If x is
a vector, the product of all of its elements is returned. If x is a
matrix, a vector is returned, each element of which is the product of the
elements in the corresponding columns of x. If x is a
table, the product function is applied to each of its members.readnum function for more information.
The random number generator may be "seeded" with the srand
function. If you don't call srand, the seed is based on the
system's clock. Note that rand and randn share the same
seed.
See also randn and srand.
readnum function for more information.
The random number generator may be "seeded" with the srand
function. If you don't call srand, the seed is based on the
system's clock. Note that rand and randn share the same
seed.
See also rand and srand.
reverse function simply reverses the order of the elements in
vector v. For example, reverse(1,2,3) returns the vector
(3,2,1).matrix and vector.select function selects one element at random from the given
vector. The result is a scalar. The vector x must have at least
one element.seq function returns a vector of consecutive integers from 1
to n. The argument n is first rounded to the nearest
integer; if it is less than 1, a zero-length vector is returned.sign function returns an array with the same size as
v. Each element of the returned array is either -1, 0, or +1,
depending on whether the corresponding element of v is negative,
zero, or positive, respectively. The argument v must be a numeric
scalar or array.sort(20,40,10,30) returns the vector
(10,20,30,40). Complex numbers are sorted primarily by real
value and secondarily by imaginary value.
See also isort, max, min, and set.
srand function is used to "seed" the random number
generator rand. (OK, they're really pseudo-random numbers.) The
seed s determines the sequence of numbers returned by rand.
If s is NULL (or if srand is never called at all), then the
seed is taken from the system's clock. See also rand.sum function sums elements of arrays. If x is a
scalar, it is returned unchanged. If x is a vector, the sum of
all of its elements is returned. If x is a matrix, a vector is
returned, each element of which is the sum of the elements in the
corresponding columns of x.surprise function returns a more-or-less random matrix with
characteristics given by the input arguments. It may be useful for
testing or timing purposes.
Each of the arguments to surprise may be a vector, in which case
one of its members is chosen at random for that characteristic. For
example, if rows is the vector 4,5,6 then the resulting
array will have either 4, 5, or 6 rows.
The rows and cols arguments specify the number of rows and columns in the array. The choice of dimensions has a lower priority than the other characteristics, so any dimensions that are incompatible with other choices are discarded. For example, consider the call
surprise (3,4; 4,5; "real"; ; "general","symmetric") |
This specifies a real matrix that is either general or symmetric. If general symmetry is chosen, then the array will have either 3 or 4 rows and either 4 or 5 columns. However, if a symmetric array is chosen, then it must be square and the result will have 4 rows and 4 columns.
On the other hand, the specified dimensions can affect the choice of the other characteristics. For example, in the call
surprise (3; 4; "real"; ; "general","symmetric") |
the result will not be symmetric, because the dimensions preclude it.
The type argument may contain any or all of the character strings
"integer", "real", or "complex". Only numeric
types are supported. The default is "real".
The density argument specifies the approximate ratio of nonzero to total array elements. For diagonal arrays, only the diagonal elements are considered. The default is 1.0.
The symmetry argument may be "general", "symmetric",
or "hermitian". The default is "general".
The other argument may be "none", "diagonal", or
"positive_definite". The default is "none".
See also hermitian.
tril function returns the lower "triangular" part of the
matrix a. The return matrix is identical to a except that,
for any i, all elements a[i;i+k+1], a[i;i+k+2],
etc. are zero. If k is NULL, it defaults to zero. Here is an
example:
> tril (magic (4); 1)
[ 9 7 0 0 ]
[ 4 14 15 0 ]
[ 16 2 3 13 ]
[ 5 11 10 8 ]
|
triu function returns the upper "triangular" part of the
matrix a. The return matrix is identical to a except that,
for any i, all elements a[i;i+k-1], a[i;i+k-2],
etc. are zero. If k is NULL, it defaults to zero. Here is an
example:
> triu (magic (4); 1)
[ 0 7 6 12 ]
[ 0 0 15 1 ]
[ 0 0 0 13 ]
[ 0 0 0 0 ]
|
unlabel function removes the labels from its argument. If
x is a vector or a matrix, it sets the labels to NULL and returns
the result. If x has any other class, it is returned
unchanged.
See also scalar and matrix.
readnum function for more
information.
If a vector or matrix is returned, it is stored in sparse form. This
may or may not be the most efficient storage for your application; use
the dense function to convert it if desired.
See also dense.
complement function returns the relative complement of
a in b; that is, the set of elements of b which are
not in a. The return value is a set, which we define as a sorted,
irredundant vector. "Sorted" means sorted by increasing value;
character strings are sorted by ASCII values and complex numbers
are sorted with the real values as the primary key and the imaginary
values as the secondary key.
See also intersection, set, and union.
intersection function returns the intersection of set a
in b; that is, the set of elements that are in both a and
b. The return value is a set, which we define as a sorted,
irredundant vector. "Sorted" means sorted by increasing value;
character strings are sorted by ASCII values and complex numbers
are sorted with the real values as the primary key and the imaginary
values as the secondary key.
See also complement, set, and union.
If x has labels, then the returned vector also has labels that
correspond to the retained elements. If x has no labels, the
labels of the returned vector give the index of the corresponding
element in x. If x has redundant elements, the element
retained is not necessarily the first. For example,
set(1,1,1).eid might be 1, 2, or 3.
See also complement, intersection, and union.
union function returns the union of sets a and b.
The return value is a set, which we define as a sorted, irredundant
vector. "Sorted" means sorted by increasing value; character strings
are sorted by ASCII values and complex numbers are sorted with the
real values as the primary key and the imaginary values as the secondary
key.
See also complement, intersection, and set.
backsub function solves the set of linear equations
A*x=b for x. (Despite the name, it does both forward and
backward substitution.) The argument afac must contain the
factorization of A returned as a table from either of the
functions factor or chol. The argument b provides
the right-hand side.This function handles both LAPACK (dense) and SuperLU (sparse) factors. It determines the appropriate solution method based on the members in the afac table.
If BCSLIB-EXT (a sparse math library produced by Boeing Computer
Services) is available on your system and Algae has been compiled
for it, then the backsub function can use it. If afac
contains factors computed by BCSLIB-EXT routines (in factor
or chol), then backsub will call the corresponding
routines to perform the back substitution.
See also chol, factor, and solve.
chol function computes the Cholesky factorization of the
symmetric, positive definite matrix m. The factors are returned
in a table, which may contain a variety of members depending on the type
and density of m. For example, if m is a real, dense
matrix, then "L" (the output of the LAPACK routine
DSYTRF routine) and "RCOND" are returned. For now, you'll
have to go digging in the code if you really want to understand the
output. The idea, though, is that functions like backsub can
make sense of them so that you don't have to.
Although its other members vary, the table returned by chol
always contains a member named "RCOND". This non-negative real
scalar is an estimate of the reciprocal of the condition number of
m. If "RCOND" is very small (that is, the condition number
is very large), then the matrix m is ill-conditioned. A warning
message is produced if "RCOND" is less than 1.0E-8.
The key difference between the chol and factor functions
is that chol requires that the matrix be symmetric and positive
definite. In that case, no pivoting is required; this generally results
in a savings in execution time over the factor function.
Unless the matrix m happens to be diagonal (or BCSLIB-EXT is
used), chol uses a dense method from LAPACK. Consequently,
factor is preferred over chol for large, sparse
systems.
If BCSLIB-EXT (a sparse math library produced by Boeing Computer
Services) is available on your system and Algae has been compiled
for it, then chol may use its routines to factor sparse matrices.
In turn, the backsub and solve functions can use these
factors to solve linear equations. It is often useful to compute the
factors of a matrix and then save them (using put) for later use.
Keep in mind, however, that if these factors were computed with the
BCSLIB-EXT routines, then the Algae version with which you
intend to use them must also support BCSLIB-EXT. The table
returned by chol contains the member HOLD iff a
BCSLIB-EXT routine was used to produce it. BCSLIB-EXT is
never used on dense matrices.
See also backsub, factor, and solve.
eig function computes eigenvalues and eigenvectors for the
standard and generalized eigenvalue problemsA * x = lambda * x |
A * x = lambda * B * x |
If a is the only matrix argument, then both lambda and x are computed for the standard problem. If a matrix b is given as the second argument, then the generalized problem is solved.
A table of options opts may be given as the last argument. If it
contains the member no_right, then the right eigenvectors x
are not computed. If it contains the member left, then the
corresponding left eigenvectors are computed.
This function returns a table containing the eigenvalues in a member
called values and, if appropriate, the (right) eigenvectors in a
member called vectors. If the left option was given and
the problem is unsymmetric, then the left eigenvectors are returned in a
member called left_vectors.
The solution method used depends on the type and symmetry of the coefficient arrays. Some of these methods entail strict requirements, such as positive definiteness, on the matrices. Each case is described below, along with the LAPACK routine called to solve it. See the LAPACK Users' Guide for more information.
For real, symmetric, standard problems, DSYEV is used. The eigenvalues are real and in ascending order; the eigenvectors are orthonormal.
For complex, hermitian, standard problems, ZHEEV is used. The eigenvalues are real and in ascending order; the eigenvectors are orthonormal.
For nonsymmetric, standard problems, DGEEV is used for real type and ZGEEV for complex type. The eigenvalues are complex, and complex conjugate pairs appear consecutively with the eigenvalue having the positive imaginary part first. The eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.
For real, symmetric, generalized problems, DSYGV is used. This method requires that a be symmetric and that b be symmetric and positive definite. The eigenvalues are real and in ascending order. The eigenvectors are b-orthonormal.
For real, nonsymmetric, generalized problems, DGGEV is used. This
method performs full balancing on a and b. The eigenvalues
are complex, and complex conjugate pairs appear consecutively with the
eigenvalue having the positive imaginary part first. Note, however,
that the eigenvalues are returned in a table; the member num
contains their numerators and member denom contains their
denominators. The quotients may easily overflow or underflow, and the
denominator (and maybe the numerator, too) may be zero. A good
beginning reference is the book Matrix Computations by G. Golub
and C. Van Loan. Each eigenvector is scaled so that the sum of the
absolute values of the real and imaginary parts of the largest component
is 1, except that for eigenvalues with numerator and denominator both
zero, the corresponding eigenvector is zero.
For complex, hermitian, generalized problems, ZHEGV is used. This method requires that a be hermitian and that b be hermitian and positive definite. The eigenvalues are real and in ascending order. The eigenvectors are b-orthonormal.
For complex, nonsymmetric, generalized problems, ZGGEV is used. The
eigenvalues and eigenvectors are complex. Note, however, that the
eigenvalues are returned in a table; the member num contains
their numerators and member denom contains their denominators.
The quotients may easily overflow or underflow, and the denominator (and
maybe the numerator, too) may be zero. A good beginning reference is
the book Matrix Computations by G. Golub and C. Van Loan. Each
eigenvector is scaled so that the sum of the absolute values of the real
and imaginary parts of the largest component is 1, except that for
eigenvalues with numerator and denominator both zero, the corresponding
eigenvector is zero.
See also iram.
equilibrate function computes row and column scalings
intended to equilibrate a matrix and reduce its condition
number. The input matrix A need not be square, but it must be
numeric and both dimensions must be nonzero. If A is hermitian or
real symmetric, it must be positive definite.
If equilibrate succeeds, it returns a table containing the
following members:
Even if the matrix A is symmetric or hermitian, the row and column scale factors may or may not be identical. If not, the resulting equilibrated matrix would no longer be symmetric or hermitian. You will have to decide whether the improvement in matrix condition is worth the loss of symmetry. The row and column scale factors will be identical for dense, real, symmetric matrices and for dense, complex, hermitian matrices (both of which must also be positive definite); otherwise, there's no guarantee.
Oftentimes, going on to apply the equilibration scale factors is not worth the effort. Roughly speaking, if amax is not close to overflow or underflow and both rowcnd and colcnd are greater than 0.1, it probably isn't worth doing.
If equilibrate is unsuccessful, it returns a scalar integer
instead of a table. If this value is positive, it indicates the row
causing the problem. If it is negative, its absolute value indicates
the column causing the problem.
To apply the equilibration, the matrix must be pre- and post-multiplied
by the scale factors. If r and c are identical, you'll
probably want to use the transform function to preserve symmetry.
Here is an example usage:
e = equilibrate (A);
if (equal (e.r; e.c))
{
Ae = transform (A; diag (e.r));
else
Ae = diag (e.r) * A * diag (e.c);
}
|
Notice, though, that this example does not check the result of
equilibrate to determine if it succeeded, it applies the scale
factors without checking to see if it's worth it, and it assumes that
A obeys the requirement regarding positive definiteness.
factor function computes a triangular factorization of the
matrix m. The factors are returned in a table, which may contain
a variety of members depending on the type, density, and definition of
m. For example, if m is a real, dense, symmetric matrix,
then "LD" and "IPIV" (the output of the LAPACK
routine DSYTRF routine) are returned. For now, you'll have to go
digging in the code if you really want to understand the output. The
idea, though, is that functions like backsub can make sense of them
so that you don't have to.For sparse systems, the SuperLU package is used. These routines reorder the columns of the matrix to increase its sparsity, and perform dynamic pivoting of the rows for numerical stability. (No advantage is taken of symmetry.) For the column ordering, the Multiple Minimum Degree (MMD) algorithm is used if the matrix is symmetric; otherwise, the Column Approximate Minimum Degree (COLAMD) algorithm is used. Currently, there is no provision to supply your own ordering.
Although its other members vary, the table returned by factor
always contains a member named "RCOND". This non-negative real
scalar is an estimate of the reciprocal of the condition number of
m. If "RCOND" is very small (that is, the condition number
is very large), then the matrix m is ill-conditioned. A warning
message is produced if "RCOND" is less than 1.0E-8.
If BCSLIB-EXT (a sparse math library produced by Boeing Computer
Services) is available on your system and Algae has been compiled
for it, then factor may use its routines to factor sparse matrices.
In turn, the backsub and solve functions can use these
factors to solve linear equations. It is often useful to compute the
factors of a matrix and then save them (using put) for later use.
Keep in mind, however, that if these factors were computed with the
BCSLIB-EXT routines, then the Algae version with which you
intend to use them must also support BCSLIB-EXT. The table
returned by factor contains the member HOLD iff a
BCSLIB-EXT routine was used to produce it. BCSLIB-EXT is
never used on dense matrices.
The key difference between the chol and factor functions
is that chol requires that the matrix be symmetric and positive
definite. In that case, no pivoting is required; this generally results
in a savings in execution time over the factor function.
See also backsub, chol, and solve.
fft function computes the complex discrete Fourier
transform of the vector (or matrix--see below) x using the fast
Fourier transform method. The result is ordered such that the
corresponding frequencies are increasing (from negative to positive).
Note that this is a different order than that returned by many FFT
routines.
If x has N elements, the transform is a complex vector with
the same length. If N is even, then element k of the
transform is equal to
sum (x @ exp (-2*i*pi*(k-N/2)*((1:N)-1)/N)) |
N is odd, then element k is equal to
sum (x @ exp (-2*i*pi*(k-(N+1)/2)*((1:N)-1)/N)) |
The method is most efficient when the number of elements in x is the product of small primes. For example, a vector with 1021 elements (a large prime number) takes roughly 100 times longer to transform than a vector with 1024 elements (a power of 2). On the other hand, 1152 elements (with prime factors of 2 and 3) is almost as fast as 1024.
If x has integer or real labels, they are taken to be the corresponding abscissae (i.e., time or distance) and are transformed accordingly to form labels for the results vector. In that case, the label values must be evenly spaced.
If x is a matrix, the transform is computed for each row or
column. By default this is done by column, but you can change it to
work by rows by supplying the row option, as described below.
The argument opts is an options table; it may contain any or all of the following members. (The values of these members are irrelevant; only their existence is signficant.)
row
estimate
measure
patient
measure, but considers a wider range of algorithms. This
may produce a "more optimal" plan, but at the expense of a longer
planning time.
exhaustive
patient, but with an even wider range of algorithms and
substantially increased planning time.
forget
forget option causes this accumulated
information to be discarded.
See also ifft.
filter function computes a digital filter from the standard
difference equationy[n] = b[1]*x[n] + b[2]*x[n-1] + b[3]*x[n-2] + ...
- a[2]*y[n-1] - a[3]*y[n-2] - ...
|
This filter is implemented using the "direct form II transposed" method. For more information, see Chapter 6 of the book "Discrete-Time Signal Processing" by Oppenheim and Schafer.
ifft function computes the complex inverse discrete Fourier
transform of the vector (or matrix--see below) x using the fast
Fourier transform method. The input is ordered such that the
corresponding frequencies are increasing (from negative to positive).
Note that this is a different order than that required by many inverse
FFT routines.
If x has N elements, the transform is a complex vector with
the same length. If N is even, then element n of the
transform is equal to
(1/N) * sum (x @ exp (2*i*pi*((1:N)-N/2)*(n-1)/N)) |
N is odd, then element n is equal to
(1/N) * sum (x @ exp (2*i*pi*((1:N)-(N+1)/2)*(n-1)/N)) |
The method is most efficient when the number of elements in x is the product of small primes. For example, a vector with 1021 elements (a prime number) takes roughly 100 times longer to transform than a vector with 1024 elements (a power of 2). On the other hand, 1152 elements (with prime factors of 2 and 3) is almost as fast as 1024.
If x has integer or real labels, they are taken to be the corresponding abscissae (i.e., frequencies) and are transformed accordingly to form labels for the results vector. In that case, the label values must be evenly spaced.
If x is a matrix, the transform is computed for each row or
column. By default this is done by column, but you can change it to
work by rows by supplying the row option.
For a description of the valid options and their effects, see the
fft function.
See also fft.
inv function computes the inverse of a matrix. This is
sometimes useful, but you should be aware that an alternative approach
using factor (or chol) and solve is usually more
efficient and accurate. (See below for more details.)
The inv function is implemented simply by calling the
solve function with an identity matrix as the second argument.
The argument opts is optional and is passed as the third argument
to solve. For example,
inv (A; {pos});
|
computes the inverse of the positive-definite matrix A.
The matrix inverse is often used to solve systems of linear equations.
If Ax=b is such a system, where A is a square matrix with
full rank, then its solution may be computed as
x = inv (A) * b; |
A more efficient and accurate approach, though, is to use the
solve function as
x = solve (A; b); |
The use of solve is appropriate even when several right-hand
sides are to be considered separately. For example, the code
x = {};
A_inv = inv (A);
for (b in members (blist))
{
x.(b) = A_inv * blist.(b);
}
|
would be better written as
x = {};
A_fac = factor (A);
for (b in members (blist))
{
x.(b) = solve (A_fac; blist.(b));
}
|
See also backsub, chol, factor, and solve.
The iram function provides almost all of the capability of
ARPACK. However, it is a relatively low-level function and requires
very careful use. Because of it's design, many simple input errors
cannot be detected by the function and may well result in incorrect
results. I hope to soon see a high-level function written, perhaps
called arnoldi, that calls iram but trades off some of its
flexibility in exchange for a simplified and less error-prone
interface.
To use iram effectively, you will need more information regarding
the various solution modes than can be included here. Please consult
the ARPACK users manual.
The most general problem solved by iram is to compute a few
eigenvalues (lambda) and eigenvectors (x) for
A * x = M * x * lambda |
where A and B are real or complex square matrices and M is hermitian (symmetric when real). This is a standard eigenvalue problem if M is the identity matrix; otherwise it is a generalized problem.
The argument n specifies the dimension of the problem. This should be an integer scalar, or at least able to be cast to one.
Rather than providing A and B directly, iram requires
input functions that perform certain operations on those matrices. This
allows numerous solution approaches (called "modes"), which are
described in detail below. The two user-supplied functions are called
as
op_func (v; params) b_func (v; params) |
The operations required of the op_func function differ depending
on the solution mode. In the simplest case (mode 1), it is just the
matrix-vector multiplication A*v. The argument params is
passed unchanged from the params argument to iram---it is
commonly a table containing the coefficient arrays A and
B.
The requirements for the b_func function also depend on the
solution mode. It is not used at all for some modes, in which case it
may be NULL. As with op_func, the params argument is passed
unchanged from the iram input.
The options argument is a table; its members specify various problem and solution options:
mode
1
A*v, and b_func should be NULL.
This mode is mostly suited to finding a few of the eigenvalues with
largest magnitude. It may be used with any problem type.
2
inv(B)*A*v (conceptually, but you
should use factor and solve rather than inv), and
b_func should return B*v. This mode may not be used with
real, symmetric problems.
3
inv(A-sigma*B)*B*v
(concepturally, but you should use factor and solve rather
than inv), and b_func should return B*v. This mode
may be used with any problem type. Convergence is enhanced to
eigenvalues that are near the value of sigma. For this mode, the
which option refers to the inverted problem, so use "LM" to
refer to the eigenvalues closest to sigma.
4
inv(A-sigma*B)*A*v (concepturally,
but you should use factor and solve rather than
inv), and b_func should return A*v. This mode may
only be used for real, symmetric problems. Convergence is enhanced to
eigenvalues that are near the value of sigma. For this mode, the
which option refers to the inverted problem, so use "LM" to
refer to the eigenvalues closest to sigma.
5
inv(A-sigma*B)*(A+sigma*B)*v
(concepturally, but you should use factor and solve rather
than inv), and b_func should return B*v. This mode
may only be used for real, symmetric problems. Convergence is enhanced
to eigenvalues that are near the value of sigma. For this mode,
the which option refers to the inverted problem, so use "LM" to
refer to the eigenvalues closest to sigma.
bmat
"I" for standard problems or "G" for
generalized problems. The default is "I".
which
"LM"). For real, symmetric problems, which must be one of
the following: "LA" (largest amplitude), "SA" (smallest
amplitude), "LM" (largest magnitude), "SM" (smallest
magnitude), or "BE" (both ends of the spectrum). For
non-symmetric or complex problems, which must be one of these:
"LM" (largest magnitude), "SM" (smallest magnitude), or
"LR" (largest real), "SA" (smallest real), "LI"
(largest imaginary), or "SI" (smallest imaginary). The default
is "LA" for symmetric problems and "LR" for the
others.
nev
n/2. For non-symmetric and complex problems, nev must be
less than n-1.
ncv
2*nev and
n. For real, non-symmetric problems, ncv must be less than
or equal to n, must exceed nev+1, and defaults to the
lesser of 2*nev+1 and n. For complex problems, ncv
must be less than or equal to n, must exceed nev, and
defaults to the lesser of 2*nev and n.
mxiter
basis
test function), then iram
returns, along with its other results, an orthonormal basis for the
associated approximate invariant subspace.
vectors
test function), then iram
returns the eigenvectors along with its other results.
tol
resid
shift
An underdetermined system is one for which A has fewer rows than columns. The rank of A must be equal to the number of rows of A. The return value is the minimum norm solution.
An overdetermined system is one for which A has more rows than columns. This function treats a square A as if it were an overdetermined system. The rank of A must be equal to the number of columns of A. The return value is the solution to the least squares problem
minimize || B - A*X || |
See also solve.
Only the leftmost non-NULL arguments are used, and at most 10 arguments are allowed. They will be converted to matrices if necessary (and possible), but in a context-free manner. A scalar becomes a 1x1 matrix, and a vector becomes a matrix with 1 row. This is in contrast to Algae's multiplication operator, which makes different conversions according to the context.
As an example, consider an NxN matrix A and an N vector
b. If you give Algae the product A'*A*b, it first
multiplies A'*A which involves an order of N^3 operations. But
if you do it as A'*(A*b), there are only an order of N^2
operations. The difference can be huge! With this function, you'd do
it as mult(A';A;b'). (Note that we have to explicitly convert
b to a matrix with the right shape.)
A*x=b for x.
Conceptually it returns inv(A)*b, except that solve is
generally much more efficient and accurate. If A has already been
factored by the factor or chol functions, then its factors
(in the form of a table returned by those functions) may be given in the
place of A. If the optional argument opts contains the
member "pos", then the Cholesky method is used to factor
A.
See also backsub, chol, and factor.
ssi function uses subspace iteration to compute the first few
eigenvalues and vectors for a real, symmetric, generalized eigenvalue
problem A*v=B*v*w. The terms w and v are the
eigenvalues and eigenvectors, respectively.Subspace iteration is an effective approach for large, sparse problems and for problems for which good estimates of the eigenvectors are available. Otherwise, it sucks.
Both A and B must be real, symmetric, non-singular matrices. The matrix V contains a set of starting vectors, one to a column. A convergence tolerance may be specified with eps; 1.0E-6 is used if eps is NULL. A table is returned, containing the members "values" and "vectors".
See also eig.
A=U*S*V',
where S is a real matrix with the same dimensions as A and
U and V are the orthogonal (unitary) left and right singular
vectors. The matrix S contains the singular values on its
diagonal and is zero elsewhere.
When A is not square, one or more rows or columns of S are
known a priori to be zero. By default, svd then returns a
"compact" form of the singular vectors, in which the corresponding
columns of U or V are not included. This behavior may be
modified by using the full option, as described below.
The argument opts is an options table; it may contain the members
novectors and full. (The values of these members are
irrelevant; only their existence is signficant.) If opts contains
the member novectors, then only the singular values are computed.
If opts contains only the member full, then the full
singular vectors are computed.
The results of svd are returned in a table. The member
sigma contains the vector of singular values, sorted in
decreasing order. Unless the novectors option is specified, the
members U and V contain the corresponding left and right
singular vectors.
This function calls the LAPACK routines DGESVD and ZGESVD.
transform function performs the linear coordinate
transformation P'*A*P, where A is a square matrix. In
fact, there are only two differences between the operator expression
P'*A*P and the function expression transform(A;P):
transform function results in the most appropriate symmetry.
For example if A and P are both real, and A is
symmetric, then the result is a symmetric matrix. With the operator
expression this would not be the case; the result would be exactly the
same (neglecting any roundoff errors) except that it would be marked
with general symmetry.
transform function is considerably faster than the operator
expression for some matrix combinations. In particular, it may be two
or three times faster when P is dense and A is large,
sparse, and either hermitian or real symmetric.f(x)=0. The ordinates A
and B must bracket a root, meaning that f(A)*f(B) must be
negative. If param is given, it is passed as the second argument
to the function f.
The tol argument specifies an error tolerance. (This argument
is required: there is no default.) There is no limit to the number of
points chosen; monte continues to pick points until its one
standard deviation error estimate is less than tol. Generally,
the accuracy increases only as the square root of the number of points
evaluated.
The Monte Carlo method provides a simple, easy approach for solving (numerically) even complicated, multi-dimensional integrals. However, it is notoriously inefficient. Each additional digit of accuracy requires 100 times more effort.
Here is a simple example, in which we compute the area of a circle with unit radius. With domain, we specify a square centered on the origin within which to integrate. The number of independent variables in the problem (2 in this case) is determined by the number of rows of domain. The function f takes as its argument a vector of the independent variables and returns the function value at that point. It will be called repeatedly, with the independent variables spread randomly within the domain. For this problem, f returns 1 if the point is inside the circle, and 0 otherwise.
domain = [-1,1; -1,1];
f = function (x) { return norm(x) < 1; };
monte (f; domain; 0.1)?
3.040
monte (f; domain; 0.01)?
3.120
monte (f; domain; 0.001)?
3.142
|
As you can see from the results, the third call to monte gave
a more accurate result, but it took 10,000 times the effort used in
the first call.
By the way, this problem could be solved faster with a one-dimensional integration. Then again, we already know the answer, don't we?
ode4 function uses fourth and fifth order Runge-Kutta
formulas to solve an initial value problem involving a system of
ordinary differential equations. The function f, called as
f( t; x ), is expected to return the first derivative of the
state vector x with respect to the independent variable t.
The arguments t0 and tf are the initial and final values of
the independent variable, and x0 is the initial state vector. The
argument y may be a function, called as y( t; x ), that
returns a vector of output values. If y is NULL, the output
vector is the state vector itself.The optional argument tol specifies the desired accuracy; we use 0.001 if tol is NULL. An initial stepsize may be suggested with h. The s argument may be used to provide better control over the accuracy--if s is not NULL, then a step is accepted if the estimated error in each state, divided by the corresponding term of s, is not greater than tol.
The return value is a matrix containing the output history. Each column contains the output from y for the particular value of the independent variable given by the column label.
n+1
elements, the equation solved (for x) is
c[1]*x^n + c[2]*x^(n-1) + ... + c[n]*x + c[n+1] = 0 |
A vector with zero elements is returned if the equation has no solution
(roots(1), for example). NULL is returned if it has an infinite
number of solutions (roots(0), for example). Otherwise, the
vector returned contains all of the solutions.
If spline is called with two arguments, a and x, then
the spline a (or a spline interpolation of a, if a is
not a spline) is evaluated at the abscissa values given by the vector
x.
As an example, use trapz to compute the definite integral of
1/x from 1 to 2.
x = linspace (1; 2; 11);
y = label (1/x; x);
trapz (y)
0.6938
log(2)-log(1)
0.6931
|
close function closes the specified file. The file name
file must be given exactly as it was when the file was opened.
(Whitespace is significant.) Any buffered data is flushed.
digits function may be used to specify the number of digits
that Algae uses when printing real and complex numbers. This
applies to the print function and to "printing statements"
(those with question mark or newline terminators). The argument n
specifies the number of digits. If n is NULL, no action is
taken.The function returns the number of digits being used. By default, Algae uses 4 digits.
Here's an example interaction, with user input preceded by the ">" prompt:
> digits () # the default
4
> acos (-1) # pi, to 4 significant digits
3.142
> digits (20); # set to 20 digits
> acos (-1) # pi, to 20 digits (but only accurate to 17)
3.1415926535897931000
|
This function simply sets the global variable $digits, so you
can do the same thing by assigning to it directly.
See also print.
fread function reads data from a binary file into a numeric
array. This is a low-level routine, and it requires that you know
precisely the format of the data that is being read. To read Algae's own
binary files (those written with put), see get.The file argument (a character scalar) is the file name. As with any file name, if its first character is an exclamation mark then the rest of the string is given to the shell as a filter. If file is NULL, then standard input is read.
The length argument (an integer scalar) gives the number of items
to read. If it's NULL, the entire file is read. Otherwise, it should
be a non-negative integer. A vector is returned that contains the data
that was read. If fread reaches the end of the file before
having read length elements, then the vector it returns will
contain fewer than length elements.
The type argument is an optional table; its member's names specify the type of data to read and, implicitly, the type of the returned array. (The values of this table's members are inconsequential; only their names are important.) These names may include any of the C language types:
char
integer.
int
integer.
float
real.
double
real.
The C language type qualifiers are also accepted in type:
long
int and double types. Both are cast
to real.short
int. Casts to integer.
signed
char and int. Casts to integer.
unsigned
char and int. An unsigned char
casts to integer, while an unsigned int casts to
real.unsigned double)
are disallowed here, too. If type is NULL or empty, double
is assumed.Besides the C language type qualifiers just described, the following special qualifiers are also accepted:
big
little
ieee
fprintf function prints formatted output to file. The
format argument is a character string that contains two types of
objects: characters and conversion specifications. Characters are
printed directly; conversion specifications cause the next argument to
be formatted and printed.Conversion specifications are introduced using the percent sign `%'. Although extensions are planned, the conversion specifications are currently identical to those for the "fprintf" function of the ANSI C language.
If the first character of file is an exclamation mark, then the
rest of file is given to the shell as a filter to which the
printed output of fprintf is sent. For example, the command
fprintf("!sort";"a\nc\nb\n");close("!sort"); sends "a", "c",
and "b", each on a separate line, to the system's sort function, which
then presumably sorts and prints it. If file is NULL, stdout is
used.
See also message, print, printf, sprintf,
and put.
fwrite function writes binary data to a file from a numeric
array. This is a low-level routine, and it requires that you specify
precisely the format of the data being written. To write an entity for
use in Algae, you'll probably want to use put instead.The file argument (a character scalar) is the file name. As with any file name, if its first character is an exclamation mark then the rest of the string is given to the shell as a filter. If file is NULL, then the data is written to standard output.
The data to be written is given in the argument v, a numeric vector.
The type argument is an optional table; its member's names specify
the type of data to write. (The values of this table's members are
inconsequential; only their names are important.) These names may
include any of the C language types and qualifiers, and a few others.
See fread for details.
message function is used to write error messages to stderr.
Its arguments are exactly like printf, except that a maximum of
10 arguments is allowed. The message written includes the program and
file names as well as the line number.
See also printf.
print function prints its argument x to a file named by
file. This output goes through the same formatting that
Algae uses to print entities to the screen. In particular, the
variable $term_width controls the width of this output, and the
$digits variable specifies the number of significant digits
printed for real and complex numbers.If file is NULL, then the output is written to stdout.
See also digits, fprintf, printf, sprintf,
and put.
printf function prints formatted output to stdout. It is
exactly the same as calling fprintf with NULL as the first
argument.
See also fprintf, message, print, put, and
sprintf.
fprintf
for details. If file is NULL, the text is read from standard
input.
When read reaches the end of file, it returns NULL.
See also get and readnum.
shape[1] elements is returned. If shape
has two elements, then a matrix is returned with shape[1] rows and
shape[2] columns.
The values are read from file, or from stdin if file is NULL
or not given. The `#' symbol is special--it and all remaining
characters on a line are ignored. If the end of file is reached before
the specified number of values is read, then the result is filled out
with zeros. After a call to readnum, the number of values actually
read is contained in the variable $read.
For example, let's say you have a file called `jnk' that contains the following:
# this line is ignored - numbers like 1.23 are skipped the first two numbers are 42 and 99. 3.14 is the next one |
In the following interaction we read the numbers from the file into a matrix:
> x = readnum (2,4; "jnk")?
[ 42.00 99.00 3.140 0.000 ]
[ 0.000 0.000 0.000 0.000 ]
> $read?
3
> readnum (3; "jnk")?
( 0.000, 0.000, 0.000 )
|
As you can see, anything that doesn't look like a number is simply
skipped. The readnum function found three numbers and filled in
the rest of the matrix with zeros. The global variable $read was
set to 3 - the number of values actually read.
Subsequent calls to readnum with the same file name simply
continue reading from the last position in the file. In the last
statement of the previous example, no more numbers were found so the
vector was filled with zeros. To read again from the start of the file,
you must first close the file with the close function.
The readnum function recognizes integers and floating point
numbers. A floating point number consists of an optionally signed
integer part, a decimal point, a fraction part, an `e' or `E',
and an optionally signed integer exponent. All of the following
examples are recognized as numbers: `1', `-1.2', `1e2',
`-1.e2', `.1e2', `-1.2e3', and `1.2e-3'.
Fortran users take note! A Fortran double precision number such as
`1.23D4' will not be read by readnum as a single number.
The character `D' is not recognized as part of a floating point
number, so readline will read this as the two numbers `1.23'
and `4'.
sprintf function writes formatted output to a character
scalar. It works the same as printf, except that the output is
returned instead of printed. For example, x=sprintf("y=%d";1)
puts the string "y=1" in the scalar x.
See also fprintf, printf, and put.
tmp_file function creates a temporary file and returns its
name. This file will be deleted automatically when Algae
exits.
By default, Algae writes temporary files in the directory
`/tmp'. You can override that with an environment variable called
TMPDIR. For example, if you use the Bourne shell and put the
lines
TMPDIR=/usr/tmp export TMPDIR |
in your `~/.profile' file, then Algae will put its temporary files in the `/usr/tmp' directory.
This function not only creates the file, but opens it for output. If you want to read from the temporary file, then you should close it first. (Of course, you'll need to arrange to get something in it, first.) Notice that temporary files are not deleted when closed, but only when Algae exits.
get function reads an Algae entity from a binary file named
file. If no argument is given, standard input is used. Use this
function to read files written with put. The file name can
specify a pipe; see fprintf for more details.
NULL is returned if an error occurs in reading the file; otherwise
get returns the entity it read.
The file name can specify a pipe; see fprintf for more
details.
getmat function reads matrices from a MATLAB binary file
named fname. It returns the matrices in a table.
The file name can specify a pipe; see fprintf for more
details.
Notice that getmat doesn't assign the matrices to the global
symbol table as MATLAB's "load" function does. To do that, you might
use something like the following:
Load = function (fname)
{
local (t; n);
t = getmat (fname);
for (n in members (t))
{
$$.(n) = t.(n);
}
};
|
The getmat function reads the MATLAB v.4 "Level-1 MAT-file"
format. MATLAB can read and write in this format, but its default for
writing is now Level-2. Using the newer format would require us to link
to their libraries, meaning that you'd need to actually have MATLAB on
that machine. We'll stick with the older format.
One unfortunate consequence of using the older file format is that
matrices cannot be written in sparse form. To get around this, convert
the matrix to coordinate form (which is dense) within MATLAB and then
change it back to sparse form in Algae with the mksparse
function.
See also get, put, and putmat.
load function reads a table from a file named fname and
assigns its members as global variables. This can be used, for example,
to restore an Algae session that was saved with the save
function.
If an error occurs in reading the file, load returns 0; otherwise
it returns 1.
See also get, put, and save.
put function writes x to a binary file named
fname (or standard output if fname is NULL). This binary
file can be read with get. The file name can specify a pipe; see
fprintf for more details.NULL is returned if an error occurs in writing the file; otherwise a 1 is returned.
fprintf for more details.
The file is not closed when putdyn finishes, so subsequent writes
to the same file will append data to the end of it unless you close it
first. On some machines, this could be used to build a complete DYNASTY
file with several different calls to putdyn. The disadvantage,
though, is that the DYNASTY file format doesn't allow this on all
machines.
putmat function writes a table of matrices t to a
MATLAB binary file named fname. The file name can specify a pipe;
see fprintf for more details.NULL is returned if an error occurs in writing the file; otherwise a 1 is returned.
NULL members in t are skipped; otherwise each member is converted to a matrix if necessary and then written to the file.
The putmat function writes in the MATLAB v.4 "Level-1 MAT-file"
format. MATLAB can read and write in this format, but its default for
writing is now Level-2. Using the newer format would require us to link
to their libraries, meaning that you'd need to actually have MATLAB on
that machine. We'll stick with the older format.
One unfortunate consequence of using the older file format is that
matrices cannot be written in sparse form. To get around this, use the
exsparse function to convert the matrix to coordinate form (which
is dense) and then change it back to sparse form within MATLAB.
See also get, getmat, and put.
save function writes all non-NULL global variables as a table
to a file named fname. This can be used, for example, to save an
Algae session that you wish to start again later using the
load function.
NULL is returned if an error occurs in writing the file; otherwise
save returns 1.
See also get, load, and put.
src function with the file name
file. In this way, you can put off reading a function's
definition until it is actually called.
If file is NULL, then autosrc uses the value of
name for the file name.
For example, consider the following code:
autosrc ("yow");
yow ();
|
The first line defines an autosrc function named yow. When
yow is called on the second line, it reads the real d