A method and system for multi-level parallel acceleration of gridding of ocean observation data

Abstract

The application discloses a kind of marine observation data gridding multi-level parallel acceleration calculation method and system.The application is based on the coarse-grained parallel implementation of climate model between MPI parallel architecture and parallel computing acceleration of ocean depth data, mainly includes the IO parallel reading of data, the calculation and distribution of climate model data, the synchronization of gridding data etc.;Based on the fine-grained parallel computing of OpenMP parallel architecture, the parallel computing of algorithm level between CPU core in node is realized, mainly includes the data gridding of each grid point in grid;For the gridding calculation of data in grid point, the distribution of data structure in memory is adjusted, and the specific calculation operation is accelerated using SIMD.The application can flexibly configure parallel scale according to specific hardware resources, realize the maximization of computing power utilization, and can provide computing power support for the construction of long-time sequence complete and reliable ocean grid data.

Description

Technical Field

[0001] This invention belongs to the interdisciplinary field of parallel computing, marine observation and data processing, and specifically relates to a multi-level parallel accelerated computing method and system for gridded marine observation data using Message Passing Interface (MPI), Open Multi-Processing (OpenMP), and Single Instruction Multiple Data (SIMD) instruction sets. Background Technology

[0002] Long-term, complete, and reliable ocean temperature, salinity, and gridded data are fundamental to climate change and marine science research. However, the quality of ocean observations varies considerably. Limited by observation methods and costs, sparseness and uneven coverage are key characteristics of in-situ oceanographic data. These features significantly restrict the development of physical oceanography theory and technology, representing a major challenge in current marine science research.

[0003] To reconstruct gridded data with complete spatiotemporal coverage using sparse data and address problems in mainstream international statistical or dynamic methods, the research team previously proposed an ensemble optimal interpolation method. This new method outperforms current international gridded optimization methods in many aspects, and has therefore gained widespread application and recognition in the past five years, becoming a core technology with significant international influence.

[0004] However, this method is based on serial computing. Taking an Intel i5-10650KF processor with 48GB RAM as an example, calculating the gridded global sea-level heat content data at a resolution of 1°×1° for January 1960 takes 3185.7 seconds, and constructing a complete gridded temperature dataset for more than 1000 months takes 34.7 days. When the resolution is adjusted to 0.5 degrees, the computation time increases exponentially, with the computation time for a single month's data increasing to 4 hours, and the computation time for the total data taking more than 4 months. For global observation data at a resolution of 1 degree, the spatial range represented by each grid point is approximately 78,700 square kilometers. It is obviously insufficient to represent such a large-scale physical change using only a single physical data point. Therefore, theoretically, the higher the resolution when gridding the data and the smaller the grid division of each point, the more accurate the physical range represented by the calculated data will be. At this point, it is clear that this time-consuming and inefficient gridded computing method based on traditional serial computing has seriously restricted the further promotion of this method, and has even become a bottleneck for constructing large-scale, high-quality, long-term, complete and reliable ocean observation gridded datasets. Summary of the Invention

[0005] The purpose of this invention is to provide a method and system for parallel acceleration of the gridding calculation process of ocean observation data using a multi-level parallel model, so as to make full use of the computing power of a computing cluster with multiple computing nodes and multiple core processors, and solve the problem of low computational efficiency of traditional serial computing programs when gridding ocean data.

[0006] This invention is built upon an MPI+OpenMP hybrid architecture, achieving improved computational efficiency through multi-level hybrid parallel acceleration of the data gridding process. Specifically, it includes coarse-grained parallel computing based on the MPI parallel architecture to accelerate parallel computation between climate model and ocean depth data, mainly encompassing parallel data I / O reading, climate model data computation and distribution, and gridded data synchronization; fine-grained parallel computing based on the OpenMP parallel architecture to achieve algorithm-level parallel computation between CPU cores within nodes, mainly involving data gridding at each grid point; and for gridding computation of data within a grid point, SIMD is used to accelerate specific computational operations by adjusting the data structure and its distribution in memory. This invention allows for flexible configuration of the parallel scale based on specific hardware resources, maximizing the utilization of computing power and providing computational support for the construction of complete and reliable long-term ocean gridded data. The invention specifically includes the following steps:

