Алгоритмы

Некоторые полезные алгоритмы

Метод Гаусса для решения систем линейных уравнений

Достатачно просто решать такие уравнения методами матричной алгебры (вычисляя детерминанты и т.д.) но чаще удобнее использовать метод Гаусса. Суть метода в последовательном исключении из уравнений переменных, домножая уравнение на некоторый коэффициент и вычитая из него другое уравнение системы так, чтобы коэффициент при одной из переменных в этом уравнении стал нулем. После таких преобразований система выглядит так:

a*x + b*y + c*z + d = 0
e*x + f*y + g = 0
i*x + j = 0

Такую систему решать значительно легче, чем исходную. Пример Gauss.pas.

Обратная польская нотация

Наиболее простая система записи математических выражений для электронной обработки - обратная польская нотация. Например, выражение "(a + b)*(c - d)/sin(e)" в этой системе записывается как "a b + c d - * e sin /". В этой системе не нужны скобки. и все выражения могут быть вычислены с помощью простого алгоритма, использующего стек. Когда при проходе строки встречается переменная, она помещается в стек, а когда встречается знак математической операции, из стека извлекается нужное количество переменных, осуществляется операция и результат снова помещается в стек. Когда вся строка пройдена, остается лишь извлечь из стека результат (если формула была записана правильно, то это будет единственный элемент стека).

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

Функции (напр. sin, cos) и скобки
Умножение и деление
Сложение и вычитание

Алгоритм перевода таков:

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

Пример Calculator.pas.

Приблеженные методы

Почти все приближенные методы основаны на разложении Тейлора:

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

Пример Tailorrw.pas

Чтобы решать системы любых уравнений, можно использовать такой метод. Пусть исходная система записана в форме f_i(x₁,x₂,...,x_n)=0, i=1..n. Используя линейное приближение, найдем k'е уточненное значение как решение такой ситемы линейных ураванений относительно x_j^k:

В случае одномерной функции эти операции выглядят так:

Некоторые полезные уравнения

- Двумерные волны на поверхности жидкости с учетом вязкости.
- Одномерный солитон на поверхности жидкости
- Одномерный солитон в плазме, на поверхности жидкости (без вязкости)
- Уравнение теплопроводности
- Уравнение Шредингера
- Уравнение Ван-Дер-Поля для лампового генератора

Моделирование голографии и интерференции

Моделирование интерференции и голографии как частного случая интерференции требует сложения большого количества когерентных полей (с разными фазами и амплитудами, но одинаковыми частотами): . Результирующая амплитуда вычисляется как
Пример Golography.zip

Движение заряженных частиц.

Все моделирование этого явления заключается в решении уравнений движения в суммарном поле всех остальных частиц. Пример Motion.zip

Сжатие изображений с использованием двумерного разложения Фурье

Как и в случае одной переменной, в двумерном случае разложение Фурье задается такой формулой:

Теперь мы имеем дело с волновыми векторами k и фаза вычисляется как f=k*r, где r = {x, y}. Исходное изображение преобразуется в набор волновых векторов. Качество отличное ;)), но ведь это только экспериментальная программа. Пример Compression.pas.

Игра "Жизнь"

Жизнь - игра:
сюжет не очень, зато графика отличная

В этой версии популярной игры есть "Мирные" (положительные числа) и "Чужие" (отрицательные числа), всех их можно "кормить", причем "Мирные" и "Чужие" едят "положительную" и "отрицательную" пищу. Клавиши "up","down","left","right" передвигают текущую позицию, клавишами +/- (без Shift) устанавливаются положительные и отрицательные "существа", клавишами "0" и "9" на месте курсора устанавливается еда (положительная либо отрицательная). Нажмите Esc чтобы начать развитие существ.
Еще два интересных примера - Жизнь в трехмерном пространстве и Жизнь на поверхности тора. Поверхность тора топологически ограничена, но границ не имеет.
Пример программы Life.pas
Следующих два примера - для Delphi и OpenGL под Delphi: Life3D.zip LifeTor.zip
Готовые программы: LifeTor.exe

Метод наименьших квадратов

Метод наименьших квадратов - простой и быстрый способ получить неизвестные параметры в функциональных зависимостях и оценить их погрешности. Пусть ожидаемая теоретическая зависимость -
y = f(x), и мы получили ряд значений(x_i, f(x_i)). Тогда величину ошибки можно оценить как сумму квадратов всех отклонений от теоретической зависимости: (x_i - x₀)², где х₀ - среднее значение х. Для достижения наилучшей точности эта ошибка должна быть минимальной. Возьмем от полученной суммы по всем параметрам производные и приравняем их к нулю - получим систему уравнений для этих параметров, решениями которой и будут наиболее вероятные их значения. В случае линейной зависимости (а практически любая зависимость может быть линеаризована) имеем:
y = a*x + b


