1.1.3.1. Feature Compatibility

In principle, the functions defined in oct-files must have the same capabilities as functions in m-files. In particular, the input and output arguments lists must be able to carry different numbers of arguments, which are moreover of different type, just as m-file functions can. This means that there must be a way of mapping a function from Octave like

function [array, real-scalar, integer] =
    func(complex-scalar, array, list, integer) 
    ## func frobs the snafu, returning all gromniz coefficients.

    -- actual function code goes here

endfunction

to C++. To this end, Octave's core interface supplies the macro (in OCT/src/defun-int.h)

octave_value_list
DEFUN_DLD(function_name,
          const octave_value_list& input_argument_list,
          [ int number_of_output_arguments ],            -- optional argument
          const char* documentation_string)

We have decorated the macro arguments with their types. Note that the first argument is the name of the function to be defined.

Of course the macro has to be defined somewhere. The easiest way to pull in all necessary definitions is to include OCT/src/oct.h. Now, our skeletal source file of the dynamic extension has the following shape:

#include <oct.h>

DEFUN_DLD(func, args, nargout,
          "func frobs the snafu, returning the gromniz coefficients.")
{
    -- actual function code goes here

    return possibly_empty_list_of_values;
}

Often functions do not depend on the actual number of output parameters. Then, the third argument to DEFUN_DLD can be omitted.

	Naming convention
	DEFUN_DLD gives the user the freedom to choose any name for input_argument_list and number_of_output_arguments, but conventionally args and nargout are used (thus reminding of the parameter names in main, which are argc and argv).

1.1.3.2. Essential types and methods

Before we start to write a low-level implementation of matpow, we look at the most essential types and methods used to handle data that flows through the Octave-interpreter interface.

As has been said above, the arguments passed to dynamically loaded functions are bundled in an octave_value_list. Results returned from a function are also passed in an octave_value_list. The default constructor, octave_value_list, creates an empty list which, when used as return value, yields an Octave function returning "void". To access the elements of a list the following methods form class octave_value_list (defined in OCT/src/oct-obj.h) are useful:

{
public octave_value_list();public octave_value& operator()(int index);  public const octave_value operator()(int index);  public const int length();}

length returns the number of elements in the list. The overloaded parenthesis operator selects a single element from the list. The index base for index is 0 and not 1 as Fortran adepts might infer. The value returned by operator() is an octave_value, i.e., any type known to Octave, e.g. integers, real matrices, complex matrices, and so on.

Knowing how to access the arguments of unspecified type, the next thing we would like to do is get their values. Octave follows a uniform naming scheme: all functions that return an object of a certain type ends in _type. Some of the more important of these methods (defined in OCT/src/ov-base.h) are:

{
public const int int_value();  public const double double_value();  public const Complex complex_value();  public const Matrix matrix_value();  public const ComplexMatrix complex_matrix_value();}

1.1.3.3. Low-level implementation

Now, we are ready to implement matpow as a dynamically loaded extension to Octave.

#include <octave/oct.h>


static bool is_even(int n);


DEFUN_DLD(matpow, args, ,
          "Return b = a^n for square matrix a, and non-negative, integral n.")
{
    const int n = args(1).int_value();  

    if (n == 0)
        return octave_value(
               DiagMatrix(args(0).rows(), args(0).columns(), 1.0) 
               );
    if (n == 1)
        return args(0);

    Matrix p(args(0).matrix_value());   
    Matrix b(p);
    int np = n - 1;
    while (np >= 1)
    {
        if (is_even(np))
        {
            np = np / 2;
        }
        else
        {
            np = (np - 1) / 2;
            b = b * p;                  
        }
        p = p * p;
    }

    return octave_value(b);            
}


bool
is_even(int n)
{
    return n % 2 == 0;
}
    

We learn from the example that the C++ closely resembles the Octave function. This is due to the clever class structure of the Octave library.

1.1.3.4. Compiling

Now that we are done with the coding, we are ready to compile and link. The Octave distribution contains the script mkoctfile, which does exactly this for us. In the simplest case it is called with the C++-source as its only argument.

$ ls
Makefile  matpow0.m  matpow1.m  matpow.cc  test_matpow.m
$ mkoctfile matpow.cc

Good, but not good enough! Presumably, we shall compile several times, so we would like to run our test suite and finally remove all files that can be re-generated from source. Enter: Makefile.

# Makefile for Octave extension matpow

.phony: all
all: matpow.oct

.phony: clean
clean:
        rm -f matpow.oct *.o

.phony: distclean
distclean: clean
        rm -f *~ octave-core core

.phony: test
test: matpow.oct
        octave --silent test_matpow.m

%.oct: %.cc
        mkoctfile $<
    

1.1.3.5. Running

If the oct-file is in the LOADPATH, it will be loaded automatically – either when requesting help on the function or when invoking it directly.

