Flow and temperature distributions in a disc type winding-Part II: Natural cooling modes

The liquid flow and temperature distributions in oil natural (ON) cooled power transformer windings are of great importance for thermal design, thermal rating and life estimation of ON power transformers. In this paper, both experimental and numerical investigations of liquid flow and temperature distributions in a disc-type winding model are conducted. The influences of representative operating conditions on liquid flow and temperature distributions are studied using Particle Image Velocimetry (PIV) and temperature measuring systems. Reverse liquid flows and consequential localised overheating at the top part of winding passes are recorded in the experiments. 2D and 3D CFD simulations are implemented and compared with the measurements. 3D CFD simulations are required for reliable predictions of temperature distribution when reverse flows occur. Dimensional analysis is adopted to interpret the measurement and CFD simulation results in order to provide simplicity, insight and universally applicable conclusions. The liquid flow and temperature distributions are found to be mainly controlled by the Richardson number (Ri). The higher Ri the more uneven the liquid flow and temperature distributions would be. At low Ri region flow and temperature distributions are relatively uniform with a low hotspot factor. Above the region, flow and temperature distributions become rapidly uneven and the hotspot factor increases steeply with the increase of Ri. Further increase of Ri results in the occurrence of reverse flows and fluctuation of the hotspot factor.


Introduction
The hottest temperature within a power transformer winding, usually referred to as the hotspot temperature, determines the severest thermal ageing rate of the winding paper insulation. As such, studying of the hotspot temperature attracts interests from the wider transformer community.
Liquid thermosiphon-driven naturally cooled transformers, henceforth referred to as ON transformers according to IEC standards [1,2], are widely used in the electric grid as they could potentially offer more reliability and require less maintenance compared to pump-driven forced and directed cooled (OD) transformers. In addition, ON transformers are especially suitable for low ambient temperature conditions, in which the liquid can be too viscous to pump.
Thermal hydraulic network models (THNM) [3][4][5] and computational fluid dynamics (CFD) models [6][7][8][9][10] were developed to predict liquid flow and temperature distributions within ON transformer windings. Since THNM and CFD share the same physical principles and CFD results are more detailed and reliable at the expense of much higher computational requirement, only the findings from CFD modelling are summarized in the following paragraphs.
Hot streaks were reported from CFD modelling as a consequence of high Prandtl number (Pr) of the liquid [11]. The hot streak dynamics lead to strong coupling of liquid flow distributions between winding passes in ON transformers [12].
The liquid flow distribution in ON transformer windings features high flow fractions in pass bottom horizontal ducts and low flow fractions in pass top horizontal ducts, originated from the hydrostatic pressure distribution in the windings [6,13]. In addition, reverse flows were reported in horizontal ducts at the top part of the winding pass [7,12]. 3D liquid flow and temperature distribution patterns were observed in 3D CFD modelling [14,15], indicating the necessity of performing 3D modelling.
Experimentally, an ON transformer model was established to investigate the total liquid flow rate in [16]. It was found that the total liquid flow rate is proportional to the square root of the total power loss in the winding and approximately proportional to the square root of the thermal head, the height difference between the radiator central point and the winding central point [16]. Experimental and 3D CFD investigations of effects of different disc-type ON transformer winding arrangements on total liquid flow rate, liquid flow and temperature distributions in the winding were conducted in [15]. 3D CFD modelling of the complete cooling loop was found inaccurate in determining the total liquid flow rate, and therefore resulted in large differences of liquid flow and temperature distributions between the CFD and the measurement results [15]. The 3D CFD model was then reduced to the winding part with measured winding inlet flow rate and liquid inlet temperature as the boundary conditions to achieve satisfactory comparisons with measurements [15]. In addition, experimental study of the influence of different liquids on the flow and temperature distributions in transformer windings in ON cooling modes was presented in [17].
In this study, liquid flow and temperature distributions in a disctype power transformer winding model in ON cooling modes are investigated systematically by experimental measurements and CFD modelling. It is the first time the influence of operational conditions (power loss non-uniformity, winding pass inlet temperature, inlet velocity and total power loss) on an ON transformer winding was investigated with direct experimental measurements of the liquid flow distribution in the winding horizontal ducts. In addition, dimensional analysis is adopted to generalise the results to provide insight of the ON transformer thermal behaviour and guidelines for ON transformer thermal design review. This is a companion paper of the study on OD cooling modes (Part I). Methodology of the study is provided in Section 2. Liquid flow and temperature distribution results from measurements and CFD models are presented and compared in dimensional forms in Section 3. The results are then generalised in dimensionless forms in Section 4, followed by conclusion in Section 5.

