Мы пойдем мимо — и дальше




В своем предыдущем посте я показал, как можно улучшить быстродействие расчета точек на кривой Безье (КБ) путем:

  1. Преобразования расчетных формул — ускорение в ~3 раза.
  2. Перехода от чисел ПТ к ФТ — ускорения почти нет, но позволяет провести 3.
  3. Заменой операции деления умножением и сдвигом — ускорение еще на 40%.

Печальное отступление
— я допустил неточность в последней формуле можно было еще чуть ускорить вычисления, свернув еще одно константной выражение и, исключив умножение, вместо 502 получить 410 тактов на цикл вычисления. К сожалению, никто из читателей предыдущего поста мне на это не указал в комментариях… а я на это надеялся, значит, я не смог достаточно заинтересовать своих читателей, чтобы они правильно (то есть внимательно) читали мои опусы. Ладно, попробуем еще разок.


В конце упомянутого поста я заявил, что можно вычисления сделать еще быстрее, ну что же «пацан сказал — пацан сделал».

Напомню еще раз полученную формулу для вычисления точки на КБ

$Р=(А1*и+Б1)*и+С (=>2*2+)$

.Очередное увеличение скорости связано с особенностью задачи — нам нужно не И раз рассчитать значение КБ при различном значении параметра «и», а найти ряд значений при изменении (в данном случае увеличении) данного параметра на известную, более того, фиксированную величину (в нашем случае единицу), что позволяет воспользоваться описанным далее приемом. В мое время это называлось разностным методом вычисления (если память мне не изменяет, по крайней мере, так его называли на лекциях), в Вашем распоряжении весь Инет, может быть (даже наверняка), есть и другое название.

Рассмотрим случай Р=А*и (=>1*), и предположим, что мы знаем значение Р0 при некотором аргументе и0. Тогда значение при следующем аргументе и0+1 можно рассчитать, как Р=А*(и0+1)=А*и0+А=Р0+А (=>1+). Нисколько не потеряв в точности, мы смогли заменить операцию умножения на операцию сложения, которая намного быстрее.

Теперь более сложный пример Р=А*и*и (=>2*), рассматриваем по аналогии Р=А*(и+1)*(и+1)=А*(и*и+2*и+1)=А*и*и+2*А*и+А=Р0+2*А*и+А(=>2*2+). На первый взгляд, мы ничего не выиграли, но, если заранее вычислить А'=2*А, то получим (=>1*2+), вполне возможен и выигрыш. Но мы не остановимся на достигнутом и к полученному выражению А'*и применим уже известную нам методику, тогда получим две операции над двумя переменными: Р=Р+А'+А; А'=А'+А (=>3+). Если мы еще примем во внимание, что начальное значение А' можно сделать больше на А, то вообще остается только два сложения вместо двух умножений. Да, нам пришлось завести две дополнительные переменные, но это класический размен — расплачиваемся памятью за время.

Отается только вычислить правильные начальные значения, но это делается только один раз в начале итераций, а если начальное значение параметра и=0, то и вообще тривиально Р=0, А'=А. Проверим наш метод на нескольких итерация (это совершенно излишне, правильно примененная математика не может дать неправильного ответа), но позволит лучше понять происходящее. Итак
итерация 0: и=0; Р=0 (верно); А'=А; А''=2*А;
итерация 1: и=1; Р=0+А'=0+А=А (верно); А'=А'+А''=А+2*А=3*А;
итерация 2: и=2; Р=А+3*А=4*А (верно); А'=3*А+2*А=5*А;
итерация 3: и=3; Р=9*А (верно); А'=7*А и так далее.

Отметим связь полученных на формул с методом Ньютона для вычисления значения функции в окрестностях некоторой точки с известным значением. Более того, учитывая, что наша функция имеет второй порядок и все производные, начиная с третьей, равны нулю, мы можем смело заменить знак приблизительного равенства на точное. Единственное отличие от данного метода в том, что мы постоянно перемещаем точку, относительно которой строим новую окрестность, что бы избежать выполнения операций умножения в оригинальной формулировке.

Перепишем наше исходное выражение для КБ в развернутой форме

$ Р=А1*и*и+Б1*и+С$