Если же параметр равен нулю по определению, то формулы становятся еще проще:
y = a*x

Такой метод можно применять для любых функциональных зависимостей. Например, такая задача возникла во время исследований свойств фракталов (см. Теория протекания). Результаты серии экспериментов - булевы величины Res_i, зависящие от некоторого вещественного параметра q_i. Необходимо определить такое критическое значение параметра q_cr, что для большинства q_i<q_cr результат - true (позитивный), а для большинства q_i>q_cr результат - false (отрицательный). Чтобы найти это значение параметра я рассмотрел функцию ошибок R_i(Res_i,q_i,q_cr):
R_i(Res_i,q_i,q_cr) = 1, если q_i<q_cr и Res_i=False или q_i>q_cr и Res_i=True
R_i(Res_i,q_i,q_cr) = 0 во всех остальных случаях.
Суммарная ошибка для некоторого параметра q_cr :
E =SR_i(Res_i,q_i,q_cr)(q_i-q_cr)²
Чтобы минимизировать ее, будем рассматривать R_i(Res_i,q_i,q_cr) как константы и используем некоторое число в качестве начального значения для q_cr. Дифференцируя функцию суммарной ошибки E по q_cr, находим первое приближение для q_cr. Это рекуррентный процесс, и j-е приближение для q_cr:
q_cr j =SR_i(Res_i,q_i,q_cr j-1)*q_i / SR_i(Res_i,q_i,q_cr j-1)
После нескольких итераций может возникнуть неопределенность "0 / 0". В этом случае итерации можно остановить, так как правильное значение уже найдено.

Тексты программ для общего случая и случая прямой пропорциональности соответственно LSM.pas LSM1.pas
Готовые программы LSM.exe LSM1.exe

Вычисление неявных функций, заданных определенными интегралами

Функции часто задаются неявно при помощи неопределенных интегралов:
g(a) = т^b1_b2 f(x,a)dx
Как b1, так и b2 могут стремится к бесконечности, либо сама функция f(x,a) может обращаться в бесконечность в некоторых точках между b1 и b2. Приближенное вычисление неопределенного интеграла - не самый быстрый способ получить значение такой неявной функции. Более того, этот способ может оказаться практически бесполезным, когда пределы интегрирования b1 и\или b2 стремятся к бесконечности. Разложение функции в ряд и аналитическое интегрирование всех членов ряда - довольно быстрый метод. Простейшие ряды - ряд Тейлора и ряд Фурье. Ряду Тейлора присущ один недостаток - его нельзя проинтегрировать, когда пределы интегрирования стремятся к бесконечности. Тогда для разложения можно использовать отрицательные степени x или любые другие функции, интегрируемые на бесконечности. Например:
f(x) = a₂/x² + a₃/x³ ... или:
f(x) = a₁e^-x + a₂e^-2x ...
и так далее. Коэффициенты a_i могут быть получены так же, как и для ряда Тейлора:
a_i = 1/i! dⁱf(x)/df₀ⁱ
где f₀ - базисная функция разложения. Коэффициенты могут быть также найдены приближенно, например, при помощи метода наименьших квадратов. Теперь интегрирование осуществляется элементарно:
т^Ґ_a f(x)dx = т^Ґ_a a₂/x² + т^Ґ_a a₃/x³ + ... = -a₂/a - 0.5*a₃/a²...
Остается лишь просуммировать полученный ряд с необходимой точностью. Очень важно, чтобы разложение по отрицательным степеням начиналось с x^-2 - это позволяет проинтегрировать ряд.
Бета-функция и гамма-функция - важные неявные функции.
Бета-функция определяется как:
B(p,q) = т¹₀x^p-1(1-x)^q-1 dx (p>0, q>0)
Гамма-функция определяется как:
G(s) = т^Ґ₀x^s-1*Exp(-x) dx (s>0)
Процедура, возвращающая коэффициенты разложения (1-x)^q в ряд Тейлора по x и процедура, вычисляющая бета-функцию, используя разложение в ряд Тейлора при интегрировании, включены в модуль Matrix.pas.

Сайт создан в системе uCoz