home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1994 January / usenetsourcesnewsgroupsinfomagicjanuary1994.iso / sources / unix / volume10 / complex-lib / cx_test.c < prev    next >
Encoding:
C/C++ Source or Header  |  1987-07-14  |  4.2 KB  |  182 lines
  1. /*
  2.     ctest -- complex arithmetic test
  3.  
  4.     last edit:    86/01/04    D A Gwyn
  5.  
  6.     SCCS ID:    @(#)cx_test.c    1.1 (modified for public version)
  7. */
  8.  
  9. #include    <stdio.h>
  10. #include    <math.h>
  11.  
  12. #include    <complex.h>
  13.  
  14. #define DEGRAD    57.2957795130823208767981548141051703324054724665642
  15.                     /* degrees per radian */
  16. #define Abs( x )    ((x) < 0 ? -(x) : (x))
  17. #define Max( a, b )    ((a) > (b) ? (a) : (b))
  18.  
  19. extern void    exit();
  20.  
  21. #define    Printf    (void)printf
  22.  
  23. #define    TOL    1.0e-10            /* tolerance for checks */
  24.  
  25. static int    errs = 0;        /* tally errors */
  26.  
  27. static void    CCheck(), RCheck();
  28. static double    RelDif();
  29.  
  30.  
  31. /*ARGSUSED*/
  32. main( argc, argv )
  33.     int    argc;
  34.     char    *argv[];
  35.     {
  36.     complex a, *bp, *cp;
  37.  
  38.     /* CxAllo test */
  39.     bp = CxAllo();
  40.     if ( bp == NULL )
  41.         {
  42.         Printf( "CxAllo failed\n" );
  43.         exit( 1 );
  44.         }
  45.  
  46.     /* CxReal, CxImag test */
  47.     CxReal( bp ) = 1.0;
  48.     CxImag( bp ) = 2.0;
  49.     RCheck( "CxReal", CxReal( bp ), 1.0 );
  50.     RCheck( "CxImag", CxImag( bp ), 2.0 );
  51.  
  52.     /* CxCons test */
  53.     cp = CxCons( &a, -3.0, -4.0);
  54.     CCheck( "CxCons 1", a, -3.0, -4.0 );
  55.     CCheck( "CxCons 2", *cp, -3.0, -4.0 );
  56.  
  57.     /* CxNeg test */
  58.     cp = CxNeg( &a );
  59.     CCheck( "CxNeg 1", a, 3.0, 4.0 );
  60.     CCheck( "CxNeg 2", *cp, 3.0, 4.0 );
  61.  
  62.     /* CxCopy test */
  63.     cp = CxCopy( bp, &a );
  64.     (void)CxCons( &a, 1.0, sqrt( 3.0 ) );
  65.     CCheck( "CxCopy 1", *bp, 3.0, 4.0 );
  66.     CCheck( "CxCopy 2", *cp, 3.0, 4.0 );
  67.  
  68.     /* CxAmpl, CxPhas test */
  69.     RCheck( "CxAmpl 1", CxAmpl( &a ), 2.0 );
  70.     RCheck( "CxPhas 1", CxPhas( &a ) * DEGRAD, 60.0 );
  71.     /* try other quadrants */
  72.     a.re = -a.re;
  73.     RCheck( "CxAmpl 2", CxAmpl( &a ), 2.0 );
  74.     RCheck( "CxPhas 2", CxPhas( &a ) * DEGRAD, 120.0 );
  75.     a.im = -a.im;
  76.     RCheck( "CxAmpl 3", CxAmpl( &a ), 2.0 );
  77.     RCheck( "CxPhas 3", CxPhas( &a ) * DEGRAD, -120.0 );
  78.     a.re = -a.re;
  79.     RCheck( "CxAmpl 4", CxAmpl( &a ), 2.0 );
  80.     RCheck( "CxPhas 4", CxPhas( &a ) * DEGRAD, -60.0 );
  81.     /* one more for good measure */
  82.     RCheck( "CxAmpl 5", CxAmpl( bp ), 5.0 );
  83.  
  84.     /* CxPhsr test */
  85.     cp = CxPhsr( &a, 100.0, -20.0 / DEGRAD );
  86.     RCheck( "CxPhsr 1", CxAmpl( &a ), 100.0 );
  87.     RCheck( "CxPhsr 2", CxPhas( &a ) * DEGRAD, -20.0 );
  88.     RCheck( "CxPhsr 3", CxAmpl( cp ), 100.0 );
  89.     RCheck( "CxPhsr 4", CxPhas( cp ) * DEGRAD, -20.0 );
  90.  
  91.     /* CxConj test */
  92.     cp = CxConj( bp );
  93.     CCheck( "CxConj 1", *bp, 3.0, -4.0 );
  94.     CCheck( "CxConj 2", *cp, 3.0, -4.0 );
  95.  
  96.     /* CxScal test */
  97.     cp = CxScal( bp, 2.0 );
  98.     CCheck( "CxScal 1", *bp, 6.0, -8.0 );
  99.     CCheck( "CxScal 2", *cp, 6.0, -8.0 );
  100.  
  101.     /* CxAdd test */
  102.     cp = CxAdd( CxCons( &a, -4.0, 11.0 ), bp );
  103.     CCheck( "CxAdd 1", a, 2.0, 3.0 );
  104.     CCheck( "CxAdd 2", *cp, 2.0, 3.0 );
  105.  
  106.     /* CxSub test */
  107.     cp = CxSub( CxCons( &a, 4.0, 7.0 ), bp );
  108.     CCheck( "CxSub 1", a, -2.0, 15.0 );
  109.     CCheck( "CxSub 2", *cp, -2.0, 15.0 );
  110.  
  111.     /* CxMul test */
  112.     cp = CxMul( CxCons( bp, -1.0, 3.0 ), CxCons( &a, 1.0, 2.0 ) );
  113.     CCheck( "CxMul 1", *bp, -7.0, 1.0 );
  114.     CCheck( "CxMul 2", *cp, -7.0, 1.0 );
  115.  
  116.     /* CxDiv test */
  117.     cp = CxDiv( bp, &a );
  118.     CCheck( "CxDiv 1", *bp, -1.0, 3.0 );
  119.     CCheck( "CxDiv 2", *cp, -1.0, 3.0 );
  120.  
  121.     /* CxSqrt and overlapping CxMul tests */
  122.     (void)CxCons( &a, -1.0, 2.0 );
  123.     cp = CxSqrt( CxMul( &a, &a ) );
  124.     CCheck( "CxSqrt 1", a, -1.0, 2.0 );
  125.     CCheck( "CxSqrt 2", *cp, -1.0, 2.0 );
  126.     (void)CxCons( &a, 3.0, 2.0 );
  127.     cp = CxSqrt( CxMul( &a, &a ) );
  128.     CCheck( "CxSqrt 3", a, 3.0, 2.0 );
  129.     CCheck( "CxSqrt 4", *cp, 3.0, 2.0 );
  130.  
  131.     /* CxFree "test" */
  132.     CxFree( bp );
  133.  
  134.     return errs;
  135.     }
  136.  
  137.  
  138. static void
  139. CCheck( s, c, r, i )            /* check complex number */
  140.     char    *s;            /* message string for failure */
  141.     complex    c;            /* complex to be checked */
  142.     double    r, i;            /* expected real, imaginary parts */
  143.     {
  144.     if ( RelDif( CxReal( &c ), r ) > TOL
  145.       || RelDif( CxImag( &c ), i ) > TOL
  146.        )    {
  147.         ++errs;
  148.         Printf( "%s; s.b. (%f,%f), was (%g,%g)\n",
  149.             s, r, i, c.re, c.im
  150.               );
  151.         }
  152.     }
  153.  
  154.  
  155. static void
  156. RCheck( s, d, r )            /* check real number */
  157.     char    *s;            /* message string for failure */
  158.     double    d;            /* real to be checked */
  159.     double    r;            /* expected value */
  160.     {
  161.     if ( RelDif( d, r ) > TOL )
  162.         {
  163.         ++errs;
  164.         Printf( "%s; s.b. %f, was %g\n", s, r, d );
  165.         }
  166.     }
  167.  
  168.  
  169. static double
  170. RelDif( a, b )            /* returns relative difference:    */
  171.     double    a, b;        /* 0.0 if exactly the same,
  172.                    otherwise ratio of difference
  173.                    to the larger of the two    */
  174.     {
  175.     double    c = Abs( a );
  176.     double    d = Abs( b );
  177.  
  178.     d = Max( c, d );
  179.  
  180.     return d == 0.0 ? 0.0 : Abs( a - b ) / d;
  181.     }
  182.