Современные численные методы в объектно-ориентированном изложении на C#

Вычислительные эксперименты в моделировании волн-убийц

Показывать лекцию целиком

Цель лекции: Продемонстрировать возможности современных вычислительных экспериментов.

На последней лекции мы рассмотрим реальные вычислительные эксперименты, которые были проведены в Институте океанологии им. П.П. Ширшова РАН совместно с Физическим институтом им. П.Н. Лебедева РАН и Университетом Аризоны (США). Были проведены масштабные вычислительные эксперименты по моделированию волн-убийц в мировом океане.

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

После получения натурных данных началось систематическое изучение этого эффекта. В лаборатории нелинейных волновых процессов Института океанологии им. П.П. Ширшова РАН проводились систематические вычислительные эксперименты по изучению волн-убийц. С помощью доказательных вычислений было показано, что волны-убийцы могут возникать в ходе нелинейной динамики идеальной жидкости со свободной поверхностью. Опишем последние вычислительные эксперименты.

Мы решали численно уравнения Эйлера, описывающие глубокую идеальную жидкость со свободной поверхностью в двумерной геометрии $$0<x<2\pi$$, $$-\infty<y<\eta(x)$$. Граничные условия на концах интервала $$x=0,2\pi$$, предполагались периодическими.

Течение предполагалось потенциальным, а жидкость несжимаемой,$$v=\nabla\phi,\quad {\rm div}\, v=0.$$ Так что потенциал подчинялся уравнению Лапласа$$\Delta\phi=0.$$ Мы осуществляли конформное отображение области, занятой жидкостью на нижнюю полуплоскость, координаты на которой $$w=u+iv$$. Отображение задается функцией $$z=z(w)$$, $$z=x+iy$$.