[0007] S1: Parallel Data Reading: The gridding method reconstructs a gridded dataset of ocean observations with complete spatial coverage using ocean observation data, ocean observation error data, Coupled Model Intercomparison Project (CMIP) model data, global ocean-land boundary topographic data, and global ocean average observation data. For large amounts of input data, parallel I / O is used to improve reading efficiency.

[0008] S1-1: Reading observation data O: Parallel I / O is implemented to read observation data through OpenMP multithreading; as the resolution increases, the amount of observation data grows exponentially. The reading speed of parallel I / O is improved by adjusting the number of threads.

[0009] S1-2: Read global ocean average observation data M;

[0010] S1-3: Reading historical CMIP model data. C: Climate model data contains various models and types of data from the past 100 years. The total data volume is large, but the storage format is a multi-grid matrix structured storage. MPI+OpenMP is used to achieve multi-process, multi-threaded hybrid parallel IO mode to realize parallel accelerated reading of historical CMIP model data from multiple climate models and multiple depths.

[0011] S1-4: Read topographic data L of global ocean and land boundaries;

[0012] S1-5: Read the average state observation error data D.

[0013] S2: Data Preparation for Coarse-grained Parallel Computation Based on MPI Multi-Process: Coarse-grained parallel computation is implemented through MPI-based multi-process parallel processing. The basic unit for gridded computation of ocean observation data is a grid matrix with a resolution of Longitude° × Latitude°; grid matrices at different depths are extended downwards to obtain underwater three-dimensional information; through the calculation of gridded underwater three-dimensional gridded data under different climatic conditions, rapid reconstruction of three-dimensional data from multiple sources in the global ocean area is achieved. The observation data grid matrix with specified climatic conditions and depths is used as the smallest basic unit for gridding. Coarse-grained observation data grid matrix parallel gridded computation at different climatic conditions and depths is implemented through MPI multi-process processing. The specific steps are as follows:

[0014] S2-1: Obtain all the observation data grid matrices for climate state and depth, denoted as total; assume the current observation data grid matrix index is index;

[0015] S2-2: If the index does not exceed the total, proceed to S2-3 to perform the gridding process of the observation data grid matrix; if the index exceeds the total, proceed to S6.

[0016] S2-3: Request for computing resources on the node where the child process resides;

[0017] S2-4: The main process calculates the background perturbation field data A, spatial coordinate matrix P, observation field Q, analysis field U, and perturbation field V based on the current observation data grid matrix. The calculation method is as follows: For grid point (x,y), if |O (x,y) |<6 and L (x,y) If the value is ≥1, then the grid point needs to be processed by data gridding calculation, and the coordinate matrix of the grid point to be calculated is recorded:

[0018]

[0019] Where P i Given the number of rows P, calculate the background disturbance field data A and record the average data of the analysis field:

[0020]

[0021] U i Let k be the number of rows in U, k = x * d4 + y. The total number of grid cells in the longitude direction is the same for each latitude, which is d4. Record the ensemble observation field data:

[0022]

[0023] Q i Let Q be the row number, and if the observation error data D is... (x,y) >0, then

[0024] V (x,y) =D (x,y) ,

[0025] otherwise,

[0026]

[0027] Where m is the number of climatological models in the CMIP model, repmat(X,r1,...,rN) specifies a scalar list r1,...,rN to describe how the copies of X are arranged in each dimension; mean(X) returns a row vector containing the mean of each column (i.e., the ensemble mean field of the background field); V i To analyze the number of rows in the perturbation field V, M i Let x be the row number of the global ocean average observation matrix M, and y be the grid index in the longitude and latitude directions, respectively.

[0028] S2-5: If the spatial coordinate matrix P is not empty, proceed to S3; otherwise, proceed to S5-10.

[0029] S3: The main process distributes computing tasks to each computing node: This is a key step based on MPI coarse-grained parallel acceleration.

[0030] S3-1: Confirm the distribution strategy for coarse-grained tasks:

[0031] There are three task distribution strategies, and one strategy can be selected based on the specific situation. Assume there are a total of TCnt computing tasks, numbered from 0 to TCnt-1, and NCnt hardware computing resources, numbered from 0 to NCnt-1. The group ID (GID) after grouping tasks is determined by the number of the hardware computing resource.

[0032] (1) Static task distribution strategy: The static task distribution strategy is to distribute all computing tasks evenly across computing nodes. The group number of each task is calculated as GID = TCnt%NCnt. The advantage is that it is simple to implement and appears to be an even distribution. The disadvantage is that it ignores the fact that in actual operation, the computing time of different computing tasks is different, which will result in wasted time for computing resources to wait.

[0033] (2) Dynamic task distribution strategy: The dynamic task distribution strategy dynamically schedules each computing task and distributes new computing tasks to idle nodes that have completed their computing tasks. The advantage is that it makes better use of hardware computing resources, but the disadvantage is that it is more complex to implement and will incur additional task scheduling overhead.

[0034] (3) Hybrid Distribution Strategy: The hybrid distribution strategy is a two-stage distribution. First, the computing tasks are grouped. Second, the grouped computing tasks are distributed according to the hardware computing resources using a static distribution strategy. Third, the computing tasks within each group are distributed a second time using a dynamic distribution strategy. This strategy combines the advantages of both static and dynamic distribution strategies and is suitable for situations where the number of computing tasks far exceeds the hardware computing resources, but its implementation is complex.

[0035] S3-2: The main process (PID == 0) distributes the observation data field, perturbation data field, spatial coordinate matrix P, observation field Q, analysis field U, and perturbation field V at a specified depth to the child process (PID != 0).

[0036] S4: Data preparation for fine-grained parallel computation based on OpenMP multithreading:

[0037] S4-1: Calculate the corresponding observation perturbation matrix Z:

[0038] Z0 = Q T -mean(V T ) T

[0039] Z1 = U T

[0040] Z (:,0) =Z0

[0041] Z (:,1) =Z1

[0042] S4-2: Calculate the scaled gain vector matrix H, the scaled set observation anomaly matrix N, and the intermediate state matrix H′, and set the localization radius:

[0043]

[0044] N = R × Z0

[0045] H = R × (O - repmat (mean (V)) T ) T ,1,m))

[0046] Where m is the number of climatological models in the CMIP model. Since the row dimension of the scaled gain vector matrix H and the scaled ensemble observation anomaly matrix N is the same, both being the number of all grid points in the observation data grid matrix that need to be calculated for data gridding, and N is a column vector with 1 column, N is extended into H according to the number of columns to form H′, whose dimension is row(H)×(col(H)+1); the localization radius is determined based on factors such as the typical spatial scale of the ocean, the type of ocean observation data processed, and the actual annual rate scale.

[0047] S4-3: Confirm the fine-grained parallel task distribution strategy.

[0048] S5: Gridding of data points in the current observation data grid matrix: Let the current grid point index be jndex;

[0049] S5-1: If jndex does not exceed the total number of grid points in the current observation data grid matrix, then proceed to S5-2, the gridding process of the observation data grid matrix; otherwise, index = index + 1, then proceed to S2-2.

[0050] S5-2: Calculate the distance matrix F:

[0051] If nx≤1,

[0052] F = abs(P - jndex * I)

[0053] otherwise

[0054]

[0055] Where I is the identity matrix constructed with the number of rows of matrix P as its dimension, nx is the maximum number of grid points in the longitude direction, jndex is the grid point index, X and Y represent the distance matrices in different directions, and X = jndex % nx * IP. (:,0) Y = (jndex / nx) * IP (:,1) ;P (:,0) P (:,1) These are matrices formed by all elements in the first and second columns of matrix P, respectively.

[0056] S5-3: Calculate the localized radius matrix J:

[0057] J = exp(-0.5 * F′ * F′)

[0058] Where F′=F / radius.

[0059] S5-4: Find the local observation matrix K:

[0060] Calculate the elements greater than 0 in J. The linear indices of these elements form a vector, which is the local observation matrix K:

[0061] K = find(J>0)

[0062] S5-5: Using the values ​​of the elements in K as row indices, extract the elements of all columns in all indexed rows of J, and update the localized radius matrix J:

[0063] J = J (K,:)

[0064] S5-6: Update the calculated background disturbance field data A:

[0065] Using SIMD instructions to accelerate computation through multi-data parallelism of H′ on S′:

[0066] S′=H′ (K,:) .*repmat(J,1,m+1)

[0067] Where H′ (K,:) This indicates that the element values ​​in K are used as row indices, and the matrix formed by the row and all columns specified by the row index in H′ is obtained. .* indicates that the corresponding elements of the matrix are multiplied one by one.

[0068] Take all the elements from columns 1 to m of S′ and assign them to S as a whole, that is, S = S′. (:,1:m)

[0069] G = (I + S) T ×S) -1

[0070] Elements A in row jndex of A (jndex,:) =A (jndex,:) ×sqrt(G); As jndex increases, the complete A can be obtained.

[0071] S5-7: Calculate the background mean data field B:

[0072] Elements of all columns in row jndex of B (jndex,:) =C (jndex,:) ×(G×S T )×S′ (:,m+1) S′ (:,m+1) This represents the elements of all rows in the (m+1)th column of S′; as jndex increases, the complete B can be obtained.

[0073] S5-8: Integrate the average data field X and the perturbation data field Y at a specified depth:

[0074] X (lat,lon) =B (jndex)

[0075] Y(index,jndex) =A (:,:)

[0076] Among them, A (:,:) This represents all elements in all rows and columns of A, where `index` is the grid index. For longitude grid indexing, This is a latitudinal grid index.

[0077] S5-9: The gridding of the grid data with index jndex is completed. The corresponding observation data grid matrix index, average data field X, and perturbation data field Y are sent to the main process.

[0078] S5-10: jndex = jndex + 1, enter S5-1.

[0079] S6: Gridded Data Aggregation and Organization

[0080] S6-1: The main process receives average data field and perturbation data field data of a specified depth sent by the child process;

[0081] S6-2: The main process aggregates and organizes the average data field and the perturbation data field based on the observation data grid matrix index to form a three-dimensional ocean observation gridded set of data with complete spatial coverage;

[0082] S7: Outputs a gridded ensemble of three-dimensional ocean observation data with complete spatial coverage.

[0083] The multi-level parallel acceleration method and system based on MPI+OpenMP proposed in this invention can make full use of the hardware resources of multiple computing nodes and multiple cores of the computing cluster to achieve multi-level hybrid parallel acceleration of the gridding process of ocean observation data, thereby improving the computing efficiency and providing computing power support for the construction of long-term complete and reliable ocean gridded data. Attached Figure Description

[0084] Figure 1 This is a schematic diagram of a multi-level parallel accelerated computing method.

[0085] Figure 2 This is a schematic diagram of coarse-grained process parallel task partitioning.

[0086] Figure 3 A schematic diagram illustrating the partitioning of fine-grained thread-parallel tasks.

[0087] Figure 4 This diagram illustrates the SIMD-based parallel computation of the observed anomaly matrix and gain vector.

[0088] Figure 5 This is a diagram illustrating a static task distribution strategy.

[0089] Figure 6 This is a diagram illustrating a dynamic task distribution strategy. Detailed Implementation

[0090] The present invention will now be described in further detail with reference to the accompanying drawings. The examples given are only for explaining the present invention and are not intended to limit the scope of the present invention.

[0091] To make the objectives, content, and advantages of the present invention clearer, the specific implementation of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments, taking global ocean heat content observation data with a resolution of 1°×1° from January 1970 as an example. The specific steps are as follows:

[0092] S1: Parallel data reading:

[0093] S1-1: Reading ocean observation data O: This data is a three-dimensional matrix with dimensions of 41×180×360. The observation data is read in parallel IO mode through OpenMP multi-threading to improve the reading speed of observation data.

[0094] S1-2: Read global ocean average observation data M;

[0095] S1-3: Reading historical CMIP mode operation data. C: This data is a four-dimensional matrix with dimensions of 40×41×180×360. It uses MPI+OpenMP to achieve multi-process and multi-threaded hybrid parallel IO, enabling parallel accelerated reading of historical CMIP mode operation grid data of multiple climate modes and multiple depths.

