Текст книги "Рассказы о математике с примерами на языках Python и C (СИ)"

Автор книги: Дмитрий Елисеев

сообщить о нарушении

Текущая страница: 4 (всего у книги 4 страниц)

15. Распределение случайных величин

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

Равномерное распределение

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

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

Нормальное распределение

Возьмем группу людей, например в 100 человек, и измерим их рост. Очевидно, что будет некоторое количество людей небольшого роста, некоторое количество высоких людей, совсем мало очень высоких, и совсем мало очень низких. Такое распределение естественно для многих объектов, не только людей, потому оно и называется нормальным.

Формула нормального распределения совпадает с формулой Гаусса:

Подбирая коэффициенты, можно получить разные виды распределения.

Касаемо роста людей, согласно сайту http://tall.life, график роста для мужчин и женщин имеет следующий вид:

Распределение Пуассона

Следующий вид распределения не менее интересен. Рассмотрим события, происходящие с некоторой известной интенсивностью независимо друг от друга, например приход покупателей в магазин. Допустим, в магазин приходит в среднем 10 покупателей в минуту. Какова вероятность, что в какой-то момент времени в магазин придет 20 покупателей?

Вероятность таких событий описывается распределением Пуассона:

График распределения имеет примерно такой вид (в нашем примере λ = 10):

Этим же распределением описываются различные случаи, от вероятности неисправностей (если 0,01% телевизоров имеют неисправность, какова вероятность что в партии из 20 штук окажется 2 неисправных телевизора), до скорости роста колоний в чашке Петри.

Вернемся к нашему примеру с 20 покупателями. В интернете можно найти таблицы значений Пуассона для λ=10. По ним можно найти, что вероятность прихода сразу 20 человек составляет 0,19%.

16. Измерение скорости света

С бытовой точки зрения, скорость света практически мгновенна. Действительно, свет за секунду может обогнуть Землю 8 раз, а за 2 секунды пролетает расстояние от Земли до Луны. Поэтому до 17-го века про реальную скорость света никто не знал. Как же ее вычислили?

Сегодня опыт по измерению скорости света можно провести даже в школе – достаточно длинного куска кабеля, генератора импульсов и осциллографа. Действительно, задержка сигнала в куске кабеля длиной 50 м, будет равна 50/300000000, или 0,016 мкс – величина которую покажет даже дешевый осциллограф с максимальной частотой 10-20 МГц. Но как же это сделали в 17-м веке, когда не было не то что осциллографов, даже до появления лампы накаливания было еще 200 лет ожидания? Помогли астрономия, геометрия, и разумеется, математика.

Говоря точнее, помогло наблюдение Юпитера и его спутников. Спутники Юпитера были открыты еще Галилеем, увидеть их может каждый, даже с балкона в небольшой телескоп. С увеличением около 300х они видны примерно так:

Период вращения спутников Юпитера невелик, и составляет примерно 2 дня. Уже в 17-м веке измерение времени было достаточно точным (маятниковые часы изобрел голландский физик и математик Гюйгенс в 1657 г.), чтобы датский астроном Олаф Ремер в 1676 году обнаружил расхождение расчетного и реального положения спутника примерно в 16 минут (величина, которую трудно не заметить даже при технологиях 17 века).

Для измерения орбит спутников Юпитера Ремер использовал момент, когда спутник входит в тень Юпитера – момент, который можно измерить довольно-таки точно. Как догадался Ремер, запаздывание во времени было связано с движением Земли по орбите.

Картинка с сайта http://www.speed-light.info:

В момент второго измерения расстояние до Юпитера больше примерно на диаметр орбиты Земли (период обращения Юпитера вокруг Солнца – 12 лет, что гораздо больше земного). Это и приводило к тому, что свет от Юпитера приходил с большим запаздыванием, чем при первом измерении. Сделав подсчеты, Ремер получил значение скорости света в 220000 км/c. В то время конечность скорости света казалась настолько невероятной, что после публикации во французской академии наук далеко не все поверили молодому ученому. Разумеется, последующие измерения подтвердили правильность метода.

