Задачи: Python, NumPy и Numba

Численные методы для физиков

Автор

Дмитрий В. Наумов

← Вернуться к курсу

Ноутбук лекции: открыть · скачать .ipynb

Задача 1. Аналитический эталон

Частица блуждает по квадратной решётке на плоскости: на каждом шаге она сдвигается на единицу в одном из четырёх направлений с вероятностью \(1/4\).

  1. Докажите, что \(\langle \mathbf r_N^2\rangle = N\), где \(\mathbf r_N\) — положение после \(N\) шагов.
  2. Найдите \(\langle x_N^2\rangle\) для блуждания на прямой, у которого длина шага равномерно распределена на \([-1,1]\).
  3. Проверьте оба ответа моделированием на NumPy: \(M=2000\) блужданий по \(N=1000\) шагов, генератор np.random.default_rng с явным seed. Сравните отклонение от эталона со статистической ошибкой.

Задача 2. Численный эксперимент

Объём единичного шара в \(d\) измерениях \[ V_d=\frac{\pi^{d/2}}{\Gamma(d/2+1)} . \]

  1. Оцените \(V_d\) методом Монте-Карло для \(d=2,\dots,10\): \(n=10^6\) точек в кубе \([-1,1]^d\), доля точек внутри шара, умноженная на \(2^d\). Используйте один массив формы \((n,d)\) и операции NumPy без циклов по точкам.
  2. Постройте относительную ошибку оценки как функцию \(d\) и сравните с \(\sqrt{(1-p)/(np)}\), где \(p=V_d/2^d\).
  3. Объясните, почему при фиксированном \(n\) точность падает с ростом \(d\). Сколько точек нужно для относительной точности 1% при \(d=10\)?

Задача 3. Исследовательская задача

Блуждание с шагами \(\pm1\) начинается в нуле и заканчивается, когда частица впервые достигает \(x=\pm L\). Среднее время выхода равно \(\langle T\rangle=L^2\).

  1. Напишите функцию, которая моделирует одно такое блуждание шаг за шагом и возвращает \(T\). Почему эту задачу трудно векторизовать в NumPy?
  2. Оцените \(\langle T\rangle\) для \(L=10, 20, 40, 80\) по \(10^4\) блужданиям и сравните с \(L^2\) с учётом статистической ошибки.
  3. Ускорьте функцию с помощью @njit. Измерьте время чистого Python и Numba (без первого вызова) и объясните, откуда берётся разница.
  4. Найдите распределение \(T\) при \(L=20\). Чем оно отличается от гауссова?