*Testing has shown the efficacy of a multi-objective algorithm for optimization of powerhouses with multiple units. This article also presents a comparison of results from the algorithm to actual operational data from the Grand Coulee facility.*

**By Robert Bass, Michael Curtiss, Lee Sheldon and Jonas Parker**

A decision support tool is available to aid owners of multi-unit hydro facilities in minimizing water flow for a given power request. The tool uses sequential quadratic programming (SQP) to minimize the objective function - flow as a function of power setpoints - while adhering to constraints pertinent to a specific powerhouse.^{1}

A case study demonstrating efficacy of the algorithm is presented, consisting of data from four days of operation at the 6,809-MW Grand Coulee facility on the Columbia River in Washington. The U.S. Department of Interior's Bureau of Reclamation (USBR) provided characteristic data from four groups of units (Units 1 through 9, 10 through 18, 19 through 21 and 22 through 24), including capability curves and nominal capacities, condensing capabilities, dispatch priorities, discharge equations, and operational notes. USBR also provided four days of operational data against which the unit selection and load sharing solutions from the algorithm could be compared.

This article discusses Type 2 optimization, SQP and constraints specific to hydropower optimization. The article concludes with a discussion of results from the Grand Coulee case study.

## Theory of Type 2 optimization

There are five types of optimization pertaining to hydroelectric power. The one focused on in this study is Type 2, collective coordination of turbines within a powerhouse.

Type 2 optimization depends on first achieving Type 1 (optimization of individual turbines within a powerhouse), which involves determining a unit's absolute optimum efficiency profile. Methods to achieve Type 1 optimization vary based on turbine type and powerhouse configuration. For Francis, Kaplan and fixed-blade propeller turbines, a table of power vs. flow relationships is assembled at increments of head. Kaplan turbines also require characterization of their blade-to-gate relationship to determine optimal blade positioning. These data are compiled in a table of efficiency vs. power output and a third column is derived indicating the rate of change of flow with respect to power. Equality of this parameter is the criteria for Type 2 optimization.

Type 2 optimization concerns coordination of unit selection and load sharing for all units within a powerhouse to achieve a total power setpoint using minimal flow. Load sharing involves selecting wicket gate openings such that the derivative of flow with respect to power output (dq/dp) is equal for all units.

Even with turbines of the same design, efficiency profiles will be slightly different due to abnormalities in manufacturing and method of operation. Wear-and-tear over years can exacerbate these differences and result in peak unit operating efficiency differences and shifts in peak efficiency points of 2% to 10%, based on one author's experience with testing Francis turbines.

Optimization algorithms take advantage of these differences by proper unit selection and load sharing, given a power setpoint and known head. Experienced operators can make these decisions with impressive results, although there is room for additional improvement through the use of an efficient optimization algorithm.

USBR provided discharge equations characterizing groups of units. These equations provide flow as a function of power setpoint and head:

Equation 1

*q(p, h) = A + Bp + Ch ^{-1} + Dp^{2} + Eh^{-2} + Fph^{-1} + Gp^{3} + Hph^{-2} + Jp^{2}h*

where:

- q is flow in cubic feet per second;

- p is power output in MW;

- h is head in feet; and

- A through J are characteristics of the turbine groups.

USBR and Bonneville Power Administration (BPA) use these discharge equations in the day-to-day operation of Grand Coulee, as well as for planning and modeling. The groups are each characterized by a unique discharge equation that is used to derive the dq/dp curves:

Equation 2

Each unit would have its own dq/dp curve, thereby taking advantage of the variation between all turbines to achieve optimal unit selection and load sharing, rather than characterizing units in groups. Nonetheless, the use of group discharge equations proved sufficient to test the algorithm's capabilities.

## Sequential quadratic programming

SQP, a non-linear technique used to minimize an objective function subject to equality and inequality constraints, has been applied to various aspects of power systems, such as distribution system design, placement of VAR compensators and more.^{2,3,4,5,6,7,8} Systems whose objective functions can be represented as sets of polynomials are good candidates for SQP optimization. If the objective function cannot be well-represented as a polynomial, SQP is not appropriate.

For Type 2 optimization, the objective function is total flow to be minimized, and p is a vector of independent variables representing the units' power setpoints. SQP determines the p vector of power setpoints such that q(p) is minimized.

Equation 3

s.t. *c _{i} (p) = 0 i = l,...,m_{c}*

