P.S. По моим прикидкам следующего члена будет достаточно (cs7 и cs8) соответственно. Вопрос, разумеется остаётся к качеству самой таблички. Я, если честно, не проверял.
P.P.S. Проверьте существующие значения cs и sn. Они могут отличатся от ряда Тейлора. Если так, то рассчитайте их сами точно.
Если Вам интересна практическая реализация то почему бы Вам не адаптировать s_sin? Вам нужна большая точность? Если так то самым простым образом можно решить эту проблему увеличением числа членов в разложении синуса и косинуса:
...
s = x + (dx + x * xx * (sn3 + xx * sn5)); // здесь
c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6)); // и здесь
SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs);
...
1) А если в процессоре нет аппаратной поддержки FP арифметики? Почему вы говорите, что
но это в 10 раз медленнее чем double*double.
Если у Вас есть какой-то конкретный процессор возьмите его и сделайте эксперимент.
2) a) Есть такая суперштука как FMA. С помощью этих комманд можно проводить точное умножение.
b) Если мантисса заполнена на половину, то результат умножения будет точным. Поэтому использование битовых масок решит проблему. Правда понадобится много комманд.
Этот результат по-прежнему не опровергает заявленных авторами posit свойств — превосходящую double точность.
В вашей статье Вы работаете с 32 битным POSIT и сравниваете его с 64 битным double? И делаете вывод, что POSIT менее точный? Очевидно, что даже если все 32 бита занять мантиссой их всё равно будет меньше, чем переменной double. Приведите, плз, ссылку, где создатели posit утверждают, что posit32 точнее чем double.
И не опровергает катастрофическое снижение точности при вычислении сильно разных порядков.
Просто по структуре переменной понятно, что она будет плохо работать на больших порядках. Зато сильно лучше около 1. Это же очевидно.
Добавлю краткое резюме. Если кто-то прочитал эту статью и подумал, что posit не подходит для реального вычисления синуса, то это абсолютно не так. Вот таким простым не оптимизированным кодом кодом можно посчитать синус в posit точнее, чем во float для приведённых автором значений.
template< typename T >
T sin_smp(T x) {
T xx = x * x;
T res = 1.0;
for(int i = 31; i >= 3; i -= 2) {
res = T(1.0) - res * xx / T(float(i * (i - 1)));
}
return x * res;
}
Данный код вы можете добавить в проект, описанный в статье и запустить. Результат такой
В своих экспериментах я не увидел качественного увеличения точности при перетасовке аргументов при сложении. Характер искажений — меняется, но незначительно, и чуть большая точность на одних входных данных компенсируется чуть меньшей на других.
Данное утверждение идёт в разрез не только с моими данными, но и данными этой статьи. Поэтому оно требует большего обоснования с Вашей стороны.
Вот в этой статье я считал синус схемой Горнера и получил вполне неплохие результаты — при том, что решал задачу «в лоб» и оптимальные коэффициенты никак не подбирал.
В чём «неплохость» Вашего результата? Точность, как видно нет, скорость тоже нет.
Если бы вы действительно внимательно прочитали мою статью, то обратили бы внимание, что там речь идёт о вычислении полинома, а не о вычислении синуса, и рассматривается не абсолютная, а накопленная погрешность вычислений, полученная при использовании различных типов данных.
Цитирую Вашу статью:
Используя Wolfram Mathematica и команду PadeApproximant[Sin[x], {x, 0, {11, 11}}]
получим вот такой вот рациональный полином для аппроксимации синуса, обеспечивающий точность double в диапазоне примерно от -2 до 2:
Т.е. сейчас вы говорите, что там просто какой-то полином. А зачем тогда вы его сравниваете с синусом? Что тогда с выводом статьи? Как абстрактный полином является реальным вычислением? Цитирую Вашу статью:
Преимущество Posit, демонстрируемое на единичных вычислениях — не более чем фокус, расплатой за которое является катастрофическое снижение точности в «тяжёлых» реальных вычислениях.
Я очень внимательно прочитал Вашу статью. Вывод из неё один. Вы нашли полином, который плохо аппроксимирует синус в double и очень плохо в posit. Никаких реальных вычислений, о которых Вы говорите в выводе статьи в самой статье не показано. Реальные вычисления приведены ниже.
Если вы помните, то мои комментарии были посвящены вычислению синуса с 30-ю знаками — то есть в 2 раза больше, чем у Вас.
Если я правильно помню, то Вы не хотите поделится с нами вашими результатами? Я не собираюсь считать в double-double. У меня нет такой цели. Всё, что я увидел из Вашей статьи, что вы не умете грамотно считать ни в double ни в posit. Добавьте в свой код вот эту функцию и вы увидите, что на некоторых значениях (например на 0.5) posit на 1 десятичный знак(!) точнее.
template< typename T >
T sin_smp(T x) {
T xx = x * x;
T res = 1.0;
for(int i = 31; i >= 3; i -= 2) {
res = T(1.0) - res * xx / T(float(i * (i - 1)));
}
return x * res;
}
Я посмотрел внимательно Вашу статью. Теперь я понял, откуда берутся Ваши заблуждения о простоте синуса, независимости точности при перетасовке аргументов и пр. Я бы не хотел никого обидеть, но ваши эксперименты сделаны очень поверхностно, к сожалению, я вынужден сказать, что безграмотно. На примере вашей же статьи я поясню своё утверждение.
Мне не удалось в точности воспроизвести ваши условия (VC++). Много лет уже не пользовался Windows.))) Компилировал GNU и clang для x86 и ARM в режиме -O0 и -O2. Результаты идентичные между системами и Вашими результатами в таблице в статье.
Вывод на экран исходной программы
**********************
* rational poly test *
**********************
x = 0.5
-------------------------
sin(x): 0.479425538604203
double: 0.479425538604203
x = 1
-------------------------
sin(x): 0.8414709848078965
double: 0.8414709848078965
x = 2
-------------------------
sin(x): 0.9092974268256817
double: 0.9092974268256816
x = 3
-------------------------
sin(x): 0.1411200080598672
double: 0.1411200080585958
x = 4
-------------------------
sin(x): -0.7568024953079282
double: -0.7568024960833886
x = 5
-------------------------
sin(x): -0.9589242746631385
double: -0.9589243758030122
Что же плохо в коде Вашей программы:
1) Использование чисел вида 1.0 и 0.5 подобных сильно уменьшает ошибку вычисления. Добавил немного «сложных» чисел, которые заполняют всю мантиссу.
2) Делать «позитивные» выводы по 5-7 точкам просто неправильно. В моих статьях я использую 16 млн точек и считаю, что надо больше.
Чтобы показать, что ваш подход даёт плохую точность надо немного изменить вашу программу program.cpp. А именно поменять входные значения для расчётов и сделать вывод немного более информативным:
Посмотрите, например на 3-ий x = 0.047619… Результаты отличаются на 2 ULP даже при таком реально малом значении x при котором классический ряд Тейлора сходится на 9 степени.
Посмотрите на x=0.99999914192480920949. Погрешность в 4 ULP!!! При значении близком к единице. Реально такой подход нуждается в глубокой переработке. Ваш расчёт за такое же количество операций умножения (а ещё плюс деление), как и ряд Тейлора даёт гораздо худшую точность.
Чтобы не было сомнения, вы можете проверить результаты в Mathematica. Например, sin(0.99999914192480920949). Видно, что это значение гораздо ближе к «внутреннему» синусу, чем к Вашему.
Уважаемый Refridgerator, я призываю Вас перепроверить Ваши результаты, например, используя framework моей программы. Я так же думаю, что было бы очень честно, если прежде чем писать комментарии Вы бы разобрались бы так же глубоко в моих статьях, как я в Вашей. А иначе Ваши комментарии просто сбивают людей с толку и направляют на ложный путь.
Спасибо за комментарий! Всё абсолютно правильно. Более того добавлю, что метод наименьших квадратов использовать нельзя. Нужно минимизировать максимальное значение отклонения. А градиент от этого посчитать ещё более невозможно. В следующих статьях я покажу, как можно аналитически аппроксимировать последний коэффициент. Последние коэффициенты, очевидно можно варьировать больше, чем последние биты. Там есть красивый математических подход.
Не хочу устраивать беспредметную дискуссию, но логарифм кажется проще. По крайней мере в glibc. После преобразований он просто считается «полиномом, очень похожим на ряд Тейлора.» (не хочу снова вступать в дискуссию о том, можно ли называть рядом Тейлора полином у которого коэффициенты различаются на несколько ULP).
Тут вопрос сложный. ИМХО диссер сейчас скорее квалификационная работа, чем настоящий результат сам по себе. Причём делается под непосредственным руководством. Аспиранты, попавшие к разным научным руководителям, находятся уже в априоре разных условиях, в независимости от уровня их знаний.
«…одной из главных причин потока научной литературы является то, что, когда исследователь достигает стадии, на которой он перестаёт видеть за деревьями лес, он слишком охотно склоняется к разрешению этой трудности путём перехода к изучению отдельных листьев».
«Ланцет», декабрь 1980 г.
Вообще есть замечательная книга «Физики шутят». Эту цитату я знаю из неё. Очень интересно, что над чем смеялись 40 лет назад и казалось шуткой, теперь норма жизни и «передовой край современной науки». Откройте например главу «Как писать научные статьи» (1966 год!). Цитирую…
Целый ряд причин (от обычной графомании до стремления
улучшить свое общественное положение) побуждает человека
писать и публиковать свои научные работы. Я не буду вдаваться
в подробности и ограничусь рассмотрением лишь четырех
главных мотивов: 1) бескорыстное стремление к распро-
странению знаний; 2) забота о собственном приоритете;
3) беспокойство за свою профессиональную репутацию;
4) стремление продвинуться по службе.
Под влиянием первой причины пишут главным образом
молодые люди, и то, по-видимому, лишь при подготовке своего
первого научного труда. Число таких авторов невелико, и для
большинства из них первая статья бывает последней.
Следовательно, первую причину нельзя ставить в один ряд с
другими более сильными мотивами, хотя забывать о ней все же
не следует.
А сейчас Светлана в своей статье уже пишет на полном серьёзе, что статьи пишутся для создания научной репутации. О «бескорыстном распространение научных знаний» в 2020-ом никто уже и не думает. Цитирую (помните, это юмроистически-сатирическая книжка):
Третья причина — забота о профессиональной репутации.
Высокой профессиональной репутации можно достичь
различными способами. Достаточно, например, сделать выдающееся
изобретение или, еще лучше, получить Нобелевскую премию,
и ваша компетентность в данном вопросе будет вне
всякого сомнения. Однако для подавляющего
большинства научных работников единственным доступным
способом является написание большого количества статей,
каждая из которых вносит в науку хотя бы небольшой вклад.
Целесообразно при этом несколько первых статей ограничить
узкой темой (например, «Соединения в волноводах»), чтобы
завоевать признание. Однако позднее автор должен
засвидетельствовать свою многосторонность, написав несколько
работ, охватывающих более широкую тему (например,
«Сверхвысокочастотные колебания»). После опубликования трех
десятков статей известность автора выйдет на насыщение и уже
не будет возрастать при дальнейшем увеличении числа печатных
работ. Тут наступает самый подходящий момент, чтобы внезапно
прекратить печататься (несколько обзорных статей не в счет) и
попытаться занять приличную руководящую должность.
Повторюсь, о чём я писал в предыдущем комментарии. Круг сам по себе порочный. А что делать я не знаю. Жду результат естественной эволюции.
Спасибо Светлане за прекраснейшую статью. Очень было интересно прочитать. Очень хорошая подача качественной информации. Как бывший учёный, я бы хотел взглянуть на проблему с другой точки зрения. Мне кажется, что проблема гораздо глубже и системнее. И происходит она от картинки 2 в статье (белоголубая со стрелочками). А именно, пока во главу угла будут ставится публикации, финансирование и научная репутация, проблемы будут только усугубляться. Поясню на примере. Наука базируется не только на журналах. Есть ещё приборы, программы, базы данных, сайты и пр. Я считаю, что труд людей которые это делают должен быть оценен. Вот и появляются публикации с gift авторами. И не обязательно gift автор это большой начальник. Может просто это человек, который сидит на каком-то научном приборе и выдаёт графики. Стоит его включать в статью или нет?
Вторая проблема современной науки — отсутствие негативный статей. Рецензирование часто формально, и, если ты уже пропихнул статью, то можешь чувствовать себя спокойно. Никогда и никогда не напишет, что статья плохая. Максимум напишут, что Ваши результаты не подтвердились и ссылку на Вас. Это порождает безнаказанность у авторов и и грантовых комиссий. И это на самом деле сидит просто в подсознании у ученых. Даже Вы в вашей статье не захотели назвать «героев» спрятав их под номерами и воздерживаясь от откровенных обвинений. Как эти люди должны почувствовать, что делают плохо? Один только «Диссернет» в поле воин. Причём в моей области знаний (физика конденсированного состояния) написание диссертаций является никому не нужным трудом. Всё читают статьи, а не диссер. В многих иностранных ВУЗах диссер это просто склеенный препринты.
Я считаю, что пока не придёт обратно настоящая научная дискуссия, как в первой половине XX (конечно на новом уровне) озвученные Вами проблемы не решить.
Спасибо за комментарий. Приводить аргумент к нужному диапазону, конечно, буду в следующих статьях. А как разберётесь, почему точность разная — пишите. Очень интересно.
Ну матчасть у вас абсолютно такая же, как в libm, которую я описал в этой статье. Но в любом случае, даже если библиотека не готова, интересно узнать, как вы использовали FPU для ускорения вычислений, какие у Вас получились точности, как вы поженили SSE и «double double». Сравнение со стандартной библиотекой и пр. Почему 4 таблицы по 64 значения. Как математически вы подобрали размер таблицы. Ну и, естественно, производительность.
Это как раз способ, который лучше избегать. Посмотрите, если не сложно, параграф «Распределительное свойство умножения в мире чисел с плавающей точкой» и в прошлой статье сравнение sin_e4 и sin_e5.
P.P.S. Проверьте существующие значения cs и sn. Они могут отличатся от ряда Тейлора. Если так, то рассчитайте их сами точно.
Если у Вас есть какой-то конкретный процессор возьмите его и сделайте эксперимент.
2) a) Есть такая суперштука как FMA. С помощью этих комманд можно проводить точное умножение.
b) Если мантисса заполнена на половину, то результат умножения будет точным. Поэтому использование битовых масок решит проблему. Правда понадобится много комманд.
sin(96889.85999999591149389743804931640625)
В вашей статье Вы работаете с 32 битным POSIT и сравниваете его с 64 битным double? И делаете вывод, что POSIT менее точный? Очевидно, что даже если все 32 бита занять мантиссой их всё равно будет меньше, чем переменной double. Приведите, плз, ссылку, где создатели posit утверждают, что posit32 точнее чем double.
Просто по структуре переменной понятно, что она будет плохо работать на больших порядках. Зато сильно лучше около 1. Это же очевидно.
Данный код вы можете добавить в проект, описанный в статье и запустить. Результат такой
Данное утверждение идёт в разрез не только с моими данными, но и данными этой статьи. Поэтому оно требует большего обоснования с Вашей стороны.
В чём «неплохость» Вашего результата? Точность, как видно нет, скорость тоже нет.
Цитирую Вашу статью:
Т.е. сейчас вы говорите, что там просто какой-то полином. А зачем тогда вы его сравниваете с синусом? Что тогда с выводом статьи? Как абстрактный полином является реальным вычислением? Цитирую Вашу статью:
Я очень внимательно прочитал Вашу статью. Вывод из неё один. Вы нашли полином, который плохо аппроксимирует синус в double и очень плохо в posit. Никаких реальных вычислений, о которых Вы говорите в выводе статьи в самой статье не показано. Реальные вычисления приведены ниже.
Если я правильно помню, то Вы не хотите поделится с нами вашими результатами? Я не собираюсь считать в double-double. У меня нет такой цели. Всё, что я увидел из Вашей статьи, что вы не умете грамотно считать ни в double ни в posit. Добавьте в свой код вот эту функцию и вы увидите, что на некоторых значениях (например на 0.5) posit на 1 десятичный знак(!) точнее.
Мне не удалось в точности воспроизвести ваши условия (VC++). Много лет уже не пользовался Windows.))) Компилировал GNU и clang для x86 и ARM в режиме -O0 и -O2. Результаты идентичные между системами и Вашими результатами в таблице в статье.
Что же плохо в коде Вашей программы:
1) Использование чисел вида 1.0 и 0.5 подобных сильно уменьшает ошибку вычисления. Добавил немного «сложных» чисел, которые заполняют всю мантиссу.
2) Делать «позитивные» выводы по 5-7 точкам просто неправильно. В моих статьях я использую 16 млн точек и считаю, что надо больше.
Чтобы показать, что ваш подход даёт плохую точность надо немного изменить вашу программу program.cpp. А именно поменять входные значения для расчётов и сделать вывод немного более информативным:
И чтобы не быть голословным — сразу привожу результаты выполнения программы.
Посмотрите, например на 3-ий x = 0.047619… Результаты отличаются на 2 ULP даже при таком реально малом значении x при котором классический ряд Тейлора сходится на 9 степени.
Посмотрите на x=0.99999914192480920949. Погрешность в 4 ULP!!! При значении близком к единице. Реально такой подход нуждается в глубокой переработке. Ваш расчёт за такое же количество операций умножения (а ещё плюс деление), как и ряд Тейлора даёт гораздо худшую точность.
Чтобы не было сомнения, вы можете проверить результаты в Mathematica. Например, sin(0.99999914192480920949). Видно, что это значение гораздо ближе к «внутреннему» синусу, чем к Вашему.
Уважаемый Refridgerator, я призываю Вас перепроверить Ваши результаты, например, используя framework моей программы. Я так же думаю, что было бы очень честно, если прежде чем писать комментарии Вы бы разобрались бы так же глубоко в моих статьях, как я в Вашей. А иначе Ваши комментарии просто сбивают людей с толку и направляют на ложный путь.
Какую функцию Вам кажется было бы посчитать сложнее всего. Напишите, я посмотрю, и может быть сделаю статью)))
«…одной из главных причин потока научной литературы является то, что, когда исследователь достигает стадии, на которой он перестаёт видеть за деревьями лес, он слишком охотно склоняется к разрешению этой трудности путём перехода к изучению отдельных листьев».
«Ланцет», декабрь 1980 г.
Вообще есть замечательная книга «Физики шутят». Эту цитату я знаю из неё. Очень интересно, что над чем смеялись 40 лет назад и казалось шуткой, теперь норма жизни и «передовой край современной науки». Откройте например главу «Как писать научные статьи» (1966 год!). Цитирую…
А сейчас Светлана в своей статье уже пишет на полном серьёзе, что статьи пишутся для создания научной репутации. О «бескорыстном распространение научных знаний» в 2020-ом никто уже и не думает. Цитирую (помните, это юмроистически-сатирическая книжка):
Повторюсь, о чём я писал в предыдущем комментарии. Круг сам по себе порочный. А что делать я не знаю. Жду результат естественной эволюции.
Вторая проблема современной науки — отсутствие негативный статей. Рецензирование часто формально, и, если ты уже пропихнул статью, то можешь чувствовать себя спокойно. Никогда и никогда не напишет, что статья плохая. Максимум напишут, что Ваши результаты не подтвердились и ссылку на Вас. Это порождает безнаказанность у авторов и и грантовых комиссий. И это на самом деле сидит просто в подсознании у ученых. Даже Вы в вашей статье не захотели назвать «героев» спрятав их под номерами и воздерживаясь от откровенных обвинений. Как эти люди должны почувствовать, что делают плохо? Один только «Диссернет» в поле воин. Причём в моей области знаний (физика конденсированного состояния) написание диссертаций является никому не нужным трудом. Всё читают статьи, а не диссер. В многих иностранных ВУЗах диссер это просто склеенный препринты.
Я считаю, что пока не придёт обратно настоящая научная дискуссия, как в первой половине XX (конечно на новом уровне) озвученные Вами проблемы не решить.