This commit is contained in:
lwj
2025-06-16 22:51:11 +08:00
parent 6baaac12c5
commit 08d293657d
4 changed files with 127 additions and 0 deletions

View File

@@ -0,0 +1,49 @@
import math
# 获取下三角矩阵的索引
# 列 行
def getIndexFromDownMatrix(col, row):
if row > col:
row, col = col, row
return (1+col)*col//2+row
# 平方根法求解
def SqrtSolve(A,b):
n = len(b)
L = [0]*(n*(n+1)//2)
for j in range(n):
for i in range(j,n):
if i == j:
L[getIndexFromDownMatrix(i,j)] = A[getIndexFromDownMatrix(i,j)]
for k in range(j):
L[getIndexFromDownMatrix(i,j)] -= L[getIndexFromDownMatrix(j,k)]**2
L[getIndexFromDownMatrix(i,j)] = L[getIndexFromDownMatrix(i,j)]**0.5
else:
L[getIndexFromDownMatrix(i,j)] = A[getIndexFromDownMatrix(i,j)]
for k in range(j):
L[getIndexFromDownMatrix(i,j)] -= L[getIndexFromDownMatrix(i,k)]*L[getIndexFromDownMatrix(j,k)]
L[getIndexFromDownMatrix(i,j)] /= L[getIndexFromDownMatrix(j,j)]
print(L)
for i in range(n):
for k in range(i):
b[i] -= L[getIndexFromDownMatrix(i,k)]*b[k]
b[i] /= L[getIndexFromDownMatrix(i,i)]
for i in range(n-1,-1,-1):
for k in range(i+1,n):
b[i] -= L[getIndexFromDownMatrix(k,i)]*b[k]
b[i] /= L[getIndexFromDownMatrix(i,i)]
return b
#把A,b换成题干的数值###########################################
if __name__ == "__main__":
##########################################################################
# 储存下三角矩阵 a11, a21, a22, a31, a32, a33 ...,按照这个格式写
A = [4,2,2,-2,-3,14]
b = [10,5,4] # b向量常数项
print("x:")
print(SqrtSolve(A,b))

View File

@@ -0,0 +1,35 @@
#模 范数
def Norm(x,v):
if len(x[0]) == 1:
if v == 1:
return sum([abs(i[0]) for i in x])
elif v == 2:
return (sum([i[0]**2 for i in x]))**0.5
elif v == float("inf"):
return max([abs(i[0]) for i in x])
else:
if v == 1:
return max([sum([abs(x[i][j]) for i in range(len(x))]) for j in range(len(x[0]))])
elif v == float("inf"):
return max([sum([abs(i) for i in x[j]]) for j in range(len(x))])
return None
# 计算矩阵的点积
def Dot(A,B):
if len(A[0]) != len(B):
return None
return [[sum([A[i][j] * B[j][k] for j in range(len(A[0]))]) for k in range(len(B[0]))] for i in range(len(A))]
#改成题干的矩阵##########################
if __name__ == "__main__":
A = [[1, 3], [-2, 4]] # 注意储存形式
x = [[1], [-1]] # 注意储存形式,单独一列的相量,每个数字都要中括号
print("Norm of x with v=1:", Norm(x, 1))
print("Norm of x with v=inf:", Norm(x, float("inf")))
print("Norm of x with v=2:", Norm(x, 2))
print("Dot product of A and x:", Dot(A, x))
print("Norm of Ax with v=2:", Norm(Dot(A, x), 2))
print("Norm of A with v=inf:", Norm(A, float("inf")))
print("Norm of A with v=1:", Norm(A, 1))

View File

@@ -0,0 +1,42 @@
import math
# 追赶法求解
def ZGsolve(A,b):
n = len(b)
beta = [0]*n
for i in range(n):
if i == 0:
beta[i] = A[i][2] / A[i][1]
else:
beta[i] = A[i][2] / (A[i][1] - A[i][0]*beta[i-1])
print("beta:")
print(beta[:-1])
for i in range(n):
if i == 0:
b[i] = b[i] / A[i][1]
else:
b[i] = (b[i] - A[i][0]*b[i-1]) / (A[i][1] - A[i][0]*beta[i-1])
print("y:")
print(b)
for i in range(n-2,-1,-1):
b[i] = b[i] - beta[i]*b[i+1]
return b
if __name__ == "__main__":
##########################################################################
#把A,b换成题干的数值###########################################
# 储存追赶法A矩阵
A = [
[0,4,-1],
[-1,4,-1],
[-1,4,-1],
[-1,4,-1],
[-1,4,0]
]
# A第一行前面有个0最后一行后面有个0
b = [100,200,200,200,100]
print("x:")
print(ZGsolve(A,b))

View File

@@ -50,6 +50,7 @@ def MullerSolve(fx,x0,x1,x2,err1,err2,N):
if __name__ == "__main__": if __name__ == "__main__":
##############################################################################################################
err1 = 1e-5 # 精度要求 P152 err1 = 1e-5 # 精度要求 P152
err2 = 1e-5 # 精度要求 P152 err2 = 1e-5 # 精度要求 P152
N = 100 # 最大迭代次数 N = 100 # 最大迭代次数