Changeset 3767


Ignore:
Timestamp:
Jan 12, 2015, 7:38:53 PM (7 years ago)
Author:
sam
Message:

lolremez: useless refactoring.

Location:
trunk/tools/lolremez
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/tools/lolremez/lolremez.cpp

    r3765 r3767  
    8282    }
    8383
    84     RemezSolver solver(degree + 1, 20);
    85     solver.Run(xmin, xmax, f, g);
     84    remez_solver solver(degree + 1, 20);
     85    solver.run(xmin, xmax, f, g);
    8686
    8787    return 0;
  • trunk/tools/lolremez/solver.cpp

    r3766 r3767  
    2727using lol::real;
    2828
    29 RemezSolver::RemezSolver(int order, int decimals)
     29remez_solver::remez_solver(int order, int decimals)
    3030  : m_order(order),
    3131    m_decimals(decimals),
     
    3434}
    3535
    36 void RemezSolver::Run(real a, real b, char const *func, char const *weight)
     36void remez_solver::run(real a, real b, char const *func, char const *weight)
    3737{
    3838    using std::printf;
     
    5050    m_epsilon = pow((real)10, (real)-(m_decimals + 2));
    5151
    52     Init();
    53 
    54     PrintPoly();
    55 
    56     real error = -1;
     52    remez_init();
     53    print_poly();
    5754
    5855    for (int n = 0; ; n++)
    5956    {
    60         real newerror = FindExtrema();
    61         printf(" -:- error at step %i: ", n);
    62         newerror.print(m_decimals);
    63         printf("\n");
    64 
    65         Step();
    66 
    67         if (error >= (real)0 && fabs(newerror - error) < error * m_epsilon)
     57        real old_error = m_error;
     58
     59        find_extrema();
     60        remez_step();
     61
     62        if (m_error >= (real)0
     63             && fabs(m_error - old_error) < m_error * m_epsilon)
    6864            break;
    69         error = newerror;
    70 
    71         PrintPoly();
    72 
    73         FindZeroes();
    74     }
    75 
    76     PrintPoly();
     65
     66        print_poly();
     67        find_zeroes();
     68    }
     69
     70    print_poly();
    7771}
    7872
     
    8175 * order N+1 and get a good initial polynomial estimate.
    8276 */
    83 void RemezSolver::Init()
     77void remez_solver::remez_init()
    8478{
    8579    /* m_order + 1 zeroes of the error function */
     
    9589    {
    9690        m_zeroes[i] = (real)(2 * i - m_order) / (real)(m_order + 1);
    97         fxn.Push(EvalFunc(m_zeroes[i]));
     91        fxn.Push(eval_func(m_zeroes[i]));
    9892    }
    9993
     
    127121 * of order N+2 to both refine the estimate and compute the error.
    128122 */
    129 void RemezSolver::Step()
     123void remez_solver::remez_step()
    130124{
    131125    /* Pick up x_i where error will be 0 and compute f(x_i) */
    132126    array<real> fxn;
    133127    for (int i = 0; i < m_order + 2; i++)
    134         fxn.Push(EvalFunc(m_control[i]));
     128        fxn.Push(eval_func(m_control[i]));
    135129
    136130    /* We build a matrix of Chebishev evaluations: row i contains the
     
    147141    for (int i = 0; i < m_order + 2; i++)
    148142    {
    149         real error = fabs(EvalWeight(m_control[i]));
     143        real error = fabs(eval_weight(m_control[i]));
    150144        system[i][m_order + 1] = (i & 1) ? error : -error;
    151145    }
     
    171165}
    172166
    173 void RemezSolver::FindZeroes()
     167void remez_solver::find_zeroes()
    174168{
    175169    m_stats_cheby = m_stats_func = m_stats_weight = 0.f;
     
    183177
    184178        a.value = m_control[i];
    185         a.error = EvalEstimate(a.value) - EvalFunc(a.value);
     179        a.error = eval_estimate(a.value) - eval_func(a.value);
    186180        b.value = m_control[i + 1];
    187         b.error = EvalEstimate(b.value) - EvalFunc(b.value);
     181        b.error = eval_estimate(b.value) - eval_func(b.value);
    188182
    189183        static real limit = ldexp((real)1, -500);
     
    200194             * out of this situation. */
    201195            c.value = newc == c.value ? (a.value + b.value) / 2 : newc;
    202             c.error = EvalEstimate(c.value) - EvalFunc(c.value);
     196            c.error = eval_estimate(c.value) - eval_func(c.value);
    203197
    204198            if (c.error == zero)
     
    215209    }
    216210
    217     printf(" -:- times for zeroes: estimate %f func %f weight %f\n",
     211    printf(" -:- timings for zeroes: estimate %f func %f weight %f\n",
    218212           m_stats_cheby, m_stats_func, m_stats_weight);
    219213}
    220214
    221215/* XXX: this is the real costly function */
    222 real RemezSolver::FindExtrema()
     216void remez_solver::find_extrema()
    223217{
    224218    using std::printf;
     
    227221     * compute the relative error, since its extrema are at slightly
    228222     * different locations than the absolute error’s. */
    229     real final = 0;
     223    m_error = 0;
    230224
    231225    m_stats_cheby = m_stats_func = m_stats_weight = 0.f;
     
    251245                {
    252246                    ++evals;
    253                     real error = fabs(EvalEstimate(c) - EvalFunc(c));
    254                     real weight = fabs(EvalWeight(c));
     247                    real error = fabs(eval_estimate(c) - eval_func(c));
     248                    real weight = fabs(eval_weight(c));
    255249                    /* if error/weight >= maxerror/maxweight */
    256250                    if (error * maxweight >= maxerror * weight)
     
    281275            {
    282276                real e = fabs(maxerror / maxweight);
    283                 if (e > final)
    284                     final = e;
     277                if (e > m_error)
     278                    m_error = e;
    285279                m_control[i] = (a + b) / 2;
    286280                break;
     
    289283    }
    290284
    291     printf(" -:- times for extrema: estimate %f func %f weight %f\n",
     285    printf(" -:- timings for extrema: estimate %f func %f weight %f\n",
    292286           m_stats_cheby, m_stats_func, m_stats_weight);
    293287    printf(" -:- calls: %d rounds, %d evals\n", rounds, evals);
    294 
    295     return final;
    296 }
    297 
    298 void RemezSolver::PrintPoly()
     288    printf(" -:- error: ");
     289    m_error.print(m_decimals);
     290    printf("\n");
     291}
     292
     293void remez_solver::print_poly()
    299294{
    300295    using std::printf;
     
    316311}
    317312
    318 real RemezSolver::EvalEstimate(real const &x)
     313real remez_solver::eval_estimate(real const &x)
    319314{
    320315    Timer t;
     
    324319}
    325320
    326 real RemezSolver::EvalFunc(real const &x)
     321real remez_solver::eval_func(real const &x)
    327322{
    328323    Timer t;
     
    332327}
    333328
    334 real RemezSolver::EvalWeight(real const &x)
     329real remez_solver::eval_weight(real const &x)
    335330{
    336331    Timer t;
  • trunk/tools/lolremez/solver.h

    r3764 r3767  
    1414
    1515//
    16 // The RemezSolver class
    17 // ---------------------
     16// The remez_solver class
     17// ----------------------
    1818//
    1919
     
    2222#include "expression.h"
    2323
    24 class RemezSolver
     24class remez_solver
    2525{
    2626public:
    27     RemezSolver(int order, int decimals);
     27    remez_solver(int order, int decimals);
    2828
    29     void Run(lol::real a, lol::real b,
     29    void run(lol::real a, lol::real b,
    3030             char const *func, char const *weight = nullptr);
    3131
    3232private:
    33     void Init();
    34     void FindZeroes();
    35     lol::real FindExtrema();
    36     void Step();
    37     void PrintPoly();
     33    void remez_init();
     34    void remez_step();
    3835
    39     lol::real EvalEstimate(lol::real const &x);
    40     lol::real EvalFunc(lol::real const &x);
    41     lol::real EvalWeight(lol::real const &x);
     36    void find_zeroes();
     37    void find_extrema();
     38
     39    void print_poly();
     40
     41    lol::real eval_estimate(lol::real const &x);
     42    lol::real eval_func(lol::real const &x);
     43    lol::real eval_weight(lol::real const &x);
    4244
    4345private:
     
    5355    lol::array<lol::real> m_control;
    5456
    55     lol::real m_k1, m_k2, m_epsilon;
     57    lol::real m_k1, m_k2, m_epsilon, m_error;
    5658
    5759    /* Statistics */
Note: See TracChangeset for help on using the changeset viewer.