Алгоритм Ремеза в экселе
Вот и до него руки дошли, причина станет ясна в следующем посте. Изучать чужие библиотеки было лениво (в том же BOOSTе сам чёрт ногу сломит), писать самому - тоже, а вот "в полуавтоматическом режиме" в Excel - вышло в пол-пинка.
Для примера возьмём приближение
)
при условии x2+y2=1.
Чтобы подобрать коэффициенты, делаем "замену" x=cos(z), y=sin(z), и по сути:
=\alpha\sin{z}+\frac{\beta}{2}\sin{2z})
Используя метод наименьших квадратов и "вручную" подбирая диапазон оптимизации, мы смогли подобрать α, β чтобы арктангенс считался с точностью не хуже 2,5 угловых минут. Позже, сделав интерполяцию по узлам Чебышёва, смог получить результат чуть лучше: 2,36 угловых минут. Интересно, можно ли улучшить этот результат, и насколько сильно.
Нам понадобится график ошибки, который приводили в прошлый раз:
Мы должны выписать все значения z, на которых достигаются минимумы и максимумы этой функции ошибки. Я его строил с шагом в 0,1 градус, поэтому значения получаются такие:
13,8° (минимум); 36,3° (максимум); 45° (конец диапазона, "в минус").
Обычно для самой первой итерации, когда ещё НИКАКОГО ПРИБЛИЖЕНИЯ НЕТ, вместо нахождения этих минимумов и максимумов просто берут узлы Чебышева, сколько нужно. В моём случае весь диапазон от -45 до 45 градусов, на нём получится 6 максимумов/минимумов, значит, взяли бы n=6.
Вышло бы 11,65°; 31,82°, 43,46° и такие же со знаком "-". Довольно близко...
И далее решаем систему линейных уравнений: мы хотим, чтобы
Неизвестными здесь являются α, β и величина ошибки E. Т.е мы ТРЕБУЕМ, чтобы функция ошибки прыгала вверх-вниз, и чтобы все её максимумы и минимумы были одинаковы!
Это система линейных уравнений. Я в Excel сначала выписал значения z1, z2, z3:
Формулы частично на русском языке - то ещё "удовольствие", ну да ладно. Главное, у меня переключение раскладки стоит по Ctrl+Shift, потому что Alt+Shift слишком уж часто "теряет фокус", куда-то на меню его уносит, бесило страшно!
Нам, понятно, это в радианах нужно, градусы это мне для удобства. Выходит так:
Далее ввожу формулы для коэффициентов системы уравнений. Одну строчку "ручками", затем "растягиваю" на строчку вниз, и единицу во второй строке заменяю на "умножение прошлого значения на -1" , и растягиваю ещё на строчку. Вот показываю формулы:
В виде значений:
Далее требуется "ловкость рук", чтобы найти обратную матрицу. Сначала выделяем пустую область 3х3 где-нибудь ниже:
И вводим формулу:
Но закрыв скобку, нажимаем Ctrl+Shift+Enter! Только в этом случае он воспримет результат этой операции в виде массива 3х3 - и распределит его по выбранной области:
Осталось умножить эту обратную матрицу на столбец (Z1,Z2,Z3), который у нас уже торчит наверху. Для этого выделяем пустой столбец в 3 ячейки, вводим в него формулу =МУМНОЖ(A10:C12;B2:B4) и снова Ctrl+Shift+Enter. Вуаля:
На этом одна итерация алгоритма Ремеза завершена!
Мы получили новые значения α и β и ещё значение E=0,000664579 радиан = 0,038 градусов. Правда, не нужно радоваться раньше времени! данное значение ещё не означает, что ошибка будет именно такова - такая она получилась в этих 3 точках, а больше ни о чём данный алгоритм "не знает". Чтобы узнать реальную ошибку, надо снова построить графики. Я это делаю "по-старинке", вот так:
Понятно: первую строку ввожу "от руки", причём ткнувшись мышкой по "альфе", тут же жму F4, что "оборачивает" её в доллары, т.е при копировании данная ячейка смещаться никуда не будет. А дальше достаточно "растянуть" вниз, пока не дойдёт до угла в 45 градусов. Вот что у нас выходит:
Очень чётко все 3 максимума/минимума совпадают по абсолютной величине: 0,038°. Можно ещё "честно" найти минимальное и максимальное значения на этой "простыне", и мы получим:
-0,038079292
0,038082253
Вот этому значению уже вполне можно доверять!
Давайте проведём ещё одну итерацию, это очень просто:
надо найти координаты новых максимумов и минимумов, это 13,72°, 36,22° и, разумеется, 45°. Подставляем их в столбец z1,z2,z3, с которого мы начали - и далее всё пересчитается автоматически.
Но "визуально" ничего не изменилось, но теперь мы получаем следующие максимум и минимум:
-0,038080136
0,038080136
Вот теперь они совпали, насколько хватает разрядности, и на самую-самую малость максимальная ошибка снизилась, аж на 0,0076 угловой секунды!
В общем, делать ещё итерации большого смысла уже нет.
Как я и ожидал, удалось добиться улучшения, но не очень большого: с 2,36 угловых минут до 2,28 угловых минут (ещё минус 3,3%) и это, похоже, предел для такой формулы.
Кстати, далеко не факт, что "стандартная реализация" алг. Ремеза мне бы помогла. Тут были заморочки с тем, что функция нечётная, поэтому я брал полный интервал, но на деле использовал точки только на его половине и ни капли не возмущался, что в нуле ошибка уходит в ноль. По-моему, если здесь подойти слишком "в лоб" - система уравнений получилась бы вырожденной, с повторяющимися значениями, и пришлось бы танцевать с бубном, соображая, как так ему это впихнуть. А ещё далеко не все принимают такой интересный набор базисных функций, как синусы. Для большинства используются многочлены, но только от одной переменной. Да и самостоятельно писать программу, которая бы позволяла всё это дело учесть - я бы зарылся ещё на недельку-другую...
Для примера возьмём приближение
при условии x2+y2=1.
Чтобы подобрать коэффициенты, делаем "замену" x=cos(z), y=sin(z), и по сути:
Используя метод наименьших квадратов и "вручную" подбирая диапазон оптимизации, мы смогли подобрать α, β чтобы арктангенс считался с точностью не хуже 2,5 угловых минут. Позже, сделав интерполяцию по узлам Чебышёва, смог получить результат чуть лучше: 2,36 угловых минут. Интересно, можно ли улучшить этот результат, и насколько сильно.
Нам понадобится график ошибки, который приводили в прошлый раз:
Мы должны выписать все значения z, на которых достигаются минимумы и максимумы этой функции ошибки. Я его строил с шагом в 0,1 градус, поэтому значения получаются такие:
13,8° (минимум); 36,3° (максимум); 45° (конец диапазона, "в минус").
Обычно для самой первой итерации, когда ещё НИКАКОГО ПРИБЛИЖЕНИЯ НЕТ, вместо нахождения этих минимумов и максимумов просто берут узлы Чебышева, сколько нужно. В моём случае весь диапазон от -45 до 45 градусов, на нём получится 6 максимумов/минимумов, значит, взяли бы n=6.
Вышло бы 11,65°; 31,82°, 43,46° и такие же со знаком "-". Довольно близко...
И далее решаем систему линейных уравнений: мы хотим, чтобы
Неизвестными здесь являются α, β и величина ошибки E. Т.е мы ТРЕБУЕМ, чтобы функция ошибки прыгала вверх-вниз, и чтобы все её максимумы и минимумы были одинаковы!
Это система линейных уравнений. Я в Excel сначала выписал значения z1, z2, z3:
Формулы частично на русском языке - то ещё "удовольствие", ну да ладно. Главное, у меня переключение раскладки стоит по Ctrl+Shift, потому что Alt+Shift слишком уж часто "теряет фокус", куда-то на меню его уносит, бесило страшно!
Нам, понятно, это в радианах нужно, градусы это мне для удобства. Выходит так:
Далее ввожу формулы для коэффициентов системы уравнений. Одну строчку "ручками", затем "растягиваю" на строчку вниз, и единицу во второй строке заменяю на "умножение прошлого значения на -1" , и растягиваю ещё на строчку. Вот показываю формулы:
В виде значений:
Далее требуется "ловкость рук", чтобы найти обратную матрицу. Сначала выделяем пустую область 3х3 где-нибудь ниже:
И вводим формулу:
Но закрыв скобку, нажимаем Ctrl+Shift+Enter! Только в этом случае он воспримет результат этой операции в виде массива 3х3 - и распределит его по выбранной области:
Осталось умножить эту обратную матрицу на столбец (Z1,Z2,Z3), который у нас уже торчит наверху. Для этого выделяем пустой столбец в 3 ячейки, вводим в него формулу =МУМНОЖ(A10:C12;B2:B4) и снова Ctrl+Shift+Enter. Вуаля:
На этом одна итерация алгоритма Ремеза завершена!
Мы получили новые значения α и β и ещё значение E=0,000664579 радиан = 0,038 градусов. Правда, не нужно радоваться раньше времени! данное значение ещё не означает, что ошибка будет именно такова - такая она получилась в этих 3 точках, а больше ни о чём данный алгоритм "не знает". Чтобы узнать реальную ошибку, надо снова построить графики. Я это делаю "по-старинке", вот так:
Понятно: первую строку ввожу "от руки", причём ткнувшись мышкой по "альфе", тут же жму F4, что "оборачивает" её в доллары, т.е при копировании данная ячейка смещаться никуда не будет. А дальше достаточно "растянуть" вниз, пока не дойдёт до угла в 45 градусов. Вот что у нас выходит:
Очень чётко все 3 максимума/минимума совпадают по абсолютной величине: 0,038°. Можно ещё "честно" найти минимальное и максимальное значения на этой "простыне", и мы получим:
-0,038079292
0,038082253
Вот этому значению уже вполне можно доверять!
Давайте проведём ещё одну итерацию, это очень просто:
надо найти координаты новых максимумов и минимумов, это 13,72°, 36,22° и, разумеется, 45°. Подставляем их в столбец z1,z2,z3, с которого мы начали - и далее всё пересчитается автоматически.
Но "визуально" ничего не изменилось, но теперь мы получаем следующие максимум и минимум:
-0,038080136
0,038080136
Вот теперь они совпали, насколько хватает разрядности, и на самую-самую малость максимальная ошибка снизилась, аж на 0,0076 угловой секунды!
В общем, делать ещё итерации большого смысла уже нет.
Как я и ожидал, удалось добиться улучшения, но не очень большого: с 2,36 угловых минут до 2,28 угловых минут (ещё минус 3,3%) и это, похоже, предел для такой формулы.
Кстати, далеко не факт, что "стандартная реализация" алг. Ремеза мне бы помогла. Тут были заморочки с тем, что функция нечётная, поэтому я брал полный интервал, но на деле использовал точки только на его половине и ни капли не возмущался, что в нуле ошибка уходит в ноль. По-моему, если здесь подойти слишком "в лоб" - система уравнений получилась бы вырожденной, с повторяющимися значениями, и пришлось бы танцевать с бубном, соображая, как так ему это впихнуть. А ещё далеко не все принимают такой интересный набор базисных функций, как синусы. Для большинства используются многочлены, но только от одной переменной. Да и самостоятельно писать программу, которая бы позволяла всё это дело учесть - я бы зарылся ещё на недельку-другую...