Ldu sparse matrix computation acceleration method, system, product, medium, and device

Abstract

This invention provides an accelerated method, system, product, medium, and device for LDU sparse matrix computation, relating to the field of high-performance numerical computation technology. It is applied to the Sunway many-core processor, which includes a main core and multiple slave cores. The method includes: discretizing the simulation region of the target computation example into an unstructured mesh to generate an lduMatrix format sparse matrix; dividing the cells and surfaces of the unstructured mesh into blocks and allocating them to each slave core, generating a communication graph containing a computation mapping table and a communication mapping table; based on the communication graph, each slave core performs two-stage parallel computation in parallel, and finally writes the complete computation results back to main memory in batches. This invention effectively overcomes the computational bottleneck of unstructured meshes on the Sunway many-core processor, significantly reduces the drag on overall performance from main memory access latency, greatly leverages the computing power of the Sunway many-core processor, and can meet the stability requirements of industrial-grade CFD simulation for long-term, large-scale computation.

Description

Technical Field

[0001] This invention relates to the field of high-performance numerical computing technology, specifically to methods, systems, products, media, and equipment for accelerating LDU sparse matrix computation. Background Technology

[0002] Computational Fluid Dynamics (CFD) is a key technology that uses computers to solve fluid dynamics equations to analyze phenomena such as flow and heat transfer. With the increasing complexity of engineering problems, unstructured meshes have become mainstream due to their strong adaptability to complex geometric boundaries. OpenFOAM, as the most popular open-source CFD software package internationally, is based on the Finite Volume Method (FVM) and written in C++, possessing excellent scalability. In the OpenFOAM solution process, solving linear equations consumes a significant amount of computation time, and sparse matrix-vector multiplication (SpMV), as a core fundamental operation, directly impacts solution efficiency. OpenFOAM has designed the lduMatrix sparse matrix storage format specifically for unstructured meshes, including diagonal coefficients, upper and lower triangular coefficients, and a topological addressing array, to efficiently handle face traversal calculations in the Finite Volume Method.

[0003] Porting OpenFOAM's native lduMatrix matrix-vector multiplication algorithm directly onto the Sunway many-core processor architecture faces significant challenges. The Sunway architecture has a unique storage hierarchy, with each slave core equipped with a limited high-speed local data storage (LDM), requiring explicit management to achieve high performance. Existing technologies generate numerous fine-grained, non-contiguous DMA access requests when handling indirect addressing caused by unstructured meshes, resulting in inefficient data transfer. Furthermore, atomic operation conflicts arising from accumulating results to the same cell during parallel computation of multiple faces are difficult to implement efficiently on the Sunway many-core architecture. In addition, the failure to fully utilize the data reuse characteristics of LDM leads to repeated readings of the same data from main memory, severely limiting computational performance.

[0004] Therefore, when the existing OpenFOAM unstructured grid sparse matrix-vector multiplication (Amul function) is ported to the domestic Sunway many-core supercomputing platform, there are technical problems such as low memory access efficiency, serious parallel write conflicts, and limited on-chip storage space. Summary of the Invention

[0005] To address the aforementioned issues, this invention proposes an LDU sparse matrix computation acceleration method, system, product, medium, and device. This effectively overcomes the computational bottleneck of unstructured grids on the Sunway many-core processor, significantly reduces the drag on overall performance from main memory access latency, greatly leverages the computing power of the Sunway many-core processor, and can meet the stability requirements of industrial-grade CFD simulation for long-term, large-scale computation.

[0006] According to some embodiments, the present invention adopts the following technical solution:

[0007] An LDU sparse matrix computation acceleration method is applied to the Sunway many-core processor, which includes a master core and multiple slave cores. The slave cores are configured with a Local Data Store (LDM), comprising:

[0008] The target computational fluid dynamics software is read into the main kernel. The simulation region of the target computational fluid dynamics software is discretized into an unstructured mesh to generate a sparse matrix in lduMatrix format. The numerical storage of the matrix is ​​decoupled from the topological structure through the internal surfaces.

[0009] Based on the topology of the unstructured mesh, the cells and faces of the sparse matrix in the lduMatrix format are divided into blocks and assigned to each slave core. By traversing the computation range of each slave core through simulation, the required remote cell indexes are collected and a mapping relationship between the global cell index and the local LDM buffer index is established, thereby generating a communication graph containing a computation mapping table and a communication mapping table.

[0010] Based on the communication graph, each slave core performs a two-stage parallel computation, including: dividing a buffer in the LDM, first traversing the internal face it is responsible for in the first stage, calculating the local contribution to neighboring cells according to the computation mapping table and storing it in the first buffer, and after synchronization, entering the second stage, pulling remote data from other slave cores LDM according to the communication mapping table and accumulating it into the second buffer, then completing the diagonal contribution and the contribution to the owner cell, and finally writing the complete computation results back to main memory in batches.

