Finding a rational approximation of a real number when both the numerator and denominator have a res-CodePudding

I'm trying to learn how to find a rational approximation of a real number when both the numerator and denominator are constrained. I've looked at many pages now, including the following two, and learned about the continued fractions method, the Farey sequence, and the Stern-Brocot tree. However, in all of these examples, either the numerator or denominator are unconstrained.

Algorithm for simplifying decimal to fractions
https://gist.github.com/mikeando/7073d62385a34a61a6f7

Here's my situation:

I'm testing mixed-signal ICs.
In one of our tests, in order to find the IC's maximum operating frequency, the clock signal going into the IC is set to 12 MHz and continually decreased until the IC is able to run a simple digital sequence.
The test platform's main clock has a range of 25 to 66 MHz, and the function that sets it takes in a double.
In the current version of the test, it's set to a constant 50.0 MHz, then a function that divides that frequency is called in a loop. The divisor is an integer that can be anywhere between 1 and 4096.

However, this results in imprecise measurements.
The devices always pass at:
50 / 5 = 10 Mhz
50 / 6 = 8.333 MHz

If possible, to get more precise measurements, I'd like to be able to change the main clock's frequency AND the clock divisor in each loop iteration. That's why I'm trying to learn how to write something like a continued fractions algorithm with constraints to both the numerator and denominator. I'm envisioning something like this:

while(dFmax > dFmin)
{
    std::pair<double, int> bestSettings = GetBestClkSettings(dFmax);
    double dFreq = bestSettings.first;
    int iDiv = bestSettings.second;

    // Set up clock and timesets        
    clkset(dFreq);
    clkdivide(iDiv);

    // Run pattern
    // ...
    
    // Get results
    // ...
    
    dFmax -= 0.1;
}

Not only have I spent hours experimenting with the code in the second link, but I also tried writing a function that uses something like binary searching just to see what would happen. I'm fully aware that this is terrible code that I can't use to achieve my goals; I just wanted to show that I have been trying.

#include <iostream>
#include <stdio.h>
#include <cmath>


// The fraction struct and the main() function were largely taken from:
// https://gist.github.com/mikeando/7073d62385a34a61a6f7

struct fraction {
    int n;
    int d;
    
    fraction()
    {
        this->n = -1;
        this->d = -1;
    }
    fraction(int n, int d)
    {
        this->n = n;
        this->d = d;
    }
    
    double asDouble()
    {
        double dReal = static_cast<double>(n) / static_cast<double>(d);
        return dReal;
    }
};



