Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/libcfa/rational.c

    r6c6455f r561f730  
    1010// Created On       : Wed Apr  6 17:54:28 2016
    1111// Last Modified By : Peter A. Buhr
    12 // Last Modified On : Mon May 15 21:29:23 2017
    13 // Update Count     : 149
     12// Last Modified On : Sun May 14 17:25:19 2017
     13// Update Count     : 131
    1414//
    1515
     
    190190// conversion
    191191
    192 forall ( otype RationalImpl | arithmetic( RationalImpl ) | { double convert( RationalImpl ); } )
    193 double widen( Rational(RationalImpl) r ) {
    194         return convert( r.numerator ) / convert( r.denominator );
    195 } // widen
    196 
    197 // http://www.ics.uci.edu/~eppstein/numth/frap.c
    198 forall ( otype RationalImpl | arithmetic( RationalImpl ) | { double convert( RationalImpl ); RationalImpl convert( double ); } )
    199 Rational(RationalImpl) narrow( double f, RationalImpl md ) {
    200         if ( md <= (RationalImpl){1} ) {                                        // maximum fractional digits too small?
    201                 return (Rational(RationalImpl)){ convert( f ), (RationalImpl){1}}; // truncate fraction
    202         } // if
    203 
    204         // continued fraction coefficients
    205         RationalImpl m00 = {1}, m11 = { 1 }, m01 = { 0 }, m10 = { 0 };
    206         RationalImpl ai, t;
    207 
    208         // find terms until denom gets too big
    209         for ( ;; ) {
    210                 ai = convert( f );
    211           if ( ! (m10 * ai + m11 <= md) ) break;
    212                 t = m00 * ai + m01;
    213                 m01 = m00;
    214                 m00 = t;
    215                 t = m10 * ai + m11;
    216                 m11 = m10;
    217                 m10 = t;
    218                 double temp = convert( ai );
    219           if ( f == temp ) break;                                                       // prevent division by zero
    220                 f = 1 / (f - temp);
    221           if ( f > (double)0x7FFFFFFF ) break;                          // representation failure
    222         } // for
    223         return (Rational(RationalImpl)){ m00, m10 };
    224 } // narrow
     192// forall ( otype RationalImpl | arithmetic( RationalImpl ) )
     193// double widen( Rational(RationalImpl) r ) {
     194//      return (double)r.numerator / (double)r.denominator;
     195// } // widen
     196
     197// // http://www.ics.uci.edu/~eppstein/numth/frap.c
     198// forall ( otype RationalImpl | arithmetic( RationalImpl ) )
     199// Rational(RationalImpl) narrow( double f, RationalImpl md ) {
     200//      if ( md <= 1 ) {                                                                        // maximum fractional digits too small?
     201//              return (Rational(RationalImpl)){ f, 1};                 // truncate fraction
     202//      } // if
     203
     204//      // continued fraction coefficients
     205//      RationalImpl m00 = 1, m11 = 1, m01 = 0, m10 = 0;
     206//      RationalImpl ai, t;
     207
     208//      // find terms until denom gets too big
     209//      for ( ;; ) {
     210//              ai = (RationalImpl)f;
     211//        if ( ! (m10 * ai + m11 <= md) ) break;
     212//              t = m00 * ai + m01;
     213//              m01 = m00;
     214//              m00 = t;
     215//              t = m10 * ai + m11;
     216//              m11 = m10;
     217//              m10 = t;
     218//              t = (double)ai;
     219//        if ( f == t ) break;                                                          // prevent division by zero
     220//        f = 1 / (f - (double)t);
     221//        if ( f > (double)0x7FFFFFFF ) break;                          // representation failure
     222//      }
     223//      return (Rational(RationalImpl)){ m00, m10 };
     224// } // narrow
    225225
    226226
Note: See TracChangeset for help on using the changeset viewer.