CN106709181A - Parallel programming and modular method based distributed hydrologic model calibration method - Google Patents

Abstract

The invention discloses a parallel programming and modular method based distributed hydrologic model calibration method and aims to improve calibration efficiency and simulation precision of a basin hydrologic model. The method includes: 1) on the basis of a parallel programming environment constructed based on MPICH, coupling a source program of the open-source distributed hydrologic model with a multi-objective optimization algorithm; 2) calibrating the hydrologic model aiming at a daily runoff sequence and a peak value sequence respectively, and establishing corresponding model frames; 3) adopting a modular method for combining runoff outputs of the models to form recombined simulation runoff. The parallel programming and modular method based distributed hydrologic model calibration method has advantages that by parallel multi-objective operation, calibration efficiency of the distributed hydrologic model can be effectively improved while calibration time is shortened; the modular method takes full advantages of a peak value model and a daily runoff model, and accordingly simulation precision can be improved remarkably.

Description

A kind of hydrological distribution model rating method based on multiple programming and agent approach

Technical field

The present invention relates to hydrological model technical field, a kind of distribution based on multiple programming and agent approach is related in particular to Formula hydrological model rating method.

Background technology

With urban sprawl and the aggravation of climate change, the generation of the natural calamity such as flood arid is more and more frequent, seriously National economy and people life property safety are threaten, accurate topotype is carried out to flood and drought process using Watershed Hydrologic Models Fit and forecast with unerring accuracy to alleviating the condition of a disaster important in inhibiting.

Hydro-Model Parameter Calibration Technology calibration is the important step of hydrological simulation, and the calibration of model parameter directly influences hydrological model Prediction performance.Hydrological distribution model is related to parameter numerous, and calibration is complicated, and ginseng not only inefficiency is adjusted manually, expends a large amount of Time and efforts, and be difficult to find that optimal parameter group, bring very big uncertainty to model result, significantly limit The application of hydrological distribution model.Therefore, the automatic calibration technology of parameter is appealed to, the calibration of hydrological distribution model will be significantly improved Efficiency and simulation precision.But, as the calculating process of hydrological distribution model is increasingly sophisticated fine, hydrological simulation is more and more smart While really, also cause the time of single simulation more and more long, even if using automatic calibration technology, time of consuming is still very It is long.

Practical Project shows to be usually present multiple conflicting object functions in hydrologic process, traditional based on monocular Target prioritization scheme, only only accounts for hydrologic process feature in a certain respect, it is impossible to fully reflect the different dynamic of Hydrology Scholarship and moral conduct is to cause simulation precision not enough.

The content of the invention

Be to solve the deficiencies in the prior art, it is an object of the invention to provide it is a kind of based on multiple programming and agent approach point Cloth hydrological model rating method, to solve hydrological distribution model, time-consuming during calibration, and efficiency is low, and precision is not high Problem.

To realize above-mentioned target, a kind of hydrological distribution model rating method based on multiple programming and agent approach, its spy Levy and be, comprise the following steps：

1) the hydrometeorological data that hydrological model needs are collected, meteorological data is interpolated into the resolution ratio of model needs, really Cover half type treats calibration parameter；

2) determine the threshold value of super sizing technique, chosen using super sizing technique (Peaks Over Threshold methods) and ground The peak value sequence of diurnal courses during studying carefully.

3) MPI (Message-Passing Interface) Parallel Programming Environment based on MPICH, by multiple-objection optimization Algorithm ε-NSGA II and hydrological model source program are coupled, with step 1) collect meteorological data, and actual measurement hydrographic data make It is input, assorted efficiency factor NSE/D and relative deviation Bias/D is received with diurnal courses sequence, and peak value sequence receives assorted efficiency Coefficient NSE/P and relative deviation Bias/P is used as object function, calibration hydrological model；

4) target Watershed Hydrologic Models are set up according to calibration result, linear reservoir is integrated using agent approach, form restructuring Inflow Sequence.

In above-mentioned technical proposal, step 2) it is specially：

(1) threshold value T is set, and the setting of threshold value need to ensure that the peak value of final selection meets independence requirement；

(2) it is chosen for runoff peak value more than the run-off of the threshold value.

Described independence is required：

(1) the number of days interval between two neighboring peak value must be more than the natural logrithm of drainage area (sq. mi) plus 5；

Intermediate flow between (2) two adjacent peaks must be less than 75% of smaller value in the two peak-datas.

Step 3) described in multi-objective optimization algorithm ε-NSGA II and hydrological model source program are coupled, specific side Method is as follows：

(1) relevant parameter of the algorithms of setting ε-NSGA II：Initial population number N, evolutionary generation gens, crossover probability, variation Probability etc., while the correlation for setting VIC models treats the number and scope of calibration parameter；

(2) MPI is initialized, what is included in acquisition MPI process groups enters number of passes cSize and process of this process in communication set Numbering pRank；

(3) Population Size is randomly generated by main program master and is the parent population P0 (i.e. N number of parameter group) of N, and passed through Processor active task is distributed to cSize subprogram (being separately operable hydrological model according to parameter group) by MPI_Bcast functions；

(4) according to result of calculation, calculate each individual fitness function in parent population P0 and (receive assorted efficiency factor and relative Deviation), carry out quick non-dominated ranking and crowding is calculated, produce A0 ε non-dominant to achieve solution；

(5) 25% injection mechanism (Injection Scheme) is used, A0 ε non-dominant in (4) achieves solution and formed 25% individuality in new parent population Pt, 75% by randomly generating, Population Size is Nt for remaining；

(6) progeny population Qt is produced by selection, intersection, mutation operator, Population Size is Nt；By progeny population Qt and father Merge the new population Rt=Pt ∪ Qt of composition for Pt, Population Size is 2Nt；

(7) processor active task is distributed to by cSize subprogram, and root by MPI_Bcast functions by main program master The calculating of quick non-dominated ranking and crowding is carried out according to operation result, according to each individual non-dominated ranking in colony and gathering around Squeeze distance selection and produce elite parent population Nt, while producing At ε non-dominant to achieve solution；

(8) the 25% injection mechanism according to (5), by the Population Size of archive population A t, increases random generation sample and is formed Population Size is the parent population Pt+1 of 4At.

(9) circulation step (6)-(8), evolutionary generation or other end conditions until reaching setting.

Optimum results according to the multi-objective Algorithms of ε-NSGA II are the diurnal courses parameter group different with the selection of peak value runoff, structure Build the model/D models specifically designed for daily stream flow and the model/P models for peak analogue.

The output result of model/D models and model/P models is carried out integrating as follows to obtain final calibration As a result Q_C：

Wherein m_model/DAnd m (Q)_model/P(Q) it is the distribution coefficient of model/D and model/P models, Q_model/A(t) and Q_model/PT () is respectively the linear reservoir of model/D and model/P models.Q_o,maxDiurnal courses data are periodically surveyed for rate most Big value, T is the threshold value of the super sizing technique described in claim 3, and E (T) is the expectation of all threshold values that claim 3 is chosen.Q_o It is measuring runoff.

By using above-mentioned technological means, beneficial effects of the present invention to be：

(1) parallel automatic rating method solves that hydrological distribution model efficiency during calibration is low, time-consuming asks Topic, is remarkably improved model efficiency, in same time, can obtain degree of accuracy calibration result higher.

(2) multi-objective optimization algorithm significantly reduces the uncertainty of model parameter, according to multiple target solution, ginseng Number determination process are more objective, it is to avoid the subjectivity of artificial judgement.

Further aspect of the present invention, by the description of detailed description below, is clearly stated that.

Brief description of the drawings

Fig. 1 is the schematic flow sheet of example of the present invention.

Fig. 2 is the false code of multiple programming.

Fig. 3 is the solution NSE/D and NSE/P of multi-objective Algorithm, for parameter determination.

