5 Oct 2021
Решение кубического уравнения
x^3 + a*x^2 + b*x + c = 0
a, b, c <- R
Всегда есть один вещественный корень.
1. С помощью замены x' = x + s можно добиться, чтобы a == 0
a + 3*s == 0 => s = -a/3
Будем считать, что уравнение имеет вид
x^3 + a*x + b = 0
Основная идея: сделаем замену x = u + v
(u + v)^3 + a*(u + v) + b = 0
u^3 + 3uv(u + v) + v^3 + a*(u + v) + b = 0
u^3 + v^3 + b + (u + v)(3uv + a) = 0
------------- ================
{ u^3 + v^3 + b = 0
{ 3uv + a = 0
Надо решить эту систему. Выразим v через u из второго уравнения
v = -a/(3*u)
и подставим в первое уравнение:
u^3 - a^3/(27*u^3) + b = 0
Сделаем замену y = u^3
y - a^3/(27y) + b = 0
y^2 + b y - a^3/27 = 0
Решаем квадратное уравнение. Обозначим
d = b^2 + 4*a^3/27
d = complex(d)**(1/2)
y = (-b + d)/2
Дома: решить уравнение
x^3 + a*x^2 + b*x + c = 0
a, b, c <- R
Решение присылайте мне на почту:
vladibor266@gmail.com
===============================================
14 Oct 2021
x^3 - x - 1 = 0
std::vector
CArray
size()
capacity()
****************..........
************************** + *
***************************............
m0
m1 = m0*d, d > 1
m2 = m1*d = m0*d^2
m3 = m2*d = m0*d^3
Как правильно выбрать d?
d = 2 в Linux
d = 1.5 в Windows
d = 1.25 в Java
d = 1.125 в MFC
Идея: надо выбрать d так, чтобы можно было использовать память,
освобожденную на предыдущих шагах.
m3 <= m0 + m1, в экстремальном случае
m3 = m0 + m1
d^3*m0 = m0 + d*m0
d^3 - d - 1 = 0
Вещественное решение x = 1.32471795724475
(Сродни константе "золотое сечение"
x = (1 + sqrt(5))/2 = 1.61803398874989
которая есть решение уравнения x^2 - x - 1 = 0).
=====================================================================
Домашние задания
1. Всегда оформляйте алгоритм в виде функции.
Функция не должна ничего печатать, она должна возвращать
результат.
2. Лучше оформлять любую программу на Питоне как модуль.
Это означает, что при загрузке модуля не должен выполняться
никакой код на внешнем уровне (вне функций).
==================================================================
Рабочие лошадки:
1) наибольший общий делитель
d = gcd(a, b)
2) алгоритм быстрого возведения в степень:
p = fastPower(a, n)
p = powermod(a, n, m)
List comprehension
(умное задание списка, математическое задание списка)
Математика:
M = { x: x <- [1, 2, ..., 100], x простое }
Python:
s = [ x for x in range(1, 101) if isPrime(x) ]
Факторизация числа (разложение на множители)
60 = 2^2 * 3 * 5
Функция возвращает список пар
[(2, 2), (3, 1), (5, 1)]
def factor(m):
'''Возвращает список пар: (простое число, его степень),
составляющих разложение m на множители'''
rho-алгоритм факторизации Полларда
Дано целое число m
Надо: найти нетривиальный делитель числа m
f: Z_m --> Z_m
f: x --> x^2 + 1 (mod m)
x -- случайное число
y = f(x)
d = gcd(x - y, m)
while d == 1:
x = f(x)
y = f(f(y))
d = gcd(x - y, m)
x, f(x), f(f(x)), ...., x_i, x_i+1 = f(x_i)
x0, x1, x2, x3, ....
Орбита элемента x: начальный отрезок непериодический и затем
входит в цикл.
x_i+k == x_i цикл длины k.
Рассматриваем эти числа как элементы Z_m (по модулю m)
Пусть p | m
Если те же самые числа x_i рассмотреть по модулю p, то цикл этой
последовательности по модулю p будет, скорее всего, меньше.
x_i == x_j (mod p)
x_i != x_j (mod m)
p | (x_i - x_j)
m не делит (x_i - x_j)
Значит, gcd(x_i - x_j, m) делится на p и не делится на m =>
это нетривиальный делитель числа m.
Реализация алгоритма на Питоне:
------------------------------------------------
import random
from math import gcd
def f(x, m):
return (x*x + 1)%m
def pollardRhoFactor(m, maxSteps = 10000000):
x = random.randrange(2, m-1)
y = f(x, m)
d = gcd(x - y, m)
step = 0
while d == 1 and step <= maxSteps:
x = f(x, m)
y = f(y, m)
y = f(y, m)
d = gcd(x - y, m)
step += 1
if d == m:
d = 1
return d
Дома: написать функцию factor с использованием rho-алгоритма Полларда
http://mech.math.msu.su/~vvb/
===============================================================
19 Oct 2021
Найти все счастливые билеты длины n для четного n.
Билет счастливый, когда сумма первых n//2 цифр равна сумме
последних n//2 цифр.
def nextTicket(t):
n = len(t)
pos = n - 1
while pos >= 0:
if t[pos] < 9:
t[pos] += 1
return True
t[pos] = 0
pos -= 1
return False
def isLucky(t):
n = len(t)
k = n // 2
s0 = 0
for i in range(0, k):
s0 += t[i]
s1 = 0
for i in range(k, n):
s1 += t[i]
'''
if s0 == s1:
return True
else:
return False
'''
return (s0 == s1)
def luckyTickets(n):
assert n > 0 and n%2 == 0
t = [0]*n
s = []
while True:
if isLucky(t):
# s.append(t) ERROR!
s.append(t.copy()) # Correction: save a copy of mutable object
if not nextTicket(t):
break
return s
Дома, задача 1 (для всех):
найти число счастливых билетов длины n (именно число, а не список!).
=====================================================================
Вернемся к теории чисел.
Схема кодирования с открытым ключом RSA
m = p*q, где p, q -- большие простые числа
m открыто, разложение m на множители секретно.
p, q 200-значные, m -- 400-значное.
Умеем создавать большие простые числа, но не умеем
раскладывать на множители большие числа.
Хотим проверить простоту числа.
Какие эксклюзивные свойства есть у простого числа?
1. Малая теорема Ферма
m простое, b != 0 (mod m) => b^(m-1) == 1 (mod m)
Греки заблуждались, считая, что обратная теорема
к Малой теореме Ферма также верна (и даже для конкретного
основания 2):
2^(m-1) == 1 (mod m) => m простое.
В начале XIX века французский П.Сарруз предъявил контрпример:
m = 341 = 31*11
2^340 = (2^10)^34 = (1024)^34 == (1)^34 == 1 (mod 341)
1024 = 341*3 + 1 == 1 (mod 341)
Американский математик Р.Кармайкл (Carmichael)
в 1910 г. нашел число m = 561, которое ведет себя
как простое число в Малой теореме Ферма:
если gcd(b, m) = 1, то b^(m-1) == 1 (mod m)
Такие числа называеются Кармайкловыми.
561, 1105, 1729, ...
Их бесконечно много.
Задача на дом для нечетных номеров:
Найти все кармайкловы числа, не превышающие n.
Китайская теорема об остатках:
m = m1*m2*...*mk, gcd(mi, mj) == 1 =>
Z_m == Z_m1 + Z_m2 + ... + Z_mk (изоморфизм колец,
сумма в правой части прямая)
x = (x1, x2, x3, ..., x_k) x_i == x (mod m_i)
105 = 3 * 5 * 7
Сколько квадратных корней из 1 в Z_105?
Z_105 = Z_3 + Z_5 + Z_7
x == (x%3, x%5, x%7)
x^2 = (x^2%3, x^2%5, x^2%7) = (1, 1, 1)
x_i - квадратный корень из 1 => x_i = +- 1
Всего 2^3 = 8 корней:
(1, 1, 1)
(-1, 1, 1)
(1, -1, 1)
(-1, -1, 1)
. . .
(-1, -1, -1)
Теорема.
m кармайклово <=>
1) m свободно от квадратов: m = p1*p2*...*p_k
2) k >= 3
3) для любого i (p_i - 1) | (m - 1)
Пример:
561 = 3 * 11 * 17,
(3 - 1)|560, (11 - 1)|560, (17 - 1)| 560
Это был негатив: Малую теорему Ферма нельзя использовать в
качестве теста простоты числа.
Но, оказывается, МТФ несложно чуть-чуть подправить, чтобы
получился (вероятностный) тест простоты.
Тест Миллера-Рабина:
для числа m может выдавать 2 ответа:
1) m составное (выдается доказательство этого!);
2) не знаю.
Причем вероятность второго ответа для составного числа <= 0.25.
Тест использует датчик случайных чисел.
Если для m мы получили 20 ответов "не знаю",
то вероятность того, что оно составное, не превышает (1/4)^20.
То есть число m, скорее всего, простое, но, к сожалению,
мы не получаем доказательства простоты.
http://mech.math.msu.su/~vvb/
Задача на дом для четных номеров:
Сгенерировать случайное простое число заданной длины в
двоичной записи:
p = generatePrime(n)
p имеет длину ровно n в двоичной записи
2^(n-1) <= p < 2^n
Простых чисел очень много ~ n/ln(n) (Теорема Чебышева)
Тест Миллера-Рабина.
Две идеи:
1) Малая теорема Ферма
2) в любом поле есть не больше двух корней из 1:
1 и -1 (поскольку корни удовлетворяют уравнению x^2 - 1 = 0,
которое имеет не больше 2-х решений по теореме Безу).
Проверяем на простоту нечетное число m.
1. Генерируем случайное основание b, 2 <= b < m-1
2. m - 1 = t * 2^s, где t -- нечетное число;
3. Вычисляем последовательность
a0 = b^t (mod m)
a1 = a0^2 = b^(2*t) (mod m)
a2 = a1^2 (mod m)
. . .
a_s = a_s-1^2 = b^(m - 1) (mod m)
каждое следующее число в этой последовательности есть
квадрат предыдущего по модулю m.
Число составное, если
1) либо a_s != 1 (mod m) (в этом случае не выполняется
Малая теорема Ферма),
2) либо в последовательности встретится фрагмент
*, *, *, 1, 1, ..., 1
где * отлична от +-1 (mod m) (в этом случае мы нашли
нетривиальный корень из 1,
которых в поле Z_m
при простом m не бывает)
Ответ "не знаю" выдается для последовательности вида либо
1, 1, ..., 1
либо
*, *, ..., *, -1, 1, ..., 1
Запись алгоритма Миллера-Рабина на псевдокоде есть на странице
http://mech.math.msu.su/~vvb/Master/NumberTheory.html
Реализация на Питоне оставлена в качестве задачи для самостоятельного
решения.
Желательно оформлять домашние задания в виде функций,
на вход которых подаются исходные данные и на выходе
функция возвращает результат (функция не должна ничего
вводить с клавиатуры или печатать на консоль).
def main():
# Выполняется, когда файл выполняется как скрипт
. . .
if __name__ == "__main__":
main()
import moduleName
====================================================
26 Oct 2021
Сайт курса:
http://mech.math.msu.su/~vvb/
или
http://158.250.33.65/~vvb/
Прошлое домашнее задание:
1) вычислить число счастливых билетов длины n (четное число).
Не список счастливых билетов, а именно число, причем программа
должна работать по крайней мере для n = 16.
Идея совсем простого алгоритма:
перебираем не номера длины n, а все номера длины k = n//2
(например, если n = 6, то вместо миллиона номеров перебираем
1000 номеров). Для каждой комбинации цифр вычисляем их сумму.
Главная идея: сумма может принимать значения от 0 до k*9.
Создаем список comb количества комбинаций с данной суммой,
размер списка равен k*9 + 1.
comb[sum] = число комбинаций цифр длины k с суммой s.
цикл для всех комбинаций t цифр длины k
| sum = сумма цифр(t)
| comb[sum] += 1
конец цикла
res = 0
цикл для всех элементов x списка comb
| res += x*x
конец цикла
ответ = res
Задача 2 (по теории чисел), 2 варианта:
1) написать функцию carmichaelNumbers(n), которая возвращает
список всех кармайкловых чисел <= n.
Кармайкловы числа -- числа, которые не являются простыми,
но ведут себя как простые в Малой теореме Ферма.
Малая теорема Ферма:
b != 0 (mod p) => b^(p-1) == 1 (mod p)
Число m называется кармайкловым, если
1) m -- составное число;
2) для всякого b, взаимно простого с m (т.е. gdc(b, m) = 1),
b^(p-1) == 1 (mod p)
Первые 3 кармайкловых числа -- это
561, 1105, 1729.
Теорема
m кармайклово <=> выполняются 3 условия:
1) m свободно от квадратов
m = p1 * p2 * . . . *p_k, p_i простые;
2) k >= 3;
3) для всякого i
(p_i - 1) | (m - 1)
Например,
561 = 3 * 11 * 17
(3 - 1) | 560, (11 - 1) | 560, (17 - 1) | 560
Небольшие кармайкловы числа можно найти с помощью следующей программы:
from math import gcd
primeList = [2, 3, 5, 7, 11, 13, 17, 19]
def pseudoPrime(m):
"""m behaves like a prime number
in Fermat's Little Theorem"""
for b in primeList:
if gcd(b, m) != 1:
continue
x = pow(b, m - 1, m)
if x != 1:
return False
return True
def isPrime(m):
"""Check whether m is a prime number"""
if m in (2, 3, 5, 7):
return True
if m <= 10 or m%2 == 0:
return False
d = 3
while d*d <= m:
if m%d == 0:
return False
d += 2
return True
def carmichael(n):
"""Compute a list of Carmichael numbers < n"""
res = []
for m in range(3, n, 2):
if pseudoPrime(m) and (not isPrime(m)):
res.append(m)
return res
Второй вариант второй задачи:
сгенерировать случайное простое число, в двоичной записи которого
ровно n цифр
2^(n-1) <= p < 2^n pi(n) ~ n/ln n
Для небольших n можно воспользоваться
наивным алгоритмом проверки простоты:
import random
def isPrime(m):
"""Check whether m is a prime number"""
if m in (2, 3, 5, 7):
return True
if m <= 10 or m%2 == 0:
return False
d = 3
while d*d <= m:
if m%d == 0:
return False
d += 2
return True
def generatePrime(n = 16):
"""Generate a random prime number
of length n in binary form"""
assert n > 2
lowerBound = 2 ** (n-1)
upperBound = lowerBound*2
while True:
m = random.randrange(lowerBound, upperBound)
m |= 1 # Make m odd
if isPrime(m):
return m
Для того, чтобы эта программа работала для больших
чисел n (например, когда нужно найти 1000-значное простое число),
следует заменить наивный тест простоты (функцию isPrime(m))
на вероятностный тест простоты Миллера--Рабина.
====================================================
Работа со списками
Умное задание списка, математическое задание списка:
list comprehension
В математике, множество простых чисел, меньших 100,
задается как
s = { x : x <- 1, 2, ..., 100, x -- простое число }
Список тех же простых чисел в Питоне задается с использованием
конструкции list comprehension:
s = [ x for x in range(1, 101) if isPrime(x) ]
Задача:
найти все пифагоровы тройки чисел (a, b, c):
a^2 + b^2 = c^2, где числа a, b взаимно простые и a <= b,
в интервале 1..100
Решение с помощью list comprehension:
[ (a, b, c) for a in range(1, 101)
for b in range(a, 101) for c in range(b, 101)
if a*a + b*b == c*c and gcd(a, b) == 1 ]
Результат:
[(3, 4, 5), (5, 12, 13), (7, 24, 25), (8, 15, 17),
(9, 40, 41), (11, 60, 61), (12, 35, 37), (13, 84, 85),
(16, 63, 65), (20, 21, 29), (28, 45, 53), (33, 56, 65),
(36, 77, 85), (39, 80, 89), (48, 55, 73), (65, 72, 97)]
2*3*5 = 30
Найти все обратимые элементы в Z_30:
[ d for d in range(1, 30) if gcd(d, 30) == 1 ]
Результат:
[1, 7, 11, 13, 17, 19, 23, 29]
Эти числа используются в алгоритме разложения числа на множители
"Колесо со спицами"
Самый простой алгоритм разложения на множители состоит в том,
что мы сначала проверяем, делится ли число на 2; если нет,
то мы пытаемся делить число на все пробные нечетные делители.
========
Идея "Колеса со спицами" сходная:
сначала делим число на 2, 3, 5. Если оно не делится ни на одно
из них, то затем делим число на пробные делители вида
30*k + d_i,
где k = 0, 1, 2, ...
для каждого k индекс i = 1, 2, ..., 8,
d <- [1, 7, 11, 13, 17, 19, 23, 29]
Это делители
2, 3, 5,
\1, 7, 11, 13, 17, 19, 23, 29, + 30
31, 37, 41, 43, 47, 49, 53, 59, + 30
61, 67, 71, 73, 77, 79, 83, 89,
. . .
(только надо пропустить делитель единицу при k = 0).
Почему используется название "Колесо со спицами":
представим колесо с длиной обода 30, у которого
из 30 спиц большинство выломано, остались только спицы
1, 7, 11, 13, 17, 19, 23, 29
Оно катится по числовой оси, на которой отмечены целочисленные точки,
и отмечает среди них пробные делители.
Дома:
1) реализовать алгоритм проверки простоты числа типа "Колесо со спицами"
с длиной обода
m = 2*3*5*7*11*13*17*19 = 9699690
Тогда список обратимых элементов кольца вычетов Z_m
[ d for d in range(1, m) if gcd(d, m) == 1 ]
состоит из 1658880 элементов (это значение функции Эйлера phi(m));
2) та же задача для разложения числа на множители.
Студенты с нечетными номерами по Журналу делают первую задачу,
с четными -- вторую.
===================================================================
2 Nov 2021
Классификация языков программирования
Традиционные -- объектно-ориентированные
Императивные -- функциональные
Императивные языки -- есть переменные, циклы, операторы языка --
это команды в повелительном наклонении.
x = x + 1
В функциональных языках нет ни переменных, ни циклов; как бы нет времени.
Вместо переменных есть обозначения функций,
n = 10
n = 20
Главными объектами в функциональных языках есть функции и
списки (списки однонаправленные).
List comprehension -- конструкция, которая пришла из функциональных языков.
[ x**2 for x in range(1, 21) ]
Пифагоровы тройки
l = [ (a, b, c) for a in range(1, 101) for b in range(a+1, 101)
for c in range(a+b+1, 101)
if gcd(a, b) == 1 and a**2 + b**2 == c**2 ]
l = []
for a in range(1, 101):
for b in range(a+1, 101):
for c in range(a+b+1, 101):
if gcd(a, b) == 1 and a**2 + b**2 == c**2:
l.append((a, b, c))
Раньше был популярен Python2, и там range(1, 10) означало список.
Например, когда мне нужно было найти все пары простых чисел-близнецов
(p, q): q = p + 2 в интервале до 1000000, то я написал
for m in range(1000000):
. . .
сразу создается огромный список! Зачем???
Был добавлен тип xrange, а в Питоне3 он вообще выкинут, а range --
это то, что называлось раньше xrange.
Один из главных принципов финкциональных языков -- lasy evaluation,
ленивые вычисления. Мы определяем объекты, но конкретные вычисления
проводим только тогда, когда это становится неизбежным.
s = [1, 3, 5, ..]
В Питоне: также можно делать бесконечные объекты, но они не помещаются целиком
в память (потенциальная бесконечность), а можно по ним итерировать
(перебирать их элементы). Имеется несколько конструкций:
есть класс итератор, у которого есть 2 метода:
конструктор
метод iter, возвращает сам себя
метод __next__() вызывается с помощью функции next(it) === it.__next__()
он возвращает значение, на которое ссылается итератор, при этом сам
итератор переходит к следующему значению структуры данных.
Если итератор уже переместился за конец структуры данных, то
при выполнении next(it) произойдет исключение типа StopIteration.
for x in s:
action(x)
это синтаксический сахар для
i = iter(s)
while True:
try:
x = next(i)
action(x)
except StopIteration:
break
Упрощенный способ создания итераторов:
1) генераторы (generators);
2) генераторные выражения (generator expressions).
Пример генератора: перебираем все целочисленные точки плоскости в первом
найти квадранте, у которых сумма координат не превышает maxNum.
def nodes(maxNum):
n = 0
while n <= maxNum:
x = n; y = 0
while x >= 0:
yield (x, y)
x -= 1; y += 1
n += 1
С помощью этого генератора решим задачу:
найти хорошее приближение числа pi в виде целочисленной дроби.
Ищем приближение с ошибкой, не превышающей единицы, деленной
на квадрат знаменателя дроби.
from math import pi
prec = 1e30 # infinity
m = (3, 1)
for p in nodes(1000):
if p[1] == 0:
continue
d = abs(pi - p[0]/p[1])
if d < 1/(p[1]**2) and d < prec:
m = p
prec = d
print(m)
print("pi is approximated as", m[0], "/", m[1])
Вот результат выполнения этой программы:
(3, 1)
(19, 6)
(22, 7)
(333, 106)
(355, 113)
pi is approximated as 355 / 113
2) генераторные выражения (generator expressions)
такое же выражение, которое используется в list comprehension,
но в круглых скобках:
s = (d for d in range(1, 30) if gcd(d, 30) == 1)
Домашнее задание
Задача:
Написать генератор для всех перестановок порядка n.
Перестановки должны перебираться в лексикографическом порядке.
[1, 2, 3]
[1, 3, 2]
[2, 1, 3]
[2, 3, 1]
[3, 1, 2]
[3, 2, 1]
Замечание:
удобно реализовать вспомогательную функцию nextPermutation(s),
на вход подается список чисел, на выходе этот список изменяется.
Функция возвращает True, если перестановка была не последней,
и False, если следующей перестановки нет.
. . . 9, 12, 10, 8, 7, 5
. . . 10, 5, 7, 8, 9, 12
Необязательная задача:
Написать 2 функции:
1) по перестановке получить ее разложение в произведение
независимых циклов;
2) обратная функция: по списку независимых циклов получить перестановку.
Выписать все перестановки в лексикографическом порядке
в виде произведений циклов.
===========================================================
9 Nov 2021
Обработка текстов и регулярные выражения
Regular Expressions, или re
Коротко о теории формальных языков
Хомский предложил определение формальных грамматик и
языков, задаваемых ими.
Sigma -- конечный алфавит, |Sigma| < inf
Для описания грамматики добавляется алфавит нетерминалов
N -- алфавит метасимволов, или нетерминалов
S -- начальный нетерминал
Контекстно-свободная грамматика задается правилами вида
A --> u
A <- N -- нетерминал, u <- (Sigma U N)* -- любая цепочка из терминалов
и нетерминалов
Задача:
текстовая строка содержит запись формулы:
(2 + 3.5)*10
Надо вычислить значение формулы.
Рассмотрим грамматику арифметических формул с маркером конца строки
(конец формулы отмечается маркером $):
S -> F $
F -> a | F + F | F - F | F * F | F / F | (F)
a + (a*a)/(a - a) $
Это плохая грамматика, поскольку одна и та же формула
может иметь несколько разных левых выводов (или, что эквивалентно,
несколько разных синтаксических деревьев вывода).
Вывод формулы a + (a*a)/(a - a) $
S => F $ => F + F $ => a + F $ => a + F/F $ => a + (F)/F $ =>
=> a + (F*F)/F $ => a + (a*F)/F $ => a + (a*a)/F $ =>
=> a + (a*a)/(F) => a + (a*a)/(F - F) $ =>
=> a + (a*a)/(a - F) $ => a + (a*a)/(a - a) $
Пример двух разный синтаксических деревьев для одной и той же формулы
(т.е. трактовка формулы неоднозначна):
a + a*a $
S S
/ \ / \
F $ F $
/|\ /|\
F + F F * F
| /|\ /|\ |
a F * F F + F a
| | | |
a a a a
Хорошая (однозначная) грамматика формул учитывает приоритеты операций:
формула есть сумма (разность) термов;
терм есть произведение (частное) множителей;
множитель есть либо число, либо любая формула, взятая в скобки.
S -> F $
F -> T | F + T | F - T
T -> M | T * M | T / M
M -> e | ( F )
Строение фраз языка описывается на двух уровнях:
1) лексика, т.е. строение отдельного слова языка;
2) синтаксис, т.е. строение фраз языка.
Компилятор всегда работает в два прохода:
1) сканер, или лексический анализатор, разбивает входной поток
букв на слова (лексемы, токены). Например, в языке Си лексемы --
это целые и вещественные константы, ключевые слова языка,
имена переменных, скобки, знаки операций;
Сканер: поток символов ---> поток лексем
2) парсер, или синтаксический анализатор, получает на вход поток лексем,
на выходе он решает нужную задачу (восстановить вывод фразы, построить
синтаксическое дерево, перевести текст с языка высокого уровня на
язык более низкого уровня и т.п.).
Программисты (вернее, ученые от computer science) научились решать
эффективно обе эти задачи.
Хомский, праволинейные грамматики
A -> uB
A -> v u, v <- Sigma* -- терминальные цепочки.
Праволинейные языки == автоматные языки.
Языки, порожденные конечными автоматами:
(картинки)
Теорема Клини
Следующие классы языков совпадают:
1) праволинейные языки
2) леволинейные языки
3) автоматные языки
4) регулярные языки
Доказательство: все импликации очевидны, кроме 3 => 4
Импликация 3 => 4 хоть и простая, но очень красивая
(доказательство мы здесь не приодим).
Регулярные языки, способ построения
1. Язык, состоящий из пустой цепочки, регулярен
()
2. Язык, состоящий из одноэлементной цепочки, регулярен
a (цепочки из одного символа)
3. L1, L2 регулярны => L1 U L2 регулярно (объединение)
(e1 | e2)
4. L1, L2 регулярны => L1 L2 регулярно
(произведение L1 L2 = { uv, u<-L1, v<-L2 } )
e1 e2
5. L регулярен => итерация Клини L* регулярна
L* = { u1 u2 ... u_k, k >= 0, u_i <- L }
Регулярные языки можно охарактеризовать как цепочки, которые можно
распознать на машине с конечной памятью. Как правило, множество
слов (лексем) любого языка образует регулярный язык.
Задание лексем с помощью конечного автомата удобно для компьютера;
но для человека удобнее задавать лексемы языка в виде регулярных выражений
(потому что это текст, а не картинка).
В любом языке программирования есть набор расширенных
регулярных выражений:
1) дипазон, набор символов
[a-zA-Z] == (a|b|c|...|z|A|B|C|...|Z)
2) необязательная последовательность
e? == (()|e)
3) повторение 1 или более раз
e+ == e e*
4) повторение от m до n раз
e{m,n}
5) \d -- любая десятичная цифра == [0-9]
\s -- white space (пробельный символ)
\w -- либо буква, либо цифра == [a-zA-Z_0-9]
^e -- цепочка e в начале строки (или в начале текста)
e$ -- цепочка e в конце строки (или в конце текста)
Дома:
1) написать регулярное выражение для поиска всех целочисленных констант
в тексте (они могут быть десятичные, восьмеричные, шестнадцатеричные
и двоичные);
2) написать регулярное выражение для поиска всех вещественных констант
в тексте (не целых констант!)
0.5, 12.301, .5, 0., 100., 1.5e+30, 10e30, 1E-100, .5e3
Константа вещественная, если либо
1) в ней есть десятичная точка, при этом либо целая, либо дробная
часть должна быть непустой, десятичный порядок может быть
или отсутствовать;
2) либо в ней нет десятичной точки, зато обязательно есть порядок.
1) [0-9]+\.[0-9]*((e|E)(\+|\-)?[0-9]+)? |
[0-9]*\.[0-9]+((e|E)(\+|\-)?[0-9]+)? |
2) [0-9]+(e|E)(\+|\-)?[0-9]+
===========================================================
small_primes = [2, 3, 5, 7, 11, 13, 17, 19]
m = 1
for p in small_primes:
m *= p
# divizors = [ d for d in range(1, m) if gcd(d, m) == 1 ]
divizors = [] # Список пробных делителей -- глобальный объект, который
# заполняется один раз при первом обращении к функции
def spoked_wheel(m):
global small_primes, m, divizors
if len(divizors) == 0:
divizors = [ d for d in range(1, m) if gcd(d, m) == 1 ]
. . .
Бывает также nonlocal
========================================
16 Nov 2021
Задача 1.
Написать регулярное выражение для поиска всех целочисленных констант
(в смысле Python3) в тексте.
0, 1203, (запрещено 0377), 0xFF, 0x00ab,
0o377 (octal), 0b11111111 (binary)
Задача 2. То же самое для вещественных (float) констант.
0.5, 123.25, 0.5e+6, 123.44E30, 12.5e-10,
.5, 0., 25.
1e10, 1e-15, 100e3
Указание:
Константа вещественная, когда выполняется одно из трех условий:
1) целая часть непустая, есть десятичная точка, дробная часть
может быть пустой, порядок необязателен;
123.e+10, 123.456E-20
2) целая часть может быть пустой, есть десятичная точка, дробная часть
непустая, порядок необязателен;
.123, .5e-10, 12.34e30
3) целая часть непустая, десятичная точка отсутствует,
порядок обязателен
1e10, 25E-6
Регулярные выражения и регулярные языки
1) L = { eps } ()
2) L = { a } a
3) L1 U L2 (e1|e2)
4) L1 L2 (e1 e2)
5) L* итерация Клини (e*)
L* = { u1 u2 ... u_k, k >= 0 }
eps <- L*
Проект "Компилятор формул"
(2 + 3)*10 - 5 RPN 2, 3, +, 10, *, 5, -
= 45
scanner
sin(1/x) ---> Func LPAR NUMBER DIV X RPAR ---> RPN
(Reverse Polish Notation)
1, x, /, sin
scanner parser интерпретатор
text ---> список лексем ---> обратная польская запись ---> значение формулы
(2 + 3) * 10 $
(NUMBER, "2", 2.0)
===========================================================
23 Nov 2021
Прошлый раз был реализован модуль scanner для языка арифм. формул
Формулы включают числа, знаки арифм. операций, круглые скобки
и имена стандартных функций
sin(1)*cos(2*x + 3) - 17*(sin(x) + 10)
Компилятор формул принимает на вход цепочки символов с маркером
конца строки
sin(1)*cos(2*x + 3) - 17*(sin(x) + 10) $
Сканер преобразует поток символов в поток лексем (или токенов).
Мы хотим реализовать парсер (синтаксический анализатор), который
формулу в инфиксной записи преобразует в формулу в обратной польской записи.
Ян Лукасиевич предложил два способа записи формул, не использующие скобки
для указания порядка операций. В обычной (инфиксной) записи в формуле
знаки оперций стоят между аргументами
(2+3)*(17-6)
В прямой польской записи знаки ставятся перед аргументами:
*, +, 2, 3, -, 17, 6
Но наиболее удобна обратная польская запись (RPN -- Reverse Polish Notation),
в которой знак операции указывается после аргументов:
2, 3, +, 17, 6, -, *
1) 2) 3) 4) 5) 6) 7)
Для вычисления значения формулы в RPN используется стековый вычислитель:
0)
2 1) 2
3 2) 2, 3
+ 3) 5
17 4) 5, 17
6 5) 5, 17, 6
- 6) 5, 11
* 7) 55
Обратная польская запись, как правило, используется в объектно-ориентированных
языках для представления промежуточного кода. Например, компилятор Java
переводит тест программы на языке Java в так называемый байткод,
который продставляет собой просто обратную польскую запись программы.
То же самое для Питона, для C# (.Net, основа которой -- язык CIL,
Common Intermediate Language), CIL компилируется "на лету" (JIT -- Just In Time
Compiling). Есть еще язык PostScript.
Наша цель -- создать класс LLParser, который преобразует поток лексем,
составляющих формулу, в обратную польскую запись.
Задание языков с помощью контектстно-свободной грамматики (Хомский).
Плохая грамматик языка арифм. формул
S -> F $
F -> e | F + F | F * F | (F)
e + e * e $
S => F $ => F + F $ => e + F $ => e + F * F $ => e + e * F $ => e + e * e $
e + (e * e)
S => F $ => F * F $ => F + F * F => e + F * F => ... => e + e * e $
(e + e) * e
Преобразуем грамматику, учитывая приоритеты операций:
S -> F $
F -> T | F + T | F - T
T -> M | T * M | T / M
M -> e | ( F )
Эта грамматика однозначная (любая выводимая цепочки имеет единственное
дерево вывода), причем она учитывает приоритет операций:
мультипликативные операции имеют более высокий приоритет, причем
операции одинакового приоритета выпоняются слева направо.
1 - 2 - 3 ==== (1 - 2) - 3, а не 1 - (2 - 3)
Добавим в наш язык функции от одного аргумента и переменную x
S -> F $
F -> T | F + T | F - T
T -> M | T * M | T / M
M -> E | Func
E -> x | number | ( F )
Func -> name ( F )
Как пишутся парсеры? Есть два популярных алгоритма разбора:
1) рекурсивный, или нисходящий разбор, или LL(1)-разбор;
L -- читаем цепочку слева направо,
L -- строим левые выводы,
1 -- учитываем 1 символ аванцепочки
2) восходящий разбор, или разбор с помощью детерминированного
конечного автомата со стеком (или магазином), или LR(1) -- разбор
L -- читаем цепочку слева направо,
R -- строим правые выводы,
1 -- учитываем 1 символ аванцепочки
Все современные компиляторы используют LR(1)-разбор, поскольку
класс языков, допускающих LR-разбор, намного шире.
Даже если язык допускает рекурсивный разбор, то его рекурсивная грамматика
намного менее естественна для человека, чем грамматика того же языка,
допускающая LR-разбор.
Но: строить вручную таблицы автомата с магазинной памятью -- занятие
не для лентяев. Поэтому умные (ленивые) люди написали программу, которая
строит такой парсер прямо по записи грамматики -- это
yacc (Yet Another Compiler Compiler), или bison -- свободная версия.
Мы напишем рекарсивный (LL(1)) парсер.
Идея рекурсивного парсера
Пусть есть правило
A -> aBCdE | FBC
Ему соответствует функции
processA(), которая из входной цепочки обрабатывает
начальный кусок, который выводится из нетерминала A.
void processA() {
skipToken(a);
processB();
processC();
skipToken(d);
processE();
}
Грамматика
S -> F $
F -> T | F + T | F - T
T -> M | T * M | T / M
M -> E | Func
E -> x | number | ( F )
Func -> name ( F )
не допускает рекурсивный разбор. Однако в ней можно избавиться от
левой рекурсии. Общий способ:
пусть для нетерминала A есть правила
A -> Au1 | Au2 | ... Au_k | w1 | w2 | ... | w_s
A => Au3 => Au1 u3 => A u u u u u => w u u u u u u u
Упрощенно:
A -> Au | w
A => Au => Auu => Auuu...u => wuuu...u
Эти правила для A можно заментить на
A -> w Atail
Atail -> eps | u Atail
Применив этот метод к нашей грамматике, получим:
S -> F $
F -> T Ftail
Ftail -> eps | + T Ftail | - T Ftail
T -> M Ttail
Ttail -> eps | * M Ttail | / M Ttail
M -> E | Func
E -> x | number | ( F )
Func -> name ( F )
Дома:
1) реализовать метод evaluate(x) парсера, используя стековый вычислитель;
2) добавить унарный минус в грамматику (операция изменения знака);
3) добавить константы pi, e
4) добавить операцию возведения в степень,
которая обозначается либо как в Питоне **, либо как в математике
и в SageMath "^"
sin(x)^2 + cos(x)^2
(последний пункт требует изменений и в сканере, и в парсере).
============================================================
7 Dec 2021
Надо сделать дома задание по проекту Scanner + Parser + Plotter
1. Элементы теории вероятности
Рассмотрим последовательность чисел, записанную в файле.
Надо вычислить
1) среднее арифметическое, или математическое ожидание случайной величины;
2) среднее геометрическое;
3) среднее гармоническое -- обратная величина к среднему из обратных
величин последовательности
s = [1, 2, 3, 4]
mean = sum[s]/len(s)
meanGeom = product[s] ** (1/len(s))
sinv = [1/x for x in s]
meanHarm = 1/(sum(sinv)/len(sinv))
Дисперсия (variance):
D = E((x - E(x))^2) -- мера уклонения случайной величины x от ее среднего
значения, размерность в метрах в квадрате.
Физический смысл имеет квадратный корень из дисперсии:
standard deviation = средне-квадратичное уклонение = sqrt(D)
Как вычислить дисперсию за 1 проход (за один просмотр последовательности)?
Это можно сделать. Просматривая последовательность, достаточно хранить
3 величины:
1) число элементов (длина просмотренной части последовательности) n;
2) сумма элементов (первый момент случайной величины) S1;
3) сумма квадратов элементов (второй момент) S2.
E = сумма(x_i)/n = S1/n
D = сумма((x_i - E)^2)/n = сумма(x_i^2)/n - 2*E*сумма(x_i)/n +
+ сумма(E^2)/n = S2/n - 2*E^2 + E^2 = S2/n - E^2 =
= S2/n - S1^2/n^2
f0 = S2/n - S1^2/n^2
Формула f0 по некоторой выборке дает оценку значения дисперсии (или
какой-то другого параметра случайной величины). Оценка называется
состоятельной, если lim E(f0) = D. Оценка называется несмещенной,
n -> inf
если E(f0) = D. Приведенная оценка для дисперсии является состоятельной,
но не является несмещенной.
Несмещенная оценка для дисперсии при n > 1:
f1 = сумма((x_i - E)^2)/(n - 1) =
= f0 * n / (n - 1) = (S2/n - S1^2/n^2) * n / (n - 1) =
= S2/(n - 1) - S1^2/(n*(n - 1))
Проблема: это несмещенная оценка для дисперсии, но когда мы извлекаем
квадратный корень из дисперсии, оценка вновь становится смещенной. Вроде
бы для sqrt(D) надо делить на (n - 1.5):
f1 = сумма((x_i - E)^2)/(n - 1.5)
sd = sqrt(f1) ???????????????
Попробуйте экспериментальным путем определить оптимальный параметр s
f1 = сумма((x_i - E)^2)/(n - s)
sd = sqrt(f1) ???????????????
при котором мы получаем наиболее близкое значение к sd при n = 10 при n = 100
для двух случайных величин:
1) равномерное распределение от -1 до 1;
2) нормальное распределение с sigma = 1.
----------------------------------------------------------
Модуль numpy -- реализация на Питоне многомерных массивов.
Вопрос: чем плохи обычные списки в Питоне?
import numpy as np
a = [[1, 2], [3, 4]]
b = np.array([[1, 2], [3, 4]])
Разница между списком в Питоне и объектом типа numpy.ndarray в том,
что второй хранит элементы матрицы в виде линейного массива (так же,
как в языке С).
b_00, b_01, b_10, b_11
numpy.array хранит только те типы данных, которые реализованы в компьютере
(или являются базовыми типами языка С): например, int64, но не int
в смысле Питона.
Реализованы векторные операции с (многомерными) матрицами:
сложение (строк) матрицы, произведение матриц и т.д.
Реализован ли в модуле numpy метод Гаусса приведения матрицы к ступенчатому
виду? По-английски он называется row echelon form, или
reduce to row echelon form. (Такой метод есть в модуле sympy.)
Домашнее задание: реализовать метод Гаусса
def Gauss(a):
""" a -- matrix as a 2-dimensional numpy.array
return row echelon form of a"""
Программистский метод Гаусса отличается от математического тем,
что
1) надо выбирать в качестве разрешающего элемента максимальный
по модулю элемент в столбце;
2) для инвариантности определителя
вместо обмена строк можно при этом менять знак одной из строк.
Обмен строк:
+-------------+
| *?????? |
| 0 *??????? |
| 000*------+ i
| 0| |
| 0| |
| 0 | |
| | | k a[[i, k]] = a[[k, i]] обмен строк i, k
| | |
| 0| |
+-------------+
j
Задачи:
1) дан набор векторов. Найти среди них максимальную линейно независимую
систему;
2) дана вырожденная матрица. Вычислить базис в пространстве решений
системы A*x = 0.