[0011] According to some embodiments, the present invention adopts the following technical solution:

[0012] The LDU sparse matrix computation acceleration system includes:

[0013] The read-in module is configured to: read the target calculation example from the computational fluid dynamics software into the main kernel, discretize the simulation region of the target calculation example into an unstructured mesh, generate a sparse matrix in lduMatrix format, and decouple the numerical storage of the matrix from the topological structure through the internal surfaces.

[0014] The generation module is configured to: divide the cells and faces of the sparse matrix in the lduMatrix format according to the topology of the unstructured mesh and assign them to each slave core; traverse the computation range of each slave core by simulation; collect the required remote cell indexes and establish the mapping relationship between the global cell index and the local LDM buffer index; and then generate a communication graph containing a computation mapping table and a communication mapping table.

[0015] The calculation module is configured to perform two-stage parallel calculations in parallel by each slave core based on the communication graph. This includes: dividing the buffer in the LDM; firstly, traversing the internal surfaces each core is responsible for in the first stage; calculating the local contribution to neighboring units according to the calculation mapping table and storing it in the first buffer; after synchronization, entering the second stage; pulling remote data from other slave core LDMs according to the communication mapping table and accumulating it in the second buffer; then completing the diagonal contribution and the contribution to the owner unit; and finally writing the complete calculation results back to main memory in batches.

[0016] According to some embodiments, the present invention adopts the following technical solution:

[0017] A computer program product includes a computer program that, when executed by a processor, implements the LDU sparse matrix computation acceleration method.

[0018] According to some embodiments, the present invention adopts the following technical solution:

[0019] A non-transitory computer-readable storage medium is provided for storing computer instructions, which, when executed by a processor, implement the LDU sparse matrix computation acceleration method.

[0020] According to some embodiments, the present invention adopts the following technical solution:

[0021] An electronic device includes a processor, a memory, and a computer program; wherein the processor is connected to the memory, the computer program is stored in the memory, and when the electronic device is running, the processor executes the computer program stored in the memory to enable the electronic device to implement the LDU sparse matrix calculation acceleration method.

[0022] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0023] The overall execution efficiency of sparse matrix-vector multiplication is significantly improved: the timing statistics of sparse matrix-vector multiplication on the entire unstructured mesh are performed using the GPTL (General Purpose Timing Library) performance analysis tool. The test results show that, under the same mesh size, the optimization method proposed in this invention achieves an overall performance speedup of approximately 2.44 times compared to the native serial / insufficiently optimized parallel method. This result proves that the "full LDM residency" and "communication graph preprocessing" strategies adopted in this invention effectively overcome the computational bottleneck of unstructured meshes on the Sunway many-core processor and significantly reduce the drag on overall performance caused by main memory access latency.

[0024] The utilization rate of the core computing part is greatly improved: when comparing the GPTL timing data of the computing part alone, the computing time of the method of this invention is significantly lower than that of the original method, and the computing core efficiency is improved by about 4.78 times. This is due to the "lock-free accumulation" mechanism designed in this invention. At the same time, the asymmetric upper triangular / lower triangular processing logic enables about 50% of the computing to be completed locally without cross-core communication. This design greatly leverages the computing power of the Sunway many-core.

[0025] It exhibits excellent linear scalability and long-term operational stability: After continuous testing with the iteration steps gradually increased from 100 to 700 (step size of 100), no performance degradation or additional overhead caused by memory leaks was observed. The speedup advantage of this invention compared to the original method remains stable, which indicates that the communication and memory management mechanism of this invention has extremely high robustness and can meet the stability requirements of industrial-grade CFD simulation for long-term, large-scale computation. Attached Figure Description

[0026] The accompanying drawings, which form part of this invention, are used to provide a further understanding of the invention. The illustrative embodiments of the invention and their descriptions are used to explain the invention and do not constitute an improper limitation of the invention.

[0027] Figure 1 This is a flowchart of the method in Example 1.

[0028] Figure 2 This is an example diagram of the surface in Example 1.

[0029] Figure 3 This is a storage hierarchy diagram of the Sunway many-core processor in Example 1.

[0030] Figure 4 This is the framework flowchart in Example 1.

[0031] Figure 5 This is a schematic diagram of the LDM space division in Example 1.

[0032] Figure 6 This is the communication graph data structure diagram in Example 1.

[0033] Figure 7 This is a diagram illustrating data communication strategies of different scales in Example 1. Detailed Implementation

[0034] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0035] It should be noted that the following detailed descriptions are exemplary and intended to provide further illustration of the invention. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0036] It should be noted that the terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the scope of exemplary embodiments according to the invention. As used herein, the singular form is intended to include the plural form as well, unless the context clearly indicates otherwise. Furthermore, it should be understood that when the terms "comprising" and / or "including" are used in this specification, they indicate the presence of features, steps, operations, devices, components, and / or combinations thereof.