Fig. 4 is diurnal courses single goal calibration result.

Fig. 5 is peak value runoff single goal calibration result.

Fig. 6 is the multiple target calibration result recombinated using agent approach.

Specific embodiment

The present invention will be further described with reference to the accompanying drawings and detailed description.It is protrusion advantage of the invention, adopts The variable infiltration capacity model VIC proposed with University of Washington is embodied as case.

As shown in figure 1, a kind of hydrological distribution model rating method based on multiple programming and agent approach of the invention, bag Include following steps：

The GIS data and Law of DEM Data in target basin are extracted and located using arcgis Reason, obtains the data such as basin DEM, flow direction and profile, and sizing grid takes 5 × 5.

Vegetation parameter, soil parameters and meteorological Parameter File that preparation model needs, file resolution is 5 × 5.

Vegetation pattern is divided into 11 classes using the global 1km soil cover classification of University of Maryland's development.Vegetation parameter bag Include vegetation structure impedance, the impedance of minimum stomata, the leaf area index of vegetation, zero-plane displacement, albedo, roughness and root area Distribution proportion in each sublevel soil etc..

Soil data is according to China 1:1000000 pedological maps cut, and soil parameters mainly includes base flow exponent, satisfies With hydraulic conductivity, soil moisture diffusion parameter, soil antecedent water content, net lattice control elevation, soil horizon thickness, average soil Earth temperature, soil pressure in bubbles, soil content, soil total body density, soil bulk density, critical soil moisture content ratio, Wilting point soil moisture content ratio, exposed soil surface roughness, accumulated snow rugosity, year-round average temperature, residual moisture content etc..Above-mentioned parameter Data can be extracted by arcgis directly to determine, and be directly used in model, typically without modification.

Variable infiltration capacity curve parameter, the non-linear base flow maximum factor, base flow Peak Flow Rate, base flow hair in soil file Raw maximum soil water storage capacity factor measure and thickness of soil are determined by parameter calibration.

Meteorological file includes daily rain amount total amount, highest day temperature, minimum day temperature, five variables of wind speed and relative humidity.

The peak value sequence of diurnal courses during choosing research using super sizing technique, the threshold value of super sizing technique is put down for many years according to basin Equal run-off determines that value is twice average annual flow, while ensuring that peak-data is more than 30 in the research phase.

Multi-objective genetic algorithm ε-NSGA II and VIC model source codes are coupled based on MPICH, concrete thought It is (as shown in Fig. 2 false codes)：

(1) relevant parameter of the algorithms of setting ε-NSGA II：Initial population number N, evolutionary generation gens, crossover probability, variation Probability, while the correlation for setting VIC models treats the number and scope of calibration parameter；

(2) MPI is initialized, what is included in acquisition MPI process groups enters number of passes cSize and process of this process in communication set Numbering pRank；

(3) Population Size is randomly generated by main program master and is the parent population P0 (i.e. N number of parameter group) of N, and passed through Processor active task is distributed to cSize subprogram (being separately operable hydrological model according to parameter group) by MPI_Bcast functions；

(4) according to result of calculation, calculate each individual fitness function in parent population P0 and (receive assorted efficiency factor and relative Deviation), carry out quick non-dominated ranking and crowding is calculated, produce A0 ε non-dominant to achieve solution；

(5) 25% injection mechanism (Injection Scheme) is used, A0 ε non-dominant in (4) achieves solution and formed 25% individuality in new parent population Pt, 75% by randomly generating, Population Size is Nt for remaining；

(6) progeny population Qt is produced by selection, intersection, mutation operator, Population Size is Nt；By progeny population Qt and father Merge the new population Rt=Pt ∪ Qt of composition for Pt, Population Size is 2Nt；

(7) processor active task is distributed to by cSize subprogram, and root by MPI_Bcast functions by main program master The calculating of quick non-dominated ranking and crowding is carried out according to operation result, according to each individual non-dominated ranking in colony and gathering around Squeeze distance selection and produce elite parent population Nt, while producing At ε non-dominant to achieve solution；