$ ls
Makefile  matpow0.m  matpow1.m  matpow.cc  matpow.o  matpow.oct  test_matpow.m
$ octave -q
octave:1> help matpow
matpow is the dynamically-linked function from the file
/home/cspiel/hsc/octave/src/matpow/matpow.oct

Return b = a^n for square matrix a, and non-negative, integral n.

octave:2> matpow([2, 0; 0, -1], 4)
ans =

  16   0
   0   1

octave:3> Ctrl-D

1.1.4.1. Input parameter checking

Why do we need parameter checking at all? A rule of thumb is to perform a complete parameter checking on

• every function that is exposed to the user (e.g. at Octave's interpreter level), i.e. for every function created with DEFUN_DLD, or

• any function that might be used in a different context than the writer intended Functions that are used only internally do not require full parameter checking.

One problem which arises frequently with parameter checking of (interfacing) functions is that the checks easily take up more space than the routine itself, thereby distracting the reader from the function's original purpose. Often, including all the checks into the main function bloats it beyond the desired maintainability. Therefore, we have to take precautions against these problems by consistently factoring out all parameter checks.

The rule of thumb here is to group logically related tests together in separate functions. The testing functions get the original functions' arguments by constant reference and return a boolean value. Constant referencing avoids any constructor calls and – in addition to that – protects the arguments against modification. In our example a single test function is enough:

static bool any_bad_argument(const octave_value_list& args);

any_bad_argument prints a detailed message, raises an Octave error, and returns true if the arguments fail any test; otherwise it returns false and we can continue processing. The only thing we have to change in matpow is to add a call to any_bad_argument and on failure to return an Octave-void value, i.e. an empty octave_value_list. The first few lines of matpow then take the following form:

DEFUN_DLD(matpow, args, ,
          "Return b = a^n for square matrix a, and non-negative, integral n.")
{
    if (any_bad_argument(args))
        return octave_value_list();

    const int n = args(1).int_value();

    -- rest of matpow is unchanged
}

As we are convinced that we have to check the input arguments, the question is how to do the checks.

	Ordering of Tests
	The correct ordering of tests is essential! Actual number of arguments Type of argument Size of argument (where applicable) Range of argument (where applicable and necessary) Inter-argument relations (if necessary) Item 1 always goes first. Items 2-4 usually repeat for each of the arguments. The final tests check the relations between the arguments, i.e., belong to Item 5.

Octave supplies the user with all necessary type- and size-testing functions. The type-tests (defined in OCT/src/ov-base.h) share the common prefix is_. Here are the most commonly used:

{
public const bool is_real_scalar();  public const bool is_complex_scalar();  public const bool is_real_matrix();  public const bool is_complex_matrix();}

To examine the sizes of different Octave objects, the following methods prove useful:

{
public const int length();}

{
public const int rows();  public const int columns();  public const int length();}

Remember that the methods available depend on the underlying type. For example, a ColumnVector only has a length (OCT/src/ov-list.h), whereas a Matrix has a number of rows and columns (OCT/src/ov-base-mat.h).

We have all the knowledge we need to write the argument testing function to augment matpow.

static bool
any_bad_argument(const octave_value_list& args)
{
    if (!args(0).is_real_matrix())
    {
        error("matpow: expecting base A (arg 1) to be a real matrix");
        return true;
    }

    if (args(0).rows() != args(0).columns())
    {
        error("matpow: expecting base A (arg 1) to be a square matrix");
        return true;
    }

    if (!args(1).is_real_scalar())
    {
        error("matpow: exponent N (arg 2) must be a real scalar.");
        return true;
    }

    if (args(1).scalar_value() < 0.0)
    {
        error("matpow: exponent N (arg 2) must be non-negative.");
        return true;
    }

    if (floor(args(1).scalar_value()) != args(1).scalar_value())
    {
        error("matpow: exponent N (arg 2) must be an integer.");
        return true;
    }

    return false;
}

Out final duty is to update the test-frame and run it. For brevity, we only list the new tests in testmp.m

###
### The following tests should trigger the error exits.
###

## number of arguments

unit_test_err("error exit, too few arguments",
              "matpow:", "matpow([1,1; 1,1])");
                                                                                                                            
unit_test_err("error exit, too few arguments",
              "matpow:", "matpow()");
                                                                                                                            
unit_test_err("error exit, too many arguments",
              "matpow:", "matpow([1,1; 1 1], 2, 1)");

## argument type and size
                                                                                                                           
unit_test_err("error exit, A not a matrix",
              "matpow:", "matpow(1, 1)");
                                                                                                                            
unit_test_err("error exit, A not a square matrix",
              "matpow:", "matpow([1 2 3; 4 5 6], 1)");
                                                                                                                            
unit_test_err("error exit, N not a real scalar(here: complex)",
              "matpow:", "matpow(a, 1+i)");
                                                                                                                            
unit_test_err("error exit, N not a real scalar(here: non-scalar)",
              "matpow:", "matpow(a, ones(2,2))");
                                                                                                                            
unit_test_err("error exit, N negative",
              "matpow:", "matpow(a, -1)");
                                                                                                                            
unit_test_err("error exit, N non-integral",
              "matpow:", "matpow(2.5)");

$ octave -q test_matpow.m
................
# of testcases attempted   16
# of expected passes       16
# of expected failures     0
# of unexpected passes     0
# of unexpected failures   0
# of unresolved testcases  0
 
Unit testing completed successfully!

1.1.4.2. Texinfo documentation strings

Our much improved matpow still carries around the poor documentation string:

"Return b = a^n for square matrix a, and non-negative, integral n."

Let us improve on this one! Octave supports documentation strings – docstrings for short – in Texinfo format. The effect on the online documentation will be small, but the appearance of printed documentation will be greatly improved.

The fundamental building block of Texinfo documentation strings is the Texinfo-macro @deftypefn, which takes two arguments: The class the function is in, and the function's signature. For functions defined with DEFUN_DLD, the class is Loadable Function. A skeletal Texinfo docstring therefore looks like this:

-*- texinfo -*- 
@deftypefn{Loadable Function} {return_values =} function_name(arguments)

@cindex index term 

Help text in Texinfo format...
@end deftypefn
    

Texinfo has several macros which control the markup. In the group of format changing commands, we note @var{variable_name}, and @code{code_piece}. The Texinfo format has been designed to generate output for online viewing with text-terminals as well as generating high-quality printed output. To these ends, Texinfo has commands which control the diversion of parts of the document into a particular output processor. Two formats are of importance: info and TeX. The first is selected with

@ifinfo
text for info only
@end ifinfo

the latter with

@iftex
@tex
text for TeX only
@end tex
@end iftex

If no specific output processor is chosen, by default, the test goes into both (or, more precisely: every) output processors. Usually, neither @ifinfo, nor @iftex appear alone, but always in pairs, as the same information must be conveyed in every output format.

-*- texinfo -*-
@deftypefn{Loadable Function} {@var{b} =} matpow(@var{a}, @var{n})

@cindex matrix power

Return matrix @var{a} raised to the @var{n}-th power.  Matrix @var{a} is
a square matrix and @var{n} a non-negative integer.
@iftex 
@tex 
$n = 0$
@end tex
@end iftex
@ifinfo 
n = 0
@end ifinfo
is explicitly allowed, returning a unit-matrix of the same size as
@var{a}.

Note: @code{matpow}  duplicates part of the functionality of the built-in
exponentiation operator
@iftex
@tex
``$\wedge$''.
@end tex
@end iftex
@ifinfo
`@code{^}'.
@end ifinfo


Example:
@example 
matpow([1, 2; 3, 4], 4)
@end example
@noindent
returns
@example
ans =

  199  290
  435  634
@end example

The algorithm takes
@iftex
@tex
$\lfloor \log_{2}(n) \rfloor$
@end tex
@end iftex
@ifinfo
floor(log2(n))
@end ifinfo
matrix multiplications to compute the result.
@end deftypefn
    

For further information about Texinfo consult the Texinfo documentation. For TeX-beginners we recommend "The Not So Short Introduction to LaTeX" by Tobias Oetiker et. al.

One thing we held back is the true appearance of a Texinfo docstring – mainly because it looks so ugly. The C++-language imposes the constraint that the docstring must be a string-constant. Moreover, because DEFUN_DLD is a macro, every line-end has to be escaped with a backslash. The backslash does not insert any whitespace and TeX separates paragraphs with empty lines, so that we have to put in new-lines as line-separators. Thus, the Texinfo docstring in source form has each line end decorated with "\n\".

DEFUN_DLD(matpow, args, ,
          "-*- texinfo -*-\n\
@deftypefn{Loadable Function} {@var{b}} = matpow(@var{a}, @var{n})\n\
\n\
@cindex matrix power\n\
\n\
Return matrix @var{a} raised to the @var{n}-th power.  Matrix @var{a} is\n\
a square matrix, and @var{n} a non-negative integer.\n\
...")

At least the formatted versions look much better. The info-version, which will be used in Octave's online help has approximately the following appearance:

b = matpow(a, n)                                        Loadable Function
Return matrix a raised to the n-th power. Matrix a is a square matrix, and
n a non-negative integer. n = 0 is explicitely allowed, returning a
unit-matrix of the same size as a.

Note: matpow duplicates part of the functionality of the built-in
exponentiation operator `^'.

Example:
matpow([1, 2; 3, 4], 4)

returns
ans =

  199  290
  435  634

The algorithm takes floor(log2(n)) matrix multiplications to compute the
result.

whereas the TeX-version will look like this: