Hydraulic network analysis

Initialization

If mass flows have been calculated in a previous time step or iteration, then the previously calculated values are used.

If no mass flows have yet been calculated for any of the hydraulic elements, all mass flows, are initialized for 2-node hydraulic elements with a duct Reynolds Number Re = 10000.

The first time the hydraulic loop is entered, the initial total pressures, P_T, of elements without defined pressure sink boundaries are set to zero, and the initial static pressures, P_S, are set to zero gauge pressure.

Increase iteration counter

The iteration counter hn is increased by one.

Flow properties calculation

Mass flow adjustment

If the mass flow through a flow resistance is entirely due to buoyancy forces, to increase the stability of the solution the change in mass flow from that of the previous iteration is multiplied by the damping factor.

Duct Reynolds number calculation

The duct Reynolds Numbers are computed, as

where:

Re_duct is the duct Reynolds Number
is the mass flow through a hydraulic flow resistance between elements I and J.
D_duct is the hydraulic diameter of the element.
A is the cross–sectional area of the element.
μ is the viscosity of the fluid.

Mass flow limits check

The mass flows are then checked to see whether they fall outside the limits of 1<Re_duct<108. If the mass flows fall outside this limit, they are appropriately adjusted.

Density calculation for gases

If the compressible flow option is not specified, the density for gases is computed using the Ideal Gas Law:

where:

ρ is the density at the element.
P_S is the static pressure.
T is the temperature of the element.
T_ABS is the absolute temperature.
R_U is the Universal Gas Constant.

T_ENV is the standard temperature at which the fluid densities are defined.
P_SENV is the standard pressure at which the fluid densities are defined.
ρ_nom is the nominal fluid density.

Density change check

In order to avoid large solution oscillations, a check is performed to see whether the change in density from that of the previous value brings the new computed density below 25% of the nominal density specified for the element. If it does, the change in density is clipped. Then, the change in density is damped by multiplying it with the dumping parameter.

For liquids the change in density = 0.

Buoyancy pressure rise calculation

The buoyancy pressure rise over each 2-node hydraulic element is computed as:

where:

Δ(P_B) is the buoyancy pressure increase over the element.
g is the acceleration of gravity.
ρ_amb is the density of the ambient air at standard temperature and pressure.
Δ(h) is the difference in elevation between the two ends of the 2–node hydraulic element.

Flow velocities calculation

Next, the flow velocities are computed:

where V is the velocity.

Dynamic pressure calculation

The dynamic pressures are computed:

where P_D is the dynamic pressure.

Static pressure recalculation

The static pressures are recomputed:

Centrifugal effects calculation

If hydraulic elements undergo rotation about an axis., Δ P_B is calculated by:

where:

ρ_CG is the density at the center of gravity (CG) of the element.
ω is the circular frequency of the centrifugal force.
H1,H2 are the perpendicular distances of the two ends of the hydraulic element from the axis of rotation.

Compressible flow calculations

If the compressible flow option is specified, the density and static pressure calculations take into account compressible effects for all elements with Mach numbers > 0.1. For these elements, the following procedure is followed.

Upstream calculations

For each two-node hydraulic element I, the upstream velocity, density, and static pressure are calculated from the upstream temperature, mass flow, and total pressure:

where

V_U is the upstream velocity.
A_U is the upstream cross-sectional area.
ρ_U is the upstream density.
T_U is the upstream temperature.
P_SU is the upstream static pressure.
P_TU is the upstream total pressure.
is the mass flow through the hydraulic element.
T_ABS is the absolute temperature.
R_U is the Universal Gas Constant.

Upstream Mach number calculation

The upstream Mach number is calculated using the ratio of the specific heats and the upstream temperature:

where:

Ma_U is the upstream Mach number.
k_ratio is the ratio of the specific heat at constant pressure to the specific heat at constant volume (fixed at 1.4).

Upstream stagnation temperature calculation

where T_SU is the upstream stagnation temperature.

Heat net input calculation

The net heat input per unit mass flow for the hydraulic element I is computed:

where:

C_p is the specific heat for the element.
Q_unit is the net heat input per unit mass flow for the hydraulic element.
GT_{SUM_I} is the sum of the conductances times connected element temperatures for element I.
Q_I is heat input into hydraulic element I.
G_{SUM_I} is the sum of the conductances of element I.
T_I is the temperature of element I.

Downstream stagnation temperature calculation

where T_SD is the downstream stagnation temperature. T_SD is only affected by the heat input into the flow.

Equivalent friction parameter calculation

The equivalent friction parameter K_{loss_eq} for the hydraulic element is computed.

For non-FANPUMP type elements as

For FANPUMP type elements:

where:

P_{T_D} is the downstream total pressures of the FANPUMP element.
K_loss is the head loss factor.
P_D is the dynamic pressure of the element at its narrower end.

Mach Number differential equation

The differential equation for Ma²(x), the square of the Mach number at a normalized distance x along the hydraulic element is solved:

where

A_D is the downstream cross-sectional area.
Ma²(x) is the square of the Mach number at a distance x in the element.

Applying the boundary condition:

where Ma_U is the upstream Mach number.

The downstream Mach number is then computed from:

where Ma_D is the downstream Mach number. The following limitation is imposed 0<Ma_D<1.

When choke flow (Ma_D∼1) occurs, a warning message is written to the verbose log file.

Downstream temperature calculation

Downstream static pressure calculation

Next, the downstream static pressure P_SD is calculated:

where:

P_SD is the downstream static pressure.
P_SU is the upstream static pressure.

Downstream velocity calculation

