Synthetic Unit Hydrographs

Least Squares Procedure

The following procedure produces estimates of the unit hydrograph (UH) by minimizing the sum of the squares of the estimation errors. Estimation errors are defined as the difference between the observation and the prediction, (Q observed - Q predicted ).

The discrete convolution equation has been rewritten in matrix form as:

[Q]=[P][U]

where [Q] is the column vector containing the N direct runoff observations; [P] is the precipitation matrix containing the volumes of M rainfall pulses; and [U] is the column vector containing the unit hydrograph (UH) ordinates. In what follows, N is the number of direct runoff ordinates, M is the number of effective rainfall pulses, and DT is the duration of the observation intervals. Thus, the instantaneous direct runoff is observed at every DT, and the rainfall pulses are volumes of rainfall of constant intensity over DT.

The following table presents the effective rainfall hyetograph and the direct runoff hydrograph for a hypothetical basin. Observations of discharge were made at 3/4-hour intervals, and the hyetograph represents averages over 3/4-hour periods.

Time(h) Average Intensity (cm/h)

0.00 - 0.75 0.2

0.75 - 1.50 1.5

1.50 - 2.25 1.0

2,25 - 3.00 0.5

The corresponding direct runoff discharge hydrograph is:

Time(h) Discharge (m3/s)

0.75 30

1.50 250

2.25 500

3.00 400

3.75 180

4.50 30

Initialize linear algebra library.

> with(linalg);

> with(student);

> N:=6; M:=4;

Define discharge hydrograph matrix.

> Q:=matrix(N,1,[30,250,500,400,180,30]);

Define rainfall hyetograph matrix.

> H1:=matrix(M,1,[0.2,1.5,1.0,0.5]);

Define time interval, DT, in seconds.

> DT:=2700;

Compute area of the basin by determining the area A that would produce an input volume (effective rainfall volume) equal to the observed outflow volume (direct runoff volume).

> A:=sum('Q[k,1]','k'=1..N)*DT/(sum('H1[l,1]','l'=1..M)*DT/100/3600);

Transform the rainfall hyetograph matrix from units of intensity (cm/h) to units of flux (m3/s).

> H:=scalarmul(H1,A/3600/100);

Define the P matrix of the discrete convolution equation. Observe that now both P and Q are in consistent units of m3 and m3/s.

> P:= matrix(6,3,[[H[1,1]*DT,0,0],[H[2,1]*DT,H[1,1]*DT,0],[H[3,1]*DT,H[2,1]*DT,H[1,1]*DT],[H[4,1]*DT,H[3,1]*DT,H[2,1]*DT],[0,H[4,1]*DT,H[3,1]*DT],[0,0,H[4,1]*DT]]);

Obtain pseudoinverse to solve the discrete convolution equation.

> PT:=transpose(P);

> PTP:=multiply(PT,P);

> PTPINV:=inverse(PTP);

> Z:=multiply(PTPINV,PT);

Obtain the Least Squares estimate of the UH ordinates, ULS.

ULS:=multiply(Z,Q);

>

Obtain the Forward Substitution estimate of the UH ordinates, UFS.

> UFS:=matrix(3,1,[Q[1,1]/P[1,1],(Q[2,1]-P[2,1]*Q[1,1]/P[1,1])/P[1,1],(Q[3,1]-P[3,1]*Q[1,1]/P[1,1]-P[2,1]*(Q[2,1]-P[2,1]*Q[1,1]/P[1,1])/P[1,1])/P[1,1]]);

> sum('ULS[m,1]','m'=1..N-M+1)*DT;

> sum('UFS[m,1]','m'=1..N-M+1)*DT;

Predict discharge for observed rainfall using ULS.

> QHATLS:=multiply(P,ULS);

Predict discharge for observed rainfall using UFS.

> QHATFS:=multiply(P,UFS);

Define new P matrix using the following observed effective rainfall hyetograph.

Time(h) Average Intensity (cm/h)

0 - 0.75 0.3

0.75 - 1.5 0.7

1.5 - 2.25 1.5

2.25 - 3.0 0.3

3.0 - 3.75 0.1

> H1NEW:=matrix(5,1,[0.3,0.7,1.5,0.3,0.1]);

> HN:=scalarmul(H1NEW,A/3600/100);

> PNEW:= matrix(7,3,[[HN[1,1]*DT,0,0],[HN[2,1]*DT,HN[1,1]*DT,0],[HN[3,1]*DT,HN[2,1]*DT,HN[1,1]*DT],[HN[4,1]*DT,HN[3,1]*DT,HN[2,1]*DT],[HN[5,1]*DT,HN[4,1]*DT,HN[3,1]*DT],[0,HN[5,1]*DT,HN[4,1]*DT],[0,0,HN[5,1]*DT]]);

Use the Least Squares estimate of the UH to predict the discharge that would be observed for that rainfall event.

> QHATLS:=multiply(PNEW,ULS);

Use the Forward Substitution estimate of the UH to predict the discharge that would be observed for that rainfall event.

> QHATFS:=multiply(PNEW,UFS);

> sum('ULS[m,1]','m'=1..N-M+1)*DT;

> sum('UFS[m,1]','m'=1..N-M+1)*DT;

>