Более точное значение было получено лишь через 200 лет, французский физик Луи Физо с помощью зубчатого колеса и двух зеркал получил значение в 312000 км/c. Расстояние между зеркалами было 8,6 км, одно зеркало было расположено в доме отца Физо недалеко от Парижа, второе зеркало было расположено на Монмартре. Физо нашел такую скорость вращения колеса, при котором луч света проходящий через зубец колеса затемнялся, что означало что запаздывание света соответствует скорости вращения колеса.

17. Можно ли своими глазами увидеть прошлое?

С предыдущим вопросом связан простой забавный факт. Можно ли лично своими глазами увидеть событие, происходящее например, миллион лет назад? Да, и это очень просто сделать – достаточно посмотреть на небо.

Как наверное читатели уже догадались, все дело в скорости света. Расстояние от Земли до Солнца составляет 150 млн. километров, свет преодолевает его за 8 минут. Таким образом, когда мы смотрим на небо, мы видим Солнце так, как оно было 8 минут назад. Смотря на Юпитер, мы видим его с запаздыванием примерно в полчаса. Расстояние до других звезд гораздо больше. Сириус мы видим таким, каким он был 8 лет назад, звезда Вега имеет расстояние в 25 световых лет, столько нужно свету чтобы долететь до Земли. Туманность Андромеды мы видим с запаздыванием в 2,5 миллиона лет – время сопоставимое с периодом жизни мамонтов и первобытных людей. Если бы там висело гигантское зеркало, в котором отражалась бы Земля, гипотетически можно было бы их увидеть, хотя такого размера телескоп конечно же физически невозможен.

А если серьезно, то данный факт запаздывания света очень помогает астрономам изучать историю Вселенной – наблюдая удаленные галактики, находящиеся на расстоянии миллиардов световых лет, мы видим их в прошлом – их свет шел до нас эти самые миллиарды лет. И таким образом, чем мощнее телескоп, тем дальше во времени он позволяет заглянуть, тем самым, позволяя видеть то время, когда и галактики и звезды еще были молодыми.

С другой стороны, в запаздывании света есть и большой минус. Даже если будут созданы космические корабли, способные долететь до другой звезды, обмениваться сообщениями с космонавтами придется с тем же запаздыванием. К примеру, радиосообщение до Сириуса дойдет до получателя через 8,6 лет, и столько же придется ждать ответа. Уже сейчас теоретически можно поговорить по телефону с астронавтами на МКС (в 2015 году британский астронавт Тим Пик ошибся номером, и удивил неизвестную женщину вопросом «Здравствуйте, это Земля?»), а вот для Марса время задержки составит около 15 минут – так что поговорить по телефону или по Скайпу с марсианской колонией было бы невозможно.

18. Сколько вольт в электросети?

Глупый вопрос, подумают многие. Каждый школьник в России знает что напряжение в сети 220 вольт (в США каждый школьник знает что напряжение в сети 110 вольт). Полезно привести такую картинку:

Кстати, в 90-е годы, когда поездки за границу только становились доступными, некоторые привозили американскую электронику, но работала она зачастую не долго, из-за того что сетевое напряжение отличается в 2 раза. А сейчас даже чуть больше, по российскому стандарту 2003 года, напряжение в сети должно составлять 230 В. Предельно допустимым отклонением считается 10%, т. е. Значения 210-250 В в принципе возможны.

Но вопрос заголовка не в этом. Будем для простоты считать напряжение равным «условным» 220 вольт. Однако подключим осциллограф к электросети, и увидим примерно такую картинку:

Что это значит? Где «наши» 220 вольт?