and *h _{i} (p) ≥ 0 i = m_{c} + l,...,m*

where:

- q(p) is the objective function;

- c(p) is the mc equality constraint; and

- h(p) is the m-mc inequality constraint;

Vectors of Lagrange multipliers ( for equality constraints and π for inequality constraints) are used to constrain the objective function, resulting in a Lagrange function - L(p,λ,π) - that is simply a constrained version of the objective function q(p).

Equation 4

The gradient of the Lagrangian when set to zero gives the conditions under which the optimum may be found.

Equation 5

An iterative quasi-Newton method is applied to solve for p, and π that uses a second-order Taylor series expansion to provide a quadratic approximation of the Lagrangian at iteration k+1:

Equation 6

s.t.

and

where:

- H(p_{k}) is the Hessian of the Lagrange function, expressed as

An approximate Hessian is updated after every iteration. Using these steps, SQP can be used to solve constrained nonlinear problems. Provided with reasonable initial guesses for the arguments of the Lagrange function, the algorithm should converge within a practical number of iterations.

## Constraints

Unit commitment and load sharing of synchronous turbines within a powerhouse are subject to constraints: minimum loading, up-margin, condensing/motoring, run and down times, and start-up/shut-down priorities. To optimize powerhouse operation, these constraints must factor into the algorithm. We categorize these constraints as non-temporal inequality, temporal inequality, equality and "other."

*Non-temporal inequality constraints*

Non-temporal inequality constraints, h(p), confine the search space for the minimized solution to q(p) by modifying the Lagrangian as shown in Equation 4. These constraints are time-independent. If a constraint has an upper and lower bound - say hmin ≤ h(p) ≤ hmax such as for minimum and maximum loading - it is redefined as two constraints ha(p) = h(p) - hmin ≥ 0 and hb(p) = hmax - h(p) ≥ 0 to match the form in Equation 3.

- Rough zones/cavitation zones constraints. Also known as "exclusion zones," these constraints define bands of operation where turbine run-out and/or cavitation is high. When operating a unit within these zones, the goal is to ramp through as quickly as possible.

Exclusion zone constraints may be addressed in two ways. First, if the best solution places a turbine within an exclusion zone, it is thrown out and the next-best solution chosen. This process iterates until the best solution is found that does not place any turbines within an exclusion zone. Second, inequality constraints defining exclusion zones for every generator are included within the definition of the Lagrangian. The optimal solution is then calculated based on these constraints. The latter approach is more computationally efficient because it restricts the search space for the solution to the objective function and does not require multiple iterations to determine a second-best solution.

- Steady-state unit constraint. When a unit is generating, its load must not be changed. A steady-state constraint is handled by setting the unit's power output, subtracting that output from the powerhouse power setpoint, and then excluding it from unit selection within the optimization algorithm.

- Minimum loading constraint. When a unit is operating, it must be running at or above a defined generation load. This is the same as an exclusion zone boundary extending from 0 MW to the given lower megawatt bound.

- Up-margin constraint. Up-margin represents the additional power that can be generated without starting another unit. Up-margin reserve may be built into the program as an upper-bound constraint. When additional power is needed, the upper-bound constraint would be shifted upward. The load-sharing portion of the algorithm would then be re-run (excluding the unit-selection portion) to determine new power setpoints for all units given the new powerhouse setpoint.

- Shared penstock constraint. If two generators share a penstock, the discharge rate of each unit can depend upon the loading of the other unit. This is accounted for in the optimization algorithm by allowing for multiple discharge equations for a single generator. When the constraint begins to affect efficiency, one or both units will be dropped. Whichever of the two is the better choice will return, individually, if it fits into the solution.

*Temporal inequality constraints*

Temporal inequality constraints are a time-limited form of the inequality constraints h(p). These constraints enter into the Lagrangian if a solution is sought within a specified time frame. The time frame is triggered by a unit shutdown or start-up. After the time frame has expired, the constraint is removed.

Four temporal inequality constraints are applied to multi-unit hydropower optimization: minimum run-time, how long a unit must be run before being shut down; maximum run-time, the maximum time a unit may be run before being shut down; minimum down-time, the minimum time a unit must stay shut down before being started; and maximum down-time, the maximum time a unit can be shut down before it must be restarted.

