Index: trunk/psLib/src/dataManip/psMatrix.c
===================================================================
--- trunk/psLib/src/dataManip/psMatrix.c	(revision 4385)
+++ trunk/psLib/src/dataManip/psMatrix.c	(revision 4388)
@@ -21,6 +21,6 @@
  *  @author Robert DeSonia, MHPCC
  *
- *  @version $Revision: 1.33 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-06-25 00:51:28 $
+ *  @version $Revision: 1.34 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-06-25 01:15:01 $
  *
  *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -185,5 +185,5 @@
 /*****************************************************************************/
 
-psImage* psMatrixLUD(psImage* outImage, psVector** outPerm, psImage* inImage)
+psImage* psMatrixLUD(psImage* out, psVector** perm, psImage* in)
 {
     psS32 signum = 0;
@@ -194,27 +194,27 @@
 
 
-    #define psMatrixLUD_EXIT {psFree(outImage); return NULL;}
-
-    // Error checks
-    PS_ASSERT_GENERAL_IMAGE_NON_NULL(inImage, psMatrixLUD_EXIT);
-    PS_CHECK_POINTERS(inImage, outImage, psMatrixLUD_EXIT);
-    PS_CHECK_DIMEN_AND_TYPE(inImage, PS_DIMEN_IMAGE, psMatrixLUD_EXIT);
-    PS_ASSERT_GENERAL_PTR_NON_NULL(outPerm, psMatrixLUD_EXIT);
-
-    outImage = psImageRecycle(outImage, inImage->numCols, inImage->numRows, inImage->type.type);
-
-    PS_CHECK_SQUARE(inImage, psMatrixLUD_EXIT); // gsl_linalg_LU_decomp would fail on non-square input.
-    PS_CHECK_SQUARE(outImage, psMatrixLUD_EXIT);
+    #define psMatrixLUD_EXIT {psFree(out); return NULL;}
+
+    // Error checks
+    PS_ASSERT_GENERAL_IMAGE_NON_NULL(in, psMatrixLUD_EXIT);
+    PS_CHECK_POINTERS(in, out, psMatrixLUD_EXIT);
+    PS_CHECK_DIMEN_AND_TYPE(in, PS_DIMEN_IMAGE, psMatrixLUD_EXIT);
+    PS_ASSERT_GENERAL_PTR_NON_NULL(perm, psMatrixLUD_EXIT);
+
+    out = psImageRecycle(out, in->numCols, in->numRows, in->type.type);
+
+    PS_CHECK_SQUARE(in, psMatrixLUD_EXIT); // gsl_linalg_LU_decomp would fail on non-square input.
+    PS_CHECK_SQUARE(out, psMatrixLUD_EXIT);
 
     // Initialize data
-    numRows = inImage->numRows;
-    numCols = inImage->numCols;
+    numRows = in->numRows;
+    numCols = in->numCols;
 
     // Initialize GSL data
     permGSL.size = numCols;
     if (sizeof(size_t) == 4) {
-        *outPerm = psVectorRecycle(*outPerm, numCols, PS_TYPE_S32);
+        *perm = psVectorRecycle(*perm, numCols, PS_TYPE_S32);
     } else if (sizeof(size_t) == 8) {
-        *outPerm = psVectorRecycle(*outPerm, numCols, PS_TYPE_S64);
+        *perm = psVectorRecycle(*perm, numCols, PS_TYPE_S64);
     } else {
         psError(PS_ERR_UNKNOWN, true,
@@ -224,10 +224,10 @@
     }
 
-    (*outPerm)->n = numCols;
-    permGSL.data = (psPtr)((*outPerm)->data.U8);
+    (*perm)->n = numCols;
+    permGSL.data = (psPtr)((*perm)->data.U8);
     lu = gsl_matrix_alloc(numRows, numCols);
 
     // Copy psImage data into GSL matrix data
-    psImageToGslMatrix(lu, inImage);
+    psImageToGslMatrix(lu, in);
 
     // Calculate LU decomposition
@@ -235,58 +235,58 @@
 
     // Copy GSL matrix data to psImage data
-    gslMatrixToPsImage(outImage, lu);
+    gslMatrixToPsImage(out, lu);
 
     // Free GSL data
     gsl_matrix_free(lu);
 
-    return outImage;
-}
-
-psVector* psMatrixLUSolve(psVector* outVector, const psImage* inImage, const psVector* inVector,
-                          const psVector* inPerm)
+    return out;
+}
+
+psVector* psMatrixLUSolve(psVector* out, const psImage* LU, const psVector* RHS,
+                          const psVector* perm)
 {
     psS32 numRows = 0;
     psS32 numCols = 0;
     gsl_matrix *lu;
-    gsl_permutation perm;
+    gsl_permutation permGSL;
     gsl_vector *b = NULL;
     gsl_vector *x = NULL;
 
-    #define LUSOLVE_CLEANUP {psFree(outVector); return NULL;}
-
-    // Error checks
-    PS_ASSERT_GENERAL_IMAGE_NON_NULL(inImage, LUSOLVE_CLEANUP);
-    PS_CHECK_DIMEN_AND_TYPE(inImage, PS_DIMEN_IMAGE, LUSOLVE_CLEANUP);
-    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(inImage, LUSOLVE_CLEANUP);
-    PS_ASSERT_GENERAL_VECTOR_NON_NULL(inVector, LUSOLVE_CLEANUP);
-    PS_CHECK_DIMEN_AND_TYPE(inVector, PS_DIMEN_VECTOR, LUSOLVE_CLEANUP);
-    PS_ASSERT_GENERAL_VECTOR_NON_NULL(inPerm, LUSOLVE_CLEANUP);
-
-    outVector = psVectorRecycle(outVector, inImage->numRows, inImage->type.type);
-
-    PS_CHECK_POINTERS(outVector, inVector, LUSOLVE_CLEANUP);
-    PS_CHECK_POINTERS(inVector, inPerm, LUSOLVE_CLEANUP);
-    PS_CHECK_POINTERS(outVector, inPerm, LUSOLVE_CLEANUP);
+    #define LUSOLVE_CLEANUP {psFree(out); return NULL;}
+
+    // Error checks
+    PS_ASSERT_GENERAL_IMAGE_NON_NULL(LU, LUSOLVE_CLEANUP);
+    PS_CHECK_DIMEN_AND_TYPE(LU, PS_DIMEN_IMAGE, LUSOLVE_CLEANUP);
+    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(LU, LUSOLVE_CLEANUP);
+    PS_ASSERT_GENERAL_VECTOR_NON_NULL(RHS, LUSOLVE_CLEANUP);
+    PS_CHECK_DIMEN_AND_TYPE(RHS, PS_DIMEN_VECTOR, LUSOLVE_CLEANUP);
+    PS_ASSERT_GENERAL_VECTOR_NON_NULL(perm, LUSOLVE_CLEANUP);
+
+    out = psVectorRecycle(out, LU->numRows, LU->type.type);
+
+    PS_CHECK_POINTERS(out, RHS, LUSOLVE_CLEANUP);
+    PS_CHECK_POINTERS(RHS, perm, LUSOLVE_CLEANUP);
+    PS_CHECK_POINTERS(out, perm, LUSOLVE_CLEANUP);
 
     // Initialize data
-    numRows = inImage->numRows;
-    numCols = inImage->numCols;
+    numRows = LU->numRows;
+    numCols = LU->numCols;
 
     // Initialize GSL data
     lu = gsl_matrix_alloc(numRows, numCols);
-    psImageToGslMatrix(lu, inImage);
-    b = gsl_vector_alloc(inVector->n);
-    psVectorToGslVector(b, inVector);
-    x = gsl_vector_alloc(inVector->n);
-
-    outVector->n = numCols;
-    perm.size = inPerm->n;
-    perm.data = (psPtr)(inPerm->data.U8);
+    psImageToGslMatrix(lu, LU);
+    b = gsl_vector_alloc(RHS->n);
+    psVectorToGslVector(b, RHS);
+    x = gsl_vector_alloc(RHS->n);
+
+    out->n = numCols;
+    permGSL.size = perm->n;
+    permGSL.data = (psPtr)(perm->data.U8);
 
     // Solve for {x} in equation: {b} = [A]{x}
-    gsl_linalg_LU_solve(lu, &perm, b, x);
+    gsl_linalg_LU_solve(lu, &permGSL, b, x);
 
     // Copy GSL vector data to psVector data
-    gslVectorToPsVector(outVector, x);
+    gslVectorToPsVector(out, x);
 
     // Free GSL data
@@ -295,5 +295,5 @@
     gsl_matrix_free(lu);
 
-    return outVector;
+    return out;
 }
 
