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.
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.
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.
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.
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 numeric scalars, vectors, or matrices.
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.
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.
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
.
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.
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 uses the Gibbs-Poole-Stockmeyer method to
reduce the profile of matrix m. "Profile" is defined as the sum
over all columns of the number of elements above the diagonal but not
above the highest non-zero element in the column. Profile reduction can
be very effective in reducing the time and memory required to factor a
matrix.
The matrix m must be symmetric in structure.
This function uses the GPSKCA subroutine written by John Lewis.
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
.
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.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.)
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")
.
The form
function differs from this function only in that it
pads with zeros rather than reusing elements of x.
See also form
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
.
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","","","")
.
The fill
function differs from this function only in that it
reuses the elements of x rather than padding with zeros.
See also fill
and zero
.
full
function converts a matrix stored in "sparse_upper"
form to the "sparse" form. Frankly, you shouldn't have to worry about
the storage form--the existence of this function shows that we haven't
met that goal yet.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
.
symmetric
.ident
function returns an identity matrix with n rows
and columns.
imax(v)
is not necessarily equal to imin(-v)
.
v must be (or be convertible to) a vector, and must not have complex type.
See also imin
, max
, min
, isort
, and sort
.
imin(v)
is not necessarily equal to imax(-v)
.
v must be (or be convertible to) a vector, and must not have complex type.
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
.
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 same as
lose(a[1];b),lose(a[2];b),...
.
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
.
Argument v must be (or be convertible to) a vector, 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.
Argument v must be (or be convertible to) a vector, and must not have complex type.
See also max
, imax
, imin
, isort
, and sort
.
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
.
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
.
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.
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. 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.
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
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.
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.
The key difference between the chol
and factor
functions
is that chol
requires that the matrix be positive definite. In
that case, no pivoting is required; this generally results in a
significant savings in execution time over the factor
function.
See also backsub
, factor
, and solve
.
eig
function computes eigenvalues and eigenvectors for the
standard and generalized eigenvalue problems
A * x = lambda * x
and
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 positive definite. The eigenvalues are real and in ascending order. The eigenvectors are b-orthonormal.
For real, nonsymmetric, generalized problems, DGEGV 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.
No capability for complex generalized problems currently exists.
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.
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 positive definite. In
that case, no pivoting is required; this generally results in a
significant 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 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))
If 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.
See also ifft
.
filter
function computes a digital filter from the standard
difference equation
y[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 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))
If 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.
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
.
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
.
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 a linear coordinate
transformation P'*A*P
, where A is a symmetric matrix. The
only real advantage to using transform
is that its result is
known by Algae to be symmetric.f(x)=0
. The ordinates A
and B must bracket a root, meaning that f(A)*f(B)
must be
negative.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.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.
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
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
.
Combinations that are disallowed in C (such as 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
function 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
.
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.
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.
put_ascii
. The file name can specify a pipe; see
fprintf
for more details.
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 Level-1 MAT-file format. Some day
MATLAB will have to define a Level-2 format that accommodates sparse
matrices. When they do, we'll have to update this 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.
This function is very crude, and will probably be removed soon. Only real, dense matrices are handled. Labels are kept, but symmetry is lost.
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 the Level-1 MAT-file format. Some day
MATLAB will have to define a Level-2 format that accommodates sparse
matrices. Until then, sparse matrices are written in dense form.
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 definition of
yow
by calling src
with the file name "yow"
. It
then calls its replacement (with the same arguments, of course).
The current implementation of autosrc
restricts the number of
arguments to 10. (It's an ugly restriction, but if you need that many
arguments perhaps you should be using a table, instead.) A zero
returned indicates success.
Suppose you have some old Fortran code that computes your chances of
being hit by a comet. You can write a C wrapper for it, link them as a
shared object, use builtin
to attach it, and, voila!, it's a
builtin function.
Unlike other function in Algae that deal with files, builtin
does
not accept the special "bang" type file names. (These are file names
that begin with an exclamation point and are given to the shell as a
command.)
It isn't hard to create shared objects for builtin
, but currently
there's no documentation for it. Perhaps you can find some code to
imitate.
cd
function changes the current default directory to that
specified by the character scalar path. If path is NULL,
the directory is changed to that given by the HOME
environment
variable.
An exception is raised if an error occurs; otherwise, cd
returns
the integer 1.
eval
function evaluates the string s as an expression
and returns its value. Variable references, if any, have global scope.
For example, consider the following interaction:
> x = 0.5; > str = "sin(x)"; > eval (str) 0.4794 > x = 1.0; > eval (str) 0.8415
Although the value of str
remains "sin(x)"
, the result of
eval(str)
changes when the value of x
changes.
Unlike its prominent place in some matrix programming languages,
eval
is rarely necessary in Algae. Most examples of its use are
rather contrived. Here's some code which evaluates an expression typed
in by the user:
x = eval (read ("!echo -n \"expression? \" 1>&2 ; head -1"))
One use for eval
that may not be quite so contrived involves
building a function within a function--something that you can't do
directly in Algae. Here's an example:
add_const = function (c) { return eval (sprintf ("function (x) { return x+%g; }"; c)); };
This could then be used as follows:
> ten_more = add_const( 10 ); > ten_more(1)? 11
Probably the most common mistake made with eval
is giving it a
statement instead of an expression. If that's really what you want to
do, use the exec
function instead.
Remember, too, that all variable references are at global scope. For example, the function
function( a; b ) { return eval( "a + b" ); }
returns the sum of the global variables a
and b
; it
completely ignores the function's arguments, as they are local
variables.
See also exec
.
exec
function causes algae
to parse and execute the
Algae statements contained in the character string s. Execution
begins at global scope, even if exec
is called from a function.
It returns a negative integer if an error occurs, and zero otherwise.
See also eval
and source
.
exit
should be numeric scalar. If
exit
is called with anything else, Algae returns 0 as its exit
status.
A call to exit
is not mandatory in your Algae programs. Algae
terminates when it reaches the end of its input, returning exit status
0. If an error causes Algae to terminate, exit status 1 is returned.
file
function returns the name of the file from which
statements are currently being executed.get_path
function returns a path list, for use in file
searching. The argument env, if it isn't NULL, gives the name of
an environment variable containing a colon-separated list of paths. The
def argument gives a default value, in case env isn't given
or doesn't exist.
See also search_path
and src
.
getenv
function returns a character scalar that contains the
value of the environment variable named by ename. If no such
variable exists, NULL is returned.
line
function returns the line number from which it is
called.provide
function with feature, the
name of that feature. This means that the facilities associated with
feature have been made available for other Algae programs.
Programs that require a particular feature can source the appropriate
file with the require
function.
The argument feature must be a character scalar. It is added to
the global variable $features
if it is not already present in
that vector. The provide
function returns feature.
require
with the name of
that feature. The use of named features simplifies the task of
determining whether required assignments have been made.
The require
function checks for the feature named feature
in the current Algae session; if it isn't present, then require
reads the source file filename with src
. If filename
is NULL
, then the value of feature is used as the file
name.
Before it finishes, the source file read by require
is expected
to call the provide
function, indicating that the desired feature
has been provided. An exception is raised if the source file cannot be
read or does not provide the feature named feature.
search_path
function searches for a file in the directories
named in path. The file name begins with the string given by
fn, with one of the suffices given by suffices.
algae
to parse and execute the file named
fname. Execution begins at global scope, even if source
is
called from a function. It returns NULL if it has a problem opening
fname, a negative integer if a parse or run time error occurs, and
zero otherwise.
The file is closed after it is read.
See also exec
and sources
.
sources
function causes algae
to parse and execute
specified files. The character scalar fname is passed through
sh(1) to ls(1), so file name generation is available and multiple names
may be given. For example,
sources( "*.A" );
sources all of the files ending with `.A' in your working
directory. Although you could use the source
function to do this
as
source( "!cat *.A" );
using sources
is preferable because it maintains the correct file
and line information for the interpreter.
The "pipe" capability normally available to Algae I/O functions is not
permitted in sources
. Command substitution works, though, as
demonstrated here:
sources( "`find $HOME -name '*.A' -print`" );
This example sources every file ending with `.A' found in any subdirectory under your home directory.
src
function finds and opens a file of Algae code, executes
the statements within it at global scope, and then closes the file.
This function, rather than source
or sources
, is intended
to be the normal way to source files directly. Other functions that may
be used to source files are require
and autosrc
.
To find the file, src
looks first for a file named
`filename.A', that is, it appends `.A' to the name given
by the argument filename. If a file with that name exists, it is
sourced; otherwise, src
looks for a file named filename
with nothing appended, and sources it if it exists.
If filename is a relative file name, such as `foo' or
`foo/bar.A',
src
searches for the file using the variable
$src_path
. It appends filename to each of the directories
listed in the character vector $src_path
(both with and without
the `.A'), and sources the first file it finds whose name matches.
The current default directory is tried only if it is specified in
$src_path
.
No directory searching is performed if filename is an absolute
file name. File names that begin with `/', `./', `../',
or `~/' are absolute. If the file name begins with `~/', the
tilde is replaced with the value of the
HOME
environment
variable.
A file name that begins with the `!' character is a pipe; the rest
of the file name is given as a command line to the shell, the standard
output of which is then read by
src
. For example,
src("!m4 foo.A")
runs the file `foo.A' through the m4
macro preprocessor
first before src
reads it.
vector is initialized from the environment variable
ALGAE_SRC_PATH
, if it exists. It may be modified in either of
the startup files `/usr/local/lib/algae/3.3.4/algae.A' or
`$HOME/.algae'. (See section Startup Files.)
The syntax of ALGAE_SRC_PATH
is the same as that of PATH
for
the shell; fields are separated by `:', and `.' is used for
the current default directory. To set ALGAE_SRC_PATH
in the Bourne
shell, type something like the following:
export ALGAE_SRC_PATH ALGAE_SRC_PATH=.:~/algae:/usr/local/lib/algae
Execution of the Algae code in filename begins at global scope.
For example, if the source file contains the single statement
x=1
, the assignment is made to the global variable x
regardless of whether src
was called from within a function or
not.
The src
function raises an exception if filename cannot be
found. It returns NULL
if it cannot be opened, a negative integer
if a parse or run-time error occurs, and 0 otherwise.
strip
function removes all file and line information from the
given function f. The parser ordinarily includes this information
so that, for example, run-time error messages can identify the line on
which they occurred. If an error occurs in a function that has been
stripped, then the error is attributed instead to the calling
expression.
This function may also be useful in conjunction with profiling. Any time spent in a call to a stripped function gets charged to the file and line from which it was called.
system("ls")
lists the files in the current directory. If the
call succeeds, the exit status is returned; otherwise a -1 is returned.
npsol
function uses a sequential quadratic programming method
to minimize a given objective function subject to given constraints.
The domain of the objective function is called the "design space", and
the coordinates of the design space are called "design
variables".
This function uses the NPSOL package of Fortran routines developed and
sold by Stanford University. Your version of Algae might not have the
npsol
builtin function; if not, it probably means that you don't
have a license for NPSOL or that there was some problem installing it on
your computer.
The objective argument may be either a function or a table. If a function, then it's assumed to take one argument, a vector of design variables, and return the value of the objective at that point. If objective is a table, it may have the following members:
objf
params
exists in objective, then it is passed as a second
argument to objf
.
objgrd
objf
.
params
objf
and objgrd
, providing a means for passing additional
parameters to those functions.
Only the objf
member is required.
The start argument is a vector that specifies the starting point.
The constraints argument is a table that may contain any of all of the following members:
side_constraints
linear_constraints
coefficients
bounds
side_constraints
, that gives the upper and lower
bounds for the linear constraints. If not given, defaults are provided
as with side_constraints
.nonlinear_constraints
values
conf
is a function returning the values and congrd
is a
function returning the gradients. The constraint gradients, if
provided, should be given as a matrix in which the rows correspond to
the constraints and the columns correspond to the design variables. As
in objectives, a member called params
may also be
included.
bounds
side_constraints
, that gives the upper and lower
bounds for the nonlinear constraints. If not given, defaults are
provided as with side_constraints
.{ hessian = "yes" }
specifies the
corresponding option in NPSOL. Use an underscore instead of blanks
between words, as in { major_print_level = 1 }
. The valid
options are:
central_difference_interval
cold_start
warm_start
cold_start
, the first working set is chosen by npsol
based on the values of the variables and constraints at the initial
point. Broadly speaking, the initial working set will include equality
constraints and bounds or inequality constraints that violate or
"nearly" satisfy their bounds within crash_tolerance
. With a
warm_start
, the user must provide the state
table. A warm
start will be advantageous if a good estimate of the initial working set
is available--for example, when npsol
is called repeatedly to
solve related problems.
crash_tolerance
cold_start
, which is the default. When making a cold start, the
QP algorithm in npsol
must select an initial working set. The
initial working set will include, if possible, bounds or general
inequality constraints that lie within crash_tolerance
of their
bounds. The default value is 0.01. If crash_tolerance
is less
than zero or greater than one, the default value is used.
derivative_level
0
1
objgrd
function in
argument objective. Constraint derivatives are not
provided.
2
congrd
function in argument constraints. Objective derivatives are not
provided.
3
objgrd
function in
argument objective, and the nonlinear constraint derivatives are
provided by the congrd
function in argument constraints.npsol
is more reliable and will usually be more efficient when
all derivatives are exact.
If derivative_level
is 0 or 2, npsol
will estimate the
objective gradients using finite differences. The computation of
finite-difference approximations usually increases the total run time,
since a call to the objective function is required for each constraint.
Furthermore, less accuracy can be attained in the solution.
If derivative_level
is 0 or 1, npsol
will approximate the
constraint gradients. One call to the constraint function is required
for each variable. At times, central differences are used rather than
forward differences, in which case twice as many calls to the objective
and constraint functions are needed. The switch to central differences
is not under the user's control.
difference_interval
npsol
determine constant elements in the objective and constraint
gradients.
feasibility_tolerance
function_precision
function_precision
should be set to
1.0E-6. In contrast, if the objective function is typically on the
order of 1.0E-4 and the first 6 digits are known to be correct, then
function_precision
should be set to 1.0E-10. The choice of
function_precision
can be quite complicated for badly scaled
problems. The default value is machine precision raised to the 0.9
power; this is appropriate for most simple functions that are computed
with full accuracy. However, when the accuracy of the computed function
values is known to be significantly worse than full precision, then
function_precision
should be large enough that npsol
will
not attempt to distinguish between function values that differ by less
than the error inherent in the calculation.
hessian
"yes"
or "no"
. (The default is
"no"
). This option controls the contents of the R
matrix
in state. The user should set hessian
to "yes"
if
the state
table returned is to be used for a subsequent warm
start.
infinite_bound_size
infinite_step_size
infinite_step_size
, the objective function is considered to be
unbounded below in the feasible region. The default is the greater of
infinite_bound_size
and 1.0E10.
iteration_limit
linear_feasibility_tolerance
nonlinear_feasibility_tolerance
npsol
, an iterative procedure is executed in order to
find a point that satisfies the linear constraints and bounds on the
variables to within linear_feasibility_tolerance
. All subsequent
iterates will satisfy the linear constraints to within the same
tolerance, unless this value is comparable to the finite difference
interval.
For nonlinear constraints, nonlinear_feasibility_tolerance
defines the largest constraint violation that is acceptable at an
optimal point. Since nonlinear constraints are generally not satisfied
until the final iterate, this value acts as a partial termination
criterion for the iterative sequence.
These tolerances should reflect the precision of the corresponding
constraints. for example, if the variables and the coefficients in the
linear constraints are of order unity, and the latter are correct to
about 6 decimal digits, it would be appropriate to specify
linear_feasibility_tolerance
as 1.0E-6.
linesearch_tolerance
list
print_level
npsol
. The levels are as follows:
0
1
5
10
20
minor_iteration_limit
minor_print_level
npsol
. The levels are as follows:
0
1
5
10
20
optimality_tolerance
optimality_tolerance
is 1.0E-6 and npsol
terminates successfully, the final value of
the objective function should have approximately 6 correct figures. The
npsol
function will terminate successfully if the iterative
sequence is judged to have converged and the final point satisfies the
first-order optimality conditions.
verify_level
objgrd
and congrd
. It may take on
one of the following values:
-1
0
1
2
3
npsol
. Its members are:
ISTATE
CLAMDA
R
npsol
function returns a table containing the following
members:
objective
solution
inform
0
optimality_tolerance
, i.e., the projected gradient and
active constraint residuals are negligible. (Success.)
1
2
3
4
major_iteration_limit
, has been reached.
6
7
9
iter
state
state
table described above.
plot
, splot
, replot
, and unplot
functions provide an interface to the gnuplot program for plotting
data. Use plot
for plotting lines in two or three dimensions;
splot
is for plotting surfaces in three dimensions.
To use these plotting functions effectively, you will probably need some
familiarity with the gnuplot commands. For example, to add a title to an
existing plot, you'd type something like replot("set title 'Good
Stuff'")
. Read the gnuplot manual or its on-line help for more
information.
The plot
function accepts as many as three arguments. If either
or both of the first two arguments are (or can be converted to)
character vectors, then their elements are sent, each on a separate
line, to gnuplot as commands. For example, plot("set term X11")
causes gnuplot to generate X11 output. Even commands that request
information from gnuplot, like plot("show term")
are acceptable
(although you might temporarily lose sight of Algae's prompt). Don't use
the "exit" or "quit" commands of gnuplot--they cause gnuplot to
exit without Algae knowing about it.
If more than one argument is given to plot
and the last one is an
integer scalar, then it is taken to be an identifier for the plot
process. This allows you to have more than one plot open at a time. If
no identifier is given, then the active plot is killed and replaced by a
new one with the same identifier. The "active" plot is the one last
referenced, or 1 if there are no plots open. The replot
function
may be used to modify an open plot.
Any other arguments given to plot
are taken to be data to be
plotted. Vectors are plotted using their element values as ordinates
and their labels for the abscissae. Matrices with 2 columns are plotted
using the second column for ordinates and the first column for
abscissae. A matrix with 3 columns is plotted as a curve in three
dimensions, with the third column specifying the ordinates. Surface
plots are obtained using the splot
function.
If the data given to plot
is a vector, but has no labels, then
the element indices are used instead of the labels.
Multiple curves may be shown on the same plot, either by giving
plot
two data arguments or (better yet) supplying the data in a
table. When data is given in a table (even if it contains only one
member), the member name is used in the plot legend. Otherwise,
plot
uses the names "y_1" and "y_2".
If the data is given in a table, any members that are tables themselves
are given to gnuplot in a single data file. This means that the same
line and point type is used. For example, if x
and y
are
vectors to be plotted, then plot({x;y})
plots the two, each
with a different line and point type and described by name in the
legend. On the other hand, plot({z={x;y}})
plots them with
the same line and point type and names the combination "z" in the
legend.
If plot
is called with no arguments, it returns a vector of the
open plot identifiers. This vector is sorted from most to least
recently referenced, so the "active" plot is first.
Here's an example that simulates the Van der Pol system and plots the results:
xdot = function( t; x ) { return x[1]*(1-x[2]^2)-x[2], x[1]; } x = ode( xdot; 0; 20; 0,.25 ); plot( "set data style lines"; { x1=x[1;]; x2=x[2;] } );
See also splot
, replot
, and unplot
.
replot
function is used to modify or redisplay plots created
with either the plot
or splot
functions. It takes at most
2 arguments. If the first argument has character type, it is passed to
gnuplot as commands.
The last argument is the optional plot identifier. If not given, the active plot is used by default.
See also plot
, splot
, and unplot
.
splot
function is used to plot surfaces in three dimensions.
Except for the form of the data, its input is identical to that of the
plot
function. The data is specifed as a matrix of ordinates.
The labels of the matrix, or the corresponding indices if the labels
don't exist, are used for the abscissae.
See also plot
, replot
, and unplot
.
umin
function performs unconstrained minimization using the
Nelder-Mead direct search method. It does not have the features or
sophistication of npsol
, but it works well for some
cases--particularly when the objective surface is not smooth.
The objective argument is a function that takes one argument, a vector of design variables, and returns the integer or real value of the objective at that point.
The start argument is an integer or real vector that specifies the starting point.
The options augument may be either NULL or a table containing convergence and other specifications. The meaningful options are:
bound
display
evals
iter
right
size
tol
verbose
The return value from umin
is a table with the following
members:
evals
inform
0
1
2
3
iter
msg
obj
sol
This code is based on nmsmax
, a MATLAB function by Nick
Higham.
See also npsol
.
unplot
function is used to terminate plots created with the
plot
or splot
functions. The integer scalar or vector
argument id specifies the identifiers of the plots to be
terminated. If no id is given, the active plot is terminated.
You can terminate all open plots with unplot(plot())
.
See also plot
, replot
, and splot
.
all
function evaluates the "truth" of its argument x
in the same way that the function test
does, except that vectors
and matrices are "true" only if all of their elements are "true".
For example, all(1,0)
returns 0 while test(1,0)
returns
1. If x has no elements, all
returns 0.
See also test
and equal
.
atof
function converts character strings to real numbers.
The argument s may be a scalar, vector, or matrix. The function
reads only up to the first unrecognized character of each string, and
ignores anything that remains. If the first character of a string is
unrecognized, then the value is taken to be zero.
char
function converts the vector v into a character
string; each element of v contributes a single character according
to its ASCII value. For example, char(65,66,67)
returns the
string "ABC"
.
If an element of v is less than 0 or greater than 255, then it
will "wrap" (modulo 256). Thus char(65)
and
char(65+256)
both return the string "A"
.
Algae's character strings are terminated with a NUL (0) character. For
the char
function, that means that if an element of v is
zero (or a multiple of 256) then the string is terminated at that
point. For example, char(65,0,66)
yields the string
"A"
.
class
function returns a character string (such as "scalar"
or "table") that describes the class of x.
a==b
would.
See also test
.
info("sort")
takes you directly to the description of the sort
function.
If possible, info
uses an X-based html browser (like netscape
or mosaic). Otherwise, a character-based html browser (like lynx)
will be tried. As a final resort, the GNU Info browser is called.
The names of the available browsers are assigned by the startup code
(see section Startup Files) as members xhtml
, html
, and
info
of the global table $programs
. To prevent
algae
from using a particular browser, simply set the
appropriate member of $programs
to a zero-length string. If
you prefer to use Info, for example, you could simply put the
line
$programs.xhtml = $programs.html = "";
in your `.algae' file.
prof
function reads an `algae.out' file produced by
algae
's execution profiler (the `-p' option), and summarizes
it by file and by line number. The infile argument specifies the
file name. outfile specifies the output file; if it's NULL,
stdout is used.
The threshold argument may be used to truncate the summary listings. The truncation occurs after enough entries have been printed to account for threshold percent of the total number of hits. If threshold is NULL, no truncation is performed.
Below is an example of prof
output, obtained with the command
prof("algae.out";;50);
.
Algae execution profile listing. 628 total hits --- by file --- hits % hits cum % file 201 32.01 32.01 lqrtest.A 127 20.22 52.23 /usr/local/lib/algae/ode.A --- by line --- hits % hits cum % line file 130 20.70 20.70 88 lqrtest.A 42 6.69 27.39 1 /usr/local/lib/algae/plot.A 40 6.37 33.76 38 lqrtest.A 21 3.34 37.10 1 /usr/local/lib/algae/ode.A 17 2.71 39.81 44 /usr/local/lib/algae/ode.A 15 2.39 42.20 42 /usr/local/lib/algae/ode.A 15 2.39 44.59 41 /usr/local/lib/algae/ode.A 14 2.23 46.82 1 /usr/local/lib/algae/ode4.A 12 1.91 48.73 1 /usr/local/lib/algae/spline.A 10 1.59 50.32 1 lqrtest.A
First, notice that the number of hits counted in this example was only 628. Four digits of precision in the statistics are clearly unjustifed with such a small sample. The first column reports the number of profiler hits counted toward each particular file or line number. The second column restates that as a percentage of the total number of hits. The third column gives the cumulative percentage.
All profiler hits that occur while algae
is parsing a file are
counted toward its first line. This is normally insignificant compared
to execution, but it explains why line 1 shows up so often in the "by
line" report in the example above.
No files are closed in prof
. Normally, this means that you must
close them yourself if prof
is to be called more than
once.
split
function takes a character scalar argument s,
splits it into tokens, and returns the tokens in a character vector.
Each token is delimited by one or more characters from w. For
example, split("/bin:/usr/bin";":")
returns the vector
( "/bin" , "/usr/bin" )
If w is NULL, the string " \t\n"
is used. Thus
split("This is a test.")
returns the vector
( "This" , "is" , "a" , "test." )
string
function converts its argument e to character
type. If e is NULL, a zero-length character scalar is returned.
For example, string(1/(1:3))
returns the vector
( "1" , "0.5" , "0.333333" )
Currently, the output format cannot be changed.
substr
function returns the substring of the character string
c, beginning at index start, with length length. The
integer start must be greater than zero; if it exceeds the length
of c, then a zero-length string is returned. If length is
NULL, then all of the remaining characters of c are returned;
otherwise, length must not be negative.
See also dice
, split
, and string
.
test
function evaluates the "truth" of e in the same
way that Algae's if
statement would, returning 1 for true and 0
for false. The following entities are false:
""
.
All others are true. Notice that if e is an array
(scalar
, vector
, or matrix
), then test
returns true if any element of e is nonzero (or has nonzero
length, for character type).
An example of the use of test
is Algae's equal
function,
which is written as
equal = function( a; b ) { return !test( a != b ); }
If a
and b
are both matrices, then a!=b
returns a
matrix that is all zeros only if every element pair is equal. In that
case, test
returns 0 and !
changes that to 1.
See also equal
and all
.
time
function returns the number of seconds of user time that
the current process has used. On most machines, its precision is not
much more than 1/10 of a second.tolower
function converts strings to lower case. Every
letter is converted as appropriate for the current locale. The argument
S may be a character scalar, vector, or matrix.
See also toupper
.
toupper
function converts strings to upper case. Every
letter is converted as appropriate for the current locale. The argument
S may be a character scalar, vector, or matrix.
See also tolower
.
what
function returns a table containing all of the global
functions. See also who
.
who
function returns a table containing all global variables
that are neither functions nor NULL. Variables with names that begin
with `$' are excluded, unless the optional argument opt
equals "$"
. See also what
.Go to the first, previous, next, last section, table of contents.