[0096] S1-4: Read topographic data L of global ocean and land boundaries;

[0097] S1-5: Read the ocean mean state observation error data D.

[0098] S2: Data Preparation for Coarse-Grained Parallel Computation Based on MPI Multi-Process: Coarse-grained parallel computing is achieved through MPI-based multi-process parallel computing. The basic unit for gridded calculation of ocean observation data is a 180×360 grid matrix; grid matrices at different depths are extended downwards to obtain underwater three-dimensional information; through the calculation of gridded underwater three-dimensional gridded data under different climatic conditions, rapid reconstruction of three-dimensional data from multiple sources in the global ocean area is achieved. Parallel gridded calculation of coarse-grained observation data grid matrices based on climatic conditions and depth is implemented through MPI multi-process. The specific steps are as follows:

[0099] S2-1: Obtain all the observation data grid matrices for climate state and depth, totaling 1640; assume the current observation data grid matrix index is index;

[0100] S2-2: If the index does not exceed the total, proceed to S2-3 to perform the gridding process of the observation data grid matrix; if the index exceeds the total, proceed to S6.

[0101] S2-3: Resource request from the node where the child process resides;

[0102] S2-4: The main process calculates the background disturbance field data A, spatial coordinate matrix P, observation field Q, analysis field U, and disturbance field V based on the current observation data grid matrix;

[0103] S2-5: If the spatial coordinate matrix P is not empty, proceed to S3; otherwise, proceed to S5-10.

[0104] S3: The main process distributes computing tasks to each computing node:

[0105] S3-1: Confirm the distribution strategy for coarse-grained tasks: In this specific embodiment, a static task distribution strategy is adopted to group and distribute coarse-grained computing tasks in parallel;

[0106] S3-2: The main process (PID == 0) distributes the observation data field, perturbation data field, spatial coordinate matrix P, observation field Q, analysis field U, and perturbation field V at a specified depth to the child process (PID != 0).

[0107] S4: Data preparation for fine-grained parallel computation based on OpenMP multithreading:

[0108] S4-1: Calculate the corresponding observation perturbation matrix Z:

[0109] Z0 = Q T -mean(V T ) T

[0110] Z1 = U T

[0111] Z (:,0) =Z0

[0112] Z (:,1) =Z1

[0113] Where, when index = 0, Z is Z 11023×2 ;

[0114] S4-2: Calculate the scaled gain vector matrix H, the scaled set observation anomaly matrix N, and the intermediate state matrix H′, and set the localization radius:

[0115]

[0116] N = R × Z0

[0117] H = R × (O - repmat (mean (V)) T ) T ,1,m))

[0118] Where m is the number of climatological models in the CMIP model; since the row dimensions of the scaled gain vector matrix H and the scaled ensemble observation anomaly matrix N are the same, which are the number of all grid points in the observation data grid matrix that need to be calculated by data gridding, and N is a column vector with 1 column, N is extended into H according to the number of columns to form an intermediate state matrix H′ with a dimension of 11203×41;

[0119] S4-3: Confirm the fine-grained parallel task distribution strategy: In this specific embodiment, a static task distribution strategy is used to group fine-grained computing tasks.

[0120] S5: Gridding of data points in the current observation data grid matrix: Let the current grid point index be jndex;

[0121] S5-1: If jndex does not exceed the total number of grid points in the current observation data grid matrix, then proceed to S5-2, the gridding process of the observation data grid matrix; otherwise, index = index + 1, then proceed to S2-2.

[0122] S5-2: Calculate the distance matrix F:

[0123] If nx≤1,

[0124] F = abs(P - jndex * I)

[0125] otherwise

[0126]

[0127] Where nx is the maximum number of grid points in the longitude direction, jndex is the grid point index, and X = jndex % nx * IP (:,0) Y = (jndex / nx) * IP (:,1) ;P (:,0) P (:,1) These are matrices formed by all elements in the first and second columns of matrix P, respectively.

[0128] S5-3: Calculate the localized radius matrix J:

[0129] J = exp(-0.5 * F′ * F′)

[0130] Where F′=F / radius.

[0131] S5-4: Find the local observation matrix K:

[0132] Calculate the elements greater than 0 in J. The linear indices of these elements form a vector, which is the local observation matrix K.

[0133] K = find(J>0)

[0134] S5-5: Using the values ​​of elements in K as row indices, extract the elements of all columns in all specified rows of J, and update the localized radius matrix J:

[0135] J = J (K,:)

[0136] S5-6: Update the calculated background disturbance field data A:

[0137] Using SIMD instructions to accelerate computation through multi-data parallelism of H′ on S′:

[0138] S′=H′ (K,:) .*repmat(J,1,m+1)

[0139] S=S′ (:,1:m)

[0140] G = (I + S) T ×S) -1

[0141] A (jndex,:) =A (jndex,:) ×sqrt(G)

[0142] S5-7: Calculate the background mean data field B:

[0143] B (jndex,:) =C (jndex,:) ×(G×S T )×S′ (:,m+1)

[0144] S5-8: Integrate the average data field X and the perturbation data field Y at a specified depth:

[0145] X (lat,lon) =B (jndex)

[0146] Y (index,jndex) =A (:,:)

[0147] Where index is the grid index, lat = jndex / 360 is the longitude grid index, and lon = jndex%360 is the latitude grid index.

[0148] S5-9: After gridding of the jndex grid is completed, send the corresponding observation data grid matrix index, average data field X, and perturbation data field Y to the main process;

[0149] S5-10: jndex = jndex + 1, enter S5-1.

[0150] S6: Gridded Data Aggregation and Organization

[0151] S6-1: The main process receives average data field and perturbation data field data of a specified depth sent by the child process;

[0152] S6-2: The main process aggregates and organizes the average data field and the perturbation data field based on the grid matrix index of the observation data to form a gridded set of three-dimensional ocean observation data with complete spatial coverage.

[0153] S7: Outputs a gridded ensemble of three-dimensional ocean observation data with complete spatial coverage.

[0154] Although specific embodiments of the invention have been disclosed for illustrative purposes to aid in understanding and implementing the invention, those skilled in the art will understand that various substitutions, variations, and modifications are possible without departing from the spirit and scope of the invention and the appended claims. Therefore, the invention should not be limited to the content disclosed in the preferred embodiments, and the scope of protection claimed by the invention is defined by the claims.

Claims

1. A multi-level parallel accelerated computation method for gridded ocean observation data, comprising the following steps: S1: Acquire ocean observation data O, ocean observation error data D, CMIP model data C, global ocean-land boundary topographic data L, and global ocean average observation data M for ocean observation data gridding assimilation. S2: Gridding of all observation data grid matrices based on MPI multi-process: Obtain all climatological and depth observation data grid matrices, the total number of which is denoted as total; Assume the index of the current observation data grid matrix is ​​index; S3: Check if the index exceeds the total: S3-1: If index does not exceed total, confirm the distribution strategy of coarse-grained tasks, and then the main process distributes the observation data field, perturbation data field, spatial coordinate matrix P, observation field Q, analysis field U and perturbation field V of the specified depth to the child process. The single observation data grid matrix with index index is gridded, and the average data field and perturbation data field data calculated by the child process are sent to the main process. index = index + 1, and proceed to step S3. S3-2: If the index exceeds the total, the gridded data is aggregated and organized: The main process receives the average data field and perturbation data field data at each specified depth sent by the child process; then, based on the observation data grid matrix index, the average data field and perturbation data field are aggregated and organized to form a three-dimensional ocean observation gridded set of data with complete spatial coverage.