[0037] Example 1

[0038] One embodiment of the present invention provides a method for accelerating LDU sparse matrix computation, applied to a Sunway many-core processor containing a master core and multiple slave cores, wherein the slave cores are configured with a local data storage (LDM), such as... Figure 1 As shown, it includes:

[0039] Step S1: Read the target case from the computational fluid dynamics software into the main kernel, discretize the simulation region of the target case into an unstructured mesh, generate a sparse matrix in lduMatrix format, and decouple the numerical storage of the matrix from the topology through the internal surfaces.

[0040] Step S2: Based on the topology of the unstructured mesh, the cells and surfaces of the sparse matrix in the lduMatrix format are divided into blocks and assigned to each slave core. By traversing the computation range of each slave core through simulation, the required remote cell indexes are collected and a mapping relationship between the global cell index and the local LDM buffer index is established, thereby generating a communication graph containing a computation mapping table and a communication mapping table.

[0041] Step S3: Based on the communication graph, each slave core performs a two-stage parallel computation, including: dividing a buffer in the LDM, first traversing the internal face it is responsible for in the first stage, calculating the local contribution to neighboring cells according to the computation mapping table and storing it in the first buffer, and after synchronization, entering the second stage, pulling remote data from other slave cores LDM according to the communication mapping table and accumulating it into the second buffer, then completing the diagonal contribution and the contribution to the owner cell, and finally writing the complete computation results back to main memory in batches.

[0042] The acceleration method in this embodiment effectively overcomes the computational bottleneck of unstructured grids on the Sunway many-core processor, significantly reduces the drag on overall performance from main memory access latency, and greatly leverages the computing power of the Sunway many-core processor. It can meet the stability requirements of industrial-grade CFD simulation for long-term, large-scale computation. To better illustrate the specific implementation process, the theoretical knowledge will be explained first:

[0043] I. Computational Fluid Dynamics and OpenFOAM

[0044] Computational Fluid Dynamics (CFD) is a technique that uses computers to solve fluid dynamics equations to analyze physical phenomena such as fluid flow, heat transfer, and chemical reactions. It typically involves discretizing the simulation domain into a finite number of grid cells for solution. As the complexity of engineering problems increases, the demand for computational accuracy also rises. Unstructured grids, due to their strong adaptability to complex geometric boundaries, have become the mainstream grid form in the field of CFD.

[0045] OpenFOAM (Open Field Operation and Manipulation) is one of the most popular open-source CFD software packages internationally. It is based on the Finite Volume Method (FVM), written in C++, and makes extensive use of object-oriented programming and template metaprogramming techniques, giving it excellent scalability. In the process of solving OpenFOAM, solving the linear equation system accounts for a large part of the computation time, and sparse matrix-vector multiplication (SpMV) is the most core and frequently called basic operation in linear solvers (such as PCG, GAMG, etc.).

[0046] II. The lduMatrix sparse matrix storage format in OpenFOAM

[0047] To efficiently process sparse matrices generated by unstructured meshes, OpenFOAM defines a unique storage format—lduMatrix (Lower-Diagonal-Upper Matrix). Unlike the common CSR (Compressed SparseRow) format, lduMatrix is ​​designed specifically for the face traversal characteristics of the finite volume method. It utilizes the principle that "two adjacent cells generate an internal face," decoupling the numerical storage of the matrix from the topological structure. The data structure of lduMatrix mainly consists of three parts:

[0048] (1) Diagonal coefficient (diagPtr) array: stores the diagonal coefficients of all grid cells, i.e. the influence coefficient of the cell on itself, and the length is equal to the total number of grid cells.

[0049] (2) Upper and lower triangular coefficients (upperPtr & lowerPtr) array: Stores the influence coefficients between the two cells of the internal face.

[0050] The two cells of each internal face are named the owner cell and the neighbor cell, respectively. The cell with the smaller index is the owner cell, and the cell with the larger index is the neighbor cell. The upperPtr stores the influence coefficient of the neighbor cell on the owner cell, while the lowerPtr stores the influence coefficient of the owner cell on the neighbor cell.

[0051] (3) Topology addressing array (lPtr&uPtr), also known as owner index array and neighbor index array: respectively record the global index of the "owner cell" and "neighbor cell" connected to each internal face, that is, record the connection relationship between the face and the cell.

[0052] For example Figure 2 Taking the face shown as an example, 1 and 8 are the indices of two adjacent cells, 5 is the internal face index generated by cell 1 and cell 8, the diagonal coefficient (diag) records the influence coefficients of cell 1 and cell 8 on themselves, the upper and lower triangular coefficients (upper & lower) record the influence coefficients of cell 8 on cell 1 and cell 1 on cell 8, respectively, the owner cell (Owner) records that cell 1 is the owner cell of face 5, and the neighbor cell (Neighbour) records that cell 8 is the neighbor cell of face 5.