Всё просто (хотя и не совсем). Ток в сети переменный – он меняет свое направление с частотой 50 раз в секунду. В отличие к примеру, от батарейки – если на ней написано 1,5 вольта, это значит что на ней действительно 1,5 вольта и направление тока не меняется. Но вернемся к розетке. Ток в нее подается не просто так, а с целью выполнения какой-либо работы. Как измерить работу переменного тока, который в разные моменты времени движется то в одну, то в другую сторону? Для этого было введено понятие действующего напряжения – величины постоянного тока, способного выполнить ту же работу (например нагрев спирали электроплитки). Напряжение, которое показывает осциллограф – называется амплитудным. Эти величины связаны простой формулой:

220 умноженное на √2, дает как раз 310 В. Разумеется, обычный тестер откалиброван в «бытовых» единицах, в режиме измерения переменного тока он покажет 220 В. А если выпрямить напряжение, например диодным мостиком, то тестер покажет как раз 310 В постоянного тока.

И еще немного про переменный ток. Откуда берется напряжение в 380 вольт? Ток от трансформатора подается по 3-м фазам: это 3 линии, напряжение в которых сдвинуто на разный угол друг относительно друга.

Картинка из Википедии:

Нулевой провод – общий. В квартиры подается напряжение с одной из фаз, значением в стандартные 220 вольт. Это напряжение называется фазным. Если же используется 3-фазная сеть целиком, то напряжение между двумя фазами, например в точках a и c на рисунке, составляет как раз 380 вольт. Это напряжение называется линейным.

Математически, оба напряжения связаны простой формулой:

Действительно, 220 * √3 = 380.

Кстати, обрыв нулевого провода в доме – серьезная неисправность, из-за чего в квартиры может быть подано линейное напряжение, составляющее те самые 380 В. Такой случай произошел лично с автором, впрочем ущерб оказался невелик, перегорели лишь настенные электронные часы и несколько блоков питания. Но при отсутствии в доме людей это может привести и к пожару, такие случаи не редкость. Так что тем, у кого в квартире старая проводка, рекомендуется установить в электрощиток устройство защиты от перенапряжения, его цена невелика, и явно дешевле ремонта в квартире.

19. Приложение 1 – Вычисления с помощью видеокарты

Еще 20 лет назад, во времена процессоров 80386, пользователям приходилось покупать математический сопроцессор, позволяющий быстрее выполнять вычисления с плавающей точкой. Сейчас такой сопроцессор покупать уже не надо – благодаря прогрессу в игровой индустрии, даже встроенная видеокарта компьютера имеет весьма неплохую вычислительную мощность. Например, даже бюджетный видеочип Intel Graphics 4600 имеет 20 вычислительных блоков, что превышает количество ядер «основного» процессора. Разумеется, каждое ядро GPU по отдельности слабее CPU, но здесь как раз тот случай, когда количество дает преимущество над качеством. Вычисления с помощью GPU сейчас очень популярны – от майнинга биткоинов до научных расчетов, диапазон ценовых решений также различен, от «бесплатной» встроенной видеокарты до NVIDIA Tesla ценой более 100 тыс. рублей. Поэтому интересно посмотреть, как же это работает.

Есть две основные библиотеки для GPU-расчетов – NVidia CUDA и OpenCL. Первая обладает большими возможностями, однако работает только с картами NVIDIA. Библиотека OpenCL работает с гораздо большим числом графических карт, поэтому мы рассмотрим именно ее.

Основной принцип GPU-расчетов – параллельность вычислений. Данные, хранящиеся в «глобальной памяти» (global & constant memory) устройства, обрабатываются модулями (каждый модуль называется «ядром»), каждый из которых работает параллельно с другими. Модуль имеет и свою собственную память для промежуточных данных (private memory). Так это выглядит в виде блок-схемы:

Таким образом, если задача может быть разбита на небольшие блоки, параллельно обрабатывающие небольшой фрагмент блока данных, такая задача может эффективно быть решена на GPU.

Рассмотрим пример: необходимо проверить, какие числа в массиве являются простыми. Массив может быть большим, например миллион элементов. Такая задача идеальна для распараллеливания: каждое число может быть проверено независимо от предыдущего.

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