fraction ApproximateFrequency(double dTargetFreqMHz, double dTol)
{
    fraction result;
    
    if (dTargetFreqMHz < (25.0 / 4096) || dTargetFreqMHz > 66.0)
    {
        return result;
    }
    else if (dTargetFreqMHz >= 25.0 && dTargetFreqMHz <= 66.0)
    {
        result.n = static_cast<int>(dTargetFreqMHz);
        result.d = 1;
        return result;
    }

    int iFrqLo = 25;
    int iFrqHi = 66;
    int iDivLo = 1;
    int iDivHi = 4096;

    int iFrqCurr = (iFrqLo   iFrqHi) / 2;
    int iDivCurr = (iDivLo   iDivHi) / 2;
    double dFreq = static_cast<double>(iFrqCurr) / static_cast<double>(iDivCurr);
    double dPrevFreq = 0;
    int iNumIters = 1;
    
    while (fabs(dTargetFreqMHz - dFreq) > dTol && fabs(dFreq - dPrevFreq) > 1e-8 && iNumIters < 25)
    {
        dPrevFreq = dFreq;
        
        if (dFreq < dTargetFreqMHz)
        {
            // The frequency needs to be increased.
            
            // The clock frequency could be increased:
            int iFrqNew = (iFrqCurr   iFrqHi) / 2;
            double dFrqIfClkInc = static_cast<double>(iFrqNew) / static_cast<double>(iDivCurr);
            double dClkIncDiff = fabs(dTargetFreqMHz - dFrqIfClkInc);
            
            // Or the divider could be decreased:
            int iDivNew = (iDivLo   iDivCurr) / 2;
            double dFrqIfDivDec = static_cast<double>(iFrqCurr) / static_cast<double>(iDivNew);
            double dDivDecDiff = fabs(dTargetFreqMHz - dFrqIfDivDec);
            
            // Find the option that produces a better result:
            if (dClkIncDiff < dDivDecDiff && iFrqNew >= 25 && iFrqNew <= 66)
            {
                iFrqCurr = iFrqNew;
            }
            else if (dDivDecDiff < dClkIncDiff && iDivNew >= 1 && iDivNew <= 4096)
            {
                iDivCurr = iDivNew;
            }
        }
        else
        {
            // The frequency needs to be decreased.
            
            // The clock frequency could be decreased:
            int iFrqNew = (iFrqLo   iFrqCurr) / 2;
            double dFrqIfClkDec = static_cast<double>(iFrqNew) / static_cast<double>(iDivCurr);
            double dClkDecDiff = fabs(dTargetFreqMHz - dFrqIfClkDec);
            
            // Or the divider could be increased:
            int iDivNew = (iDivCurr   iDivHi) / 2;
            double dFrqIfDivInc = static_cast<double>(iFrqCurr) / static_cast<double>(iDivNew);
            double dDivIncDiff = fabs(dTargetFreqMHz - dFrqIfDivInc);
            
            // Find the option that produces a better result:
            if (dClkDecDiff < dDivIncDiff && iFrqNew >= 25 && iFrqNew <= 66)
            {
                iFrqCurr = iFrqNew;
            }
            else if (dDivIncDiff < dClkDecDiff && iDivNew >= 1 && iDivNew <= 4096)
            {
                iDivCurr = iDivNew;
            }
        }
        
        // See the frequency attainable with the current settings
        dFreq = static_cast<double>(iFrqCurr) / static_cast<double>(iDivCurr);
        std::cout << "prev = " << dPrevFreq << ", current = " << dFreq << std::endl;
        iNumIters  ;
    }
    
    result.n = iFrqCurr;
    result.d = iDivCurr;
    return result;
}



int main(int argc, char* argv[])
{
    double dTargetFreqMHz = 20.0;
    std::cout << "Target: " << dTargetFreqMHz << "\n\n";
    double dTol = 0.05;
    fraction mcf = ApproximateFrequency(dTargetFreqMHz, dTol);
    printf("tol=%f, n/d = %d/%d = %f (err=%f)\n", dTol, mcf.n, mcf.d, mcf.asDouble(), mcf.asDouble()-dTargetFreqMHz);
}

Any advice or hints would be greatly appreciated. Thank you in advance.

CodePudding user response：

Since your range is so constrained, you could just brute force it. There's only 172,032 possible combinations of numerator and denominator to check. The algorithm could be made more efficient by iterating from 25 to 66 and calculating the closest two denominators, in which case you only have to check 84 possibilities:

fraction ApproximateFrequency(double dTargetFreqMHz, double dTol)
{
    fraction result;

    if (dTargetFreqMHz < (25.0 / 4096) || dTargetFreqMHz > 66.0)
    {
        return result;
    }
    else if (dTargetFreqMHz >= 25.0 && dTargetFreqMHz <= 66.0)
    {
        result.n = static_cast<int>(dTargetFreqMHz);
        result.d = 1;
        return result;
    }

    double smallestError = 66.0;
    int closestNum = 0;
    int closestDenom = 0;
    for (int num = 25; num <= 66; num  )
    {
        int denom = floor((double)num / dTargetFreqMHz);
        if (denom >= 1 && denom <= 4096)
        {
            double freq = (double)num / double(denom);
            double err = fabs(dTargetFreqMHz - freq);

            if (err < smallestError)
            {
                closestNum = num;
                closestDenom = denom;
                smallestError = err;
            }
            if (denom <= 4095)
            {
                freq = (double)num / double(denom   1);
                err = fabs(dTargetFreqMHz - freq);
                if (err < smallestError)
                {
                    closestNum = num;
                    closestDenom = denom   1;
                    smallestError = err;
                }
            }
        }
    }

    result.n = closestNum;
    result.d = closestDenom;

    return result;
}

the dTol parameter isn't used so you could get rid of it.