[0053] Based on the above data structure, in the existing lduMatrix sparse matrix-vector multiplication, the physical quantities of each cell are updated using its own influence coefficient, the influence coefficient of the owner cell on its neighboring cells, and the influence coefficient of the neighboring cells on the owner cell. It should be noted that in computational fluid dynamics (CFD), the updated physical quantities include the velocity describing flow and the temperature describing heat transfer. Both sparse matrix-vector multiplication and the lduMatrix scheme are used to solve the high-level partial differential equations describing fluid flow and heat transfer. The formula for updating the physical quantities is:

[0054]

[0055] in, This is the unit currently awaiting update. Representation unit Updated physical quantity values, To express one's own contribution, Representation unit The coefficient of influence on itself Representation unit Physical quantity values ​​before the update; This indicates that all owner units have access to the current unit. The contribution, where n is the unit As the number of faces of neighboring units, For unit As the face of the neighboring unit, For noodles Owner unit, for Owner unit to neighbor unit Influence coefficient, For owner unit The physical quantity values; m represents the unit. The number of faces of the owner unit. This indicates all neighboring cells in relation to the current cell. Contribution For unit As the face of the owner unit, For noodles The neighboring unit, for Neighboring unit to owner unit Influence coefficient, For neighboring units The physical quantity value.

[0056] It is important to note that a numbering rule based on cell index order is used for the numbering of faces in the discrete mesh: that is, the faces associated with each cell are numbered sequentially in ascending order of cell number. Based on this numbering method, the related index array (lPtr) in the matrix structure constructed from face adjacency relationships exhibits good orderliness and local continuity, thereby improving the locality of data access and reducing memory access overhead. The optimization proposed in the algorithm utilizes this orderliness feature under this numbering rule to rearrange and optimize the data access and communication process.

[0057] III. Storage Characteristics of the Sunway Many-Core Processor Architecture

[0058] like Figure 3 As shown, the Sunway many-core processor consists of a main core and a processing core group composed of several slave cores (CPEs). The main core storage hierarchy includes registers, instruction cache, data cache, a shared L2 cache for data and instructions, and main memory. The slave core storage hierarchy includes registers, on-chip data storage space (configurable as a software-managed local data storage LDM or a hardware-managed L1 data cache), instruction cache, and main memory. Each slave core has a high-speed local data storage space LDM, which has a limited capacity and is a scarce resource of the processor.

[0059] LDM space offers fast memory access, but to achieve high performance, users need to explicitly manage LDM usage, considering how to use LDM during algorithm design. Slave cores can initiate DMA operations to transfer data between LDM and main memory, improving the efficiency of slave core access to main memory. Slave cores can also use RMA operations to transfer data between local LDM and other slave core LDMs. RMA operations enable efficient data exchange between slave cores, alleviating the problem of limited LDM capacity.

[0060] IV. LDU Sparse Matrix Computation Acceleration Framework on the Sunway Many-Core Processor

[0061] like Figure 4 As shown, the main process of this framework is as follows:

[0062] After obtaining the sparse matrix in lduMatrix format to be calculated, firstly, the original physical quantity values ​​of all units and the data structure of lduMatrix are stored in main memory. Then, the main core constructs a communication graph before the slave cores perform parallel computation, transforming irregular memory access and computation dependencies into a deterministic instruction sequence that the slave cores can directly execute by looking up tables. Secondly, a two-stage parallel mode based on the full resident state of the intermediate result LDM is designed and adopted in the parallel strategy of the slave cores to efficiently utilize the data patterns of the lduMatrix format and completely solve the write conflict problem. Finally, the complete computation results resident on LDM are written back to main memory in batches.

[0063] Based on the above process, a specific example is provided:

[0064] 1. LDM physical memory environment preparation:

[0065] Define a usable LDM buffer for each slave core, set its size to 90KB, manually manage the memory of this buffer, and divide the LDM buffer into sections as follows: Figure 5 The three logical regions shown are:

[0066] [0 ~ size1]: local_ApsiPtr (cannot be modified, other kernels need to retrieve data), which is the first buffer area, used to store the local calculation results generated in the first calculation stage;

[0067] [size1 ~ size2]: local_final_ApsiPtr (stores the final complete result), i.e., the second buffer area, used to store the final calculation result;

[0068] [size2 ~ end] : rma_buffer (RMA temporary buffer), which is the third buffer area, used as an RMA temporary buffer to store remote data pulled back from other kernel LDMs.

[0069] 2. Create a data structure for communication from the core network – a communication graph:

[0070] Communication graph data structure types and composition, such as Figure 6 As shown, the top-level data container is MyFormat, which contains the following fields:

[0071] (1) apsi_size (LDM buffer allocation table): an integer array of size 64, apsi_size[i] records the size of the LDM buffer required by the i-th slave core in the first computation stage, that is, the size of the first buffer area.