where V_D is the downstream velocity.

Downstream density calculation

where ρ_D is the downstream density.

Equivalent compressible heat input calculation

In order to compute the downstream temperature using the standard temperature computational routines Q_CD is computed:

where Q_CD is an equivalent compressible heat input value which is applied to the immediately downstream 1-node hydraulic element to account for the compressible heating effects. Q_CD is summed with all the other heat inputs into the downstream element.

Downstream total pressure calculation

where P_TD is the downstream total pressure.

Flow resistance calculation

where RES_IJ is the value of the hydraulic flow resistance.

Friction heating effects calculation

The friction heating effects on the wall are computed. Due to friction, under adiabatic conditions the wall temperature will approximately equal the temperature of the stagnation temperature of the fluid. To model this effect, first the difference between the static and dynamic temperatures is computed for the 2-node hydraulic element I:

where ΔT_I is the difference between the temperature and stagnation temperature of the fluid at the center of the element.

Next, for each wall element connected to the fluid element with convective conductance an additional heat input Q_{W_I} is computed.

where Q_{W_I} is the corrective heat input applied element I to account for the friction heating effects. To maintain energy conservation, Q_{W_I} is also subtracted from Q_I.

Flow resistances calculation

Next, if compressibility effects are not taken into account, flow resistances RES_IJ are computed as described in Hydraulic flow resistances for incompressible flows.

Clipping on flow resistances

If the fractional change between two subsequent iterations in the value of a hydraulic flow resistance is greater than the hydraulic iteration damping parameter, the change in the flow resistance is appropriately clipped. This is done to avoid solution oscillations.

Mass flow and flow resistance table interpolation

Next, table interpolations are performed. Table interpolation may be used to multiply computed hydraulic flow resistances (see Hydraulic flow resistances for incompressible flows), and to calculate FANPUMP mass flows and pressure rises.

Each 2–node FANPUMP element must reference a table, which is used to specify velocities, mass flows, or total pressure rise as a function of mass or volume flow through the element.

Velocities and volume flows are transformed into mass flows using the fluid density and FANPUMP cross–sectional area values. These in turn are transformed into a mass flow couple: a negative mass flow at the inlet and a positive mass flow at the outlet 1-node elements of the FANPUMP element. The FANPUMP element itself has an infinite hydraulic flow resistance. The total pressure at the center of the FANPUMP element is estimated to be the average value of the total pressures at the end elements.

If total pressure rise as a function of mass flow is specified for a FANPUMP element, first the mass flows in and out of the end elements are computed from the current estimates of the flow resistances from the ends of the FANPUMP element and the total pressures of the elements they are joined to. Then, the total pressure rise value over the FANPUMP element is interpolated from the table using the average of these mass flow values. The total pressure rise value at the CG of the FANPUMP element is estimated to be the average of the total pressures values of the end elements. Last, the total pressures of the end elements are then computed to be P_{T_CG}+ΔP_T/2, and P_{T_CG}−ΔP_T/2 the end elements are defined to be pressure sink elements.

Using the user–written subroutine

Next, the user–written subroutine is entered. At this time, you can alter hydraulic flow resistances and total pressures for sink elements.

Buoyancy-driven mass flow calculation for each branch

Next, the total buoyancy pressures, the total flow branch flow resistances, and the buoyancy–driven mass flow are computed for each flow branch.

where:

is the computed mass flow through the flow branch due to buoyancy effects in the flow branch.
ΔP_Btot is the sum of the buoyancy pressures in the flow branch computed by summing the buoyancy pressures of each element over the branch.
RES_tot is the total hydraulic flow resistance of the flow branch.

Once is computed, it is modeled as a mass flow couple: positive and negative mass flow boundary conditions are added at the beginning and end 1-node elements of the flow branch.

Calculation of total pressures and mass flows in branches

Once the mass flow inputs and total pressure values for the pressure sink elements are defined, the matrix is simplified by reducing it to flow branches. Here the flow branches are as defined above, with the exception that the endpoints of a FANPUMP element are considered to consist of the start of a new flow branch. The purpose of the simplification is to reduce the size of the matrix to be solved. The simplification consists of summing the flow resistances in a flow branch.

Next, the total pressures at the endpoints of each flow branch are computed from the mass balance equation, i.e. the sum of the mass flows into each hydraulic element must be zero.

hence:

where:

IE,JE are 1-node elements at the ends of a flow branch.
is the impressed mass flow into element I, computed from FANPUMP elements and from couples.
P_{T_IE}, P_{T_JE} are the total pressure values at elements IE and JE where
is the mass flow through hydraulic flow resistance RES_IE,JE.

This matrix equation is solved by a conjugate–gradient sparse matrix inversion routine.

Mass flows through elements calculation

The total pressures for each element I in each branch are calculated from the computed total pressures P_{T_IE} and P_{T_JE} . Starting from the end IE, the first total is computed:

and so on.

Next, the mass flows are computed:

Convergence check

If no clipping has been performed on any conductances, the maximum total pressure change for any element from that of the previous iteration is evaluated and is compared with the convergence parameter for the pressure/flow solution, the maximum allowable total pressure change for convergence to have occurred. The convergence parameter for the pressure/flow solution defaults to 1% of the difference between the maximum and minimum total pressure values in the model, but may be overridden.

If the solution has not converged, the sequence of steps from increasing iteration counter onwards is repeated.

If the solution has not converged after the maximum number of iterations, a warning message is printed to the verbose log file. The maximum number of iterations defaults to 100, however this may be overridden.

Calculation completed

Total pressure/mass flow calculations completed.