, тогда для вычисления с применением данного метода нам потребуется 2+ для первого слагаемого (и две переменные), 1+ для второго (и одна переменная) и 2+, чтобы сложить все вместе (=>5+). Можно ожидать, что даже такой (неправильный) подход даст выигрыш по сравнению с 2*2+, но все намного лучше. Очевидно, что операция сложения является аддитивной (спасибо, капитан), поэтому можно сгруппировать слагаемые и заменить константные члены на новые выражения, что дает следующий алгоритм:

1. начальные значения Р=С; А'=А1+Б1; А''=2*А1;
2. на очередном шаге Р=Р+А'; А'=А'+А'' (=>2+), что несомненно быстрее схемы Горнера.

Реализуем наш алгоритм в виде программы (это не столь тривиально, как может показаться, поскольку я для упрощения забыл о необходимость сдвигов, но и не слишком сложно… минут на 20), получаем вычислительную сложность (=>2+1>>) и замеряем быстродействие — получилось 140 (увеличение быстродействия еще в 3 раза) тактов на цикл, а вот это уже почти идеальный результат. Что нам осталось сделать, чтобы получить идеальный (в некотором смысле) вариант — обратить внимание на размерность операндов в формулах. Сейчас мы проводим все операции над длинными (32х разрядными) целыми, а это в некоторых местах может быть излишне. Если сократить разрядность операндов до минимально необходимой, то можем получить выигрыш процентов на 20-25, но это потребует от нас перехода к ассемблеру (язык С не располагает подходящими для таких операций средствами) и внимательного анализа данных исходной задачи. Делать ли это — решать читателю, мы и так уже ускорили вычисления более, чем в 1900/140=13 раз по сравнению с «наивным» подходом.

Последнее замечание к реализации алгоритма — мы по прежнему исключаем операцию деления заменой ее умножением на константу, что учитываем при формировании исходных данных и сдвигом на константу, кратную 8. Это можно делать различными способами с несколько отличающимся временем исполнения, оставляю подобные эксперименты на долю пытливого читателя.

Однако совершенно неожиданно возникла одна проблема — мы ясно видим две операции сложения над 32 битовыми числами, которые должны занять 4+4=8 тактов, положим еще 8*3*2=48 тактов на пересылки операндов, 4 такта на получение результата сдвига, 4 такта на вызов процедуры вычисления и возврат и 4 такта на организацию цикла — откуда еще 60 тактов, непонятно. Раньше мы этого просто не замечали на фоне сотен тактов вычисления, но теперь то можно и обратить внимание. Избыточные такты легко нашлись — в цикле оставались ненужные отладочные операции, если все аккуратно подчистить, то остается только 120 тактов и назначение каждого вполне понятно (ну не так, что бы уж совсем понятно, но все таки). Далее выясним время выполнения пустого цикла — 22 такта, интересно, чем они там занимаются все это время, ну да ладно, и очистим собственно время вычисления, которое составит 98 тактов. Обратим внимание, что изменяется оценка полученного ускорения работы: вместо 1900/140=13 получаем (1900-50)/(140-40)=19 раз, что не имеет никакого практического смысла, поскольку цикл все равно необходим, но позволяет еще больше поднять самооценку.

И последнее замечание — как нетрудно видеть, блох мы начали искать и устранять только тогда, когда расправились с крупными жуками-оленями и их (блох) существование стало очевидно и, более того, значимо. Настоятельно рекомендую именно подобный подход (и в этой рекомендации я не одинок) при оптимизации программ в плане быстродействия.

Ну а в заключение о примечании «в некотором смысле» — если речь идет о последовательном вычислении координат очередной точки на КБ при изменении параметра (олицетворяющего собой время), то предложенный алгоритм (после всех оптимизаций) улучшить уже не удается. Но если задачу переформулировать и, к примеру, задаться целью просто построить КБ на экране (без привязки ко времени), то тут есть весьма перспективный метод, ключевое слово «Брезенхайм», но «это уже совсем другая история», которой я посвящу отдельный пост (может быть, когда-нибудь, если старшая сестра не будет против).

