Разложение векторов по базису
Федеральное государственное образовательное бюджетное учреждение
Высшего образования
«ФИНАН
«ФИНАНСОВЫЙ УНИВЕРСИТЕТ ПРИ ПРАВИТЕЛЬСТВЕ
РОССИЙСКОЙ ФЕДЕРАЦИИ»
Департамент анализа данных, принятия решений и финансовых технологий
С.А.Зададаев
Расширение double-арифметики(R S tudio)
Учебно-методические рекомендации для проведения
семинара №30 по компьютерному практикуму
Для бакалавров направления 38.03.01 «Экономика»
Электронное издание
Москва 2018
Расширение double-арифметики
В данной части мы познакомимся с вычислениями в семантике обыкновенных дробей, которые позволяют сохранить точность результатовпри оперировании с целыми или рациональными числами.
Алгебраически точное решение матричных уравнений*
В языке Rразработана специальная библиотека gmp, которая позволяет оперировать не с десятичным форматом чисел, а с алгебраическим.
Поясним смысл сказанного, предварительно загрузив библиотеку gmpчерез меню Tools/InstallPackages… и активировав ее кодом:
library(gmp) # Библиотека расширенных типов целых и рациональных чисел
Теперь попробуем продемонстрировать разницу между обычной doubleарифметикой as . numericи расширенной с помощью специальных типов as . bigzи as . bigq.
Задание 1.Вычислить с максимальной точностью число 
.
Решение.Подключим максимальную выводимую на экран точность doubleарифметики:
|
|
|
options(digits = 22)# Выводим на экран консоли максимальное число разрядов
и вычислим число :
2^200
с результатом, записанным с помощью digits= 22 символов:
>options(digits = 22)>2^200[1] 1.6069380442589903e+60 Это означает, что число мы представляем себе примерно, как:
–
число из 61 цифры, первые 17 из которых мы видим, как значащие.
А остальные 44 цифры, вместо которых стоят нули, чему равны? Ответить на это вопрос здесь невозможно, т.к. мы использовали обычную doubleарифметику, и памяти, выделяемой на такие числа, недостаточно для получения более точного ответа.
Попробуем теперь воспользоваться специальным типом данных bigz, который выделяет под такие вычисления практически неограниченную память в рамках физически доступной на каждом конкретном компьютере (bigz – дословно «большие целые»):
as.bigz(2^200) # Вычислить выражение со всеми значащими цифрами
as.bigz(2)^200 # То же самое
с абсолютно точным результатом:
>as.bigz(2^200) # Вычислить выражение со всеми значащими цифрами BigInteger ('bigz') :[1] 1606938044258990275541962092341162602522202993782792835301376вкотором хорошо видны все 61 значащие цифры числа!
Далее прокомментируем идею абсолютной точности рациональных чисел в алгебраической форме, т.е. чисел, которые можно представить в виде несократимой дроби 
.
|
|
|
Потеря точности при вычислении таких дробей связана с возникновением периодических десятичных дробей, которые неизбежно приходится где-то обрывать. Однако, если результаты, включая промежуточные, представлять в виде обыкновенных дробей, то точность будет сохраняться – эту идеюкак раз и используетв своих алгоритмах специальный тип данных bigq(bigq – дословно «большие рациональные»).
Задание 2. Для матрицы A
найти обратную матрицу с максимальной точностью.
Решение.Считаем матрицу из Excel (см. Рис.1)
Рис.1
# В excel сформировать таблицу чисел и скопировать в буфер обмена
Data<- read.table("clipboard",h=FALSE,dec=",",sep = "\t") # Чтение из буфера обмена excel-формата
A<-data.matrix(Data) # Объявить таблицу чисел матрицей в R
A # Посмотреть результат
и для сравнения укажем, что обычная максимальная точность double-арифметикидала бы такой ответ для обратной матрицы :
solve(A) # Найти обратную матрицу в десятичном представлении
с результатом:
>solve(A) # Найти обратную матрицу в десятичном представлении [,1] [,2] [,3]V1 0.27272727272727276 0.18181818181818188 -0.09090909090909087V2 0.18181818181818188 0.45454545454545464 0.27272727272727282V3 1.00000000000000022 1.00000000000000022 1.00000000000000022В случае представления результата через обыкновенные дроби типа bigq имеем:
|
|
|
solve(as.bigq(A)) # Найти обратную матрицу в обыкновенных дробях
solve.bigq(A) # Tоже самое
срезультатом:
>solve(as.bigq(A)) # Найти обратную матрицу в обыкновенных дробяхBig Rational ('bigq') 3 x 3 matrix:[,1] [,2] [,3] [1,] 3/11 2/11 -1/11[2,] 2/11 5/11 3/11 [3,] 1 1 1 Проверим, насколько точной получится единичная матрица 
, если перемножить 
:
A%*%solve(as.bigq(A)) # Проверка точности
с абсолютно точным результатом:
>A%*%solve(as.bigq(A)) # Проверка точностиBig Rational ('bigq') 3 x 3 matrix:[,1] [,2] [,3][1,] 1 0 0 [2,] 0 1 0 [3,] 0 0 1
Задание 3. Решить систему уравнений
в обыкновенных дробях.
Решение. Матрица коэффициентов системы совпадает с матрицей A из предыдущего примера
,
а столбец правой части 
.
Для сравнения приводим два решения: приближенное и точное:
B<- 1:3 # Задаем вектор B(1, 2, 3)
solve(A, B) # Находим решение системы AX=B в десятичном представлении
solve.bigq(A, B) # Находим решение системы AX=Bв обыкновенных дробях
A%*%solve.bigq(A, B) # Проверяем второй результат (должен получиться столбец B)
Результаты приведены на рис. 2
Рис.2
Разложение векторов по базису
Теперь мы готовы найти точное разложение векторов по базису для случая, если координаты всех векторов заданы в виде целых чисел или посредством обыкновенных дробей (рациональные числа).
|
|
|
Задание 4.Разложить по базису
a) вектор 
b) вектор 
Решение.Для того, чтобы найти разложение вектора 
по базису 
, 
, 
, необходимо найти коэффициенты 
и 
в уравнении:
.
Если теперь вместо каждого вектора подставить столбец его координат, то мы получим алгебраическую систему уравнений 
с матрицами,составленными из столбцов координат соответствующих векторов. А именно:
, 
и 
.
Поступая аналогично предыдущим вычислениям, образуем в Excelрасширенную матрицу W и считываем ее в матрицы A и B:
# Работа с буфером обмена
Data<- read.table("clipboard", h=FALSE, dec=",", sep = "\t") # Чтение из буфера обмена excel-формата
W <- data.matrix(Data) # Объявить таблицу чисел матрицей в R
A <- W[ ,1:3]; A # Считываем матрицу A системы
B <- W[ ,4]; B # Считываем правый столбец системы
с результатом, представленным на рисунке 3
Рис. 3
Далее остается запустить операцию решения системы:
solve(A, B) # Приближенно решаем систему AX=B
solve.bigq(A, B) # Точно решаем систему AX=B
с результатом для обоих случаев:
>solve(A, B) # Приближенно решаем систему AX=B
V1 V2 V3
3.0000000000000000 6.0000000000000000 -4.9999999999999991
>solve.bigq(A, B) # Точно решаем систему AX=B
Big Rational ('bigq') 3 x 1 matrix:
[,1]
[1,] 3
[2,] 6
[3,] -5
Таким образом, искомое разложение по базису принимает вид:
Аналогично, для вектора получим:
B <- c(6, -3.5, 14.5); B
solve(A, B) # Приближенно решаем систему AX=B
solve.bigq(A, B) # Точно решаем систему AX=B
> B <- c(6, -7/2, 29/2); B[1] 6.0 -3.5 14.5
>solve(A, B) # Приближенно решаем систему AX=B
V1 V2
5.0000000000000000e-01 5.0000000000000000e-01
V3
-2.0438670584833814e-17
>solve.bigq(A, B) # Точно решаем систему AX=B
Big Rational ('bigq') 3 x 1 matrix:
[,1]
[1,] 1/2
[2,] 1/2
[3,] 0
с ответом:
Заключение
Заметим, что если в координатах вектора указать обыкновенную дробь, например, для вектора 
запустив строку кода:
b<- c(1/3, 2/7, 4/5)
b
>b<- c(1/3, 2/7, 4/5);
>b
[1] 0.3333333 0.2857143 0.8000000
то последующее решение solve.bigq(A, b) уже не будет точным, т.к. дроби 1/3 и 2/7успеют потерять точность после округления в десятичном представлении double-арифметики.
Для сохранения точности необходимо такой вектор 
изначально заявлять как bigq:
b<- c(as.bigq(1,3),as.bigq(2,7),as.bigq(4,5))
b
>b<- c(as.bigq(1,3),as.bigq(2,7),as.bigq(4,5))
>b
Big Rational ('bigq') object of length 3:
[1] 1/3 2/7 4/5
Если же координаты вектора заданы как иррациональные числа, например, 
, то ответомsolve.bigq(A, b) станут обыкновенные дроби, наиболее близкие к истинным иррациональным результатам. То есть сам по себе ответ в виде обыкновенных дробей еще не говорит об абсолютной точности! Проверяйте результаты прямыми подстановками.
Замечание .Важно знать, что функция solve.bigqиногда отказывается работать с вполне хорошими невырожденными матрицами, ошибочно объявляя их вырожденными. В такие моменты всегда непосредственно проверяйте равенство нулю определителя det(A). В случае ошибки остается использовать лишь приближенное решение с помощью функции solve, которая с практической точки зрения ничем не уступаетsolve.bigq.
Впрочем, нет сомнений, что со временем подобная ошибка будет устранена авторами библиотеки.
Дата добавления: 2021-07-19; просмотров: 171; Мы поможем в написании вашей работы! |
Мы поможем в написании ваших работ!