数值分析第一次作业

韩旭 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$