home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD2.mdf / c / library / dos / math / cephes / remes / remes.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-11-17  |  2.4 KB  |  149 lines
  1. /* remes.c
  2.  *
  3.  * This is an interactive program that computes least maximum polynomial
  4.  * and rational approximations.
  5.  * Last rev: 21 Sept 1986
  6.  */
  7. #include "remes.h"
  8.  
  9. main()
  10. {
  11. int i;
  12. int chgflg; /* flags changes to default values */
  13. double x;
  14. char s[40];
  15. char *sp;
  16.  
  17.  
  18. printf( "\nRational Approximation by Remes Algorithm\n\n" );
  19.  
  20. START:
  21.  
  22. /* Get operator commands. */
  23. remesa();
  24.  
  25. /* Jump to operator intervention point. */
  26. goto showg;
  27.  
  28.  
  29. LOOP:
  30.  
  31. iter += 1;
  32. printf( "Iteration %d\n", iter );
  33.  
  34. if( search != 0 )
  35.     {
  36. /* Search for error maxima. */
  37.     remess();
  38.     goto solveq;
  39.     }
  40.  
  41. showg:
  42. /* Display old values of guesses */
  43. /* and let user change them if desired */
  44. chgflg = 0;
  45. /* First get the step size if rational form */
  46. if( (d > 0) && (search != 0) )
  47.     { /* there is a denominator polynomial */
  48.     printf( "delta = %.4E ? ", delta );
  49.     gets( s );
  50. /* If input is not a null line, */
  51. /* then decode the number. */
  52.     if( s[0] != '\0' )
  53.         {
  54.         chgflg = 1;
  55.         sscanf( s, "%lf", &delta );
  56.         }
  57.     for( i=0; i<=neq; i++ )
  58.         step[i] *= delta;
  59.     }
  60.  
  61. /* Read in guesses for locations of solution. */
  62. for( i=0; i<neq; i++ )
  63.     {
  64.     printf( "x[%d] = %.4E ? ", i, xx[i] );
  65.     gets( s );
  66.     if( s[0] != '\0' )
  67.         {
  68.         chgflg = 1;
  69.         sscanf( s, "%lf", &xx[i] );
  70.         if( (d == 0) && ((config & ZER) == 0) )
  71.             mm[i] = xx[i];
  72.         }
  73.     }
  74.  
  75. if( (d > 0) || ((config & ZER) != 0) )
  76.     {
  77.     for( i=0; i <=neq; i++ )
  78.         {
  79.         printf( "peak[%d] = %.4E ? ", i, mm[i] );
  80.         gets( s );
  81.         if( s[0] != '\0' )
  82.             {
  83.             chgflg = 1;
  84.             sscanf( s, "%lf", &mm[i] );
  85.             }
  86.         }
  87.     }
  88.  
  89. /* If there were any changes to the default values */
  90. /* then reinitialize the step size array. */
  91. if( chgflg )
  92.     stpini();
  93.  
  94. solveq:
  95.  
  96. /* Solve equations. */
  97. remese();
  98.  
  99. goto whtnxt;
  100.  
  101.  
  102. ptabl:
  103.  
  104. /* Display the results */
  105. remesp();
  106.  
  107.  
  108. whtnxt:
  109. /* Test for convergence criteria. */
  110. if( (delta < 1.0e-15) || (spread < 1.0e-15) )
  111.     askitr = iter;
  112. if( askitr > iter )
  113.     goto LOOP;
  114.  
  115. /* Ask what to do next */
  116. printf(
  117. "Enter #, p(rint), w(rite), g(uess), x(it), t(runc), n(one) ?"
  118. );
  119. /* Get command line from operator */
  120. sp = &s[0];
  121. gets( sp );
  122.  
  123. if( *sp == 'w' )
  124.     { /* Write results to file */
  125.     remesw();
  126.     goto whtnxt;
  127.     }
  128. if( *sp == 'g' ) /* Modify the guesses */
  129.     goto showg;
  130. if( *sp == 'p' ) /* Display results */
  131.     goto ptabl;
  132. if( (*sp >= '1') && (*sp <= '9') )
  133.     { /* Numeric input is iteration count */
  134.     sscanf( sp, "%d", &askitr );
  135.     askitr += iter;
  136.     goto LOOP;
  137.     }
  138. if( *sp == 't' )
  139.     {
  140.     config |= TRUNC;
  141.     goto whtnxt;
  142.     }
  143. /* Close files and exit */
  144. if( *sp == 'x' )
  145.     exit(0);
  146. else
  147.     goto START;
  148. }
  149.