Динамические уравнения формулируются для переменных Дьяченко$$R=\frac{1}{z'_w}\quad V=i\frac{\partial \Phi}{\partial z}$$ и имеют вид$$\begin{gathered} R_t(u,t)=i(UR_u-U_uR), \\ V_t(u,t)=i(UV_u-B_uR) + g(R-1), \\ U=P(VR^*+RV^*), \\ B=P(VV^*), \end{gathered}$$ где $$P$$ - оператор проектирования на нижнюю полуплоскость $$P=\frac{1}{2}(1+iH)$$, $$H$$ - аналог оператора Гильберта для периодического случая$$H[f](y)=\frac{1}{2\pi}v.p.\int_0^{2\pi}\frac{f(u')}{\tan(\frac{u'-u}{2})}du'$$ Система 22.1 в настоящее время широко используется.

В наших экспериментах начальные условия определялись как ансамбль бегущих в одну сторону волн со средним значением волнового числа $$K_0=25$$.

Мы предполагали, что начальное возмущение поверхности задается суммой гармоник со случайными фазами$$\eta_0(x)=\sum\limits_{-\frac{1}{2}K_{max}}^{\frac{1}{2}K_{max}}\phi(k-k_0)\cos(kx-\xi_k)$$ Здесь $$K_{max}$$ - полное число спектральных мод, $$\xi_k$$ - случайная величина, равномерно распределенная на интервале $$-\frac{1}{2}K_{max}<k<\frac{1}{2}K_{max}$$.

Начальные значения поля скоростей предполагались связанными с 22.2 формулами линейной теории. Конформное преобразование осуществлялось при помощи итерационного алгоритма, предложенного А.И.Дьяченко.

Функция $$\phi(k)$$ определялась по формуле$$\phi(k)=\left\{% \begin{array}{ll} \delta_k, |k|>K_w; \\ \kappa\exp{(-\alpha k^2)}+\delta_k, |k|\le K_w. \\ \end{array}% \right.$$ Здесь $$\delta_k$$ - независимые случайные параметры равномерно распределенные на интервале $$-\frac{1}{2}K_{max}>k>\frac{1}{2}K_{max}$$.

Число $$1\le K_w\le 10$$ определяло спектральную ширину, $$\kappa$$, $$\alpha$$ - "внутренние" параметры спектра, определенные так, чтобы "внешние" параметры - средняя крутизна $$\mu$$$$\mu^2=\frac{1}{2\pi}\int\limits_0^{2\pi}\eta^2_xdx$$ и дисперсия $$D$$$$D=\left(\int\limits_{-K_w}^{K_w}k^2e^{-\alpha k^2}dk\right)\left(\int\limits_{-K_w}^{K_w}e^{-\alpha k^2}dk\right)^{-1}$$ принимали заданные значения. Далее, мы вычисляем точные значения полной энергии $$E$$ и следили за тем, чтобы вклад в нее случайного шума составлял не более трех процентов. Было проделано 5000 "элементарных" экспериментов. В каждом эксперименте время менялось в интервале $$0<t<200$$, что соответствовало приблизительно 500 периодам волн. Если происходило обрушение волн, счет прекращался досрочно. В расчетах полное число гармоник было $$K_{max}=2048$$ или $$K_{max}=4096$$ в зависимости от полной энергии, которая менялась в пределах $$1.5\cdot 10^{-4}\le E\le 4\cdot 10^{-4}$$.

Регистрация волн-убийц производилась следующим образом. После окончания элементарного эксперимента рассчитывалась величина $$\nu$$$$\nu=\frac{\max \eta(x,t)}{\langle|\eta|\rangle}$$ Здесь максимум в числителе берется по координате и по времени за интервал $$0<t<T$$,$$\langle|\eta|\rangle=\frac{1}{T}\int\limits_0^T\max\limits_{x\in(0,2\pi)}|\eta(x,t)|dt.$$

"Волна-убийца" фиксировалась если параметр $$\nu$$ превышал критическое значение $$\nu=1.8$$. Данное определение количественно лишь не существенно отличается от общепринятого, когда считается, что волны-убийцы вдвое превышают существенную высоту (significant wave height) Требовалось также, чтобы локальная крутизна волны $$|\eta_x|$$ превышала критическое значение $$\max\limits_{0<x<2\pi}|\eta_x|\le0.3$$. Это требование вызвано очевидными физическими соображениями и является весьма существенным.

Результаты экспериментов приведены в таблице 22.1. По горизонтали отложены значения дисперсии, по вертикальной - значения квадрат крутизны. Число "активных" мод начального условия для каждого эксперимента также показано. Из наших данных следует, что даже для волн довольно умеренной крутизны ( $$\mu^2\simeq 2.06\cdot10^{-3}$$, $$\mu\simeq 0.045$$ ) образование экстремальной волны за столь короткий отрезок времени как 500 периодов (при периоде 10 секунд это меньше полутора часов) есть весьма вероятное событие даже, если спектральная ширина по волновым числам сравнима с несущим волновым числом. Собственно, этот эксперимент и подчеркивает "обыденность" экстремальных волн. На рисунке 22.1 приведен профиль начальной волны со средней крутизной - $$\mu^2=2.56\cdot10^{-3}$$ и дисперсией $$D=4$$. На рисунке 22.2 показан профиль волны-убийцы для этой волны. Время образование волны-убийцы - $$t=67.2$$, параметр $$\nu=2.13$$ максимальная крутизна - $$0.558$$. На рисунке 22.3 приведена плотность импульса в момент образования этой волны-убийцы.

(рис 22.2) Профиль начальной волны(рис 22.1) Профиль волны-убийцы(рис 22.3) Плотность импульса в момент образования волны-убийцы

Интересно, что вероятность возникновения экстремальных волн, рассматриваемая как функция от средней крутизны при заданной дисперсии имеет максимум при весьма умеренных крутизнах ( $$\mu^2=2.0\cdot10^{-3}$$ ) а затем убывает при увеличении крутизны. Этот факт объясняется увеличением силы конкурирующего эффекта - обрушения волн. Использованная нами схема счета позволяет вести эксперимент только до первого обрушения. Мы полагаем, что при использовании более совершенных методик, зависимость вероятности возникновения экстремальных волн от крутизны остается монотонной.

На рисунке 22.4 приведены частоты возникновения волн-убийц в зависимости от дисперсии: $$\mu^2=1.54\cdot10^{-3}$$ - сплошная линия, $$\mu^2=2.56\cdot10^{-3}$$ - "тире", $$\mu^2=2.06\cdot10^{-3}$$ - точечная линия, $$\mu^2=3.08\cdot10^{-3}$$ - "точка-тире".

(рис 22.4) Частоты возникновения волн-убийц$$\begin{tabular}{|c|c|c|c|c|} \multicolumn{5}{r}{\bf Таблица 22.1}\\ \hline \textsl{$D$} \textsl{$\mu^2=$} \textsl{$\mu^2=$} \textsl{$\mu^2=$} \textsl{$\mu^2=$} \rule{0em}{1.1em}\\ \textsl{$1.54\cdot10^{-3}$} \textsl{$2.06\cdot10^{-3}$} \textsl{$2.56\cdot10^{-3}$} \textsl{$3.08\cdot10^{-3}$}\rule{0em}{1.1em}\\ \hline $\begin{array}{c} D=0.07 \\ K_w=1 \end{array}$ $0.141$ $0.638$ $0.828$ $0.849$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=2 \\ K_w=1 \end{array}$ $0.152$ $0.457$ $0.616$ $0.554$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=4 \\ K_w=2 \end{array}$ $0.011$ $0.231$ $0.346$ $0.272$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=6 \\ K_w=3 \end{array}$ $0.000$ $0.192$ $0.305$ $0.246$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=8 \\ K_w=4 \end{array}$ $0.011$ $0.154$ $0.280$ $0.195$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=10 \\ K_w=5 \end{array}$ $0.022$ $0.125$ $0.247$ $0.186$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=12 \\ K_w=6 \end{array}$ $0.010$ $0.173$ $0.256$ $0.172$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=14 \\ K_w=7 \end{array}$ $0.000$ $0.058$ $0.216$ $0.170$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=16 \\ K_w=8 \end{array}$ $0.000$ $0.136$ $0.208$ $0.151$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=18 \\ K_w=9 \end{array}$ $0.000$ $0.118$ $0.219$ $0.134$ \rule{0em}{1.1em} \\ \hline $\begin{array}{c} D=20 \\ K_w=10 \end{array}$ $0.034$ $0.127$ $0.206$ $0.099$ \rule{0em}{1.1em} \\ \hline \end{tabular}$$

Ключевые термины

Волна-убийца - поверхностная волна экстремальной амплитуды в океане.

Уравнение Эйлера - основное уравнение, описывающее динамику идеальной жидкости.

Краткие итоги: Дано описание вычислительных экспериментов. Приведены результаты о статистических характеристик волн-убийц.

Практическое занятие "Матричные игры"

Цель занятия

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

Практическая задача

Рассмотрим матричные игры на примере детской игры "камень-ножницы-бумага". Это игра двух игроков с нулевой суммой. Смысл игры состоит в том, что два игрока одновременно выбирают одну из трех стратегий: камень, ножницы или бумага. Один игрок выигрывает у другого если

  • он выбрал камень, а противник - ножницы
  • он выбрал ножницы, а противник - бумага
  • он выбрал бумага, а противник - камень
  • Если оба игрока выбрали одинаковую стратегию, то фиксируется ничья.

    Теперь мы построим математическую модель данной игры. Множество стратегий состоит из трех стратегий$$S=\{\mbox{камень},\ \mbox{ножницы},\ \mbox{бумага}\}.$$ Для простоты мы занумеруем в этом порядке эти стратегии. Поскольку мы имеем антагонистическую игру, то мы имеет платежную матрицу $$3\times3$$. Эта матрица будет указывать выигрыш или проигрыш первого игрока. Будем считать, что при выигрыше противник платит первому игроку единицу, а при проигрыше первый игрок платит второму единицу. Вот платежная матрица нашей игры$$A=\left(% \begin{array}{ccc} 0 -1 1 \\ 1 0 -1 \\ -1 1 0 \\ \end{array}% \right)$$ В этой игре также нет состояния равновесия, но есть решение в смешанных стратегиях одинаковое для обоих игроков:$$X=\left(\frac{1}{3},\frac{1}{3},\frac{1}{3}\right),$$ $$Y=\left(\frac{1}{3},\frac{1}{3},\frac{1}{3}\right).$$ Цена этой игры равна нулю.

    Для проверки проведем вычислительные эксперименты с нашей игрой. Построим класс который с помощью которого мы будем проводить вычислительные эксперименты.

    $$\begin{verbatim} class TGame { protected double[,] A; protected int m = 0, n = 0; Random rnd; public TGame() { rnd = new Random(); } public double Calc(double[] X, double[] Y, int Count) { double res = 0; int i, j; for (int k = 0; k < Count; k++) { i = Release(X, m); j = Release(Y, n); res += GetAij(i, j); } return res / (double)Count; } \end{verbatim}$$ $$\begin{verbatim} public int Release(double[] Z, int N) { double p = rnd.NextDouble(); double a = 0; for (int i = 1; i <= N; i++) { a += Z[i]; if (p <= a) { return i; } } return N; } public double GetAij(int i, int j) { return A[i, j]; } public double GetC(double[] X, double[] Y) { double res = 0; int i, j; for (i = 1; i <= m; i++) { for (j = 1; j <= n; j++) { res += A[i, j] * X[i] * Y[j]; } } return res; } } \end{verbatim}$$ $$\begin{verbatim} class TGame3 : TGame { public TGame3() : base() { m = 3; n = 3; A = new double[4, 4]; A[1, 1] = 0; A[1, 2] = -1; A[1, 3] = 1; A[2, 1] = 1; A[2, 2] = 0; A[2, 3] = -1; A[3, 1] = -1; A[3, 2] = 1; A[3, 3] = 0; } } \end{verbatim}$$

    Запускать нашу игру будем следующим образом.

    $$\begin{verbatim} double[] X; double[] Y; double l3 = 1 / 3; TGame3 Game3 = new TGame3(); X = new double[4] { 0, l3, l3, l3 }; Y = new double[4] { 0, l3, l3, l3 }; Console.WriteLine("Res = {0}", Game3.Calc(X, Y, 1000000)); \end{verbatim}$$

    В данном случае мы используем оптимальные смешанные стратегии для обоих игроков. Вот результат рассчитанной цены игры и результат розыгрыша $$1000000$$ игр.

    $$Res = 0$$

    Лабораторная работа "Матричные игры"

    Цель занятия

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

    Сценарий лабораторной работы

  • Написать программу для моделирования матричных игр в интерактивном режиме для двух лиц.
  • Задать различные матричные игры, имеющие оптимальные смешанные стратегии.
  • Рассчитать для этих матричных игр оптимальные смешанные стратегии.
  • Разыграть серию игр с реальными участниками.
  • Проанализировать эффективность использования оптимальных смешанных стратегий.
  • Указания

    Игры двух игроков с нулевой суммой удобно записывать в матричном виде. Пусть множество стратегий первого игрока равно $$n>1$$, а второго - $$m>1$$, тогда запишем в виде матрицы значения функции выигрышей$$A=\left(% \begin{array}{cccc} a_{11} a_{12} \dots a_{1n} \\ a_{21} a_{22} \dots a_{2n} \\ \dots \dots \dots \dots \\ a_{m1} a_{m2} \dots a_{mn} \\ \end{array}% \right)$$ Игра в этом случае состоит в том, что первый игрок выбирает строку, а второй игрок (одновременно!) выбирает столбец. Число, стоящее на пересечении выбранных строки и столбца, означает выигрыш первого игрока и проигрыш второго игрока.

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

    При смешанном расширении понятия матричной игры, игроки выбирают свои смешанные стратегии: первый игрок$$X=(x_1,\dots,x_m),\ x_i\ge0,\ \sum\limits_{i=1}^mx_i=1,$$ $$Y=(y_1,\dots,y_m),\ y_i\ge0,\ \sum\limits_{i=1}^my_i=1,$$

    Вернуться к учебному плану