SAFPWR home page

Kinetic momentum balances for coolant flow calculations

Differential form

Let us consider an element of linear hydraulic circuit, of transverse area sw and length dl oriented in the direction of the flow w.
Let us suppose, for the moment, that the element is homogeneous and that the 2 phases move at the same celerity c.

The differential form of the kinetic momentum conservation equation (Euler equation) applied to the element and projected on the "neutral fiber" of the circuit is a scalar differential equation (1), where
dz is the elevation increase,
dp is the pressure variation
From (2), the temporal acceleration term becomes (3).
The friction term dpf is generally represented (4) as fraction df of the kinetic energy par kg c2/(2 vm)
To simplify notations, let us introduce (5) the variables a and y so that the "kinetic balances", integrated over the part e to s of the circuit simplifies to (6).

Integration on the node

Picture euler_2 represents the node, delimited by the sections e and s.

For the sake of simplification, the transverse injection flow wy and the internal "expansion source" s= (m1-m2)/δs effects are supposed to be concentrated in a thin, constant transverse surface area sw, entry volume, delimited by the sections e and ê of the node, so that oulet values ws and vms of the entry volume will apply to the whole remaining volume between the sections ê and s.

It should be remembered that, in the course of the (m,h,b)-balances, the pressure, as well as the reference entry flows(wel, wea, wed for primary, wsan for secondary), are mean values, constant over δs .
The last m balances do, circuit_name have provided the nodal expansion sources s=(m1-m2)/dt.
The kinetic balances, enabled by the key-word euler recuperates those sources as such and aim now at updating the reference inlet flow of the flow branch containing the node .
Thus w1 actually represents the w assumed in the course of the last mass balances and w2 the searched new flows to be updated by means of the Euler balances.

(9) details the calculation of the spatial acceleration term.
The primary configuration allows (pic_3) arranging the nodes as a set of non intersecting branches having common terminal points O and A : the branches l=(1,l9), the branch abc and the branch d.
The pressure drop in the outlet o is neglected.
A is not a node but a "flow connection point"
The branch through the pressu is not retained for flow calculation because the flow wso outgoing o towards the expansion branch X of the pressu system is determined by closure of total primary system volume.
Likewise, the flow wespr entering the spray is imposed by the pressu regulation and, therefore, considered here as a fixed value.
In (9), for two successive nodes in the branch, the outlet term as2 ys of the first node eliminates whith the inlet term -ae2 ye of the second, and if the branches are joined as a loop, like the recirculation loop of the sg those terms disappear, as well as the Δp terms.
Define a nodal average an from (10)
and a composite nodal friction + spatial acceleration pressure drop factor f (11).
f depends on the node geometry and the distribution of singular and distributed friction factors.
It is calculated outside of the program. If requested, the f dependence with ws could be represented by a fitted correction.
With notations (13:16), (17) represents the simplified form of the pressure loss pe - ps across any node.

Temporal discretization

Summing on the nodes li of loop_l, using the time discretization (18), and taking the flow-dependant terms at eos to insure stability, we obtain
(19) the kinetic balance relation for loop_l.
The vm are the last calculated values by the m,h-balances.
The loop gravity pressure increment apzl is calculated (20) from pump head azpl [Accretion,Pressure,due to elevation Z,for loop L], and node geometrical elevations.
As azpl is updated at bos (step, euler) it is advisable to put euler as the first keyword of the working sequence, in which case the vm values as well of the "sources sli" are those of preceding step.
(21) is the kinetic relation for branch abc.

Note that, although the outlet o does belong to a branch, we have added a gravity contribution apzo for it, in order to insure that the sum of the az on closed loop abcod is zero.


The relation (22) for the dome branch had to be simplified because of the way the "cold dome" mixing effect is modeled: the pressure loss due to friction and spatial acceleration are represented by a single term using the upstream values of w and wm.
With (18), the loop_l kinetic balance is (19).
Note that to insure stability, the quadratic flow terms are taken at eos. The density effects (vm2) are those resulting from the last mass balances.
Linearization
In order to linearize the equations (w2)2 is developed (23) about its current value w whose first guess is the bso w1, and the quadratic term is neglected.
Continuity relations
The nodal outlet flows w2li are evaluated (24, 25) in terms of the reference branch inlet w2el and the nodal "sources" sli obtained from last m,h-balances.
An external continuity relation (26) links the reference branch flows.