Governing equations and dimensional analysis
As for the OD cooling modes, the fluid flow and heat transfer in ON cooling modes are governed by the conservation laws for mass, momentum and energy. The difference between the two types of cooling modes is that, in ON cooling modes, the flow is driven by buoyancy force and OD flow by pumps. When Boussinesq approximation is adopted, the governing equations for planar 2D problem in nondimensional forms for the fluid domain are shown in equation (1)-(3). The detailed deduction of these equations is given in [7,18]. (2) The Richardson number (Ri), Reynolds number (Re), Grashof number (Gr) and Pr are defined as where the del operator ∇ * incorporates a characteristic length D h to make it dimensionless (∇ = ∇ * D h ), with the characteristic length chosen to be the hydraulic diameter of the pass inlet vertical duct, D h = 4A/L. The coordinates are normalised against the hydraulic diameter, velocity normalised against winding pass inlet velocity and the pressure against twice the dynamic pressure at the winding pass inlet. The dimensionless temperature is defined as . The definition of dimensionless temperature takes the same format as in defining the hotspot factor in IEC standard [1,2]. Therefore, the highest dimensionless temperature is the hotspot factor (H).
hs to aw to bo (8) In ON cooling modes, the fluid flow and heat transfer are coupled by the buoyance force term ( *

Ri T
· ) in the momentum equation. When Ri is much smaller than 1, the buoyance force is negligible and the cooling is  in a forced convection (OD) regime; when Ri is close to or larger than 1, the cooling is in a combined or natural convection (ON) regime. From the above dimensionless governing equations, one can conclude that for a fixed winding geometry in ON cooling modes, the dimensionless flow distribution or volumetric flow fraction distribution and the dimensionless temperature distribution are dictated by Ri, Pr and Re at the winding pass inlet. Ri is influential because the flow is driven by buoyancy force. In addition, Re is more influential than Pr because Re has the same role as Pr in the energy equation but Re affects flow in the momentum equation.

Experimental setup
The experimental setup is the same as the one presented in the OD companion paper (Part I), where more detailed geometric descriptions are provided. The winding model was designed to model a section of a disc-type winding between two sets of spacers/sticks. The winding model was made of 3 passes with 10 plates/11 horizontal ducts per pass. The plates are numbered from bottom to top from 1 to 30. Horizontal cooling duct height was fixed to 4 mm and vertical/axial cooling duct width on both sides was fixed to 10 mm. The plates have a rectangular shape with the dimensions of 100 × 104 × 10 mm, each embedded with two resistive cartridge heaters and two K-type thermocouples calibrated to accuracy of ± 1°C. Liquid flow rates within the third pass horizontal cooling ducts were measured using a Particle Image Velocimetry (PIV) system. The ducts in pass 3 were numbered from 1 to 11 from pass bottom to top.
The total liquid flow rate in ON cooling modes is of vital importance in determining flow and temperature distributions in the winding. It is dictated by the power loss in the winding, the thermal head and the hydraulic resistance of the complete cooling loop. In this study, flow and temperature distributions in the winding are the focus. Therefore, a pump is used to control the total liquid flow rate in a reasonable range to represent ON cooling modes. The total liquid flow rate in ON cooling modes is typically one order of magnitude lower than the corresponding OD cooling mode [19]. In addition, the Ri's of ON cooling modes are usually larger than 1. For the winding test rig studied, the winding pass inlet velocity was set to range from 0.017 m/s to 0.033 m/s.

CFD simulations
CFD models for the experimental winding model were built using commercial software COMSOL Multiphysics, in which the governing equations for mass, momentum and energy conservation for laminar flow conditions were solved directly using finite element methods. Since the winding model was symmetric, only half of the winding was simulated in 3D models. The same mesh as shown in the OD companion paper (Part I) was adopted to guarantee mesh independent results. The same oil properties were used as shown in the companion paper (part I). The workstation used for the CFD simulation has 768 GB RAM and 2 Intel E6 CPUs. A typical 2D simulation requires tens of GB RAM and computation time of several hours, while a 3D simulation requires hundreds of GB RAM and computation time of several days to about 2 weeks.

Experimental results
In this section, some representative experimental results of temperature and liquid flow distributions in the winding model are shown under experimental parametric sweeps of different operational parameters using mineral oil as the coolant. The ambient temperature in the ON tests ranged from 22 to 24°C. It was taken that a steady state was reached when the temperature change in 1 h was less than ± 1°C. A steady state was reached in 4 h from the start of the test. The desired inlet temperature was obtained by slightly altering the external heater to counter any fluctuations in ambient temperature during the test such that the fluctuation of the inlet temperature was within ± 1°C.

Non-uniform losses (Q factor)
The power loss distribution was tested to see its effect on liquid flow distribution in ON cooling modes. The power injection in the plates was divided into two parts: the bottom 25 plates were subject to uniform power injection and the topmost 5 plates were subject to higher power injection to simulate eddy current losses at the top of the winding. The ratio of the highest power loss to the average power loss is defined as the Q factor. The power loss distributions are shown in Table 1. Tests were conducted under pass inlet velocity V in = 0.021 m/s and pass 1 inlet temperature T in ≈ 55°C. Results of flow and temperature distributions are shown in Fig. 1.
The elevated losses in the topmost 5 plates increased the temperatures of these plates but did not significantly influence the flow distribution. Therefore, uniform power loss cases still generate representative liquid flow distribution results. In the following tests, all the power injections are set to be uniform.

Pass inlet velocity (V in ) and inlet temperature (T in )
Tests were conducted under fixed pass 1 inlet temperature of T in = 60°C and uniform plate power loss of P loss = 1206 W m 2 on each plate surface. Pass inlet velocity V in was varied from 0.017 m/s to 0.033 m/s to cover a typical ON flow range. Liquid flow and temperature distributions are shown in Fig. 2. As V in is reduced, the flow distribution becomes more uneven causing more liquid to go into pass bottom cooling ducts and less liquid to flow in pass top cooling ducts. This eventually causes liquid stagnation to occur in duct11 under V in = 0.025 m/s and in duct 10 under V in = 0.021 m/s and then liquid reverse flow in ducts 8 and 9 under V in = 0.017 m/s. Flow distribution directly impacts temperature distribution. An increase in top plates temperatures is noticed as V in is reduced. Liquid reverse flow at The tests were repeated under T in = 32°C and the results are shown in Fig. 3. In comparison to Fig. 2, similar conclusions can be drawn. However, lower pass 1 inlet temperature causes relatively less distorted flow distributions, and hence the occurrence of reverse flow and its downward shift are delayed.

Inlet temperature (T in )
Tests were conducted under fixed pass inlet velocity of V in = 0.021 m/s and uniform loss distribution in each plate of P loss = 1206 W/m 2 . T in was controlled using an external heating unit connected in series to the pipe feeding the winding model and varied from 30°C to 64°C. The flow and temperature rise over the inlet temperature distributions are shown in Fig. 4. It can be observed that with the increase in T in the flow in duct 10 is gradually reversed. Because the flow magnitudes do not change dramatically, the plate temperature rises over the pass inlet temperature are not significantly different for these cases.

Power losses (P loss )
Tests were conducted under fixed pass 1 inlet temperature of T in = 32°C and inlet velocity of V in = 0.025 m/s. Power loss of each plate was varied from 804 W m 2 to 2010 W m 2 . Liquid flow and temperature distributions in the third pass are shown in Fig. 5. For the lowest power loss, the flow and temperature distributions are the evenest. As P loss increases, more liquid flows into pass bottom cooling ducts and less liquid flows into pass top cooling ducts, almost reaching a stagnation state in duct 11 under P loss = 1608 W/m 2 and eventually in duct 9 and 10 under P loss = 2010 W/m 2 . Correspondingly, temperature distribution becomes more uneven with the increase in P loss . The increase in hotspot temperature corresponds to two factors. The first one is the intrinsic effect of higher P loss leading to higher temperatures. The second one is the more uneven liquid flow distribution caused by higher P loss as pass top ducts receive lower flow rates in comparison to pass bottom cooling ducts.

Comparisons with CFD results
The measured flow and temperature results were compared with 2D and 3D CFD simulation results. The comparisons for two representative cases are shown in this section, one without reverse flows and the other with reverse flows.

Without liquid reverse flows
The case without reverse flows is with pass 1 inlet velocity of = V in 0.033 m/s, pass 1 inlet temperature of = T in 60°C and power loss of P loss = 1206 W m 2 . In the 2D and 3D CFD simulations, V in and T in were set as the pass 1 inlet boundary conditions; P loss = 1206 W m 2 corresponds to 15 W power injection for each cartridge heater; the heat transfer coefficient on the winding housing outer surface was set to be 10 W (K·m ) 2 to consider air natural convection and radiation to the surrounding. It is worth mentioning that because the housing is 12 mm thick and the thermal conductivity of the polycarbonate housing is about 0.2 W (K·m), heat transfer coefficient being 10 W (K·m ) 2 on the outer surface has limited influence on the hotspot temperature in the winding model. When the housing is treated as adiabatic, the hotspot temperature increases less than 1°C.
The comparisons for flow and temperature distributions in pass 3 (top pass) are shown in Fig. 6. Flow distribution from 2D CFD shows more uneven flow pattern than PIV measurements and 3D CFD results with lower flow rate at the top cooling duct. The overall temperature distribution comparisons show good agreement with temperature differences being within approximately 1°C. 2D CFD temperature results show subtle zigzag pattern, arising from more uneven flow distribution and hot streak dynamics. On the other hand, zigzag temperature pattern from measurements is mainly due to geometric deviation and measurement uncertainty.

