IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
May 1, 2007, 6:20:06 PM (19 years ago)
Author:
gusciora
Message:

Added tabular file of all psLib functions in Chapter 6 of the SDRS (Data
Manipulation). Most of the test file changes in this check-in are fairly
cosmetic.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psLib/test/math/tap_psSparse.c

    r13123 r13124  
    1 #include <stdio.h>
     1    #include <stdio.h>
    22#include <string.h>
    33#include <pslib.h>
     
    66#include "pstap.h"
    77
    8 int main (void)
     8int main(void)
    99{
    10     plan_tests(26);
     10    psLogSetFormat("HLNM");
     11    psLogSetLevel(PS_LOG_INFO);
     12    plan_tests(40);
     13
     14
     15    // test psSparseAlloc()
     16    {
     17        psMemId id = psMemGetId();
     18        psSparse *matrix = psSparseAlloc(10, 20);
     19        ok(matrix != NULL, "psSparse successfully allocated");
     20        skip_start(matrix == NULL, 12, "skipping tests because psSparseAlloc() returned NULL");
     21        ok(matrix->Aij != NULL, "psSparseAlloc() set ->Aij correctly");
     22        ok(matrix->Aij->n == 0, "psSparseAlloc() set ->Aij->n correctly");
     23        ok(matrix->Si != NULL, "psSparseAlloc() set ->Si correctly");
     24        ok(matrix->Si->n == 0, "psSparseAlloc() set ->Si->n correctly");
     25        ok(matrix->Sj != NULL, "psSparseAlloc() set ->Sj correctly");
     26        ok(matrix->Sj->n == 0, "psSparseAlloc() set ->Sj->n correctly");
     27        ok(matrix->Bfj != NULL, "psSparseAlloc() set ->Bfj correctly");
     28        ok(matrix->Bfj->n == 10, "psSparseAlloc() set ->Bfj->n correctly");
     29        ok(matrix->Qii != NULL, "psSparseAlloc() set ->Qii correctly");
     30        ok(matrix->Qii->n == 10, "psSparseAlloc() set ->Qii->n correctly");
     31        ok(matrix->Nelem == 0, "psSparseAlloc() set ->Nelem correctly");
     32        ok(matrix->Nrows == 10, "psSparseAlloc() set ->Nrows correctly");
     33        skip_end();
     34        psFree(matrix);
     35        ok(!psMemCheckLeaks(id, NULL, NULL, false), "no memory leaks");
     36    }
     37
    1138
    1239    // test psSparseSolve for a simple normal example matrix
    1340    {
    1441        psMemId id = psMemGetId();
    15 
    16         //solve a normalized matrix equation with psSparseSolve
    17 
    1842        // the basic equation is Ax = b
    1943
     
    2246        // Ax = b for x using psSparseSolve.  compare with the input values for x.
    2347
    24         psSparse *matrix = psSparseAlloc (100, 100);
     48        psSparse *matrix = psSparseAlloc(100, 100);
    2549        ok(matrix != NULL, "psSparse successfully allocated");
    2650        skip_start(matrix == NULL, 5, "Skipping tests because psSparseAlloc() failed");
    2751
    28         for (int i = 0; i < 100; i++) {
    29             psSparseMatrixElement (matrix, i, i, 1.0);
     52        for(int i = 0; i < 100; i++) {
     53            psSparseMatrixElement(matrix, i, i, 1.0);
    3054            if (i + 1 < 100) {
    31                 psSparseMatrixElement (matrix, i + 1, i, 0.1);
     55                psSparseMatrixElement(matrix, i + 1, i, 0.1);
    3256            }
    3357        }
     
    3761        psSparseResort(matrix);
    3862
    39         psVector *xRef = psVectorAlloc (100, PS_TYPE_F32);
    40         for (int i = 0; i < 100; i++)
    41         {
     63        psVector *xRef = psVectorAlloc(100, PS_TYPE_F32);
     64        for (int i = 0; i < 100; i++) {
    4265            xRef->data.F32[i] = 1.0;
    4366        }
     
    6083        float dS = 0;
    6184        float dS2 = 0;
    62         for (int i = 0; i < 100; i++)
    63         {
     85        for (int i = 0; i < 100; i++) {
    6486            float dX = xRef->data.F32[i] - xFit->data.F32[i];
    65             // fprintf (stderr, "%5.3f %5.3f %5.3f %5.3f\n", bVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
     87            // fprintf(stderr, "%5.3f %5.3f %5.3f %5.3f\n", bVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
    6688            dS2 += PS_SQR(dX);
    6789            dS += dX;
    6890        }
    6991        dS /= 100.0;
    70         dS2 = sqrt (dS2/100.0 - dS*dS);
    71 
    72         is_float_tol (dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
    73 
    74         psFree (matrix);
    75         psFree (xRef);
    76         psFree (bVec);
    77         psFree (xFit);
    78 
    79         skip_end();
    80         ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     92        dS2 = sqrt(dS2/100.0 - dS*dS);
     93
     94        is_float_tol(dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
     95
     96        psFree(matrix);
     97        psFree(xRef);
     98        psFree(bVec);
     99        psFree(xFit);
     100
     101        skip_end();
     102        ok(!psMemCheckLeaks(id, NULL, NULL, false), "no memory leaks");
    81103    }
    82104
     
    84106    {
    85107        psMemId id = psMemGetId();
    86 
    87         //solve a non-normalized matrix equation with psSparseSolve
    88108        // the basic equation is Ax = b
     109
    89110        // create a matrix A with diagonals of 1 and a small number of off diagonal elements.
    90111        // construct a vector x, construct the corresponding vector b by multiplication. solve
    91112        // Ax = b for x using psSparseSolve.  compare with the input values for x.
    92113
    93         psSparse *matrix = psSparseAlloc (100, 100);
     114        psSparse *matrix = psSparseAlloc(100, 100);
    94115        ok(matrix != NULL, "psSparse successfully allocated");
    95116        skip_start(matrix == NULL, 5, "Skipping tests because psSparseAlloc() failed");
    96117
    97118        for (int i = 0; i < 100; i++) {
    98             psSparseMatrixElement (matrix, i, i, 5.0);
     119            psSparseMatrixElement(matrix, i, i, 5.0);
    99120            if (i + 1 < 100) {
    100                 psSparseMatrixElement (matrix, i + 1, i, 0.1);
     121                psSparseMatrixElement(matrix, i + 1, i, 0.1);
    101122            }
    102123        }
     
    106127        psSparseResort(matrix);
    107128
    108         psVector *xRef = psVectorAlloc (100, PS_TYPE_F32);
    109         for (int i = 0; i < 100; i++)
    110         {
     129        psVector *xRef = psVectorAlloc(100, PS_TYPE_F32);
     130        for (int i = 0; i < 100; i++) {
    111131            xRef->data.F32[i] = 1.0;
    112132        }
     
    123143        constraint.paramDelta = 1e8;
    124144
    125         // solve for normalization terms (need include local sky?)
     145        // solve for normalization terms(need include local sky?)
    126146        psVector *xFit = psSparseSolve(NULL, constraint, matrix, 4);
    127147
     
    129149        float dS = 0;
    130150        float dS2 = 0;
    131         for (int i = 0; i < 100; i++)
    132         {
     151        for (int i = 0; i < 100; i++) {
    133152            float dX = xRef->data.F32[i] - xFit->data.F32[i];
    134             // fprintf (stderr, "%5.3f %5.3f %5.3f %5.3f\n", bVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
     153            // fprintf(stderr, "%5.3f %5.3f %5.3f %5.3f\n", bVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
    135154            dS2 += PS_SQR(dX);
    136155            dS += dX;
    137156        }
    138157        dS /= 100.0;
    139         dS2 = sqrt (dS2/100.0 - dS*dS);
    140 
    141         is_float_tol (dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
    142 
    143         psFree (matrix);
    144         psFree (xRef);
    145         psFree (bVec);
    146         psFree (xFit);
    147 
    148         skip_end();
    149         ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     158        dS2 = sqrt(dS2/100.0 - dS*dS);
     159
     160        is_float_tol(dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
     161
     162        psFree(matrix);
     163        psFree(xRef);
     164        psFree(bVec);
     165        psFree(xFit);
     166
     167        skip_end();
     168        ok(!psMemCheckLeaks(id, NULL, NULL, false), "no memory leaks");
    150169    }
    151170
     
    164183    {
    165184        psMemId id = psMemGetId();
    166 
    167         //solve a simple, small matrix equation
    168 
    169         // the basic equation (Ax = b) is:
     185        // the basic equation(Ax = b) is:
    170186        // |S B'||x| = |f|
    171187        // |B Q ||y| = |g|
     
    180196
    181197        // construct the sparse matrix
    182         psSparse *matrix = psSparseAlloc (Nrows, Nrows);
     198        psSparse *matrix = psSparseAlloc(Nrows, Nrows);
    183199        ok(matrix != NULL, "psSparse successfully allocated");
    184200        skip_start(matrix == NULL, 5, "Skipping tests because psSparseAlloc() failed");
     
    190206
    191207        // border region has a width of 1:
    192         psSparseBorder *border = psSparseBorderAlloc (matrix, Nborder);
     208        psSparseBorder *border = psSparseBorderAlloc(matrix, Nborder);
    193209
    194210        // construct the B component:
    195         psSparseBorderElementB (border, 0, 0, 0.1);
    196         psSparseBorderElementB (border, 1, 0, 0.2);
     211        psSparseBorderElementB(border, 0, 0, 0.1);
     212        psSparseBorderElementB(border, 1, 0, 0.2);
    197213
    198214        // construct the T component:
    199         psSparseBorderElementT (border, 0, 0, 0.5);
     215        psSparseBorderElementT(border, 0, 0, 0.5);
    200216
    201217        // construct the X and Y vectors:
    202         psVector *xRef = psVectorAlloc (Nrows, PS_TYPE_F32);
     218        psVector *xRef = psVectorAlloc(Nrows, PS_TYPE_F32);
    203219        xRef->data.F32[0] = 1.0;
    204220        xRef->data.F32[1] = 1.0;
    205         psVector *yRef = psVectorAlloc (Nborder, PS_TYPE_F32);
     221        psVector *yRef = psVectorAlloc(Nborder, PS_TYPE_F32);
    206222        yRef->data.F32[0] = 0.5;
    207223
     
    210226
    211227        // test the support functions: LowerProduct
    212         fVec = psSparseBorderLowerProduct (NULL, border, xRef);
    213         is_float (fVec->n, 1.0, "f dimen: %d", fVec->n);
    214         is_float (fVec->data.F32[0], 0.3, "f(0): %f", fVec->data.F32[0]);
    215         psFree (fVec);
     228        fVec = psSparseBorderLowerProduct(NULL, border, xRef);
     229        is_float(fVec->n, 1.0, "f dimen: %d", fVec->n);
     230        is_float(fVec->data.F32[0], 0.3, "f(0): %f", fVec->data.F32[0]);
     231        psFree(fVec);
    216232
    217233        // test the support functions: Upper Product
    218         fVec = psSparseBorderUpperProduct (NULL, border, yRef);
    219         is_float (fVec->n, 2.0, "f dimen: %d", fVec->n);
    220         is_float (fVec->data.F32[0], 0.05, "f(0): %f", fVec->data.F32[0]);
    221         is_float (fVec->data.F32[1], 0.10, "f(0): %f", fVec->data.F32[1]);
    222         psFree (fVec);
     234        fVec = psSparseBorderUpperProduct(NULL, border, yRef);
     235        is_float(fVec->n, 2.0, "f dimen: %d", fVec->n);
     236        is_float(fVec->data.F32[0], 0.05, "f(0): %f", fVec->data.F32[0]);
     237        is_float(fVec->data.F32[1], 0.10, "f(0): %f", fVec->data.F32[1]);
     238        psFree(fVec);
    223239
    224240        // test the support functions: Square Product
    225         fVec = psSparseBorderSquareProduct (NULL, border, yRef);
    226         is_float (fVec->n, 1.0, "f dimen: %d", fVec->n);
    227         is_float (fVec->data.F32[0], 0.25, "f(0): %f", fVec->data.F32[0]);
    228         psFree (fVec);
     241        fVec = psSparseBorderSquareProduct(NULL, border, yRef);
     242        is_float(fVec->n, 1.0, "f dimen: %d", fVec->n);
     243        is_float(fVec->data.F32[0], 0.25, "f(0): %f", fVec->data.F32[0]);
     244        psFree(fVec);
    229245
    230246        fVec = NULL;
    231247        psVector *gVec = NULL;
    232         psSparseBorderMultiply (&fVec, &gVec, border, xRef, yRef);
    233         is_float (fVec->data.F32[0], 1.15, "f(0): %f", fVec->data.F32[0]);
    234         is_float (fVec->data.F32[1], 1.20, "f(1): %f", fVec->data.F32[1]);
    235         is_float (gVec->data.F32[0], 0.55, "g(0): %f", gVec->data.F32[0]);
     248        psSparseBorderMultiply(&fVec, &gVec, border, xRef, yRef);
     249        is_float(fVec->data.F32[0], 1.15, "f(0): %f", fVec->data.F32[0]);
     250        is_float(fVec->data.F32[1], 1.20, "f(1): %f", fVec->data.F32[1]);
     251        is_float(gVec->data.F32[0], 0.55, "g(0): %f", gVec->data.F32[0]);
    236252
    237253        // supply the fVec and gVec data to the border
     
    253269        psVector *yFit = NULL;
    254270        psSparseBorderSolve(&xFit, &yFit, constraint, border, 4);
    255         is_float_tol (xFit->data.F32[0], 1.0, 1e-4, "f(0): %f", xFit->data.F32[0]);
    256         is_float_tol (xFit->data.F32[1], 1.0, 1e-4, "f(1): %f", xFit->data.F32[1]);
    257         is_float_tol (yFit->data.F32[0], 0.5, 1e-4, "g(0): %f", yFit->data.F32[0]);
    258 
    259         psFree (xFit);
    260         psFree (yFit);
    261         psFree (fVec);
    262         psFree (gVec);
    263         psFree (xRef);
    264         psFree (yRef);
    265         psFree (border);
    266         psFree (matrix);
    267         skip_end();
    268         ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
     271        is_float_tol(xFit->data.F32[0], 1.0, 1e-4, "f(0): %f", xFit->data.F32[0]);
     272        is_float_tol(xFit->data.F32[1], 1.0, 1e-4, "f(1): %f", xFit->data.F32[1]);
     273        is_float_tol(yFit->data.F32[0], 0.5, 1e-4, "g(0): %f", yFit->data.F32[0]);
     274
     275        psFree(xFit);
     276        psFree(yFit);
     277        psFree(fVec);
     278        psFree(gVec);
     279        psFree(xRef);
     280        psFree(yRef);
     281        psFree(border);
     282        psFree(matrix);
     283        skip_end();
     284        ok(!psMemCheckLeaks(id, NULL, NULL, false), "no memory leaks");
    269285    }
    270286
     
    272288    {
    273289        psMemId id = psMemGetId();
    274 
    275         // solve a simple, small matrix equation
    276 
    277         // the basic equation (Ax = b) is:
     290        // the basic equation(Ax = b) is:
    278291        // |S B'||x| = |f|
    279292        // |B Q ||y| = |g|
     
    288301
    289302        // construct the sparse matrix
    290         psSparse *matrix = psSparseAlloc (Nrows, Nrows);
     303        psSparse *matrix = psSparseAlloc(Nrows, Nrows);
    291304        ok(matrix != NULL, "psSparse successfully allocated");
    292305        skip_start(matrix == NULL, 5, "Skipping tests because psSparseAlloc() failed");
     
    301314
    302315        // border region has a width of 1:
    303         psSparseBorder *border = psSparseBorderAlloc (matrix, Nborder);
     316        psSparseBorder *border = psSparseBorderAlloc(matrix, Nborder);
    304317
    305318        // construct the B component:
    306         for (int i = 0; i < Nrows; i++)
    307         {
     319        for (int i = 0; i < Nrows; i++) {
    308320            for (int j = 0; j < Nborder; j++) {
    309                 psSparseBorderElementB (border, i, j, 0.1);
     321                psSparseBorderElementB(border, i, j, 0.1);
    310322            }
    311323        }
     
    313325
    314326        // construct the Q component:
    315         for (int i = 0; i < Nborder; i++)
    316         {
     327        for (int i = 0; i < Nborder; i++) {
    317328            for (int j = 0; j < Nborder; j++) {
    318                 psSparseBorderElementT (border, i, j, i+j+2);
     329                psSparseBorderElementT(border, i, j, i+j+2);
    319330            }
    320331        }
    321332
    322333        // construct the X and Y vectors:
    323         psVector *xRef = psVectorAlloc (Nrows, PS_TYPE_F32);
    324         for (int i = 0; i < Nrows; i++)
    325         {
     334        psVector *xRef = psVectorAlloc(Nrows, PS_TYPE_F32);
     335        for (int i = 0; i < Nrows; i++) {
    326336            xRef->data.F32[i] = 1.0;
    327337        }
    328         psVector *yRef = psVectorAlloc (Nborder, PS_TYPE_F32);
    329         for (int i = 0; i < Nborder; i++)
    330         {
     338        psVector *yRef = psVectorAlloc(Nborder, PS_TYPE_F32);
     339        for (int i = 0; i < Nborder; i++) {
    331340            yRef->data.F32[i] = 1.0;
    332341        }
     
    335344        psVector *fVec = NULL;
    336345        psVector *gVec = NULL;
    337         psSparseBorderMultiply (&fVec, &gVec, border, xRef, yRef);
     346        psSparseBorderMultiply(&fVec, &gVec, border, xRef, yRef);
    338347
    339348        // supply the fVec and gVec data to the border
     
    351360        constraint.paramDelta = 1e8;
    352361
    353         // solve for normalization terms (need include local sky?)
     362        // solve for normalization terms(need include local sky?)
    354363        psVector *xFit = NULL;
    355364        psVector *yFit = NULL;
     
    359368        float dS = 0;
    360369        float dS2 = 0;
    361         for (int i = 0; i < Nrows; i++)
    362         {
     370        for (int i = 0; i < Nrows; i++) {
    363371            float dX = xRef->data.F32[i] - xFit->data.F32[i];
    364             // fprintf (stderr, "%5.3f %5.3f %5.3f %5.3f\n", fVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
     372            // fprintf(stderr, "%5.3f %5.3f %5.3f %5.3f\n", fVec->data.F32[i], xRef->data.F32[i], xFit->data.F32[i], dX);
    365373            dS2 += PS_SQR(dX);
    366374            dS += dX;
    367375        }
    368376        dS /= Nrows;
    369         dS2 = sqrt (dS2/Nrows - dS*dS);
    370         is_float_tol (dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
     377        dS2 = sqrt(dS2/Nrows - dS*dS);
     378        is_float_tol(dS2, 0.0, 1e-4, "scatter: %.20f", dS2);
    371379
    372380        // measure stdev between yFit and yRef
    373381        float dY = yRef->data.F32[0] - yFit->data.F32[0];
    374         // fprintf (stderr, "%5.3f %5.3f %5.3f %5.3f\n", gVec->data.F32[0], yRef->data.F32[0], yFit->data.F32[0], dS);
    375         is_float_tol (dY, 0.0, 2e-4, "scatter: %.20f", dY);
    376 
    377         psFree (xRef);
    378         psFree (yRef);
    379         psFree (xFit);
    380         psFree (yFit);
    381         psFree (fVec);
    382         psFree (gVec);
    383         psFree (matrix);
    384         psFree (border);
    385 
    386         skip_end();
    387         ok(!psMemCheckLeaks (id, NULL, NULL, false), "no memory leaks");
    388     }
    389 
    390     return exit_status();
     382        // fprintf(stderr, "%5.3f %5.3f %5.3f %5.3f\n", gVec->data.F32[0], yRef->data.F32[0], yFit->data.F32[0], dS);
     383        is_float_tol(dY, 0.0, 2e-4, "scatter: %.20f", dY);
     384
     385        psFree(xRef);
     386        psFree(yRef);
     387        psFree(xFit);
     388        psFree(yFit);
     389        psFree(fVec);
     390        psFree(gVec);
     391        psFree(matrix);
     392        psFree(border);
     393
     394        skip_end();
     395        ok(!psMemCheckLeaks(id, NULL, NULL, false), "no memory leaks");
     396    }
    391397}
Note: See TracChangeset for help on using the changeset viewer.