[0072] (2) g2l_face_idxs (computation mapping table): a pointer array of size 64. Each element in g2l_face_idxs is a pointer that points to a mapping array required by the slave core in the first computation stage. For example, g2l_face_idxs[i][j] means that the influence of the Owner Cell corresponding to the j-th face processed by slave core i in the first computation stage on the Neighbour Cell (i.e., the product of the influence coefficient of the owner cell on the neighbor cell and the physical quantity value stored in the main memory) is stored in the g2l_face_idxs[i][j] position of the slave core's local LDM.

[0073] (3) lPtr_start_idxs (local neighbor iteration index): an integer pointer. Since the face list is sorted according to the owner cell index corresponding to the face, this array records the continuous start and end positions of all faces in the face list with these cells as owner cells when calculating the cells that belong to the current core. This allows the slave core to directly read the corresponding coefficients and corresponding physical quantity values ​​to calculate the cells currently managed by the slave core without performing complex cross-core communication.

[0074] (4) segment_size indicates the total number of segments, i.e. the size of the segment array below.

[0075] (5) segments (communication map): a pointer to the Segment array, which records all the communication task sequences that the kernel needs to execute in the second computation stage.

[0076] A Segment is a data structure used to manage all RMA requests for the same logical block. It represents a list of all communication tasks to be performed by the core and contains two fields:

[0077] rmas_size: Integer, records how many independent RmaInfo requests the current segment contains.

[0078] rmas: A pointer to RmaInfo, recording the RmaInfo requests contained in the current segment. RmaInfo is the smallest unit of information used for remote memory access, and its fields include:

[0079] 1) src (source core identifier): Integer, indicating which slave core's LDM the remote data to be pulled is currently residing in.

[0080] 2) offset (remote offset): Integer, indicating the relative starting position of the remote data to be fetched in the target slave kernel LDM buffer.

[0081] 3) size (data length): Integer, representing the number of consecutive data elements in this batch. This field directly determines whether subsequent communication uses the GAS command (when size < 8) or the RMA engine (when size ≥ 8).

[0082] 4) map (local mapping table): an integer pointer that points to an array. Since the remote data fetched back logically belongs to different units, this array records the specific index position of the local result vector where the k-th element should be accumulated. This local mapping table abstracts the data fetching operation into a standard memory address offset operation.

[0083] In this context, map[k] records the physical quantity of the kth number obtained in this communication that needs to be accumulated in the map[k] unit. Each slave core is responsible for the unit from start to end, but the value of the start unit is stored in the local_final_ApsiPtr[0] position of the local LDM. Therefore, the kth number obtained should be accumulated in local_final_ApsiPtr[map[k]–start].

[0084] An RmaInfo instance represents a communication. For example, if an RmaInfo instance contains the data {src=1, offset=5, size=3, map}, this indicates that the slave core possessing this RmaInfo instance needs to perform such a communication.

[0085] Starting from local_ApsiPtr[5] of kernel 1, three elements are retrieved consecutively and remotely pulled into the local LDM. The 0th number retrieved is added to the position corresponding to map[0]-start, and the 1st number is added to the position corresponding to map[1]-start. This process is repeated to complete one communication cycle.

[0086] The MyFormat data structure allows for two main advantages. First, it enables the preprocessing of irregular, unstructured grids into regular data flow graphs suitable for multi-core parallel processing. Second, during the iterative computation of OpenFOAM, the grid topology changes multiple times. MyFormat can cache the data flow graphs processed at different levels of grid topology, avoiding an expensive reconstruction every time the grid topology changes.

[0087] When dealing with a target computational instance, constructing a communication graph is not necessary. In this embodiment, a communication graph is constructed and parallel computation from the kernel is performed based on the communication graph. This method is effective for large-scale grids, but requires separate processing for small-scale grids.

[0088] Specifically, after acquiring standard data for the target example, the mesh size is determined first because subsequent code needs to build a complex communication graph, which will have a large initialization overhead. If the mesh is small, the overhead of building the graph is even greater than that of computing directly on the main core. Therefore, for small meshes (total number of faces less than 10,000), the native serial code of OpenFOAM is executed directly.

[0089] For large-scale meshes, the cache is checked first. If the mesh topology has not been processed before, the communication graph construction process is initiated. Otherwise, the previously constructed communication graph is retrieved directly from the cache. That is, the communication graph needs to be constructed only after the mesh topology changes. Otherwise, the communication graph construction step on the main core is omitted to reduce computational overhead.

[0090] After entering the communication graph construction process, both faces and cells are divided into blocks and evenly distributed to all slave cores. In this embodiment, there are 64 slave cores. Based on the communication graph data structure type and composition designed above, a format object to be filled is created, simulating the traversal of the range managed by each slave core, including the allocated face range and cell range, specifically:

[0091] Step 1: Based on the neighboring units of the face, construct the computational mapping table g2l_face_idxs and the communication mapping table segments of MyFormat.

[0092] Calculate the range of faces currently being processed by the kernel and analyze which cells the kernel is involved in. Here, only the neighbor cells corresponding to all faces are collected, not the owner cells, because the indexes of the owner cells are ordered and correspond to different processing flows. That is, in the second step below, a set v is used to store all the collected cells. After collection, the cells are sorted in ascending order according to the indexes of the neighbor cells.

[0093] The cells here use global cell indices. Based on the sorted cell indices, a local index of the local buffer is constructed. A mapping relationship is established between the global cell index and the local buffer index. Then, based on this mapping relationship, the calculation mapping table g2l_face_idxs required for the execution phase is established. That is, the relationship between the i-th face processed by the kernel and the index of the corresponding Neighbour cell in the local buffer is established, and this relationship is stored in the MyFormat structure.

[0094] Next, the data flow is analyzed. The cells (sorted v) corresponding to all faces processed by the current core are traversed. The segment number (i.e., the slave core number) to which this cell belongs is calculated. A greedy strategy is used to find a continuous segment of faces where all corresponding cells belong to the same target segment. This allows multiple small data points managed by the same slave core to be merged into a single RMA request, improving bandwidth utilization. A map is constructed by recording the global cell index corresponding to each data point, packaged into an RmaInfo, and this request is added to the Segment of the target segment, completing the initial construction of the communication graph.

[0095] Step 2: Based on the owner unit of the face, construct the local neighbor iteration index lPtr_start_idxs of MyFormat.

[0096] Utilizing the ordered pattern of lPtr in OpenFOAM, calculate the segment index of lPtr and save it to the current MyFormat.

[0097] Since the face list is sorted according to the Owner cell index corresponding to the face, the continuous start and end positions of all faces in the face list with the cell currently managed by the slave core as the Owner cell are calculated to form the local neighbor iteration index lPtr_start_idxs, which is saved to the current MyFormat. This allows the slave core to directly read the values ​​of the corresponding coefficients and source vectors to calculate the cells currently managed by the slave core without the need for complex cross-core communication.

[0098] At this point, the unstructured and irregular sparse matrix multiplication has been forcibly decomposed and reassembled into two sets of mapping tables:

[0099] The computation mapping table (g2l_face_idxs) tells each slave core how to compute its local tasks in the first computation phase, how to save the results, and other information.

[0100] Communication mapping table (Segment): tells each core which other cores it needs to retrieve data from after the first stage of computation, and where the retrieved data should be stored, etc.

[0101] Through this "preprocess once, use multiple times" strategy, the Sunway supercomputer can perform high-efficiency computations during the actual sparse matrix-vector multiplication stage, eliminating the overhead of conditional judgments and dynamic searches and resolving the serious write conflict problem.

[0102] 3. Two-stage calculation from the core

[0103] By designing and adopting a two-stage parallel mode based on the intermediate result LDM full resident mode from the core parallel strategy, the data patterns of the lduMatrix format are efficiently utilized, and the write conflict problem is completely solved. The specific steps are as follows:

[0104] 1) First calculation stage

[0105] By embedding assembly code to read the slave core ID register, the computational mapping table g2l_face_idxs[slave core ID][face ID] and buffer size of the current slave core are retrieved from the pre-processed structure. A first buffer region local_ApsiPtr of the corresponding size is allocated and initialized in the initially declared slave core LDM local buffer. The mapping relationships included in the current slave core are sequentially retrieved from the computational mapping table g2l_face_idxs. Each face handled by the current slave core is traversed, and the influence of the owner unit of this face on its neighboring units is calculated and accumulated into the local LDM. Finally, the local contribution of all faces in the current slave core to their specific neighboring units is obtained. This process does not require any locks, and the writes are made to the first buffer region in the slave core LDM, making it extremely fast.

[0106] 2) Synchronization

[0107] In the second computation phase, each slave core needs to pull the corresponding local contribution from the LDM of other slave cores. Therefore, all slave cores need to be stopped and wait until all slave cores have completed the task of the first computation phase.

[0108] 3) Second calculation stage

[0109] In the second computation phase, each slave core's computation revolves around managing the units it owns, specifically:

[0110] First, it is necessary to capture the local data scattered in other slave cores LDM and store it in the RMA temporary buffer; specifically, for each slave core, the pre-processed communication mapping table is retrieved, and data with an offset of offset is retrieved from the slave core with src number according to each RmaInfo, for a total of size.

[0111] The data retrieval method described here addresses the significant differences in communication granularity in sparse matrix multiplication, such as... Figure 7As shown, this embodiment designs two adaptive communication strategies based on the size, making full use of the characteristics of the Sunway architecture to significantly reduce the fine-grained communication latency commonly found in unstructured meshes.