Because of the temporal nature of these constraints, any optimization algorithm employing them would have to be run on a continuous loop and incorporate an internal clock. As a stand-alone decision support tool, the optimization algorithm does not account for temporal inequality constraints. However, it would be straightforward to modify the software to incorporate these constraints. This could be done by monitoring the on/off status tag of each turbine via hooks to the supervisory control and data acquisition system, which would then be used to include or remove the appropriate constraints within the Lagrangian.

*Equality constraints*

Equality constraints, c(p), confine the search space for the minimized solution to q(p). These constraints either set specific required values that modify the Lagrangian or result in a unit being removed from unit selection consideration.

- Must-run constraint. A constraint must be made if the unit should not be shut down, such as a pump that is electrically tied to a specific generator. While pumping, that generator must remain generating. A second example is a unit or group of units armed for gen drop, which is a mechanism that can be used by a power marketing agency to quickly shed generation if a problem occurs with the transmission system.

There are two possibilities for implementing this constraint. First, the unit can be "fixed" to generate a specific output and excluded from the optimization matrix using an equality constraint. Second, if specific output is not required but the unit must stay on, the unit would be a required part of the solution but its load sharing would be optimized. For the first example, a minimum inequality constraint would be set to ensure sufficient power was available to meet the needs of the pump.

- Motoring/condensing constraints. Units may be motored (kept spinning but not producing power) or condensed (kept online to provide reactive power support) instead of shut down. This could be done for operational or voltage support issues. Motoring may be done to limit cycling on breakers, for example. Condensing is often used to buck VARs from the system during light loading periods.

To implement these constraints, the real power consumption of the unit when it is in a condensing or motoring state is added to the powerhouse setpoint such that additional power from adjacent units is available to meet the power requirements set by dispatch and real power consumption of the condensing units.

- Ancillary services constraints. These include spinning reserve and voltage/frequency reserve. In most cases, spinning reserve is the same as up-margin. There are a few exceptions where condensing units might contribute to spinning reserve. Operation of a plant to meet voltage/frequency control requirements is typically driven by a unit's governor. An example would be to run a unit at speed-no-load only for VAR support and not power generation.

These constraints can all be handled by fixing generation levels of selected units outside the optimization matrix. Units that are in condensing service would simply be noted as such and not enter the optimization matrix.

*Others*

Two other constraints of the optimization algorithm do not fall within the equality or inequality definitions pertinent to Equation 3. These are priorities start-up and priorities shut-down, which are definitions of priorities that should be used when determining an order in which to stop or start units. Examples are trying to start Unit A before Unit B because of a bus issue or unit location or trying to shut down Unit B before Unit A.

These constraints are handled by withholding a unit from being added/removed before its paired unit.

## Results

Through our case study of Grand Coulee, we have demonstrated potential flow savings of 1% to 2%. A decrease in flow within this range may be possible if operators use this SQP-based algorithm as a decision support tool. Further, the algorithm calculated results for Grand Coulee in less than one second, allowing operators to make unit selection and load sharing adjustments without delay. Adjustments to unit selection and/or load sharing at Grand Coulee were made roughly every 10 seconds.

Figure 1 on page 38 shows flow as a function of time for Grand Coulee on Aug. 11, 2007 (called the "actuals") in comparison to the flow that would have resulted if the algorithm were used as a decision support tool. For this data set, the algorithm would have resulted in the use of 0.92% less water. The difference in flows between the actuals and the algorithm has an energy value of about 734 MWhQ (water potential energy) over the day.

Table 1 on page 44 presents results from the four sets of data against which the algorithm was verified, which included information on megawatt request, up-margin, forebay and tailbay elevations (head), individual generator status (online or offline) and generator megawatt output (condensing mode provides a negative number), updated about every ten seconds. Using the megawatt request and head setpoints, the algorithm was set to determine unit selection and load sharing for each time interval, thereby producing its own generator status profiles and megawatt output values. Using the group discharge equations, we compare the flow resulting from the actuals to that resulting from the algorithm set-points. Table 1 presents these results in terms of energy per day of flow (in MWhQ), from which we derive an estimate of the electricity "savings" (in MWhe) based on the moment-by-moment efficiency of the powerhouse. Because the power setpoint in both cases is identical, we know there are no actual electrical energy savings. Rather, translation of the conserved flow to MWhE allows us to place a momentary price signal on the water savings, provided a reasonable price per MWhE, such as the Mid-Columbia on-peak spot price.

