After repeatedly running out of paper, Cabi finally decided that relying on physical paper was fundamentally unsafe.
A sheet of paper may be lost, torn apart, given to suffering students during a lab, or destroyed by an unexpected earthquake.
Sometimes, Cabi falls down and causes a huge earthquake, which may destroy his remaining paper. Therefore, Cabi began designing an earthquake-resistant storage facility exclusively for his remaining paper.
Unfortunately, Cabi had no idea how to determine whether the structure would survive an earthquake.
To analyze the building, engineers model the displacement of a solid object using the equations of elastodynamics: , where
After studying for the equation for several minutes, Cabi successfully understood none of it ._. .
Fortunately, a mysterious engineer had already transformed the entire problem into a system of linear equations: .
The matrix is large, sparse, symmetric, and positive-definite.
Cabi first tried to solve the system using the Conjugate Gradient algorithm that you provided for him. However, the convergence was not fast enough to guarantee the safety of his paper before the next time he fell down.
He therefore decided to use a more powerful method: the Preconditioned Conjugate Gradient algorithm, also known as PCG.
Unfortunately, implementing it required paper to calculate the preconditioner.
Cabi had no paper.
Therefore, this problem now belongs to you.
Now, your task is to help Cabi solve the problem. Write a program to solve the linear system using the Preconditioned Conjugate Gradient (PCG) algorithm.
Unlike Gaussian Elimination, CG is an iterative method perfect for large, sparse, and symmetric positive-definite matrices. However, there is an even stronger version of CG called Preconditioned Conjugate Gradient (PCG), which uses a preconditioner to improve the convergence of the algorithm.
Here is how it works:
Prepare a starting vector called , the matrix and the vector .
Initiate variables: let
Iterate the below for
As a kind TA, Cabi has already written a basic structure for you, but some of the functions are still missing. Your task is to complete those functions.
N:The dimension of the linear system. This is also the size of the vectors and the matrixstep: the number of iterations you need to perform for the PCG algorithm.d: a list of length . d[i] represents the diagonal element of matrix at row (i.e.,).pos: a list of lists. pos[i] contains the column indices of the non-zero off-diagonal(非對角線) elements in row of the matrix .val: a list of lists. val[i] contains the values of the non-zero off-diagonal elements in row , corresponding to pos[i] of the matrix .b: the right-hand side of the equation of length .x0: the initial guess vector of length .Suppose we have a matrix :
The inputs for this matrix will look like this:
4 4 5 (this is d)
1 1 -1 (Row 0: 1 off-diagonal element at column 1, value is -1)
2 0 2 -1 -2 (Row 1: 2 off-diagonal elements at column 0 and 2, values are -1 and -2)
1 1 -1 (Row 2: 1 off-diagonal element at column 1, value is -1)
And for such it is
To help you build the PCG algorithm, you need to implement the following helper functions:
dot(v1, v2): Returns the dot product of two vectors v1 and v2(i.e.,,)vec_add(v1, v2): Returns a new vector which is the element-wise addition of v1 and v2.vec_sub(v1, v2): Returns a new vector which is the element-wise subtraction of v1 and v2.vec_const_mul(c, v): Returns a new vector where each element of v is multiplied by a constant c.spmul(d, pos, val, v): Performs sparse matrix-vector multiplication and returns a new vectorPCG(A_d, A_pos, A_val, b, x0, step): The main Preconditioned Conjugate Gradient loop. Use the formulas provided above and your helper functions to find
def spmul(d:list,pos:list,val:list,v:list)->list:
#TODO
pass
def dot(v1:list,v2:list)->float:
#TODO
pass
def vec_add(v1:list,v2:list)->list:
#TODO
pass
def vec_sub(v1:list,v2:list)->list:
#TODO
pass
def vec_const_mul(c:float,v:list)->list:
#TODO
pass
def PCG(A_d:list,A_pos:list,A_val:list,b:list,x0:list,step:int)->list:
#TODO
pass
import sys
input = sys.stdin.readline
if __name__=='__main__':
n,step=map(int,input().split())
b=list(map(int,input().split()))
x0=list(map(int,input().split()))
d=list(map(int,input().split()))
pos=[ [] for _ in range(n) ]
val=[ [] for _ in range(n) ]
for i in range(n):
tmp_list=list(map(int,input().split()))
m=tmp_list[0]
if m ==0 : continue
tmp=tmp_list[1:]
pos[i]=tmp[0:m]
val[i]=tmp[m:2*m]
ans_x = PCG(d,pos,val,b,x0,step)
for x in ans_x:
print(f'{x:.5f}',end=' ')
the first line contains two integers:
N step
The line contains integers, representing the vector b.
The line contains integers, representing the initial guess vector x0.
The line contains integers, representing d(the diagonal elements of ).
The next lines describe the off-diagonal non-zero elements for row to row . Each line follows this format:
M_i pos_i_0 pos_i_1 ... pos_i_{M-1} val_i_0 val_i_1 ... val_i_{M-1}
If the line will just be 0.
Output the final vector .
Print the elements of the vector separated by a single space. Please format each number to decimal places.