Fast Hypotenuse Algorithm for Embedded Processor?

Если результат не должен быть особенно точным, вы можете получить сырой аппроксимация довольно проста:

Возьмите абсолютные значения a и b и поменяйте местами, если необходимо, так, чтобы у вас было a <= b . Затем:

h = ((sqrt(2) - 1) * a) + b

Чтобы интуитивно понять, как это работает, рассмотрим способ, которым пологая наклонная линия отображается на пиксельном дисплее (например, с использованием алгоритма Брезенхема). Это выглядит примерно так:

+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
| | | | | | | | | | | | | | | | |*|*|*|    ^
+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+    |
| | | | | | | | | | | | |*|*|*|*| | | |    |
+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+    |
| | | | | | | | |*|*|*|*| | | | | | | | a pixels
+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+    |
| | | | |*|*|*|*| | | | | | | | | | | |    |
+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+    |
|*|*|*|*| | | | | | | | | | | | | | | |    v
+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
 <-------------- b pixels ----------->

Для каждого шага в направлении b следующий пиксель, который должен быть нанесен, находится либо сразу вправо, либо на один пиксель вверх и вправо.

Идеальная линия от одного конца до другого может быть аппроксимирована путем, соединяющим центр каждого пикселя с центром соседнего. Это серия из a сегментов длиной sqrt (2) и b-a сегментов длиной 1 (если в качестве единицы измерения используется пиксель). Отсюда вышеприведенная формула.

Это явно дает точный ответ для a == 0 и a == b ; но дает завышенную оценку промежуточных значений.

Ошибка зависит от отношения b / a ; максимальная ошибка возникает, когда b = (1 + sqrt (2)) * a и оказывается 2 / sqrt (2 + sqrt (2)) , или около 8,24% выше истинного значения. Это не очень хорошо, но если этого достаточно для вашего приложения, то преимущество этого метода состоит в том, что он простой и быстрый.(Умножение на константу можно записать как последовательность сдвигов и сложений.)

Одна возможность выглядит так:

#include <math.h>

/* Iterations   Accuracy
 *  2          6.5 digits
 *  3           20 digits
 *  4           62 digits
 * assuming a numeric type able to maintain that degree of accuracy in
 * the individual operations.
 */
#define ITER 3

double dist(double P, double Q) {
/* A reasonably robust method of calculating `sqrt(P*P + Q*Q)'
 *
 * Transliterated from _More Programming Pearls, Confessions of a Coder_
 * by Jon Bentley, pg. 156.
 */

    double R;
    int i;

    P = fabs(P);
    Q = fabs(Q);

    if (P<Q) {
        R = P;
        P = Q;
        Q = R;
    }

/* The book has this as:
 *  if P = 0.0 return Q; # in AWK
 * However, this makes no sense to me - we've just insured that P>=Q, so
 * P==0 only if Q==0;  OTOH, if Q==0, then distance == P...
 */
    if ( Q == 0.0 )
        return P;

    for (i=0;i<ITER;i++) {
        R = Q / P;
        R = R * R;
        R = R / (4.0 + R);
        P = P + 2.0 * R * P;
        Q = Q * R;
    }
    return P;
}

Это все еще делает пару делений и четыре кратных за итерацию, но вам редко требуется больше трех итераций (а двух часто бывает достаточно) на ввод. По крайней мере, с большинством процессоров, которые я видел, это обычно будет быстрее, чем sqrt сам по себе.

На данный момент он написан для double s, но, если предположить, что вы реализовали основные операции, преобразование его для работы с фиксированной точкой не должно быть очень сложным.

Возможно, вы могли бы использовать некоторые из библиотек ассемблера Elm Chans и адаптировать функцию ihypot к вашему ATtiny. Вам нужно будет заменить MUL и, возможно (я не проверял) какие-то другие инструкции.

Рассмотрите возможность использования методов CORDIC. У доктора Добба есть статья и связанный с ней источник библиотеки здесь . Квадратный корень, умножение и деление рассматриваются в конце статьи.

Вы можете начать с переоценки, нужен ли вам вообще sqrt. Много раз вы вычисляете гипотенузу только для того, чтобы сравнить ее с другим значением - если вы возведете в квадрат значение, с которым сравниваете, вы можете полностью исключить квадратный корень.