#### Abstract

By focusing on wind-rain two-way coupling algorithm, simulation iterations of wind field and raindrops in the world highest cooling tower (210m) in northwest China were carried out using continuous phase and discrete phase models based on CFD numerical simulation. Firstly, influence laws of 9 wind velocity-rainfall intensity combinations on wind-induced rainfall, raindrop additional force, and equivalent pressure coefficient on internal and external surface of the tower body were discussed. On this basis, speed flow line, turbulence energy strength, raindrop running speed, and track on the tower body in the wind-rain coupling field were disclosed. Finally, qualitative and quantitative contrastive analyses on wind pressure, rain pressure, and equivalent pressure coefficient on internal and external surfaces of the tower body were conducted under different working conditions. Thus, the most unfavorable wind-rain combination was identified. Calculation formulas of equivalent internal and external pressure coefficients of super-large cooling towers were fitted from nonlinear least square method. Research results demonstrate that the 3D effect of equivalent internal and external pressure coefficients with considerations to wind-rain two-way coupling is more prominent. Particularly, there is strong transition on the windward region of the external surface and leeside region at bottom of internal surface. The quantity of caught raindrops on the structural surface is negatively related to wind velocity but is positively related to rainfall intensity. Rain load and rainfall coefficients on the external surface are significantly higher than those on the internal surface. Equivalent internal pressure coefficient has a sharp reduction on the leeside region under different working conditions. Besides, equivalent internal pressure coefficient of different meridians decreases with the increase of height. The maximum and minimum are -0.574 and -0.282, respectively. The proposed equivalent internal and external pressure coefficients of super-large cooling tower can predict wind load under extreme climate conditions accurately.

#### 1. Introduction

With universal uses of high-capacity and high-parameter power units, cooling tower which is one of the important structures of fire/nuclear power plants is developing toward (super-) large types [1, 2]. Compared to conventional cooling tower, super-large cooling tower is a typical wind sensitive structure. Wind load is the control load of internal force design of tower body. Existing structural specifications of cooling towers all neglect impacts of rainfall on tower body. Nevertheless, structures have to bear combined effect of strong wind and rainstorm under extreme climatic conditions. Due to collaborative driving by wind force and gravity, the movement track of raindrops will deviate. Some will impact on external surface of the cooling tower strongly at the leading edge, while some will impact on inner wall strongly through the open tower top, resulting in significant changes of aerodynamic force distribution on surface. Meanwhile, rainstorm will deteriorate turbulence effect of turbulence wind where there is no rain. Under this circumstance, airflows close to tower body and in the tower are complicated, which changes movement track, additional force, and internal and external pressure of raindrops. For this reason, it has important theoretical significance and engineering value to study wind load distribution characteristics on internal and external surfaces of super-large cooling tower under different wind velocities and rainfall intensities as well as the corresponding mechanism of action.

Researches concerning wind load of super-large cooling towers mainly focus on stochastic behaviors of single wind load [3], extreme wind pressure distribution [4], static interference effects [5], and dynamic interference effect [6]. Chinese code [7] offers the mean wind pressure coefficient distribution curve on external surface of hyperbolic cooling tower and internal pressure coefficient of single tower (-0.5), but it neglects variations of internal attractive force along height and meridians completely. Germany code [8] gives the mean pressure coefficient distribution pattern on external surface according to external surface roughness but still overlooks effects of 3D characteristics and rainfall on internal pressure coefficient. However, Shen et al. [9] simulated wind load on internal surface of cooling towers under cold-hot air self-circulating system and external wind field by using the realizable k-*ε* turbulence model and multiphase flow model. They found that internal pressure which was produced in external wind field changes dramatically with height and attitude. Zou et al. [10] analyzed 3D effects of wind pressure distribution in internal surface of tower body under different ventilation rates (100% and 30%) by rigid pressure test, finding nonuniform distribution of wind pressure on internal surface along height and circumferential directions. Associated results could solve values of wind pressures on internal and external surfaces of conventional cooling towers. Nevertheless, there are hardly researches about aerodynamic performances on internal and external surfaces of super-large cooling towers under wind-rain combined effect, not to mention qualitative and quantitative comparisons of internal and external pressures under different wind velocity-rainfall intensity combinations.

In addition, some researches [11–19] concerning wind-rain combined effect mainly focus on low-rise buildings, bridge cable, power transmission tower, and wind turbine. They only consider effects of wind-induced rainfall but overlook reaction of rain to wind. Currently, there are rare researches that have discussed whether wind-rain coupling affects aerodynamic performances of internal and external surface of super-large cooling towers under poor climate conditions.

In this study, by focusing on wind-rain two-way coupling algorithm, simulation iterations of wind field and raindrops in the world highest cooling tower (210m) in northwest China were carried out using continuous phase and discrete phase models based on CFD numerical simulation. On this basis, mechanism of action of internal and external pressures of the tower body under wind-rain combined effect was studied. Influence laws of different wind velocity-rainfall intensity combinations on wind-induced rainfall, raindrop additional force, and equivalent pressure coefficient on internal and external surface of the tower body were discussed. Finally, the most unfavorable wind-rain combination was identified. The calculation formulas of equivalent internal and external pressure coefficients of super-large cooling towers were fitted based on nonlinear least square method.

#### 2. Wind-Rain Two-Way Coupling Algorithm

##### 2.1. Rainfall Intensity

Classifications of rainfall intensity are shown in Table 1 [20–22]. Two rainfall intensity levels have different sampling times. Measurement results of the same rainfall differ significantly. Structure checking based on hourly rainfall intensity is safer than that based on daily rain intensity. Hourly rainfall intensity can intuitively reflect influences of instantaneous rainfall intensity on structural performances under extreme climate conditions which are mostly concerned in engineering.

In addition, the main purpose of this paper is to compare the influence of different wind speed-rain intensity combinations on wind load characteristics and action mechanism on internal and external surfaces of super-large cooling towers. Considering the possibility that wind speed and rainfall both reach a certain intensity under wind-rain combined effects, this paper chooses hourly rainfall intensity, rather than the rainfall intensity for minor durations with a very low probability of occurrence, and adopts the data based on the downpour grades under extreme conditions which may occur in China.

##### 2.2. Raindrop Spectral Distribution

Rainfall spectrum approximately obeys negative exponential distribution. Common models [23] have Best Spectrum, Marshall-Palmer Spectrum, and Gamma Spectrum. In this study, Marshall-Palmer Spectrum was used:where* D*_{p} is the diameter of raindrop (mm),* n*(*D*_{p}) is concentration spectrum of raindrops with different diameters,* N*_{0} is concentration (=8,000), and *λ* is a scale parameter.* λ *can be expressed aswhere

*R*is rainfall intensity (mm/h).

##### 2.3. Final Velocity of Raindrops

When the falling height is over 20m, raindrops of almost all diameters can reach the final velocity [24]. Raindrops of diameter smaller than 2mm fall approximately as balls, but larger raindrops will suffer large deformations due to air resistance. The vertical final velocity of raindrops can be calculated by the following equation [25, 26]:where* v*(*D*_{p}) is the vertical final velocity of raindrops with a diameter of* D*_{p} (m/s).

##### 2.4. Wind-Rain Two-Way Coupling Algorithm

During downpours, the volume fraction of raindrops in air is far smaller than 10%. Therefore, raindrop can be simulated by DPM model. Simulated raindrops are added to the continuous phase as the second phase after stabilization of wind field computing for wind-rain two-way coupling [27, 28]. The motion equilibrium equation of raindrops in wind field iswhere is discrete phase particle velocity, is continuous-phase fluid velocity, (-) is the drag force of unit particle mass, *ρ*_{p} and *ρ* are particle and fluid densities, and is interaction force between discrete phase and continuous phase.where *μ* is viscosity coefficient of fluid,* D*_{p} is particle diameter, and* Re* is relative Reynolds number.

After discrete phase of raindrops is considered, the basic control equation of wind continuous phase iswhere* S*_{m} is the second discrete phase mass which is added to the continuous phase,* p* is pressure, is tensor of stress, and is gravity. Among them, tensor of stress can be expressed aswhere* I* is unit tensor and the second term in the right is volume expansion.

##### 2.5. Wall Collision Equation

The impacting process of raindrops onto tower walls complies with the law of conservation of momentum. The key to solve impact force lies in impact time. Evaporation, splash, and rupture of raindrops during the impacting process are neglected. It is believed the raindrop-structure interaction follows Newton’s second law. According to theorem of momentum,where* f*(*t*) is impact force vector of single raindrop (N) and* v* is velocity vector of raindrop.

The impact force of raindrop onto the structure in unit time (*F*(*τ*)) is

If the falling raindrops are viewed as balls approximately, it has

Since raindrop diameter is generally lower than 6mm and the horizontal final velocity before the impact is relatively high, for calculation simplification, the impact time *τ* is determined:

Then, the impact force of raindrops to structure can be simplified as

#### 3. Brief Introduction to the Project and Working Condition Setting

##### 3.1. Brief Introduction to the Project

The super-large reinforced concrete hyperbola indirect air-cooling tower with natural ventilation built is 210m high. Elevation of the throat is 157.5m and elevation of the inlet mouth is 32.5m. Inner diameters of the throat and inlet mouth are 110m and 159m, respectively. The zero diameter is 180m. The tower body employs exponential thickening. The minimum thickness (0.37m) lies in the throat section and the maximum thickness (2.0m) is at the lower ring beam. The tower body is connected with circumferential plate foundation by 52 pairs of X-shaped columns. X-shaped columns have rectangular section (1.8m×1.2m). The circumferential foundation is a cast-in-place reinforced concrete structure: 12.0m (width)×2.5m (height). Main structural size and picture of the cooling tower are shown in Table 2.

##### 3.2. Working Conditions

This cooling tower locates in the B-type landform. With considerations to conventional working state of louver of cooling tower, the open effect of louver was considered according to 30% of ventilation rate [1, 2]. Influences of combinations of three wind velocities and three rain intensities on aerodynamic performance of the cooling tower were discussed. Weak wind, moderate wind, and strong wind were divided according to the maximum wind velocities in the return period of 10 years, 50 years, and 100 years. Rainfall intensity was classified into weak downpour, moderate downpour, and heavy downpour. Therefore, a total of 9 working conditions could be set (Figure 1).

#### 4. Numerical Simulation of Wind-Rain Two-Way Coupling

##### 4.1. Establish the Wind-Rain Field Model

To ensure the cooling tower in the rainfall region and full development of wake flow, dimension of calculation domain was set 3000m (wind length)1500m (wind width)600m (height). The bottom tower center was chosen as the origin of the coordinate system and the following wind direction is consistent with X-axis. To consider both computing efficiency and accuracy, the cooling tower was divided into local and periphery wind-rain fields during gridding. The local wind-rain field contained the cooling tower model and was meshed by unstructured grid. The periphery wind-rain field had a regular shape and was meshed by high-quality structured grid. The minimum grid size in the core region was 0.2m. The total number of grids in the overall model was over 18 millions and grid mass was higher than 0.40 (it is required to be higher than 0.1 in order to prevent negative volume). Both number and quality of grids can meet the calculation requirements. Overall gridding and local gridding are shown in Figure 2.

**(a) Overall gridding**

**(b) Local gridding**

Inlet of the calculation domain is the boundary of velocity inlet and the outlet is boundary of pressure outlet. Two side walls and top surface use symmetry boundaries and the cooling tower and ground adopt wall boundary. The overlapping face of local and periphery calculation domains are set the interface. The calculation domain of the wind-rain field and its boundary conditions are shown in Figure 3.

##### 4.2. Wind-Rain Field Coupling Calculation

Wind-rain two-way coupling simulation theory and discrete phase tracking method are applied to this super-large Reynolds number structure, which proposes extremely high requirements on memory capacity of computer. Therefore, numerical calculation in this study is accomplished based on large compute server in the wind turbine aerodynamics high-performance computing center of the project group (Figure 4). The server processer used Intel(R) Xeon(R) CPU E5-2650 v3 @ 2.30GHz 2.30GHz (2 processor) and the installed memory was 256GB. The 64-bit operating system was employed.

**(a) Calculation center**

**(b) Calculation server**

Numerical calculation used 3D single precision and separated solver. Flow velocity of fluid field was used as the absolute velocity. Air model was equivalent to the ideal incompressible fluid and the turbulence model chose the k-*ω* shear-stress (SST) control equation. Inlet of the calculation domain used the wind profile model with an exponential index of 0.15. Wind velocities at 10m over the ground were set as three standard wind velocities in Section 3.2. The flow field was solved through coupling between seed and pressure by using the SIMPLEC algorithm. The solving formulation of convective term was the second order. During the calculation process, grid slope correction was set in order to increase hybrid grid computing effect. Calculation residue of the control equation was set as 1×1. Finally, the wind field was initialized for iterative computations. Simulation curves and theoretical curves of mean wind velocity and turbulence profile are shown in Figure 5. Results indicate that simulated results of mean wind velocity and turbulence profile conform to theoretical values, indicating that simulation standards of wind field meet engineering requirements.

**(a) v0=20m/s**

**(b) v0=23.7m/s**

**(c) v0=25.3m/s**

After stabilization of wind field solving, the discrete phase was added for iterative operation of wind-rain field coupling. Here, raindrops with six different diameters in the range of 1.0~6.0mm were used to simulate rainfall of continuous diameter distribution (Table 3). Number and volume occupancy of raindrops with different diameters were determined by Marshall-Palmer spectrum in Section 2.2. Raindrops were “plane” released at the horizontal velocity of 0. Driven by wind, they kept consistent movement velocity with horizontal wind velocity in the location. The vertical release velocity was -5m/s. The collaborative effects of gravity and resistance make raindrops reach the final velocity in (3) to the maximum height range.

##### 4.3. Validity Verification

The mean wind pressure coefficients in the throat area on external surface of the cooling tower under three wind velocities were compared with existing Chinese and foreign relevant codes [7, 8] and measurement curves [29] (Figure 6). According to analysis, the angle between negative pressure extreme point and separation point on the mean wind pressure distribution curve in the throat area conforms to that of measurement curve under all three wind velocities. The negative pressure value in the leeside region is slightly higher than fire engineering codes and measured values but is lower than Germany code. Other numerical values of angles basically agree with codes and measurement results. To sum up, numerical simulation results are valid to some extent.

#### 5. Contrastive Analysis of Results

##### 5.1. Analysis of Wind Field

Three-dimensional wind velocity flow lines and vorticity distribution of the cooling tower under three standard wind velocities before raindrops are added are shown in Figures 7 and 8, respectively. Clearly, the incoming flow runs through the tower body and splits in the windward side. It is accelerated along two external side walls of the tower body and forms different sizes of vortexes on the leeside surface. Some airflow enters into the tower body through louver. Meanwhile, inflows at tower top skim over and change partial rising airflow directions. Airflow generates reflux and forms a large-sized vortex at the throat area. With the increase of wind speed, airflow runs at an accelerated velocity, which intensifies vortex falloff, dense distribution of wind velocity flow lines, and interaction between rising airflows and incoming flows at the tower top. Turbulence energy intensity is positively related to wind velocity, especially in middle and lower parts of the leeside surface. Peak mainly lies in top of windward region and throat area of the leeside region.

**(a) v0=20m/s**

**(b) v0=23.7m/s**

**(c) v0=25.3m/s**

**(a) v0=20m/s**

**(b) v0=23.7m/s**

**(c) v0=25.3m/s**

Three-dimensional distribution patterns of external and internal pressure coefficients under 9 working conditions are shown in Figures 9 and 10.

**(a) Working condition 1**

**(b) Working condition 2**

**(c) Working condition 3**

**(d) Working condition 4**

**(e) Working condition 5**

**(f) Working condition 6**

**(g) Working condition 7**

**(h) Working condition 8**

**(i) Working condition 9**

**(a) Working condition 1**

**(b) Working condition 2**

**(c) Working condition 3**

**(d) Working condition 4**

**(e) Working condition 5**

**(f) Working condition 6**

**(g) Working condition 7**

**(h) Working condition 8**

**(i) Working condition 9**

Distribution pattern of external and internal pressure coefficients under all working conditions is basically consistent with numerical values. They are all highly symmetric around the wind axis.

External pressure coefficients on different cross-sections of tower body under all working conditions all present a V-shaped variation trend with the expansion of circumferential angle and become stable close to the leeside region. Numerical values of external pressure coefficient are relatively stable close to windward region and leeside region except for gentle fluctuation in the crosswind region.

Influenced by airflow which runs through louver, the internal pressure coefficient at tower bottom jumps strongly along circumferential direction. The numerical value at the leeside surface is significantly lower than those in other regions. With the gradual growth of internal airflow from tower bottom, the circumferential random impact of tower body from middle and upper airflows is relatively weak. Therefore, the circumferential numerical value of internal pressure coefficient at middle and upper parts of the tower body basically remains constant.

##### 5.2. Analysis of Rain Field

Raindrop movement tracking based on resultant velocity of particles in wind-rain field under 9 working conditions is shown in Figure 11. Raindrop density was processed by equal proportional coarsening. The following can be seen from Figure 11:

**(a) Working condition 1**

**(b) Working condition 2**

**(c) Working condition 3**

**(d) Working condition 4**

**(e) Working condition 5**

**(f) Working condition 6**

**(g) Working condition 7**

**(h) Working condition 8**

**(i) Working condition 9**

Raindrops on near wall of the cooling tower show complicated movement patterns. Some raindrops impact on the windward surface of the cooling tower under wind actions, while the remaining raindrops separate at two structural sides along the airflow. Only few adhere onto side walls of the cooling tower and most enter into the wake flow with winds.

Raindrops at upper front end of the cooling tower run through the outlet of the cooling tower and enter into the inside region and impact on middle and upper leeside regions of internal surface at a great velocity. Moreover, the stronger rainfall intensity is, the more raindrops caught by the internal wall will be.

With the increase of wind velocity, horizontal acting force of raindrops increases dramatically, which makes raindrops move along the wind direction at accelerating speeds. Tremendous raindrops skim over the tower top and move to rear of the tower body. As a result, number of raindrops which enter into the cooling tower drop sharply.

Impacts of number of raindrops with different diameters, velocity, and occupancy of speed on external and internal surfaces of the cooling tower under 9 working conditions are shown in Figures 12 and 13.

**(a) Number of raindrops**

**(b) Occupancy of speed**

**(c) Horizontal velocity**

**(a) Number of raindrops**

**(b) Occupancy of speed**

**(c) Horizontal velocity**

Diameter of raindrops caught by internal and external surfaces of the cooling tower mainly ranges within 3~6mm. Raindrops with 5mm diameter occupy the highest proportion. This is because, driven by the same wind intensity, movement speed of smaller raindrops increases quickly. Before reaching the height range of the cooling tower, raindrops have skimmed over the tower body in horizontal direction and entered into the wake flow region.

Number of raindrops caught by internal and external surfaces decreases with the increase of wind velocity. With the increase of rainfall intensity, number of raindrops on external surface is far higher than (about 10 times) that on the internal surface.

Occupancy of horizontal speed of raindrops which are caught by external surface shows basically the same distribution laws. It increases firstly and then decreases with the increase of final horizontal velocity. The mean horizontal velocity of raindrops which are caught by internal and external surfaces is significantly lower than the standard wind velocity. Moreover, the horizontal impact velocity of raindrops is positively related to diameter of raindrops.

Rain loads on internal and external surfaces of the cooling tower under 9 working conditions are calculated from (14). Meanwhile, the rain-wind load ratios at different heights of internal and external surfaces of the cooling tower are given. Results are shown in Figure 14 and Tables 4–6.

Under all working conditions, rain load on external surface of the cooling tower decreases firstly and then increases with the increase of height, reaching the valley at 0.23~0.31H and the peak at tower bottom or tower top.

Under all working conditions, rain load on internal surface of the cooling tower increases with the increase of height. There is no raindrop at places below 0.69H and most rain loads are concentrated in the height range of 0.90~1.0H. The ratio is about 95%.

Rain-wind load ratio at different heights is very small. The highest rain load is only 3.655 of wind load, which occurs in the height range of 0.90~1.0H under working condition 3. Besides, rain load on external surface is significantly higher than that on internal surface under all working conditions.

Rain load in internal and external surface under different wind velocities climbs up significantly with the increase of rainfall intensity. Given the fixed rainfall intensity, rain loads on internal surface as well as middle and lower parts of external surface decline as wind velocity increases. In middle and lower parts of the cooling tower, number of raindrops influences rain load more significantly than wind velocity. However, accelerating movement of raindrops on middle and upper parts of the cooling tower increases rain load dramatically.

To display location and number of raindrops at typical heights and angles as well as corresponding pressure coefficients more clearly, 3D distributions of raindrops and rain-induced pressure coefficients on internal and external surfaces under all working conditions are shown in Figures 15 and 16.

**(a) Working condition 1**

**(b) Working condition 2**

**(c) Working condition 3**

**(d) Working condition 4**

**(e) Working condition 5**

**(f) Working condition 6**

**(g) Working condition 7**

**(h) Working condition 8**

**(i) Working condition 9**

**(a) Working condition 1**

**(b) Working condition 2**

**(c) Working condition 3**

**(d) Working condition 4**

**(e) Working condition 5**

**(f) Working condition 6**

**(g) Working condition 7**

**(h) Working condition 8**

**(i) Working condition 9**

Under all working conditions, raindrop impact positions mainly concentrate on windward region of external surface and leeside region of internal surface. Driven by airflow vortexes, there are few raindrops adhered on leeside region of external surface and windward region of internal surface.

There are more raindrops caught on internal and external surface of the cooling tower under working condition 3 than the remaining working conditions. Number of caught raindrops drops sharply with the increase of wind velocity and increases gradually with the increase of rainfall intensity. Meanwhile, number of caught raindrops on external surface is far higher than that on internal surface.

Under all working conditions, rain-induced external pressure coefficient mainly concentrates on 60° range of two sides of the windward surface, and it is basically 0 in the remaining ranges. The maximum value is 0.184, which is observed in the range of 0.15~0.23H under working condition 3.

Under all working conditions, rain-induced internal pressure coefficient mainly concentrates in the meridian range of 0.9H~1.0H and 90° range of two sides of the leeside region. It is basically 0 in the remaining ranges.

##### 5.3. Analysis of Equivalent Pressure Coefficient

For quantitative comparison of wind and rain-induced pressure distribution on the cooling tower under different working conditions, the equivalent internal and external pressure coefficients are defined. Specific calculation steps are as follows: ① rain loads at different monitoring points of internal and external surfaces are converted into rain-induced pressure. ② Calculate the ratios between rain pressure at monitoring points and wind pressure at the corresponding reference heights, that is, rain-induced internal and external pressure coefficients. ③ Calculate the sum of vectors of rain-induced coefficient and wind pressure coefficient, that is, equivalent internal and external pressure coefficients under wind-rain collaborative effects.

Comparison curves of equivalent internal and external pressure coefficients of 6 internal and external typical sections under different working conditions are exhibited in Figures 17 and 18.

**(a) 0.27H (56m)**

**(b) 0.51H (104m)**

**(c) 0.73H (152m)**

**(d) 0.80H (168m)**

**(e) 0.87H(184m)**

**(f) 0.95H (200m)**

**(a) 0.27H (56m)**

**(b) 0.51H (104m)**

**(c) 0.73H (152m)**

**(d) 0.80H (168m)**

**(e) 0.87H (184m)**

**(f) 0.95H (200m)**

Distribution laws and numerical values of equivalent external pressure coefficients of the section at the same height are consistent under different working conditions except for small differences between windward region and leeside region. The external pressure coefficient decreases firstly, increases, decreases, and finally stabilizes from the windward region to leeside region.

