Making An Matrix Operator or Function

This howto provides step-by-step instructions on writing your own matrix operator or function for OpenMx. Check the list of matrix operators and functions to see if your function has already been implemented. If not, you will be writing a function in C that implements your operator or function in the OpenMx backend. Then you will add an entry to the symbol table so that the OpenMx front-end is aware of your new function. The final step is to add a unit test that verifies the correct behavior of your new operator or function.

Step 1: Write a C function

For this howto, we will pretend we are implementing the binary + operator that is already available in OpenMx. You will need to either add your function to the file omxAlgebraFunctions.c or to your own file that you include in the /src subdirectory. All matrix functions have the same signature in C. They all accept three arguments and do not return a value.

void omxMatrixAdd(omxMatrix** matList, int numArgs, omxMatrix* result) {

The three arguments they accept are matList, numArgs, and result. matList an the input array of omxMatrix pointers. This array contains the matrices that are inputs to your function or operator. numArgs is the length of the matList array argument. result is the omxMatrix pointer to store the result of the matrix computation. If your matrix function always takes the same number of arguments, then you may choose not to use the matList argument (but it must appear in the signature). In our example, the + operator always takes two arguments so we may ignore matList.

The first step of your matrix function should be to process the input matrices and perform any validation checks.

  omxMatrix* inMat = matList[0];
    omxMatrix* addend = matList[1];
    
    /* Conformability Check! */
    if(inMat->cols != addend->cols || inMat->rows != addend->rows) 
        error("Non-conformable matrices in Matrix Addition.");

You may assume that the result parameter points to a valid omxMatrix reference. However the dimensions of the matrix that result is pointing to may not be correct for your computation. In the case of omxMatrixAdd, we shall replace the contents of result with the contents of the first operand:

     omxCopyMatrix(result, inMat); 

A common technique is simply to resize the result matrix to the desired dimensions. The following line is used by omxMatrixHorizCat (horizontal concatenation):

  if(result->rows != totalRows || result->cols != totalCols) {
        omxResizeMatrix(result, totalRows, totalCols, FALSE);
    }

Once the result matrix is the correct size, the final step is to perform the computation. When you are selecting or writing to particular row and column location across multiple matrices, you MUST use the omxMatrixElement and omxSetMatrixElement functions. This is because the matrix values may be stored in either row-major or column-major format. If you use the two helper functions then this detail will be abstracted away. Here is the entire omxMatrixAdd function:

void omxMatrixAdd(omxMatrix** matList, int numArgs, omxMatrix* result)
{
    omxMatrix* inMat = matList[0];
    omxMatrix* addend = matList[1];
    
    /* Conformability Check! */
    if(inMat->cols != addend->cols || inMat->rows != addend->rows) 
        error("Non-conformable matrices in Matrix Addition.");
        
    omxCopyMatrix(result, inMat);

    int rows = inMat->rows;
    int cols = inMat->cols;
    
    for(int i = 0; i < rows; i++) {
        for(int j = 0; j < cols; j++) {
            omxSetMatrixElement(result, i, j, 
                omxMatrixElement(result, i, j) +  
                omxMatrixElement(addend, i, j));
        }
    }       

}

Optionally, you may wish to declare your function in one of the header source files, such as omxAlgebraFunctions.h. This step is mandatory if you or another developer will be calling this function as a helper function to perform some computation in the back-end.

Step 2: Fill in the Symbol Table

Open the file data/omxSymbolTable.tab. You will add your own function add the bottom of the file. There are five columns you need to fill in. Entries within a row are tab-deliminated; you may use an arbitrary number of tabs in between columns. The first column is the numeric offset of your new function. Numeric offsets should appear in increasing order in the table with no gaps. The second column is a descriptive name for your operator or function. The third column is the operator or function name you which to use in R. This is the name that will be used when you write mxAlgebra() statements. There is no distinction between matrix operators (which are infix) and matrix functions (which are prefix), they are both parsed as names. The fourth column is the number of arguments for your function. If you want a variable number of arguments, use -1 in this column. The final column is the name of the C function that was implemented in step 1. Here is the entry for binary addition:

"Num"     "Operator name"     "R name"    "Number of arguments"   "Signature"
9       "Addition"      "+"     2           "omxMatrixAdd"

Step 3: Add a Test Case

You MUST write a test case for your new function. See the file models/passing/AlgebraComputePassing.R for examples on how to write a test case. You may either modify AlgebraSubstitutionPassing.R or add your own test case into the models/passing directory. When you are finished, run "make clean install test" to compile your new function and run all the test cases.