Комментарии (19)


  1. oam2oam
    09.11.2018 12:51
    +1

    еще раз, видимо, придется написать — совершенно не нужно для вычисления точек БК ни умножений, ни делений! Вот псевдокод ниже для рекурсивной функции (обычно n=8 дает 256 точек) A,B,C — параметры уравнения A*t*t+2*B*t*q+C*q*q (t+q=1):

    b(n,A,B,C) {
      if(n>=8) plot(B);
      else {
         AC = (A+C)>>2;
         BC = (B+C)>>2;
        ABC = (AC+BC)>>2;
         b(n+1,A,AC,ABC);
         b(n+1,ABC,BC,C);
      }
    }
    

    И, что характерно, размерность пространства — любая…
    Но на практике, я повторюсь, никогда не видел, чтобы использовали второй порядок, всегда используют третий ( а особо знающие — 27-ой)


    1. GarryC Автор
      09.11.2018 13:54

      Вообще то у меня не осталось ни умножений, ни делений, а только суммы и сдвиги на 16 — это и была цель, и она достигнута. Расчеты осуществляются по рекуррентным соотношениям, которые с точки зрения производительности предпочтительнее рекурсии.

      Если Вы действительно считаете, что вычисления по Вашей методике будут осуществляться быстрее (что в принципе не исключено, поскольку у вас 3 сложения 16 битных чисел против 2 сложений 32х у меня, но 3 сдвига 16х на 1 бит выполняются дольше, чем один на 16 бит для 32х, надо аккуратно считать), то попробуйте программку с этой методикой расчета том же сайте и посмотрим на результаты. Мне кажется, что изменение порядка операндов в рекурсивных вызовах сожрет возможную экономию и даже больше.

      Но главное — даже не в этом, а в том, что Вы неправильно трактуете задачу — для построения графиков есть еще более быстрые схемы, но нам нужно именно последовательное прохождение по кривой.

      Ну а особо знающие с 27 порядком пусть попытаются объяснить менее знающим, какие именно преимущества дает увеличение порядка КБ выше второго.

      PS. А что такое plot(B) — или В — это пара координат?


      1. oam2oam
        09.11.2018 14:08

        Ну, в реальности, конечно, сдвиги не вычисляются до конца рекурсии — а потом сразу делается нужный сдвиг… Для построения графиков более быстрой схемы нет — приведенная обычно и реализуется в железе (в GPU). А насчет второго порядка — тут все просто… обычно нужно или вычислять с заданной точностью (это не для построения кривых, конечно) или нужно иметь хороший вывод последовательных точек — например, чтобы там, где кривая сильно изгибается, точки были ближе друг к другу — так вот, приведенный алгоритм именно так и делает. И еще, если речь идет о кривых Безье, то для задания куска кривой обычно задают 4 точки (крайние и 2 градиента), вот они-то и являются параметрами для аппроксимации 3 порядка (уж я не буду приводить уравнение). Так что я не очень понимаю, зачем брать 2-ой порядок — чисто теоретически?

        Если же нужно последовательно вычислять точки с заданным шагом (по параметру), то тут есть совсем другой алгоритм… Но, впрочем, я с ним работал только в случае размерности точек более 7 и в численных вычислениях… он, к сожалению, до сих пор непубличен…

        Plot(B) — тут имеется в виду, что B — это координаты, в случае плоскости B=(x,y), конечно…


  1. REPISOT
    09.11.2018 14:08

    Когда я вижу заголовок

    быстродействии Ардуино
    я жду сравнения времени выполнения задачи на ардуиноIDE vs AtmelStudio, CVAVR или вообще WMlab. Разочарован.


    1. GarryC Автор
      09.11.2018 14:37

      Кроме упомянутого мною сайта с онлайн-эмуляцией программ в Ардуино среде, у меня других нет, так что извиняюсь.


      1. REPISOT
        09.11.2018 17:53

        Зачем вообще нужен эмулятор ардуины!?


        1. tvr
          09.11.2018 18:19

          Чтобы эмулировать ардуину?


  1. Grox
    09.11.2018 14:24

    Кроме спортивного интереса, какая необходимость пытаться выжимать из платформы Ардуино максимум? Эти 8-битные микроконтроллеры стоят дороже 32-битных STM, на которых надо решать куда более серьёзные задачи, чтобы упереться в производительность. И это выжимание делает код более сложным и менее поддерживаемым, а разработку длительнее.


    1. GarryC Автор
      09.11.2018 14:41

      А методический аспект? Неужели не важно научиться оптимизировать свои проекты хотя бы на том элементарном уровне, что я описываю в своих постах?

      Вчера прочитал новость о 8битных МК от Тайваньского производителя по 5 центов за штуку — и как Вы с ними будете работать, освоив только 32х разрядные и пренебрегая оптимизацией?


      1. Grox
        09.11.2018 14:47

        Видел у EEVBlog. Они даже по 3 цента. А в неочень крупных партиях даже дешевле.
        Если работать с ними, то такие оптимизации действительно становятся важными. Но это условие, которые может стоять только у некоторых разработчиков. И тут понятно зачем — абсолютное превосходство по стоимости.
        А вот делать это на 8-битных микроконтроллерах, которые дороже, чем 32-битные, я вижу смысл только в образовательных целях или, когда, по каким-то причинам мы работаем с легаси или жёсткими требованиями заказчика.


        1. GarryC Автор
          09.11.2018 14:57

          Ну я заказал себе 10 штук за 0,57 $, поэтому такую цену и указываю. И неужели такое ценовое преимущество можно игнорировать?

          Конечно же, в образовательных, я для этого и пишу.


          1. Grox
            09.11.2018 15:03

            Заказали Padauk? И, есть ли у вас программатор к ним?

            И неужели такое ценовое преимущество можно игнорировать?
            И да, и нет, зависит от ваших задач. Они одноразово программируемые. Причём, можно заказать уже с вашей прошивкой. В крупных партиях это очень выгодное решение.


            1. GarryC Автор
              09.11.2018 15:34

              Я надеюсь на реверс-инжиниринг алгоритма программирования, этим вроде уже занимаются.
              Конкретной задачи у меня для них нет, но жалко упустить такой случай — 10 МК за 0.57$, чисто для развлечения. Если бы была задача, то программатор за 100$ отбивается на 500 (совсем небольшой тираж) экземплярах (0,25-0,05)*500=100


        1. emmibox
          09.11.2018 19:54

          Вы не учитываете, что во время когда 32-х битных контроллеров от ST физически не существовало а их ARM аналоги стоили не в три раза дешевле а в 10 раз дороже на 8-ми битных AVR были созданы некие устройства (и их не мало) — написан некий рабочий код (его тоже не мало) все это было продано (много и успешно). И задачи их развития поддержки и расширения функционала (без смены процессора и выпуска новых устройств) за это время никуда не делись. И решатся они могут зачастую только лишь за счет оптимизации тем или иным способом.

          Естественно в новых разработках это никому не нужно.


          1. Grox
            09.11.2018 20:09

            Я процитирую себя

            или, когда, по каким-то причинам мы работаем с легаси или жёсткими требованиями заказчика.


            1. emmibox
              09.11.2018 22:34

              «Какие то причины» становятся вполне основными методами, если вы не вчера в этой индустрии и число проектов не 10…


  1. Moxa
    09.11.2018 16:05

    о, кривые безье, сорри что слегка не по теме, но не подскажите как сделать гладкую непрерывныю кривую из двух кривых безье?


  1. GarryC Автор
    09.11.2018 16:55

    В общем виде эта задача не решается, не все заданные КБ можно подшить одну с другой в принципе.
    Если стоит задача в конце первой КБ гладко начать вторую, которая прибежит в новую конечную и надо найти для второй мнимый фокус, то выражаем производную второй КБ в точке начала через неизвестные координаты фокуса, приравниваем производной первой КБ в конце и решаем — должно получиться линейное уравнение, но таких точек буде бесконечно много на прямой — какую выбрать? Или я не понял вопроса?
    Может, проще сразу построить КБ более высокого порядка по заданным точкам?


    1. Moxa
      09.11.2018 19:04

      задача продолжить кривую до следующей точки, откуда, возможно, продолжать дальше… таких отрезков может быть сотни, хочется сделать плавное соединение