See also this document
It is frequently required to simulate fluid-flow and heat-transfer processes in and around solids which are, partly as a consequence of the flow, subject to thermal and mechanical stresses.
Often, indeed, it is the stresses which are of major concern, while the fluid and heat flows are of only secondary interest.
(b) Practical occurrence
Engineering examples of fluid/heat/stress interactions include:
It has been customary for two computer codes to be used for the solution of such problems, one for the fluid flow and the other for the stresses
Iterative interaction between the two codes is then employed, often with considerable awkwardness.
This interfacing is rendered unnecessary by the solid-stress option of PHOENICS, which permits fluid flow, heat flow and solid deformation, and the interactions between them, all to be simulated at the same time.
The method of doing so exploits the similarity between the equations governing velocity (in fluids) and those governing displacement (in solids).
In order to demonstrate what this method can achieve, the results from PHOENICS library case s400 will be shown. This concerns the calculation of fluid-flow, conjugate heat transfer and solid stress in the situation illustrated by the sketch below.
cooling | air
|
| V |/////// hot radiating wall ///////////|
| ----------------------------------------
| duct -----> exit
|------------- -------------
|// steel ///| cavity |/// steel /|
|-------------------------------------------
|////////////// aluminium /////////////////|
|-------------------------------------------
The task is to discover what stresses arise in the metal blocks as a result of the non-uniformities of temperature.
Further details are:
(b) Vector and contour plots
Vectors
The velocity vectors displayed in Fig.1 reveal the nature of the air-flow pattern.
They are calculated at the same time as the displacement vectors shown additionally Fig.2; but, of course, the two sets of vectors have different scales, and indeed dimensions.
The solids are supposed to be confined by a stiff-walled box, but are allowed to slide relative to its walls. This is why the displacement vectors are vertical near the confining-box walls.
They are however not allowed to slide relative to each other; this is what causes the concentrations of stress at their contact surfaces.
The stresses in the x- (horizontal) and y- (vertical) directions are displayed in Fig.3 and Fig.4 respectively.
They have been deduced from the strains shown in
Fig.5 for the x-direction , and in
Fig.6 for the y-direction,
The strains have been deduced from the displacements by
differentiation.
They are displayed in
Fig. 7 for the x-direction
and
Fig. 8 for the y-direction.
However, their implications are less visually dramatic,
Fig. 9 displays contours of temperature in the air and the solid, of which the most noticeable feature is that the air is heated by contact with the 80-degree-Celsius top wall and by the metal blocks which have been receiving heat by radiation from above.
Temperature differences within the high-conductivity solids are too small to be discerned visually.
It is interesting to compare Fig. 9 with
Fig. 10,
This displays
the distribution within the air space of the "radiation temperature",
T3,
which is computed by IMMERSOL, and which is defined as the temperature
which would be taken up by a probe which was affected only
by radiation.
Obviously, and understandably, T3 and TEM1, have very different values.
The solid temperature influences the stresses and strains, of course, primarily through the agency of the temperature-dependent thermal- expansion distribution. However, its variations with position are too slight to be revealed by a contour diagram, as inspection of Fig. 11 will reveal.
The IMMERSOL model, of which the solution of the T3 equation is the major feature, enables the radiant heat fluxes in the coordinate directions to be established by post-processing.
The results are displayed in
Fig. 12 for the
x-direction. and by
Fig. 13 for the
y-direction.
The shapes of these contour diagrams, if studied and interpreted in physical terms, will be found to be plausible.
A crucial feature of the IMMERSOL model is its use of the distribution of the "distance between the walls", WGAP. This quantity, which has an easily-understood meaning when the walls are near, and nearly parallel, is computed from the solution of the "LTLS" equation.
The distributions of these two quantities are shown by
Fig. 14
for the former, and by
Fig. 15
for the latter.
It will be seen that WGAP has a uniform value in the region of between the top of the duct and the tops of the upper metal slabs, between which the actual distance is 0.008 meters.
Further, it has approximately twice this value near the convex corners; and it becomes zero in the concave corners.
The flow field was calculated by means of the LVEL turbulence model, which makes use of the wall-distance field which is also derived from the LTLS distribution.
The contours of WDIS are displayed in
Fig. 16.
which exhibits the expected maximum of 0.004 between the parallel
horizontal walls, and a somewhat greater value near the cavity.
where the true distance from the wall depends on the direction in
which it is measured.
LVEL, like IMMERSOL, is a "heuristic" model, by which is meant that it is incapable of rigorous justification, but is nonetheless useful.
WDIS is calculated once for all, at the start of the computation. From it, and from the developing velocity distribution, the evolving distribution of ENUT, the effective (turbulent) viscosity is derived.
The resulting contours of ENUT are shown in Fig. 17.
Since the laminar viscosity is of the order of 1.e-5 m**2/s, it is evident that turbulence raises the effective value, far from the walls, by an order of magnitude.
The similarities referred to in section 1 (d) are here decribed for only one cartesian direction; but they prevail for all three directions.
[del**2]* U + [d/dx]* [ D*C1 - Te*C3 ] + Fx*C2 = 0
[del**2]* u - [d/dx]* [ p*c1 ] + fx*c2 = 0 ,
where
Notes:
[del**2]* U + [d/dx]* [ D*C1 - Te*C3 ] + Fx*C2 = 0
and
fx * c2 = Fx * C2
The strains (ie extensions ex, ey and ez) are obtained from differentiation of the computed displacements.
Thus: ex = [d/dx]* u, ey = [d/dx]* v, ez = [d/dx]* w .
Then the corresponding normal stresses, sx, sy, sz, and shear stresses tauxy, tauyz, tauzx, are obtained from the strains by way of equations such as:
sx = {YM / (1 - PR**2)} * {ex + PR*ey}
tauxy = {YM / (1 - PR**2)} * {0.5 * (1 - PR)*gamxy}
where:
PHOENICS employs (a variant of) the "SIMPLE" algorithm of Patankar & Spalding (1972) for computing velocities from pressures, under a mass-conservation constraint.
Its essential features are:
All that it is necessary to do, in order to modify this algorithm for solving for displacements simultaneously, is to treat the dilatation B as the mass source error and to ensure that p obeys the above linear relation to it.
Therefore a CFD code based on SIMPLE can be made to solve the displacement equations by:
The "staggered grid" used as the default in PHOENICS proves to be extremely convenient for solid-displacement analysis; for the velocities and displacements are stored at exactly the right places in relation to p.
The following sketch illustrates a combined mechanical-and-thermal-stress situation, in which the fluid flow plays a part. It is the focus of attention in a series of case to be found in the Input-File Library of the solid- stress option,
------------------------------------------------------
load V V <--------- ^ y
|-----------| |
|///////////| cold-air jet |
|// heated /| |
|// block //| |
|///////////| <------z-----O
----------------------------------------------
The task is to compute displacements, and thence strains and stresses, in the heated block, subject to a variety of constraints as follows:
PHOENICS solves the fluid-flow, displacement and thermal-energy equations simultaneously, being enabled to do so easily because:
The user needs only to set STRA=T in the Q1 file, and to provide appropriate boundary conditions for displacement, etc; then PHOENICS will test, cell-by-cell, whether a solid or a fluid is present, and will solve whichever equation is appropriate there.
This is done whenever the solid-stress option is present.
Only linear stress-strain relations are allowed for; but the Young's Modulus and the thermal- expansion coefficient can be temperature-dependent.
Graphically-displayed results for this case are as follows:
Fig.18 shows simultaneously the displacement vectors and the velocity vectors; the second shows the temperature field which caused the thermal expansion.
The arrows represent velocity vectors in the upper region and displacement vectors in the block, and Fig.19 shows the computed temperature distribution which caused (most of) the displacements.
A library of input-files is supplied with the solid-stress option of the PHOENICS.
It has not been common to choose the displacements as the dependent variables in numerical-solution procedures. However, this has been done by:
Their numerical procedure differ from that used here, which was that of
Reports describing the early work include:
