11/01/2015 15:17
А как посчитать?
![[personal profile]](https://www.dreamwidth.org/img/silk/identity/user.png)
Вот кусок из библиотеки программ NOVAS, используемой в USNO для вычисления высокоточных эфемерид небесных тел.
Кусок интерполирует барицентрические координаты и скорости планеты, записанные в файле. Обратите внимание: интерполяция выполняется по формуле Лагранжа. Не только заголовок, но и сам текст об этом говорит.
Программу писал студент 1-го курса? Или сотрудники USNO в принципе не слышали про численные метода, поэтому и взяли первый, на который наткнулись в справочнике?
На самом деле, это частный пример, который показывает общую ситуацию. Не скрою, курс методов вычислений в СПбГУ был моим любимым. Любимым, потому что все его содержание, и даже больше, я знал задолго до 4-го курса, когда он читался. И уже тогда мне не нравился подход учебников к изложению предмета. Применительно к решением дифференциальных уравнений и методам Рунге-Кутты я об этом уже когда-то писал. А сейчас, сочиняя бибилиотеку программ для решения уравнений интегральных, скажу применительно к решению систем линейных алгебраических уравнений.
Что описывают учебники? Правильно, метод Гаусса. Некоторые еще рассказывают про схему Жордана1. Иван Петрович Мысовских в лекциях, конспект которых у меня сохранился, добавил сюда метод ортогонализации строк. Но дело даже не в том, что описывают мало методов. Дело в том, как и для чего их описывают. Вот, допустим, есть три метода. Чем один лучше другого? Молчок! Если один явно лучший, зачем тратить время на другие? Загадка! Максимум, посчитают количество операций для реализации метода. Все!
Но, дорогие мои, при вычислении на компьютере не менее важным, чем количество операций, вопросом является погрешность. Мало найти остаточный член, например, квадратурной формулы. Нужно уметь оценить реальную погрешность реального вычислительного процесса, который выполняется (неожиданно, да?) с конечным числом разрядов. Решать методом Гаусса систему 10000x10000 может быть и не так долго на современных ЭВМ, но совершенно бесполезно из-за нарастающей погрешности округления. Кто и в каком учебнике показал это? Кто предложил лучший метод?
На самом деле, эти вопросы исследованы и решены. И замечательная книга Дж. Уилкинсона "Алгебраическая проблема собственных значения" даже переведена на русский язык. Я считаю эту книгу эталоном того, как следует писать книги по численным методам. Все погрешности выявлены, проанализированы, и даны советы по их минимизации для разных типов машинных арифметик. Если вы хотите знать, как решить систему уравнений или найти собственные значения матрицы - вы узнаете это их книги Уилкинсона. Более того, автор не ограничился математическим анализом задачи. Совместно с Райншем он издал сборник конретных программ для ЭВМ, тоже переведенный, хотя и безобразно, на русский язык. Программы из этого сборника широко разошлись по миру и включены во многие известные пакеты прикладных программ, например SSP, IMSL, NAG.
Плохо одно. Книга Уилкинсона - не учебник. Хороший студент прочтет ее (если узнает) со вниманием и будет использовать вс работе. Такой студент никогда не станет писать программу интерполяции по формуле Лагранжа. Но много ли их, хороших студентов? Средний выпускник, оказавшись потом в среднем НИИ, вспомнит, когда придется решать реальную задачу, про метод Гаусса и найдет в справочнике простую программку на Бейсике, которую и использует для решения. А то, что решение оказалось неправильным, хотя в программе ошибок нет, выяснится, когда очередной "Протон" под строгим взглядом Путина пробурит очередную скважину на пути исследования мирового поземелья.
1 Схему эту придумал Вильгельм Йордан, немецкое написание фамилии которого - Jordan - неправильно транслитерировано как Жордан. Некоторые, впрочем, думают, что он итальянец, и пишут: "Метод Жордано".
* PERFORM LAGRANGIAN INTERPOLATION FOR POSITION AND VELOCITY AT
* EPOCH TJD
60 TI = TJD - XJD(LMIDDL) + LMIDDL
DO 63 J = 1, 3
POS(J) = 0.D0
VEL(J) = 0.D0
DO 62 L = 1, INTPTS
AK = ASTART + DFLOAT(L)
P = 1.D0
DO 61 I = 1, INTPTS
IF ( I .EQ. L ) GO TO 61
AI = ASTART + DFLOAT(I)
P = P * (TI-AI) / (AK-AI)
61 CONTINUE
POS(J) = POS(J) + P * BPOS(LSTART+L,J)
VEL(J) = VEL(J) + P * BVEL(LSTART+L,J)
62 CONTINUE
63 CONTINUE
GO TO 99
Кусок интерполирует барицентрические координаты и скорости планеты, записанные в файле. Обратите внимание: интерполяция выполняется по формуле Лагранжа. Не только заголовок, но и сам текст об этом говорит.
Программу писал студент 1-го курса? Или сотрудники USNO в принципе не слышали про численные метода, поэтому и взяли первый, на который наткнулись в справочнике?
На самом деле, это частный пример, который показывает общую ситуацию. Не скрою, курс методов вычислений в СПбГУ был моим любимым. Любимым, потому что все его содержание, и даже больше, я знал задолго до 4-го курса, когда он читался. И уже тогда мне не нравился подход учебников к изложению предмета. Применительно к решением дифференциальных уравнений и методам Рунге-Кутты я об этом уже когда-то писал. А сейчас, сочиняя бибилиотеку программ для решения уравнений интегральных, скажу применительно к решению систем линейных алгебраических уравнений.
Что описывают учебники? Правильно, метод Гаусса. Некоторые еще рассказывают про схему Жордана1. Иван Петрович Мысовских в лекциях, конспект которых у меня сохранился, добавил сюда метод ортогонализации строк. Но дело даже не в том, что описывают мало методов. Дело в том, как и для чего их описывают. Вот, допустим, есть три метода. Чем один лучше другого? Молчок! Если один явно лучший, зачем тратить время на другие? Загадка! Максимум, посчитают количество операций для реализации метода. Все!
Но, дорогие мои, при вычислении на компьютере не менее важным, чем количество операций, вопросом является погрешность. Мало найти остаточный член, например, квадратурной формулы. Нужно уметь оценить реальную погрешность реального вычислительного процесса, который выполняется (неожиданно, да?) с конечным числом разрядов. Решать методом Гаусса систему 10000x10000 может быть и не так долго на современных ЭВМ, но совершенно бесполезно из-за нарастающей погрешности округления. Кто и в каком учебнике показал это? Кто предложил лучший метод?
На самом деле, эти вопросы исследованы и решены. И замечательная книга Дж. Уилкинсона "Алгебраическая проблема собственных значения" даже переведена на русский язык. Я считаю эту книгу эталоном того, как следует писать книги по численным методам. Все погрешности выявлены, проанализированы, и даны советы по их минимизации для разных типов машинных арифметик. Если вы хотите знать, как решить систему уравнений или найти собственные значения матрицы - вы узнаете это их книги Уилкинсона. Более того, автор не ограничился математическим анализом задачи. Совместно с Райншем он издал сборник конретных программ для ЭВМ, тоже переведенный, хотя и безобразно, на русский язык. Программы из этого сборника широко разошлись по миру и включены во многие известные пакеты прикладных программ, например SSP, IMSL, NAG.
Плохо одно. Книга Уилкинсона - не учебник. Хороший студент прочтет ее (если узнает) со вниманием и будет использовать вс работе. Такой студент никогда не станет писать программу интерполяции по формуле Лагранжа. Но много ли их, хороших студентов? Средний выпускник, оказавшись потом в среднем НИИ, вспомнит, когда придется решать реальную задачу, про метод Гаусса и найдет в справочнике простую программку на Бейсике, которую и использует для решения. А то, что решение оказалось неправильным, хотя в программе ошибок нет, выяснится, когда очередной "Протон" под строгим взглядом Путина пробурит очередную скважину на пути исследования мирового поземелья.
1 Схему эту придумал Вильгельм Йордан, немецкое написание фамилии которого - Jordan - неправильно транслитерировано как Жордан. Некоторые, впрочем, думают, что он итальянец, и пишут: "Метод Жордано".
no subject
no subject
На моих глазах один аспирант - выпускник Казанского университета - искал собственные значения матрицы 18-го порядка построением характеристического многочлена и отысканием его корней. То, что полученные таким путем числа скорее всего не имеют отношения к действительности, просто не пришло ему в голову.
no subject
no subject
Почему же не имеют? Поставить, скажем, "precision 300;" - и вперёд! Меняя этот параметр, программист может легко определить требуемую точность, не имея никакого представления о сколь-либо "продвинутых" численных методах. Наверное, здесь вполне хватит 30-40 знаков после запятой, если случай не совсем патологический (скажем, среди матричных элементов нет частичных сумм рядов Фурье и т.п.).
no subject
no subject
no subject
no subject
1) "Чтобы избежать катастрофического влияния вычислительной погрешности, применяют метод Гаусса с выбором главного элемента" (С. 257). До этого - ни слова о том, что погрешность может быть катастрофической.
2) Простой прием для оценки погрешности решения на с. 258. Характеристика полезности этой оценки от авторов: "...отличается не в очень большое число раз".
3) На с. 262 упоминается, что в ходе вычислений мера обусловленности матрицы не должна меняться сильно.
Все!
Основной упор сделан на подсчет числа операций. О том, что при плохой обусловленности матрицы (и том что такое обусловленность) любой из изложенных методов даст абсурдный результат - ни слова.
А между прочит, БЖК - основной вузовский учебник по численным методам.