(8) the 25% injection mechanism according to (5), by the Population Size of archive population A t, increases random generation sample and is formed Population Size is the parent population Pt+1 of 4At.

(9) circulation step (6)-(8), evolutionary generation or other end conditions until reaching setting.

The above-mentioned meteorology of input, soil, vegetation data and measuring runoff data, according to following object function come calibration hydrology mould Type：

N and M is the quantity of diurnal courses and peak value runoff, Q in formula_obsAnd Q_simRespectively survey daily flow and simulating sun stream Amount, Q_obs,PAnd Q_sim,PRespectively survey and simulated peak.

The solution of the multi-objective optimization algorithm according to Fig. 3 is the diurnal courses parameter different with the selection of peak value runoff Group, builds the VIC/D models specifically designed for daily stream flow and the VIC/P models for peak analogue.

The output result of VIC/D models and VIC/P models is carried out as follows to integrate to obtain final calibration result Q_C：

Wherein m_VIC/DAnd m (Q)_VIC/P(Q) it is the distribution coefficient of VIC/D and VIC/P models, Q_VIC/D(t) and Q_VIC/P(t) point Not Wei VIC/D and VIC/P models output.Q_o,maxTo survey the maximum of diurnal courses data, T is the threshold value of super sizing technique, E (T) be choose all threshold values expectation.

Using above method calibration result as Figure 4-Figure 6, it can be seen that using agent approach integrate Inflow Sequence without By little with measuring runoff deviation in runoff high or low runoff plot.

The above is only implemented to example of the invention, is not intended to limit the invention, for the research of this area For person, the present invention can have various modifications and variations.It is every claim of the invention restriction in the range of, that is done is any Modification, equivalent, improvement etc., all should be within protection scope of the present invention.

Claims (6)

1. a kind of hydrological distribution model rating method based on multiple programming and agent approach, it is characterised in that including following step Suddenly：

1) the hydrometeorological data that hydrological model needs are collected, meteorological data is interpolated into the resolution ratio of model needs, determine mould Type treats calibration parameter；

2) threshold value of super sizing technique is determined, the peak value sequence of diurnal courses during choosing research using super sizing technique；

3) Parallel Programming Environment based on MPICH, coupling is carried out by multi-objective optimization algorithm ε-NSGA II and hydrological model source program Close, with step 1) collect meteorological data, and actual measurement hydrographic data as input, with the assorted efficiency factor of receiving of diurnal courses sequence NSE/D and relative deviation Bias/D, and peak value sequence receives assorted efficiency factor NSE/P and relative deviation Bias/P as target Function, calibration hydrological model；

4) target Watershed Hydrologic Models are set up according to calibration result, linear reservoir is integrated using agent approach, form the runoff of restructuring Sequence.

2. the hydrological distribution model rating method based on multiple programming and agent approach according to claim 1, its feature Be, step 2) in the super sizing technique that uses, detailed process is：

(1) threshold value T is set, and the setting of threshold value need to ensure that the peak value of final selection meets independence requirement；

(2) it is chosen for runoff peak value more than the run-off of the threshold value.

3. the hydrological distribution model rating method based on multiple programming and agent approach according to claim 2, its feature It is that described independence requirement is specifically included：

(1) the number of days interval between two neighboring peak value must be more than the natural logrithm of drainage area (sq. mi) plus 5；

Intermediate flow between (2) two adjacent peaks must be less than 75% of smaller value in the two peak-datas.

4. the hydrological distribution model rating method based on multiple programming and agent approach according to claim 1, its feature Be, step 3) described in multi-objective optimization algorithm ε-NSGA II and hydrological model source program are coupled, specific method It is as follows：

(1) relevant parameter of the algorithms of setting ε-NSGA II：Initial population number N, evolutionary generation gens, crossover probability, mutation probability, The correlation for setting VIC models simultaneously treats the number and scope of calibration parameter；

