Не совсем. Во-первых, тут либо ты находишь 2 синуса одновременно, либо синус и косинус. Забегу вперед, код быстрого, но неточного синуса svml_d_sin4_core_avx2.S. Во-вторых, присмотритесь к коду вычисления обоих функций. Он существенно разный. Значит полной парализации не получится.
Да. Вы правы. То, что это происходит на последней итерации это реально важно. Тут проще будет понять на примере. Представьте, что у Вас десятичное плавающее число с 4-мя знаками точности. Например x*(1+y). y=1.234E-2, x=5.678E-2, 1=1.000E0. Попробуйте провести операции с данными числами, не забывая округлять.
Вы абсолютно правы. Не означает. Для «cредней» точности служит параметр STD. Хотя, безусловно, статистику можно сделать лучше. Лучше статистика будет чуть дальше.
Спасибо за вопрос. Ответ достаточно простой. Потому для double(float64) просто недостижима точность quad (float128) потому, что в переменной double меньше значащих знаков. Всё, что мы можем сделать это попытаться как можно ближе подобраться к теоретическому пределу в 0.5ULP. И данное сравнение показывает кто всё-таки ближе. Написанная функция, или встроенная.
Начну со конца. Ну напишите, наконец-то полином Чебышева, который будет сходится быстрее, а расчётная сложность будет меньше. Уже написана программа для тестирования. Вам реально там остаётся добавить одну (!) функцию. Ну покажите нам всем это! Это же, я надеюсь, не секрет какой-то?
По первому вопросу DrSmile вы реально думаете, что читатели (или я) не понимают, что небольшая вариация коэффициентов в разложении может уменьшить общую ошибку? Это абсолютно банальное утверждение по вашему в статье поменяет что? Она поменяет схему оптимизации скорости? Она поменяет вывод что, сумму надо вести с мелких элементов? С другой стороны, если вы всё прекрасно знаете, то вы не можете не понимать, что количество и значение этих коэффициентов зависит от пределов вычисления синуса. Для [0, pi/2] одна, для [0, pi/256] совершенно другая. И что в самом начале это давать просто нерационально. Откройте код и посмотрите следующий метод, который я буду описывать.
Я поставил цель при написании этой статьи не самопиар, а планомерное объяснение шаблонов работы с числами с плавающей точкой и не только. Я правда призываю Вас не тратить ни свои ни мои силы на написание этих бесполезных комментариев, а взять и написать статью лучше. Получить много кармы от радостных подписчиков. Я, правда, считаю, что подбор этих коэффициентов выходит за рамки данного цикла. Как я написал в PS я только покажу подбор последнего члена. Это реально можно сделать аналитически красиво.
К сожалению я вижу в Ваших комментариях только желание принизить данную работу. У каждого, конечно, есть право на своё мнение и желания. У меня тоже есть право больше не отвечать на ваши комментарии, чем я и воспользуюсь.
Во-первых, вы путаете понятия математического разложения функции, и схемы её вычисления. Если вам так будет проще, то в строке 126 ряд Тейлора вычисляется по схеме Горнера. Абсолютно так же как мой пример sin_e4. Если сомневаетесь, то раскройте скобки в выражении и приведите получившийся полином к каноническому виду. Вы увидите, что это не что иное как запись разложения синуса в ряд Тейлора.
Во-вторых, я не заявляю, что любой полином выглядит как ряд Тейлора. Я говорю, что какие бы вы полиномы не использовали, Вам не удастся достичь вычислительной сложности (меньше сложений/умножений), лучше чем у ряда Тейлора при такой же «математической» (все операции проводятся абсолютно точно) точности при приведении конечного выражения к каноническому виду. Считаете по-другому — так приведите контрпример.
P.S. Есть один способ улучшения точности последнего члена разложения. О нём я поговорю позже.
Простыми словами, разработчики glibc очень умные люди и не стали бы они делать плохо.
Чуть более сложно: Функция fsin некоторых случаях даже точнее за счёт использования расширенных регистров. Но в некоторых случаях она даёт неприемлемый результат, и в добавок ко всему, сильно медленнее. Я очень подозреваю, что данная функция оставлена только для совместимости с ранними программами. Косвенно на это указано во втором параграфе документа.
P.S. Статья обновлена. Можно сравнивать.
Спасибо за замечательный вопрос. Конечно инструкция fsin есть и для аргумента double, но
1) У неё есть одна неприятная особенность, связанная с точностью intel fsin. Вызывать её просто так нельзя (см. график и «Conclusion» в статье). Те, кто считают CORDIS новым, современным методом, прочитайте, пожалуйста первый параграф документа по ссылке.
2) Она медленнее. На моём железе (i5 mobile) ~ в 4 раза. Через некоторое время я отредактирую исходную статью (полный текст программы) и добавлю туда ассемблерный вызов синуса. Вы сможете у себя проверит разницу в производительности.
Спасибо Артём за ваш ответ В свою очередь я тоже, конечно, не спорю, что для других функций другие методы могут быть лучше. Рад что мы нашли общий язык. Сами понимаете, что не все полезут в код этой библиотеки по вашей ссылке и не все разберутся, что там на самом деле не минмакс полином. Поэтому мне нужно было отдельно это отметить. В остальном я прочитал ваши комментарии с большим интересом и принял к сведению. Замечание, что ряд Тейлора это не универсальная вещь, конечно, полезно.
При вызове функции из библиотеки, она, конечно, исполняется в машинных кодах. Но этот год генерируется из С кода s_sin.c, так и из ассемблера svml_d_sin4_core_avx2.S. Безусловно это обсудим в следующих статьях.
Это практическая работа. Научности тут не будет. И не потому, что я не умею так писать. Умею. Я учёный, у меня больше 20 работ в иностранных журналах. Этим циклом я бы хотел дать читателям возможность руками пощупать все подводные камни. Если вам кажется, что задачи слишком простые, может в следующих статьях вы найдёте для себя что-то более сложное.
Ход мыслей у Вас очень правильный. Но утверждение 2: «происходит суммирование более-менее сопоставимых величин.» не верно. Величины 1.0 и y, абсолютно так же несапоставимы как и x и x * y, если x не 0 разумеется. Попробуете докрутить?
Ответ на Ваш вопрос не нельзя назвать простым. Ни одна функция в стандартной библиотеке не гарантирует абсолютной точности 0.5ULP. Это, как абсолютный ноль, величина не достижимая. Ну по крайней мере для сложных функций. Поэтому точность функции это всегда компромисс между её скоростью и самой точностью. Главный мой интерес, как программиста — разобраться, как устроена библиотека, повторить, а может быть в чём-то превзойти. Например, при такой же точности сделать быстрее.
Жаль Вас разочаровывать, но таких исследований в данном цикле статей не будет. Впрочем, метод, который я буду обсуждать в следующей статье может быть адаптирован для любой заданной точности.
Текущим моим комментарием я не хочу отвечать на негативные отзывы. Я хочу обратится к людям, которым статья понравилась, и, прочитав огромное количество комментариев (которые уже превышают размер статьи) могут усомниться в том, что я написал. Если вы начали сомневаться подумайте о следующих вещах.
1) «Такая точность не нужна» в различных вариациях. Безусловно она нужна не всегда. Но, если бы она никогда не была нужна, зачем вообще тогда придумали числа double? Зачем в С они по умолчанию (напр константа 1.0 это double, а 1.0f float (32-bit)). sin(x) это double, а sinf float (32-bit). Зачем nvidia добавляет в свои топовые карты поддержку double, если это ни кому не нужно?
2) «Ряд Тейлора безнадёжно устарел» в различных вариациях. Утверждение, конечно, сильное но совершенно не обоснованное. Никто (кроме ArtemKaravaev) не написал код лучше, либо не дал ссылку. А ссылка от Артёма приводит к полиному, который выглядит как ряд Тейлора и считается как ряд Тейлора. От того, что его назвать по-другому не изменится ничего.
3) «Статья ничего нового не рассказывает». Может быть, но если всё до нельзя банально, то почему никто из комментаторов не ответил на вопрос в конце? Почему никто не объяснил в двух словах почему у sin_e4 и sin_e5 результаты разные? И так, чтобы всем стало понятно.
Хочу сказать спасибо и обратиться к тем, кому понравилась эта статья. Если у Вас есть желание меня поддержать, 10 минут времени и базовые знания программирования, то я буду очень благодарен, если вы перейдёте по ссылке от Артёма (sinq_kernel.c) и напишите в комментариях является ли формула в строке 126 одинаковой с второй формулой в моей статье или нет.
По первому вопросу DrSmile вы реально думаете, что читатели (или я) не понимают, что небольшая вариация коэффициентов в разложении может уменьшить общую ошибку? Это абсолютно банальное утверждение по вашему в статье поменяет что? Она поменяет схему оптимизации скорости? Она поменяет вывод что, сумму надо вести с мелких элементов? С другой стороны, если вы всё прекрасно знаете, то вы не можете не понимать, что количество и значение этих коэффициентов зависит от пределов вычисления синуса. Для [0, pi/2] одна, для [0, pi/256] совершенно другая. И что в самом начале это давать просто нерационально. Откройте код и посмотрите следующий метод, который я буду описывать.
Я поставил цель при написании этой статьи не самопиар, а планомерное объяснение шаблонов работы с числами с плавающей точкой и не только. Я правда призываю Вас не тратить ни свои ни мои силы на написание этих бесполезных комментариев, а взять и написать статью лучше. Получить много кармы от радостных подписчиков. Я, правда, считаю, что подбор этих коэффициентов выходит за рамки данного цикла. Как я написал в PS я только покажу подбор последнего члена. Это реально можно сделать аналитически красиво.
К сожалению я вижу в Ваших комментариях только желание принизить данную работу. У каждого, конечно, есть право на своё мнение и желания. У меня тоже есть право больше не отвечать на ваши комментарии, чем я и воспользуюсь.
Во-вторых, я не заявляю, что любой полином выглядит как ряд Тейлора. Я говорю, что какие бы вы полиномы не использовали, Вам не удастся достичь вычислительной сложности (меньше сложений/умножений), лучше чем у ряда Тейлора при такой же «математической» (все операции проводятся абсолютно точно) точности при приведении конечного выражения к каноническому виду. Считаете по-другому — так приведите контрпример.
P.S. Есть один способ улучшения точности последнего члена разложения. О нём я поговорю позже.
Чуть более сложно: Функция fsin некоторых случаях даже точнее за счёт использования расширенных регистров. Но в некоторых случаях она даёт неприемлемый результат, и в добавок ко всему, сильно медленнее. Я очень подозреваю, что данная функция оставлена только для совместимости с ранними программами. Косвенно на это указано во втором параграфе документа.
P.S. Статья обновлена. Можно сравнивать.
1) У неё есть одна неприятная особенность, связанная с точностью intel fsin. Вызывать её просто так нельзя (см. график и «Conclusion» в статье). Те, кто считают CORDIS новым, современным методом, прочитайте, пожалуйста первый параграф документа по ссылке.
2) Она медленнее. На моём железе (i5 mobile) ~ в 4 раза. Через некоторое время я отредактирую исходную статью (полный текст программы) и добавлю туда ассемблерный вызов синуса. Вы сможете у себя проверит разницу в производительности.
1) «Такая точность не нужна» в различных вариациях. Безусловно она нужна не всегда. Но, если бы она никогда не была нужна, зачем вообще тогда придумали числа double? Зачем в С они по умолчанию (напр константа 1.0 это double, а 1.0f float (32-bit)). sin(x) это double, а sinf float (32-bit). Зачем nvidia добавляет в свои топовые карты поддержку double, если это ни кому не нужно?
2) «Ряд Тейлора безнадёжно устарел» в различных вариациях. Утверждение, конечно, сильное но совершенно не обоснованное. Никто (кроме ArtemKaravaev) не написал код лучше, либо не дал ссылку. А ссылка от Артёма приводит к полиному, который выглядит как ряд Тейлора и считается как ряд Тейлора. От того, что его назвать по-другому не изменится ничего.
3) «Статья ничего нового не рассказывает». Может быть, но если всё до нельзя банально, то почему никто из комментаторов не ответил на вопрос в конце? Почему никто не объяснил в двух словах почему у sin_e4 и sin_e5 результаты разные? И так, чтобы всем стало понятно.
Хочу сказать спасибо и обратиться к тем, кому понравилась эта статья. Если у Вас есть желание меня поддержать, 10 минут времени и базовые знания программирования, то я буду очень благодарен, если вы перейдёте по ссылке от Артёма (sinq_kernel.c) и напишите в комментариях является ли формула в строке 126 одинаковой с второй формулой в моей статье или нет.