With liquid reverse flows
The case with reverse flow is with pass 1 inlet velocity of = V in 0.017 m/s, pass 1 inlet temperature = T in 60°C and plate power loss of P loss = 1206 W m 2 . The same CFD simulation settings as in the case without liquid reverse flow were adopted.
The comparisons of liquid flow and temperature distributions are shown in Fig. 7. The flow distributions show good agreement with each other. The measured temperature distribution and the 3D CFD temperature distribution share the same pattern with maximum temperature difference being 3.6°C at plate 29. However, the temperature distributions show prominent differences between 2D CFD and the other two at the top part of the pass where reverse flows occur. The temperature contours of 2D and 3D CFD simulations are shown in Fig. 8.
2D CFD flow rates in duct 8-11 are closer to PIV measured results than their 3D counterparts. However, 2D temperatures of plate 28-30 are about 10°C higher than the measurements and 3D results. The reason for the big temperature differences is that a slight difference in the flow rate of the critical duct that is about to suffer reverse flow (duct 7 in this case) can bring big differences in hot streak dynamics and therefore big differences in the temperature distribution.
In duct 7, the PIV measurement and 3D CFD showed the average velocity being close to 0, while 2D CFD showed a noticeable average  velocity of 1.32 mm/s, as boxed in Fig. 7. To trace out the ins and outs of flow in duct 7, streamlines that only go through duct 7 entrance are shown in Fig. 9(a) and (c) for 2D and 3D CFD respectively. The 2D flow in duct 7 splits into two parts. The first part enters and then exits shortly from duct 7 entrance and eventually passes through the topmost duct. The second part with an average velocity being 1.32 mm/s passes through duct 7 and then is sucked into duct 8 becoming a reverse flow.
On the other hand, the 3D CFD flow in duct 7, as shown in Fig. 9(c), behaves like the first part of the 2D flow without passing through the duct. In duct 8, the liquid flow is reversed. To trace out the ins and outs of flow in duct 8, streamlines that only go through duct 8 exit are shown in Fig. 9(b) and (d) for 2D and 3D CFD respectively. The 2D CFD flow in duct 8 is mainly from duct 7, while the 3D flow in duct 8 is mainly from duct 6. 3D CFD oil temperature at the exit of duct 6 is about 10°C lower than the 2D CFD oil temperature at the exit of duct 7 and herein lies the reason why plate 27 and 28 temperatures from 2D CFD are approximately 5°C and 10°C higher than the 3D CFD and measurement counterparts, though 2D flow rates in duct 8 and duct 11 are closer to the PIV measurements than the 3D CFD counterparts. For the 2D CFD results, because plate 27 and 28 temperatures are changed and so are the consequential hot streak dynamics, the temperatures of plate 29 and 30 are approximately 14°C and 17°C higher than the 3D CFD     counterparts.
For cases with reverse flows, 2D and 3D CFD comparisons show repeatedly the pattern of comparable flow distribution and comparable temperatures at the lower part of the pass but significant temperature deviations at the higher part of the pass where reverse flows occur and 2D CFD overestimates the plate temperatures.
For ON cases without reverse flows and the flows are not close to reverse flow conditions, both 2D and 3D CFD flow and temperature results are comparable with measurements because moderate differences in flow distribution do not alter the hot streak dynamics significantly. On the other hand, for ON cases with reverse flows, 3D CFD results are more reliable than 2D CFD results because ON liquid flow pattern is not exactly 2D due to the side wall effect [14,20] and a slight difference in the flow rate of the critical duct that is about to suffer reverse flow can alter the hot streak dynamics significantly. When reverse flows occur, 2D CFD models overestimate temperatures in the reverse flow region and 3D CFD simulations are needed.

Dimensionless representation
According to the governing equations shown in (1)-(3), the dimensionless flow and temperature distributions in the winding are controlled by Re at the winding pass inlet, Ri and Pr. The Re's and Ri's for all the tested cases are shown in Fig. 10. These cases can be grouped into two clustered Re groups as shown in Table 2. In each group the relative variation of Re is small and can therefore be treated as invariant. There is an outlier with Re being 155, Ri 0.97 and Pr 52.
Pr is less influential than Re and Ri. In addition, for the two groups of tested cases, the product of Re and Pr, which is the reciprocal of the coefficient in Eq. (3), is quite similar. Therefore, the dimensionless flow (flow fraction) and dimensionless temperature distributions can be treated as only controlled by Re and Ri. Group 1 has a wider range of Ri and is therefore taken as an example to show how Ri controls flow fraction and dimensionless temperature distributions.

Flow fraction distribution
The liquid flow fractions in the horizontal cooling ducts are converted from the PIV measurements and the evolution of the flow fraction distribution with Ri for tests in Group 1 is shown in Fig. 11. At low Ri (close to 1), the flow distribution is quasi-uniform. With the increase in Ri, i.e. the increase in buoyancy force weight over inertia force, the flow distribution gets more and more uneven with higher flow fractions at the bottom part of the pass and lower flow fractions at the top part of the pass. When Ri is further increased, reverse flows start to occur from the penultimate duct and shift downwards to the middle part of the pass.
3D CFD simulations were carried out with T in = 36°C, v in = 0.025 m/s and P loss ranging from 402 W/m 2 to 1809 W/m 2 . Because of the increase in power injection, Re increases from 54 to 82, Ri from 0.55 to 3.9 and Pr decreases from 108 to 74. The Re range and Pr range of this group of 3D CFD simulations are similar to those of the tested cases in Group 1. Therefore, the 3D CFD dimensionless results are comparable with those of the tested cases in Group 1. It is worth mentioning that the power loss in the 3D CFD sweep does not increase to the level that brings Ri to 6 because the purpose of the 3D CFD sweep is to show the flow and temperature distribution variation trends with Ri up to the condition where reverse flows have already occurred. In addition, it becomes increasingly computationally costly to get the CFD model converged when reverse flows start to shift downwards. The flow fraction distribution evolution with Ri from the 3D CFD simulations is shown in Fig. 12. It shares the same trend shown by the PIV measurements.
The ON flow distribution pattern is dictated by the hydrostatic pressure distribution in the two vertical ducts [6] and the hot streak dynamics [12]. To elaborate on this point, a pass with its vertical ducts divided into liquid columns, each next to a horizontal duct, is shown in Fig. 13. Assumptions and denotations are made as the following: (1) The hydrostatic pressure at the entrance or exit of a horizontal duct is the pressure at the bottom of the adjacent liquid column. The difference of the hydrostatic pressure (dp) generated by the two liquid columns on the two sides of a horizontal duct is 11 11' 11 11' When there are no reverse flows, the temperature of the column in Outlet vertical duct Therefore, the hydrostatic pressure gradient over the horizontal duct accumulates from the top of the pass to the bottom and herein lies the reason why ON horizontal duct flow fraction increases from pass top to bottom.
However, when the total liquid flow rate decreases or the power loss increases, i.e. Ri increases, the liquid flow fractions passing through the bottom cooling ducts can be so high that the temperature of the liquid reaching the top part of the outlet vertical duct is lower than that of the liquid at the top part of the inlet vertical duct. Therefore, the hydrostatic pressure gradient at the top part of the pass is reversed. However, the reverse flow cannot occur in the topmost horizontal duct because that will contradict the upward flow currents. The reverse flow occurs first in the penultimate horizontal duct and then returns through the topmost duct. When Ri further increases, the reverse flow in the penultimate duct will be further heated so that the hot liquid reversely flowed to the inlet vertical duct can result in the accumulative hydrostatic pressure on the entrance side of the horizontal duct beneath being lower than the hydrostatic pressure on the exit side, leading to reverse flows extending and shifting downwards. The reverse flows eventually pass through the top horizontal ducts and therefore cause the gradual recovery of the flow direction and magnitude for the pass top horizontal ducts, as shown in Fig. 11 and Fig. 12. When the hot streak dynamics are fully considered in scenarios with reverse flows, the situation becomes more complicated as pass to pass coupling is more pronounced [12]. However, the overall flow pattern remains the same.

Dimensionless temperature distribution
The evolution of the liquid flow distribution brings corresponding temperature evolution, as shown in Fig. 14 from measurements and Fig. 15 from 3D CFD simulations. The uniform flow at low Ri leads to uniform temperature distribution. Uneven flow distribution at high Ri causes uneven temperature distribution with top part of the pass being overheated. There is a relatively narrow region of Ri, which is 2-3 for the group 1 tested cases, where flow and temperature distributions become rapidly uneven with the increase of Ri. When reverse flows occur and the flow reversed region shifts downwards, the hotspot temperature fluctuates, and its position shifts accordingly. It is worth noting that in practice top part of the winding experiences extra power losses due to eddy current losses and the oil starvation at the top part of the pass can therefore result in an even higher temperature rise.

Hotspot factor
The highest dimensionless temperature is in fact the hotspot factor (H), which characterises the uniformity of temperature distribution in the winding [7]. The variations of hotspot factor with Ri for the tested  cases and 3D CFD simulations are shown in Fig. 16. Two more sets of 3D CFD power loss parametric sweeps with modified viscosity and pass inlet temperature were conducted to cover a wider range of Re.
The CFD results share the same H-Ri variation pattern. The curve can be divided into 3 regions. The first region is a flat region with low Ri, and it indicates uniform temperature distribution. The smaller than 1 hotspot factor is because the winding model investigated has no paper insulation wrapped on the winding plates. The second region shows a steep rise of H against Ri, caused by local overheating at the top part of the winding due to increasingly uneven flow distribution and reverse flows. The third regions features a fluctuation of H against Ri, arising from the shift of the reverse flow region. In addition, higher Re shifts the curve to the left, meaning the uneven flow distribution and reverse flows occur at a lower Ri. The experimental results follow the same trend as shown by the CFD parametric sweeps.
For practical ON transformer thermal design, the region with high Ri that brings the fluctuation of H is of less interest because it needs to be avoided. On the other hand, the region where H starts to increase steeply with increasing Ri is critical in the thermal design. Adequate total liquid flow rate is needed to guarantee a low Ri and a relatively uniform temperature distribution to avoid localised overheating at the top part of the winding.

Conclusion
Liquid flow and temperature distributions in an ON transformer winding model have been investigated experimentally by PIV and temperature measurements for representative operational conditions of different liquid inlet velocities, inlet temperatures and winding power losses. The ON uneven liquid flow and reverse flow pattern is elucidated by a detailed hydrostatic pressure gradient analysis in the two vertical cooling ducts.
Comparisons of flow and temperature distributions are made between experimental tests and CFD simulations. Before the occurrence of reverse flows, 2D and 3D CFD results are comparable with measurements. When reverse flows occur or is about to occur, 2D CFD flow distribution results are still similar to 3D CFD results obtained from the symmetry plane and PIV measurements, but 2D CFD temperature results are significantly higher than those of 3D CFD and measurements at the reverse flow region. This is arising from a subtle change of hot streak dynamics in the horizontal duct that is about to have reverse flows. 3D CFD results are comparable with measurements.
From dimensional analysis, the dimensionless liquid flow (flow fraction in horizontal ducts) and dimensionless temperature distributions are found to be mainly controlled by Ri. The dimensionless flow and temperature distributions get increasingly uneven with the increase in Ri. There is a relatively narrow region of Ri where flow and temperature distributions get rapidly uneven with the increase in Ri. When Ri further increases, reverse flows and localised overheating occur at the top part of the pass and shift downwards. These flow and temperature distribution trends are generically applicable to any other disctype winding geometries.
Flow fraction and dimensionless temperature distribution results from 3D CFD simulations match well with those from measurements. The comparisons of the highest dimensionless temperature, the hotspot factor (H), show a clear H-Ri relationship and the H-Ri curve shifts left with the increase in Re. The H-Ri curve features a flat and small H region at low Ri region, then a steep rise of H with increasing Ri and followed by a fluctuation of H at high Ri region.
Comparing OD cooling modes studied in the companion paper (Part I) and ON cooling modes studied in this paper (Part II), one can find differences in the dimensional and dimensionless parameters. For dimensional parameters, OD cooling modes have higher total liquid flow rate and higher power loss than ON cooling modes. For dimensionless parameters, OD cooling modes have high Re and Ri usually much smaller than 1, whereas ON cooling modes have lower Re and Ri usually larger than 1. The flow and temperature distributions in OD cooling   modes are mainly controlled by Re with high Re causing reverse flows at the bottom part of the pass, indicating that it is not necessarily the higher the total liquid flow rate the better the cooling performance. The flow and temperature distributions in ON cooling modes are mainly controlled by Ri with high Ri causing reverse flows at the top part of the pass, indicating the importance of adequate total liquid flow rate for ON transformer thermal design.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.