@@ -383,5 +383,5 @@
 }
 
-psImage* psMatrixMultiply(psImage* outImage, psImage* inImage1, psImage* inImage2)
+psImage* psMatrixMultiply(psImage* out, psImage* in1, psImage* in2)
 {
     gsl_matrix *m1 = NULL;
@@ -389,33 +389,33 @@
     gsl_matrix *m3 = NULL;
 
-    #define MULTIPLY_CLEANUP { psFree(outImage); return NULL; }
-
-    // Error checks
-    PS_ASSERT_GENERAL_IMAGE_NON_NULL(inImage1, MULTIPLY_CLEANUP);
-    PS_ASSERT_GENERAL_IMAGE_NON_NULL(inImage2, MULTIPLY_CLEANUP);
-    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(inImage1, MULTIPLY_CLEANUP);
-    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(inImage2, MULTIPLY_CLEANUP);
-    PS_CHECK_DIMEN_AND_TYPE(inImage1, PS_DIMEN_IMAGE, MULTIPLY_CLEANUP);
-    PS_CHECK_DIMEN_AND_TYPE(inImage2, PS_DIMEN_IMAGE, MULTIPLY_CLEANUP);
-    PS_CHECK_POINTERS(inImage1, outImage, MULTIPLY_CLEANUP);
-    PS_CHECK_POINTERS(inImage1, inImage2, MULTIPLY_CLEANUP);
-
-    if (inImage1->numRows != inImage2->numCols) {
-        psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Invalid operation: number of rows of inImage1 != number of cols of inImage2.");
+    #define MULTIPLY_CLEANUP { psFree(out); return NULL; }
+
+    // Error checks
+    PS_ASSERT_GENERAL_IMAGE_NON_NULL(in1, MULTIPLY_CLEANUP);
+    PS_ASSERT_GENERAL_IMAGE_NON_NULL(in2, MULTIPLY_CLEANUP);
+    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(in1, MULTIPLY_CLEANUP);
+    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(in2, MULTIPLY_CLEANUP);
+    PS_CHECK_DIMEN_AND_TYPE(in1, PS_DIMEN_IMAGE, MULTIPLY_CLEANUP);
+    PS_CHECK_DIMEN_AND_TYPE(in2, PS_DIMEN_IMAGE, MULTIPLY_CLEANUP);
+    PS_CHECK_POINTERS(in1, out, MULTIPLY_CLEANUP);
+    PS_CHECK_POINTERS(in1, in2, MULTIPLY_CLEANUP);
+
+    if (in1->numRows != in2->numCols) {
+        psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Invalid operation: number of rows of in1 != number of cols of in2.");
         MULTIPLY_CLEANUP;
     }
-    if (inImage1->type.type != inImage2->type.type) {
-        psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Invalid operation: data types of inImage1 and inImage2 must match.");
+    if (in1->type.type != in2->type.type) {
+        psError(PS_ERR_BAD_PARAMETER_SIZE, true, "Invalid operation: data types of in1 and in2 must match.");
         MULTIPLY_CLEANUP;
     }
 
-    outImage = psImageRecycle(outImage, inImage2->numCols, inImage1->numRows, inImage2->type.type);
+    out = psImageRecycle(out, in2->numCols, in1->numRows, in2->type.type);
 
     // Initialize GSL data
-    m1 = gsl_matrix_alloc(inImage1->numRows, inImage1->numCols);
-    psImageToGslMatrix(m1, inImage1);
-    m2 = gsl_matrix_alloc(inImage2->numRows, inImage2->numCols);
-    psImageToGslMatrix(m2, inImage2);
-    m3 = gsl_matrix_alloc(outImage->numRows, outImage->numCols);
+    m1 = gsl_matrix_alloc(in1->numRows, in1->numCols);
+    psImageToGslMatrix(m1, in1);
+    m2 = gsl_matrix_alloc(in2->numRows, in2->numCols);
+    psImageToGslMatrix(m2, in2);
+    m3 = gsl_matrix_alloc(out->numRows, out->numCols);
 
     // Perform multiplication
@@ -424,5 +424,5 @@
 
     // Copy GSL matrix data to psImage data
-    gslMatrixToPsImage(outImage, m3);
+    gslMatrixToPsImage(out, m3);
 
     // Free GSL structs
@@ -431,8 +431,8 @@
     gsl_matrix_free(m3);
 
-    return outImage;
-}
-
-psImage* psMatrixTranspose(psImage* outImage, const psImage* inImage)
+    return out;
+}
+
+psImage* psMatrixTranspose(psImage* out, const psImage* in)
 {
     psU32 i = 0;
@@ -443,18 +443,18 @@
     psS32 numColsOut = 0;
 
-    #define TRANSPOSE_CLEANUP { psFree(outImage); return NULL; }
-    // Error checks
-    PS_ASSERT_GENERAL_IMAGE_NON_NULL(inImage, TRANSPOSE_CLEANUP);
-    PS_CHECK_DIMEN_AND_TYPE(inImage, PS_DIMEN_IMAGE, TRANSPOSE_CLEANUP);
-    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(inImage, TRANSPOSE_CLEANUP);
-    PS_CHECK_POINTERS(inImage, outImage, TRANSPOSE_CLEANUP);
-
-    outImage = psImageRecycle(outImage, inImage->numCols, inImage->numRows, inImage->type.type);
+    #define TRANSPOSE_CLEANUP { psFree(out); return NULL; }
+    // Error checks
+    PS_ASSERT_GENERAL_IMAGE_NON_NULL(in, TRANSPOSE_CLEANUP);
+    PS_CHECK_DIMEN_AND_TYPE(in, PS_DIMEN_IMAGE, TRANSPOSE_CLEANUP);
+    PS_ASSERT_GENERAL_IMAGE_NON_EMPTY(in, TRANSPOSE_CLEANUP);
+    PS_CHECK_POINTERS(in, out, TRANSPOSE_CLEANUP);
+
+    out = psImageRecycle(out, in->numCols, in->numRows, in->type.type);
 
     // Initialize data
-    numRowsIn = inImage->numRows;
-    numColsIn = inImage->numCols;
-    numRowsOut = outImage->numRows;
-    numColsOut = outImage->numCols;
+    numRowsIn = in->numRows;
+    numColsIn = in->numCols;
+    numRowsOut = out->numRows;
+    numColsOut = out->numCols;
 
     if(numRowsIn!=numColsOut && numRowsOut!=numColsIn) {
@@ -463,8 +463,8 @@
     }
 
-    if(outImage->type.type == PS_TYPE_F32) {
+    if(out->type.type == PS_TYPE_F32) {
         for(i=0; i<numRowsOut; i++) {
             for(j=0; j<numColsOut; j++) {
-                outImage->data.F32[i][j] = inImage->data.F32[j][i];
+                out->data.F32[i][j] = in->data.F32[j][i];
             }
         }
@@ -472,10 +472,10 @@
         for(i=0; i<numRowsOut; i++) {
             for(j=0; j<numColsOut; j++) {
-                outImage->data.F64[i][j] = inImage->data.F64[j][i];
+                out->data.F64[i][j] = in->data.F64[j][i];
             }
         }
     }
 
-    return outImage;
+    return out;
 }