2. The method according to claim 1, characterized in that, The method for assimilating all grid points of a single observation data grid matrix based on OpenMP multithreading is as follows: S4-1: Calculate the observed perturbation matrix Z: Z0 = Q T -mean(V T ) T Z1 = U T ; S4-2: Calculate the scaled gain vector matrix H, the scaled set of observation anomaly matrix N, and the intermediate state matrix H. ′ And set the localization radius; where, N=R×Z0, H=R×(O-repmat(mean(V T ) T ,1,m); m is the number of climatological models in the CMIP model; the gain vector matrix H and the ensemble observation anomaly matrix N have the same row dimension, which is the number of all grid points in the observation data grid matrix that need to be calculated by data gridding. N is a column vector with 1 column. N is extended into H according to the number of columns to form the intermediate state matrix H. ′ The localization radius (radius) is determined based on the typical spatial scale of the ocean, the types of ocean observation data processed, and the actual annual rate scale. S4-3: Group and distribute fine-grained computing tasks according to the fine-grained parallel task distribution strategy, and then perform parallel gridded computation on each grid point one by one.

3. The method according to claim 2, wherein the method for parallel gridded computation of grid points is as follows: S5-1: If the current grid point index jndex does not exceed the total number of grid points in the current observation data grid matrix, then proceed to S5-2) the gridding process of the observation data grid matrix; otherwise, index = index + 1, and return to step S3. S5-2: Calculate the distance matrix F: If nx ≤ 1, then F = abs(P - jndex * I); otherwise in, I is an identity matrix constructed with the number of rows of the spatial coordinate matrix P as its dimension, nx is the maximum number of grid points in the longitude direction, jndex is the grid point index, and X = jndex % nx * IP. (:,0) Y = (jndex / nx) * IP (:,1) ;P (:,0) P (:,1) These are matrices formed by all elements in the first and second columns of matrix P, respectively. S5-3: Calculate the localized radius matrix J = exp(-0.5*F ′ *F ′ ); where F ′ =F / radius; S5-4: Find the local observation matrix K: Combine the linear indices of the elements greater than 0 in J to form the local observation matrix K; S5-5: Using the values ​​of the elements in K as indices, extract the elements of all columns in the row specified by the index in J, and update the localized radius matrix J; S5-6: Update the calculated background perturbation field data A: Implement using SIMD instructions and H ′ Parallel accelerated computation yields S ′ =H ( ′ K,:) .*repmat(J,1,m+1); then take S ′ All elements in columns 1 to m are assigned to S as a whole, i.e., S = S ( ′ :,1:m) ; Calculate G = (I + S) T ×S) -1 The elements of A in all columns of row jndex in A (jndex,:) =A (jndex,:) ×sqrt(G); H ( ′ K,:) Indicate that the element value in K is used as the row index, and H is taken. ′ The matrix consisting of the row specified by the row index and all columns, where .* indicates that the corresponding elements of the matrix are multiplied one by one; S5-7: Calculate the background average data field B: Calculate the elements of all columns in the jndex row of B. (jndex,:) =C (jndex,:) ×(G×S T )×S ( ′ :,m+1) S ( ′ :,m+1) S represents ′ Elements in all rows of the (m+1)th column; S5-8: Integrate the average data field X and the perturbation data field Y at a specified depth: X (lat,lon) =B (jndex) Y (index,jndex) =A (:,:) Among them, A (:,:) This represents all elements in all rows and columns of A, where index is the grid index, lat is the longitude grid index, and lon is the latitude grid index. S5-9: The gridding of the grid data with index jndex is completed. The corresponding observation data grid matrix index, average data field X, and perturbation data field Y are sent to the main process. Then jndex = jndex + 1, and proceed to step S5-1.

4. The method according to claim 2, characterized in that, The fine-grained parallel task distribution strategy can be a static task distribution strategy, a dynamic task distribution strategy, or a hybrid task distribution strategy.

5. The method according to claim 1, characterized in that, The coarse-grained parallel task distribution strategy can be a static task distribution strategy, a dynamic task distribution strategy, or a hybrid task distribution strategy.

6. The method according to claim 4 or 5, characterized in that, The static task distribution strategy is as follows: all computing tasks are evenly distributed across computing nodes, and the group number of each task is GID = TCnt%NCnt; TCnt is the total number of computing tasks, NCnt is the total number of hardware computing resources, and the number of the hardware computing resources is the group number after task grouping.

7. The method according to claim 4 or 5, characterized in that, The dynamic distribution strategy is as follows: dynamically schedule each computing task and distribute new computing tasks to idle nodes that have completed their computing tasks.

8. The method according to claim 1, characterized in that, In step S1, parallel I / O is used to acquire ocean observation data O, ocean observation error data D, CMIP model data C, global ocean-land boundary topographic data L, and global ocean average observation data M for ocean observation data grid assimilation.

Images