Under the same working condition, there is a small difference in numerical values of equivalent external pressure coefficients on 6 typical sections. However, they are all highly symmetric around the wind axis. The maximum negative pressure presents an inverted V-shaped variation with the increase of height. The negative pressure close to the throat area is about -1.5.

Three-dimensional distribution effect of internal pressure coefficient is remarkable. Under the same working condition, equivalent internal pressure coefficients on 6 typical sections are different, but they are all symmetric around the wind axis. The circumferential distribution laws of equivalent internal pressure coefficients are basically consistent under different working conditions.

Equivalent internal pressure coefficients of all sections decrease in the leeside region of the cooling tower (which is in the blue dotted framework). This is mainly because airflow runs through the louver and impacts on internal surface of the leeside region. Influenced by the incoming wind, raindrops will impact on upper internal surface of the cooling tower. Moreover, raindrops influence internal pressure coefficient more than the incoming wind.

Equivalent internal and external pressure coefficient curves of four typical meridians (0°, 75°, 120°, and 180°) of the tower body are shown in Figures 19 and 20.

**(a) 0°**

**(b) 75°**

**(c) 120°**

**(d) 180°**

**(a) 0°**

**(b) 75°**

**(c) 120°**

**(d) 180°**

The equivalent external pressure coefficient of the 0° meridian is basically stable, while those of the 75° and 180° meridians increase firstly and then decrease with the increase of height. The equivalent external pressure coefficient of the 120° meridian remains basically constant with the increase of height, but the numerical value decreases above the throat area.

Distribution trends of internal pressure coefficients of all four meridians are basically consistent. They all decrease gradually with the increase of height. The maximum numerical value is -0.574, which is at inlet height of 75° meridian. The minimum numerical value is -0.282, which is at the tower top of 180° meridian.

Compared to other meridians, the equivalent internal pressure coefficient of internal surface bottom of the 180° meridian decreases sharply, while the equivalent internal pressure coefficient of the tower top declines slightly. This is mainly because bottom equivalent internal pressure coefficient is influenced by impacts of airflow that runs through the louver on internal wall of the cooling tower, whereas the top equivalent internal pressure coefficient is caused by raindrop impact load on internal surface of the leeside region.

##### 5.4. Two-Dimensional Fitting Surface of Equivalent Pressure Coefficient

For the convenience of engineering designers to get equivalent internal and external pressure coefficients of super-large cooling towers, the calculation formula of equivalent external pressure coefficient is fitted based on the principle of nonlinear least square method by using meridian height and circumferential angle as the objective function. The cooling tower is divided into N_{1} sections uniformly along the circumferential direction and divided into N_{2} sections uniformly along the meridian direction. Let N=N_{1}N_{2}. The specific definition of the formula iswhere ** I** is the N1 matrix which covers only elements valuing 1, is the N1 matrix that circulates by N

_{2}times that using N

_{1}angles as the circulation unit,

**is the N1 matrix that circulates by N**

*Z*_{1}times that using N

_{2}same heights as the circulation unit, is multiplying corresponding elements of the matrix, is dividing corresponding elements of the matrix, is the n power of corresponding matrix elements, exp () is the matrix of exponent that uses each element of the matrix in the return brackets e as the base number, and is the N1 matrix that uses wind pressure coefficient of N

_{1}circumferential angles as the unit and changes by N

_{2}times along the meridian height. a

_{i}, b

_{i}……o

_{i}(

*i*=1, 2,…9) are fitting coefficients. The fitting coefficients under the most unfavorable working condition (condition 3) are shown in Table 7.

The comparison curves between two-dimensional actual values and fitting values of equivalent external pressure coefficients under the most unfavorable working condition (condition 3) are shown in Figure 21. Numerical values of scattered points are actual equivalent external pressure coefficients of the cooling tower. The corresponding numerical values of the hook surface are equivalent external pressure coefficients which are gained by two-dimensional fitting formula. It can be seen from Figure 21 that equivalent external pressure coefficient of the cooling tower has obvious two-dimensional features along the meridian and circumferential directions. Viewed from overall distribution of fitting hook surface, the variation laws along meridian and circumferential directions basically conform to distribution of actual equivalent external pressure coefficient except for slight difference (<10%) in the region of maximum negative pressure. The comparison results demonstrate that the proposed two-dimensional fitting formula can be used as the valuing reference of equivalent external pressure coefficient of super-large cooling towers.

