Программа курса:
-- язык программирования C/C++;
Ассемблер -- символическая запись команд компьютера в форме,
удобной для человека.
add ax, cx
Главный небостаток Ассемблера -- зависимость от (архитектуры)
компьютера.
Языки высокого уровня -- независимые от компьютера.
Транслятор (компилятор) -- программа, которая переводит запись другой
программы с языка высокого уровня на язык машинных команд.
Самым удачным языком высокого уровня оказался Fortran (Formula Translation)
y = sin(alpha)*r
Algol-60 -- более структурированный язык (понятия блока)
begin
любой код
end
Проблема Алгола -- был придуман "теоретиками", не очень хорошо
знакомыми с практикой программирования. Алгол никогда не был полностью
реализован в соотвествии со стандартом, использовались языки типа
"Алгамс" -- частичная реализация Алгола-60, зависящая от компьютера
и компилятора.
Никлаус Вирт -- на основе Алгола придумал язык Pascal.
Проблема с Pascal -- изначально вся программа должна была быть написана
в виде одного файла (одного текста). Реальные программные продукты
состоят из тысяч файлов и сотен тысяч или даже миллионов строк, их
невозможно написать и отладить в виде одного файла.
В Pascal была добавлена конструкция Unit.
Но это сделал не Н.Вирт, в результате получились разные версии Pascal.
Сам Н.Вирт придумал как развитие Паскаля язык Modula-2.
Затем Вирт в процессе написания компилятора для Modula-2
решил выкинуть из языка все не очень нужные конструкции, которые
препятствуют его эффективной реализации, и получился язык Oberon.
Мы не будем использовать Паскаль и его производные языки
Modula-2, Oberon, Delphy. Паскаль остался как язык для обучения в
школе и некоторых ВУЗах, он на практике почти не используется.
Как появился язык C? (Видимо, были языки A и B.)
Изначально он был предназначен для замены Ассемблера для написания
операционной системы и компиляторов. Язык C создавали не теоретики,
а практики, хорошо знакомые с операционными системами и Ассемблером.
Была написана операционная система Unix полностью на C
grep, lex, yacc (Yet Another Compiler Compile).
Также Internet возник как сетевая архитектура для компьютеров Unix
(середина 1970-х годов). Internet стала очень популярна после
1989 г., когда в Церне был начат проект WWW (Всемирная паутина).
Операционные системы:
Unix-подобные системы: Linux, Free BSD, ...
Apple-Mac: Mac OS-X (сейчас просто OS-X)
----
MS Windows 1992 (Xenix -- Unix от Microsoft 1980-е)
До Win NT (1992) была Window 3.1 (16-разрядная система
без виртуальной памяти).
Как выяснилось, язык C вполне заменяет и язык Fortran, т.е. C/C++
можно использовать и для научных и инженерных задач.
В 1985 г. Бьярне Страуструп создал объектно-ориентированный язык
C++ на основе C.
Современное состояние:
есть C/C++, на которых пишутся быстрые и надежные программы,
операционные системы, утилиты, драйверы устройств и т.п.
Есть чисто объектно-ориентированные языки: Java, C# (как часть .Net
от Microsoft), Python3. Однако модули для Python пишутся на C++.
===========
Всё сказанное -- оправдание того, что мы выбираем для изучения
язык C/C++.
Я буду писать программы на языке C++ (используя компилятор C++),
оставаясь по сути в рамках языка C (пока не используя классы,
потоковый ввод-вывод, ссылки и др. конструкции из C++, которых нет в C).
Изучение языка
От вас требуется:
-- компьютер с ОС MS Windows (желательно >= 10), Apple c OS-X
или Linux;
-- компилятор gcc/g++;
-- хороший текстовый редактор, в MS Windows вместо системного Notepad
надо использовать свободный редактор Notepad++;
-- отладчик gdb (GNU Debugger);
-- в MS Windows надо установить систему CygWin, которя реализует
все утилиты Unix по Windows -- в частности, оболочка bash
(командная строка вместо cmd и PowerShell от MS Windows),
gcc/g++ и make.
OS-X: gcc/g++, текстовый редактор TextWrangler, терминал
ПЕРЕРЫВ 10 минут до 14:30
GNU -- Gnu is Not Unix
========================================
Язык C
Основная книга по C:
Керниган, Ритчи
Язык программирования C (издание 1989 г. или более позднее, ANSI C, или C-90).
По C++:
Б. Страуструп
Язык программирования C++ (3-е издание 1998 г. или более позднее).
Керниган, Ритчи придумали способизучения языка:
сразу написать минимальную программу, которая что-то делает:
печатает на консоли "Здравствуй, Мир!".
#include <stdio.h>
int main() {
printf("Hello, World!\n");
return 0;
}
Кнструкция цикла while:
while (условие продолжения цикла) {
// Тело цикла
//...
}
Кнструкция цикла for:
for (n = 1; n <= 20; ++n) {
//...
}
Цикл for в языке C сводится к while:
n = 1;
while (n <= 20) {
//...
++n;
}
while и for -- циклы с ПРЕДУСЛОВИЕМ:
условие проверяется ПЕРЕД каждой итерацией цикла.
В C есть еще цикл с постусловием:
do {
// Тело цикла
. . .
} while (условие);
Рекомендую никогда такой цикл не использовать!
Если нужно выйти из середины цикла, то можно использовать
оператор break
while (true) {
. . .
if (условие завершения)
break;
. . .
}
Если надо пропустить итерацию, то можно использовать
оператор continue
i = 10;
for (int j = 1; j <= 20; ++j) {
if (j == i)
continue;
// тело цикла
}
тело цикла будет выполняться для j = 1, 2, ..., 9, 11, 12, ..., 20.
Пример 3:
решение квадратного уравнения.
Файл sqeq.cpp
#include <stdio.h>
#include <math.h>
int main() {
// Решаем квадратное уравнение
// a*x^2 + b*x + c = 0
double a, b, c;
printf("a, b, c:\n");
if (scanf("%lg%lg%lg", &a, &b, &c) < 3) {
printf("Error input.\n");
return (-1);
}
if (fabs(a) <= 0.) {
printf("Coeff. a must be nonzero.\n");
return (-1);
}
double d = b*b - 4.*a*c;
if (d < 0.) {
printf("Not roots.\n");
return 0;
}
d = sqrt(d);
double x1 = (-b + d)/(2.*a);
double x2 = (-b - d)/(2.*a);
printf("x1 = %g, x2 = %g\n", x1, x2);
return 0;
}
Сайт курса:
http://mech.math.msu.su/~vvb/
http://mech.math.msu.su/~vvb/PhysChem/
Дома:
для пользователей MS Windows
установить CygWin + Notepad++
(в CygWin установить g++, make, gdb).
Решить задачу: напечатать все пары простых чисел близнецов
(p, q): p, q простые, q == p + 2 и p, q <= n
в интервале от 2 до n.
Программа должна ввести n (диапазон чисел) и
напечатать все пары простых чисел-близнецов.
Результат пошлите мне на vladibor2016@yandex.ru
===================================================
15 Sep 2021
Оформление текста программы на C/C++.
1. Вложение структурных операторов показывается с помощью отступов
(identation):
while (условие) {
********
********
********
*** Тело цикла -- сдвинуто на 4 пробела вправо
*** Не должно быть табуляций!!!
********
********
********
********
}
*********
*********
*********
if (условие) {
******
******
******
******
******
for (.......) {
***************
***************
***************
***************
***************
***************
}
*************
*****************
}
if (....) {
*******;
}
bool isPrime(int m) {
******************
******************
******************
******************
******************
******************
******************
******************
}
======== Другой стиль (обычно используется фирмой Microsoft)
bool isPrime(int m)
{
*********
**********
while (******)
{
************
************
************
************
************
}
************
}
Другие языки
while (**********)
**************
**************
**************
**************
endwhile
if (**********)
**************
**************
**************
**************
endif
if (**********)
**************
**************
**************
**************
fi
-----------------------------------------------
Мы в прошлый раз написали программу, которая решает
квадратное уравнение по школьной формуле. Оказывается,
что в разных ситуациях часто возникают также кубические
уравнения.
Например, уравнение
x^3 - x - 1 = 0
(Напомним, что решением уравнения x^2 - x - 1 = 0
является "золотое сечение" (1 + sqrt(5))/2.)
Как егo решить?
Выведем формулу Кардано
(На самом деле, кубическое уравнение имеет вид
a*x^3 + b*x^2 + c*x + d = 0
Сначала разделим на a:
x^3 + b'*x^2 + c'*x + d' = 0, где b' = b/a
Сделаем замену x = y + s и подберем s так, чтобы коэффициент при
x^2 обнулился
(y + s)^3 + b*(y + s)^2 + b*(y + s) + d = 0
y^3 + 3*s*y^2 + . . . + b*y^2 + ...
y^2 + (3*s + b)*y^2 + ...
Положим s = (-b)/3. Уравнение примет вид
x^3 + a*x + b = 0
К такому виду сводится любое кубическое уравнение.)
Идея: сделаем замену
x = u + v
где u, v -- независимые независимые переменные.
x^3 - x - 1 = 0
(u+v)^3 - (u+v) - 1 = 0
u^3 + 3*u^2*v + 3*u*v^2 + v^3 - (u+v) - 1 = 0
u^3 + v^3 - 1 + 3*u*v*(u + v) - (u+v) = 0
u^3 + v^3 - 1 + (u+v)*(3*u*v - 1) = 0
------------- =================
Чтобы найти один корень, достаточно решить систему
{ u^3 + v^3 - 1 = 0
{ (u+v)*(3*u*v - 1) = 0
Сократим второе уравнение на (u+v), получим систему:
u^3 + v^3 - 1 = 0
3*u*v - 1 = 0
Из второго уравнения выразим v:
v = 1/(3*u)
Подставим в первое уравнение:
u^3 + [1/(3*u)]^3 - 1 = 0
u^3 + 1/(27*u^3) - 1 = 0
Обозначим y = u^3
y + 1/(27*y) - 1 = 0
y^2 - y + 1/27 = 0
Достаточно найти одно решение:
y = (1 + sqrt(1 - 4*(1/27)^2))/2
u = y^(1/3)
v = 1/(3*u)
x = u + v
Проблема: далеко не любое уравнение можно решить точно.
Есть формула Кардано для кубического уравнения, есть формула
для уравнения четвертой степени (сводится к формуле Кардано
и к биквадратному уравнению).
Вопрос: а для уравнений 5, 6, ... степеней, есть ли формулы?
Алгебра как наука началась с результата Абеля и Галуа о том,
что не существует формулы для решения общего уравнения 5-й степени.
x^5 + a*x^3 + b*x^2 + c*x + d = 0
Для такого "общего" уравнения не существует формулы (Абель).
Галуа сформулировал критерий существования формулы для решения уравнения
любой степени в радикалах: группа Галуа должна быть разрешима.
Как быть, если нужно решить уравнение общего вида?
Нам всегда нужно приближенное значение корня. Его для любой непрерывной
функции можно найти простым способом половинного деления отрезка.
ПЕРЕРЫВ до 14:25
Нахождение корня методом деления пополам:
файл root.cpp
#include <stdio.h>
#include <math.h>
#include <assert.h>
double f(double x);
double root(double (*f)(double x), double a, double b);
const double EPS = 1e-15;
int main() {
while (true) {
double a, b;
printf("a, b:\n");
if (scanf("%lg%lg", &a, &b) < 2)
break;
if (f(a)*f(b) > 0.) {
printf("No roots on [a, b]\n");
} else if (a >= b) {
printf("a must be greater than b\n");
} else {
double x = root(&f, a, b);
printf("Root of function: %.15g\n", x);
}
}
return 0;
}
double f(double x) {
// return sin(x);
return x*x*x - x - 1.;
}
double root(double (*f)(double x), double a, double b) {
assert(a <= b);
int sa, sb;
double fa = f(a);
if (fa > 0) {
sa = 1;
} else if (fa < 0) {
sa = (-1);
} else {
return a;
}
double fb = f(b);
if (fb > 0) {
sb = 1;
} else if (fb < 0) {
sb = (-1);
} else {
return b;
}
assert(sa*sb < 0);
while (fabs(b - a) > EPS) {
double c = (a + b)/2.;
if (c <= a || c >= b)
break;
double fc = f(c);
int sc;
if (fc > 0) {
sc = 1;
} else if (fc < 0) {
sc = (-1);
} else {
return c;
}
if (sa*sc < 0) {
// Select the left part of [a, b]: b <-- c
b = c;
fb = fc;
sb = sc;
} else {
// Select the right part of [a, b]: a <-- c
a = c;
fa = fc;
sa = sc;
}
} // end while
return (a + b)/2.;
}
Дома:
тем, кто не сделал первую программу -- сделайте ее!
те, кто ее сделал, могут либо не делать ничего,
либо сделать еще одну задачу:
найдите все совершенные числа в диапазоне от 2 до n.
Целое число называется совершенным, если оно равно сумме
своих собственных делителей.
6 = 1 + 2 + 3;
28 = 4*7 = 1 + 2 + 4 + 7 + 14;
Все должны сдать по первой теме хотя бы одну из этих двух задач
(либо обе).
=================================================
22 Sep 2021
На "чистом" языке C уже почти никто не пишет программы,
все используют язык (и компилятор) С++. При этом нормальные
C-программы целиком компилируются в С++. При этом C++
запрещает некоторые патологические преобразования, например,
указателей, более естественный контроль ошибок и т.п.
Мы пока не используем таких черт C++, как
1) классы,
2) полиморфизм,
3) переопределение операторов для классов,
4) template, ...
остаемся в основном в рамках С.
Единственная вещь из C++, которую мы используем -- это захват
и освобождение динамической памяти: new и delete[].
Текст программы очень важен!!!
Если программа работает правильно, но ее текст плохой,
то такая программа никуда не годится. Результат работы программиста --
это именно текст программы.
К сожалению, язык C (и за ним C++) не накладывает никаких ограничений
на текст программы.
Неформальные соглашения по оформлению текста C-программы.
1. Не должно быть табуляций в тексте.
2. Вложенность конструкций определяется по отступам в началах строк.
Общепринятый размер отступа -- 4 пробела.
3. После фигурных скобок нельзя ничего писать, кроме комментариев.
if (a == b) { c += b; ++d; }
Вместо этого пишем
if (a == b) {
c += b; ++d;
}
4. Существуют 3 стиля расстановки фигурных скобок:
Unix-стиль (Java, Python, ...)
while (услови) {
// Тело цикла
**************
**************
**************
for (....) {
*************
*************
*************
*************
}
**************
**************
}
*************
*************
*************
MS Windows-стиль
while (условие)
{
// Тело цикла
**************
**************
**************
for (....)
{
*************
*************
*************
*************
}
**************
**************
}
*************
*************
*************
Третий стиль (прошу, не используйте его!!!)
while (условие)
{
// Тело цикла
**************
**************
**************
for (....)
{
*************
*************
*************
*************
}
**************
**************
}
*************
*************
*************
===================================
while (условие)
if (.....) {
********
*********
} else
*****
Плохо! Лучше:
while (условие) {
if (.....) {
********
*********
} else {
*****
}
}
5. Избегайте длинных строк!!!
Во всех текстах, которые мы читаем (книги, газеты, ...),
длина строк должна быть не больше 60 символов.
В программировании допускаюся строки не длиннее 72 символов.
std::sout << "Пара простых чисел-близнецов: ("
<< p << ", " << q << ")" << std::endl;
6. Прототипы всех функций должны быть описаны в начале файла
(или, лучше, в h-файле). Например,
#include <stdio.h>
#include <assert.h>
bool isPrime(int m);
int main() {
int n;
printf("n:\n");
if (scanf("%d", &n) < 1 || n < 1) {
printf("Incorrect upper bound.\n");
return -1;
}
for (int m = 2; m <= n; ++m) {
if (isPrime(m)) {
printf("%d\n", m);
}
}
return 0;
}
bool isPrime(int m) {
if (m <= 1)
return false;
if (m == 2)
return true;
if (m%2 == 0)
return false;
assert(m >= 3 && m%2 != 0);
int d = 3;
while (d <= m/d) {
if (m%d == 0)
return false;
d += 2; // d = d + 2;
}
return true;
}
ПЕРЕРЫВ 10 минут до 14:40
Массивы и указатели в языке C/C++
=================================
Массив (таблица) -- набор ячеек памяти с последовательными
адресами. Размер массива в нормальном С (С90) и в С++ описывается
заранее и вычисляется статически (по тексту программы, а не на стадии
выполнения).
int f(int x, double y) {
int a[100];
...
}
Элементы массива доступны по индексам, которые начинаются с нуля
(индекс элемента == смещение относительно начала массива)
a[0], a[1], a[2], ..., a[99]
они идут подряд в памяти по возрастанию индексов.
Элементы массива можно читать, или в массив можно записывать
новые значения.
int x = a[12]; int y = 2*a[10] + 1;
a[5] = x + 3;
Массив очень важен, потому что память компьютера представляет
собой массив байтов, или ее можно рассматривать как массив слов.
Байт -- 8 двоичны разрядов.
Слово -- это либо 4 байта, адрес младшего байта делится на 4
(в 32-разрядной архитектуре);
либо 8 байтов, адрес младщего делится на 8 (в 64-разрядной).
Массив -- структура данных произвольного доступа (Random Access Memory).
Это означает, что любой элемент массива читается за одно действие
(вернее, за фиксированное число операций, не зависящее от индекса).
Отметим, что память компьютера используется по разному. Есть 3 вида памяти
(с точки зрения программы, физически любая память одинакова):
1) стековая память -- захватывается в стеке в момент входа в функцию;
2) статическая память -- выделяется до начала работы работы
функции main();
3) динамическая память (куча) -- резервуар памяти (очень большой
по объему), в которой можно захватить кусок нужного размера,
поработать с ним и затем освободить его, когда работа закончена.
Размер стека программы по умолчанию 4M
0| |
4| |
8| |
12| |
16| |
...|========|
| |
|--------|<--- SP
| |
| |
|--------|<--- SP (Stack Pointer)
| top |
| |
| |
2^32-4| |
+--------+
Мораль:
описывать большие массивы внутри функции нельзя!!!
Как быть?
1. Можно описать массив как статический -- определить его как глобальный
вне функции.
2. Статических переменных желательно по возможности избегать --
с ними возникают трудности при параллельном программировании
(например, в оконных системах).
3. Большие массивы захватыватся в динамической памяти (куче).
Указатели -- это переменные (или выражения), которые содержат адреса
объектов.
int n = 10, k = 100;
int *p;
p = &n;
// После этого *p -- это тот объект, на который указывает p
// Операция * называется разыменование (dereference)
*p = 77; // то же самое, что n = 77
p = &k;
n = *p + 50; // n == 150
Такого рода указатели, содержащие адреса объектов, были в разных языках
программирования, но ни в одном, кроме C, не додумались до масштабирования
указателей! Именно в этом главное достоинство языка С, которое обеспечило
его полную победу над всеми другими языками.
p содержит адрес целочисленной переменной, который делится на 4.
Пусть, допустим, p содержит 1000000.
++p;
Чему равно содержимое p?
p будет указывать на следующее число, т.е. содержать адрес 1000004
Какие операции с указателями возможны?
1. Увеличить ++, уменьшить --
2. К указателю можно прибавить или отнять целое число.
p = (int*) 1000000;
int n = (int)(p + 4); // n = 1000016
3. Разность двух указателей -- целое число
int n = p - q;
(содержимое p - содержимое q)/sizeof(n)
4. Сумма двух указателей БЕССМЫСЛЕННА!
5. Имя массива есть указатель на его начальный элемент
int a[100];
int *p = a; // p = &(a[0]);
// после этого выражения a[i] и *(p + i) означают одно и то же!
6. Указатель можно использовать как массив:
p[i] === *(p + i)
Захват сегмента размера n в динамической памяти:
int *p = new int[n];
После окончания работы с сегментом памяти надо его освободить:
delete[] p;
Пример: вычислить первые n простых чисел. В предыдущей программе
для проверки простоты числа p мы использовали пробные деления на все нечетные
числа до квадратного корня из p. В новой программе мы проверяем делимость
числа только на простые числа, найденные на предыдущих шагах алгоритма,
также не превышающие квадратного корня из p. Найденные простые числа хранятся
в массиве размера n, который захватывается в динамической памяти.
#include <stdio.h>
int main() {
while (true) {
int n; // Number of primes to find
printf("n:\n");
if (scanf("%d", &n) < 1 || n <= 0)
break;
int *primes = new int[n];
primes[0] = 2;
int numPrimes = 1;
int p = 3;
while (numPrimes < n) {
bool pIsPrime = true;
int k = 1;
while (pIsPrime && k < numPrimes && primes[k] <= p/primes[k]) {
if (p % primes[k] == 0) {
pIsPrime = false;
} else {
++k;
}
}
if (pIsPrime) { // if p is prime, then
primes[numPrimes] = p; // store p in the array
++numPrimes;
}
p += 2;
}
// Print the prime numbers found
for (int i = 0; i < n; ++i) {
if (i > 0) {
// Delimeter
if (i%5 == 0) { // Print 5 numbers in a line
printf("\n");
} else {
printf(" ");
}
}
printf("%d", primes[i]);
}
printf("\n");
delete[] primes;
}
return 0;
}
Дома:
1) теорема Лагранжа: любое целое положительное число можно представить
в виде суммы не более чем 4-х квадратов.
Написать программу, которая для данного числа находит его
представление в виде МИНИМАЛЬНОГО числа квадратов;
2) номер четной длины n называется счастливым,
если сумма первых его n/2 цифр совпадает с суммой последних
n/2 цифр.
Найти количество счастливых билетов длины n, где n четное, n <= 16.
Пример счастливого билета: 812731
====================================
29 Sep 2021
Домашние задания:
оформление текста
отступ равен уровню вложенности * 4
каждая открывающая фигурная скобка увеличивает уровень вложенности на 1,
закрывающая уменьшает на 1.
После фигурной скобки запрещено писать всё, кроме комментариев.
Строки не должны быть длиннее 72 символов.
Табуляций быть не должно! (Настройти редактор так, чтобы он заменял
табуляции на пробелы с отступом в 4 позиции.)
Не используйте российскую локализацию!
3.1415926
3,1415926
2, 3, 5.7, 37.65,
2; 3; 5,7;
csv-файлы (Comma Separated Values).
4,256,128,376
setlocale(LC_ALL, "C");
===================================
Есть неприятное расширение языка C99
(в стандарте С++ этого расширения нет!):
VLA -- Variable Length Arrays
int f() {
int n;
. . .
// double a[n]; // Так писать нельзя! VLAs are forbidden!
double* a = new double[n];
. . .
delete[] a;
return res;
}
Прошлый раз разобрали задачу: найти первые n простых чисел
(путем пробных делений на простые числа до квадратного корня,
найденные на предыдущих шагах; эти числа мы сохраняем в массиве).
Есть более быстрый способ нахождения простых чисел, не проевосходящих n:
решето Эратосфена.
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26
- - - = - = -- == -- == -- == -- == -- == \\ --
2 3 5 7 11 13 17 19 23
Числа 0, 1 вычеркнуты как не простые.
Двойка -- первое простое число.
Вычеркиваем все четные числа, которые делятся на 2, начиная с 2^2.
Потом находим первое невычеркнутое число -- это 3, второе простое число.
Вычеркиваем все числа, кратные 3, начиная с 3^2 = 9, с шагом 3.
Находим первое невычеркнутое число -- это 5, третье простое число.
Вычеркиваем все числа, кратные 5, начиная с 5^2 = 25, с шагом 5.
. . .
Заканчиваем, когда квадрат очередного (невычеркнутого)
простого числа больше, чем n.
#include <stdio.h>
int main() {
while (true) {
int n;
printf("n:\n");
if (scanf("%d", &n) < 1 || n <= 2) {
break;
}
bool* prime = new bool[n + 1];
// Инициализация
for (int i = 0; i <= n; ++i) {
prime[i] = true;
}
prime[0] = false;
prime[1] = false;
int p = 2;
int numPrimes = 1;
while (p*p <= n) { // Пока есть что вычеркивать
for (int i = p*p; i <= n; i += p) {
prime[i] = false; // Вычеркиваем числа, кратные p,
// начиная с p^2
}
// Ищем первое невычеркнутое число
++p;
// while (prime[p] == false) // ТАК НЕ НАДО!!!
while (p <= n && !prime[p]) {
++p;
}
if (p > n)
break;
++numPrimes;
}
printf("Number of primes <= n: %d\n", numPrimes);
printf("Prime numbers:\n");
int k = 0;
for (int i = 0; i <= n+1; ++i) {
if (prime[i]) {
// Печатаем разделитель
if (k > 0) {
if (k%5 == 0) {
printf("\n");
} else {
printf(" ");
}
}
printf("%d", i);
++k;
}
}
printf("\n");
}
return 0;
}
Задача про счастливые билеты
Примитивная программа, находящая число счастливых билетов длины 6,
состоит из 6 вложенных циклов.
// Compute the number of lucky tickets of length 6
#include <stdio.h>
bool lucky(const int *number);
int main() {
int number[6];
int n = 0; // Number of lucky tickets
for (int i0 = 0; i0 < 10; ++i0) {
number[0] = i0;
for (int i1 = 0; i1 < 10; ++i1) {
number[1] = i1;
for (int i2 = 0; i2 < 10; ++i2) {
number[2] = i2;
for (int i3 = 0; i3 < 10; ++i3) {
number[3] = i3;
for (int i4 = 0; i4 < 10; ++i4) {
number[4] = i4;
for (int i5 = 0; i5 < 10; ++i5) {
bool lucky(const int *number);number[5] = i5;
if (lucky(number)) {
printf(
"%d%d%d%d%d%d\n",
i0, i1, i2, i3, i4, i5
);
++n;
}
}
}
}
}
}
}
printf("Number of lucky tickets: %d\n", n);
printf("Probability of lucky ticket: %g\n", double(n)/1000000.);
return 0;
}
bool lucky(const int *number) {
int s0 = number[0] + number[1] + number[2];
int s1 = number[3] + number[4] + number[5];
return (s0 == s1);
}
Как реализовать переменное число вложенных циклов?
Пусть k -- число вложенных циклов.
Идея: имеем массив контрольных переменных длины k.
Пусть, например, k = 6.
Допустим, переменные имеют значения
1, 3, 2, 7, 9, 9
Как найти следующий номер?
Надо к этому набору как бы прибавить 1.
Идем справа налево, к последней цифре прибавляем 1.
Если цифра < 9, то на этом заканчиваем.
Если 9, то меняем ее на 0 и переносим 1 в следующий разряд
(предыдущая по порядку цифра):
.
1, 3, 2, 7, 9, 9
|
V
. перенос обозначаем точкой над цифрой
1, 3, 2, 7, 9, 0
|
V
.
1, 3, 2, 7, 0, 0
|
V
1, 3, 2, 8, 0, 0
bool nextDigits(int* number, int k) {
int pos = k - 1;
while (pos >= 0) {
if (number[pos] < 9) {
++number[pos];
return true;
}
assert(number[pos] == 9);
number[pos] = 0;
--pos;
}
return false;
}
int -- 32 разряда
long long -- 64 разряда
формат %lld
================================
6 Oct 2021
Задача о теореме Лагранжа:
любое натуральное число представимо в виде суммы
не более чем четырех квадратов.
Требовалось найти минимальное представление:
20 = 16 + 4 = (4^2 + 2^2)
18 = 9 + 9 = 3^2 + 3^2
18 = 16 + 1 + 1
19 = 9 + 9 + 1
19 = 16 + 1 + 1 + 1
int n = sqrt(m);
int n = 1.5; // n == 1
int k = -2.9; // k == -2
int n = int(sqrt(m)) sqrt(49) = 6.999999999999
int(6.999999999999) == 6
вместо
int n = sqrt(m);
надо писать либо
int n = int(round(sqrt(m)));
либо
int k = int(floor(sqrt(m))); // ceil(x) >= x
// floor(x) <= x
2451
sum = 6 C[s] -- число комбинаци их 2-х цифр
с суммой s
Число счастливых билетов =
C[0]^2 + C[1]^2 + ... + C[18]^2
00 0 0211
01 1
02 2
...
09 9
10 1
11 2
12 3
...
19
Почему надо писать
++i;
а не
i++;
add (R0)++, R1 Автоинкрементный режим адресации:
содержимое регистра увеличивается ПОСЛЕ
выполнения чтения элемента памяти
add --(R2), R3 Автодекрементный режим адресации
содержимое регистра уменьшается ПЕРЕД
выполнения чтения элемента памяти
int n = 5;
int k = n++ - 3; // Запоминается старое значение переменной n
// в некоторой вспомогательтной переменной,
// значение переменной n увеличивается,
// после чего в выражении используется
// старое значение n, сохраненное во
// вспомогательной переменной.
// n == 6, k == 2
int n = 5;
int k = (--n)*2;
// n == 4, k == 8
"Настоящие" программисты пишут фрагменты типа
int a[100];
...
int *p = a;
int *e = a + 100;
int sum = 0;
while (p != e) {
sum += *p++;
}
Я всегда выбираю самое тупое, но при этом понятное решение:
int a[100];
...
int *p = a;
int *e = a + 100;
int sum = 0;
while (p != e) {
sum += *p;
++p;
}
Или совсем тупо:
int a[100];
...
int sum = 0;
for (int i = 0; i < 100; ++i) {
sum += a[i];
}
Строку
setlocale(LC_ALL, "RUS");
не нужно использовать!!!
В российской локализации дробная часть вещественного числа
отделяется запятой, а не точкой, что противоречит стандартам,
принятым в программировании:
дробная часть отделяется точкой
3.1415926
а запятая используется для перечисления элементов списка:
1.5, 37, 28.96, 300, 200, 123.76
Лучше либо ничего не писать, либо использовать строку
setlocale(LC_ALL, "C");
"C" -- это минимальная локализация, которую понимают во всех странах.
==================================================================
13 Oct 2021
Требования к текстам программ:
1. Правильные отступы
int f(int x) {
***********
*************
while (...) {
***********
***********
***********
***********
if (...) {
**************
**************
**************
**************
}
*********
}
*****************
*****************
*****************
return res;
}
int main() {
**********
***********
cout << "Bye-bye!" << endl;
return 0;
}
=========================================
Другой стиль:
int f(int x)
{
***********
*************
while (...)
{
***********
***********
***********
***********
if (...)
{
**************
**************
**************
**************
}
*********
}
*****************
*****************
*****************
return res;
}
======================================
Утечка памяти:
программа захватывает память с помощью оператора new,
но не освобождает ее с помощью delete[]. Это приводит к утечку памяти,
что, в свою очередь, рано или поздно приведет к зависанию компьютера.
Есть утилита valgrind, которая отслеживает утечку памяти.
int n;
. . .
double *p = new double[n];
. . .
// Используем указатель как массив p[i]
// По завершению освобождаем память
delete[] p;
=============================================
cout << "Hello!";
==============================================
Зачем-то в теорема Лагранжа о четырех квадратах
использовали переменные типа double
double n1, n2, n3, n4, m;
. . .
if (n1*n1 + n2*n2 + n3*n3 + n4*n4 == m) {
. . .
}
Так писать нельзя!
Дело в том, что операции с вещественными производятся приближенно,
и надо с ними быть очень осторожным. Например,
(a + b) + c может не равняться a + (b + c)
Или, например, 1./3. + 1./4. != 7./12. в компьютере (разница небольшая,
но есть). Для сравнения чисел мы проверяем, что из разность по абсолютной
величине не превосходит eps
if (fabs(x - y) <= 1e-12) {
// Считаем, что числа x и y равны
. . .
}
======================================================
for (int i=0;i<100;i++)
ну хотя бы так:
for (int i = 0; i < 100; i++) {
. . .
}
Я ОЧЕНЬ советую никогда не писать i++, так все пишут,
но это ПЛОХО! Надо писать ++i
x*y + n++ - 96
Для запоминания старого значения n нужна дополнительная переменная.
А здесь она не нужна:
x*y + ++n - 96
В Питоне операции n++ вообще нет! Зато есть префиксная операция
n += 1 (что соответствует операция языка C
++n;
n += 1;
x*y + (n += 1) - 96
==================================================
Новая тема: вычисление функций на последовательнстях.
Дано некоторое множество X (например, числа) и последовательность w
элементов этого множества:
w = { x0, x1, x2, ..., x_n-1 }
Обозначим череp X* множество конечных последовательностей
элементов их X, включая пустую последовательность lambda <- X*
Рассматриваем некоторую функцию на множестве таких последовательностей:
F: X* --> Y
y = f(w) для w <- X*
Мы хотим предложить общую схему для вычисления подобных функций.
Эта схема основана на понятии индуктивной функции.
Функция F называется индуктивной, если при добавлении к последовательности
w еще одного одного элемента x в конец w новое значение функции F
можно вычислить, зная только старое значение F(w) и добавленный элемент x.
Более формально: функция
F: X* --> Y
индуктивна, если существует функция G
G: Y x X --> Y
такая, что для всякой последовательности w <- X* и всякого x <-X
F(w&x) = G(F(w), x)
-------------------------------------
w = { x0, x1, x2, ..., x_n-1 }
y = F(w)
Приписываем справа x к w
w' = { x0, x1, x2, ..., x_n-1, x }
Тогда
F(w') = G(y, x)
Большинство функций НЕ индуктивно, поэтому индуктивные функции --
это как алмазы в гравии (их очень мало, но они очень ценные).
Примеры индуктивных функций:
1) сумма элементов последовательности
S(w & x) = S(w) + x, в качестве G берется сумма двух элементов:
G(y, x) = y + x;
2) произведение элементов последовательности
P(w & x) = P(w)*x
3) максимум из элементов последовательности w
F(w & x) = max(F(w), x)
4) минимум, аналогично;
5) значение многочлена p(x) в фиксированной точке x = t,
коэффициенты многочлена заданы в последовательности
по убыванию степеней.
p(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
последовательность коэффициентов многочлена
w = { a0, a1, ..., a_n }
значение функции
F(w) = p(t)
Схема вычисления индуктивной функции: читаем последовательность
слева направо (от начала к концу), в любо момент храним вычисленное
значение функции на прочитанном отрезке последовательности. В конце
у нас получится значение функции на все последовательности.
Очень важно: мы начинаем со значения функции на ПУСТОЙ ПОСЛЕДОВАТЕЛЬЬНОСТИ
lambda.
y = F(lambda)
цикл пока в последовательности w есть непрочитанные элементы
| x = прочесть очередной элемент последовательности w
| y = G(y, x)
конец цикла
ответ = y
Для применения этой схемы надо
1) понять, чему равно значение функции на пустой последовательности.
Сумма(lambda) = 0 Сумма { x } = 0 + x
Произведение(lambda) = 1
0 и 1 -- это нейтральные элементы для операций + и *
Элемент e называется нейтральным для операции $,
если для любого x
e $ x = x
(другое название -- "левая единица")
Что является нейтральным элементов для операции max?
max(e, x) = x для любого x
Это должно быть число, которое заведомо меньше или равно
всех других чисел. Это элемент минус бесконечность, или -inf
(infinity)
Если брать целые числа, то существует минимальное целое число,
эта константа описана в файле limits.h
#iclude <limits.h>
она называется INT_MIN
d + d = 0 d = -d для d = -2147483648
INT_MAX = 2147483647 *.vob -- видео-файлы MPEG-2, размер 2GB
Для типа long long (64 разряда)
LLONG_MAX 9223372036854775807LL
Минимальные и максимальные вещественные константы описаны в файле
#include <float.h>
но он не всегда представлен в системе. Там константа
DBL_MIN -- это минимальное ПОЛОЖИТЕЛЬНОЕ число, а вовсе не
минус бесконечность.
Минус бесконечность тоже есть, но ее очень опасно использовать!
Но можно вместо минус бесконечности просто использовать очень
большое по абсолютной величине отрицательное число.
Например, я использую -1е30, это число
-1000000000000000000000000000000.
Напишем программу, которая читает последовательность чисел из файла
и вычисляет какую-нибудь индуктивную функцию, например, сумму
элементов последовательности. Последовательность записана
в файле "input.txt" в текущей директории.
Программа в файле sum.cpp
Разберем программу вычисления значения многочлена в точке x=t,
коэффициенты заданы в последовательности по убыванию степеней.
В школе это называется схемой Горнера.
p(x) = 2x^3 - 3x^2 + 5x - 4 =
= ((2*x - 3)*x + 5)*x - 4
p_n(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
Припишем к последовательности еще один коеффициент:
p_n+1(x) = a0*x^(n+1) + a1*x^n + ... + a_n-1*x^2 + a_n*x + a_n+1 =
= p_n(x)*x + a_n+1
p_n+1(x) = p_n(x)*x + a_n+1
Напишем программу polvalue.cpp, которая вычисляет значение
многочлена в заданной точке. Коэффициенты многочлена заданы по
убыванию степеней. В файле "input.txt" записана сначала точка t,
в которой мы вычисляем значение многочлена, и затем последовательность
коэффициентов многочлена по убыванию степеней.
=======================================================
20 Oct 2021
Индуктивные функции и их индуктивные расширения.
Общая схема вычисления значения функции на последовательности
Индуктивная функция
X -- некоторое множество
X* -- множество конечных последовательностей элементов из X
w = { x0, x1, ..., x_k }, k >= 0
lambda <- X* -- пустая последовательность
F: X* --> Y -- функция на множестве последовательностей
Функция F называется индуктивной, если при добавлении еще одного
элемента x в конец последовательности w новое значение функции
y' = F(w & x) можно вычислить, зная только старое значение функции
y = F(w) и добавленный элемент x. Это означает, что существует
отображение
G: Y x X --> Y (отображение прямого произведения Y на X в Y)
такое, что
F(w & x) = G(F(w), x)
Отображение G обычно называют действием алфавита X на множестве Y.
Примеры индуктивных функций:
1) сумма элементов числовой последовательности;
2) произведение элементов;
3) максимальный элемент;
4) минимальный элемент;
5) значение многочлена p(x) в точке x = t, коэффициенты
многочлена заданы в последовательности по убыванию степеней.
p(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
w = { a0, a1, a2, ..., a_n }
F(w) = p(x)| x = t
В школе схему вычисления этой функции называют схемой Горнера.
Схема вычисления индуктивной функции:
y = F(lambda)
цикл пока в последовательности есть непрочитанные элементы
| прочесть очередной элемент последовательности в x
| y = G(y, x)
конец цикла
ответ = y
В этой схеме для каждой конкретной функции F надо понять,
чему равно значением функции на пустой последовательности F(lambda)
и как вычислить новое значение функции F(w & x) при добавдении
еще одного элемента в конец w (т.е. чему равно отображение G).
1) сумма: F(lambda) = 0, G(y, x) = y + x
2) произведение: F(lambda) = 1, G(y, x) = y*x
3) максимум: F(lambda) = -inf, G(y, x) = max(y, x)
0 + x = x
1 * x = x
max(-inf, x) = x
F(lambda) -- это левый нейтральный элемент e
для бинарной операции o
e o x = x для любого x
4) минимум: F(lambda) = +inf, G(y, x) = min(y, x)
5) значение многочл p(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
ена:
F(lambda) = 0, G(y, a) = y*t + a
Но индуктивность функции -- это приятное исключение.
Большинство реальных функций на последовательностях не индуктивны.
Пример: среднее арифметическое значение элементов последовательноси.
Покажем, что эта ф-ция не индуктивна.
w1 = { 1, 2 } F(w1) = 1.5
w2 = { 1.5 } F(w2) = 1.5
Добавим к обеим последовательностям элемент x = 3. Если бы функция
F была бы индуктивна, то новые значения F на обеих последовательностях
совпадали бы. Но
w1 & 3 = { 1, 2, 3} F(w1 & x) = 2
w2 & 3 = { 1.5, 3 } F(w2 & x) = 2.25
F(w1 & x) != F(w2 & x)
Значит, функция не индуктивна.
Другие примеры неиндуктивных функций:
число максимальных элементов целочисленной последовательности;
число локальных максимумов вещественной последовательности;
значение многочлена в точке, коэффициенты заданы по возрастанию
степеней:
p(x) = a0 + a1*x + a2*x^2 + ... + a_n*x^n
Как быть? Как же вычислить значение неиндуктивной функции?
Правильный ответ: свести вычисление неиндуктивной функции к
вычислению некой иньуктивной функции, зная значение
которой, можно легко вычислить значение исходной
функции.
Такая функция называется индуктивным расширением исходной неиндуктивной
функции. Формальное определение:
Пусть дана функция F: X* --> Y (предположим, что она не индуктивная)
Функция Phi: X* --> Z называется индуктивным расширение
функции F, если:
1) Phi -- индуктивная функция,
2) зная значение Phi, можно легко вычислить значение F.
На математической языке, существует отображение
P: Z --> Y (оно обычно называется проекцией)
такое, что композиция P o Phi = F
(это означает, что для любой последовательности w
F(w) = P(Phi(w)).
Как найти индуктивное расширение? Это всегда немножко творческая задача,
но, неформально, надо просто понять, что мешает вычислить новое
значение функции F(w & x) при добавлении элемента x к последовательности w,
зная только старое значение y = F(w) и добавленный элемент x?
Какой информации не хватает, чтобы вычислить новое значение?
Просто надо добавить эту информацию к значению функции F.
Пример: среднее арифметическое элементов последовательности.
Для вычисления нового значения нам не хватает длины
последовательности.
Удобно хранить сумму элементов последовательности и число элементов
(длину) последовательности.
F(w) = среднее арифметическое значение элементов(w)
Phi(w) = (сумма элементов w, длина w)
P(s, n) = s/n
Алгоритм вычисления среднего арифметического:
s = 0.
n = 0
цикл пока в w есть непрочитанные элементы
| прочесть очередной элемент в x
| s += x
| ++n
конец цикла
ответ = s/n
Среднее геометрическое: корень n-й степени из произведения элементов
последовательности, где n -- длина последовательности.
Со средним геометрическим мы сталкиваемся в жизни:
логарифмическая шкала в регуляторах громкости;
хроматическая гамма в музыке, ...
Среднее гармоническое:
обратная величина к среднему арифметическому обратных величин.
w = { x0, x1, x2, ..., x_n-1 }
y = (1/x0 + 1/x1 + ... + 1/x_n-1)/n
F(w) = 1/y
Где мы в жизни встречаемся со средним гармоническим?
Едем на автомобиле, путь разбит на n равных отрезков.
На i-м отрезке мы едем со скоростью v_i.
Тогда средняя скорость на всём пути есть среднее гармоническое
от скоростей на отрезках.
Пример: путь 100 км разбит на 2 участка по 50 км.
на первом участке мы едем со скорость 50 км/ч,
на втором -- со скоростью 100 км/ч.
Какова средняя скорость?
средняя скорость = путь / время
путь = 100 км
время = 50 км / 50 км/ч + 50 км / 100 км/ч = 1 + 1/2
средняя скорость = 100 / 1.5 = 66.666 км/ч
А среднее арифметическое 75 км/ч
Зависимость расхода бензина от скорости квадратичная,
в основном мощность тратится на сопротивление воздуха.
А сила сопротивления воздуха пропорциональна квадрату скорости.
Работа = сила * расстояние.
p(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
Вопрос: я проехал со скоростью 100 км/ч за время 1 час израсходовал
10 литров бензина. Следующий час я ехал со скоростью 200 км/ч.
Сколько бензина я израсходовал за второй час? Ответ: 80 литров
(сопротивление воздуха возрастает в 4 раза, расстояние в 2 раза =>
работа в единицу времени, т.е. мощность, возрастает в 8 раз).
Это очень хорошо было инвестно в начале XX века,
когда на автомобилях устанавливали рекорды скорости:
для того, чтобы удвоить рекорд скорости, надо было увеличить
мощность двигателя в 8 раз (кубическая зависимость).
--------------------------------------------------------------
ПЕРЕРЫВ 10 минут до 14:20
Давайте напишем программу, вычисляющую 3 средних значения
для последовательности положительных вещественных чисел.
Последовательность записана в файле "input.txt"
(в текстовом виде).
Программа в файле "meanVal.cpp".
#include <stdio.h>
#include <math.h>
bool meanValues(
FILE *f,
double* meanAriphm,
double* meanGeom,
double* meanHarm
);
int main() {
FILE* f = fopen("input.txt", "rt");
if (f == NULL) {
perror("Cannot open an input file");
return (-1);
}
double meanAriphm, meanGeom, meanHarm;
bool res = meanValues(f, &meanAriphm, &meanGeom, &meanHarm);
fclose(f);
if (!res) {
printf("Error contents of input file.");
return (-1);
}
printf(
"Mean ariphmetic = %g, geometric = %g, harmonic = %g\n",
meanAriphm, meanGeom, meanHarm
);
return 0;
}
bool meanValues(
FILE *f,
double* meanAriphm,
double* meanGeom,
double* meanHarm
) {
double sum = 0.;
double product = 1.;
double sumInv = 0.;
double n = 0.;
while (true) {
double x;
if (fscanf(f, "%lg", &x) < 1)
break;
if (x <= 0.) // Incorrect sequence
return false;
sum += x;
product *= x;
sumInv += 1./x;
n += 1.;
}
if (n <= 0.) // Empty sequence
return false;
*meanAriphm = sum/n;
*meanGeom = pow(product, 1./n);
// *meanHarm = 1./(sumInv/n)
*meanHarm = n/(sumInv);
return true;
}
Очень важный (изящный!) пример:
Вычислить значение производной многочлена в точке x = t.
Коэффициенты многочлена заданы по убыванию степеней.
p(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
w = { a0, a1, ..., a_n }
F(w) = p'(x) | x = t
Как вычислить производную при добавлении еще одного коэффициента
к последовательности коэффициентов многочлена?
p_n(x) = a0*x^n + a1*x^(n-1) + ... + a_n-1*x + a_n
p_n+1(x) = a0*x^(n+1) + a1*x^n + ... + a_n-1*x^2 + a_n*x + a_n+1 =
= p_n(x)*x + a_n+1
Дифференцируем это равенство и получаем выражением для производной:
p_n+1(x) = p_n(x)*x + a_n+1
p'_n+1(x) = p'_n(x)*x + p_n(x)
#include <stdio.h>
#include <math.h>
bool polDeriv(FILE* f, double *deriv);
int main() {
FILE* f = fopen("input.txt", "rt");
if (f == NULL) {
perror("Cannot open an input file");
return (-1);
}
double d;
bool res = polDeriv(f, &d);
fclose(f);
if (!res) {
printf("Incorrect file format.\n");
return (-1);
}
printf("Derivative of polynomial = %g\n", d);
return 0;
}
bool polDeriv(FILE* f, double *deriv) {
double t;
if (fscanf(f, "%lg", &t) < 1) {
return false; // Incorrect format of input file
}
double p = 0.;
double dp = 0.;
while (true) {
double a;
if (fscanf(f, "%lg", &a) < 1)
break;
dp = dp*t + p;
p = p*t + a;
}
*deriv = dp;
return true;
}
Проверка
Многочлен x^3 - x - 1 вычисляем в точке x = 5.
Содержимое файла input.txt
5
1 0 -1 -1
/home/vvb/PhysChem>./polderiv
Derivative of polynomial = 74
Проверка:
Производная
3*x^2 - 1
3*5^2 - 1 = 3*25 - 1 = 74
Ответ правильный.
Замечание:
для данной неиндуктивной функции может быть много индуктивных
расширений. Однако на множестве этих функций можно ввести частичный
порядок, относительно которого существует наименьший элемент,
т.е. минимальное индуктивное расширение.
На практике всегда минимальное индуктивное расширение
оказывается самым удобным! Искусство программиста состоит в том,
чтобы найти именно минимальное индуктивное расширение
(бритва Окама: не хранить ненужной, излишней информации).
Пример:
Найти количество строго возрастающих участков
последовательности целых чисел.
(Считаем, что строго возрастающий участок последовательности
состоит не менее чем из двух элементов.)
Надо хранить:
1) последний элемент;
2) нужно ли хранить предпоследний элемент???
1 1 2 3 4 5 5 5 4 3 4 7 9
--------- -------
Первое (плохое!) решение: храним
1) число строго возрастающих участков
2) храним два последних элемента
Это не минимальное индуктивное раширение: лишняя информация!
w1 = -10 1
w2 = -5 1
В частности, будут трудности для определиния значения ф-ции
на пустой последовательности.
Хорошее решение (минимальное инд. расширение): храним
1) число строго возрастающих участков;
2) последний элемент;
3) логическую переменную: true, если последний
элемент есть конец строго возрастающего участка, и
false в противном случае.
Значение на пустой последовательности:
1) число строго возрастающих участков = 0
2) последний элемент = -inf
3) false, поскольку посл. элемент не является концом возр. участка
n = 0
lastElem = -inf
assending = false
цикл пока есть непрочитанные элементы
| прочесть очередной элемент в x
| если x > lastElem
| | если не assending
| | | n += 1
| | | assending = true
| | конец если
| иначе
| | assending = false
| конец если
| lastElem = x
конец цикла
ответ = n
------------------------------------------
Очень важный для практики пример:
Для последовательности вещественных чисел найти
среднее квадратическое отклонение от среднего арифметического
(дисперсию случайной величины)
Характеристики случайной величины X:
математическое ожидание E(X) (Expectation)
мера разброса, дисперсия (Dispersion, или Variance)
случайной величины: E(|X - E(X)|^2)
измеряется в метрах в квадрате.
Средне-квадратичное уклонение -- квадратный корень из дисперсии,
выражается в метрах (standard deviation, sd).
E(X) = (x0 + x1 + ... + x_n-1)/n
Var(X) = ( (x0 - E(X))^2 + (x1 - E(X))^2 + ... + (x_n-1 - E(X))^2 ) / n
(содержательная оценка)
В несмещенной оценке надо делить на n-1.
Надо построить индуктивное расширение этой функции:
w = { x0, x1, x2, ..., x_n-1 }
Надо вычислить дисперсию
d = ( Sum_i (x_i - M)^2 )/n, где M = (Sum_i x_i)/n
Ответ: достаточно хранить первый момент S1 случайной величины,
второй момент S2 и длина выборки n
S1 = x0 + x1 + ... + x_n-1 M = S1/n
S2 = x0^2 + x1^2 + ... + x_n-1^2
( Sum_i (x_i - M)^2 )/n = 1/n*( Sum x_i^2 - 2*M*Sum x_i + Sum M^2 ) =
= S2/n + 1/n*(-2*M^2*n + n*M^2) = S2/n - S1^2/n^2
D = S2/n - S1^2/n^2 (состоятельная оценка оценка)
Несмещенная оценка
D = (n/(n - 1)) * (S2/n - S1^2/n^2) =
= S2/(n-1) - S1^2/((n-1)*n)
sd = sqrt(D)
==================================
3 Nov 2021
Работа с массивами
Предполагаем, что данные уже прочитаны (из файла, из сети
и т.д.). Надо либо
1) вычислить некоторую функцию, зависящую от массива
(например, найти число различных элементов в массиве,
множество элементов, значение многочлена или его производной
в заданной точке и т.п.);
2) преобразовать каким-либо образом массив (например,
сортировка -- упорядочить элементы массива по возрастанию).
Каждый из студентов должен сделать 4 задачи на эту тему.
Рассмотрим задачи
-- циклический сдвиг элементов массива вправо (задача 4);
-- сортировка (задача 12);
Упорядочить элементы массива по возрастанию
a[0] <= a[1] <= a[2] <= ... <= a[n-1]
Это известная очень важная задача.
Несложно доказать, что любой алгоритм сортировки на самом плохом
для него входе выполнит не меньше чем log_2(n!) сравнений.
log_2(n!) ~ n*log_2(n)
Формула Стирлинга:
n! ~ n^n e^(-n) sqrt(2*pi*n) при n --> infinity
log_2(n!) ~ log_2( n^n e^(-n) sqrt(2*pi*n) ) =
n*log_2(n) - n log_2(e) + 1/2*(1 + log_2(pi) + log_2(n))
============================================
o(n*log_2(n))
~ n*log_2(n)
Есть оптимальные алгоритмы сортировки, работающие за время O(n*log(n))
log_2(1000000) = log_2(1000^2) = 2*log_2(1000) ~= 2*log_2(2^10) = 20
1000 ~= 2^10 = 1024
Пузырьковая сортировка работает за время O(n^2)
Для n = 1000000 она будет работать в 50000 раз дольше
(если сортировка слиянием работает за время 1 секунда, то пузырьковая
будет работать за время 50000 секунд = 13.9 часов).
В заданиях по массивам время работы не имеет значения.
Перерыв 10 минут до 14:10
Разберем несколько других задач, которые сводятся к сортировке.
Задача 5. Определить количество различных значений в массиве целых чисел.
Могут быть разные решения. Например, для каждого элемента проверяем,
встречался ли он раньше. Если нет, то увеличиваем счетчик числа элементов.
int numElements(const int* a, int n);
. . .
int numElements(const int* a, int n) {
int num = 0;
for (int i = 0; i < n; ++i) {
bool duplicate = false;
for (int j = 0; j < i; ++j) {
if (a[j] == a[i]) {
duplicate = true;
break;
}
}
if (!duplicate) {
++num;
}
}
return num;
}
Мне такая реализация ОЧЕНЬ НЕ НРАВИТСЯ, потому что она работает
за время O(n^2). Хорошее решение состоим сначала в применении
алгоритма сортировки и затем решение этой задачи для уже упорядоченного
массива. Хотя пузырьковая сортировка также работает за время O(n^2),
но ведь можно применить какой-нибудь приличный алгоритм сортировки,
и тогда мы получим время работы O(n*log(n)).
int numElems(
int *a,
int n
) {
bubbleSort(a, n);
if (n == 0)
return 0;
int num = 1;
for (int i = 1; i < n; ++i) {
if (a[i] > a[i - 1]) {
++num;
}
}
return num;
}
Задача 7.
Удалить из массива целых чисел все отрицательные значения,
а оставшиеся неотрицательные сдвинуть (уплотнить)
к началу массива с сохранением их исходного порядка.
"Освободившиеся" элементы в конце массива заполнить нулями.
Идея: выполняем "сортировку", считая, что любой неотрицательный элемент
строго меньше, чем любой отрицательный, а любые два неотрицательных
элемента равны, так же как и равны любые два отрицательных элемента.
То есть надо применить пузырьковую сортировку (или любую стабильную
сортировку), изменив лишь критерий сравнения элементов.
void task7(int *a, int n) {
bool wasInversions = true;
int pass = 0;
while (wasInversions) {
wasInversions = false;
for (int i = 0; i < n - 1 - pass; ++i) {
// if (a[i + 1] < a[i]) {
if (a[i] < 0 && a[i + 1] >= 0) {
wasInversions = true;
// Swap the elements a[i] and a[i + 1]
int tmp = a[i];
a[i] = a[i + 1];
a[i + 1] = tmp;
}
}
++pass;
}
// Annihilate all negative elements in the end of array
for (int i = n - 1; i >= 0; --i) {
if (a[i] < 0) {
a[i] = 0;
} else {
break;
}
}
}
К сортировке сводится задача 9:
Назовем массив целых чисел "плотным", если множество значений
элементов этого массива целиком заполняет некоторый отрезком [p,q]
в множестве целых чисел, т.е. для любого целого x найдется элемент a[i],
равный x. Массив, состоящий из одинаковых элементов, равных x,
считается плотным, и его отрезок значений есть [x,x].
Определить, является ли массив целых чисел "плотным".
Результат является корректным при наличии в массиве хотя бы одного числа.
Ответ: два числа p,q, если массив является "плотным",
и одно число -1, если нет.
Задачи 14, 15 также сводятся к сортировке аналогично задаче 7.
Нетривиальная задача 20, решение несложное, если правильно
применить рекурсию.
Очень интересна задача 17:
Для заданного натурального числа n заполнить массив целых чисел
размера n+1 биномиальными коэффициентами C_n^k, k=0, 1, ..., n.
Для вычисления биномиальных коэффициентов использовать
треугольник Паскаля.
1
1 1
1 2 1
1 3 3 1
1 4 6 4 1
1 5 10 10 5 1
1 6 15 20 15 6 1
1 7 21 35 35 21 7 1
1 8 28 56 70 56 28 8 1
. . .
=============================================
10 Nov 2021
Представление объектов в компьютере:
-- целые числа,
-- вещественные числа
-- тексты.
Целые числа
тип char 8 битов (1 байт)
short 16 битов
int == long 32 бита
long long 64 бита
модификаторы типа: signed или unsigned
числа со знаком, числа без знака
Тип char -- это не символ, а целое число!!!
небольшого размера
если signed (например, в процессорах Intel)
то -128..127
d = -128 d + d = 256 == 0 (mod 2^8)
d = -d
Поскольку все языки программирования придумывались
в США, а американцам наплевать на остальной мир,
поэтому им было достаточно набора ASCII-символов
American Standard Code of Information Interchange
Набор включает большие и маленькие латинские буквы,
знаки препинания, некоторые знаки операций --
всего 128 символов (использующие 7 младших битов).
Unicode -- взяли 2 байта (65536 символов).
В настоящее время в Unicode более миллиона символов!
(21 бит)
UTF-8 совместима с ASCII, не зависит от
порядка байтов в слове процессора, одинакова во
всех странах, во всех операционных системах,
во всех компьютерах.
Вещественные числа
float -- атавизм старых компьютеров. Современные компьютеры
не работают с типом float!!!
Они работают с типом double, преобразуя при необходимости
тип float к типу double.
Тип float занимает 4 байта,
точность примерно 7 значащих цифр.
double (двойная точность) -- 8 байтов (64 бита, в Intel
внутреннее представление 80 битов), примерно 15 значащих
цифр.
Linux: gcc компилятор, long double == 80 битов
MS Windows, MS Visual C++ long double == double
bool
логический тип, 1 байт
true -- любое ненулевое значение
false -- 0
bool b = (0 == 0); // Значение равно 1
bool c = (1 == 0); // Значение равно 0
if (b == true) { // Грубейшая ошибка!!!
// Надо писать:
if (b) {
...
}
Надо четко понимать:
НИКАКИХ ЦЕЛЫХ ЧИСЕЛ В КОМПЬЮТЕРЕ НЕТ!!!
При всех операциях с целыми числами от результата остается
только 32 младших бита для типа int (8 битов для char,
16 для short, 64 для long long). То есть от результата мы всегда
берем остаток при делении на 2^32 = 4294967296.
Кольцо вычетов по модулю m.
Зафиксируем целое число m > 0
x ~ y <==> m | (x - y)
def
Всё множесто целых чисел разбивается на классы эквивалентных
элементов. Это фактор-множество по данному отношению эквивалентности
называется кольцом вычетов по модулю m и обозначается Z_m.
Например, при m = 5 имеется 5 классов эквивалентности
{ ..., -15, -10, -5, 0, 5, 10, 15, 20, 25, ... }
{ ..., -14, -9, -4, 1, 6, 11, 16, 21, 26, ... }
{ ..., -13, -8, -3, 2, 7, 12, 17, 22, 27, ... }
{ ..., -12, -7, -2, 3, 8, 13, 18, 23, 28, ... }
{ ..., -11, -6, -1, 4, 9, 14, 19, 24, 29, ... }
Это элементы кольца Z_5 вычетов по модулю 5
Операции в фактор-множестве определяются через представителей:
[ x1 ] + [ x2 ] = [ x1 + x2 ]
Система остатков: из каждого класса можно выбрать по одному
представителю, это называется система остатков. Обычно в Z_m
рассматриваеют две системы остатков:
1) неотрицательная система
0, 1, 2, ..., m-1
Например, для Z_5 это система
{ 0, 1, 2, 3, 4 }
2 + 4 == 1 (mod 5)
2) "симметричная" система остатков (различаются случаи
четного и нечетного m)
Четное m:
{ -m/2, -m/2 + 1, ..., -1, 0, 1, 2, m/2 - 1 }
Например, для m = 8
{ -4, -3, -2, -1, 0, 1, 2, 3 }
Нечетное m, действительно симметрия:
{ -(m-1)/2, -(m-1)/2 + 1, ..., -1, 0, 1, 2, ..., (m-1)/2 }
Например, для m = 5
{ -2, -1, 0, 1, 2 }
"Целые" числа в компьютере -- это элементы кольца вычетов Z_m по модулю m,
где m = 2^8 = 256 для типа char,
m = 2^16 = 65536 для типа short,
m = 2^32 = 4294967296 для типа int,
m = 2^64 = 18446744073709551616 для типа long long.
Тип unsigned char -- это когда выбираем неотрицательную систему остатков:
{ 0, 1, 2, ..., 255 }
Тип signed char -- это когда выбираем симметричную систему остатков:
{ -128, -127, ..., -2, -1, 0, 1, 2, 3, ..., 127 }
Аналогично все остальные "целые" типы.
Что нужно понимать:
никакого "переполнения" при работе с "целыми" числами не бывает!
Просто операции выполняются в Z_m, а не в Z.
Для типа signed char есть элемент d != 0 такой, что d + d == 0 (mod 256)
(-128) + (-128) = 0
Иди 126 + 4 = (128) + 2 == (-128) + 2 = -126
126 + 4 = -126 (mod 256)
1000000 == -128
1000000
--------
1|0000000
ПЕРЕРЫВ 10 минут до 14:25
В современных компьютерах "числа" представляются в двоичном виде,
используется либо симмеричная система остатков (для "чисел" со знаком),
либо неотрицательная система остатков.
На самом деле, все операции выполняются независимо от системы остатков,
просто разная интерпретация результатов. (Число типа char, равное 255,
можно трактовать как -1.)
Как получить доичный код числа -x?
Это число, которое, будучи прибавленным к x, даст 0.
x = 01101101
Алгоритм:
1) инвертируем x
2) прибавляем 1 к младшему разряду с учетом переносов
x = 01101101
invert 10010010
+1 10010011
Проверка
x 01101101
+ -x 10010011
--------------------
1|00000000
Это иногда (неграмотно!) называют представлением отрицательных чисел
в дополнительном коде
1 00000001
-1 11111111
2 00000010
-2 11111110
Отметим, что все "отрицательные" числа имеют старший бит, равный 1
(разряды в программировании нумеруются СПРАВА НАЛЕВО, начиная с нуля).
Почему используется двоичная система счисления?
Понятно, что 0 или 1 удобно представлять уровнями напряжения,
или элементами, которые могут находиться в двух состояниях.
Казалось бы, что десятичная система экономичнее. Но
сколько элементов надо, чтобы представить все числа от 0 до 999
в десятичной системе счисления?
Надо 3*10 = 30 элементов.
В двоичной системе 2^10 = 1024 надо 10 разрядов, каждый разряд
представляется 2-мя элементами, итого 2*10 = 20 элементов.
Двоичная система счисления намного экономичнее десятичной!
Но, оказывается, есть и более экономичная система счисления:
троичная.
Рассуждение: пусть надо записать все числа в диапазоне 0..N-1
В системе счисления с основанием b надо log_b(N) разрядов,
всего элементов b*(число разрядов) = b*log_b(N)
Найдем минимум этой функции b
f(b) = b*log_b(N)
f(b) = b*ln(N)/ln(b)
Для нахождения минимума вычислим ноль производной.
f'(b) = ln(N) (b/ln(b))'
(b/ln(b))' = b' * (1/ln(b)) + b * (1/ln(b))' =
= 1/ln(b) - b * (-1/b)/(ln(b))^2 =
= (ln(b) - 1) / (ln(b))^2
Получаем f'(b) = 0 при ln(b) = 1, т.е. b = e ~ 2.718281828459
Число 3 ближе всего к e.
Поэтому самая лучшая система счисления -- троичная.
В СССР в 60-е годы был спроектирован и реально изготовлен компьютер
"Сетунь", который работал в троичной системе счисления.
Возьмем N = 2^11 = 2048 число двоичных разрядов 11, элементов
2*11 = 22
Тогда 3^7 = 2187 число троичных разрядов 7, элементов
3*7 = 21
В троичной системе можно использовать 3 цифры: 0, 1, 2
Но замечательно то, можно использовать сбалансированную систему цифр:
-1, 0, 1
Обычно вместо -1 пишут единицу, надчеркнутую сверху:
_
1
или это заменяют буквой T
Представление чисел в сбалансированной троичной системе:
0 0
1 1
2 1T
3 10
4 11
5 1TT = 1*3^2 + (-1)*3 + (-1)
6 1T0
7 1T1
8 10T
9 100
10 101
Сбалансированная система хороша тем, что для изменения знака числа
достаточно изменить знак каждой цифры:
0 0
1 1 -1 T
2 1T -2 T1
3 10 -3 T0
4 11 -4 TT
5 1TT -5 T11
6 1T0 -6 T10
7 1T1 -7 T1T
8 10T -8 T01
9 100 -9 T00
10 101 -10 T0T
EC ЭВМ -- угробила советсткую и российскую электронику.
========================================================
Двоичная система не очень удобна для человека,
а десятичная система плохо согласуется с двоичной.
Но
1) можно использовать восьмеричную систему счисления.
Если каждую цифру (0..7) восьмеричного представления числа
записать 3-мя двоичными битами, то мы получим двоичную запись.
В 60-е, 70-е годы в осноном использовалась восьмеричная система
для записи машинных команд, чисел и т.п.
В языке C/C++ восьмеричные константы начинаются с нуля:
0377 == 255
Но -- восьмеричные константы плохо сочетаюся с байтами
поскольку 8 не делится на 3.
Поэтому в современных языках чаще используется
16-ричная система счисления. Цифры:
0 0000
1 0001
2 0010
3 0011
4 0100
5 0101
6 0110
7 0111
8 1000
9 1001
A 1010
B 1011
C 1100
D 1101
E 1110
F 1111
Как записать число в системе счисления с основанием b?
Младшая цифра -- это остаток при делении на b
Отнимаем младшую цифру и делим на b
Предпоследняя цифра -- это остаток при делении на b
Отнимаем предпоследнюю цифру и делим на b
. . .
Записываем цифры справа налево, пока не получим 0.
Для примера, переведем 1234 в 16-ричную форму.
1234 % 16 = 2
2
(1234 - 2)/16 = 77
77 % 16 = 13 -- цифра D
(77 - 13)/16 = 4
Ответ: 1234 == 0x4D2 = 0b010011010010
Содержимое одного байта (8 двоичных разрядов) записывается как
две шестнадцатеричные цифры.
Помним, что -1 == 0xFF для типа char
-1 == 0xFFFF short
-1 == 0xFFFFFFFF int
-1 == 0xFFFFFFFFFFFFFFFF long long
По теме "Двоичное представление целых чисел"
надо сделать ровно одну задачу: ту, номер которой
совпадает с вашим по журналу по модулю 6. Например,
если номер студента 15, то он делает задачу 3 == 15 (mod 6).
Как сделать обратную задачу?
Пусть есть текстовое представление числа в b-ичной системе счисления.
. . .
Отложим на следующий раз (идея: использовать схему Горнера
для вычисления значения многочлена по убыванию).
Как получить значение 16-ричной цифры?
int hexNumber(char digit) {
if ('0' <= digit && digit <= '9') {
return (digit - '0');
} else if ('A' <= digit && digit <= 'F') {
return (digit - 'A');
} else if ('a' <= digit && digit <= 'f') {
return (digit - 'a');
}
assert(false);
return 0;
}
0x1a23f = (((1*16 + 10)*16 + 2)*16 + 3)*16 + 15
======================================================
17 Nov 2021
Представление целых чисел в компьютере
Современные компьютеры работают в двоичной системе счисления.
Число от 0 до 999 в 10-тичной системе счисления предствляются
тремя цифрами, каждая цифра принимает 10 значений =>
всего надо 30 элементов.
То же самое в 2-ичной системе счисления:
числа от 0 до 1023 в 2-исчной системе представляются 10
двоичными числами, каждая цифра принимает 2 значения =>
всего надо 20 элементов.
Самая экономичная система счисления -- 3-ичная.
(Теоретически самая экономичная система счисления с основанием e = 2.718.)
Сравним 2-ичную и 3-ичную системы:
для записис чисел в диапазоне 0..2047 = 2^11 - 1 надо 11 2-ичных разрядов =>
всего 11*2 = 22 элемента.
Для записи чисел в диапазоне 0..2186 = 3^7 - 1 надо 7 троичных разрядов
всего 7*3 = 21 элемент.
Но главное преимущество троичной системы счисления состоит в ее
симметричности. Бывает обычная троичная система с цифрами 0, 1, 2
и "сбалансированная" троичная система с цифрами -1, 0, 1.
Цифра -1 обычно обозначается единицей с надчеркиванием
_
1
или буквой T (эта буква похожа на 1 с надчеркиванием).
Преимущество в том, что нам не нужет знак минус для отрицательных чисел:
для изменения знака числа достаточно каждую цифру заменить на противоположную.
Арифметика получается достаточно интересной, вместо 8-разрядный двоичных
байтов обычно используются 6-разрядные троичные; в конструкции компьютера
используюся электронные элементы, принимающие 3 значения.
Задачи на представление целых чисел.
Зафиксируем систему счисления (2-ичную, 3-ичную сбалансированную,
8-ричную, 10-ичную, 16-ричную). Есть две задачи:
1) дана текстовая запись числа, получить само число;
2) обратная задача: дано число, получить его текстовую запись.
Во всех случаях эти задачи делаются одинаково.
Задача 1)
d0 d1 d2 ... d_n -- цифры числа, b -- основание иситемы счисления.
Тогда само число вычисляется по схеме Горнера как значение многочлена
с коэффициентами по убыванию степеней в точке x=b.
n = d0*b^n + d1*b^(n-1) + ... + d_n-1*b + d_n =
= (...((d0*b + d1)*b + d2)*b + ...)*b + d_n
Задача 2)
последняя цифра числа n вычисляется как остаток от деления
n на b. Важно отметить, что можно рассматривать неотрицательный
остаток в диапазоне 0..b-1, но можно рассмотривать и симметричный остаток,
когда b нечетное:
диапазон -(b-1)/2, -(b-1)/2+1, ..., -1, 0, 1, ..., (b-1)/2
В частности, для сбалансированной троичной системы счисления
остаток принимает значения -1, 0, 1
Как вычислить такой остаток?
Пусть n > 0.
Используем следующий фрагмент кода на C++:
int n;
. . .
int r = n%3; int q = n/3;
assert(0 <= r && r < 3);
assert(n == q*3 + r);
if (r > 1) {
r -= 3;
++q;
}
assert(-1 <= r && r <= 1);
assert(n == q*3 + r);
Вычислив последнюю цифру числа n как остаток или сбалансированный
остаток r при делении n на b, мы число n заменяем на частное от
деления q, где n = q*b + r.
Затем вычисляем предпоследнюю цифру, и так далее.
Пример: записать число 20 в сбалансированной троичной системе
счисления.
20 = 6*3 + 2 = (6+1)*3 + (2-3) = 7*3 + (-1)
Последняя цифра -1 записывается как T
Число 20 заменяем на 7 и находим предпоследнюю цифру.
7 = 2*3 + 1
Предпоследняя цифра 1.
Число 7 заменяем на 2 и находим третью с конца цифру
2 = 0*3 + 2 = 1*3 + (-1)
Третья с конца цифра -1 записывается как T.
Число 2 заменяем на 1 и находим четвертую с конца цифру
1 = 0*3 + 1
Четвертая с конца цифра 1.
Записываем цифры в обратном порядке и получаем запись числа 20
в сбалансированной троичной системе:
1T1T
Проверка:
1*3^3 + (-1)*3^2 + 1*3 + (-1) =
= 27 + (-9) + 3 + (-1) = 20
ПЕРЕРЫВ 10 минут до 14:35
Продолжим.
Переделаем программу с троичной системы на шестнадцатеричную.
Представление вещественных чисел
В истории программирования использовались две формы представления
вещественных чисел:
1) в форме с фиксированной десятичной точкой;
2) в плавающей форме.
Пусть, к примеру, у нас есть дешевый компьютер, но для расчетов нам нужны
вещестенные числа. Тогда, например, можно представляеть вещественные
числа с точностью до 0.001, используя целые числа и зафиксировав
положение десятичной точки. Например, число 10.25 представляется
как
10250 === 10.250
Число 0.03 представляется как
30 === 0.030
Как выполнять операции с такими числами?
Сложение, вычитание -- просто складываем целые числа
1.5 + 32.25
1500
+ 32250
-----
33750 === 33.750
Целое число m представляет вещественное число x = m/1000
Как умножить x*y?
x = m1/1000
y = m2/1000
x*y = m1*m2/(1000^2) = (m1*m2/1000)/1000
============
Целое число
Таким образом, число x*y представляется целым числом
m1*m2/1000
Пример: 1.5*10.333 x = 1.5 = 1500/1000 m1 = 1500
y = 10.333 = 10333/1000 m2 = 10333
x*y = (1500*10333/1000) / 1000 =
= 15499 / 1000
x*y = 15.499500000000001 представляется числом 15499 = 1500*10333//1000
Деление
x = m1/1000
y = m2/1000
x/y = m1/m2 = (m1*1000/m2) / 1000
Число x/y представляется целым числом m1*1000/m2
Проблема:
если мы примерно знаем диапазон используемых чисел, то
числа с фиксированной точкой вполне нас удовлетворяют.
В случае инжекторного двигателя автомобиля вполне достаточно точности
0.001 и не очень больших чисел, так что дешевый процессор,
не поддерживающий вещественную арифметику на аппаратном уровне,
впрлне нас устраивает. Но
допустим, мы заранее не знаем, какие задачи будет решать компьютер,
мы можем работать как с размерами галактик, так и с размерами атомов.
Точность задается не минимальным числом (как 0.001 в нашем примере),
а количеством значащих цифр (в современных компьютерах 15-16
значащих цифр). Задаются отдельно цифры числа и отдельно его
порядок. Число представляется в современных компьютерах
в так называемой плавающей форме (float):
x = s * m * 2^e
где
s -- знак числа, плюс или минус s <- {-1, 1}
m -- мантисса, 1 <= m < 2
(примерно 15 десятичных или 53 двоичные цифры)
e -- порядок (exponent), принимает значения от -1023 до 1024
2^1023 ~ 10^307
123000000000000000000000000000000000000000000000000000000000000.0
0.000000000000000000000000000000000000000000000000000000000000123
Рассмотрим примеры, чему равны мантисса и порядок.
100 = 100*2^0 = 50*2^1 = 25*2^2 = 12.5*2^3 = 6.25*2^4 =
= 3.125 * 2^5 = 1.5625 * 2^6
Для числа 100 знак равен +1, мантисса равна 1.5625, порядок равен 6
Как это записывается в двоичном коде?
В C/C++ есть числа типа float (атавизм) и числа типа double.
Число float занимает 4 байта = 32 двоичных разряда.
Как они распределяются?
1 бит -- знак числа
8 битов -- порядок
23 разряда -- мантисса
Для типа double, всего 8 байтов или 64 двоичные разряда
1 бит -- знак числа
11 битов -- порядок
52 разряда -- мантисса
Для отрицательных чисел знак равен 1, для неотрицательных 0
Порядок хранится в так называемом смещенном виде:
для float хранится смещенный порядок e'
e' = e + 127
для double
e' = e + 1023
То есть для чисел double максимальный по абсолютной величине
отрицательный порядок равен -1023 (при e' = 0)
11 порядка означает, что максимальное значение
смещенного порядка
e' = 2^11 - 1 = 2047
=> максимальное значение порядка e = e' - 1023 = 1024
максимальное число имеет значение порядка 2^1024 ~ 10^308
Грубо округляя, диапазон чисел от 10^(-300) до 10^300.
Сколько значащих цифр в мантиссе?
Для типа double хранятся 52 разряда мантиссы.
1 <= m < 2
Двоичная запись m имеет вид
1.1001010100111000010100...000101
Целая часть мантиссы -- это всегда 1, поэтому она вообще
не хранится в памяти. Хранятся только двоичные знаки дробной
части, их 52 знака, или 15.65 десятичных знаков (между 15 и 16
десятичными цифрами). Грубо, точность -- примерно 15 десятичных
цифр.
1 + eps = 1
Максимальное такое число eps называется машинным эпсилоном,
или машинным нулем.
x = 1.1000101 * 2^2
y = 1.0010000 * 2^(-1)
При сложении выравниваем порядки:
y = 1.0010000 * 2^(-1) = 0.1001000 * 2^0 =
= 0.0100100 * 2^1 = 0.0010010 * 2^2
x+y 1.1000101 * 2^2
0.0010010 * 2^2
----------------------
1.1010111 * 2^2
У машинного eps мантисса при выравнивании уйдет за пределы разрядной сетки
=> это примерно 2^(-53)
Для вещественных чисел не выполняется заколн ассоциативности сложения:
(a + b) + c = a + (b + c)
Для чисел double
eps = 2.**(-53)
(1 + eps) + eps == 1
1 + (eps + eps) > 1
>>> eps = 2**(-53)
>>> eps
1.1102230246251565e-16
>>> 1 + eps
1.0
>>> 1. + eps == 1.
True
>>> 1. + 2.*eps == 1.
False
>>> (1 + eps) + eps != 1 + (eps + eps)
True
>>> (1 + eps) + eps == 1 + (eps + eps)
False
>>> (1 + eps) + eps
1.0
>>> 1 + (eps + eps)
1.0000000000000002
===============================================
24 Nov 2021
Конструкция инварианта цикла
Схема итерации
Дано некоторое множество X и в нем подмножество Q
Q есть подмножество X
Надо найти точку, принадлежащую этому подмножеству.
о многих случаях это трудно сделать сразу, поэтому
применяется схема итерации. Берем некоторую начальную точку x <- X
и строим некоторое отображение
T: X\Q --> X
такое, что оно в некотором смысле должно нас приближать к Q
x = x0 // Начальное приближени
цикл пока не Q(x)
| x = T(x)
конец цикла
утв: Q(x)
ответ = x
Пример: найти хорошее приближение числа pi в виде дроби m/n,
где m, n -- целые числа, причем ошибка приближения не превышает квадрата
знаменателя дроби
| pi - m/n | <= 1/n^2 <= eps
Возьмем eps = 1/1000000 и найдем хорошее приближение числа pi.
Используем схему итерации. Отображение T будет перебирать
все целочисленные точки первого квадранта:
(0, 0)
(1, 0), (0, 1),
(2, 0, (1, 1), (0, 2),
(3, 0), (2, 1), (1, 2), (0, 3),
. . .
#include <stdio.h>
#include <math.h>
void t(int x[2]);
bool goodAppriximation(int x[2]);
const double EPS = 1e-6;
int main() {
int x[2];
x[0] = 0; x[1] = 0; // (0, 0) is an initial point
while (!(
goodAppriximation(x) &&
fabs(double(x[0])/double(x[1]) - M_PI) <= EPS
)) {
t(x);
}
printf("Approximation of pi is %d/%d\n", x[0], x[1]);
printf(
"Precision is %g\n",
fabs(double(x[0])/double(x[1]) - M_PI)
);
return 0;
}
void t(int x[2]) {
if (x[0] == 0) {
x[0] = x[1] + 1;
x[1] = 0;
} else {
--(x[0]);
++(x[1]);
}
}
bool goodApproximation(int x[2]) {
double a = double(x[0])/double(x[1]);
double error = fabs(M_PI - a);
return (error <= 1./double(x[1]));
}
Схема инварианта цикла
Есть конфигурационное пространство задачи X
В множестве X есть два подмножества I < X, Q < X
I называется инвариантом цикла,
Q называется условием завершения цикла
Задача состоит в том, чтобы найти точку x,
принадлежащую пересечению подмножеств I и Q
Для этого сначала выбираем любую точку, принадлежащую I
x0 <- I, или на языке предикатов выполняется I(x)
Строится отображение
T: I\Q --> I
для которого I является инвариантом: T(I) есть подмножество в I
Схема решения задачи:
x = x0
утв: I(x)
цикл пока не Q(x)
| инвариат: I(x)
| x = T(x)
конец цикла
утв: I(x) и Q(x)
Сама по себе эта схема особой ценности не имеет, нужно уметь
ее применять в конкретных задачах.
Пример 1. Вычислить наибольший общий делитель двух целых чисел
(алгоритм Евклида).
Дано: целые числа m, n <- Z
Надо: вычислить gcd(m, n) Greatest Common Divizor НОД
X = { (a, b), a, b <- Z }
Инвариант: gcd(a, b) = gcd(m, n)
Условие завершения: b = 0, тогда gcd(a, 0) = a
Отображение, сохраняющее инвариант и приближающее нас к
выполнению условия завершения:
T: (a, b) --> (b, r), где r -- остаток от деления a на b.
Почему T сохранаяет инвариант:
Теорема: множества общих делителей пар (a, b) и (b, r) совпадают.
Доказательство
a = q*b + r, |r| < |b|
Пусть d | a и d | b
r = a - q*b, следовательно, r делится на d
Обратно. Пусть d | b, d | r => d | q*b + r
Теорема доказана.
Напишем функцию, вычисляющую gcd(m, n)
int gcd(int m, int n) {
int a = m; int b = n;
// assert: gcd(a, b) == gcd(m, n)
while (b != 0) {
int r = a%b;
a = b; b = r;
}
assert(b == 0);
return a;
}
Еще один замечательный алгоритм в теории чисел:
быстрое возведение в степень
Дано: число a (любое, можно даже матрица, многочлен, слово и т.п.)
целое неотрицательное число n
Надо: вычислить a^n
Применим схему инварианта цикла
X = { (b, k, p) }
Инвариант цикла:
b^k * p = a^n
Начальное значение:
b = a, k = n, p = 1
Условие завершения: k == 0, тогда из инварианта следует, что p = a^n
Отображение, сохраняющее инвариант и уменьшающее k
T( (b, k, p) ) = { (b*b, k/2, p), если k четное
{ (b, k-1, p*b), если k нечетное
int fastPower(int a, int n) {
assert(n >= 0);
int b = a;
int k = n;
int p = 1;
while (k > 0) {
if (k%2 == 0) {
k /= 2;
b *= b;
} else {
k -= 1;
p *= b;
}
}
return p;
}
ПЕРЕРЫВ до 14:30
Применение схемы инварианта цикла для приближенного вычисления
логарифма
Дано: вещественные числа a > 1, x, eps -- требуемая точность
Надо: вычислить логарифм числа x по основанию a
y = log_a(x)
с точностью не хуже eps
Разложение в ряд
ln(1+x) = x - x^2/2 + x^3/3 - x^4/x + ...
Ряд сходится при |x| < 1
Но с помощью схемы инварианта цикла мы напишем программу,
которую может написать любой школьник, не знающий разложения
логарифма в ряд, которая вычисляет значение логарифма для любого
x > 0.
Нам надо придумать инвариант цикла (творческая задача, нет общего
решения на все случаи жизни).
y = log_a(x) это означает по определению логарифма, что
a^y = x
Построим инвариант:
a^y z^t = x = const
Начальные значения:
y = 0, z = x, t = 1
Условие завершения:
для приближенного выполнения a^y = x
нужно, чтобы z^t ~= 1
Достаточно, чтобы z было бы не очень большим и не
очень маленьким и |t| <= eps:
1/a <= z <= a и |t| <= eps
Отображение в цикле, сохранающее инвариант:
{ ( y+t, z/a, t ), если z > a
T( (y, z, t) ) = { ( y-t, z*a, t ), если z < 1/a
{ ( y, z^2, t/2), если 1/a <= z <= a
double loga(double x, double a, double eps) {
double y = 0.;
double z = x;
double t = 1.;
while (z < 1./a || z > a || fabs(t) > eps) {
if (z > a) {
z /= a;
y += t;
} else if (z < 1./a) {
z *= a;
y -= t;
} else {
z *= z;
t /= 2.;
}
}
return y;
}
Вернемся к предыдущей теме: представление вещественных чисел
в компьютере.
Число типа double представляется в виде
x = s * 2^e * m
s = +-1 знак числа (1 бит в доичном представлении)
-1023 <= e <= 1024 порядок числа (11 битов в двоичном представлении)
1 <= m < 2 мантисса числа (52 бита в двоичном представлении)
всего 64 бита.
Особенности:
1) в двоичном коде хранится смещенный порядок
e' = e + 1023
целое число e' неотрицательно;
2) мантисса в двоичной записи имеет вид
1.101double loga(double x, double a, double eps) {
z double loga(double x, double a, double eps) {
z double loga(double x, double a, double eps) {
double y = 0.;
double z = x;
double t = 1.;
while (z < 1./a || z > a || fabs(t) > eps) {
if (z > a) {
z /= a;
y += t;
} else if (z < 1./a) {
z *= a;
y -= t;
} else {
z *= z;
t /= 2.;
}
} double loga(double x, double a, double eps) {
double y = 0.;
double z = x;
double t = 1.;
while (z < 1./a || z > a || fabs(t) > eps) {
if (z > a) {
z /= a;
y += t;
} else if (z < 1./a) {
z *= a;
y -= t;
} else {
z *= z;
t /= 2.;
}double loga(double x, double a, double eps) {
double y = 0.;
double z = x;
double t = 1.;
while (z < 1./a || z > a || fabs(t) > eps) {
if (z > a) {
z /= a;
y += t;
} else if (z < 1./a) {
z *= a;
y -= t;
} else {
z *= z;
t /= 2.;
}
}
return y;
}
return y;
return y;
*= z;
t /= 2.;
}
}
return y;
*= a;
y -= t;
} else {
z *= z;
t /= 2.;
}
}
return y;
1000101101110100...010
целая часть мантиссы всегда равна 1, поэтому этот бит не
хранится в памяти! Хранятся только биты дробной части мантиссы.
Как представляется число 1?
1 = + 2^0 * 1.000000.....000000000
e' = 1023
Напишем программу, которая получает биты двоичного представления
числа типа double
Прошлый раз мы указали, что есть текие числа y > 0., что
1. + y == 1.
При сложении происходит выравнивание порядков:
x1 + x2
Если порядок x1 больше, чем порядок x2, то
мы увеличиваем порядок x2 на 1, а мантиссу x2 дели пополам
т.е. сдвигаем мантиссу на 1 разряд вправо, повторяем это
до тех пор, пока порядки не выровняются.
1. = + 2^0 * 1.
eps = 2^-53 = + 2^(-53) * 1.
При выравнивании порядков мантисса второго числа 53 сдвигается
вправо, в результате она становится равной нулю. В результате
получаем
1. + eps == 1.
При этом 1. + eps*2. > 1.
Такое число называется машинным эпсилоном (или машинным нулем),
оно устанавливает предельную точность вычислений при работе
с double.
Получаем следующие странности:
1) не выполняется закон ассоциативности сложения:
(1. + eps) + eps == 1.
1. + (eps + eps) > 1.
(1. + eps) + eps != 1. + (eps + eps)
2) задача: указать массив из 10 элементов такой,
что, когда мы суммируем его слева направо,
а потом справа налево, то результат отличается
не меньше чем на 1000.
[1., eps, eps, ..., eps]
Суммируем слева направо, получаем s0 = 1.
Суммируем справа налево, получаем s1 = 1. + 9*eps > 1.
err = |s1 - s0|
Если мы все числа в массиве умножим на r, то ошибка
умножится на r.
Cуммировать можно только числа примерно одного порядка.
При неправильной вычислительной схеме вместо ответа получается
совершеннейший мусор.
===================================
1 Dec 2021
Замечания по задачам
1. Внимательнее читайте условие и требования к задачам!
1) Запрещено призводить вычисления в функции main. Цель --
научить оформлять алгоритмы в виде функций. Функция не
имеет проава что-либо печатать, ее задача -- вычислить
результат и вернуть его вызываемую программу.
Что с этим результатом будут дальше делать -- не ее забота.
2) Тексты программ, написанные без отступов или с грубыми
нарушениями отступов, не принимаются.
3) Последовательности, массивы -- целых или вещественных чисел,
внимательнее читайте условия! Например, оперции
x == y или x != y
для переменных
double x, y;
ЗАПРЕЩЕНЫ!!!
Вместо сравнения
x == y
надо сравнивать абсолютную величину разности чисел с eps:
fabs(x - y) <= 1e-12
4) Мы пока используем библиотеку стандартных функций языка C
printf, scanf, ...
Лучше не использовать (пока!) потоковый ввод-вывод из STL
языка C++
std::cout << "Hello!" << std::endl;
std::ostream& operator<<(std::ostream& s, const char* str) {
. . .
return s;
}
Для работы с файлами мы используем функции ANSI-библиотеки
FILE *f = fopen("input.txt", "r");
FILE *g = fopen("output.txt", "w");
fscanf(f, "%lg", &x);
fprintf(g, "x = %g\n", x);
5) Вместо функций языка C для работы с динамической памятью
double *a = (double *) malloc(100*sizeof(double));
. . .
free(a);
лучше использовать встроенные операторы new и delete языка C++:
double *a = new double[100];
. . .
delete[] a;
=======================================================================
Тема лекции: вычисление стандартных функций путем суммирования
степенных рядов и другими методами.
Схема инварианта цикла
Инвариант -- это некоторое статическое соотношение между
изменяющимися величинами. Понятия инварианта -- это ключевое
понятие в таких науках, как математика и физика.
Например, в физике главный инвариант -- это энергия системы.
Закон сохранения энергии формулируется как инвариант от
переменных, описывающих систему (например, координаты материальных
точек и их импульсы), аналогично закон сохранения импульса,
момента импульса и т.д.
В применении к вычислению логарифма.
Дано: x > 0, a > 1
Надо: вычислить приблизительно y = log_a(x)
т.е. найти y:
a^y ~= x
Добаим переменные z, t и сформулируем инвариант:
a^y z^t = x
Начальные условия подбираются так, чтобы инвариант выполнялся:
y = 0, z = x, t = 1
Условие завершения цикла:
z^t ~= 1, тогда ответ будет в переменной y
Для этого достаточно, чтобы
1/a <= z <= a и |t| <= eps
Действие (итерация), которое в конце-концов приведет нас
к ответу (уменьшающее t):
если z < 1/a
| z = z*a; y = y - t
иначе если z > a
| z = z/a; y = y + y
иначе
| t = t/2; z = z*z
конец если
========================================================
Еще один очень важный алгоритм: вычисление квадратного корня
из положительного вещественного числа (или вычисление корня
n-й степени), этот алгоритм был известен тысячи лет назад в
древнем Вавилоне. Сформулируем его в более современной
формулировке: метод итераций Ньютона.
Дано: число a > 0.
Надо: вычислить квадратный корень из a с точностью eps.
Ищет корень как решение уравнения
x^2 - a = 0
f(x) = x^2 - a
f'(x) = 2*x
Наришем функцию
double squareRoot(a, eps = 1e-12);
Берем начальное приближение
x0 = 1 + a/2
Легко проверить, что x0 > sqrt(a)
x0^2 = (1 + a/2)^2 = 1 + a + a^2/4 > a
Проводим касательную к графику функции f(x) в точке x0
и находим точку x1 пересечения кастельной с осью абсцисс
Точка графика (x0, y0) где y0 = x0^2 - a
dx = x0 - x1 = y0/f'(x0) = (x0^2 - a)/(2*x0) =
= (x0 - a/x0)/2
x1 = x0 - dx = x0 - (x0 - a/x0)/2 = (2*x0 - x0 + a/x0)/2 =
= (x0 + a/x0)/2
x_i ---> x_i+1 = (x_i + a/x_i)/2
ПЕРЕРЫВ 10 минут до 14:20
Задача: найти такие же приближения для корня n-й степени из a
f(x) = x^n - a
f'(x) = n*x^(n-1)
x0 = 1 + a/n
dx = x0 - x1 = (x0^n - a)/(n*x0^(n-1)) = (x0 - a/x0^(n-1))/n
x1 = x0 - dx = x0 - (x0 - a/x0^(n-1))/n =
= (n*x0 - x0 + a/x0^(n-1))/n =
= ((n-1)*x0 + a/x0^(n-1))/n
Суммирование рядов и вычисление стандартных функций
Пример 1. Вычисление sin(x)
sin(x) = x - x^3/3! + x^5/5! - x^7/7! + ...
Грубейшая ошибка -- вычислять факториал при суммировании ряда!
Надо лишь по предыдущему члену ряда (уже вычисленному на предыдущем шаге)
уметь вычислять следующий член ряда.
Напишем программу, вычисляющую sin(x) с помощью суммирования ряда
храним n и член ряда a
Формула для вычисления следующего члена ряда в случае sin(x):
a_next = -a_prev*x*x/((n+1)*(n+2)
/home/vvb/PhysChem>./mySin
x: 1
my sin(x) = 0.841470984809
standard sin(x) = 0.841470984808
x: 4
my sin(x) = -0.756802495308
standard sin(x) = -0.756802495308
x: 20
my sin(x) = 0.912945252923
standard sin(x) = 0.912945250728
x: 40
my sin(x) = -1.815509922049
standard sin(x) = 0.745113160479
x: 50
my sin(x) = -39110.458480393136
standard sin(x) = -0.262374853704
Для x = 50 суммирование ряда дает совершенно неверный ответ,
хотя математика утврждает, что ряд для sin(x) сходится всюду.
ПОЧЕМУ???
Ответ: из-за феномена машинного eps
1. + eps = 1. (при сложении происходит выравнивание порядков
и число eps исчезает)
Анагично, когда к большому числу прибавляем число, намного меньшее,
то большое число не изменяется (когда к горе прибавляем песчинку,
то гора от этого не изменяется).
Вывод: степенные ряды можно суммировать только для малых х,
желательно даже для |x| < 1 (тогда член ряда быстро уменьшается).
Для функции sin(x) можно воспользоваться периодичностью:
sin(x + k*2*pi) = sin(x)
Таким образом, можно аргумент привести к интервалу
-pi <= x <= pi
(Даже, используя формулы приведения
sin(pi/2 + x) = -cos(x)
sin(pi - x) = sin(x)
привести x к интервалу -pi/4 <= x <= pi/4, заменяя при необходимости
ряд для sin рядом для cos
cos(x) = 1 - x^2/2! + x^4/4! - x^6/6! + ...
).
Задача (очень важная!)
Вычислить y = arctg(x)
Замечание 1: в России мы пишем arctg(x), а мира общепринято обозначение
atan(x)
Замечание 2: никогда не используйте на практике фунции
arcsin(x), arccos(x)
(стандартные обозначения asin(x), acos(x))
Эти функции ОЧЕНЬ плохо себя ведут для значений
x, близких к 1 или к -1 (поскольку производная в
этих точках стремится к бесконечности).
Конкретный пример:
найти расстояние между 2-мя точками на поверхности Земли,
точки заданы своими координатами (широта, долгота)
Расстояние на сфере -- длина дуги большого круга.
По координатам точек несложно вычислить радиус-векторы v1 и v2,
проведенные из центра Земли к эти точкам.
Тогда расстояние равно углу в радианах, умноженному
на радиус Земли
d = угол(v1, v2)*R
R = 6400000 метров
Как вычислить угол?
НЕВЕРНОЕ решение:
угол = arccos(скалярное произведение(v1, v2)/R^2)
Скалярное произведение близко к 1 (угол маленький),
функция arccos даст очень большую погрешность.
ПРАВИЛЬНОЕ решение: забыть о существовании функций arcsin, arccos,
всегда пользоваться функцией arctg
(более того, программисты в основном используют
функций alpha = atan(y, x)
Как вычислить arctg(x)?
arctg(x) = x - x^3/3 + x^5/5 - x^7/y + ...
Cходится только для значений
-1 < x <= 1
Для x = 1 имеем arct(1) = pi/4
pi/4 = x - 1/3 + 1/5 - 1/7 + 1/9 - ...
Казалось бы,так можно вычислить pi. НЕЛЬЗЯ! Это ряд сходится
фантастически медленно.
Как посчитать arctg(x) для |x| > 1?
Идея та же, что и в случае sin(x):
используя какие-то свойства функции arctg(x), заменить значение
x меньшим значением.
Дано: x больше по модулю число, пусть x > 0
Надо: вычислить y = arctg(x)
tg(y) = x
Найдем z: arctg(z) = y/2
Тогда |z| < 1 и ряд сходится.
z = tg(y/2)
Применим формулу двойного угла
tg(2a) = sin(2a)/cos(2a) = 2sin(a)cos(a)/(cos(a)^2 - sin(a)^2) =
= 2 tg(a)/(1 - tg(a)^2)
tg(y) = 2 tg(y/2) / (1 - tg(y/2)^2)
x = 2*z/(1 - z^2)
Получили уравнение для z, выразим z через x
x*(1 - z^2) = 2*z
x*z^2 + 2*z - x = 0
z = (-2 +- sqrt(4 + 4*x^2))/(2*x)
z = (-2 + sqrt(4 + 4*x^2))/(2*x) = (-1 + sqrt(1 + x^2))/x
z = (-1 + sqrt(1 + x*x))/x
y/2 = arctg(z)
y = 2*arctg(z)
Покажем, что z < 1 для любого x > 0
z = (-1 + sqrt(1 + x*x))/x < (-1 + sqrt(1 + 2*x + x*x))/x =
= (-1 + x + 1)/x = 1
exp(x + n) = e^n * e(x)
Вычисление ряда для экспоненты сводится к суммированию ряда
для |x| < 1
x^a = ?
x = e^(ln(x))
x^a = e^(a*ln(x))
Вычисление степени сводится к вычислению экспоненты и вычислению
логарифма.
ln(1 + x) = x - x^2/2 + x^3/3 - x^4/4 + x^5/5 - ...
Сходится при -1 < x <= 1
=============================================
8 Dec 2021
Требования к задачам.
Последовательности:
1) запрещено производить вычисления внутри функции main.
Задача должна решаться с помощью отдельной функции, которой
передается указатель на открытый файл. Прототип этой функций
должен быть описан в начале файла с программой. Файл открывается в
функции main и затем она вызывает функцию, решающую нашу задачу.
После вызова этой функции и получения от нее ответа функция main
закрывает файл и печатает ответ.
#include <stdio.h>
#include <stdlib.h>
double meanValue(FILE* f);
int main() {
FILE* f = fopen("input.txt", "rt");
if (f == NULL) {
perror("Cannot open an input file");
return (-1);
}
double res = meanValue(f);
fclose(f);
printf("Mean arithmetic = %f\n", res);
return 0;
}
double meanValue(FILE* f) {
double s = 0.;
int n = 0;
while (true) {
double x;
if (fscanf(f, "%lg", &x) < 1)
break;
s += x;
++n;
}
return s/double(n);
}
Задачи по массивам.
Требования:
1) функция main не должна содержать никаких вычислений.
Ваша главная задача -- написать функцию, которой передается
адрес начала массива и его длина. Функция должна решить требуемую
задачу и вернуть результат. Задача может состоять в изменении
самого массива либо в вычислении некоторой функции от элементов
массива, либо в формировании другого массива и т.п.;
2) прототип этой функции обязательно должен быть описан
в начале текста программы;
сама функция может быть реализована где угодно, но естественнее ее
написать после main;
3) VLA (Variable Length Arrays) запрещены!
void f(int n, ...) {
double a[2*n + 1]; // Создается в стеке => должен быть небольшого
. . . // размера. В современном C++ запрещено!
Пример решения задачи по массивам
Переставить элементы массива вещественных чисел
по убыванию их модулей (сортировка по абсолютной величине).
#include <stdio.h>
#include <math.h>
void sortAbs(double *a, int n);
int main() {
FILE *f = fopen("input.txt", "rt");
if (f == NULL) {
perror("Cannot open an input file");
return (-1);
}
int n;
if (fscanf(f, "%d", &n) < 1) {
perror("Incorrect format of input file");
fclose(f);
return (-1);
}
double* a = new double[n];
for (int i = 0; i < n; ++i) {
if (fscanf(f, "%lg", a + i) < 1) {
fprintf(stderr, "Incorrect format of input file\n");
fclose(f);
return (-1);
}
}
fclose(f);
sortAbs(a, n);
FILE* g = fopen("output.txt", "wt");
if (g == NULL) {
perror("Could not open an output file");
return (-1);
}
for (int i = 0; i < n; ++i) {
// Write 5 number in a row
if (i > 0) {
// Write a delimiter
if (i%5 == 0) {
fprintf(g, "\n");
} else {
fprintf(g, " ");
}
}
fprintf(g, "%g", a[i]);
}
fprintf(g, "\n"); // LF must be at the end of any text file!!!
fclose(g);
return 0;
}
void sortAbs(double *a, int n) {
// Using the primitive bubble sort
bool sorted = false;
int pass = 0;
while (!sorted) {
sorted = true;
for (int i = 0; i < n - pass; ++i) {
if (fabs(a[i]) < fabs(a[i+1])) {
sorted = false;
// Swap these 2 elements
double tmp = a[i];
a[i] = a[i+1];
a[i+1] = tmp;
}
}
++pass;
}
}
ПЕРЕРЫВ 10 минут до 14:25
Задачи на вычисление стандартных функций
Вычислить arctg(x).
Вычислить arcsin(x).
Вычислить arccos(x).
Вычислить ln(x).
Ё
Вычислить x
Идея: ряды Тейлора либо сходятся не для всех значений х
(например, для arctg(x) или ln(x)), либо для больших x
из-за округлений при вычислении дают неверный ответ.
Пример: при суммировании ряда для sin(50) получаем ответ,
равный по абсолютной величине в несколько тысяч.
Поэтому надо суммировать ряд только для небольших знаячений x,
жедательно |x| < 1 (тогда, если коэффициенты с некоторого момента <= 1,
ряд сходится быстрее геометричестреской прогрессии
s = 1 + x + x^2 + x^3 + ... при |x| < 1
сходится, и сумма равна 1/(1 - x).
Ряд для arctg(x)
arctg(x) = x - x^3/3 + x^5/5 - x^7/7 + x^9/9 - ...
Сходится для -1 < x <= 1 (для x = 1 сходится КАТАСТРОФИЧЕСКИ медленно!).
Можно считать только для |x| < 1.
Для сведения к небольшому аргументу мы используем формулу
(она выводится из формулы для выражения tg(x) через tg(x/2)):
arctg(x) = 2*arctg(y)
y = (sqrt(1 + x*x) - 1)/x
Для любых x справедливо неравенство |y| < 1
Достаточно уметь суммировать ряд для |y| < 1 и уметь вычислять sqrt(x)
Квадратный корень вычисляется с помощью итераций Ньютона (формулу
для sqrt(a)
x_n+1 = (x_n + a/x_n)/2
проходять в школе, она была известна несколько тысяч лет назад
в древнем Вавилоне).
Как вычислит arcsin(x) или arccos(x)?
Тоже можно разложить в ряд, но лучше свести это вычисление
к вычислению arctg(x), потому что ряд для arcsin или arccos
сходится безобразно медленно и неточно для значений x, близких
к +-1. (Производная в этих точках стремится к бесконечности,
касательная к графику в этих точках вертикальна.)
Вообще говоря, функции arcsin или arccos нельзя точно вычислить!
Поэтому в реальных приложения для точного вычисления угла alpha надо знать
и sin(alpha), и cos(alpha), и вычислять угол, использую arctg(x)
Как свести вычисление arcsin(x) к arctg(z)?
Имеем
y = arcsin(x)
sin(y) = x => cos(y) = sqrt(1 - sin(x)^2) = sqrt(1 - x^2) =>
=> tg(y) = sin(y)/cos(y) = x/sqrt(1 - x^2)
y = arctg(x/sqrt(1 - x^2))
arcsin(x) = arctg(x/sqrt(1 - x^2))
Получается, что мы свели вычисление arcsin к arctg, значит,
значит 2 и 3 сложнее, чем задача 1.
1. Вычислить arctg(x).
2. Вычислить arcsin(x).
3. Вычислить arccos(x).
4. Вычислить ln(x).
Ё
5. Вычислить x
Поэтому в задачах 2 и 3 можно считать, что функцию arctg(x) мы уже
умеем вычислять (т.е. можно пользоваться стандартной функций atan
из математической библиотеки), т.е. вам надо написать лишь функцию
sqrt и сведение arcsin к arctg.
Задача 4. Вычислить ln(x)
Мы рассматривали решение этой задачи применением схемы
построения цикла с помощью инварианта. Альтернатива -- сосчитать
сумму ряда.
ln(1 + x) = x - x^2/2 + x^3/3 - x^4/4 + x^5/5 - x^6/6 + ...
сходится для 0 < x <= 1 (для x=1 сходится КАТАСТРОФИЧЕСКИ медленно).
Но можно легко свести задачу вычисления ln(x) для произвольного
x > 0 к вычислению ln(1 + y) при -0.8 < y < 0.8, используя простое
равенство:
ln(a*b) = ln(a) + ln(b)
ln(x) = ln(x/e) + 1 применяем, когда x > 1.8
ln(x) = ln(x*e) - 1 применяем, когда x < 0.2
Таким образом, загоним x в интервал 0.2 <= x <= 1.8
Проложим y = x - 1
и суммируем ряд для 1 + y
Задача 5.
Ё
Вычислить x
Сводится к вычислению экспоненты и логарифма. Имеем x = e^ln(x)
x^Ё = (e^ln(x))^Ё = e^(Ё*ln(x)) = exp(Ё*ln(x))
Как вычислит exp(x)?
exp(x) = 1 + x + x^2/2! + x^3/3! + x^4/4! + ...
Теоретически ряд сходится для любых x, но практически, как и в случае
sin(x), его можно суммировать только для небольших по абсолютной
величине значений x. Как свести задачу к небольшим x?
e^(y + k) = e^y * e^k, где k -- целое число.
Пусть x = y + k, k = int(x), |y| < 1.
exp(x) = exp(y) * e^k (если k < 0, то вычисляем (1/e)^(-k)).
Получается, что задача 5 сложнее задачи 4. Поэтому при ее решении
можно считать, что вычисление ln(x) уже реализовано, т.е.
можно использовать библиотечную функцию log(x).
Напишем вычисление exp(x)
#include <stdio.h>
#include <math.h>
#include <assert.h>
const double E = 2.7182818284590452354;
double fastPower(double x, int y);
double myExp(double x, double eps = 1e-15);
int main() {
while (true) {
double x;
printf("x: ");
if (scanf("%lg", &x) < 1)
break;
double y = myExp(x);
printf("My exp(x) = %.12f\n", y);
printf("Standard exp(x) = %.12f\n", exp(x));
}
return 0;
}
double myExp(double x, double eps) {
int k = int(x);
double y = x - double(k);
assert(fabs(y) <= 1);
// assert(x == y + k);
// exp(x) = exp(y)*e^k
double factor = fastPower(E, k);
// Compute the sum pf power series
// exp(y) = 1 + y + y^2/2! + y^3/3! + ...
double s = 0.;
double a = 1.;
double n = 0.;102_2021.txt
while (fabs(a) > eps) {
s += a;
n += 1.;
a *= y/n;
}
return s*factor;
}
double fastPower(double x, int y) {
// x^y
if (y < 0) {
y = (-y);
x = 1./x;
}
double p = 1.;
double b = x;
int k = y;
assert(k >= 0);
// Invariant: p*b^k = x^y
while (k > 0) {
if (k%2 == 0) {
k /= 2;
b *= b;
} else {
k -= 1;
p *= b;
}
}
return p;102_2021.txt
}
===================================
15 Dec 2021
Решения задач
Общие требования:
1) функция main не должна выполнять содержательных действий.
Только ввод исходных данных (или открытие файла),
вызов функции, решающей нужную задачу, и печать ответа
или запись ответа в выходной файл. Решение всегда должно
быть оформлено в виде функции;
2) функция не должна ничего сама печатать, ее задача -- вычислить
и вернуть нужное значение или преобразовать массив и т.п.
НИКАКИХ ПЕЧАТЕЙ ВНУТРИ ФУНКЦИИ БЫТЬ НЕ ДОЛЖНО!
3) ужасно! Некоторые студенты используют VLAs -- Variable Length Arrays.
Это неудачное расширение языка С (введены в стандарте C99,
не поддерживаются в современном стандарте языка C++).
int f(int n, double *a) {
. . .
double a[2*n + 1]; // Variable Length Array
. . .
}
Массив a создается в стеке.
4) оформление текста:
-- длина строки текста не должна превышать 70 символов;
std::cout << "Результат работы моей программы: "
<< "a = " << a[0] << ", " << a[1]
<< std::endl;
-- любую печать надо завершать переводом строки:
printf("x = %g\n", x);
-- ОБЯЗАТЕЛЬНО последний символ в текстовом файле -- это
символ перевода строки \n