... 
x, y, z = None, None, None 
RV, x, y, z = HECobject.method( x, y, z)
```
**Script 4.** Code sample 4. Passing the arguments by reference to the method of the HEC-RAS object playing the role of a function with one return value.

The main codes presented in Scripts 5–8 consist of the main elements necessary to access capabilities of the HEC-RAS application in Python code. There are also additional codes written in support.py. They are used for support of the main scripts. The support.py file contains the subroutines for data loading from ASCII files, for handling files and directory names, etc. All these additional subroutines are subject to the coder and may be written in several ways. Hence, they are not discussed here. The additional codes are loaded to the particular scripts, if there is such a need.

### **3. Results and Discussion**

### *3.1. Case Study—Model of a River Reach*

Two HEC-RAS models, called model A and model B, were used to test the presented concepts. Both of them represented a short reach of the Warta river in the central part of Poland. The reach is about 10 km long. The average width and depth of the main channel are 40 m and 1.5 m, respectively. There are two bridges along the reach. The estimated mean flow in the analyzed channel was about 47.45 m3/s. This value was estimated on the basis of the hydrological data collected at the Sieradz gauge station during 1970–2013. The observed flows for the same period vary from the minimum of 15 m3/s to the maximum of 408 m3/s [63].

In both models, the number of cross-sections was 35. The average distance between cross-sections was about 300 m. The differences between the model A and model B include the geometry and configuration of flow conditions. Model A was based on the full geometry, including bridges. The flow conditions were steady with basic discharge of 120 m3/s for the whole reach. The downstream

boundary condition was the simple Normal Depth with bottom slope as a parameter. This model was used for basic simulation in Example 1. It was also applied in Example 2, for the test of calibration.

In model B, the geometry was simplified. The bridges were removed, because the model was used for simulation of sediment transport with different sediment samples. Hence, the configuration of flow conditions was based on a quasi-unsteady flow module. The upstream boundary condition was Flow Series with a constant discharge of 269 m3/s. In the downstream boundary the Normal Depth was applied once again. The slope equals 0.2‰. The sediment transport intensity was calculated on the basis of the Meyer-Peter and Müller (MPM) formula. The simulation horizon was set to 4 months, with a 6-h time step. The time step guarantees necessary numerical stability and accuracy. The simulated period was rather short. In practice, the changes of the river bed due to sediment deposition and erosion are observed in longer periods, e.g., years. Such a configuration is chosen to avoid too long computational time in this illustrative example.

### *3.2. Example 1—Basic Simulation*

The first example presented here was composed in a way similar to the introductory example discussed by Goodell [20,21] for the reach of Beaver Creek. In this paper, model A described above was used. The main scheme of computations is shown in Figure 1. The code of this example is presented as Script 5.

```
1. import win32com.client 
2. import os, numpy 
3. from support import sc1_ShowNodes 
4. 
5. # Initiate the RAS Controller class 
6. hec = win32com.client.Dispatch("RAS503.HECRASController") 
7. hec.ShowRas() # show HEC-RAS window 
8. # full filename of the RAS project 
9. RASProject = os.path.join(os.getcwd(),r'test01\test01.prj') 
10. hec.Project_Open(RASProject) # opening HEC-RAS 
11. # to be populated: number and list of messages, blocking mode 
12. NMsg,TabMsg,block = None,None,True 
13. # computations of the current plan 
14. v1,NMsg,TabMsg,v2 = hec.Compute_CurrentPlan(NMsg,TabMsg,block) 
15. # ID numbers of the river and the reach 
16. RivID,RchID = 1,1 
17. # to be populated: number of nodes, list of RS and node types 
18. NNod,TabRS,TabNTyp = None,None,None 
19. # reading project nodes: cross-sections, bridges, culverts, etc. 
20. v1,v2,NNod,TabRS,TabNTyp = \ 
21. hec.Geometry_GetNodes(RivID,RchID,NNod,TabRS,TabNTyp) 
22. # ID of output variables: WSE, ave velocity 
23. WSE_id,AvVel_id = 2,23 
24. TabWSE = numpy.empty([NNod],dtype=float) # NumPy array for WSE 
25. TabVel = numpy.empty([NNod],dtype=float) # NumPy array for velocities
26. for i in range(0,NNod): # reading over nodes 
27. if TabNTyp[i] == "": # simple cross-section 
28. # reading single water surface elevation 
29. TabWSE[i],v1,v2,v3,v4,v5,v6 = \ 
30. hec.Output_NodeOutput(RivID,RchID,i+1,0,1,WSE_id) 
31. # reading single velocity 
32. TabVel[i],v1,v2,v3,v4,v5,v6 = \ 
33. hec.Output_NodeOutput(RivID,RchID,i+1,0,1,AvVel_id) 
34. hec.QuitRas() # close HEC-RAS 
35. del hec # delete HEC-RAS controller 
36. 
37. sc1_ShowNodes( NNod, TabRS, TabNTyp, TabWSE, TabVel )
```
**Script 5.** Example 1. Calculation of a single water profile.

**Figure 1.** Main elements of computations in Example 1.

As can be seen in the first block of the diagram in Figure 1, the computations have to start with initialization of the variable, which is the HECRASController object. In the brackets, the method used for this purpose is indicated. In the second step, the HEC-RAS project is open. The name of the project is the necessary input. Then the computations are run and the results are processed in two steps. The first is reading information about computational nodes. It is done with one of the "geometric" routines of HECRASController. The second step is reading the output of computations. Then the results may be saved, displayed, or applied in another simulation.

The program presented in Script 5 includes four main steps reflecting the sequence presented in Figure 1. These are: (a) Opening HEC-RAS and loading the project (lines 9–10), (b) running the simulation (lines 12–14), (c) reading the output variables (lines 16–35), and (d) processing the output variables (line 37), e.g., displaying them on the screen. At the beginning, the necessary modules and subroutines are loaded. They are win32com.client for handling COM libraries [34,35], os for management of files and directories and numpy for numerical data. Additionally, the subroutine sc1\_ShowNodes is loaded from the script support.py. It is used for displaying results.

The first step of the program is initialization of the HECRASController variable hec (line 6). The Dispatch function of the win32com.client module is applied. The HEC-RAS version used is 5.0.3. When the hec variable is created, the ShowRas and Project\_Open functions may be used to open the HEC-RAS windows and load the project. Detailed descriptions of these and other functions used in the paper were provided by Goodell [20]. The next two methods of HECRASCcontroller used are Compute\_CurrentPlan (line 14) and Geometry\_GetNodes (line 21). The first function runs the current computational plan. The second one reads tables of River Stations, TabRS, and node types, TabNTyp. Both of them need initialization of proper arguments (lines 12 and 18). The approach presented in Code Samples 2-4 is applied. The river and reach ID numbers are also set in line 16. The variables v1, v2, and in general vX, where X is some digit, are used when the return value does not have to be stored.

The next step is reading of the output variables. It starts at line 23, with setting of the ID numbers for water stages and velocities, WSE\_id and AvVel\_id. Then the NumPy arrays, TabWSE, and TabVel, for storing the results are prepared (lines 24–25). The results are read in a loop over nodes (line 26–33), with calls to the function Output\_NodeOutput. The values of water stages and velocities in regular cross-sections are accessed. In such nodes, the node type is an empty string (line 27).

Finally, the HEC-RAS windows is closed with QuitRas. Then the hec variable may be deleted. The last element of the code is the display of results. It is performed by sc1\_ShowNode from the support module. A screenshot with the run of Example 1 and display of results is shown in Figure 2.

**Figure 2.** Run of Example 1 and display of results.

### *3.3. Example 2—Calibration of Roughness Coefficients*

The problem of model calibration was considered. Model A described above was used here. The roughness coefficients were sought. The computations were performed in the steady mode. The calculated water surface profile should fit the imposed "observed" values. In fact, the "observed" values were calculated earlier for known roughness coefficients. Very simple distribution of their values was assumed. For the main channel along the whole reach, the value 0.025 was applied. For the left and right floodplains, one value of 0.04 was set in the same way. The purpose of the program was to find the set of optimal coefficients. It was assumed that only three values were necessary. Hence, the dimension of the problem was relatively small, but the generality of this solution is not lost.

The "observed" values were determined for each cross-section along the reach. However, they were applied with different frequencies in the search algorithm. It means that only some values were selected, and the simulated water surface was fitted only in the selected cross-sections. Hence, the objective function is as follows:

$$F(\mathbf{x}) = \sqrt{\frac{1}{N} \sum\_{i=1}^{N} \left( WSE\_i^{(c)} - WSE\_i^{(o)} \right)^2} \tag{5}$$

In the above formula, *N* is number of observations, *WSE*(*c*) *<sup>i</sup>* is the computed water surface elevation, and *WSE*(*o*) *<sup>i</sup>* is the "observed" one determined in the same *i*-th cross-section. This function should be minimized to find the optimal set of roughness coefficients *x*, where *x* = [*nLOB*, *nCh*, *nROB*] for the left floodplain, the main channel, and the right floodplain, respectively.

The algorithm for calculation of the objective function is presented in Figure 3. It starts with processing of the input parameters. At this step, the optimization variable is transformed into computational roughness coefficients. Then the coefficients are set in all model cross-sections. This step requires access to the global variable representing the object of the HECRASController. The next elements of the algorithm are similar to the operations presented in Figure 1. The computations start and the results are read. In both these steps, the HECRASController object is processed. To read the results properly, identification of the nodes where referenced or "observed" values of the water surfaces are located is necessary. This information is taken from another global variable. Then the calculated and referenced water surface elevations are used to determine the final value of the objective function according to Equation (5). The reference water surfaces are also stored in the global variable.

**Figure 3.** Algorithm for calculation of objective function values.

The main elements of the computations in this example are shown in Figure 4. The computational process starts with initialization of the HECRASController object. Then the HEC-RAS variable is used for opening of the HEC-RAS project, setting the current computational plan and running preliminary computations. These steps are necessary to compare the information about the reference nodes and values with the structure of the computational model. For this purpose, the project name and the reference values must be loaded. After these steps, the object of the HECRASController is ready and the numbers of computational nodes for comparisons with reference nodes are known. After loading of the initial estimate, which is the set of preliminary roughness coefficients, the optimization process starts. It uses the objective function presented in Figure 3.

**Figure 4.** Main computational elements of Example 2.

The starting point is an arbitrarily chosen set of roughness coefficients. The values chosen were 0.08 for the left floodplain, 0.014 for the main channel, and 0.08 for the right floodplain.

In Script 6, such modules as win32com.client, os, math, numpy and scipy.optimize are necessary. They are loaded at the beginning of the script. The function for reading the "observed" values is also loaded from the support module.

```
1. import win32com.client, os, math 
2. import numpy as np 
3. from scipy import optimize as opt 
4. from support import sc2_ReadObs 
5. 
6. def ObjFun(x): # objective function 
7. global hec,RivID,RchID,WSE_ID,ProfID,RivName,RchName,nRS,RScmp,RStyp 
8. global Nobs,iRSobs,WSEobs 
9. nLOB,nCh,nROB = x[0],x[1],x[2] 
10. for i in range(0,nRS): 
11. ErrMsg = None # list of error messages 
12. v0,v1,v2,v3,v4,v5,v0,ErrMsg = \ 
13. hec.Geometry_SetMann_LChR(RivName,RchName,RScmp[i],nLOB,nCh,nROB,ErrMsg)
14. NMsg,ListMsg,block = None,None,True # no. and list of messages, blocking
15. v1,NMsg,ListMsg,v2 = hec.Compute_CurrentPlan(NMsg,ListMsg,block) 
16. TotSum = 0.0 # total sum of square errors 
17. for i in range(0,Nobs): 
18. wse,v1,v2,v3,v4,v5,v6 = \ 
19. hec.Output_NodeOutput(RivID,RchID,iRSobs[i],100,ProfID,WSE_ID) 
20. TotSum += (wse - WSEobs[i])**2 # single square error added 
21. TotSum = math.sqrt( TotSum / Nobs ) # final value 
22. return TotSum 
23. 
24. def PokazIter(xk): # iteration control function 
25. print 'nLOB=%15.6f, nCh=%15.6f, nROB=%15.6f, ' % \ 
26. (xk[0],xk[1],xk[2]), 
27. print ' funk=%15.6f' % (ObjFun(xk)) 
28. return 0 
29. 
30. # observed: number, RSes as strings, RSes as floats, WSEs 
31. Nobs,lRSobs,dRSobs,WSEobs = sc2_ReadObs('test01','ObsH','ObsH1x.txt') 
32. # init HEC-RAS Controller 
33. hec = win32com.client.Dispatch('RAS503.HECRASController') 
34. RivID,RchID = 1,1 # ID of river and reach 
35. ProfID,WSE_ID = 1,2 # ID of profile, ID of WSE variable 
36. projekt = os.path.join(os.getcwd(),r'test01\test01.prj') # project file 
37. hec.Project_Open(projekt) # opening HEC-RAS 
38. hec.ShowRas() # show HEC-RAS window 
39. # setting current computational plan 
40. test1 = hec.Plan_SetCurrent('plan_basic') 
41. # river and reach names 
42. RivName,RchName = 'Warta','Dop_Jeziorsko' 
43. # number and list of messages, blocking mode 
44. NMsg,TabMsg,block = None,None,True 
45. # computations of the current plan 
46. v1,NMsg,TabMsg,test2 = hec.Compute_CurrentPlan(NMsg,TabMsg,block) 
47. # numbers of nodes with observations 
48. iRSobs = np.empty([Nobs],dtype=int) 
49. for i in range(0,Nobs): 
50. iRSobs[i],v1,v2,v3 = hec.Output_GetNode(RivID,RchID,lRSobs[i]) 
51. # computational: number, RS and types 
52. nRS,RScmp,RStyp = None,None,None 
53. v1,v2,nRS,RScmp,RStyp = hec.Output_GetNodes(RivID,RchID,nRS,RScmp,RStyp) 
54. 
55. print "\nIteration process: Nelder-Mead simplex" 
56. x0 = np.array([0.08,0.014,0.08]) # initial solution 
57. Xopt = opt.fmin(ObjFun,x0,callback=PokazIter) # optimization 
58. 
59. hec.QuitRas() # close HEC-RAS 
60. del hec # delete HEC-RAS controller
```
**Script 6.** Example 2. Calibration of HEC-RAS steady flow model.

Script 6 includes two functions. The first one plays the role of an objective function and is called ObjFun. It is an argument of the optimization algorithm. Its structure is discussed below. The second one, PokazIter, is used by the same algorithm to display results of the single iteration.

The main body of the script starts with the reading of "observed" values in line 31. They are stored as lists of River Stations, lRSobs and dRSobs, and a list of water surface elevations, WSEobs. Then the HEC-RAS variable, hec, is created (line 33). Before the project is opened (line 37), the HEC-RAS window is displayed (line 38) and the current plan is set (line 40), the project constants are set (line 34–35). They are the river and reach IDs, RivID and RchID, profile number, ProfID, and ID of the output variable storing water surfaces, WSE\_ID. The names of the river and the reach, RivName and RchName (line 42), are also hard-coded in the script, though they could be read using the methods of the HEC-RAS controller.

The first computations are run in line 46. The necessary variables are prepared earlier (line 44). The only purpose of these preliminary computations is to check the numbers of nodes, at which the "observed" values are determined. The function Output\_GetNode is used for this purpose. The important argument of the function is River Station, inserted as a string, lRSobs[i]. The numbers of the "observed" nodes are stored as NumPy array of integers, iRSobs.

The computational nodes, RScmp, and their types, RStyp, are also read from the results of the preliminary computations.

The optimization process with Nelder-Mead simplex starts at line 57. A detailed description of this function may be found at SciPy.org [28]. Before this the initial guess, x0, is prepared as a NumPy array. The main argument of the optimization algorithm is the objective function, ObjFun. It is defined as a separate function (lines 6–22).

The basic argument of the objective function is the set of searched parameters x. It is a 3-element NumPy array. At the very beginning, the parameters x are assigned to the variables, representing roughness coefficients nLOB, nCh, nROB. The roughness coefficients are set for any cross-section of the HEC-RAS project in a loop over nodes with the function Geometry\_SetMann\_LChR. The subsequent steps are running of the computations (line 15) and reading the results (line 19) in the loop over nodes (lines 17–20). These processes are performed with the HECRASController functions Compute\_CurrentPlan and Output\_NodeOutput. Both of them were used in the previous example. The difference between the computed, wse, and "observed", WSEobs[i], water surfaces in a single node is calculated in line 20. The final value of the objective function is determined beyond the loop, in line 21. When the optimization process stops, the HEC-RAS window is closed and the controller variable is deleted.

The results obtained in Example 2 are shown in Figure 5. Part (a) presents the convergence of the channel roughness during the run. Part (b) shows changes of the objective function (5) during the computations.

**Figure 5.** Results of Example 2: (**a**) convergence of channel roughness coefficient, (**b**) convergence of objective function.

### *3.4. Example 3—Control of Sediment Simulation*

The purpose of this example is to show that sediment simulation in HEC-RAS may be run several times for different sediment samples. Model B described above is used in this case. The samples are changed automatically by the Python code. The HEC-RAS sediment data file has the XML format. Hence, the xml.etree.ElementTree module [36] is necessary to read and handle data. The results are then read from HDF files with the methods available in the h5py module [59]. The code in this example is very simple. However, the methods presented here may be extended and applied in a more sophisticated code, e.g., calibration of the sediment routing algorithm, analysis of sediment simulation uncertainty, or assessment of flow regime changes caused by sediments.

The third example is too long to present all lines of code in one script. Hence, it is split into Scripts 7 and 8. There are four files with data and results managed in this script. The first is the file storing the tested sediment samples. The format of this file is arbitrary and it will not be discussed here. The only requirement is that the loaded sediment samples should have a form similar to those stored in HEC-RAS. Three other files are more specific. They are the XML file with HEC-RAS sediment data, the HDF file with geometric data, and the HDF file with the results.

```
1. import win32com.client 
2. import os, h5py, numpy 
3. import xml.etree.ElementTree as ET 
4. from support import sc3_LoadPrb,sc3_SaveRes 
5. 
6. # loading samples 
7. pnames,sampA,sampB,sampC,sampD = \ 
8. sc3_LoadPrb('test_sedi','samples','samples.txt') 
9. 
10. # init HEC-RAS Controller 
11. hec = win32com.client.Dispatch('RAS503.HECRASController') 
12. projekt = os.path.join(os.getcwd(),r'test_sedi\test_sedi.prj') # project name
13. hec.ShowRas() # show HEC-RAS window 
14. # loading River Stations from geometry 
15. hdfname = os.path.join(os.getcwd(),r'test_sedi\test_sedi.g03.hdf') 
16. dane = h5py.File(hdfname,'r') # opening of HDF file 
17. # link to the River Stations in HDF file 
18. ListRS = dane.get('Geometry').get('Cross Sections').get('River Stations') 
19. GeomRS = numpy.empty([len(ListRS)],dtype='S15') # NumPy array of RSes 
20. for i in range(0,len(ListRS)): 
21. GeomRS[i] = ListRS[i] 
22. dane.close() 
23. NXS = len(GeomRS) # number of cross-sections
```
**Script 7.** Example 3. Automatic control of sediment transport simulation—part 1.

An example of the sediment data file is shown in Figure 6a. It is presented in the XML Viewer. It has a structure similar to a tree. The groups include subgroups and so on. A description of the XML format can be found in the document prepared by the Python Software Foundation [36]. The group with bed gradation is opened in the XML Viewer. The name of the file and the full hierarchy to access these data are shown in Table 1. Particular gradations are stored as the group attribute. This attribute is a dictionary, where the class symbols 'GC1', 'GC2', etc. reflect the gradation in the format '%Finer'. Every XML file is an ASCII file, and it may be opened in Windows Notepad.

The HDF files may be opened in HDFViewer, as presented in Figure 6b. These are simulation results files in binary format. It is not possible to use Notepad to open them, but the HDF also has a structure similar to a tree with successive branches. Two elements are read from the opened HDF file. These are (a) the invert elevations and (b) the water surface elevations. The first is shown in Figure 6b. The hierarchy to access these data is also presented in Table 1.

```
24. 
25. # test of loaded samples 
26. sname = os.path.join(os.getcwd(),r'test_sedi\test_sedi.s01') 
27. ii = 0 
28. for sample in [sampA, sampB, sampC, sampD]: 
29. plik = ET.parse(sname) # opening XML file with sediment data 
30. dane = plik.getroot() # access to bed gradations 
31. grad = dane.find('Bed_Gradation').find('Sample').find('Gradation') 
32. for elm in sample: # assignment of tested sample 
33. grad.attrib[elm] = sample[elm] 
34. plik.write(os.path.join(os.getcwd(),r'test_sedi\test_sedi.s01')) 
35. del plik 
36. hec.Project_Open(projekt) # open HEC-RAS project 
37. test1 = hec.Plan_SetCurrent('plan sedi00') # setting sediment plan 
38. NMsg,TabMsg,block = None,None,True # no. and list of messages, blocking 
39. # computations of the current plan 
40. v1,NMsg,TabMsg,v2 = hec.Compute_CurrentPlan(NMsg,TabMsg,block) 
41. del NMsg, TabMsg, block 
42. hec.Project_Close() # close HEC-RAS project 
43. # full path to the HDF file results of the computations 
44. hdfname = os.path.join(os.getcwd(),r'test_sedi\test_sedi.p03.hdf') 
45. dane = h5py.File(hdfname,'r') # opening of HDF file 
46. # accessing bed elevations 
47. XSbed = dane.get('Results').get('Sediment').get('Output Blocks') 
48. XSbed = XSbed.get('Sediment').get('Sediment Time Series') 
49. XSbed = XSbed.get('Cross Sections').get('Invert Elevation') 
50. NTS = len(XSbed) 
51. # accessing WSE 
52. XSwse = dane.get('Results').get('Sediment').get('Output Blocks') 
53. XSwse = XSwse.get('Sediment').get('Sediment Time Series') 
54. XSwse = XSwse.get('Cross Sections').get('Water Surface') 
55. # NumPy arrays for bed elevations and WSE 
56. InitBed = numpy.empty([NXS],dtype=float) # initial bed and WSE 
57. InitWSE = numpy.empty([NXS],dtype=float) 
58. LastBed = numpy.empty([NXS],dtype=float) # final bed and WSE 
59. LastWSE = numpy.empty([NXS],dtype=float) 
60. for i in range(0,NXS): # assigment of bed and WSE form HDF 
61. InitBed[i] = XSbed[0][i] # to NumPy array 
62. InitWSE[i] = XSwse[0][i] 
63. LastBed[i] = XSbed[NTS-1][i] 
64. LastWSE[i] = XSwse[NTS-1][i] 
65. dane.close() 
66. del dane 
67. sc3_SaveRes(pnames[ii],sample,GeomRS,InitBed,InitWSE,LastBed,LastWSE) 
68. ii += 1 
69. 
70. hec.QuitRas() # close HEC-RAS 
71. del hec # delete HEC-RAS controller
```
**Script 8.** Example 3. Automatic control of sediment transport simulation—part 2.

**Table 1.** Loading data and results from XML and HDF files: type of data, type of file, access hierarchy.


**Figure 6.** View of data and result files used in Example 3: (**a**) Sediment data file in XML Viewer, (**b**) sediment results in HDFViewer.

The flow chart of the main computations in this example is shown in Figure 7. The main element is a loop over sediment samples. Before the loop starts, the HECRASController object has to be initialized and the name of the project has to be loaded. If the HEC-RAS variable is ready, the information about computational nodes may be read from the HDF5 file storing geometry data. Before the loop starts, the sediment samples have to be loaded to the script. The first step in the loop is access to the XML file with sediment data. The grain fractions are written and the XML file is closed. Then the project can be opened, because any change of the sediment data in the XML file is noticed by the HEC-RAS project

only in the phase of the data loading. Further changes of the data in files stored on the disk do not affect the project. Then the computations start and the results are read from the HDF5 file. At the end of each step, the bed elevations and positions of the computational nodes are written to the separate ASCII file, which enables preparation of the bed profiles.

**Figure 7.** Main computational elements of Example 3.

The program presented in Scripts 7 and 8 starts with import of necessary modules (lines 1–4). The modules win32com.client, os, numpy are also applied in the previous examples. The new elements are xml.etree.ElementTree and h5py. The first one is loaded as ET, and it is used for handling the XML documents. The second one is used for reading of data and results from HDF files. Two functions are loaded from the support module. They are sc3\_LoadPrb for loading sediment samples, and sc3\_SaveRes for saving results in ASCII files. For the sake of simplicity, no local functions are defined.

The first command is reading of the sediment samples for tests (lines 7–8). The samples are stored as dictionaries sampA, sampB, sampC, and sampD. Their items are composed of sediment class denotation in the format 'GC#' and weighted gradation as '%Finer'. The keys of the dictionary should be the same as those used in the HEC-RAS and XML files with sediment data. The first element of the return list, pnames, is the table with samples' names.

The HECRASController starts at line 11. The HEC-RAS window is opened (line 13) and the full path of the project file is set, but the project is not loaded. Every change in the data requires reloading of the project. Hence, the project is loaded later.

The HDF file including geometry data is opened first (line 16). The link to the River Stations is created (line 18) with the get method of the h5py module. Then the NumPy array for storing these data is prepared (line 19). The data are assigned as value-by-value in a loop over nodes. Then the first HDF file is closed (line 22). The important element taken from the River Station array is the number of cross-sections NXS (line 23). The rest of the commands are shown in Script 8.

The longest part of the script is the loop over tested samples (lines 28–68). At the beginning of each step, the XML file with sediment data is opened (line 29). Then the link to bed gradation grad is created. In the internal loop (lines 32–33) over items of the sample, the bed gradations in the XML file are replaced with bed gradations of the current sample (line 33). Then the file is saved (line 34) and closed by deleting the variable representing this file. At this moment the XML file with sediment

data is ready. The HEC-RAS project is opened (line 36) and the computational plan is set properly (line 37). After the current plan is computed (line 40), the project is closed (line 42). The last file of the current iteration is accessed (line 45). This is the HDF file with results of the computations (Table 1). Two elements are read from this file: Bed elevations and water surface elevations. These results are stored as list of successive time steps. Each time step is a list of values stored for each cross-section. The links to the bed and water surface elevations, XSbed and XSwse, are created in lines 47–49 and 52–54. The number of time steps is assigned to the variable NTS.

The NumPy arrays for storing the initial and final bed and surface elevations are created in lines 56–59. The results are assigned to them in another internal loop (lines 60–64) over cross-sections. Later, the HDF file is closed (lines 65–66). Finally, the results are saved in a file specific for each iteration, whose name is the name of the sample. At the end of the script, HEC-RAS is closed and the variable hec is deleted.

The sediment simulation running from the Python script is shown in Figure 8. The display of typical output for such simulation is also shown in the screenshot. The sieve curves of tested sediment samples are plotted in Figure 9a. The variability of assumed sieve curves is typical for a lowland river. The differences in the content of bed sediments are relatively small. However, more complex examples may also be tested in the same way. Figure 9b shows the results of the simulation as bed changes in the selected cross-sections. Considering the time horizon of the simulation, which is 4 months, the bed changes from −20 cm to +60 cm seem to be quite significant. The reason is the relatively huge discharge of 269 m3/s, kept constant during the whole simulation. The differences between the changes of bed elevations vary from a few to 10 cm. If we compare this result with the range of changes, which is about 80 cm, it transpires that the differences between consecutive simulations are up to 12.5% for a very short simulation period. This stresses the importance of the sediment sample variability for the sensitivity of the sediment transport model.


**Figure 8.** Run of sediment simulation and typical display of sediment output.

**Figure 9.** Sediment simulations: (**a**) Sieve curves of tested samples, (**b**) bed changes in selected River Stations.

### **4. Conclusions**

The examples presented in this paper explain the techniques necessary to control HEC-RAS computations with Python. The main element of the presented methods is the application of Python module Pywin32 for handling COM libraries. This module was used in all three examples. The comparison of the code developed for Example 1, with basic examples discussed by Goodell [20,21] and Leon and Goodell [25], illustrates the ease of Python application for handling HECRASController. The use of Python is not more difficult than VBA or MATLAB. Additionally, Example 1 presents a more general approach for application of Python for control of HEC-RAS computations than described in previous papers [40,41]. This example explains how to control any version of HEC-RAS, including the future developments. It also opens doors for further extension of HECRASController capabilities. Possible further extensions were shown in Examples 2 and 3.

The second example presented the application of HECRASController and Python specific modules, for the automatic calibration of the HEC-RAS model. The roughness coefficients were determined by the optimization method from the SciPy module. This approach was similar to that presented by Leon and Goodell [25], where the Genetic Algorithm from the MATLAB toolbox was applied for control of gate operation. Example 2 was also a more advanced approach for determination of roughness coefficients than that presented in the AHYDRA package [24]. Further development of the techniques presented in Example 2, could focus on application of a faster and/or more sophisticated optimization method than the Nelder-Mead simplex. It is also possible to develop a more general Python module for automatic calibration of HEC-RAS, on the basis of the concept presented here. Example 2 may also be improved by preparation of a user optimization algorithm. It could be implemented in Fortran, which is very fast and aimed at numerical computations, and then it could be translated to Python with the F2PY or ctypes modules.

Example 2 presented the linkage of HEC-RAS computations with numerical routines available in specific Python modules, NumPy and SciPy. The number of useful routines in SciPy is greater than those presented. Other extensions and HEC-RAS applications are also possible. The elements available in NumPy are also very useful for storing and handling of computational results. These features of Python are undoubted advantages of Python over standard VBA.

Example 3 illustrated an extension of HECRASController capabilities to control sediment simulation. To the best of the author's knowledge, such an approach or a similar one has not been tested yet. These examples presented a real advantage of Python, owing to the access to very complex data formats, XML and HDF5. The Python modules enabling such access are ElementTree and H5PY, respectively. There are several applications of the presented approach. The most basic is sensitivity analysis of sediment transport simulation. Another possible use is calibration of the sediment model on the basis of the historical data. These two application areas are discussed in the text. Other applications are also possible in the area of water management. Examples may be the

impact of sediment deposition on flood hazard, prediction of morphodynamic changes in the river after regulation, etc. The uncertainty of sediment transport modeling may also be estimated with the techniques used in Example 3.

The illustrated capabilities of Python are not the only advantages of this scripting language over specific VBA applications. Another is access to geoprocessing tools, such as ArcGIS and QGIS. It makes integration of Python with HEC-RAS more necessary. Not only are VBA capabilities too poor to compete with Python in this area; more dedicated numerical software such as MATLAB does not have such broad capabilities for geoprocessing.

Application of Python for control of hydrodynamic simulations with HEC-RAS opens new opportunities in such areas as water management, river engineering, and flood risk management. Good opportunities are available by access to more different data formats from simple text files to more complex XML and HDF5. A huge advantage is the linking of HEC-RAS with other software, such as numerical routines or geoprocessing. More advanced future applications are also possible.

**Software Availability:** https://sourceforge.net/projects/hec-ras-and-python/.

**Author Contributions:** Everything, T. Dysarz.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The author declares no conflict of interest.

### **References**


© 2018 by the author. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article*