*Verification of adherence to constraints*

The USBR performance data describe unit operation under non-temporal inequality and equality constraints, including up-margin, shut-down and start-up priorities, motoring/condensing, minimum loading and rough zones. With the exception of rough zones, all of these constraints were incorporated into the algorithm (temporal constraints were not incorporated). The algorithm's adherence to these constraints was verified by analyzing the results for constraint violations, as discussed below.

- Up-margin. The optimization program can be set to provide up-margin for the powerhouse. To validate, the output data files were examined to ensure that the up-margin was not violated. The sum of megawatt output and the up-margin request should never exceed the powerhouse nameplate capacity minus all offline units. For all four data sets, this was found to be the case.

- Shut-down and start-up priorities. At Grand Coulee, two specific units should be started first and stopped last. The program was modified to ensure these units adhered to these priorities, meaning they were removed from the unit selection part of the optimization algorithm. To validate, data output files were examined to ensure these units were the first started and last stopped. For all four data sets, this was found to be the case.

- Minimum loading. At Grand Coulee, Units 1 through 21 must not be set to less than 50 MW, and Units 22 through 24 must not be set to less than 25 MW. The program allows an operator to input minimum loading constraints for every generator. To validate, the algorithm output for each generator was checked to ensure power did not drop below the minimum loading threshold (see Figure 2 on page 38). For all four data sets, the minimum loadings were met.

- Motoring/condensing. Units 19 through 24 should not be shut down but should be made available for motoring/condensing, in which case the megawatt losses incurred must be accounted for. Algorithm output was checked to ensure the motoring/condensing powers were subtracted from total powerhouse output when these units were not part of a unit selection solution. Figure 2 demonstrates this constraint was met on Aug. 11, 2007. When taken offline, the units are set to motor/condense, so they consume power (negative power). For all four data sets, this was found to be the case.

*Verification of optimization*

Multi-unit optimization may be understood initially in terms of d/dp, the rate of change of efficiency with respect to power. If this derivative equals zero, the tangent to a curve of efficiency versus power is parallel to the power axis. The point of intersection between the tangent and curve is the location of peak efficiency. However, if two identical units are operating, their combined efficiency is maximized if both are operated such that the slopes of the derivatives are equal and not necessarily zero. An infinitesimal flow rate cannot be taken from one machine and given to another in any way that would result in increased combined efficiency. If the units are of different size (or different efficiency profiles), an increase of efficiency on the larger machine produces a greater increase in combined efficiency than the same increase on the smaller machine. In other words, the value of the derivative must be weighted based on size.

Equation 7 gives the basic definition of fluid power, pHP, in units of water horsepower at a constant head and specific weight of water, .

Equation 7

If the algorithm recommends optimized unit selection and load sharing, the dq/dp values for every online unit within the powerhouse will be equivalent. For the Grand Coulee case study, the algorithm results in equivalent dq/dp values for all units at all setpoints, as demonstrated in Figures 3 and 4.

Figure 3 on page 40 compares the difference between the dq/dp results from the algorithm and those from the actuals at every setpoint. The degree of optimization is visible when one overlays the dq/dp curves of several units. When overlaid, and if optimized, the dq/dp curves should all overlap. The overlay of dq/dp from the actuals shows there are variations in dq/dp over time between the various units (see top). While the actuals efficiency performance is very good, there is room for improvement. The bottom overlay shows nearly identical dq/dp curves at all times, indicating the algorithm is behaving according to theory.

Figure 4 on page 42 shows the derivative of the group discharge equations for the four turbine groups (solid lines). Overlapping these curves are four loci of dq/dp values derived from the algorithm results. The curves and data share the same head of 312 feet; recall that the units' discharge equations and therefore the dq/dp curves are functions of head. Each locus derives from a different power setpoint. Note all the data points within each locus share a dq/dp value, indicating that the algorithm has optimally allocated load sharing between the units and that all online units are adhering to their respective discharge equation curves.

## Conclusion

The Grand Coulee case study demonstrates that powerhouse efficiency may be improved if operations are complimented with an SQP-based decision support tool. However, these gains were achieved using group discharge equations rather than unit-specific efficiency curves, suggesting further improvements may be possible if all units are individually characterized such that variations between units may be exploited to maximize optimization. Such characterization could be done via absolute index testing of every unit to ascertain individual absolute performance characteristics.9,10 Nonetheless, the use of group discharge equations proved to be sufficient to suggest possible operational gains of 1% to 2%. Estimating the value of the saved flow energy in terms of electrical energy, and grossly extrapolating based on the four data sets presented in the case study, the about 627 MWh per day saved amount to an annual savings of about $6.9 million, as presented in Table 1 on page 44.