1. Написать код микроядра (kernel):

Этот код будет запускаться непосредственно на графических процессорах видеокарты. Код пишется на языке C. В данном примере мы для упрощения храним код прямо в виде строки в программе.

const char *KernelSource = «n»

«__kernel void primes(                      n»

«  __global unsigned int* input,            n»

«  __global unsigned int* output)           n»

«{                                          n»

«  unsigned int i = get_global_id(0);       n»

« //printf(»Task-%d\n", i);              n"

« output[i] = 0;                            n»

«  unsigned int val = input[i];             n»

«  for(unsigned int p=2; p<=val/2; p++) {   n»

«    if (val % p == 0)                      n»

«       return;                             n»

« }                                         n»

« output[i] = 1;                            n»

«}                                          n»

«n»;

Суть кода проста. Массив input хранит числа, которые нужно проверить, функция get_global_id возвращает индекс задачи, которую выполняет данное ядро. Мы берем число с нужным индексом, проверяем его на простоту, и записываем 0 или 1 в зависимости от результата, в массив output.

2. Инициализировать подготовку вычислений:

int gpu = 1;

clGetDeviceIDs(NULL, gpu ? CL_DEVICE_TYPE_GPU : CL_DEVICE_TYPE_CPU, 1, &device_id, NULL);

cl_context context = clCreateContext(0, 1, &device_id, NULL, NULL, &err); cl_command_queue commands = clCreateCommandQueue(context, device_id, 0, &err);

На этом этапе можно выбрать где будут производиться вычисления, на основном процессоре или на GPU. Для отладки удобнее основной процессор, окончательные расчеты быстрее на GPU.

3. Подготовить данные:

#define DATA_SIZE 1024

cl_uint *data = (cl_uint*)malloc(sizeof(cl_uint) * DATA_SIZE);

cl_uint *results = (cl_uint*)malloc(sizeof(cl_uint) * DATA_SIZE);

4. Загрузить данные и программу из основной памяти в GPU:

cl_program program = clCreateProgramWithSource(context, 1, (const char **) & KernelSource, NULL, &err);

clBuildProgram(program, 0, NULL, NULL, NULL, NULL);

cl_kernel kernel = clCreateKernel(program, «primes», &err);

cl_mem output = clCreateBuffer(context, CL_MEM_WRITE_ONLY, sizeof(cl_uint) * count, NULL, NULL);

clEnqueueWriteBuffer(commands, input, CL_TRUE, 0, sizeof(cl_uint) * count, data, 0, NULL, NULL);

clSetKernelArg(kernel, 0, sizeof(cl_mem), &output);

clGetKernelWorkGroupInfo(kernel, device_id, CL_KERNEL_WORK_GROUP_SIZE, sizeof(local), &local, NULL);

5. Запустить вычисления на GPU и дождаться их завершения:

global = DATA_SIZE;

clEnqueueNDRangeKernel(commands, kernel, 1, NULL, &global, &local, 0, NULL, NULL);

clFinish(commands);

6. Загрузить результаты обратно из GPU в основную память:

clEnqueueReadBuffer( commands, output, CL_TRUE, 0, sizeof(cl_uint) * count, results, 0, NULL, NULL );

7. Освободить данные:

free(data);

free(results);

clReleaseMemObject(input);

clReleaseMemObject(output);

clReleaseProgram(program);

clReleaseKernel(kernel);

clReleaseCommandQueue(commands);

clReleaseContext(context);

Как можно видеть, процесс довольно-таки громоздкий, но оно того стоит. Для примера, проверка простоты 250000 чисел заняла на процессоре Core i5 около 6 секунд. И всего лишь 0,5 секунд заняло выполнение вышеприведенного кода на встроенной видеокарте. Для дешевого нетбука с процессором Intel Atom этот же код выполнялся 34 секунды на основном процессоре, и 6 секунд на GPU. Т. е. разница весьма прилична.

