Preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication
By introducing the improved preprocessing pipeline co-gradient algorithm PIPEPCG in the parallel ocean program POP2, the problem of global reduction operation in the classic PCG algorithm becomes a communication bottleneck is solved, and more efficient parallel performance and simulation speed are achieved.
Patent Information
- Application Number
- CN202411627592.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2024-11-14
- Publication Date
- 2025-03-07
AI Technical Summary
In parallel ocean program POP2, the classic preprocessing conjugate gradient algorithm PCG becomes a communication bottleneck due to global reduction operations, affecting the simulation speed and accuracy.
An improved preprocessing pipeline conjugate gradient algorithm PIPEPCG is proposed. Only one global reduction operation is performed in each iteration through internal product merge, and the blocking global reduction function is replaced with non-blocking. The global reduction operation is performed overlapping with matrix vector multiplication, boundary update operation and preprocessing operation.
It significantly reduces global communication, improves parallel performance, accelerates the simulation speed of POP2 models, and provides reference for the optimization of other climate models.
Smart Images
Figure CN119576417A_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the technical field of high performance computing, and in particular relates to a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication. Background Art
[0002] High-resolution global climate models have become key tools for studying the complexity of the climate system and predicting future climate change trends. The Common Earth System Model (CESM) is one of the most widely used global climate models, including several modules such as the atmosphere, land, ocean, sea ice, and land ice. It is mainly used to study the past, present, and future climate conditions of the Earth. The parallel ocean program POP2 is an important component of CESM, responsible for simulating the physical processes of the ocean and its interactions with other climate system components such as the atmosphere, land, sea ice, and land ice. For example, in the positive pressure mode of POP2, the main task is to solve a set of elliptic equations about sea level height (SSH).
[0003] In POP2, the classic preconditioned conjugate gradient algorithm PCG is used by default, combined with a diagonal preprocessor to solve linear systems; the PCG algorithm needs to perform two inner product operations in each iteration. As the number of processes increases, the global reduction operation will become the main communication bottleneck, affecting the efficiency and accuracy of the solution, and thus affecting the user experience of the global climate model. Summary of the invention
[0004] In order to overcome the shortcomings of the above-mentioned prior art, the present invention provides a preprocessing pipeline conjugate gradient algorithm optimization method that reduces global communication, significantly improves and implements a variant of the PCG algorithm, namely the PIPEPCG algorithm, and applies the algorithm to the POP2 model to improve the simulation speed of the POP2 model in the CESM framework; compared with the classic PCG algorithm, the PIPEPCG algorithm significantly reduces the global reduction operations, and further improves the parallel performance in combination with the EVP preprocessor; the new algorithm not only accelerates the simulation speed of the POP2 model, but also provides a reference for the optimization of other climate models.
[0005] To achieve the above objectives, one or more embodiments of the present invention provide the following technical solutions:
[0006] A first aspect of the present invention provides a preprocessing pipeline conjugate gradient algorithm optimization method that reduces global communication.
[0007] A pre-processed pipeline conjugate gradient algorithm optimization method with reduced global communication for solving sea level height in parallel ocean programs includes:
[0008] Obtain the time series of the seabed depth and sea level height of the location to be predicted, construct the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving;
[0009] The improved preprocessing pipeline conjugate gradient algorithm is used to iteratively solve each sub-equation. Inner product merging is used to perform only one global reduction operation in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, boundary update operation and preprocessing operation.
[0010] The solutions of all sub-equations are combined to obtain the sea level height at the location to be predicted.
[0011] Furthermore, the specific steps of constructing the linear equation to be solved are:
[0012] Construct implicit elliptic equations for sea level height;
[0013] Based on the time series of seabed depth and sea level height at the location to be predicted, a nine-point template is used to discretize the implicit elliptic equation on a two-dimensional grid to obtain the linear equation to be solved.
[0014] Furthermore, the division into multiple sub-equations is to divide the global domain of the linear equation into multiple blocks and assign them to various processes. Each process is responsible for processing the evolution of the grid points in its assigned block, while maintaining the boundary area for data exchange and update with neighbor processes.
[0015] Furthermore, the iterative solution is to solve the sub-equation Ax=b, where A is the coefficient matrix and b is the right-hand side vector. The specific steps are:
[0016] Initialization, including setting the coefficient matrix, right-hand vector, calculating the initial guess solution, initial residual vector, preprocessed residual vector, and matrix-vector product result vector, and defining the update_halo() function for exchanging and updating boundary area values in parallel computing;
[0017] After initialization, an iterative loop is performed to perform local inner product calculations, global reduction operations, calculation operations, calculation of directional coefficients and step sizes, update of search direction vectors and related vectors, update of solution vectors and residual vectors, residual replacement operations, and judgment of convergence conditions until the convergence conditions are met, the iteration is terminated, and the current solution is output.
[0018] Furthermore, the local inner product calculation and global reduction operation are to merge two inner product operations without dependency into one inner product operation through the inner product, so that only one global reduction operation needs to be performed in each iteration;
[0019] The two merged inner product operations are: the inner product of the residual vector and the preprocessed residual vector, and the inner product of the matrix-vector product result vector and the preprocessed residual vector.
[0020] Furthermore, the global reduction operation and calculation operation are specifically:
[0021] The global reduction operation uses a non-blocking global reduction function to aggregate the local results calculated by each process to obtain a global result, including a global result after the two inner products are aggregated separately;
[0022] The calculation operation is based on a row preprocessing operation and a boundary update operation to perform matrix-vector multiplication.
[0023] Furthermore, the preprocessing operation adopts an error vector propagation block based preprocessor to accelerate the convergence of the algorithm;
[0024] The boundary update operation updates the data exchanged between the current process and the adjacent process after calculating a new matrix-vector product result vector through matrix-vector multiplication.
[0025] A second aspect of the present invention provides a preprocessing pipeline conjugate gradient algorithm optimization system that reduces global communication.
[0026] A preprocessing pipeline conjugate gradient algorithm optimization system with reduced global communication is used to solve the sea level height in parallel ocean programs, including equation building modules, parallel solution modules and solution combination modules:
[0027] The equation building module is configured to: obtain the time series of the seabed depth and sea level height of the location to be predicted, build the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving;
[0028] The parallel solution module is configured as follows: using the improved preprocessing pipeline conjugate gradient algorithm, in the iterative solution process of each sub-equation, the inner product merging is used to perform only one global reduction operation in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, the boundary update operation and the preprocessing operation;
[0029] The solution combination module is configured to combine the solutions of all sub-equations to obtain the sea level height of the position to be predicted.
[0030] The third aspect of the present invention provides a computer-readable storage medium having a program stored thereon, which, when executed by a processor, implements the steps in a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as described in the first aspect of the present invention.
[0031] The fourth aspect of the present invention provides an electronic device, comprising a memory, a processor, and a program stored in the memory and executable on the processor, wherein when the processor executes the program, the steps in a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as described in the first aspect of the present invention are implemented.
[0032] One or more of the above technical solutions have the following beneficial effects:
[0033] The improved PIPEPCG algorithm developed and implemented by the present invention is based on the PCG algorithm. By adopting pipeline technology to reorder the algorithm process, conditions are created for optimization methods such as inner product operation merging and communication and calculation overlap. The present invention significantly improves the PIPEPCG algorithm. By merging inner products, only one global reduction operation needs to be performed in each iteration. In addition, the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that these operations can be overlapped with matrix-vector multiplication, boundary update operation and preprocessing operation. The position of boundary update is also optimized, and the algorithm performance is further improved by using methods such as residual replacement operation and vector merging. At the same time, the original diagonal preprocessor is replaced by the EVP preprocessor, thereby further improving the parallel performance of the algorithm.
[0034] Advantages of additional aspects of the present invention will be given in part in the following description, and in part will become obvious from the following description, or will be learned through practice of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS
[0035] The accompanying drawings in the specification, which constitute a part of the present invention, are used to provide a further understanding of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute improper limitations on the present invention.
[0036] Figure 1 This is a flow chart of the method of the first embodiment.
[0037] Figure 2 This is a schematic diagram of vector operation merging in the first embodiment.
[0038] Figure 3 This is a comparison chart of different solvers in the positive pressure mode of the 1° POP model of the first embodiment.
[0039] Figure 4 This is a comparison chart of different solvers in the positive pressure mode of the 0.1° POP model of the first embodiment.
[0040] Figure 5 This is a diagram showing the time consumption of different solvers in the global reduction operation and boundary update operation of the first embodiment. DETAILED DESCRIPTION
[0041] It should be noted that the following detailed descriptions are illustrative and are intended to provide further explanation of the present application. Unless otherwise specified, all technical and scientific terms used in the present invention have the same meanings as those commonly understood by those skilled in the art to which the present application belongs.
[0042] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit the exemplary embodiments according to the present application. As used herein, unless the context clearly indicates otherwise, the singular form is also intended to include the plural form. In addition, it should be understood that when the terms "comprise" and / or "include" are used in this specification, it indicates the presence of features, steps, operations, devices, components and / or combinations thereof.
[0043] Embodiment 1
[0044] In one embodiment of the present disclosure, a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication is provided, which is used to solve the sea level height in a parallel ocean program, including the following steps:
[0045] Step S1: Obtain the time series of the seabed depth and sea level height of the location to be predicted, construct the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving.
[0046] Furthermore, the specific steps of constructing the linear equation to be solved are:
[0047] Construct implicit elliptic equations for sea level height;
[0048] Based on the time series of seabed depth and sea level height at the location to be predicted, a nine-point template is used to discretize the implicit elliptic equation on a two-dimensional grid to obtain the linear equation to be solved.
[0049] Furthermore, the division into multiple sub-equations is to divide the global domain of the linear equation into multiple blocks and assign them to various processes. Each process is responsible for processing the evolution of the grid points in its assigned block, while maintaining the boundary area for data exchange and update with neighbor processes.
[0050] Step S2: Using the improved preprocessing pipeline conjugate gradient algorithm, during the iterative solution of each sub-equation, inner product merging is performed so that only one global reduction operation is performed in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, boundary update operation and preprocessing operation.
[0051] Furthermore, the iterative solution is to solve the sub-equation Ax=b, where A is the coefficient matrix and b is the right-hand side vector. The specific steps are:
[0052] Initialization, including setting the coefficient matrix, right-hand vector, calculating the initial guess solution, initial residual vector, preprocessed residual vector, and matrix-vector product result vector, and defining the update_halo() function for exchanging and updating boundary area values in parallel computing;
[0053] After initialization, an iterative loop is performed to perform local inner product calculations, global reduction operations, calculation operations, calculation of directional coefficients and step sizes, update of search direction vectors and related vectors, update of solution vectors and residual vectors, residual replacement operations, and judgment of convergence conditions until the convergence conditions are met, the iteration is terminated, and the current solution is output.
[0054] Furthermore, the local inner product calculation and global reduction operation are to merge two inner product operations without dependency into one inner product operation through the inner product, so that only one global reduction operation needs to be performed in each iteration;
[0055] The two merged inner product operations are: the inner product of the residual vector and the preprocessed residual vector, and the inner product of the matrix-vector product result vector and the preprocessed residual vector.
[0056] Furthermore, the global reduction operation and calculation operation are specifically:
[0057] The global reduction operation uses a non-blocking global reduction function to aggregate the local results calculated by each process to obtain a global result, including a global result after the two inner products are aggregated separately;
[0058] The calculation operation is based on a row preprocessing operation and a boundary update operation to perform matrix-vector multiplication.
[0059] Furthermore, the preprocessing operation adopts an error vector propagation block based preprocessor to accelerate the convergence of the algorithm;
[0060] The boundary update operation updates the data exchanged between the current process and the adjacent process after calculating a new matrix-vector product result vector through matrix-vector multiplication.
[0061] Step S3: Combine the solutions of all sub-equations to obtain the sea level height of the location to be predicted.
[0062] The following is a detailed description of the implementation process of a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication in this embodiment.
[0063] This embodiment improves the preprocessing pipeline conjugate gradient algorithm, further improves parallel performance by reducing global reduction operations and combining the EVP preprocessor, thereby accelerating the simulation speed of the POP2 model.
[0064] The Preprocessed Pipelined Conjugate Gradient (PIPEPCG) algorithm is based on the PCG algorithm and reorders the algorithm flow by adopting pipeline technology, thereby creating conditions for optimization methods such as inner product operation merging and communication and calculation overlap; this embodiment significantly improves the PIPEPCG algorithm, and through inner product merging, only one global reduction operation needs to be performed in each iteration; in addition, this embodiment replaces the originally blocked global reduction function with a non-blocking global reduction function, so that these operations can be overlapped with matrix-vector multiplication, boundary update operation and preprocessing operation; this embodiment also optimizes the location of boundary update, and uses methods such as residual replacement operation and vector merging to further improve the algorithm performance; at the same time, this embodiment replaces the original diagonal preprocessor with the EVP preprocessor, thereby further improving the parallel performance of the algorithm.
[0065] The error vector propagation (EVP) block preprocessor has been proven to be an efficient parallel preprocessor. The EVP method mainly consists of two steps: preprocessing and solving. In the preprocessing stage, the influence coefficient matrix and its inverse matrix need to be calculated, and these calculations only need to be performed once at the beginning.
[0066] Taking the positive pressure mode of POP2 as an example, the main task is to solve an elliptic equation system about sea level height (SSH). POP2 simplifies the elliptic equation system into a linear system A global x global =b global to perform the calculation.
[0067] The implicit elliptic equation in POP2 can be expressed as:
[0068]
[0069] Where H represents the seafloor depth, Φ represents the function of the time step τ, and η n-1 , η n , η n+1 is the sea level height at the n-1th, nth, and n+1th time steps, and Some function representing the previous state.
[0070] Discretize implicit elliptic equations on a two-dimensional grid using a nine-point template in POP2. and are the symmetry coefficients between the grid point (i, j) and its neighboring grid points; these coefficients are determined by the ocean depth H, the time step τ, and the grid length.
[0071] The nine-point format for defining grid points (i, j) is:
[0072]
[0073] In the global region, the nine-point template can be transformed into a linear system A globalxglobal =b global , matrix A global The coefficient A 0 , A N , A E and A NE A block tridiagonal matrix.
[0074] To achieve large-scale parallel computing, POP2 divides the global area into multiple blocks and assigns them to various processes. Each process is responsible for handling the evolution of the grid points in its assigned block, while maintaining the boundary area to exchange and update data with neighboring processes to ensure the accuracy and continuity of the calculation.
[0075] The global domain size is N×N. Assume that the global region of calculation is a two-dimensional grid with a size of N×N, which means that the region is discretized into grid points with N rows and N columns, and each grid point represents a position in the calculation domain.
[0076] The global domain is divided into m×m blocks. For parallel processing, the global grid is divided into m×m small blocks, and the size of each block is (n×n), where n=N / m, which means that each process only needs to process the grid points in a specific small block.
[0077] Block allocation, each parallel process handles one or more of these small blocks and is responsible for calculating the evolution of all grid points within the block (such as changes in physical quantities such as temperature and pressure), that is, the subtask of solving the sub-equation Ax=b.
[0078] In a single process, the linear equation Ax = b is solved using the improved preprocessing pipeline conjugate gradient algorithm, such as Figure 1 As shown, the specific steps are:
[0079] A. Getting Started:
[0080] The subtasks are started in parallel.
[0081] B. Initialization:
[0082] Initialize the following parameters: the coefficient matrix A is a sparse tridiagonal matrix, the EVP preprocessor M is used to accelerate the convergence of the algorithm, x is the initial solution, and b represents the right-hand side vector of the linear equation system Ax=b.
[0083] Calculate r0=b-Ax0, where r0 represents the initial residual vector, reflecting the gap between the initial solution and the actual solution.
[0084] The update_halo() function represents a boundary update operation, which is used to exchange and update the values of the boundary area in parallel computing to ensure consistency between different computing units; the boundary update is because different small blocks belong to different processes, and the boundary area needs to exchange data with adjacent small blocks (that is, blocks allocated to other processes) to ensure that the data is consistent when the next calculation is performed; this data exchange is usually implemented through MPI communication to update the values of each small block boundary.
[0085] Calculate u0 = M -1 r0, where u0 is the preprocessed residual vector used to speed up the convergence of the PIPEPCG algorithm.
[0086] Calculate w0=Au0, where w0 represents the product of the coefficient matrix A and the vector u0, which serves as the basis for subsequent calculations.
[0087] C. Enter the iteration loop:
[0088] while i≤i max do, i represents the current number of iterations, which is a counter that starts from 0 and increases by 1 for each iteration; i max It represents the maximum number of iterations, which is a pre-set parameter used to limit the maximum number of iterations of the algorithm to prevent it from falling into an infinite loop when convergence fails.
[0089] D. Local inner product calculation:
[0090] The two inner product operations without dependencies are merged into one inner product operation through inner products, so that only one global reduction operation needs to be performed in each iteration; wherein, the two merged inner product operations are: the inner product of the residual vector and the preprocessed residual vector, and the inner product of the matrix-vector product result vector and the preprocessed residual vector, specifically:
[0091] calculate in Represents the residual vector r in the i-th iteration on a single process i And the preprocessed residual vector u i The inner product of .
[0092] calculate in represents the vector w in the i-th iteration on a single process i And the preprocessed residual vector u i The inner product of .
[0093] Here, the improved PIPEPCG algorithm achieves optimization by merging inner product operations. In this step, there is no dependency between the two inner product operations, so they can be merged into one inner product operation without affecting the correctness of the algorithm. Each time an inner product operation is performed, a global reduction operation is required. In large-scale high-performance computing applications, the global reduction operation is one of the key bottlenecks affecting application performance. Compared with the traditional PCG algorithm, the improved PIPEPCG algorithm in this embodiment significantly improves the performance of the solver by reducing one global reduction operation, especially in the case of a high number of cores, this advantage will be more obvious.
[0094] E. Global reduction operation:
[0095] The global reduction operation uses a non-blocking global reduction function to aggregate the local results calculated by each process to obtain a global result, including the global result after the two inner products are aggregated separately, which is expressed by the formula:
[0096]
[0097] Among them, the MPI_IALLREDUCE() function indicates the use of non-blocking global reduction operations. The execution time of non-blocking global reduction operations can overlap with the execution time of preprocessing operations, boundary update operations, and sparse matrix-vector multiplication (SPMV). i Indicates the local results calculated by each process through the global reduction operation The global result after aggregation, δ i Indicates the local results calculated by each process through the global reduction operation The global results after aggregation.
[0098] Here, the improved PIPEPCG algorithm achieves the overlap of the execution time of the global reduction operation and other computing tasks.
[0099] In the classic PCG algorithm, due to the data dependency between operations, a blocking global reduction function (MPI_ALLREDUCE) is used; in order to better utilize computing resources and reduce processor waiting time, this embodiment uses a non-blocking global reduction function (MPI_IALLREDUCE) in step E of the PIPEPCG algorithm; in the improved PIPEPCG algorithm, the variable γ obtained by the global reduction operation in step E is i and δ i With vector m i and vector n i There is no data dependency between them; therefore, the execution time of the global reduction operation can overlap with the execution time of the preprocessing operation, the boundary update operation, and the sparse matrix-vector multiplication (SPMV).
[0100] F. Calculation operations:
[0101] The calculation operation is based on the row preprocessing operation and the boundary update operation, and performs matrix-vector multiplication, specifically:
[0102] Calculate m i =M -1 w i , where m i represents the preprocessed w in the i-th iteration on a single process i vector.
[0103] Call the update_halo() function to perform boundary update operations, exchanging and updating the values of the boundary area in parallel computing.
[0104] Calculate n i =Am i , where n i Represents the coefficient matrix A and vector m in the i-th iteration on a single process i The product result of .
[0105] Finally, the MPI_WAIT() function is called to block the calling process until the related non-blocking global reduction communication operation is completed; this ensures that the related communication operation has been completed and the data has been correctly transmitted when the subsequent code is executed.
[0106] G. Calculate the directional coefficient and step size:
[0107] The code is as follows:
[0108] if i>0then
[0109] β i =γ i / γ i-1 ; α i =γ i / (δ i -β i γ i / α i-1 );
[0110] else
[0111] β i =0; α i =γ i / δ i ;
[0112] end if
[0113] Among them, β iRepresents the key coefficient for updating the search direction in the PIPEPCG algorithm, which ensures the efficiency of each step of optimization and the final convergence of the algorithm; α i It means to ensure that each update of the solution vector x i When, along the search direction p i Make the best progress so that the solution vector x i At each step we gradually approach the final solution.
[0114] H. Update the search direction vector and correlation vector:
[0115] z i =n i +β i z i-1 ,q i =m i +β i q i-1 、s i =w i +β i s i-1 、p i =u i +β i p i-1 , where z, q, s
[0116] It is an intermediate variable used to maintain the conjugation between vectors, ensuring the rapid convergence of the algorithm; p is the search direction vector in the conjugate gradient method; by updating the vector p, the algorithm can get closer to the solution of the equation in each iteration.
[0117] I. Update the solution vector and residual vector:
[0118] x i+1 =x i +α i p i 、r i+1 =r i -α i s i 、u i+1 =u i -α i q i 、w i+1 =w i -α i z i , where the solution vector x is updated by i+1 , the algorithm can gradually approach the true solution of the linear equation Ax = b; the residual vector r i+1 Represents the current solution x i+1 By updating the residual, the algorithm can evaluate the accuracy of the current solution and decide whether further iterations are needed; i+1is the residual vector used to optimize the search direction after preprocessing to ensure that the algorithm converges faster; vector w i+1 Used to calculate the direction and step size in the next iteration, it is an important part of maintaining the conjugation of the vector.
[0119] Here, the improved PIPEPCG algorithm of this embodiment realizes the optimization method of vector operation merging, such as Figure 2 As shown on the left, in the original implementation of the PIPECG algorithm, each vector operation function needs to be called separately to complete the update of the solution vector and the residual vector AXPY; however, Figure 2 As shown on the right, by merging vector operations into one loop, the vectors already loaded into the cache can be reused, reducing the number of reads and writes to the main memory. At the same time, merging multiple vector operations also reduces the overhead of loop control and improves overall efficiency.
[0120] J. Residual substitution operation and convergence judgment conditions:
[0121] if i%n c = =0 then
[0122] r i+1 = b-Ax i+1 ;update_halo(r i+1 );
[0123] u i+1 =M -1 r i+1 ;update_halo(u i+1 );
[0124] w i+1 =Au i+1 ;
[0125] convergence_check(r i+1 );
[0126] end if
[0127] In the PIPEPCG algorithm, in order to reduce the accumulated numerical errors, a residual replacement strategy is adopted. Specifically, after every n c After iterations (that is, when i is n c multiples of), the algorithm will recalculate the residual vector r i+1 , preprocessed residual vector u i+1 and the matrix-vector product vector w i+1 .
[0128] In addition, the update_halo() function is executed to update the boundary values of these vectors to ensure that the data between processes remains consistent in a parallel computing environment.
[0129] The convergence_check() function is used to check the new residual vector r i Is it small enough, that is, whether the current solution is close enough to the true solution of the equation; if the convergence condition is met, the algorithm can be terminated, otherwise continue to iterate.
[0130] Here, in order to improve the maximum achievable accuracy of the PIPEPCG algorithm, this embodiment introduces a residual replacement method.
[0131] Solution vector x i and the residual vector r i The update method is: i+1 =x i +α i p i ,r i+1 =r i -α i s i ; Among them, the update of the solution vector x does not directly affect the residual vector r i , which means that in the solution vector x i The error generated during the update process will not be directly reflected in the residual vector r i middle.
[0132] As the number of iterations increases, the updated residual r i and the true residual b-Ax i+1 The deviations between them will gradually accumulate; by regularly replacing the updated residuals with the true residuals, the accumulated errors can be corrected to ensure the accuracy of the residuals.
[0133] Specifically, every certain number of iterations, the true residual r is calculated using i (i.e. b-Ax i+1 ) to replace the residual vector r accumulated during the iteration i In the improved PIPEPCG algorithm, this operation occurs in step J and involves calculating r i+1 = b-Ax i+1 、u i+1 =M -1 r i+1 and w i+1 =Au i+1 And other steps.
[0134] The frequency of residual replacement in low-resolution grids and high-resolution grids is different. Compared with low-resolution grids, the longitude-latitude spacing ratio of high-resolution POP grids is closer to 1, so the condition number of the coefficient matrix is smaller, so there is no need to perform residual replacement operations as frequently as low-resolution grids.
[0135] K. End of loop iteration
[0136] In the above process, the improved PIPEPCG algorithm also optimizes the location of the boundary update operation:
[0137] In the improved PIPEPCG algorithm, M represents the EVP preprocessor; after each EVP preprocessing operation, a boundary update operation is required; at the same time, due to the numerical instability of the boundary area when solving x, a boundary update operation is also required; therefore, compared with the PCG algorithm, the improved PIPEPCG algorithm performs an additional boundary update operation in each iteration.
[0138] This embodiment verifies the effect of the improved algorithm through experiments, specifically:
[0139] 1. Experimental environment: Shanhe supercomputer uses 2.7GHz Intel(R) Xeon(R) Gold6258R processor and is connected through 100GBps InfiniBand network.
[0140] 2. Experimental case:
[0141] To focus on the performance evaluation of POP, the “G NORMAL YEAR” component set of CESM was used; in this component set, the ocean and sea ice models are dynamically simulated, while the atmospheric model is data-driven.
[0142] The two most commonly used POP horizontal grid resolutions were adopted: 1° (360×240) and 0.1° (3600×2400).
[0143] 3. Experimental parameters:
[0144] Matrix A global The dimension is N×N, where N is 3600x2400.
[0145] Right side vector b global The dimension is N.
[0146] The dimension of the initial solution vector is N.
[0147] Maximum number of iterations i max is 2000.
[0148] The convergence criterion is 10-13 .
[0149] The number of steps n for residual replacement interval c is 20.
[0150] At low resolution (1°), 56, 112, 224, 448, and 896 kernels were used for verification, respectively.
[0151] At high resolution (0.1°), 2000, 4000, 6000, 8000, and 10000 kernels were used for verification.
[0152] 4. Experimental results:
[0153] Table 1: Differences in sea surface height (SSH) between PCG and PIPEPCG
[0154]
[0155] In order to ensure that the P IPEPCG algorithm does not cause inaccuracy in the POP model, this embodiment conducted experiments using 1° POP on the Shanhe supercomputer. Table 1 shows the difference in sea level height (SSH) calculated using the PCG algorithm and the P IPEPCG algorithm. Compared with the maximum absolute SSH, the average SSH error of the P IPEPCG algorithm is very small, showing a higher accuracy.
[0156] In the experimental comparison, it was found that the larger errors mainly occurred near the coastline. Due to the drastic changes in the boundaries, these areas are prone to instability during numerical calculations. The difference between the P IPEPCG algorithm and the PCG algorithm in the calculation of sea level height (SSH) mainly increases with the extension of the simulation period. This is due to the continuous accumulation of turbulence effects throughout the simulation process, rather than a problem of the solver's own accuracy.
[0157] When using 1°POP on the Shanhe supercomputer, the execution time of different solvers in positive pressure mode was compared. The comparison results are shown in Figure 3 As shown in the figure, compared with the PCG solver using the default diagonal preprocessor (PCG+Diagonal), the PIPEPCG algorithm using the EVP preprocessor (PIPEPCG+EVP) performs better on all core numbers; at the highest number of cores (896), the execution time of the PIPEPCG algorithm per simulation cycle is 0.999 seconds, which is 1.45 times faster than that of the PCG solver with the diagonal preprocessor.
[0158] When using 0.1° POP on the Shanhe supercomputer, the execution time of different solvers in positive pressure mode was compared. The comparison results are shown in Figure 4As shown in the figure; with the increase of the number of processes, the execution time of each solver shows an increasing trend; compared with the PCG solver using the default diagonal preprocessor (PCG+Diagonal), the PIPEPCG algorithm using the EVP preprocessor (PIPEPCG+EVP) performs better on all core numbers; at 10,000 cores, the performance of the PIPEPCG algorithm is 1.98 times higher than that of the PCG solver with the diagonal preprocessor. The use of the new algorithm significantly improves the performance of POP.
[0159] In order to explain the reason for the performance improvement in more detail, this embodiment Figure 5 The time consumption of each solver in global reduction operation and boundary update operation is shown in the figure. Figure 5 The left figure shows that the trend of the global reduction execution time of each solver is basically consistent with the trend of the overall execution time, which indicates that the global reduction operation occupies the largest proportion of the execution time under a large number of cores; the performance of the PIPEPCG algorithm is better than that of the PCG algorithm because they only need to perform one global reduction operation in each iteration. Reducing the number of global reduction operations can effectively reduce the sensitivity of POP to operating system noise; at the same time, the use of the EVP preprocessor can effectively reduce the number of iterations required by the solver, thereby further reducing the number of global reduction operations; Figure 5 The figure on the right shows that due to the use of the EVP preprocessor, the P IPEPCG algorithm performs one more boundary update operation in each iteration compared to the PCG algorithm with the default diagonal preprocessor; under a large number of cores, the computational cost of the positive pressure solver is negligible compared to the global reduction operation and the boundary update operation; therefore, the additional calculation required by the EVP preprocessor has little impact on performance.
[0160] Embodiment 2
[0161] In one embodiment of the present disclosure, a preprocessing pipeline conjugate gradient algorithm optimization system with reduced global communication is provided, which is used to solve the sea level height in a parallel ocean program, including an equation building module, a parallel solving module and a solution combination module:
[0162] The equation building module is configured to: obtain the time series of the seabed depth and sea level height of the location to be predicted, build the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving;
[0163] The parallel solution module is configured as follows: using the improved preprocessing pipeline conjugate gradient algorithm, in the iterative solution process of each sub-equation, the inner product merging is used to perform only one global reduction operation in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, the boundary update operation and the preprocessing operation;
[0164] The solution combination module is configured to combine the solutions of all sub-equations to obtain the sea level height of the position to be predicted.
[0165] Embodiment 3
[0166] The purpose of this embodiment is to provide a computer-readable storage medium.
[0167] A computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps in a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as described in Embodiment 1 of the present disclosure.
[0168] Embodiment 4
[0169] The purpose of this embodiment is to provide an electronic device.
[0170] An electronic device includes a memory, a processor, and a program stored in the memory and executable on the processor. When the processor executes the program, the steps in a preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as described in Embodiment 1 of the present disclosure are implemented.
[0171] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.
Claims
1. A preprocessing pipeline conjugate gradient algorithm optimization method to reduce global communication, characterized in that: Used to solve sea level height in parallel ocean programs, including: Obtain the time series of the seabed depth and sea level height of the location to be predicted, construct the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving; The improved preprocessing pipeline conjugate gradient algorithm is used to iteratively solve each sub-equation. Inner product merging is used to perform only one global reduction operation in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, boundary update operation and preprocessing operation. The solutions of all sub-equations are combined to obtain the sea level height at the location to be predicted.
2. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 1, characterized in that: The specific steps of constructing the linear equation to be solved are: Construct implicit elliptic equations for sea level height; Based on the time series of seabed depth and sea level height at the location to be predicted, a nine-point template is used to discretize the implicit elliptic equation on a two-dimensional grid to obtain the linear equation to be solved.
3. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 1, characterized in that: The division into multiple sub-equations is to divide the global domain of the linear equation into multiple blocks and assign them to various processes. Each process is responsible for processing the evolution process of the grid points in its assigned block, while maintaining the boundary area to exchange and update data with neighbor processes.
4. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 1, characterized in that: The iterative solution is to solve the sub-equation Ax=b, where A is the coefficient matrix and b is the right-hand side vector. The specific steps are: Initialization, including setting the coefficient matrix, right-hand vector, calculating the initial guess solution, initial residual vector, preprocessed residual vector, and matrix-vector product result vector, and defining the update_halo() function for exchanging and updating boundary area values in parallel computing; After initialization, an iterative loop is performed to perform local inner product calculations, global reduction operations, calculation operations, calculation of directional coefficients and step sizes, update of search direction vectors and related vectors, update of solution vectors and residual vectors, residual replacement operations, and judgment of convergence conditions until the convergence conditions are met, the iteration is terminated, and the current solution is output.
5. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 4, characterized in that: The local inner product calculation and global reduction operation are to merge two inner product operations without dependency into one inner product operation through the inner product, so that only one global reduction operation needs to be performed in each iteration; The two merged inner product operations are: the inner product of the residual vector and the preprocessed residual vector, and the inner product of the matrix-vector product result vector and the preprocessed residual vector.
6. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 4, characterized in that: The global reduction operation and calculation operation are specifically: The global reduction operation uses a non-blocking global reduction function to aggregate the local results calculated by each process to obtain a global result, including a global result after the two inner products are aggregated separately; The computing operation is based on a row preprocessing operation, a boundary update operation, and performs a matrix-vector multiplication.
7. A preprocessing pipeline conjugate gradient algorithm optimization method for reducing global communication as claimed in claim 6, characterized in that: The preprocessing operation adopts an error vector propagation block based preprocessor to accelerate the convergence of the algorithm; The boundary update operation updates the data exchanged between the current process and the adjacent process after calculating a new matrix-vector product result vector through matrix-vector multiplication.
8. A preprocessing pipeline conjugate gradient algorithm optimization system with reduced global communication, characterized in that: Used to solve sea level height in parallel ocean programs, including equation building modules, parallel solving modules and solution combining modules: The equation building module is configured to: obtain the time series of the seabed depth and sea level height of the location to be predicted, build the linear equation to be solved, and divide it into multiple sub-equations, which are assigned to different processes on the computing nodes for parallel solving; The parallel solution module is configured as follows: using the improved preprocessing pipeline conjugate gradient algorithm, in the iterative solution process of each sub-equation, the inner product merging is used to perform only one global reduction operation in each iteration, and the originally blocked global reduction function is replaced with a non-blocking global reduction function, so that the global reduction operation is overlapped with the matrix-vector multiplication, the boundary update operation and the preprocessing operation; The solution combination module is configured to combine the solutions of all sub-equations to obtain the sea level height of the position to be predicted.
9. An electronic device, comprising: a memory for non-transitory storage of computer readable instructions; a processor for executing the computer-readable instructions; Wherein, when the computer-readable instructions are executed by the processor, the method according to any one of claims 1 to 7 is executed.
10. A storage medium, characterized in that: The computer-readable instructions are non-transitorily stored, wherein when the computer-readable instructions are executed by a computer, the method of any one of claims 1 to 7 is performed.