(2) MPI is initialized, the process of number of passes cSize and this process in communication set of entering included in MPI process groups is obtained and is numbered pRank；

(3) the parent population P0 (i.e. N number of parameter group) that Population Size is N is randomly generated by main program master, and by MPI_ Processor active task is distributed to cSize subprogram (being separately operable hydrological model according to parameter group) by Bcast functions；

(4) according to result of calculation, calculate each individual fitness function in parent population P0 and (receive assorted efficiency factor and relatively inclined Difference), carry out quick non-dominated ranking and crowding is calculated, produce A0 ε non-dominant to achieve solution；

(5) 25% injection mechanism (Injection Scheme) is used, A0 ε non-dominant in (4) achieves solution and form new 25% individuality in parent population Pt, 75% by randomly generating, Population Size is Nt for remaining；

(6) progeny population Qt is produced by selection, intersection, mutation operator, Population Size is Nt；By progeny population Qt and parent Pt Merge the new population Rt=Pt ∪ Qt of composition, Population Size is 2Nt；

(7) processor active task is distributed to by cSize subprogram by MPI_Bcast functions by main program master, and according to fortune Calculating result carries out the calculating of quick non-dominated ranking and crowding, according to each individual non-dominated ranking in colony and it is crowded away from Elite parent population Nt is produced from selection, while producing At ε non-dominant to achieve solution；

(8) the 25% injection mechanism according to (5), by the Population Size of archive population A t, increases random generation sample and forms population Size is the parent population Pt+1 of 4At.

(9) circulation step (6)-(8), evolutionary generation or other end conditions until reaching setting.

5. the hydrological distribution model rating method based on multiple programming and agent approach according to claim 1, its feature Be, step 4) constructed by target Watershed Hydrologic Models include specifically designed for diurnal courses series model model/D models and The model/P models of the peak analogue chosen for super sizing technique, two group models are respectively configured not according to multiple target solution Same parameter group.

6. the hydrological distribution model rating method based on multiple programming and agent approach according to claim 1, its feature Be, step 4) described in the specific algorithm of agent approach be：

$Q_C\left(t\right)=\frac{m_{\mathrm{mod}el/D}\·Q_{\mathrm{mod}el/D}\left(t\right)+m_{\mathrm{mod}el/P}\·Q_{\mathrm{mod}el/p}\left(t\right)}{m_{\mathrm{mod}el/D}+m_{\mathrm{mod}el/P}}$

$m_{\mathrm{mod}el/D}=\left\{\begin{array}{c}1,if{Q}_o/Q_{o,max}<\mathrm{\&gamma;}\\ 0.5,if\mathrm{\&gamma;}\&le;Q_o/Q_{o,\max }<\mathrm{\&delta;}\\ 0,if{Q}_o/Q_{o,\max }\&GreaterEqual;\mathrm{\&delta;}\end{array}\right.$

$m_{\mathrm{mod}el/P}=\left\{\begin{array}{c}0,if{Q}_o/Q_{o,max}<\mathrm{\&gamma;}\\ 0.5,if\mathrm{\&gamma;}\&le;Q_o/Q_{o,\max }<\mathrm{\&delta;}\\ 1,if{Q}_o/Q_{o,max}\&GreaterEqual;\mathrm{\&delta;}\end{array}\right.$

$\mathrm{\&gamma;}=\frac{T}{Q_{o,max}}$

$\mathrm{\&delta;}=\frac{E\left(T\right)}{Q_{o,max}}$

Wherein m_model/DAnd m_model/PIt is the distribution coefficient of model/D and model/P models, Q_model/D(t) and Q_model/P(t) difference It is the linear reservoir of model/D and model/P models, Q_o,maxThe maximum of diurnal courses data is periodically surveyed for rate, T is overdetermination amount The threshold value of method, E (T) is the expectation of the peak value that all super sizing techniques are chosen, Q_oIt is measuring runoff.