Разумеется, еще раз стоит повторить, что «игра стоит свеч» лишь в том случае, если задача хорошо распараллеливается на небольшие блоки, в таком случае выигрыш будет заметен.

Владельцы видеокарт NVIDIA (особенно игровых и достаточно мощных) могут также посмотреть в сторону библиотеки NVIDIA CUDA, расчеты с ее помощью должны быть еще быстрее.

20. Приложение 2 – Быстродействие языка Python

Язык Python очень удобен своей краткостью и лаконичностью, возможностью использования большого количества сторонних библиотек. Однако, один из его минусов, который может быть ключевым для математических расчетов – это быстродействие. Python это интерпретатор, он не создает exe-файл, что разумеется, сказывается на скорости выполнения программы.

Рассмотрим простой пример: рассчитаем сумму квадратов чисел от 1 до 1000000. Также выведем время выполнения программы.

Программа на языке Python выглядит так:

import time

start_time = time.time()

s = 0

for x in range(1,1000001):

    s += x * x

print(«Sum={}, T={}s».format(s, time.time() – start_time))

Результаты работы:

Sum = 333333833333500000, T = 0.47s

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

import time

start_time = time.time()

l = range(1000001)

s = sum(x * x for x in l)

print(«Sum = {}, T = {}s».format(s, time.time() – start_time))

Результаты работы:

Sum = 333333833333500000, T = 0.32s

Быстрее, но лишь чуть-чуть. К тому же, данный код хранит весь массив в памяти, что неудобно.

И наконец, призываем «тяжелую артиллерию»: напишем программу на языке C. Код выглядит так:

#include

#include

int main()

{

  clock_t start = clock();

  unsigned long long int sum = 0, i;

  for(i=1; i<1000001; i++) {

    sum += i*i;

  }

  clock_t end = clock();

  printf(«Sum = %llu, T = %fs», sum, (float)(end – start)/CLOCKS_PER_SEC);

  return 0;

}

Как можно видеть, он ненамного сложнее python-версии. Перед запуском программы, ее надо скомпилировать, выполнив команду C:GCCbingcc.exe "Appendix-2 - speedTest.c" -o"Appendix-2 – speedTest". Результат очевиден: T = 0,007 секунд. И еще чуть-чуть: добавляем флаг оптимизации по скорости, выполнив команду C:GCCbingcc.exe «Appendix-2 – speedTest.c» -o"Appendix-2 – speedTest" -O3. Результат: 0,0035 секунд, разница в быстродействии более 100 раз!

Увы, в более сложных задачах такого прироста реально не бывает (в последнем примере очень короткий код, который видимо полностью помещается в кеш-памяти процессора), но на некоторое улучшение быстродействия можно рассчитывать. Хотя переписывание программы – это крайний случай, сначала целесообразно поискать стандартные библиотеки, которые возможно уже решают данную задачу. К примеру, следующий код на языке Python, вычисляет сумму элементов массива за 0.1 с:

a = range(1000001)

s = 0

for x in a:

    s += x

print(s)

Можно использовать встроенную функцию sum:

a = range(1000001)

s = sum(a)

print(s)

Данный код выполняется за 0,02 секунды, т. е. в 5 раз быстрее первого варианта. Но разумеется, если заранее известно, что задача состоит в обработке большого набора чисел (например поиск простых чисел или магических квадратов), то может быть более целесообразным сразу писать программу на С или С++, в принципе это не намного сложнее, а работать программа будет быстрее.

Заключение

На этом данная книга закончена, хотя надеюсь, что не навсегда – по возможности и по мере появления новых идей, новые главы будут дописываться. Автор надеется, что хоть немного удалось познакомить читателей с увлекательным миром математики и программирования.

Продолжение следует.

Обо всех найденных неточностях или дополнениях просьба писать на электронную почту [email protected]. Наличие новой версии можно проверить на странице http://dmitryelj.spb.ru/math.htm.