Similarly, the calculation formula of equivalent internal pressure coefficient of super-large cooling towers is fitted:where* Z* is meridian height of the tower body, is circumferential angle of tower body, is fitting value of equivalent internal pressure coefficient that uses* Z* and as the target function, and a_{i}, b_{i}……o_{i} (*i*=1, 2,…9) are fitting coefficients. The fitting coefficients under the most unfavorable working condition (condition 3) are shown in Table 8.

Comparison between two-dimensional actual values and fitting values of equivalent internal pressure coefficients of cooling tower under the most unfavorable working condition (condition 3) is shown in Figure 22. Numerical values of scattering point are actual equivalent internal coefficients of the cooling tower. Corresponding numerical values of the hook surface are equivalent internal pressure coefficients which are fitted by the two-dimensional fitting formula. It can be seen from Figure 22 that equivalent internal pressure coefficient of the cooling tower decreases gradually with the increase of height, but the numerical value in the bottom leeside region decreases significantly. According to overall distribution of fitting hook surface, the variation laws of meridian direction and circumferential direction comply basically with actual equivalent internal pressure coefficients. The numerical values envelop real numerical values and the error rate is smaller than 10%. The comparison reveals that the proposed two-dimensional fitting formula can be used as the valuing reference for equivalent internal pressure coefficient of super-large cooling towers.

#### 6. Conclusions

(i)Wind field and wind pressure coefficients of cooling towers have distinct three-dimensional effect before the rain field is added. Horizontal force of raindrops which are added subsequently increases with the growth of wind velocity. Tremendous raindrops skim over the tower top and separate at two sides of the tower body and finally enter into the wake flow region. Only few raindrops adhere on external surface and enter into the tower and impact on the internal wall with airflow.(ii)Number of raindrops which are caught by internal and external surfaces is negatively correlated with wind velocity but is positively related with rainfall intensity. The number of raindrops caught by internal and external surfaces under working condition 3 (wind velocity=20m/s and rainfall intensity=200mm/h) is higher than those under the remaining working conditions. Moreover, raindrops with 5mm diameter occupy the highest proportion on internal and external surfaces. Horizontal velocity of raindrops on external surface is mainly in the range of 3~10m/s and the horizontal velocity of raindrops on internal surface is mainly in the range of 1~3m/s. Number of caught raindrops on external surface is far higher than that on internal surface.(iii)Under different working conditions, raindrop impact positions mainly concentrate in the range of 60° at two sides of the external windward surface and the range of 90° at two sides of internal leeside region near the top of cooling tower. Rain load and rain-induced pressure coefficient of external surface are significantly higher than those of internal surface. The maximum rain-induced external pressure coefficient is 0.184, which is observed near the bottom of cooling tower under working condition 3. The maximum rain-induced internal pressure coefficient is 0.0038, which is observed near the top of cooling tower under working condition 3.(iv)Under different working conditions, equivalent external pressure coefficient at section of the same height differs slightly in the windward region and leeside region. The maximum negative pressure values and negative pressure values of the leeside region increase firstly and then decrease with the increase of height. The maximum negative pressure in the throat area is about -1.5. Equivalent internal pressure coefficient under different working conditions drops dramatically in leeside region. Equivalent internal pressure coefficient of different meridians decreases gradually with the increase of height. The maximum and minimum are -0.574 and -0.282, respectively.(v)The two-dimensional fitting formulas of equivalent internal and external pressure coefficients which are proposed based on nonlinear least square method can predict load value under similar extreme climate conditions. The maximum error rate is controlled lower than 10%.

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest.

#### Acknowledgments

This project is jointly supported by the National Natural Science Foundation (U1733129; 51761165022; 51208254), Jiangsu Province Outstanding Natural Science Foundation (BK20160083), China Postdoctoral Science Foundation (2013M530255; 1202006B), and Jiangsu Qinglan Project.