To realize such gains, a decision support tool must be able to return results in significantly less time than the mean hold time of the head and power setpoints. The algorithm returns results for the 24-unit Grand Coulee facility in less than one second. Given that adjustments in unit selection and load sharing in the Grand Coulee datasets occur roughly every 10 seconds, the algorithm is sufficiently fast to serve as a decision support tool.

The authors would next like to verify the algorithm against historical data from a hydro facility that has characterized discharge equations for all units. Beyond that, a step toward realizing this tool as a means for real-time optimization would be to use it as a stand-alone decision support tool. With a properly designed graphical user interface, an operator could refer to this tool for recommendations regarding unit selection and load sharing, without the algorithm directly interfacing with the supervisory control and data acquisition system, thereby leaving the operator in full control of the plant. If proven successful for this application, a next step could be to establish hooks from the algorithm directly into the SCADA system of a hydro plant such that unit selection and load sharing results would automatically actuate the facility's turbines.

## Acknowledgments

The authors thank Toby Steves and his colleagues at Reclamation for their generous support, specifically in providing operational information and performance data for Grand Coulee. This work was partially supported by BPA through a Technology Innovation Opportunity grant, project number 1918-1556.

## Notes

^{1}Curtiss, Michael, Jonas Parker and Parker Scoggins, "Optimizing Operation of Multi-Unit Powerhouses: A New Method," *Hydro Review*, Volume 31, No. 5, July 2012, pages 88-98.

^{2}Computational Methods for Electric Power Systems, 2nd Edition, CRC Press, Boca Raton, Fla., 2009.

^{3}Finardi, E., and E. d. Silva, "Solving the Hydro Unit Commitment Problem via Dual Decomposition and Sequential Quadratic Programming," IEEE Transactions on Power Systems, Volume 21, No. 2, May 2006, pages 835-844.

^{4}Al Hajri, M., and M. El-Hawary, "Optimal Distribution Generation Sizing via Fast Sequential Quadratic Programming," Proceedings of 2007 Large Engineering Systems Conference on Power Engineering, IEEE, New York, N.Y., 2007.

^{5}Martins, L., S. Soares and A. Azevedo, "A Nonlinear Model for the Long-Term Hydro-Thermal Generation Scheduling Problem over Multiple Areas with Transmission Constraints," Proceedings of IEEE/PES Power Systems Conference and Exposition, IEEE, New York, N.Y., 2009.

^{6}Sivasubramani, S., and K. Swarup, "Sequential Quadratic Programming Based Differential Evolution Algorithm for Optimal Power Flow Problem," IET Generation, Transmission & Distribution, Volume 5, No. 11, November 2011, pages 1149-1154.

^{7}Pour, M., and A. Khodabakhshian, "A New Robust Load Frequency Control Design using Sequential Quadratic Programming (SQP) Method," Proceedings of 15th IEEE Mediterranean Electrotechnical Conference, IEEE, New York, N.Y., 2010.

^{8}Khandani, F., S. Soleymani and B. Mozafari, "Optimal Placement of SVC to Improve Voltage Profile using Hybrid Genetics Algorithm and Sequential Quadratic Programming," Proceedings of 16th Conference on Electrical Power Distribution Networks, IEEE, New York, N.Y., 2011.

^{9}Sheldon, Lee, "Optimizing the Generating Efficiency of Entire Powerhouses," Proceedings of HydroVision International 2008, PennWell Corporation, Tulsa, Okla., 2008.

^{10}Sheldon, Lee, "Method to Develop a Family of Cam Curves from a Single Index Test," HRW-*Hydro Review* Worldwide, Volume 20, No., 2, March-April 2012, pages 40-45.

**Robert Bass is an associate professor of electrical engineering at Portland State University. Michael Curtiss is a distribution engineer with Pacific Power. Lee Sheldon is an adjunct professor of hydropower engineering at the Oregon Institute of Technology. Jonas Parker is an engineer with the U.S. Army Corp of Engineers' Hydroelectric Design Center.**