数值分析第一次作业
韩旭 2020141210183
2022-09-23
题目1.
$(1.1)_{10}=(1.00011001100)_2$
用IEEE754单精度浮点数规则的表达形式为0011100011001100
最小正正规范数的表达式为0000100000000000
结果为$(-1)^0\times2^{(0001)_2-7}\times(1.00000000000)_2=1\times2^{-6}$
最大正正规范数的表达式为0111011111111111
结果为$(-1)^0\times2^{(1110)_2-7}\times(1.11111111111)_2=1.99951171875\times2^7$
题目2.
观察到在Python中
print(0.1+0,2==0.3)
>False
这就是64位计算机中无法严格地将二进制数和十进制数一一对应,导致将$0.1$从二进制转换为十进制时,只能截取其二进制表示小数点后前53位,使其无限趋近于$0.1$。我们将其打印出来发现,$0.1$在十进制小数点后第16位出现了不确定小数。因此可以估计机器精度就在$10^{-16}$左右。于是编写程序,囿于电脑内存的限制,只能精确机器精度到小数点后第24位,实现程序如下:
# 机器精度的另外一种定义是满足float(1+e)>1的最小的e
import numpy as np
# 定义函数以小数点形式输出浮点数
def float_to_str(f):
float_string = repr(f)
if 'e' in float_string:
digits, exp = float_string.split('e')
digits = digits.replace('.', '').replace('-', '')
exp = int(exp)
zero_padding = '0' * (abs(int(exp)) - 1) # minus 1 for decimal point in the sci notation
sign = '-' if f < 0 else ''
if exp > 0:
float_string = '{}{}{}.0'.format(sign, digits, zero_padding)
else:
float_string = '{}0.{}{}'.format(sign, zero_padding, digits)
return float_string
# 囿于电脑内存的限制,只能精确机器精度到小数点后第24位
for e in np.arange(10**(-24), 10**(-15), 10**(-24)):
if float(1+e)>1:
eps=e
break
else:
continue
print(eps)
print(float_to_str(eps))
输出答案为:
>1.11022303e-16
>0.000000000000000111022303
题目3.
import math
for x in [10**3, 10**5, 10**7, 10**9, 10**11]:
# print(x-math.sqrt(x*x-1))
ans = math.log(x-math.sqrt(x*x-1))
print(ans)
>-7.600902209454717
>-12.206073762186564
>-16.805431370234086
>ValueError: math domain error
改进算法,将分子有理化,求函数值$f(x)=-ln(x+\sqrt{x^2+1})$
for x in [10**3, 10**5, 10**7, 10**9, 10**11]:
# print(x+math.sqrt(x*x-1))
ans = -math.log(x+math.sqrt(x*x-1))
print(ans)
>-7.600902209541989
>-12.206072645505174
>-16.811242831518264
>-21.416413017506358
>-26.021583203494448
题目4.
import numpy as np
import math
def g(x):
ans = (math.exp(1 + float(x)) - math.exp(1)) / float(x)
return ans
def cal_absloss(a, b):
return abs(a - b)
list1, list2 = [], []
for h in np.arange(1e-16, 1, 1e-16):
# 取h从1e-16直到1,计算cal_absloss
# 实际运行中可将区间(1e-16,1)分成(1e-16, 1e-11),(1e-12, 1e-7),(1e-8, 1e-3),(1e-4, 1)提高运行效率
loss = cal_absloss(g(h), math.exp(1))
list1.append(h)
list2.append(loss)
minloss=min(list2)
dict1 = dict(zip(list2, list1))
print('loss = ', minloss)
print('h = ', dict1[minloss])
>loss = 9.694911540236717e-12
>h = 8.6857e-11
题目5.
import torch
matrix1 = torch.tensor([[10, -1],
[-1, 1]])
b1 = torch.tensor([[8],
[1]])
x0 = torch.tensor([[0],
[0]])
# print(x0.shape)
x1 = x0+(b1-torch.mm(matrix1, x0))
# print(x1, x1.shape)
x2 = x1+(b1-torch.mm(matrix1, x1))
x3 = x2+(b1-torch.mm(matrix1, x2))
x4 = x3+(b1-torch.mm(matrix1, x3))
x5 = x4+(b1-torch.mm(matrix1, x4))
x6 = x5+(b1-torch.mm(matrix1, x5))
x7 = x6+(b1-torch.mm(matrix1, x6))
x8 = x7+(b1-torch.mm(matrix1, x7))
x9 = x8+(b1-torch.mm(matrix1, x8))
x10 = x9+(b1-torch.mm(matrix1, x9))
print(x10, x10.shape)
>tensor([[-3035438298],
[ 333206829]]) torch.Size([2, 1])
import numpy as np
matrix2 = np.array([[1, -0.1],
[-1, 1]])
b2 = np.array([[0.8],
[1]])
x0 = np.array([[0],
[0]])
# print(x0.shape)
x1 = x0+(b2-np.dot(matrix2, x0))
# print(x1, x1.shape)
x2 = x1+(b2-np.dot(matrix2, x1))
x3 = x2+(b2-np.dot(matrix2, x2))
x4 = x3+(b2-np.dot(matrix2, x3))
x5 = x4+(b2-np.dot(matrix2, x4))
x6 = x5+(b2-np.dot(matrix2, x5))
x7 = x6+(b2-np.dot(matrix2, x6))
x8 = x7+(b2-np.dot(matrix2, x7))
x9 = x8+(b2-np.dot(matrix2, x8))
x10 = x9+(b2-np.dot(matrix2, x9))
print(x10, x10.shape)
>[[0.99999]
[1.99998]] (2, 1)
对于第一个线性方程组,可以得到精确解$x_{10}=[-3035438298, 333206829]^T$。
对于第二个线性方程组,计算过程中会出现无限循环小数,在进制转换时会产生方法误差,将小数点后16位截断,得到解$x_{10}=[0.99999, 1.99998]^T$