We thus obtain a linear system of l9+3 equations in searched flows w2al, w2ea, w2ed, and the common apOA ≡ pA -pO whose coefficients depend on the w1li, vm2li, sli, wli
After resolution, wli is updated as w2li, and the calculation is repeated until convergence of the reference flows, governed by the convergence criteria (?).
This method provides an exact solution of the original system of implicit, non-linear, kinetic balance equations.

Initializations

Distribution of flow in the primary hydraulic network result from the initialization (ini, down-comer, bottom,..) of the primary components m,h,b-balances where the reference flows are imposed as the volumic (wvpl) or massic (wel) for the loops, and from wed or xwld for the dome bypass.

The module (ini, euler) can thus compute (from (19, 20, 21) with the transient contributions w2li - w1li dropped) the pressure loss apOA contributed by the various branches.
These apOA generally do not coincide because the input flows are design or measured values, whilst the hydraulic friction data are best calculated.
In order to insure consistency, the friction factors cfrli, ... are operated by multipliers xcfrli, .. calculated in (ini, euler).
Here, the loop_1 azPA has been taken as reference and the other branches are adjusted. It would possibly be more sensible to impose apOA and adjust all friction factors on it.
For the recirculation loop in the sg, the friction multiplier is adjusted on the recirculation ratio.
Application test
Euler4_dat is a simple test for the euler module.
It is a single loop case (l9= 1)
with thermal power f2c generated by a heat-generator (oqc= 1)
The initial flow is imposed (owvp= t) as volumetric flow (Lstl1/wvpl).
Pump head (Lstl1/Itp_sec/azpl) is cut to zero from 20 s to simulate flow slowing down to natural circulation.
Steam generator is simulated through an uniform secondary temperature field (Lstl1/input_t1cn)
After flow slowdown, natural circulation will be activated by (heat_gen/&Lstc/Itp_sec/f2c) generating a power of 2*10e6 from 400 s.
Althoug the sg is not explicitly present, the group name S00gvg must appear in the read group in order to allocate some sg variables.
omqgv and epsqgv control the t1cn initialization iterations.
ini, call:2

the primary components are initialized as usual starting from a with a cold enthalpy hsb= 1.195179e6, corresponding to a temperature of 545 K, which should be imposed as t1cn to the secondary because the system is initially isothermal (f2c=0).
Nevertheless ini, loop_1 adjusts t1cn to insure equilibrium initial state, so that any sensible guess value may be entered as input.
As indicated above, euler is last in the ini, call:2 sequence but is first in the working sequence S4.

Results
The execution listing safpwr.lis allows following the initialization iterations of the secondary sg temperature t1cn.
Starting from a guess of 500 K, it takes 10 iterations to converge to the correct value of 545 K. As the primary system receives no power, the power lost by the primary water in the sg tubes (zl%qgv) converges to 0. As already mentioned, the iterations must be over-relaxed to insure rapid convergence.
Follow on the listing some results of ini,euler: a%apoa= pA - pO, and the branches friction multipliers. The loop_1 multiplier bl%xcv is 1 because that loop is the reference one.
The branch abc multiplier xcv_a is close to 1 which indicates the friction parameters of the loop and the abc branch are consistent with the input primary flow and pump head. On the other hand, the dome branch multiplier xcvd is far from 1, which means that the dome branch friction factor is not consistent with the imposed flow by-passed towards the dome!
On chart 1, it is first of all observed that the system remains stationary until the azpl step shutoff at 20 s.
Chart 2 depicts the branch flows slowdown after the pump head cutoff.
In fact, the flow reduction is not so steep because it does not account for the pump rotor inertia, which could possibly be entered though a azpl vs time curve.
The effect of reducing the time step dsec from 5 to 1 s is observable on the companion curves (limited to 300 s): large dsec have always the effect of damping the flow variations.
This is because, in the implicit treatment, the quadratic friction effect is taken at eos, hence under-estimated.
Note also that the abnormal transient of wed around sec= 30 is a time-step effect which disappear on the dsec=1 s case.
From sec=300, wed becomes negative so that its absolute value is plotted thereafter to allow logarithmic scale representation.
Chart 1 shows also the pressure head apoa ≡ pA - pO which pushes the flows through abc and apoa becomes also negative from sec=300 and represented as its opposite value.
However, wea remains always positive
From sec=300, a step power q2c is delivered by the heat-generator.
Graph1 compares q2c/3 (3 loops!) with the power qvgl1 lost in the sg.
It takes about 100 s to reach a new natural thermo-siphon circulation regime. Observe that apoa and wed remain negative!