[0112] The first approach, for extremely small data (size less than 8): Based on the characteristics of the Sunway architecture, a global address space pointer is constructed. The starting address of the remote data is obtained by adding offset to `local_ApsiPtr`. Bit 45 of this address is set to 1, and bits 20 to 25 are filled with the target slave core's ID, i.e., `src`. This pointer is used for a total of `size` consecutive memory accesses, reading one number at a time. The hardware intercepts this address, converts it into a data packet for the on-chip network, and retrieves the data from the corresponding slave core. The advantage is extremely fast startup with almost no software overhead; the disadvantage is low bandwidth, as a network trip is required for each data read. Therefore, for small-scale data, the overhead of fetching data via the global address space (GAS) is much smaller than that of RMA.

[0113] The second approach, for large-scale data (size greater than or equal to 8), firstly, the slave node where the remote data is located is determined based on the source slave core number src, and the starting address of the data in the slave core is located by combining the offset offset. Then, through an rma_get operation, the entire continuous data block of length size is moved to the local buffer rma_buffer. Then, the data block is traversed, and according to the mapping relationship map, each data element in the buffer is accumulated to the corresponding position in the local result array local_final_ApsiPtr.

[0114] Then, calculate the diagonal contribution and the impact on the owner unit:

[0115] Each kernel calculates the diagonal contribution of the cells it manages and adds it to local_final_ApsiPtr. Since the face list in OpenFOAM is sorted by the owner cell index, it directly traverses all faces whose owner cells are all managed cells, calculates the impact on the owner cells, and adds it to the final result local_final_ApsiPtr.

[0116] Finally, the local_final_ApsiPtr in the LDM buffer contains the complete calculation results. Each slave core writes it back to the corresponding location in main memory in batches via DMA, completing the efficient parallel computation of sparse matrix-vector multiplication in lduMatrix format.

[0117] Example 2

[0118] One embodiment of the present invention provides an LDU sparse matrix computation acceleration system, comprising:

[0119] The read-in module is configured to: read the target calculation example from the computational fluid dynamics software into the main kernel, discretize the simulation region of the target calculation example into an unstructured mesh, generate a sparse matrix in lduMatrix format, and decouple the numerical storage of the matrix from the topological structure through the internal surfaces.

[0120] The generation module is configured to: divide the cells and faces of the sparse matrix in the lduMatrix format according to the topology of the unstructured mesh and assign them to each slave core; traverse the computation range of each slave core by simulation; collect the required remote cell indexes and establish the mapping relationship between the global cell index and the local LDM buffer index; and then generate a communication graph containing a computation mapping table and a communication mapping table.

[0121] The calculation module is configured to perform two-stage parallel calculations in parallel by each slave core based on the communication graph. This includes: dividing the buffer in the LDM; firstly, traversing the internal surfaces each core is responsible for in the first stage; calculating the local contribution to neighboring units according to the calculation mapping table and storing it in the first buffer; after synchronization, entering the second stage; pulling remote data from other slave core LDMs according to the communication mapping table and accumulating it in the second buffer; then completing the diagonal contribution and the contribution to the owner unit; and finally writing the complete calculation results back to main memory in batches.

[0122] Example 3

[0123] One embodiment of the present invention provides a computer program product, including a computer program that, when executed by a processor, implements the LDU sparse matrix calculation acceleration method.

[0124] Example 4

[0125] In one embodiment of the present invention, a non-transitory computer-readable storage medium is provided for storing computer instructions, which, when executed by a processor, implement the LDU sparse matrix calculation acceleration method.

[0126] Example 5

[0127] One embodiment of the present invention provides an electronic device, including: a processor, a memory, and a computer program; wherein the processor is connected to the memory, the computer program is stored in the memory, and when the electronic device is running, the processor executes the computer program stored in the memory to enable the electronic device to implement the LDU sparse matrix calculation acceleration method.

[0128] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0129] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0130] While the specific embodiments of the present invention have been described above in conjunction with the accompanying drawings, this is not intended to limit the scope of protection of the present invention. Those skilled in the art should understand that various modifications or variations that can be made by those skilled in the art without creative effort based on the technical solutions of the present invention are still within the scope of protection of the present invention.

Claims

1. A method for accelerating LDU sparse matrix computation, characterized in that, Applied to the Sunway many-core processor, which includes a master core and multiple slave cores, wherein the slave cores are configured with a local data storage (LDM), including: The target computational fluid dynamics software is read into the main kernel. The simulation region of the target computational fluid dynamics software is discretized into an unstructured mesh to generate a sparse matrix in lduMatrix format. The numerical storage of the matrix is ​​decoupled from the topological structure through the internal surfaces. Based on the topology of the unstructured mesh, the cells and faces of the lduMatrix sparse matrix are divided into blocks and assigned to each slave core. By simulating the traversal of the computation range handled by each slave core, the required remote cell indices are collected, and a mapping relationship between global cell indices and local LDM buffer subscripts is established. This generates a communication graph containing a computation mapping table and a communication mapping table. The computation mapping table records the mapping relationship between each face processed by each slave core in the first computation stage and its local buffer. The communication mapping table records the complete sequence of communication tasks that all slave cores need to execute in the second computation stage. The communication mapping table contains several Segment structures, and each Segment structure contains several RmaInfo structures. The RmaInfo structure includes at least: source core identifier src, remote offset offset, data length size, and a mapping table map indicating the local accumulation position of the retrieved data. Based on the communication graph, each slave core performs a two-stage parallel computation, including: dividing a buffer in the LDM, first traversing the internal face it is responsible for in the first computation stage, calculating the local contribution to neighboring cells according to the computation mapping table and storing it in the first buffer, and after synchronization, entering the second computation stage, pulling remote data from other slave cores LDM according to the communication mapping table and accumulating it into the second buffer, then completing the diagonal contribution and the contribution to the owner cell, and finally writing the complete computation results back to main memory in batches; Specifically, dividing the buffer in the LDM involves dividing the LDM space of each slave core into three logical regions: The first buffer area is used to store the local calculation results generated in the first calculation stage; The second buffer area is used to store the final calculation results; The third buffer area serves as a temporary buffer for RMA to store remote data pulled back from other kernel LDMs; The remote data refers to the local computation results generated during the first computation phase, which are stored in the local LDM by other kernels.

2. The LDU sparse matrix computation acceleration method of claim 1, wherein, The lduMatrix format sparse matrix contains: A diagonal coefficient array used to record its own influence coefficient; The upper triangular coefficient array and the lower triangular coefficient array are used to record the influence coefficients between the two sides of the internal surface; Owner index array and neighbor index array are used to record the topological relationship between internal faces and cells.

3. The LDU sparse matrix computation acceleration method of claim 1, wherein, The communication graph is represented by a custom data structure MyFormat, which includes a computation mapping table g2l_face_idxs and a communication mapping table segments, as well as an apsi_size array for recording the size of the LDM buffer required by each slave core in the first computation stage.

4. The LDU sparse matrix computation acceleration method of claim 1, wherein, The first calculation stage includes: The current slave core's ID register is read by embedding assembly instructions. Based on the slave core ID, the corresponding calculation mapping table and buffer size information are retrieved from the communication graph. The calculation mapping table is traversed sequentially to complete the calculation and accumulate it into the first buffer. No atomic operations or locking mechanisms are required throughout the process.

5. The LDU sparse matrix computation acceleration method of claim 1, wherein, The batch write-back to main memory includes: Each slave core writes the complete calculation results in the second LDM buffer back to the corresponding global vector location in main memory in one go via DMA, completing a complete sparse matrix-vector multiplication operation.

6. An LDU sparse matrix computation acceleration system, characterized by, The method for accelerating LDU sparse matrix computation as described in any one of claims 1-5 includes: The read-in module is configured to: read the target calculation example from the computational fluid dynamics software into the main kernel, discretize the simulation region of the target calculation example into an unstructured mesh, generate a sparse matrix in lduMatrix format, and decouple the numerical storage of the matrix from the topological structure through the internal surfaces. The generation module is configured to: divide the cells and faces of the sparse matrix in the lduMatrix format according to the topology of the unstructured mesh and assign them to each slave core; traverse the computation range of each slave core by simulation; collect the required remote cell indexes and establish the mapping relationship between the global cell index and the local LDM buffer index; and then generate a communication graph containing a computation mapping table and a communication mapping table. The computation module is configured to perform two-stage parallel computations in parallel by each slave core based on the communication graph. This includes: dividing a buffer in the LDM; firstly, traversing the internal surfaces each core is responsible for in the first computation stage; calculating the local contribution to neighboring units according to the computation mapping table and storing it in the first buffer; then, entering the second computation stage after synchronization; pulling remote data from other slave cores LDM according to the communication mapping table and accumulating it in the second buffer; then completing the diagonal contribution and the contribution to the owner unit; and finally writing the complete computation results back to main memory in batches.

7. A computer program product comprising a computer program, characterized in that, When the computer program is executed by the processor, it implements the LDU sparse matrix calculation acceleration method according to any one of claims 1-5.

8. A non-transitory computer-readable storage medium, comprising: The non-transitory computer-readable storage medium is used to store computer instructions, which, when executed by a processor, implement the LDU sparse matrix computation acceleration method as described in any one of claims 1-5.

9. An electronic device, comprising: include: The device includes a processor, a memory, and a computer program; wherein the processor is connected to the memory, the computer program is stored in the memory, and when the electronic device is running, the processor executes the computer program stored in the memory to enable the electronic device to perform the LDU sparse matrix computation acceleration method as described in any one of claims 1-5.

Images