- Email: [email protected]
Offshore drilling rig ice accretion modeling including a surficial brine film
Offshore Drilling Rig Ice Accretion Modeling Including A Surficial Brine Film Ivar Horjen PII: DOI: Reference:
S0165-232X(15)00165-2 doi: 10.1016/j.coldregions.2015.07.006 COLTEC 2138
To appear in:
Cold Regions Science and Technology
Received date: Revised date: Accepted date:
12 May 2014 13 July 2015 20 July 2015
Please cite this article as: Horjen, Ivar, Offshore Drilling Rig Ice Accretion Modeling Including A Surficial Brine Film, Cold Regions Science and Technology (2015), doi: 10.1016/j.coldregions.2015.07.006
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT 1
OFFSHORE DRILLING RIG ICE ACCRETION MODELING INCLUDING A SURFICIAL BRINE FILM
IP
T
Ivar Horjen1, dr. techn., Trondheim, July, 2015
SC R
ABSTRACT
This article presents an extension of the two-dimensional icing model ICEMOD2. In contrast
NU
to this model the new model (called ICEMOD2.1) considers the inertia terms in the momentum equation. The derivation of the model is presented and the model is applied to two
MA
typical structures of an offshore drilling rig which may be exposed to spray icing: a corner column and a truss. A mean impact-generated spray mass flux formula (called the “normal”
D
value) based on droplet kinematics and spray measurements on the Norwegian semi-
TE
submersible rig “Treasure Scout” has been implemented in the model. The model results
CE P
showed satisfactory agreement with some actual icing registration data from the drilling rigs “SEDCO 708” and “SEDCO 709”. Additionally sensitivity tests are presented for a selected section of a rig column of diameter 5.0 m varying air temperature, wind speed and spray
AC
frequency. As expected ice accretion increases with decreasing air temperature. With the chosen formula for the spray mass flux and wind speed (17-25 m/s) ice accretion also increases with increasing wind speed. The available field results of Treasure Scout enables us only to construct a time-average spray mass flux formula independent of spray frequency. The single spray mass flux is then in the model a derived quantity obtained by dividing the timeaverage value by the product of spray frequency and spray duration. Based on these presumptions the sensitivity tests show that spray frequency variations has only a minor effect upon accreted ice load when normal spray mass flux is assumed (corresponding here to 5-10 g m-2 s-1). If however the spray mass flux is put 25 times as large as the normal value a 1
Mail address: [email protected]
ACCEPTED MANUSCRIPT 2
substantially increase in model ice load with increasing spray frequency is obtained. Probably this is due to a large increase in brine film thickness (a factor of about 10). A
T
possible theoretical explanation is given in the text Finally the significance of brine film
IP
motion in icing modeling has been discussed. It turns out that for cylinder icing ICEMOD2.1
SC R
and ICEMOD2 give the same ice load and thickness distribution. The only noticeable difference in the two model results is the brine film horizontal velocity, but only for an increased spray mass flux (of order
). For normal spray mass flux a
NU
simplified time-dependent model version with brine film velocity put equal to zero (Model01)
MA
gives similar ice accretion results as ICEMOD2.1. Even a stationary version of the latter model (Model02) gives adequate results for normal spray conditions. For the increased spray mass flux some brine may however be pushed over to the lee side by wind stress, especially if
TE
D
the cylinder is slim such as a truss, and only the complete versions of ICEMOD2.1 or ICEMOD2 can handle this situation. Model01 will still give correct ice thickness results on
CE P
the major part of the upwind side of the cylinder while Model02 only give some icing close to the cylinder shoulders (which means that this model cannot approximate periodic icing in this
AC
case).
KEY WORDS: icing modeling, rig icing, momentum equation, brine film motion
NOMENCLATURE R = instantaneous local icing intensity (kg m-2 s-1) Rw = instantaneous local water catch rate (kg m-2 s-1) = time-average spray mass flux “far” from the icing object (kg m-2 s-1) = spray mass flux of a single spray “far” from the icing object (kg m-2 s-1) Q = unit heat flux (W m-2) U = wind or spray velocity (m/s)
ACCEPTED MANUSCRIPT 3
T = wave period (s) S = salinity (‰)
T
D = initial cylinder diameter (m)
IP
Vh = free horizontal spray velocity (m/s)
P = pressure (
SC R
N = spray frequency (1/s) )
C = specific heat capacity at constant pressure (
)
NU
=mean brine film velocity with components in s and ς direction
MA
t = time (s) f = humidity
)
D
g = abs. value of acceleration of free fall (
TE
l0 = specific latent heat of fusion of pure ice (J/kg) le = specific latent heat of evaporation (J/kg)
CE P
= sand grain roughness height (m)
Λ = liquid water content (LWC) (kg m-3 ) )
AC
θ = temperature (
ρ = density (kg m-3)
δ=brine film thickness (m) μ = dynamic viscosity (kg m-1 s-1) τ = stress (
)
λ = wave length (m) σ = interfacial distribution coefficient χ=mass diffusivity of diluted salts (
)
ὐ=volume fraction = mass flux of vapor to/from the brine film surface (kg m-2 s-1)
ACCEPTED MANUSCRIPT 4
K, N = non-dimensional functions in the brine film momentum equations. E = non-dimensional extrapolation factor (defined in Eq. 50)
IP
T
Subscripts
SC R
C = collision-generated i = ice b = brine film mean value in normal direction
NU
a = air
MA
e = evaporation or condensation w = sea surface
D
s = significant or s-direction
TE
sp = spray
10 = 10 m above msl (mean sea level)
CE P
0 = wave crest (unless otherwise stated)
AC
50 = median volume value
ACCEPTED MANUSCRIPT
CE P
TE
SEDCO DRILLING RIG (from Google)
D
MA
NU
SC R
IP
T
5
1. INTRODUCTION
AC
In the 1980’s a numerical marine icing package called” ICEMOD” (Horjen and Vefsnmo, 1987, Horjen,1990) was developed at the Norwegian Hydrotechnical Laboratory (NHL). The work was going on for about ten years. This model was one-dimensional in the sense that the brine film movement was restricted to only one space dimension along the icing surface. For example, for a vertical rig column icing intensity could only be correctly calculated at the stagnation line where the brine film movement is solely downwards affected only by gravity. In order to calculate ice mass on the whole column a specified shape of the ice profile must then be postulated, for example that the ice thickness in the wind direction is constant along the whole upwind side of the column. Maximum radial ice thickness will then be at the stagnation line. This model was for example used for analyzing several drilling rig icing cases
ACCEPTED MANUSCRIPT 6
in North American waters (Brown and Horjen, 1989, Brown et al, 1988). ICEMOD was in this report compared with the Canadian model RIGICE. Ice thickness was estimated for three
T
of the platforms (SEDCO 708, SEDCO 709 and Rowan Gorilla) based on photographs.
IP
Neither ICEMOD nor RIGICE gave very convincing results for all these cases.
SC R
Recently an extension of ICEMOD to two-dimensional brine film motion was presented and applied for icing on cylinders of diameters 5-50 cm placed on a vessel (Horjen, 2013). The
NU
brine film velocity components were here found directly from the brine film momentum equations with inertia terms neglected. An interesting feature of the new model (called
MA
ICEMOD2) is the possibility of depression of the ice profile at the stagnation line. If the freezing fraction is high, i.e. close to 1, maximum ice thickness will be at the stagnation line
D
where the spray mass flux is also a maximum. For sufficient wet icing, however, the ice
TE
thickness profile on the upwind side of the cylinder is to a greater extent controlled by the heat supply to the brine film from the spray cloud and convective and evaporative heat loss.
CE P
These heat fluxes vary along the cylinder circumference. In the present article the model is extended a step further by including the equation of motion
AC
of the brine film. This version of the model will be called ICEMOD2.1. In order to simplify the model somewhat the brine film flow is assumed to be laminar. The importance, if any, of including brine film motion in a spray icing model was not discussed in Horjen (2013); in this article we will discuss this question in some detail in section 4.3. As in ICEMOD2 we will adopt the cylindrical approximation when calculating the local wind field (from potential flow theory), collection efficiency, local heat transfer coefficient and brine film movement. The model cylinder diameter increases with time assuming the accreted ice distributes uniformly on the front side of the cylinder (Horjen, 1990). Model results will be compared with observed ice thickness profiles on selected columns and trusses
ACCEPTED MANUSCRIPT 7
on an offshore drilling platform. Ice thickness growth rate depends on the icing intensity and the ice density close to the icing surface (
). Both these parameters are functions of
T
the brine film salinity. Additionally the ice density is also a function of the volume fraction , simply
IP
of air in the ice. In earlier presentations of ICEMOD/ICEMOD2 we have put
SC R
due to the lack of information of this quantity. Recently from a series of field measurements in the harbor area of Longyearbyen, Svalbard, air entrapment was for the major part of the tests (24) found to vary between 3% and 7 % on a vertical structure (Kulyakhtin et al, 2013).
NU
Icing conditions during these experiments were more or less different from field conditions on
MA
a drilling rig; the median droplet diameter was for example much smaller, in the range from 90 μm to 300 μm. Ryerson and Gow (2000) have reported a few calculated air volumes in ice
D
accreted on various part of the USCGC Midgett. The four values from ice samples collected
TE
on the bulkhead varied between 11.9 and 25.7 %. This is much higher than the results of Kulyakhtin et al, 2013. Possibly this may partly be due to some brine drainage before the
CE P
density measurements took place. For all the model results to be presented here we have simply assumed a constant air porosity of 5%. Calculated ice thickness then becomes 5.3%
AC
higher than for air-free ice.
There are no restrictions on the sign of R. If the heat input from the impinging spray cloud is very large melting of formerly accreted ice may occur, giving a negative icing intensity (a test is included in the program for this case to check if there is actually any ice to melt; if not the icing intensity is put equal to zero). This theoretically structured model is supported by several empirical parameterizations, for example the convective heat transfer coefficient, the duration of a single spray, the spray frequency, the size of the median volume droplets and the mean spray mass flux height distribution of impact-generated spray. These sub-models are important if we want to obtain icing results comparable with field results and must hence be supported by considerable
ACCEPTED MANUSCRIPT 8
scientific field experience. The main intension of this article is however to present the theoretical background of the extended model and show that calculated ice thicknesses on
T
cylinder structures are not too unrealistic compared with field results, even if the sub-models
SC R
IP
may be profitably scrutinized in future research.
NU
2. MATHEMATICS OF THE EXTENDED ICING MODEL The two-dimensional model may be used to calculate icing on various cylindrical and planar
MA
elements on both vessels and offshore rigs. The problem is two-dimensional when the icing structure is vertical or heeling, for example the bulkhead and mast of a vessel and the column,
TE
icing on rig columns and trusses.
D
diagonal truss and side wall of a rig. In this article we will however limit ourselves to study
CE P
Generally the geometry of this problem is shown in Fig.1 where the heeling angle σ is 90 for columns and typically between 50 and 70 for trusses. Icing takes place between lines parallel to the cylinder axis through the points
AC
the line through the lower point
. Origin of the coordinate system is at
at the bottom of the cylinder. If the overall collection
efficiency is equal to 1 (which we assume in this article for impact-generated droplets) will be on the y-axis.
ACCEPTED MANUSCRIPT
MA
NU
SC R
IP
T
9
TE
D
Fig. 1. Sketch of the geometry of the icing configuration showing the model’s coordinate system. On the left there is a cross-sectional view of the column and truss. On the right there
CE P
is a side-view of the cylindrical substrate.
AC
The integrated conservation equations of brine film mass, salt and momentum in s and ς direction have in fact been derived many years ago (Horjen and Vefsnmo, 1986). Since then no effort has however been made to solve the complete set of differential equations, which are: Brine mass conservation equation:
(1)
Salt mass conservation equations:
(2)
ACCEPTED MANUSCRIPT 10
Linear momentum conservation equation in s-direction:
IP
Linear momentum conservation equation in ς-direction:
T
(3)
SC R
(4)
(brine mass per unit area) and
(salt mass per unit area). Index b means
, for example salinity, across the brine film layer:
MA
mean value of a function
NU
is the mass exchange due to evaporation (negative) or condensation (positive),
D
(5)
TE
Both brine film salinity and velocity components are defined by this integral. The brine film
CE P
density and viscosity are assumed to be constant across the brine film layer. Eq. 3 and Eq. 4 are slightly modified from the original set of equations of the 1986-report due to a change in the assumed w-velocity profile across the brine film layer and some small corrections.
AC
Formerly the” icing equation” were defined to be Eq. 1 and Eq. 2; we now define the icing equations to be the whole set of the four equations 1-4 with the four dependent variables: X, Y, and
. Using the brine film temperature/salinity equilibrium conditions from Assur
(1958) and Cox and Weeks (1983,1986) Eq. 1 and Eq. 2 may be combined to give an explicit expression of the icing intensity as a function of the brine film salinity (Horjen,2013):
(6)
where
ACCEPTED MANUSCRIPT 11
is the net external heat flux (most important being heat convection, evaporation and heat is the conducive heat flux to the
T
exchange between impinging spray and brine film) and
IP
accretion( to be disregarded in this work). The heat fluxes are defined positive when heat is
SC R
transported to the brine film. 2.1 The brine film sub-model
in Eq. 3 and Eq. 4 are the brine film velocity componenets in the s and ς
NU
and
MA
direction assuming steady brine film motion ,i.e. setting the inertia terms in the momentum equations equal to zero (Horjen and Vefsnmo, 1987):
TE
D
(7)
where
AC
(9)
CE P
(8)
In Eq. 3 and Eq. 4 K and N are functions of X:
(10)
.
(11) The somewhat peculiar form of the expressions of K and N are due to the chosen second order polynomials in the n-coordinate for the velocity components (to be discussed below). Note that the expression of
(Eq. 7) is slightly different from the corresponding expression
given in Horjen (2013). In the latter article an exponential w-velocity profile is chosen and we
ACCEPTED MANUSCRIPT 12
also assumed that
, i.e. relatively low spray mass flux and/or thin brine film
layer.
direction while P is the sum of the spray and wind pressure (
T
) acting on the brine film in the s-
IP
is the sum of spray and wind stress (
) acting normal to
SC R
the brine film. g is the acceleration due to gravity (defined positive). From Horjen (1990) we have:
NU
(12)
is the shear stress exerted by the air at the brine film surface:
CE P
(13) where
is the spray bulk velocity component in the s and n
TE
direction respectively and
and
D
assume potential flow.
MA
is the wind speed above the surface boundary layer in s-direction. We will in this work
is the air/brine film friction coefficient.
AC
Using the Blasius method the theoretical local friction coefficient assuming potential flow becomes (see Schlichting, 1955):
(14)
where
(15) (16) (υ in radians)
(Reynolds number of the cylinder)
ACCEPTED MANUSCRIPT 13
The function
is a theoretical result based on the assumption of a laminar boundary layer
flow and potential, irrotational flow outside the boundary layer. The potential flow is
, giving
. The flow is called subcritical if the boundary layer is
IP
setting
T
expanded in a power series of five terms. The (theoretical) point of instability is found by
SC R
laminar up to the point of separation (as assumed above) and supercritical if the boundary layer undergoes transition and becomes turbulent before separation takes place. The transition between the two boundary layer flow regimes around a smooth cylinder takes place for a (various values are found in the literature, depending probably
NU
Reynolds number of order
MA
on test conditions such as free stream turbulent intensity). For a smooth cylinder in subcritical flow the use of Pohlhausen’s method, based on measured pressure data, results in boundary (Sobey, 2000), i.e. almost 30 degrees lower than the
D
layer separation just above
TE
theoretical result of Schlichting. In the supercritical flow regime separation is delayed up to (Incropera and deWitt, 1985) with transition for some angle
CE P
Reynolds number is very large (
. If the
) the flow becomes transcritical which means
that the boundary layer is turbulent on the major part of the cylinder circumference. For the
AC
cylinder icing cases of the SEDCO platforms the Reynolds number is larger than about , i.e. the flow is supercritical. If the cylinder is rough some part of the upwind boundary layer will be turbulent provided the cylinder Reynolds number is larger than a threshold value. Makkonen (1985) has for example found that for a cylinder of diameter D=0.15 m, Reynolds number sand grain roughness
and equivalent
the boundary layer will be turbulent for
, i.e.
almost the whole upwind side of the cylinder. For the same cylinder it was found that the whole upwind side boundary layer is laminar for a Reynolds number of
. In the
quoted work the local friction coefficient for a turbulent boundary as a function of the equivalent sand grain roughness was used (Kays and Crawford, 1980). But since
is a priori
ACCEPTED MANUSCRIPT 14
not known for rig icing cases only the laminar expression Eq.14 has for the time being been implemented in the icing model. Moreover,
for an icing surface will certainly depend on
T
both the geometry of the structure and environmental parameters. We also assume for
IP
convenience that boundary layer separation occurs at the point of theoretical instability for all
SC R
values of the Reynolds number. The local heat transfer coefficient used in the model is however based on turbulent flow theory (see section 2.2).
NU
A short outline of the derivation of Eq. 1 and 2 is given in Horjen (2013). The full derivation of Eq. 3 and 4 is rather elaborating and only a short outline will be given here as well. The
MA
starting point is the three boundary layer scalar equations of linear momentum:
D
(17)
TE
(18)
where
CE P
(19)
is the local brine film velocity with components in a local curvilinear
AC
coordinate system
, the n-axis lying along the outward normal to the icing front
(Fig.2). Index f is used for local brine film properties. Assuming circular cross section the components of acceleration of free fall are given by: (20)
ACCEPTED MANUSCRIPT
NU
SC R
IP
T
15
(21)
D
MA
Fig.2 Detail of the ice and brine film layer with control volumes V1 and V2.
TE
At the inner film boundary the brine film velocity tangent to the phase interface must be zero:
CE P
(22)
(23)
AC
Additionally the following inner boundary condition is obtained from Eq. 17 and Eq. 21:
The outer film boundary conditions may be derived using the so-called control volume method: we imagine a geometric volume V1, i.e. fixed in space, covering the outer brine film boundary. Conservation of linear momentum for the control volume letting
such that
the control volume contains the interface air/brine film then gives us three boundary conditions at (24) (25)
(see Horjen, 1990, for details):
ACCEPTED MANUSCRIPT 16
(26) Any brine film motion in normal direction has been neglected.
is the s-component of the
IP
T
shear stress in the brine film:
SC R
(27)
The control volume method was also used to define the boundary conditions for the enthalpy
NU
and salt diffusion equations from which Eq. 1 and Eq. 2 were derived. Conservation of for example mass of diluted salts in the control volumes of Fig. 2 for
and
in such
MA
a way that the interfaces air/brine film and ice/brine film is inside the control volumes results in the following two boundary conditions for the brine film salinity:
TE
D
(28)
CE P
(29)
The next step is to integrate Eq. 17-19 across the brine film layer in n-direction using the boundary conditions Eq. 21 and Eq. 23-25. After some manipulations we arrive at the
(30)
AC
following results:
(31)
where
is the gradient operator in tangential direction:
(32)
To get any further we must now specify the velocity component profiles across the brine film layer and we assume, somewhat arbitrary, a second order profile for both the u-variation and
ACCEPTED MANUSCRIPT 17
the w-variation (for a thin water film moving vertically downwards with no spray input the approximate solution of Navier-Stokes equation is also a second order velocity profile). For
T
atmospheric icing conditions Myers and Charpin (2014) has made a more thorough analysis
IP
of water film movement resulting in a second order polynomial velocity profile. In their case
SC R
water film velocity is a function of spray mass flux, heat fluxes from the icing surface and icing surface geometry. For sea spray icing the brine film velocity will also depend on brine salinity (since the four icing equations are linked together). Using now the boundary
NU
conditions Eq. 21, 22, 24 and 25 the velocity profiles to a sufficient degree of accuracy
MA
becomes: (33)
and
are the brine film velocity components at the brine film surface and
TE
where
D
(34)
AC
(35)
CE P
defined by:
is
For a vertical cylinder
is identical zero. Provided the brine film is thin, for example 0.1
mm and less, this term will be rather small for a heeling cylinder (truss). For simplicity we will in this work generally put Since now
and
. Eq. 32 and Eq. 33 may alternatively be written:
(36) (37) Putting Eq. 35 and Eq.36 in Eq. 29 and Eq. 30 we finally obtain, after some algebra, Eq. 3 and Eq. 4. Thus, ICEMOD2.1 features a rather complete hydrodynamic model of momentum
ACCEPTED MANUSCRIPT 18
conservation that includes molecular heat and salt transport in a laminar brine film layer adjacent to the icing surface (Fig. 2). Brine film salinity plays a key role in this model since
T
the icing intensity is an algebraic equation of this parameter (Eq. 6). Close to the inner and
IP
outer boundary of the brine film layer there may be large salinity gradients (Eq. 27 and Eq.
SC R
28), but since the liquid layer for most of the icing conditions studied in this article is rather thin (typically less than 0.1 mm) it is justified to use the mean salinity value normal to the icing front to calculate the icing intensity. If the mean spray mass flux is very large (examples
NU
are given in section 4.3) the brine film thickness may be as large as about 1.0 mm, and a
MA
detailed study of heat and salt transport within the brine film layer will probably be an improvement of the icing model. On the other hand an increased complexity of the model will make the computation time unacceptable long. A thick brine film means that it becomes
TE
D
more important to include the inertia terms in the momentum equations, which in itself increases the computation time significantly. Moreover, it is not certain that a more complex
CE P
model will bring about more reliable icing results, at least for field conditions. The interfacial distribution coefficient is set to a constant value
( Horjen, 2013). But
AC
it must be kept in mind that there is no physical reason that this parameter should actually be a constant. The value chosen is just a mean value of 29 Japanese, Finish and Russian observations ranging from 0.18 to 0.5 with an average deviation from the mean of 0.068. If a layer of saline water is freezing it may be shown that the interfacial distribution coefficient depends strongly on the ice growth rate (Weeks and Lofgren, 1967). If for example the growth rate increases from
to
the interfacial distribution coefficient
increases from 0.29 to 0.73 (from Weeks, 2010). In a spray icing model of Kulyakhtin, 2014, including heat conduction an expression of the ratio of ice salinity and brine film salinity at icing surface is presented which depends on a combination of the constant value used in ICEMOD2 and the value
from the theory of Weeks and Lofgren (Kulyakhtin, 2014):
ACCEPTED MANUSCRIPT 19
(38)
where
is the icing intensity caused by heat conduction to the ice only and
is the icing
IP
T
intensity caused by external heat fluxes. It is interesting that this equation indicates an
SC R
“effective interfacial coefficient” depending on icing conditions. In this work heat conduction is neglected and we have put total icing intensity equal to
. If however heat conduction is
to be included in ICEMOD2.1 (a challenge for a further improvement of the model) this term
from the total icing intensity (Eq. 6).
MA
separate the “conductive” icing intensity
NU
will be a function of brine film salinity and time. For periodic spray it is hence not possible to
Other important input parameters to the icing model include the heat transfer coefficient and
TE
2.2 The Nusselt number
D
various sea spray parameters. This will be discussed in the next two sections.
CE P
In the article of Horjen, 2013, a formula originally proposed by Lozowski et al., 1979, was used for the upwind part of cylinders:
AC
(39)
is the stagnation line Nusselt number and
is the the angle from the stagnation line in
degrees (positive in counterclockwise direction). This expression is based on experimental results of Achenbach (1977) for
. The cylinders examined in Horjen, 2013, had
diameters from 0.05 m to 0.5 m and the Reynolds number was in the range from
to
. Hence Eq. 38 is appropriate in this case. In this study the dimensions are much larger and the Reynolds number for the vertical columns of D=9.15 m was in the range from to
for the SEDCO 708 case and in the range
to
for the
SEDCO 709 case. This is outside the definition range of Eq. 38 and another Nusselt number formula should hence be applied. Based on computational fluid dynamics using the program
ACCEPTED MANUSCRIPT 20
FLUENT (ANSYS, 2009) the following formula has been proposed by Kulyakhtin and Tsarau (2014) for offshore rig columns:
SC R
IP
T
(40)
Excellent agreement (
) was found comparing this formula with the computed
results for a cylinder of diameter 90 m and wind speeds of 3, 9, 12, 15 and 22 m/s at Z=10 m.
NU
The Reynolds number was hence in the range 1.8 107-1.3 108, i.e. transcritical flow. The formula above has been adopted in ICEMOD2.1 for rig icing cases and completes the set of
from Eq. 39.
D
Eq. 38 and at
from
MA
equations contained in the model. Note that maximum Nusselt number is at
TE
2.3 Sea spray parameters for an offshore rig
CE P
For simplicity the spray bulk velocity components in Eq. 12 are calculated assuming that all the droplets in the spray cloud (both impact-generated and wind-generated) are moving in straight lines and in the same direction as the wind velocity. The spray bulk velocity
(41)
AC
components are then given by:
is the free horizontal spray velocity component normal to the cylinder axis just before impact. In the model we use the approximation
where
is the free
logarithmic wind speed at the height Z above the mean sea level. For the local water catch rate calculation we also assume that the droplets of impact-generated spray follow straight lines, i.e. the local collection efficiency is
. Assuming for
example the standard condition of section 4.2 the overall collection efficiency of median
ACCEPTED MANUSCRIPT 21
volume droplets at Z=6.75 m becomes 0.92 (at this level
from Eq.47 below).
For the much smaller droplets in wind-generated spray the droplet trajectories near the
T
cylinder are calculated from the equation of droplet motion (Horjen, 2013).
IP
Based on the field results from the drilling rig Treasure Scout operating in Norwegian waters
SC R
(Jørgensen 1984,1985) and droplet kinematics the following time-averaged spray mass flux
(42)
MA
and A and B are given by
D
where
NU
formula was presented in Horjen and Vefsnmo (1985):
TE
(non-dimensional)
CE P
is the liquid water content (
of a single spray at the wave crest and
is the
duration of a single spray. Obviously the value of A will depend on various parameters such
AC
as wind speed, sea state and platform response to incoming waves. For simplicity we will however in the model use the average of four values of A (varying between 5.16 and 6.44) obtain by adjusting the theoretical formula to field measurements: . Note that the constant A has been modified somewhat from previous works in order to agree more closely with experimental data. The water flux was registered mainly on various locations on the rig columns. The mean spray mass flux was estimated by measuring the mass of water in a spray collector after a time interval with approximate constant wind speed (one to several hours), and divide this mass with time. The spray mass flux calculated in that way will have a contribution from both impact- and windgenerated spray, but since the wind speed was not extremely large (from 5 to 9 Beaufort) we
ACCEPTED MANUSCRIPT 22
will assume that the spray mass flux given by Eq. 41 is solely impact-generated. Fig. 3 shows the vertical distribution of model spray mass flux calculated from Eq.41 for four wind speeds.
T
Note that the vertical axis is non-dimensional.
IP
When the wind speed is below 15 m/s wind-generated spray mass flux is neglected; above 15
SC R
m/s a formula given in Horjen, 1990, 2013, which is based on various measurements from the wave crest up to 7 m above the msl (Wu et al. (1984), Monaham (1968), Lai and Shemdin
NU
(1974) and Preobrazhenski (1973)) is used. Fig. 4 shows the vertical distribution of windgenerated spray mass flux for four wind speeds assuming Norwegian waters wave conditions ,
MA
(significant wave height given by Eq. 44 below). Note that above the level ξ = 1, i.e.
the spray mass flux is almost constant. This is due to the findings that droplets below a critical
D
size, depending on wind speed, stay in suspension independent of height above the sea surface
TE
(Wu et all, 1984). Comparing Fig. 3 and Fig. 4 we see that the wind-generated spray contribution to the total spray mass flux is rather small, at least for wind speeds up to 30 m/s.
CE P
For example, at the level
the mean impact-generated spray mass flux for
is about ten times as large as the wind-generated spray mass flux for
.
AC
For most rig icing cases below the main deck the omission of wind-generated spray in the model makes almost no difference, not only because the mass flux is relatively low but also because the majority of the droplets follows more or less the streamlines around the object. At great heights, for example on the derrick, the theoretical wind-generated spray mass flux may be higher than the mean impact-generated spray mass flux. The validity of the formula used for the wind-generated spray flux is however questionable above 7 m above msl. Using the spray volume concentration distribution function of Jones and Andreas (2012) Kulyakhtin and Tsarau (2014) did also conclude that wind-generated spray is unlikely to create significant contribution to icing on a rig (they tested another Norwegian model called MARICE).
ACCEPTED MANUSCRIPT 23
Due to wave washing ice may be removed mechanically under the level
of maximal run
up. According to Mitten (1994) this level may be approximated by the following formulas
IP
T
when waves are colliding with a column of diameter D:
SC R
(43)
where k is the wave number. From linear wave theory we have:
NU
(44)
MA
U10=10 m/s
IMPACT-GENERATED SPRAY MASS FLUX
D
U10=20 m/s U10=25 m/s
CE P
TE
3
2
1
0
U10=15 m/s
AC
Non-dimensional height ξ
4
0
1
2
3
4
5
Spray mass flux G (g/m2s)
Fig. 3. Impact-generated mean sea spray mass flux based on droplet kinematics and measurements on the oil rig Treasure Scout operating in Norwegian waters.
ACCEPTED MANUSCRIPT 24
U10=15 m/s
WIND-GENERATED SPRAY
U10=20 m/s
2
T
U10=25 m/s
IP
U10=30 m/s
SC R
1
0.5
0 0
100
200
300
NU
Non-dimensional height ξ
1.5
400
500
600
MA
Spray mass flux (10-3g/m2 s)
TE
D
Fig. 4. Wind-generated sea spray mass flux in Norwegian waters.
CE P
Spray frequency is an important parameter in the model. For vessel icing it has been estimated that in the mean every second significant wave produces spray ( Zakrzewski, 1987). The
AC
spray producing mechanism is however quite different on a semi-submersible oil platform. As part of the Offshore Icing project a few spray frequency registrations were made on the semisubmersible platform Treasure Scout (Jørgensen, 1984, 1985) operating in Norwegian waters. For this area the following empirical relations for significant wave height and period are found, based on numerous observations (Jørgensen, 1985 (wave height), Horjen, 1990 (wave period)): (45) (46)
ACCEPTED MANUSCRIPT 25
The 10 results for the 1984 registrations and the 4 results for the 1985 registrations are shown in Fig.5 together with the theoretical results based on the assumption that every and every
SC R
SPRAY FRQUENCY ON TREASURE SCOUT 8 7
N=1/Ts N=0.5/Ts Obs 1984 Obs 1985
NU
6
MA
5 4 3
D
Spray frquency (1/min)
IP
T
second significant wave is producing spray.
2 5
10
TE
0
15
20
25
CE P
Wind speed U10 (m/s)
Fig. 5. Impact-generated sea spray frequency N based on measurements on the oil rig
AC
Treasure Scout and theoretical results assuming that every significant spray generates spray (upper line) or every second significant spray generates spray (lower line).
There is large scatter in observed spray frequencies on Treasure Scout and no obvious relation between N and
can be read from the data in Fig. 5. A majority of the data seems however
to fall between the curves of collision frequencies of every significant and every second significant wave. In the numerical examples to follow we have chosen to assume that every significant wave produces spray on an offshore rig. A sensitivity analysis has been presented in section 4.2 showing the effect of varying the spray frequency in the model.
ACCEPTED MANUSCRIPT 26
The duration of each spray is calculated from the assumption that a single spray is generated along a distance equal to half the significant wave length (Horjen, 2013). This gives
IP
is the significant wave length.
SC R
where
T
(47)
With the assumptions made above it follows that A in Eq. 41 is mathematically a true
NU
constant, i.e. independent of wind speed, only if
. In the model the spray mass
flux of a single spray is a derived quantity obtained by dividing Eq. 41 by
. A better
MA
approach would obviously be to use directly field measurements of either liquid water content or spray mass flux of a single spray (total mass flux during a time period divided by number
TE
D
of sprays) at various heights and wind speeds.
The distribution of droplet sizes for various wind speeds and heights above mean sea level
CE P
(msl) were measured on the rig Treasure Scout as well. From these results the following relation was obtained, assuming the Best distribution (Best, 1950, 1951) of droplet sizes
(48)
AC
(Horjen and Vefsnmo, 1985):
where a = 1.31 mm and b = 0.81 mms. In order to calculate the heat exchange between the droplet cloud and the structure we must know the bulk spray temperature upon impact, which depends on the size distribution of spray droplets and the various droplets flight times
. The
calculation of impact-generated spray cooling during flight is based on the assumption that all the droplets in the spray cloud have the same size as the median volume droplet. We also assume that half of the spray cloud is moving upwards before impact without undergoing any cooling at all. The bulk impact-generated spray temperature is then the mean value of descending
droplet temperature (using a method given in Horjen, 1990) and the sea
ACCEPTED MANUSCRIPT 27
surface temperature. For simplicity we assume that this temperature does not change due to the possible contribution from wind-generated spray. Wind-generated spray is however taken
T
account of between impact-generated spray productions. For simplicity we again assume that
IP
half of the spray cloud is moving upwards without undergoing any cooling while all the
SC R
droplets in the descending part of the spray cloud attain the air temperature before impact ( partly freezing of supercooled droplets at the moment of collision is not considered in this analysis). The bulk wind-generated spray temperature upon impact is then the mean value of
NU
the sea surface and the air temperature. Probably a more realistic value is only slightly higher than the air temperature, but since the icing contribution from wind-generated spray is
MA
normally insignificant we don’t need to calculate the heat exchange due to this type of spray
D
very exact.
TE
It is reasonable to assume that a characteristic time the median volume droplet in impactgenerated spray is moving in the air is some function of the sea state. As a preliminary
CE P
hypothesis we will assume that the droplet flight time is equal to the horizontal flight time over a distance corresponding to a quarter of the significant wave length. Approximately a
(49)
AC
characteristic flight time then becomes:
This is of course a crude approximation since in reality the various droplets in a spray cloud have been in the air for quite different times before they collide with the structure.
3. THE NUMERICAL METHODS OF THE EXTENDED MODEL A finite difference method is used to solve the icing equations. The method chosen here is an extension of the two-step Lax-Wendroff scheme (see for example Lax and Wendroff, 1964)
ACCEPTED MANUSCRIPT 28
as proposed by Richtmyer (1963). The method is simple enough to be run on a personal computer. The icing surface is divided into nodes (I,J) representing the s and ζ-position. We and space steps
in s and ζ-direction
and
T
must choose appropriate time step
IP
respectively. In the original Lax-Wendroff and Richtmyer numerical scheme there are mid-
SC R
point evaluations. These have been avoided in the ICEMOD2.1 scheme by using the following approximations for an arbitrary quantity A:
NU
(49)
MA
Similar approximations are used for the J-variations. The same method was used in Horjen (2013) to solve two-dimensional icing on cylinders (D = 0.05-0.50m) using only the two first
TE
D
icing equations (Eq. 1 and Eq. 2). Stable solutions were then obtained using an angle step of and, for most of the cases studied, a time step of 0.016 s. With four equations the
CE P
numerical scheme becomes twice as large. For stable solutions the angle step for some of the cases may also now be as large as
, but the time step must at least be one order smaller
AC
than the time step mentioned above. It is hence obvious that the computer time using ICEMOD2.1 increases substantially compared with ICEMOD2. For the large column of the SEDCO platform (D = 9.15 m) the space steps were ( put
) and
. The mesh length ratio
and (D = 1,22 m) and
is hence 0.70. For the trusses we
which give the following model space steps: ,
,
(D = 1,52 m). The mesh length
ratio is now 1.11 and 1.09 respectively. Values of the space steps were for most of the cases studied here chosen in such a way that the mesh length ratio did not become very much different from 1. The effect of changing space steps for a particular icing case is discussed in section 4.2. Stable solutions for the whole icing periods of the SEDCO rigs were obtained for
ACCEPTED MANUSCRIPT 29
the large column when the time step was 0.001 s and for the truss when the time step was 0.0005 s.
T
The equations of droplet motion are non-linear, ordinary differential equation. A different
IP
numerical procedure must then be used to obtain a solution. In this article a fourth order
SC R
Runge-Kutta method, combined with the Adams-Moulton method, is applied (see Cheney and Kincaid , 1985).
NU
In the next chapter we will first analyze some real icing cases on two drilling rigs. Sensitivity tests will then be presented to examine the effects of some of the input parameters of the icing
D
MA
model. Finally the influence of brine film motion upon icing calculations has been discussed.
CE P
4.1 Two rig icing cases
TE
4. NUMERICAL RESULTS AND DISCUSSION
Unfortunately there are only a few field results which are suitable for model comparison. As
AC
part of a joint Canadian/Norwegian icing project seven of the best documented rig icing events between 1980 and 1990 in North American waters were studied (Brown and Horjen,1989). But even for these cases the measurements are not very exact; estimates of ice thicknesses are based on photos and observations from the deck. We will now compare ICEMOD2.1 results with two of the seven documented icing measurements: the semi-submersible offshore rigs SEDCO 708 and SEDCO 709. These two events were also presented in Brown et all, 1988. Weather data for the two icing periods are shown in an appendix. For simplicity we have in the model assumed that spray and wind come from the same direction during the whole icing period, normal to the cylinder axis (as in Fig. 1). In order to reduce computer time, model calculations are made for the first
ACCEPTED MANUSCRIPT 30
spray periods and ice thickness and accreted mass are then extrapolated to the time ) by multiplying the results with the factor
T
.
IP
(50)
SC R
SEDCO 708, January 4-6, 1983
This icing event was documented by Minsk (1984) when the rig was operating in the Gulf of
NU
Alaska at 56.3 N and 161.8 W. Ice thickness profiles on a windward side of a column and diagonal truss were estimated by David Minsk from photos (the report contains some of these
MA
photos). The photos were taken between 1310 LST (Local Standard Time) and 1450 LST on January 6 and model input data were hence terminated at 1500 LST that day. Totally 60
D
hourly input data set were used starting at 0400 LST on January 4. The following
TE
meteorological and oceanographic registrations were made: air, dew-point and sea surface
CE P
temperature, sea surface salinity, air pressure (mb), significant wave height and period, wind direction, visibility and type and intensity of precipitation. The last four registrations were not used in model simulations. Calculation for a column and a truss were made from Z=12 m
AC
down to the nearest integer of
.
1. Large column D=9.15 m Time evolution of horizontal thickness profiles at Z=5.0 m for observation 36, 39…57, 60 are shown in Fig. 6. We observe that maximum ice thickness for the last observation is about 54 from the stagnation line. Fig. 7 shows the time evolution of vertical ice thickness distribution along this angle assuming the air entrapment is 5%. The model ice density for each separate 1h icing period for the 25 observations from obs.36 to obs.60 varied from 905
to 923
. This may be compared with the measured ice densities of accreted ice on the bulkhead and the deck of a vessel reported by Ryerson and Gow (2000). They found that the
ACCEPTED MANUSCRIPT 31
density was in the range 690-920 kg/m3. Assuming no air is entrapped in the ice and brine drainage has not taken place the density will normally be above 950 kg/m3. Due to the
T
uncertainty of the air volume it would be better test the model against accreted mass instead
SC R
IP
of ice thickness. But this is of course difficult for a single structure on a rig.
SEDCO 708, VERTICAL COLUMN, D=9.15 m, ICEMOD2.1
Ice thickness (cm)
8
MA
6
NU
10
-72
-54
-36
-18
2
0 0
CE P
-90
TE
D
4
18
36
54
72
90
Angle from stagnation line (deg)
AC
Fig. 6. Model ice thickness cumulative horizontal profiles for obs. 36,39,42…..57,60 (3h between each observation) at the height Z=5.0 m.
ACCEPTED MANUSCRIPT 32
SEDCO 708, VERTICAL COLUMN, D=9.15 m, ICEMOD2.1 12
T IP
10 9
SC R
8 7 6 5 4 1
2
3
4
5
6
7
8
9
10
MA
0
NU
Height above msl (m)
11
Ice thickness 54o from stagnation line (cm)
TE
D
Fig. 7. Model ice thickness cumulative vertical profiles -54⁰ from the stagnation line for obs.
CE P
36,39,42…..57,60 (3h between each observation).
In Fig.8 model results just after observation 60 are compared with the estimated results from
AC
Minsk (1984) for the large column. We observe that between 5 and 8 m above msl model ice thickness is 4-6 cm larger than reported values. But this is not too bad considering the rather long icing period and the uncertainty of several of the input parameters, for example spray mean mass flux, spray frequency, duration of one spray, spray flight time and heat transfer coefficient. The discrepancy between model and observed results documented in Fig. 8 may partly be explained by a wrong height distribution of the spray mass flux. We must however expect that model ice thickness at a certain location is larger than observed ice thickness due to the assumption of constant wind direction. Ideally the icing program should also be run for the whole icing period of 60 hours and not only for 10 spray periods (typically 1-1.5 min) for
ACCEPTED MANUSCRIPT 33
each observation. This is however too time-consuming for a PC. The field measurements are also not very exact and some ice may have fallen off during the icing period.
IP
12
ICEMOD2.1, obs 60 Estimated
SC R
11 10 9 8
NU
7 6 5
MA
Height above msl (m)
T
SEDCO 708, VERTICAL COLUMN, D=9.15 m
4 3 2 2
3
4
D
1
5
6
Ice thickness (cm)
7
8
9
10
TE
0
CE P
Fig. 8. Comparison between maximum model ice thickness at each level and estimated (from photos) maximum thickness just after observation no 60.
AC
2. Diagonal (70 ) truss
Ice thickness profiles at Z=4.2 m in the AB-plane (Fig.1) for obs. 33,36,…57,60 are shown in Fig.9. We now observe that maximum ice thickness for the last observation is 36-54 from the stagnation line. Fig.10 shows the time evolution of vertical ice thickness distribution along the angle 36 for the same set of observations. In Fig.11 model results for
just after
observation 60 are compared with the estimated results (also from Minsk, 1984) for a diagonal truss. Since the type of diagonal was not specified model results for both the 70 (D = 1.22 m) and the 50 (D = 1.52 m) truss are included in Fig. 11. The differences between the model results are very small. We observe that above 4.5 m the model thickness is larger than the estimated thickness (0 - 6 cm).
ACCEPTED MANUSCRIPT 34
SEDCO 708, 70o DIAGONAL TRUSS , D=1.22 m, ICEMOD2.1 14
T
Ice thickness (cm)
12
IP
10
SC R
8 6
NU
4 2
-72
-54
-36
MA
0
-90
-18
0
18
36
54
72
90
D
Angle from stagnation line (deg)
TE
Fig. 9. Ice thickness cumulative profiles in AB-plane for obs.36,39,….57,60 at Z=4.2 m.
CE P
SEDCO 708, 70o DIAGONAL TRUSS , D=1.22 m, ICEMOD2.1
10
AC
Height above msl (m)
12
8
6
4
2 0
2
4
6
8
10
12
14
Ice thickness (cm)
Fig. 10. Model ice thickness cumulative vertical profiles for obs. 36,39,42…..57,60 at φ= -36
ACCEPTED MANUSCRIPT 35
SEDCO 708, TRUSSES D=1.22 m AND D=1.52 m
D=1.22 m
D=1.52 m
12
Estimated
11
T IP
9 8
SC R
7 6 5 4
3 2 3
6
9
12
15
18
21
MA
0
NU
Height above msl (m)
10
Ice thickness (cm)
D
Fig. 11. Vertical profile of model maximum ice thickness on D=1.22 m (υ=70 ) and D=1.52
CE P
TE
m (υ=50 ) trusses and estimated (from photos) maximum ice thickness just after obs. 60.
AC
SEDCO 709, January 24-26, 1986 This icing event was documented by the drilling platform ice and weather observer (MacLaren Planseach,1986) when the rig was operating on the Scotian Shelf at 43.7 N and 59.9 W. Ice thickness profiles on a windward column and diagonal truss were estimated from photos. A number of photos were taken at 1800 GMT on January 25 when maximum ice build-up was reported. Model calculations are based on 25 input data sets from 1400 GMT on January 24 to 1600 GMT on January 25. Ice accretion before and after this time interval is negligible. Calculations for a column and a truss were made from the nearest integer to to Z=8.0 m. 1. Large column, D = 9.15 m
up
ACCEPTED MANUSCRIPT 36
Time evolution of horizontal thickness profiles at Z=2.5 m for observation 3, 5…23,25 are shown in Fig.12 We observe that maximum ice thickness is now about 36 from the
T
stagnation line. Fig.13 shows the time evolution of vertical ice thickness distribution along
IP
this angle (obs. 1 has here been added to show the effect of wave washing). In Fig.14 model
SC R
results just after observation 25 are compared with the estimated results from photos. We now observe that model ice thickness is larger than estimated values for
(0-6.5 m).
NU
SEDCO 709, LARGE COLUMN, D=9.15 m, Z=2.5 m, ICEMOD 2.1 12
Ice thickness (cm)
MA
10
-72
-54
AC
-90
CE P
TE
D
8 6
4 2 0
-36
-18
0
18
36
54
72
Angle from stagnation line (deg)
Fig. 12. Model ice thickness cumulative horizontal profiles for obs. 3,5,7…..23,25 (2h between each observation).
90
ACCEPTED MANUSCRIPT 37
SEDCO 709, LARGE COLUMN, D=9.15 m, ICEMOD2.1 5
T
4.5
IP SC R
3.5
3 2.5 2
NU
Height above msl (m)
4
1 0
2
4
MA
1.5
6
8
10
12
D
Ice thickness (cm)
TE
Fig. 13. Model ice thickness cumulative vertical profiles -36⁰ from the stagnation line for obs.
CE P
1,3,5,7…..23,25 (2h between each observation).
SEDCO 709, LARGE COLUMN, D=9.15 m
7
Height above msl
6
ICEMOD2.1, obs 25 Estimated
AC
8
5 4 3 2 1 0 0
1
2
3
4
5
6
7
Ice thickness (cm)
8
9
10
11
12
ACCEPTED MANUSCRIPT 38
Fig. 14. Comparison between maximum model ice thickness and estimated (from photos)
IP
T
maximum ice thickness just after observation 25.
SC R
2. Diagonal (70 ) truss
Ice thickness profiles at Z=2.0 m in the AB-plane (Fig.1) for obs. 3,5,…23,25 are shown in Fig. 15. We now observe that maximum ice thickness is about 36 from the stagnation line. for obs.
NU
Fig.16 shows the time evolution of vertical ice thickness distribution for
MA
3,5….23,25. In this case only the maximum thickness was documented from photos: 15.0 cm at the height 2.0 m above msl. Maximum height of ice accretion was reported to be 7.0 m.
TE
D
SEDCO 709, 70⁰ DIAGONAL TRUSS, D=1.22 m, Z=2.0 m, ICEMOD2.1 16
12 10 8
AC
Ice thickness (cm)
CE P
14
6 4 2 0
-90
-72
-54
-36
-18
0
18
36
54
72
Angle from stagnation line (deg)
Fig. 15. Model ice thickness cumulative profiles in the AB-plane for obs. 3,5,7,…..23,25
90
ACCEPTED MANUSCRIPT 39
SEDCO 709, 70⁰ DIAGONAL TRUSS, D=1.22 m,
T IP
4
SC R
Height above msl (m)
5
3
NU
2
0
2
MA
1 4
6
8
10
12
14
16
D
Ice thickness (cm)
CE P
TE
Fig. 16. Model ice thickness cumulative vertical profiles for obs. 1,3,5,7…23,25 at φ= -36 .
Icing event
AC
A summary of the two icing events presented above is shown in Table 1. Icing object
Corner column,
Max. thickness
Height of max.
Max. height of
(cm)
thickness (m)
icing (m)
Model
Model
Model
Obs.
Obs.
Obs.
9.7
7.5
5.0
4.0
~12.0
9.0
13.1
18.0*)
4.0-4.4
3.5*)
~12.0
9.0
13.0
18.0*)
3.6-3.8
3.5*)
~12.0
9.0
SEDCO D=9.15 m 708
70⁰ Truss, D=1.22 m 50⁰ Truss, D=1.52 m
ACCEPTED MANUSCRIPT SEDCO
Corner column,
709
D=9.15 m
10.0
2.5
3.5
~5.5
7.0
15.3
15.0
2.0
2.0
~5.4
7.0
15.2
-
2.0
~5.4
-
T
70⁰ Truss,
11.5
IP
40
D=1.22 m
-
SC R
50⁰ Truss, D=1.52 m
Table 1. Comparison of observed and predicted ice accretion thickness data. *) Type of
MA
NU
truss for ice thickness observations was not specified.
Although this cannot be compared with direct observations a summary of model ice mass on a
TE
D
column and a truss and corresponding “overall freezing fraction” (ratio of total accreted ice mass and total mass of spray hitting the cylinder during all the observations, excluding for
Icing event
CE P
each observation icing and water catch rate below the level Icing object
Corner column,
Overall freezing fraction
28451
0.28
70⁰ Truss, D=1.22 m
7832
0.28
50⁰ Truss, D=1.52 m
11170
0.27
16618
0.66
70⁰ Truss, D=1.22 m
3595
0.65
50⁰ Truss, D=1.52 m
5089
0.65
AC
SEDCO 708
Ice load (kg)
) is given in Table 2.
D=9.15 m
SEDCO 709 Corner column, D=9.15 m
Table 2. Predicted final ice mass and overall freezing fraction.
ACCEPTED MANUSCRIPT 41
4.2 Sensitivity analysis Sensitivity analyses will now be presented for three of the icing model parameters: air
T
temperature, wind speed at 10 m height and spray frequency. The constant and the variable
Constant parameters Variable parameters ) sensitivity
NU
1: Air temperature (
) sensitivity
(cont. spray)
D
3: Spray frequency (N) sensitivity
MA
2: Wind speed (
SC R
Case number
IP
parameters for each test are given in Table 3.
TE
Table 3. Parameter values to be used in the sensitivity tests. Bold types in the right
CE P
column is standard condition.
AC
We will consider icing on a column of diameter 5.0 m extending between the levels Z = 6.0 m and Z = 7.5 m. The following parameters are kept constant: sea surface salinity
,
relative humidity f = 0.8, atmospheric pressure p = 1000 mb and significant wave height Hs = 5.0 m. The significant will be calculated from a formula given in Horjen (1990): (51) This formula is based on 240 wave data (
collected on the rigs Bow drill 1 and 2
and SEDCO 709 operating on the Grand Banks area off the coast of Newfoundland. With the prescribed value of the significant wave height the spray period becomes 7.4 s. The angle and height step is 4.5 and 0.25 m respectively and the time step is 0.0005s. Model calculations
ACCEPTED MANUSCRIPT 42
are again made for the first
spray periods and ice thickness and accreted mass are
then extrapolated to 24 hours (
) by multiplying the results by the factor E
T
(Eq. 50). Wave washing of accreted ice is now neglected. The overall freezing fraction is then
IP
given by:
SC R
(52)
) and the icing area of the
NU
The integration is over the period of icing calculations ( cylinder.
MA
Vertical profiles at the stagnation line for Case no 1 (air temperature sensitivity) is shown in Fig.17. As expected more ice accretes when the air temperature decreases. This is mainly due
D
to the increase in convective and evaporative heat loss and a decrease in arriving spray
to
CE P
cylinder varies from
TE
temperature (the bulk impact-generated spray temperature upon impact at the top of the for Case no 1). An inspection of Fig. 17 gives an
approximate maximum sensitivity of 6.58 cm per 4 the range of -12
to -8
, i.e. 1.6 cm/ , for air temperature in
. Fig. 18 (vertical profiles) and Fig. 19 (horizontal profiles) show
AC
the effect of varying the wind speed. In our case ice accretion increases with increasing wind speed (but this is not generally valid). At the cylinder base the ice thickness for this particular test increases approximately 0.5 cm when the wind speed increases with 1 m/s in the range 17-25m/s. Although maximum Nusselt number is located at
(Eq. 39) maximum
ice accretion at the cylinder base is from Fig.19 about 30⁰ from the stagnation line. This is due to several mechanisms, for example a relative high heat input from the spray cloud (to be discussed in more detail in the next section) and the interaction between icing intensity and brine film salinity. Note that there are some small bumps around the locations of maximum heat transfer coefficient (but this may just be a coincidence or due to numerical noise).
ACCEPTED MANUSCRIPT 43
-8 deg
CASE 1: AIR TEMPERATURE SENSITIVITY
-10 deg 7.5
-12 deg
T IP
7
SC R
6.75 6.5 6.25 6 0
4
8
NU
Height above msl (m)
7.25
12
20
MA
Ice thickness (cm)
16
TE
D
Fig. 17. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 1.
17 m/s
CASE 2: WIND SPEED SENSITIVITY
CE P
7.5 7.25 7
25 m/s
AC
Height above msl (m)
21 m/s
6.75 6.5 6.25
6 0
3
6
9
12
15
18
Ice thickness (cm)
Fig. 18. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 2.
ACCEPTED MANUSCRIPT 44
17 m/s
CASE 2: WIND SPEED SENSITIVITY, Z=6.0 m
21 m/s 24
T
25 m/s
IP
Ice thickness (cm)
18
SC R
12
0 -90
-60
-30
0
NU
6
30
60
90
MA
Angle from stagnation line (⁰)
TE
D
Fig. 19. Horizontal ice thickness profile at the stagnation line extrapolated to 24h, Case no 2.
CE P
From the discussion in section 2.4 the observations of spray frequency on an offshore rig is rather deficient. A sensitivity analysis varying this parameter is hence included assuming
AC
standard condition. Somewhat surprisingly at the first sight it turns out that for the three frequencies studied,
( continuous spray, i.e. single
sprays that follow immediately after each other) the total 24 h accreted mass of ice are more or less equal: 1346 kg, 1286 kg and 1463 kg respectively. Overall freezing fraction was 76%, 71% and 83% respectively. The stagnation line ice thickness profiles (Fig. 20) are also almost the same for the three cases except at the upper level of the cylinder; here the ice thickness for continuous spray is much larger than for the other two cases. An attempt to explain this result is made below: At the top
of the cylinder one of the boundary conditions of the model is X=0. It
ACCEPTED MANUSCRIPT 45
then follows that the ice thickness at some point on the cylinder top is given by (Horjen,2013):
IP
) is the ice density close to the icing surface. From this equation it follows that
SC R
where
T
(53)
for periodic spray the ice will grow only during spray impingement (normally
). For
the Case 3 test the spray period was 8.30 s and 16.60 s for spray frequencies
. This
NU
respectively. The single spray duration was from Eq. 46
and
MA
means that for the major part of the icing period there is no growth in ice thickness at . During continuous spray, however, ice will accumulate during the whole icing (during spray
D
period with an icing intensity at the stagnation point of
and
TE
impingement the icing intensity for the high and low spray frequency is
CE P
respectively).
Contineous spray N=1/Ts
CASE 3: SPRAY FREQUENCY SENSITIVITY
N=0.5/Ts
AC
7.5
Height above msl (m)
7.25
7
6.75 6.5 6.25 6 0
4
8
12
16
Ice thickness (cm)
. Fig. 20. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 3. .
ACCEPTED MANUSCRIPT 46
From Fig. 20 it follows that the time integral of R extrapolated to the time
is almost
independent of spray frequency. An expression of the extrapolated integral at the stagnation ) of spray
T
line may be obtained from Eq.1 for periodic spray and an integral number (
IP
periods (formally using Fourier analysis for the local water catch rate):
is the heat loss by evaporation,
NU
where
SC R
(54)
is the extrapolated time,
does not change with spray
MA
defined by Eq.50. In the model the mean spray mass flux
and E is
frequency (since A in Eq. 41 is put equal to a constant). The first term on the right hand side
D
of Eq. 54 is hence independent of spray frequency. If the local brine film thickness after
TE
sprays is rather thin (let’s say smaller than 0.1 mm) this term will be the dominant one. For
CE P
Case 3 the brine film thickness after 10 spray periods at the mean height Z=6.75 m is 38μm for spray generated by both every significant wave and every second significant wave. For the
AC
continuous spray case the brine film thickness is literally zero. From Eq.54 it is likely that the accreted ice mass depends more on spray frequency for a thick brine film. To test this hypothesis we will run Case no 3 once more, but now with a mean spray mass flux 25 times larger than the value given by Eq. 41. Stagnation line ice thickness profiles are shown in Fig. 21 (for the continuous spray case there is no ice deposit here). Ice thickness now decreases about 50% when the spray frequency decreases to half the original value. Total accreted ice mass decreases 17%. At Z=6.75 m the brine film thickness after 10 spray periods is about 360 μm for the high frequency case ( the low frequency case (
).
) and about 460 μm for
ACCEPTED MANUSCRIPT 47
N=1/Ts
CASE 3: SPRAY FREQUENCY SENSITIVITY, 25*NSF
N=0.5/Ts
7.5
T IP
7
SC R
6.75 6.5 6.25 6 0
2
4
6
NU
Height above msl (m)
7.25
8
10
12
14
MA
Ice thickness (cm)
D
Fig. 21. Vertical ice thickness profile at the stagnation line extrapolated to 24h when the
CE P
TE
spray mass flux is 25 times the “normal” value (NSF) given by Eq.41 .
As mentioned before a more realistic spray model would be the use of the spray mass flux
AC
or the LWC of a single spray as input to the icing model. The time-average impact-generated spray mass flux is then generally given by (compare Eq. 41): (55) where
is a function depending on the sea state, height above msl and type of
marine structure. Assuming that the spray frequency is proportional to the meeting frequency of significant waves (actually this is most probable for a fixed structure) it then follows that the first term on the right hand side of Eq. 54 is independent of spray frequency only if where generally
. Expressing for example the significant
wave height and period by Eq. 44 and Eq. 45 this condition gives a LWC near the wave crest
ACCEPTED MANUSCRIPT 48
proportional to
. This is a rather special condition which we cannot expect to be valid.
Nevertheless the assumption made in the model that the mean spray mass flux does not
T
change with spray frequency has proven to give reasonable approximations of ice growth rate
IP
on some marine structures.
SC R
Details of the time evolution of icing intensity at the mean height and stagnation line of the cylinder are shown in Fig.22. We observe that the variations of icing intensity are quite
NU
different for the three spray frequency cases, but the time-average value of the icing intensity is more or less equal. Note that the left hand side of Eq. 54 may be written
, which
MA
again shows that total ice load per unit area for a specified icing period do not vary much with
D
spray frequency.
N=1/Ts
TE
ICING INTENSITY VARIATIONS FOR CASE 3, Z=6.75 m
CE P
Continuous spray
1.9
AC
Icing intensity (g/m 2 s)
2.4
N=0.5/Ts
1.4
0.9
0.4 0
10
20
30
40
50
Time (s)
Fig. 22. Detail of icing intensity time variations at the stagnation line and mean height of the cylinder. Six and three spray cycles are considered for the two periodic spray cases. The spray-on periods are the horizontal or monotonically increasing parts of the curves.
ACCEPTED MANUSCRIPT 49
For the continuous spray case at the level Z = 6.75 m the stagnation point spray mass flux is and the icing intensity is
(see Fig.22). The , which means that the assumption is justified
IP
that droplet salinity upon impact is equal to the sea surface salinity
T
median droplet diameter is from Eq. 47
SC R
(Horjen,1990). Neglecting the convective transport of the brine film and mass loss by evaporation the theoretical brine film salinity is simply given by (see section 4.3):
NU
(56)
close to the ICEMOD2.1 result of
, which is
MA
Putting now the parameters given above in this equation we find .
D
There are also other model parameters which are based on very limited or none experimental
TE
documentations, for example the droplet flight time. Tests have been done making the flight
CE P
time equal to half the wave length and the eighth part of the wave length. Total ice load for the standard condition then becomes 1351 kg and 1343 kg respectively, i.e. the variations with flight time seems to be small for normal sea spray flux (NSP) given by Eq. 41. The bulk
2
AC
impact-generated spray temperatures upon impact at the top of the cylinder were
and
for the two flight times considered.
Finally some diagrams showing the time evolution of icing intensity, brine film salinity, brine film thickness and cumulative accreted ice mass per unit area at the stagnation line and mean height Z=6.75 m of the cylinder are presented in Fig. 23-26 for Case no 1. The diagrams show the variation of the four quantities during the first 10 or 20 spray cycles for wind speed 21 m/s and air temperatures -8.0
and -12.0
. The spray period is 8.30 s (from Eq.45) and, as
mentioned above, the duration of a single spray is 2.56 s. Overall freezing fraction for -8.0 and -12.0
is 67% and 83% respectively.
ACCEPTED MANUSCRIPT 50
From Fig. 24 we see that after some few spray cycles the brine film salinity more or less fluctuates around 54
before any drainage has taken place. This result may be compared with the
T
about 18
when the air temperature is -8 . Corresponding ice salinities will be
IP
vessel icing investigations of Ryerson and Gow (2000). They found that the salinity of ice
7.5
SC R
accreted on the bulkhead varied between 11.5‰ and 14.7‰. The air temperature was from to -8.5 . Time evolution of the brine film thickness is shown in Fig.25. This is the only
example with 20 spray periods. The larger number of spray periods is chosen in order to show
NU
that the brine film thickness normally stabilize after a sufficient time period. In this case we
43μm for
MA
see that after some time the brine film thickness seems to fluctuate around 63μm and and
respectively. at the
D
Fig. 26 show the time evolution of the cumulative ice mass accretion per unit area
TE
stagnation point at Z=6.75 m (integration of the R-function in Fig.23). Due to the harmonic
CE P
wave nature of the icing intensity there is a periodic-linear increase with time of the accumulated ice mass. In this example
is practically proportional with time for the whole
time period considered (10 spray periods) which is actually a condition for the extrapolation
AC
method to give approximate correct answers.
ACCEPTED MANUSCRIPT 51
-8 deg
ICING INTENSITY VARIATION FOR CASE 1, Z=6.75 m
-12 deg
T IP
2
SC R
1.5 1 0.5
0 0
10
20
30
40
50
60
70
80
90
MA
Time (s)
NU
Icing intensity (g/m 2s)
2.5
D
Fig. 23. Time history of icing intensity at the stagnation line and mean height for the first 10
CE P
TE
spray periods for a =-8⁰C and a =-12⁰C.
-8 deg
BRINE FILM SALINITY VARIATION FOR CASE 1, Z=6.75 m 120
AC
Brine film salinity (‰)
105
-12 deg
90 75 60 45 30
0
10
20
30
40
50
60
70
80
90
Time (s)
Fig. 24. Time history of brine film salinity at the stagnation line and mean height for the first 10 spray periods for a =-8⁰C and a =-12⁰C.
ACCEPTED MANUSCRIPT 52
-8 deg
BRINE FILM THICKNESS VARIATION FOR CASE 1, Z=6.75 m
-12 deg
T
70
IP
60 50 40
SC R
Brine film thickness (μm)
80
30 20
0 0
20
40
60
80
NU
10 100
140
160
180
MA
Time (s)
120
D
Fig. 25. Time history of brine film thickness at the stagnation line and mean height for the
TE
first 20 spray periods for a =-8⁰C and a =-12⁰C
CE P
-8 deg
CUMULATIVE ICE MASS ACCRETION FOR CASE 1, Z=6.75 m 160
AC
Accreted ice mass (g/m2 )
140
-12 deg
120 100
80 60 40 20 0 0
15
30
45
60
75
90
Time (s)
Fig. 26. Time history of total accreted ice mass per unit area at the stagnation line and mean height for the first 10 spray periods for a =-8⁰C and a =-12⁰C.
ACCEPTED MANUSCRIPT 53
We will now study in more detail the mutual variation of icing intensity and brine film salinity. This is shown in Fig. 27 for the three last spray cycles when the wind speed is 21 m/s vary periodically synchronized with the
T
and the air temperature is -12⁰C. Both R and
IP
model’s spray on/off. The time period from A to B corresponds to the spray-on condition and
SC R
the time period from C to D corresponds to the model’s spray-off condition. Icing intensity is lowest during spray impingement, between for example A and B. This may be explained by the relatively large heat and mass input from the spray cloud. The temperature of the
. With the assumptions made in the model he bulk spray temperature upon impact then
MA
1.34
droplets just before impact at the upper level Z = 7.5 m of the cylinder is
NU
descending
becomes
. The overall freezing fraction for the whole cylinder
D
(Eq. 52) for the 10 spray periods is about 83% , but during spray impingement from for
TE
example A to B the local freezing fraction is only 18-20%. Brine film salinity has the lowest value at the end of a single spray impact, which is due to both mixing with sea water of lower
CE P
salinity and low icing intensity.
Just after the spray has arrived icing intensity is gradually increasing from the lowest value at
AC
A to a slightly higher value at B. During the same period the brine film salinity is decreasing. Icing intensity is gradually decreasing during the no-spray period from the highest value at C to a lower value at D. Now the brine film salinity is increasing. This reflects the interesting respond brine film salinity changes have upon icing intensity and vice versa. The reason for the abrupt changes in icing intensity immediately before A and immediately after B is the discontinuous nature of the spray mass flux.
ACCEPTED MANUSCRIPT 54
R
ICING INTENSITY AND BRINE FILM SALINITY, 3 SPRAY PERIODS
Sb 2.2
90
C
1.8
70 60
D
1.6
50 40
B
1.4
30 20
A
NU
1.2
55
60
65
70
MA
1
75
Brine film salinity (‰)
IP
T
80
SC R
Icing intensity (g/m2 s)
2
10 0 80
85
D
Time (s)
TE
Fig. 27. Time history of icing intensity and brine film salinity at the stagnation line and
CE P
Z=6.75 m. Detail of Fig. 24 and Fig.25 for the three last spray periods for a =-12⁰C. Spray
AC
on: A to B. Spray off: C to D.
For the sensitivity tests presented above the angle and lateral space steps were
and
giving a mesh length ratio of 1.27. A question arise: how sensitive is the model
for angle/space steps changes? To examine this question we shall run the standard condition test again, but now with four other sets of angle/space steps. The results are summarized in Table 4. We observe that the differences in results using the various sets of angle/space steps are small; compared with test no 1 (sensitivity test condition) maximum change in total accreted ice load is only 5.1 % (test no 4). A disadvantage of using very small angle/space steps is the increase in computer time, running for example test 4 takes for example about 5 times as long time as running test 3.
ACCEPTED MANUSCRIPT
Space
Total ice
“Mean” stag.
“Mean” front
Max. ice
Overall
number
steps
load
line thickness
thickness
thickness
freezing fraction
(cm)
4.5,
0.25
1346
14.4
2
9.0 , 0.25
1349
14.5
3
18.0, 0.25
1330
14.6
4
4,5, 0.125
1278
14.3
5
4.5,
1375
0.30
MA
1
(cm)
14.3
(cm)
(%)
10.8
19.5
76.2
10.7
19.4
76.5
10.1
18.6
75.9
10.8
19.7
77.9
10.7
19.4
75.7
SC R
(kg)
NU
(⁰,m)
IP
Test
T
55
D
Table 4. Comparison of standard condition results, extrapolated to 24h, for four different sets
TE
of space steps. “Mean” thickness is defined to be the mean thickness for the upper, middle is the overall freezing fraction defined by Eq. 52.
CE P
and lower section of the cylinder and
AC
4.3 How important is the inclusion of brine film movement in icing modeling? In this section we will take a closer look on the effect of brine film movement and also find out if the convective acceleration terms have any significance at all. We will consider icing on two vertical cylinders, the first one with the same diameter as the corner columns and the second one with the same diameter as the 70⁰ trusses of the SEDCO platforms. The following parameters are common for the two tests:
Air temperature:
Sea surface temperature:
Sea surface salinity:
ACCEPTED MANUSCRIPT 56
Relative humidity:
Significant wave height:
Air pressure:
Air porosity:
IP
T
. Icing will be calculated for
SC R
The angle step chosen is
. An
extended angle range is used to account for the possibility of brine flowing over to the lee side
NU
of the cylinder at high wind speeds. On the lee side only wind stress and gravity is acting on the brine film. Note that for a laminar boundary flow we can use the Blasius formula (Eq.14)
MA
up to 110⁰ from the stagnation line before instabilities occurs. The Nusselt number equation (Eq. 39) is also valid on some part of the lee side, up to 130 from the stagnation line. For
D
icing calculations at the cylinder shoulders and beyond it was necessary to choose a relatively
TE
small angle step. The time step chosen is the same as for the sensitivity tests (5
s).
CE P
Wave period is again calculated from Eq. 51. Like the sensitivity tests wave washing of accreted ice is neglected. Model calculations are also now made for the first 10 spray periods and ice thickness and accreted mass are then extrapolated to an icing period of 24 hour by
AC
multiplying this result by the factor E (Eq. 50). CASE no 1: D = 9.15 m, Z = 5.0-7.4 m,
, Wind speed:
At the level Z=5.00 m the instantaneous local water catch rate at the stagnation point during spray impingement is
approximately. There is no contribution from wind-
generated spray (overall collection efficiency equal to zero). The bulk spray temperature just before impact at the top of the cylinder is about
. It turns out that both ICEMOD2 and
ICEMOD2.1 give an overall freezing fraction of 37.1%. There is also no difference in total ice load and “mean” (defined under Table 4) windward side ice thickness: 2743 kg and 7.8 cm respectively. The observed wind speeds for the SEDCO708 and SEDCO709 icing cases were
ACCEPTED MANUSCRIPT 57
all less than 25 m/s. It is therefore reason to believe that the results given in section 4.1 would be the same using ICEMOD2.
T
We will now investigate a case with much larger mean water catch rate by multiplying the
IP
normal sea spray flux Eq.41 (NSF) by a factor of 50. This is not unrealistic since other data
SC R
sources of spray measurements indicates spray mass fluxes 10-1000 as large as the formula values of this work (Forest et al, 2005). The instantaneous local water catch rate at the . But even now
NU
stagnation point at Z=5.0 m during spray impingement is now
the two models give almost the same total ice load: 3755 kg and 3649 kg using ICEMOD2.1
MA
and ICEMOD2 respectively. There is however one parameter which is quite different for the two icing models for large spray mass fluxes: the horizontal brine film velocity
. This is
the ICEMOD2 horizontal velocity is about 28% higher than the ICEMOD2.1 result.
TE
υ
D
shown in Fig.28 (horizontal profiles at the cylinder base at the end of spray no 10). For
CE P
The vertical brine film velocity is however almost equal for the two models (Fig.29). Due to the increased brine film thickness at the cylinder shoulders (Fig.30) maximum absolute vertical velocity occurs at these locations (compare Eq. 8). We observe that this value is not
AC
much different for the two components. It has been common practice in icing modeling to only consider the gravity driven component of the brine mass flux (Kulyakhtin and Tsarau, 2014). Although laminar lubrication theory is an idealization this work shows that for a relatively high wind speed it is important to retain both velocity components.
ACCEPTED MANUSCRIPT 58
ICEMOD 2.1
MODEL COMPARISON CASE 1
ICEMOD 2
2400
800 0 -800
-108 -90 -72 -54 -36 -18
0
18
36
-1600
-2400
72
90 108
NU
-3200
54
IP
T
1600
SC R
Brine film velocity ub (mm/s)
3200
MA
Angle from stagnation line (⁰)
Fig. 28. Horizontal profiles of brine film mean lateral velocity across the brine film (defined
D
by Eq.5) at the cylinder base for ICEMOD2.1 and ICEMOD2. Model results refer to the end
CE P
TE
of spray no 10. Spray mass flux is 50 times as large as the normal value (NSF, Eq. 41).
ICEMOD 2.1
-108 -90
-72
AC
Brine film velocity (mm/s)
MODEL COMPARISON CASE 1 VERTICAL BRINE FILM VELOCITY, Z=5.0 m, 50*NSF
-54
-36
ICEMOD 2
0 -18
0
18
36
54
72
90
108
-400
-800
-1200
-1600
Angle from stagnation line (⁰)
Fig. 29. Horizontal profiles of brine film mean vertical velocity at the cylinder base for ICEMOD2.1 and ICEMOD2 at the end of spray no 10. Spray mass flux is 50 times as large as the normal value.
ACCEPTED MANUSCRIPT 59
MODEL COMPARISON CASE 1 BRINE FILM THICKNESS, Z=5.0 m, 50*NSF
ICEMOD 2.1 ICEMOD 2
Brine film thickness ( μm)
1000
IP
T
800
SC R
600 400
0 -108 -90
-72
-54
-36
-18
0
NU
200
18
36
54
72
90
108
MA
Angle from stagnation line (⁰)
D
Fig. 30. Brine film thickness profile at the base of the test cylinder for ICEMOD2.1 and
TE
ICEMOD2 at the end of spray no 10. Spray mass flux is 50 times as large as the normal
CE P
value.
Looking at Fig.28 and Fig.29 one may suspect that the brine film velocity has no effect at all on icing calculations. Perhaps we may as well put
? The icing equations
(57)
AC
then reduce to two simple ordinary differential equations:
Periodic spray,
A further simplification is the assumption of continuous spray. The brine film salinity is then constant in time and we just have to solve two algebraic equations (see Horjen,2013):
(58)
Continuous spray,
ACCEPTED MANUSCRIPT 60
We must here assume that the conductive heat flux to the accreted ice is independent of time, but as mentioned before we have generally put
in this work. Several other icing
T
models resemble this model, for example the old model of Stallabrass (1980) and RIGICE
IP
(Brown and Horjen, 1989). It turns out that for normal sea spray mass flux ICEMOD2.1 and
SC R
the periodic spray model with no brine film movement (called Model0l) results in no difference in the distribution of accreted ice. The continuous spray model with no brine film movement (called Model02) does however depict a somewhat different ice thickness
NU
distribution (Fig. 31 and Fig. 32). We observe that Model02 gives more ice deposit at the
MA
upper part (which again may be explained by Eq. 53) and less ice deposit at the lower part of the cylinder than the other two models. This results in a total ice load for the three models
D
differing not significantly. Maximum ice thickness for all the models is at
ICEMOD 2.1+ MODEL 01 MODEL 02
12 10 8 6
AC
CE P
TE
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, Z=6.20 m, NSF
Ice thickness (cm)
.
4 2 0
-108 -90
-72
-54
-36
-18
0
18
36
54
72
90
108
Angle from stagnation line (⁰)
Fig. 31. Ice thickness horizontal profiles at the mean height for three model versions: ICEMOD2.1, Model0l and Model02. Normal spray mass flux is assumed (Eq. 41).
ACCEPTED MANUSCRIPT 61
ICEMOD 2.1+ MODEL 01
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, ϕ=0⁰, NSF
MODEL 02
T
7.4
IP SC R
6.6 6.2 5.8 5.4
5 0
2
4
NU
Height above msl (m)
7
6
8
MA
Ice thickness (cm)
D
Fig. 32. Ice thickness vertical profiles for Case no 1 at the stagnation line for three model
CE P
TE
versions. Normal spray mass flux is assumed (Eq. 41).
Although the simplified models Model01 and Model02 in many cases approximates icing on a
AC
cylinder quite well these models bring about some logical problems about continuity principles: what happens with the unfrozen brine? Since there is no motion of the brine film we must conclude that unless brine is continuously blown off by wind action the brine film thickness will increase without limit when the freezing fraction is less than 1. The growth of brine film thickness in time for the three models ICEMOD2.1, Model01 and Model02 is shown in Fig. 33. We observe that Model02 gives a linear growth in time while Model01 gives a stepwise growth following the same trend as Model02. ICEMOD2.1 on the other hand gives a brine film thickness which after some spray cycles seems to stabilize and more or less fluctuate around a constant value. Fig. 34 shows the time history of brine film salinity for the same three models. As expected Model02 results in a constant salinity (40.3‰
freezing
ACCEPTED MANUSCRIPT 62
temperature -2.14 ) all the time while the other two models give fluctuating salinity values. Note that for periodic spray the salinity decreases during periods of spray-on and increases
SC R
IP
T
during periods of spray-off.
MODEL COMPARISON CASE 1 BRINE FILM THICKNESS, ϕ=0⁰, NSF
240
MODEL 01
NU
MODEL 02
MA
160 120 80
D
40 0 15
CE P
0
TE
Brine film thickness (μ m)
200
ICEMOD 2.1
30
45
60
75
90
Time (s)
AC
Fig. 33. Model comparison of time history of brine film thickness at the stagnation line and mean height for the first 10 spray periods. Normal sea spray mass flux.
ACCEPTED MANUSCRIPT 63
MODEL COMPARISON CASE 1. BRINE FILM SALINITY, ϕ=0O, NSF
ICEMOD 2.1 MODEL 01
44
MODEL 02
T IP
40
SC R
Salinity (‰)
42
38
34 0
15
30
45
60
75
90
MA
Time (s)
NU
36
D
Fig. 34. Model comparison of time history of brine film salinity at the stagnation line and
CE P
TE
mean cylinder height for the first 10 spray periods. Normal sea spray mass flux.
We will now make the same comparison as above for a very large spray mass flux case, 50
AC
times NSF. The result is shown in Fig.35. ICEMOD2.1 and Model01 give the same profile for . A small amount of brine is pushed over to the lee side which of course cannot be calculated by Model0l. Model02 gives quite a different result with no icing for
.
The constant mean spray mass flux is in this case too large on the major part of the cylinder upwind side to give any icing. For ICEMOD2.1 and Model0l almost all icing takes place during the spray-off part of the spray period. A summary of the results for both normal and very large spray mass flux for Case no 1 is shown in Table 5.
ACCEPTED MANUSCRIPT 64
ICEMOD 2.1
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, Z=6.2 m, 50*NSF
MODEL 01
18
MODEL 02
T
16
IP
Ice thickness (cm)
14 12
SC R
10 8 6 4
NU
2 0 -72
-54
-36
-18
0
18
36
54
72
90
108
MA
-108 -90
Angle from the stagnation line (⁰)
TE
D
Fig.35. Ice thickness horizontal profiles at the mean height for three model versions. A very
CASE no 1
CE P
large spray mass flux is assumed (50*NSF).
ICEMOD2.1,
AC
Test number (
Model0l
, Periodic spray) (
M
Model02
, Periodic spray) (
M
, Cont. spray)
M
(kg)
(cm)
(%)
(kg)
(cm)
(%)
(kg)
(cm)
(%)
2743
7.78
37.1
2744
7.79
37.1
2903
7.74
39.5
2 ( 50*NSF) 3755
11.4
1.02
3447
10.6
0.93
295.7
0.753
0.08
1 (NSF)
Table 5. Summary of Case no 1 icing results for ICEMOD2.1 and two simplified models. =”mean” front thickness defined under Table 4,
CASE no 2: D =1.20 m, Z =6.0-6.48 m,
defined by Eq. 52.
, Wind speed:
ACCEPTED MANUSCRIPT 65
At the level Z=6.00 m the instantaneous local water catch rate at the stagnation point during spray impingement is approximately
.The water catch rate has now a small
T
contribution from wind-generated spray (overall collection efficiency 4%). The bulk spray . Both
IP
temperature just before impact at the top of the cylinder is also now about
SC R
ICEMOD2 and ICEMOD2.1 give an overall freezing fraction of 51.4% and there is no difference in total ice load and “mean” windward side ice thickness: 93.0 kg and 10.3 cm.
NU
We will now study the effect of a very large spray mass flux, but this time we multiply the NSF by a factor of 25. The instantaneous local water catch rate at the stagnation point at . Any significant difference
MA
Z=6.0 m during spray impingement then becomes
in ice load calculations for the two two-dimensional models is still not obtained, 129 kg
D
(ICEMOD2.1) and 131 kg (ICEMOD2). The single main difference between the two models
TE
is again the horizontal brine film velocity (Fig.36, some numerical noise occurs). Now the is about 62% lower than the ICEMOD2.1 result.
CE P
ICEMOD2 horizontal velocity at υ
The vertical brine film velocity was however almost equal at all angels for ICEMOD2.1 and ICEMOD2 (not shown here). Comparisons of ice thickness horizontal profiles for the three
AC
models ICEMOD2.1, Model0l and Model02 are shown in Fig.37. ICEMOD2.1 and Model01 give the same profile for
. More brine is now transported to the cylinder lee side
resulting in greater differences in ice loads between the three models compared with Case no 1 (Table 6). The strange Model02 ice profiles (only icing near the cylinder shoulders) is obviously again due to the constant large heat input from the spray cloud. For such cases this model must be discarded.
ACCEPTED MANUSCRIPT 66
MODEL COMPARISON CASE 2 BRINE FILM VELOCITY, Z=6.00 m, 25*NSF
ICEMOD 2.1 ICEMOD 2
400
T
300
IP
200 100 -100 -108 -90 -72 -54 -36 -18
0
SC R
0 18
36
-200 -300
-400
54
72
90 108
NU
Brine film velocity ub (mm/s)
500
-500
MA
Angle from stagnation line (⁰)
Fig. 36. Horizontal profiles of brine film mean lateral velocity at the cylinder base for
TE
D
ICEMOD2.1 and ICEMOD2. A very large spray mass flux is assumed (25*NSF).
ICEMOD 2.1
AC
Ice thickness (cm)
CE P
MODEL COMPARISON CASE 2 ICE THICKNESS PROFILES, Z=6.24 m, 25*NSF
MODEL 01 MODEL 02
25 20 15 10
5 0 -108 -90 -72 -54 -36 -18
0
18
36
54
72
90
108
Angle from stagnation line (⁰)
Fig.37. Ice thickness horizontal profiles for Case no 2 at Z=6.24 m for three model versions. A very large spray mass flux is assumed (25*NSF).
ACCEPTED MANUSCRIPT 67
ICEMOD2.1,
Test number (
Model0l
M
, Periodic spray) (
, Cont. spray)
T
, Periodic spray) (
Model02
M
M
IP
CASE no 2
(cm)
(%)
(kg)
(cm)
(%)
(kg)
(cm)
(%)
92.95
10.3
51.4
78.09
8.74
43.2
88.09
8.80
48.9
2 (25*NSF) 128.6
12.7
2.85
98.69
11.5
2.19
21.47
2.18
0.48
NU
1 ( NSF)
SC R
(kg)
MA
Table 6. Summary of Case no 2 icing results for ICEMOD2.1 and two simplified models.
D
5. CONCLUSIONS AND RECOMMENDATIONS
TE
In the old ice model ICEMOD only one space coordinate was considered. This implied that if
CE P
you wanted to calculate ice accretion along the whole circumference of for example a vertical rig column you had to specify in advance the shape the ice profile based on icing calculation at the stagnation line only. This is no longer necessary in the new two-dimensional models
AC
ICEMOD2 and ICEMOD2.1. For relative wet icing (freezing fraction much lower than 1) the shape of the ice profile is governed by mainly three factors: spray mass flux, heat balance of the brine film layer and movement of the brine film. Typical model ice profile is then like a saddle as shown in Fig. 6 and Fig. 9 for SEDCO 708 and Fig. 12 and Fig. 15 for SEDCO 709. For the vertical columns the model horizontal profiles are as expected symmetric about the stagnation line. For a truss we should expect that for very wet icing most ice accretes at the lower side of the column. This is due to the affect of gravity on the brine film movement. For the SEDCO 708 model result (Fig. 9) we observe that ice thickness is slightly higher on the lower side of the 70⁰ truss. For the SEDCO 709 case (Fig. 15) this tendency is not observed.
ACCEPTED MANUSCRIPT 68
We must however remember that the saddle form of the ice profile is here solely a theoretical result; it is impossible from the photos included in the report of Minsk, 1984, to make such a
T
conclusion. Unfortunately no photos were presented in the first theoretical analysis of the
IP
SEDCO icing cases 27 years ago (Brown and Horjen, 1988). At that time we certainly had
SC R
access to all photo material available. If these photos still exist it would of course be most interesting to make a comparison, if possible, with theoretical results. But most likely, due for example to changes in wind direction during the icing periods, the field ice profiles had not
NU
the nice symmetric form depicted in Fig. 6, 9, 12 and 15.
MA
We observe from Table 1 that model maximum ice thicknesses is up to 29% higher or lower than estimated values while the heights of maximum ice thickness are not very much different
D
from the estimated values from photos. Further we see from Fig. 8 and Fig. 11 that above
TE
Z=4.5-5.0 m maximum model ice thickness is larger than estimated thickness for the column and trusses of SEDCO 708, largest difference being about 6.0 cm for the column and about
CE P
6.4 cm for the truss. From Fig. 14 we see that model maximum ice thickness for the column of SEDCO 709 is larger than estimated thickness when
, the largest difference being
AC
about 6.5 cm. Comparing model and estimated ice thickness at specified heights show relative large differences. Assuming however that model thickness height distribution is displaced 1-2 m downwards many of the model result and the estimated results will more or less coincide; implying that there is some chance that total model and field ice load on a single column or truss will be comparable. This has not been analyzed in this work. For the sensitivity tests a constant significant wave height of 5.0 m is assumed while the wave period is calculated from an empirical relation based on field data from the Grand Bank region (not far from the icing location of SEDCO 709). From Fig. 17 we may conclude that for the base of the cylinder (Z=6.0 m) and an icing period of 24 h the ice thickness increase due to air temperature decrease is about 1.6 cm/
. When the wind speed increases by 1m/s
ACCEPTED MANUSCRIPT 69
the ice thickness increases about 0.5 cm at the same level (Fig. 18). Note however that these conclusions are only valid for the specified test conditions of the sensitivity tests. Finally the
T
effect of varying the spray frequency was investigated. Although not very realistic we have
IP
assumed that the time-average spray mass flux formula used in the icing model is independent
SC R
of spray frequencies (spray frequency was not registered simultaneously during the collection of spray mass flux data on Treasure Scout). The single spray mass flux is then a derived quantity obtained by dividing the time-average value by the product
. Reducing now the
NU
spray frequency to half the value used for rig icing modeling of this work does only have a
MA
minor effect upon ice thickness distribution (Fig.20) and total accreted ice mass (reduced by 4.5 %) for a specified icing period. For continuous spray (implying a stationary model for the
D
icing intensity) the changes are also rather small; total ice load then increases 8.7%. This
TE
conclusion is however not valid for sea spray mass flux 25 times larger than the normal value given by Eq.41 (Fig. 21). Now the total ice load and stagnation line ice thickness
CE P
increases with increasing spray frequency. Based on a theoretical consideration given in the text this is probably due to a large increase in brine film thickness (the example used in this
AC
work gives a factor of about 10). Details of the icing intensity, brine film salinity and brine film thickness variations at and mean cylinder height (Z=7.2 m) are shown in Fig. 23-26 for the initial 10 or 20 spray periods when the air temperature is -8.0
and -12.0 . All these quantities vary periodically
in time. It is interesting to note from Fig. 23 that maximum icing intensity occurs during the “no spray” portion of the on-off spraying cycle. After the first 2-3 spray cycles the brine film salinity seems to stabilize and fluctuate around a certain value. The brine film thickness does also stabilize and fluctuate around some value, but not before a longer time period has elapsed; in our examples 14-15 spray cycles (Fig. 25). As shown in Fig. 26 there is a periodiclinear increase with time of the total accumulated ice mass per unit area.
ACCEPTED MANUSCRIPT 70
The main goal of this article was the presentation of an extended icing model, ICEMOD2.1, based on four differential equations. In the beginning of this work it was expected that the
T
new model would give noticeable different and hopefully better result than the first two-
and Case 2: D = 1.20 m,
) we have
SC R
(Case 1: D = 9.15 m,
IP
dimensional model ICEMOD2. Based on the results of section 4.3 for two rig icing cases
however arrived at the following more or less unexpected conclusions:
NU
1. Concerning total ice load and thickness distribution ICEMOD2 gives the same result as ICEMOD2.1, either for the normal sea spray mass flux (Eq. 41) or a very large sea
or 25 (applied for Case 2).
MA
spray mass flux obtained by multiplying Eq. 41 by a factor of 50 (applied for Case 1)
D
2. Only one calculated parameter turns out to be noticeably different for the two models:
TE
the brine film lateral velocity for a very large sea spray mass flux (see above). Depending on icing conditions ICEMOD2.1 gives either larger (Case 2) or lower
CE P
(Case 1) maximum lateral velocity
than ICEMOD2.
3. For a very large spray mass flux the longitudinal velocity
is largest near the
AC
cylinder shoulders where we have maximum brine film thickness. 4. For NSF an icing model with no brine film movement (called Model0l) gives the same ice load and thickness distribution results as ICEMOD2.1. 5. For NSF a stationary model with no brine film movement (called Model02) gives a satisfactory total ice load approximation to ICEMOD2.1 results. There is however now a small difference in ice thickness distribution; Model02 stagnation line ice thickness being larger than ICEMOD2.1 and Model01 results at the upper section of the cylinder (a simple theoretical proof is given) and smaller at the lower section. 6. For a very large spray mass flux (in this analysis it was of order some brine may move to the lee side of the cylinder due to the action
ACCEPTED MANUSCRIPT 71
of wind stress; especially if the cylinder is slim such as for example a truss. On the lee side icing can be calculated only by using ICEMOD2 or ICEMOD2.1. For large
T
structures (rig columns) Model0l will still approximate total ice load quite well. The
(Case 2).
SC R
(Case 1) and
IP
horizontal thickness distribution of Model0l and ICEMOD2.1 are almost equal for
7. For a very large spray mass flux Model02 is useless since a large continuous mass flux will prevent icing on a major part of the upwind side of the cylinder. Using this model
NU
ice accretion is restricted to a small area near the cylinder shoulders, from about for Case 2,
MA
from the stagnation line for Case 1 and from about and this is obviously not correct.
D
With the latest extension of ICEMOD2 the main goal of making a relatively complete,
TE
theoretical structured model with a significant complement of supporting empirical parameterizations and empirical sub-models is obtained. The use of ICEMOD2.1 instead of
CE P
ICEMOD2 is however mainly of academic interest (brine film lateral velocity). In many cases a sufficient approximation of ice accretion on cylinders can be obtained by just solving two
AC
simple algebraic equations (Eq. 58). But it is after all satisfying that the complete set of icing equations in fact is solvable and gives realistic results. Improving some of the input parameters (for example based on turbulence theory) and considering other geometries may perhaps show that ICEMOD2.1 nevertheless is favorable compared with other more simple models. Further work on icing model should concentrate on studying in more detail some of the input parameters such as local heat transfer coefficient, local single spray mass flux (or alternatively liquid water content), spray frequency, single spray duration, droplet flight time and the interfacial distribution coefficient (is a constant value always a good approximation?). In most icing models heat conduction in the accreted ice close to the brine film/ice interface is
ACCEPTED MANUSCRIPT 72
neglected. A common argument for doing so has been the high liquid porosity of accreted ice near the icing front. But this is not necessarily true and a study of heat transport in the ice and
T
structure is a challenge for future research on spray icing. The problem is complicated for a
IP
time-dependent model due to the fluctuation of temperature at the icing surface.
SC R
A suggestion for future work is also an investigation of the influence of precipitation on sea spray icing. Snow showers were for example reported for several of the observations on the
NU
two SEDCO platforms during the spray icing periods. A combined spray+snow icing model was in fact introduced in the old one-dimensional version of ICEMOD (Horjen, 1990). Based
MA
on ICEMOD2 this model may be extended to a two-dimensional version. A mathematical model can however under no circumstances reproduce the nature completely.
D
A more complicated icing model will not necessarily give more reliable results for field
TE
conditions. Spray impingement on a rig structure is for example more or less a stochastic
CE P
process. The only way to really test the validity of an icing model is to compare model tests
AC
with controlled laboratory experiments.
References
Achenbach, E.G. (1977): “The effect of surface roughness on the cross flow of air”, J. of Heat and Mass Transfer, Vol. 20, pp 359-369 ANSYS ( 2009): ANSYS FLUENT 12.0 Theory guide, ANSYS Inc., Cannonsburg, 816 pp. Assur,A. (1958): “Composition of sea ice and its tensile strength”, in Arctic Sea Ice, U.S. National Academic Sciences – National Research Council Pub. 598, pp. 106-138
ACCEPTED MANUSCRIPT 73
Best, A.C. (1950): “The size distribution of raindrops”, Quart. J. Roy. Met. Soc.,Vol.77, pp.16-36
T
Best, A.C. (1951): “Drop-size distribution of clouds and fog”, Quart. J. Roy. Met. Soc., Vol.
IP
77, pp. 418-426
SC R
Brown, R. D. and Horjen, I. (1989): “Evaluation of state of the art drilling platform icing models”, Canadian Climate Centre, Report No. 89-10, Atmospheric Environment Service,
NU
Ontario, 80 pp.
MA
Brown, R. D., Horjen, I. , Jørgensen, T. and Roebber, P.(1988): “Evaluation of state- of- theart drilling platform ice accretion models”, In 4th Int. Workshop on Atmospheric icing of
D
Structures, Paris, pp. 208-213
TE
Cheney, W. and Kincaid, D. (1985): “Numerical mathematics and computing”, Brooks/Cole
CE P
Publishing Company, Monterey, California. Cox, G.F.N. and Weeks, W.F. (1983): “Equations for determining the gas and brine volumes
AC
of sea-ice samples”, Journal of Glaciology, 29(102), pp. 306-316 Cox, G.F.N. and Weeks, W.F. (1986): “Changes in the salinity and porosity of sea-ice samples during shipping and storage”, Journal of Glaciology, 32(112), pp. 371-375 Forest,T., Lozowski, E. and Gagnon, R. (2005): “Estimating marine icing on offshore structures using RIGICE04”, International Workshop of Atmospheric Icing on Structures (IWAIS), Montreal. Horjen, I. and Vefsnmo, S. (1985): “A kinematic and thermodynamic analysis of sea spray (in Norwegian), Offshore Icing – Phase II, Norwegian Hydrotechnical Laboratory (NHL). Report STF60 F85014
ACCEPTED MANUSCRIPT 74
Horjen, I. and Vefsnmo, S. (1986): “Calibration of ICEMOD - Extension to a time-dependent model”, Offshore Icing – Phase I/II, NHL- report STF60 F86040, 73 pp.
T
Horjen, I. and Vefsnmo, S. (1987): “Time-dependent sea spray icing on ships and drilling rigs
IP
– a theoretical analysis”, Offshore Icing – Phase IV, NHL- report STF60 F87130, 87 pp.
SC R
Horjen, I. (1990): “Numerical modeling of time-dependent marine icing, anti-icing and deicing”, doctor theses, Trondheim, Norway
NU
Horjen, I. (2013): “Numerical modeling of two-dimensional sea spray icing on vessel-
MA
mounted cylinders”, Cold Regions Science and Technology, 93, p. 20-35 Incropera, F. and De Witt, D.P. (1985): “Fundamentals of heat and mass transfer”, second ed.,
D
Wiley, New York
TE
Jones, K.F. and Andreas, E.L. (2012): “Sea spray concentration and the icing of fixed
144
CE P
offshore structures”, Quarterly Journal of the Royal Meteorological Society, 138(662):131-
AC
Jørgensen, T.S. (1984): “Sea spray measurements on the drilling rig “Treasure Scout” .Results from tests April 1983-February 1984” (in Norwegian), Offshore Icing-Phase II, OTTERreport STF88 F84038
Jørgensen, T.S. (1985): “Sea spray characteristics on a semi-submersible drilling rig”, Offshore Icing-Phase II, NHL-report STF60 F85015 Kays W.M. and Crawford M.E.(1980): “Convective heat and mass transfer”, McGraw-Hill, New York, 420 pp.
ACCEPTED MANUSCRIPT 75
Kulyakhtin, A., Kulyakhtin, S. and Løset, S. (2013):”Measurements of thermodynamic properties of ice created by frozen spray”, Proc. 23. International Offshore and Polar
T
Engineering, Anchorage, Alaska, p.1104-1111
IP
Kulyakhtin, A. and Tsarau, A. (2014):”A time-dependent model of marine icing with
SC R
application of computational fluid dynamics Measurements”, Cold Regions Science and Technology, 104-105, p. 33-44
NU
Kulyakhtin, A. (2014):”Numerical Modelling and Experiments on Sea Spray Icing”, doctor
MA
theses, Trondheim, Norway
Lai, R. J. and Shemdin, O.H. (1974): “Laboratory study of the generation of spray over
D
water”, J. Geophysical Res., 79(21), pp. 3055-3063
TE
Lax P.D. and Wendroff, B.(1964): “Difference schemes for hyperbolic equations with high
CE P
order of accuracy”, Comm. Pure Appl. Math., Vol. 17, pp. 381-398 Lozowski, E.P., Stallabrass, J.P. and Hearty, P.F. (1979): “The icing of an unheated nonrotating cylinder in liquid water droplet-ice crystal clouds”, Mech. Eng. Report LTR-LT-96,
AC
National Research Council, Canada, 109 pp. Makkonen, L. (1985): “Heat transfer and icing of a rough cylinder”, Cold Regions Science and Technology, vol. 10, pp.105-116 MacLaren Plansearch (1986): “Final report of the offshore icing measurement project, Scotian Shelf. Contractors final report prepared for Atmospheric Environment Service”, Bedford, Nova Scotia.
ACCEPTED MANUSCRIPT 76
Myers, T.G. and Charpin, J.P. (2004): A mathematical model for atmospheric ice accretion and water flow on a cold surface”, International Journal of Heat and Mass Transfer, 47(25):
T
5483-5500
IP
Minsk, L.D. (1984): “Ice observation program on the semisubmersible drilling vessel SEDCO
SC R
708”, CRREL Special Report 84-2, 14pp
Mitten, P. (1994): “Measurement and modeling of spray icing on offshore structures”, Final
NU
report, Atmospheric Environmental Service of Canada, Contract no. 07SE.KM169-8-7439
MA
Monaham, E.C. (1968): “Sea spray as a function of low elevation wind speed”, J. Geophysical Res., 73(4), pp. 1127-1137
D
Preobrazhenskii, L.Y. (1973): “Estimate of the content of spray drops in the near-water layer
TE
of the atmosphere”, Fluid Mechanics – Soviet Research, 2 (2): pp. 95-100
CE P
Richtmyer, R.D. (1963): “A survey of difference methods for non-steady fluid dynamics, N.C.A.R. Tech. notes 63-2
AC
Ryerson, C.R. and Gow, A.J. (2000):”Ship superstructure icing: Crystalline and physical properties”, ERDC/CRREL (Cold Regions Research and Engineering Laboratory) report TR00-11
Schlichting, H. (1955): ”Boundary layer theory”, McGraw-Hill, First English edition, 535 pp. Sobey, I. J. (2000): “Introduction to interactive boundary layer theory”, Oxford University Press, ISBN 0-19-850675-9. Stallabrass, R. (1980): “Trawler icing. A compilation of work done at N.R.C.”, Mech. Eng. Report MD-56, No 19327, Ottawa, Canada
ACCEPTED MANUSCRIPT 77
Weeks, W. and Lofgren, G. (1967): “The effective solute distribution coefficient during the freezing of NACL solutions”, In: Physics of snow and Ice, Institute of Low Temperature
T
Science, Hokkaido University, Sapporo, Japan, pp.579-597.
IP
Weeks, W. (2010): “On sea ice”, University of Alaska Press, Fairbanks, AK, 664 pp.
SC R
Wu, J., Murray, J.J. and Lai, R.J. (1984): “Production and distribution of sea spray”, Journal of Geophysical Res., 89 (C5), p. 8163-8169
NU
Zakrzewski, P. (1987): “Splashing a ship with collision-generated spray”, Cold Regions
MA
Science and Technology, 14, p. 65-83
TE
D
APPENDIX:
icing events
CE P
Weather parameters used as input to the SEDCO 708 and SEDCO 709
AC
The data in each line below is from left to right: Date (year, month, day, time), air temperature ( ), dew-point temperature ( ), air pressure (mb), sea surface temperature ( ) , wind speed (knot), significant wave height (m), significant wave period (s) and sea surface salinity (‰). For SEDCO 708 the time is LST while for SEDCO 709 the time is GMT. SEDCO 708 ICING EVENT (60 observations) 83010404, -1.7,
-3.3,
999.9,
3.8,
18.8,
1.5,
5.7, 32.5
83010405, -2.2,
-3.3,
999.9,
3.8,
19.1,
1.5,
5.2, 32.5
83010406, -2.2,
-3.3, 1000.5,
3.8,
23.4,
1.7,
5.2, 32.5
ACCEPTED MANUSCRIPT 78
-3.3, 1000.7,
3.8,
20.6,
1.6,
5.6, 32.5
83010408, -2.8,
-3.9, 1000.9,
3.8,
22.4,
1.5,
5.6, 32.5
83010409, -2.8,
-3.9, 1001.1,
3.8,
21.7,
1.6,
5.4, 32.5
83010410, -2.8,
-3.9, 1001.3,
3.8,
20.5,
83010411, -2.8,
-3.9, 1001.6,
3.8,
20.5,
83010412, -2.8,
-3.9, 1001.5,
3.8,
83010413, -3.9,
-5.0, 1001.1,
83010414, -3.3,
-4.4, 1000.2,
83010415, -2.8,
-4.4, 1000.4,
83010416, -1.7,
-3.3, 1000.4,
83010417, -2.2,
T
83010407, -2.2,
5.8, 32.5
1.6,
5.9, 32.5
21.4,
1.7,
6.1, 32.5
3.8,
17.3,
1.7,
6.1, 32.5
3.8,
17.7,
1.6,
5.8, 32.5
3.8,
14.9,
1.5,
5.9, 32.5
3.8,
19.3,
1.5,
6.1, 32.5
-3.3, 1000.1,
3.8,
20.6,
1.5,
5.7, 32.5
83010418, -2.3,
-3.4, 1000.4,
3.8,
18.5,
1.5,
5.5, 32.5
83010419, -2.3,
-3.4, 1000.4,
3.8,
16.2,
1.3,
5.4, 32.5
83010420, -2.4,
-3.5, 1000.4,
3.8,
15.9,
1.3,
5.3, 32.5
IP
SC R
NU
MA
D
TE
83010421, -2.4,
-3.5, 1000.1,
3.7,
14.7,
1.4,
5.5, 32.5
83010422, -2.5,
-3.6, 1000.2,
3.7,
16.9,
1.3,
5.3, 32.5
AC
CE P
1.6,
83010423, -2.5,
-3.6,
999.9,
3.7,
17.5,
1.3,
5.4, 32.5
83010500, -2.6,
-3.7,
999.6,
3.7,
20.3,
1.3,
5.2, 32.5
83010501, -2.6,
-3.7,
998.5,
3.8,
25.6,
1.5,
5.1, 32.5
83010502, -2.7,
-3.8,
998.0,
3.8,
30.0,
1.8,
5.5, 32.5
83010503, -2.7,
-3.8,
997.4,
3.8,
32.2,
2.2,
5.9, 32.5
83010504, -2.8,
-3.9,
997.6,
3.8,
31.6,
2.5,
6.1, 32.5
83010505, -2.8,
-3.9,
997.3,
3.8,
32.0,
2.3,
6.0, 32.5
83010506, -2.8,
-3.9,
996.7,
3.8,
35.2,
2.6,
6.4, 32.5
83010507, -2.2,
-3.3,
996.3,
3.7,
37.1,
2.8,
6.4, 32.5
83010508, -2.2,
-3.3,
996.3,
3.7,
38.9,
3.0,
6.6, 32.5
ACCEPTED MANUSCRIPT 79
-3.3,
996.4,
3.7,
39.6,
3.1,
7.2, 32.5
83010510, -2.2,
-3.3,
996.7,
3.7,
38.2,
3.5,
7.5, 32.5
83010511, -2.7,
-3.3,
996.9,
3.7,
38.8,
3.9,
8.0, 32.5
83010512, -2.8,
-3.4,
996.6,
3.7,
38.3,
83010513, -2.8,
-3.5,
996.0,
3.7,
37.9,
83010514, -2.8,
-3.6,
995.6,
3.7,
83010515, -2.8,
-3.7,
995.6,
83010516, -2.8,
-3.8,
995.8,
83010517, -2.8,
-3.9,
996.1,
83010518, -2.9,
-4.0,
996.5,
83010519, -2.9,
-4.1,
83010520, -3.0,
-4.2,
83010521, -3.1,
-4.3,
83010522, -3.2,
-4.4,
T
83010509, -2.2,
7.9, 32.5
4.2,
8.4, 32.5
39.3,
3.8,
7.9, 32.5
3.7,
35.2,
3.8,
8.1, 32.5
3.7,
33.4,
4.4,
8.4, 32.5
3.7,
34.8,
4.0,
8.0, 32.5
3.7,
33.7,
4.2,
7.9, 32.5
996.8,
3.7,
32.4,
4.4,
8.2, 32.5
997.3,
3.7,
31.7,
3.8,
8.0, 32.5
997.7,
3.7,
31.9,
3.8,
7.7, 32.5
998.0,
3.6,
30.3,
4.3,
8.0, 32.5
IP
SC R
NU
MA
D
TE
83010523, -3.3,
-4.4,
998.2,
3.6,
27.5,
4.1,
8.3, 32.5
83010600, -3.3,
-4.4,
998.6,
3.6,
29.1,
4.4,
8.6, 32.5
AC
CE P
3.5,
83010601, -3.3,
-4.4,
998.8,
3.6,
29.8,
4.1,
8.5, 32.5
83010602, -4.4,
-5.5,
998.9,
3.6,
30.4,
3.6,
8.2, 32.3
83010603, -5.6,
-6.7,
999.2,
3.7,
36.3,
3.8,
8.3, 32.3
83010604, -6.7,
-7.8,
999.2,
3.7,
38.1,
4.8,
9.0, 32.5
83010605, -6.7,
-7.8,
999.1,
3.7,
39.3,
5.2,
9.0, 32.5
83010606, -7.2,
-8.3,
999.1,
3.7,
38.8,
5.4,
9.3, 32.5
83010607, -8.3,
-9.4,
999.2,
3.7,
40.0,
5.3,
9.0, 32.5
83010608, -8.3,
-9.4,
999.2,
3.7,
39.4,
5.0,
8.8, 32.5
83010609, -8.3,
-9.4,
999.2,
3.7,
38.7,
4.4,
8.5, 32.5
83010610, -7.8,
-8.9,
999.0,
3.7,
38.6,
4.8,
8.9, 32.5
ACCEPTED MANUSCRIPT 80
-9.4,
998.8,
3.7,
38.7,
4.9,
8.7, 32.5
83010612, -8.3,
-9.4,
998.4,
3.7,
40.3,
5.2,
8.8, 32.5
83010613, -8.3,
-9.4,
997.9,
3.7,
39.1,
4.9,
9.2, 32.5
83010614, -8.3,
-9.4,
997.0,
3.7,
39.2,
83010615, -8.3,
-9.4,
996.3,
3.7,
39.1,
86012415, -5.0,
-8.0, 1028.0,
86012416, -5.0,
-8.0, 1027.5,
4.9,
9.0, 32.5
IP
SC R
-8.0, 1027.4,
1.9,
21.0,
1.4,
4.9, 34.0
1.8,
24.0,
1.4,
4.5, 34.0
1.8,
26.0,
1.5,
4.4, 34.0
86012417, -7.0, -10.0, 1027.4,
1.7,
27.0,
1.7,
4.4, 34.0
86012418, -7.0, -11.0, 1027.8,
1.7,
29.0,
1.7,
5.9, 34.0
86012419, -7.0, -11.0, 1028.2,
2.1,
28.0,
1.6,
4.4, 34.0
86012420, -7.0, -11.0, 1029.5,
1.9,
32.0,
1.8,
4.5, 34.0
86012421, -8.0, -11.0, 1030.9,
2.0,
35.0,
1.8,
4.4, 34.0
86012422, -8.0, -11.0, 1032.0,
2.0,
37.0,
2.1,
4.7, 34.0
86012423, -9.0, -11.0, 1033.6,
1.8,
34.0,
2.0,
4.5, 34.0
86012500, -9.0, -12.0, 1034.1,
1.6,
36.0,
2.0,
4.5, 34.0
86012501, -9.0, -12.0, 1034.9,
1.7,
32.0,
2.0,
4.5, 34.0
86012502,-10.0, -15.0, 1035.8,
1.8,
30.0,
2.1,
4.7, 34.0
86012503,-10.0, -15.0, 1035.8,
1.8,
30.0,
1.8,
4.5, 34.0
86012504,-10.0, -15.0, 1036.9,
2.1,
29.0,
1.8,
4.6, 34.0
86012505,-10.0, -14.0, 1037.3,
1.4,
27.0,
1.8,
4.8, 34.0
86012506,-10.0, -14.0, 1039.3,
0.8,
25.0,
1.7,
4.7, 34.0
86012507,-10.0, -14.0, 1040.4,
1.0,
21.0,
1.7,
4.7, 34.0
AC
CE P
TE
D
MA
86012414, -5.0,
8.9, 32.5
4.6,
NU
SEDCO 709 ICING EVENT (25 observations)
T
83010611, -8.3,
ACCEPTED MANUSCRIPT 81
1.4,
23.0,
1.8,
4.8, 34.0
86012509,-10.0, -16.0, 1041.8,
0.9,
21.0,
1.8,
4.9, 34.0
86012510,-10.0, -16.0, 1042.2,
0.5,
20.0,
1.6,
4.7, 34.0
86012511, -9.0, -14.0, 1043.1,
0.2,
24.0,
86012512, -9.0, -14.0, 1043.8,
1.1,
24.0,
86012513, -9.0, -14.0, 1044.9,
1.0,
86012516, -9.0, -14.0, 1044.8,
1.2,
T
86012508,-10.0, -14.0, 1041.2,
4.5, 34.0
1.4,
4.5, 34.0
22.0,
1.4,
4.5, 34.0
19.0,
1.1,
4.1, 34.0
IP
SC R
NU
1.6,
AC
CE P
TE
D
MA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ACCEPTED MANUSCRIPT 82
CE P
TE
D
MA
NU
SC R
IP
T
An extension of the two-dim. sea spray icing model ICEMOD2 is presented. Model results are compared with field results from two semi-submersible oil rigs. Sensitivity analysis of the effect of varying three input parameters is presented. The effect of including brine film motion in icing modeling has been discussed. The main difference between the new and old ICEMOD2 is the brine film velocity.
AC
S0165-232X(15)00165-2 doi: 10.1016/j.coldregions.2015.07.006 COLTEC 2138
To appear in:
Cold Regions Science and Technology
Received date: Revised date: Accepted date:
12 May 2014 13 July 2015 20 July 2015
Please cite this article as: Horjen, Ivar, Offshore Drilling Rig Ice Accretion Modeling Including A Surficial Brine Film, Cold Regions Science and Technology (2015), doi: 10.1016/j.coldregions.2015.07.006
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT 1
OFFSHORE DRILLING RIG ICE ACCRETION MODELING INCLUDING A SURFICIAL BRINE FILM
IP
T
Ivar Horjen1, dr. techn., Trondheim, July, 2015
SC R
ABSTRACT
This article presents an extension of the two-dimensional icing model ICEMOD2. In contrast
NU
to this model the new model (called ICEMOD2.1) considers the inertia terms in the momentum equation. The derivation of the model is presented and the model is applied to two
MA
typical structures of an offshore drilling rig which may be exposed to spray icing: a corner column and a truss. A mean impact-generated spray mass flux formula (called the “normal”
D
value) based on droplet kinematics and spray measurements on the Norwegian semi-
TE
submersible rig “Treasure Scout” has been implemented in the model. The model results
CE P
showed satisfactory agreement with some actual icing registration data from the drilling rigs “SEDCO 708” and “SEDCO 709”. Additionally sensitivity tests are presented for a selected section of a rig column of diameter 5.0 m varying air temperature, wind speed and spray
AC
frequency. As expected ice accretion increases with decreasing air temperature. With the chosen formula for the spray mass flux and wind speed (17-25 m/s) ice accretion also increases with increasing wind speed. The available field results of Treasure Scout enables us only to construct a time-average spray mass flux formula independent of spray frequency. The single spray mass flux is then in the model a derived quantity obtained by dividing the timeaverage value by the product of spray frequency and spray duration. Based on these presumptions the sensitivity tests show that spray frequency variations has only a minor effect upon accreted ice load when normal spray mass flux is assumed (corresponding here to 5-10 g m-2 s-1). If however the spray mass flux is put 25 times as large as the normal value a 1
Mail address: [email protected]
ACCEPTED MANUSCRIPT 2
substantially increase in model ice load with increasing spray frequency is obtained. Probably this is due to a large increase in brine film thickness (a factor of about 10). A
T
possible theoretical explanation is given in the text Finally the significance of brine film
IP
motion in icing modeling has been discussed. It turns out that for cylinder icing ICEMOD2.1
SC R
and ICEMOD2 give the same ice load and thickness distribution. The only noticeable difference in the two model results is the brine film horizontal velocity, but only for an increased spray mass flux (of order
). For normal spray mass flux a
NU
simplified time-dependent model version with brine film velocity put equal to zero (Model01)
MA
gives similar ice accretion results as ICEMOD2.1. Even a stationary version of the latter model (Model02) gives adequate results for normal spray conditions. For the increased spray mass flux some brine may however be pushed over to the lee side by wind stress, especially if
TE
D
the cylinder is slim such as a truss, and only the complete versions of ICEMOD2.1 or ICEMOD2 can handle this situation. Model01 will still give correct ice thickness results on
CE P
the major part of the upwind side of the cylinder while Model02 only give some icing close to the cylinder shoulders (which means that this model cannot approximate periodic icing in this
AC
case).
KEY WORDS: icing modeling, rig icing, momentum equation, brine film motion
NOMENCLATURE R = instantaneous local icing intensity (kg m-2 s-1) Rw = instantaneous local water catch rate (kg m-2 s-1) = time-average spray mass flux “far” from the icing object (kg m-2 s-1) = spray mass flux of a single spray “far” from the icing object (kg m-2 s-1) Q = unit heat flux (W m-2) U = wind or spray velocity (m/s)
ACCEPTED MANUSCRIPT 3
T = wave period (s) S = salinity (‰)
T
D = initial cylinder diameter (m)
IP
Vh = free horizontal spray velocity (m/s)
P = pressure (
SC R
N = spray frequency (1/s) )
C = specific heat capacity at constant pressure (
)
NU
=mean brine film velocity with components in s and ς direction
MA
t = time (s) f = humidity
)
D
g = abs. value of acceleration of free fall (
TE
l0 = specific latent heat of fusion of pure ice (J/kg) le = specific latent heat of evaporation (J/kg)
CE P
= sand grain roughness height (m)
Λ = liquid water content (LWC) (kg m-3 ) )
AC
θ = temperature (
ρ = density (kg m-3)
δ=brine film thickness (m) μ = dynamic viscosity (kg m-1 s-1) τ = stress (
)
λ = wave length (m) σ = interfacial distribution coefficient χ=mass diffusivity of diluted salts (
)
ὐ=volume fraction = mass flux of vapor to/from the brine film surface (kg m-2 s-1)
ACCEPTED MANUSCRIPT 4
K, N = non-dimensional functions in the brine film momentum equations. E = non-dimensional extrapolation factor (defined in Eq. 50)
IP
T
Subscripts
SC R
C = collision-generated i = ice b = brine film mean value in normal direction
NU
a = air
MA
e = evaporation or condensation w = sea surface
D
s = significant or s-direction
TE
sp = spray
10 = 10 m above msl (mean sea level)
CE P
0 = wave crest (unless otherwise stated)
AC
50 = median volume value
ACCEPTED MANUSCRIPT
CE P
TE
SEDCO DRILLING RIG (from Google)
D
MA
NU
SC R
IP
T
5
1. INTRODUCTION
AC
In the 1980’s a numerical marine icing package called” ICEMOD” (Horjen and Vefsnmo, 1987, Horjen,1990) was developed at the Norwegian Hydrotechnical Laboratory (NHL). The work was going on for about ten years. This model was one-dimensional in the sense that the brine film movement was restricted to only one space dimension along the icing surface. For example, for a vertical rig column icing intensity could only be correctly calculated at the stagnation line where the brine film movement is solely downwards affected only by gravity. In order to calculate ice mass on the whole column a specified shape of the ice profile must then be postulated, for example that the ice thickness in the wind direction is constant along the whole upwind side of the column. Maximum radial ice thickness will then be at the stagnation line. This model was for example used for analyzing several drilling rig icing cases
ACCEPTED MANUSCRIPT 6
in North American waters (Brown and Horjen, 1989, Brown et al, 1988). ICEMOD was in this report compared with the Canadian model RIGICE. Ice thickness was estimated for three
T
of the platforms (SEDCO 708, SEDCO 709 and Rowan Gorilla) based on photographs.
IP
Neither ICEMOD nor RIGICE gave very convincing results for all these cases.
SC R
Recently an extension of ICEMOD to two-dimensional brine film motion was presented and applied for icing on cylinders of diameters 5-50 cm placed on a vessel (Horjen, 2013). The
NU
brine film velocity components were here found directly from the brine film momentum equations with inertia terms neglected. An interesting feature of the new model (called
MA
ICEMOD2) is the possibility of depression of the ice profile at the stagnation line. If the freezing fraction is high, i.e. close to 1, maximum ice thickness will be at the stagnation line
D
where the spray mass flux is also a maximum. For sufficient wet icing, however, the ice
TE
thickness profile on the upwind side of the cylinder is to a greater extent controlled by the heat supply to the brine film from the spray cloud and convective and evaporative heat loss.
CE P
These heat fluxes vary along the cylinder circumference. In the present article the model is extended a step further by including the equation of motion
AC
of the brine film. This version of the model will be called ICEMOD2.1. In order to simplify the model somewhat the brine film flow is assumed to be laminar. The importance, if any, of including brine film motion in a spray icing model was not discussed in Horjen (2013); in this article we will discuss this question in some detail in section 4.3. As in ICEMOD2 we will adopt the cylindrical approximation when calculating the local wind field (from potential flow theory), collection efficiency, local heat transfer coefficient and brine film movement. The model cylinder diameter increases with time assuming the accreted ice distributes uniformly on the front side of the cylinder (Horjen, 1990). Model results will be compared with observed ice thickness profiles on selected columns and trusses
ACCEPTED MANUSCRIPT 7
on an offshore drilling platform. Ice thickness growth rate depends on the icing intensity and the ice density close to the icing surface (
). Both these parameters are functions of
T
the brine film salinity. Additionally the ice density is also a function of the volume fraction , simply
IP
of air in the ice. In earlier presentations of ICEMOD/ICEMOD2 we have put
SC R
due to the lack of information of this quantity. Recently from a series of field measurements in the harbor area of Longyearbyen, Svalbard, air entrapment was for the major part of the tests (24) found to vary between 3% and 7 % on a vertical structure (Kulyakhtin et al, 2013).
NU
Icing conditions during these experiments were more or less different from field conditions on
MA
a drilling rig; the median droplet diameter was for example much smaller, in the range from 90 μm to 300 μm. Ryerson and Gow (2000) have reported a few calculated air volumes in ice
D
accreted on various part of the USCGC Midgett. The four values from ice samples collected
TE
on the bulkhead varied between 11.9 and 25.7 %. This is much higher than the results of Kulyakhtin et al, 2013. Possibly this may partly be due to some brine drainage before the
CE P
density measurements took place. For all the model results to be presented here we have simply assumed a constant air porosity of 5%. Calculated ice thickness then becomes 5.3%
AC
higher than for air-free ice.
There are no restrictions on the sign of R. If the heat input from the impinging spray cloud is very large melting of formerly accreted ice may occur, giving a negative icing intensity (a test is included in the program for this case to check if there is actually any ice to melt; if not the icing intensity is put equal to zero). This theoretically structured model is supported by several empirical parameterizations, for example the convective heat transfer coefficient, the duration of a single spray, the spray frequency, the size of the median volume droplets and the mean spray mass flux height distribution of impact-generated spray. These sub-models are important if we want to obtain icing results comparable with field results and must hence be supported by considerable
ACCEPTED MANUSCRIPT 8
scientific field experience. The main intension of this article is however to present the theoretical background of the extended model and show that calculated ice thicknesses on
T
cylinder structures are not too unrealistic compared with field results, even if the sub-models
SC R
IP
may be profitably scrutinized in future research.
NU
2. MATHEMATICS OF THE EXTENDED ICING MODEL The two-dimensional model may be used to calculate icing on various cylindrical and planar
MA
elements on both vessels and offshore rigs. The problem is two-dimensional when the icing structure is vertical or heeling, for example the bulkhead and mast of a vessel and the column,
TE
icing on rig columns and trusses.
D
diagonal truss and side wall of a rig. In this article we will however limit ourselves to study
CE P
Generally the geometry of this problem is shown in Fig.1 where the heeling angle σ is 90 for columns and typically between 50 and 70 for trusses. Icing takes place between lines parallel to the cylinder axis through the points
AC
the line through the lower point
. Origin of the coordinate system is at
at the bottom of the cylinder. If the overall collection
efficiency is equal to 1 (which we assume in this article for impact-generated droplets) will be on the y-axis.
ACCEPTED MANUSCRIPT
MA
NU
SC R
IP
T
9
TE
D
Fig. 1. Sketch of the geometry of the icing configuration showing the model’s coordinate system. On the left there is a cross-sectional view of the column and truss. On the right there
CE P
is a side-view of the cylindrical substrate.
AC
The integrated conservation equations of brine film mass, salt and momentum in s and ς direction have in fact been derived many years ago (Horjen and Vefsnmo, 1986). Since then no effort has however been made to solve the complete set of differential equations, which are: Brine mass conservation equation:
(1)
Salt mass conservation equations:
(2)
ACCEPTED MANUSCRIPT 10
Linear momentum conservation equation in s-direction:
IP
Linear momentum conservation equation in ς-direction:
T
(3)
SC R
(4)
(brine mass per unit area) and
(salt mass per unit area). Index b means
, for example salinity, across the brine film layer:
MA
mean value of a function
NU
is the mass exchange due to evaporation (negative) or condensation (positive),
D
(5)
TE
Both brine film salinity and velocity components are defined by this integral. The brine film
CE P
density and viscosity are assumed to be constant across the brine film layer. Eq. 3 and Eq. 4 are slightly modified from the original set of equations of the 1986-report due to a change in the assumed w-velocity profile across the brine film layer and some small corrections.
AC
Formerly the” icing equation” were defined to be Eq. 1 and Eq. 2; we now define the icing equations to be the whole set of the four equations 1-4 with the four dependent variables: X, Y, and
. Using the brine film temperature/salinity equilibrium conditions from Assur
(1958) and Cox and Weeks (1983,1986) Eq. 1 and Eq. 2 may be combined to give an explicit expression of the icing intensity as a function of the brine film salinity (Horjen,2013):
(6)
where
ACCEPTED MANUSCRIPT 11
is the net external heat flux (most important being heat convection, evaporation and heat is the conducive heat flux to the
T
exchange between impinging spray and brine film) and
IP
accretion( to be disregarded in this work). The heat fluxes are defined positive when heat is
SC R
transported to the brine film. 2.1 The brine film sub-model
in Eq. 3 and Eq. 4 are the brine film velocity componenets in the s and ς
NU
and
MA
direction assuming steady brine film motion ,i.e. setting the inertia terms in the momentum equations equal to zero (Horjen and Vefsnmo, 1987):
TE
D
(7)
where
AC
(9)
CE P
(8)
In Eq. 3 and Eq. 4 K and N are functions of X:
(10)
.
(11) The somewhat peculiar form of the expressions of K and N are due to the chosen second order polynomials in the n-coordinate for the velocity components (to be discussed below). Note that the expression of
(Eq. 7) is slightly different from the corresponding expression
given in Horjen (2013). In the latter article an exponential w-velocity profile is chosen and we
ACCEPTED MANUSCRIPT 12
also assumed that
, i.e. relatively low spray mass flux and/or thin brine film
layer.
direction while P is the sum of the spray and wind pressure (
T
) acting on the brine film in the s-
IP
is the sum of spray and wind stress (
) acting normal to
SC R
the brine film. g is the acceleration due to gravity (defined positive). From Horjen (1990) we have:
NU
(12)
is the shear stress exerted by the air at the brine film surface:
CE P
(13) where
is the spray bulk velocity component in the s and n
TE
direction respectively and
and
D
assume potential flow.
MA
is the wind speed above the surface boundary layer in s-direction. We will in this work
is the air/brine film friction coefficient.
AC
Using the Blasius method the theoretical local friction coefficient assuming potential flow becomes (see Schlichting, 1955):
(14)
where
(15) (16) (υ in radians)
(Reynolds number of the cylinder)
ACCEPTED MANUSCRIPT 13
The function
is a theoretical result based on the assumption of a laminar boundary layer
flow and potential, irrotational flow outside the boundary layer. The potential flow is
, giving
. The flow is called subcritical if the boundary layer is
IP
setting
T
expanded in a power series of five terms. The (theoretical) point of instability is found by
SC R
laminar up to the point of separation (as assumed above) and supercritical if the boundary layer undergoes transition and becomes turbulent before separation takes place. The transition between the two boundary layer flow regimes around a smooth cylinder takes place for a (various values are found in the literature, depending probably
NU
Reynolds number of order
MA
on test conditions such as free stream turbulent intensity). For a smooth cylinder in subcritical flow the use of Pohlhausen’s method, based on measured pressure data, results in boundary (Sobey, 2000), i.e. almost 30 degrees lower than the
D
layer separation just above
TE
theoretical result of Schlichting. In the supercritical flow regime separation is delayed up to (Incropera and deWitt, 1985) with transition for some angle
CE P
Reynolds number is very large (
. If the
) the flow becomes transcritical which means
that the boundary layer is turbulent on the major part of the cylinder circumference. For the
AC
cylinder icing cases of the SEDCO platforms the Reynolds number is larger than about , i.e. the flow is supercritical. If the cylinder is rough some part of the upwind boundary layer will be turbulent provided the cylinder Reynolds number is larger than a threshold value. Makkonen (1985) has for example found that for a cylinder of diameter D=0.15 m, Reynolds number sand grain roughness
and equivalent
the boundary layer will be turbulent for
, i.e.
almost the whole upwind side of the cylinder. For the same cylinder it was found that the whole upwind side boundary layer is laminar for a Reynolds number of
. In the
quoted work the local friction coefficient for a turbulent boundary as a function of the equivalent sand grain roughness was used (Kays and Crawford, 1980). But since
is a priori
ACCEPTED MANUSCRIPT 14
not known for rig icing cases only the laminar expression Eq.14 has for the time being been implemented in the icing model. Moreover,
for an icing surface will certainly depend on
T
both the geometry of the structure and environmental parameters. We also assume for
IP
convenience that boundary layer separation occurs at the point of theoretical instability for all
SC R
values of the Reynolds number. The local heat transfer coefficient used in the model is however based on turbulent flow theory (see section 2.2).
NU
A short outline of the derivation of Eq. 1 and 2 is given in Horjen (2013). The full derivation of Eq. 3 and 4 is rather elaborating and only a short outline will be given here as well. The
MA
starting point is the three boundary layer scalar equations of linear momentum:
D
(17)
TE
(18)
where
CE P
(19)
is the local brine film velocity with components in a local curvilinear
AC
coordinate system
, the n-axis lying along the outward normal to the icing front
(Fig.2). Index f is used for local brine film properties. Assuming circular cross section the components of acceleration of free fall are given by: (20)
ACCEPTED MANUSCRIPT
NU
SC R
IP
T
15
(21)
D
MA
Fig.2 Detail of the ice and brine film layer with control volumes V1 and V2.
TE
At the inner film boundary the brine film velocity tangent to the phase interface must be zero:
CE P
(22)
(23)
AC
Additionally the following inner boundary condition is obtained from Eq. 17 and Eq. 21:
The outer film boundary conditions may be derived using the so-called control volume method: we imagine a geometric volume V1, i.e. fixed in space, covering the outer brine film boundary. Conservation of linear momentum for the control volume letting
such that
the control volume contains the interface air/brine film then gives us three boundary conditions at (24) (25)
(see Horjen, 1990, for details):
ACCEPTED MANUSCRIPT 16
(26) Any brine film motion in normal direction has been neglected.
is the s-component of the
IP
T
shear stress in the brine film:
SC R
(27)
The control volume method was also used to define the boundary conditions for the enthalpy
NU
and salt diffusion equations from which Eq. 1 and Eq. 2 were derived. Conservation of for example mass of diluted salts in the control volumes of Fig. 2 for
and
in such
MA
a way that the interfaces air/brine film and ice/brine film is inside the control volumes results in the following two boundary conditions for the brine film salinity:
TE
D
(28)
CE P
(29)
The next step is to integrate Eq. 17-19 across the brine film layer in n-direction using the boundary conditions Eq. 21 and Eq. 23-25. After some manipulations we arrive at the
(30)
AC
following results:
(31)
where
is the gradient operator in tangential direction:
(32)
To get any further we must now specify the velocity component profiles across the brine film layer and we assume, somewhat arbitrary, a second order profile for both the u-variation and
ACCEPTED MANUSCRIPT 17
the w-variation (for a thin water film moving vertically downwards with no spray input the approximate solution of Navier-Stokes equation is also a second order velocity profile). For
T
atmospheric icing conditions Myers and Charpin (2014) has made a more thorough analysis
IP
of water film movement resulting in a second order polynomial velocity profile. In their case
SC R
water film velocity is a function of spray mass flux, heat fluxes from the icing surface and icing surface geometry. For sea spray icing the brine film velocity will also depend on brine salinity (since the four icing equations are linked together). Using now the boundary
NU
conditions Eq. 21, 22, 24 and 25 the velocity profiles to a sufficient degree of accuracy
MA
becomes: (33)
and
are the brine film velocity components at the brine film surface and
TE
where
D
(34)
AC
(35)
CE P
defined by:
is
For a vertical cylinder
is identical zero. Provided the brine film is thin, for example 0.1
mm and less, this term will be rather small for a heeling cylinder (truss). For simplicity we will in this work generally put Since now
and
. Eq. 32 and Eq. 33 may alternatively be written:
(36) (37) Putting Eq. 35 and Eq.36 in Eq. 29 and Eq. 30 we finally obtain, after some algebra, Eq. 3 and Eq. 4. Thus, ICEMOD2.1 features a rather complete hydrodynamic model of momentum
ACCEPTED MANUSCRIPT 18
conservation that includes molecular heat and salt transport in a laminar brine film layer adjacent to the icing surface (Fig. 2). Brine film salinity plays a key role in this model since
T
the icing intensity is an algebraic equation of this parameter (Eq. 6). Close to the inner and
IP
outer boundary of the brine film layer there may be large salinity gradients (Eq. 27 and Eq.
SC R
28), but since the liquid layer for most of the icing conditions studied in this article is rather thin (typically less than 0.1 mm) it is justified to use the mean salinity value normal to the icing front to calculate the icing intensity. If the mean spray mass flux is very large (examples
NU
are given in section 4.3) the brine film thickness may be as large as about 1.0 mm, and a
MA
detailed study of heat and salt transport within the brine film layer will probably be an improvement of the icing model. On the other hand an increased complexity of the model will make the computation time unacceptable long. A thick brine film means that it becomes
TE
D
more important to include the inertia terms in the momentum equations, which in itself increases the computation time significantly. Moreover, it is not certain that a more complex
CE P
model will bring about more reliable icing results, at least for field conditions. The interfacial distribution coefficient is set to a constant value
( Horjen, 2013). But
AC
it must be kept in mind that there is no physical reason that this parameter should actually be a constant. The value chosen is just a mean value of 29 Japanese, Finish and Russian observations ranging from 0.18 to 0.5 with an average deviation from the mean of 0.068. If a layer of saline water is freezing it may be shown that the interfacial distribution coefficient depends strongly on the ice growth rate (Weeks and Lofgren, 1967). If for example the growth rate increases from
to
the interfacial distribution coefficient
increases from 0.29 to 0.73 (from Weeks, 2010). In a spray icing model of Kulyakhtin, 2014, including heat conduction an expression of the ratio of ice salinity and brine film salinity at icing surface is presented which depends on a combination of the constant value used in ICEMOD2 and the value
from the theory of Weeks and Lofgren (Kulyakhtin, 2014):
ACCEPTED MANUSCRIPT 19
(38)
where
is the icing intensity caused by heat conduction to the ice only and
is the icing
IP
T
intensity caused by external heat fluxes. It is interesting that this equation indicates an
SC R
“effective interfacial coefficient” depending on icing conditions. In this work heat conduction is neglected and we have put total icing intensity equal to
. If however heat conduction is
to be included in ICEMOD2.1 (a challenge for a further improvement of the model) this term
from the total icing intensity (Eq. 6).
MA
separate the “conductive” icing intensity
NU
will be a function of brine film salinity and time. For periodic spray it is hence not possible to
Other important input parameters to the icing model include the heat transfer coefficient and
TE
2.2 The Nusselt number
D
various sea spray parameters. This will be discussed in the next two sections.
CE P
In the article of Horjen, 2013, a formula originally proposed by Lozowski et al., 1979, was used for the upwind part of cylinders:
AC
(39)
is the stagnation line Nusselt number and
is the the angle from the stagnation line in
degrees (positive in counterclockwise direction). This expression is based on experimental results of Achenbach (1977) for
. The cylinders examined in Horjen, 2013, had
diameters from 0.05 m to 0.5 m and the Reynolds number was in the range from
to
. Hence Eq. 38 is appropriate in this case. In this study the dimensions are much larger and the Reynolds number for the vertical columns of D=9.15 m was in the range from to
for the SEDCO 708 case and in the range
to
for the
SEDCO 709 case. This is outside the definition range of Eq. 38 and another Nusselt number formula should hence be applied. Based on computational fluid dynamics using the program
ACCEPTED MANUSCRIPT 20
FLUENT (ANSYS, 2009) the following formula has been proposed by Kulyakhtin and Tsarau (2014) for offshore rig columns:
SC R
IP
T
(40)
Excellent agreement (
) was found comparing this formula with the computed
results for a cylinder of diameter 90 m and wind speeds of 3, 9, 12, 15 and 22 m/s at Z=10 m.
NU
The Reynolds number was hence in the range 1.8 107-1.3 108, i.e. transcritical flow. The formula above has been adopted in ICEMOD2.1 for rig icing cases and completes the set of
from Eq. 39.
D
Eq. 38 and at
from
MA
equations contained in the model. Note that maximum Nusselt number is at
TE
2.3 Sea spray parameters for an offshore rig
CE P
For simplicity the spray bulk velocity components in Eq. 12 are calculated assuming that all the droplets in the spray cloud (both impact-generated and wind-generated) are moving in straight lines and in the same direction as the wind velocity. The spray bulk velocity
(41)
AC
components are then given by:
is the free horizontal spray velocity component normal to the cylinder axis just before impact. In the model we use the approximation
where
is the free
logarithmic wind speed at the height Z above the mean sea level. For the local water catch rate calculation we also assume that the droplets of impact-generated spray follow straight lines, i.e. the local collection efficiency is
. Assuming for
example the standard condition of section 4.2 the overall collection efficiency of median
ACCEPTED MANUSCRIPT 21
volume droplets at Z=6.75 m becomes 0.92 (at this level
from Eq.47 below).
For the much smaller droplets in wind-generated spray the droplet trajectories near the
T
cylinder are calculated from the equation of droplet motion (Horjen, 2013).
IP
Based on the field results from the drilling rig Treasure Scout operating in Norwegian waters
SC R
(Jørgensen 1984,1985) and droplet kinematics the following time-averaged spray mass flux
(42)
MA
and A and B are given by
D
where
NU
formula was presented in Horjen and Vefsnmo (1985):
TE
(non-dimensional)
CE P
is the liquid water content (
of a single spray at the wave crest and
is the
duration of a single spray. Obviously the value of A will depend on various parameters such
AC
as wind speed, sea state and platform response to incoming waves. For simplicity we will however in the model use the average of four values of A (varying between 5.16 and 6.44) obtain by adjusting the theoretical formula to field measurements: . Note that the constant A has been modified somewhat from previous works in order to agree more closely with experimental data. The water flux was registered mainly on various locations on the rig columns. The mean spray mass flux was estimated by measuring the mass of water in a spray collector after a time interval with approximate constant wind speed (one to several hours), and divide this mass with time. The spray mass flux calculated in that way will have a contribution from both impact- and windgenerated spray, but since the wind speed was not extremely large (from 5 to 9 Beaufort) we
ACCEPTED MANUSCRIPT 22
will assume that the spray mass flux given by Eq. 41 is solely impact-generated. Fig. 3 shows the vertical distribution of model spray mass flux calculated from Eq.41 for four wind speeds.
T
Note that the vertical axis is non-dimensional.
IP
When the wind speed is below 15 m/s wind-generated spray mass flux is neglected; above 15
SC R
m/s a formula given in Horjen, 1990, 2013, which is based on various measurements from the wave crest up to 7 m above the msl (Wu et al. (1984), Monaham (1968), Lai and Shemdin
NU
(1974) and Preobrazhenski (1973)) is used. Fig. 4 shows the vertical distribution of windgenerated spray mass flux for four wind speeds assuming Norwegian waters wave conditions ,
MA
(significant wave height given by Eq. 44 below). Note that above the level ξ = 1, i.e.
the spray mass flux is almost constant. This is due to the findings that droplets below a critical
D
size, depending on wind speed, stay in suspension independent of height above the sea surface
TE
(Wu et all, 1984). Comparing Fig. 3 and Fig. 4 we see that the wind-generated spray contribution to the total spray mass flux is rather small, at least for wind speeds up to 30 m/s.
CE P
For example, at the level
the mean impact-generated spray mass flux for
is about ten times as large as the wind-generated spray mass flux for
.
AC
For most rig icing cases below the main deck the omission of wind-generated spray in the model makes almost no difference, not only because the mass flux is relatively low but also because the majority of the droplets follows more or less the streamlines around the object. At great heights, for example on the derrick, the theoretical wind-generated spray mass flux may be higher than the mean impact-generated spray mass flux. The validity of the formula used for the wind-generated spray flux is however questionable above 7 m above msl. Using the spray volume concentration distribution function of Jones and Andreas (2012) Kulyakhtin and Tsarau (2014) did also conclude that wind-generated spray is unlikely to create significant contribution to icing on a rig (they tested another Norwegian model called MARICE).
ACCEPTED MANUSCRIPT 23
Due to wave washing ice may be removed mechanically under the level
of maximal run
up. According to Mitten (1994) this level may be approximated by the following formulas
IP
T
when waves are colliding with a column of diameter D:
SC R
(43)
where k is the wave number. From linear wave theory we have:
NU
(44)
MA
U10=10 m/s
IMPACT-GENERATED SPRAY MASS FLUX
D
U10=20 m/s U10=25 m/s
CE P
TE
3
2
1
0
U10=15 m/s
AC
Non-dimensional height ξ
4
0
1
2
3
4
5
Spray mass flux G (g/m2s)
Fig. 3. Impact-generated mean sea spray mass flux based on droplet kinematics and measurements on the oil rig Treasure Scout operating in Norwegian waters.
ACCEPTED MANUSCRIPT 24
U10=15 m/s
WIND-GENERATED SPRAY
U10=20 m/s
2
T
U10=25 m/s
IP
U10=30 m/s
SC R
1
0.5
0 0
100
200
300
NU
Non-dimensional height ξ
1.5
400
500
600
MA
Spray mass flux (10-3g/m2 s)
TE
D
Fig. 4. Wind-generated sea spray mass flux in Norwegian waters.
CE P
Spray frequency is an important parameter in the model. For vessel icing it has been estimated that in the mean every second significant wave produces spray ( Zakrzewski, 1987). The
AC
spray producing mechanism is however quite different on a semi-submersible oil platform. As part of the Offshore Icing project a few spray frequency registrations were made on the semisubmersible platform Treasure Scout (Jørgensen, 1984, 1985) operating in Norwegian waters. For this area the following empirical relations for significant wave height and period are found, based on numerous observations (Jørgensen, 1985 (wave height), Horjen, 1990 (wave period)): (45) (46)
ACCEPTED MANUSCRIPT 25
The 10 results for the 1984 registrations and the 4 results for the 1985 registrations are shown in Fig.5 together with the theoretical results based on the assumption that every and every
SC R
SPRAY FRQUENCY ON TREASURE SCOUT 8 7
N=1/Ts N=0.5/Ts Obs 1984 Obs 1985
NU
6
MA
5 4 3
D
Spray frquency (1/min)
IP
T
second significant wave is producing spray.
2 5
10
TE
0
15
20
25
CE P
Wind speed U10 (m/s)
Fig. 5. Impact-generated sea spray frequency N based on measurements on the oil rig
AC
Treasure Scout and theoretical results assuming that every significant spray generates spray (upper line) or every second significant spray generates spray (lower line).
There is large scatter in observed spray frequencies on Treasure Scout and no obvious relation between N and
can be read from the data in Fig. 5. A majority of the data seems however
to fall between the curves of collision frequencies of every significant and every second significant wave. In the numerical examples to follow we have chosen to assume that every significant wave produces spray on an offshore rig. A sensitivity analysis has been presented in section 4.2 showing the effect of varying the spray frequency in the model.
ACCEPTED MANUSCRIPT 26
The duration of each spray is calculated from the assumption that a single spray is generated along a distance equal to half the significant wave length (Horjen, 2013). This gives
IP
is the significant wave length.
SC R
where
T
(47)
With the assumptions made above it follows that A in Eq. 41 is mathematically a true
NU
constant, i.e. independent of wind speed, only if
. In the model the spray mass
flux of a single spray is a derived quantity obtained by dividing Eq. 41 by
. A better
MA
approach would obviously be to use directly field measurements of either liquid water content or spray mass flux of a single spray (total mass flux during a time period divided by number
TE
D
of sprays) at various heights and wind speeds.
The distribution of droplet sizes for various wind speeds and heights above mean sea level
CE P
(msl) were measured on the rig Treasure Scout as well. From these results the following relation was obtained, assuming the Best distribution (Best, 1950, 1951) of droplet sizes
(48)
AC
(Horjen and Vefsnmo, 1985):
where a = 1.31 mm and b = 0.81 mms. In order to calculate the heat exchange between the droplet cloud and the structure we must know the bulk spray temperature upon impact, which depends on the size distribution of spray droplets and the various droplets flight times
. The
calculation of impact-generated spray cooling during flight is based on the assumption that all the droplets in the spray cloud have the same size as the median volume droplet. We also assume that half of the spray cloud is moving upwards before impact without undergoing any cooling at all. The bulk impact-generated spray temperature is then the mean value of descending
droplet temperature (using a method given in Horjen, 1990) and the sea
ACCEPTED MANUSCRIPT 27
surface temperature. For simplicity we assume that this temperature does not change due to the possible contribution from wind-generated spray. Wind-generated spray is however taken
T
account of between impact-generated spray productions. For simplicity we again assume that
IP
half of the spray cloud is moving upwards without undergoing any cooling while all the
SC R
droplets in the descending part of the spray cloud attain the air temperature before impact ( partly freezing of supercooled droplets at the moment of collision is not considered in this analysis). The bulk wind-generated spray temperature upon impact is then the mean value of
NU
the sea surface and the air temperature. Probably a more realistic value is only slightly higher than the air temperature, but since the icing contribution from wind-generated spray is
MA
normally insignificant we don’t need to calculate the heat exchange due to this type of spray
D
very exact.
TE
It is reasonable to assume that a characteristic time the median volume droplet in impactgenerated spray is moving in the air is some function of the sea state. As a preliminary
CE P
hypothesis we will assume that the droplet flight time is equal to the horizontal flight time over a distance corresponding to a quarter of the significant wave length. Approximately a
(49)
AC
characteristic flight time then becomes:
This is of course a crude approximation since in reality the various droplets in a spray cloud have been in the air for quite different times before they collide with the structure.
3. THE NUMERICAL METHODS OF THE EXTENDED MODEL A finite difference method is used to solve the icing equations. The method chosen here is an extension of the two-step Lax-Wendroff scheme (see for example Lax and Wendroff, 1964)
ACCEPTED MANUSCRIPT 28
as proposed by Richtmyer (1963). The method is simple enough to be run on a personal computer. The icing surface is divided into nodes (I,J) representing the s and ζ-position. We and space steps
in s and ζ-direction
and
T
must choose appropriate time step
IP
respectively. In the original Lax-Wendroff and Richtmyer numerical scheme there are mid-
SC R
point evaluations. These have been avoided in the ICEMOD2.1 scheme by using the following approximations for an arbitrary quantity A:
NU
(49)
MA
Similar approximations are used for the J-variations. The same method was used in Horjen (2013) to solve two-dimensional icing on cylinders (D = 0.05-0.50m) using only the two first
TE
D
icing equations (Eq. 1 and Eq. 2). Stable solutions were then obtained using an angle step of and, for most of the cases studied, a time step of 0.016 s. With four equations the
CE P
numerical scheme becomes twice as large. For stable solutions the angle step for some of the cases may also now be as large as
, but the time step must at least be one order smaller
AC
than the time step mentioned above. It is hence obvious that the computer time using ICEMOD2.1 increases substantially compared with ICEMOD2. For the large column of the SEDCO platform (D = 9.15 m) the space steps were ( put
) and
. The mesh length ratio
and (D = 1,22 m) and
is hence 0.70. For the trusses we
which give the following model space steps: ,
,
(D = 1,52 m). The mesh length
ratio is now 1.11 and 1.09 respectively. Values of the space steps were for most of the cases studied here chosen in such a way that the mesh length ratio did not become very much different from 1. The effect of changing space steps for a particular icing case is discussed in section 4.2. Stable solutions for the whole icing periods of the SEDCO rigs were obtained for
ACCEPTED MANUSCRIPT 29
the large column when the time step was 0.001 s and for the truss when the time step was 0.0005 s.
T
The equations of droplet motion are non-linear, ordinary differential equation. A different
IP
numerical procedure must then be used to obtain a solution. In this article a fourth order
SC R
Runge-Kutta method, combined with the Adams-Moulton method, is applied (see Cheney and Kincaid , 1985).
NU
In the next chapter we will first analyze some real icing cases on two drilling rigs. Sensitivity tests will then be presented to examine the effects of some of the input parameters of the icing
D
MA
model. Finally the influence of brine film motion upon icing calculations has been discussed.
CE P
4.1 Two rig icing cases
TE
4. NUMERICAL RESULTS AND DISCUSSION
Unfortunately there are only a few field results which are suitable for model comparison. As
AC
part of a joint Canadian/Norwegian icing project seven of the best documented rig icing events between 1980 and 1990 in North American waters were studied (Brown and Horjen,1989). But even for these cases the measurements are not very exact; estimates of ice thicknesses are based on photos and observations from the deck. We will now compare ICEMOD2.1 results with two of the seven documented icing measurements: the semi-submersible offshore rigs SEDCO 708 and SEDCO 709. These two events were also presented in Brown et all, 1988. Weather data for the two icing periods are shown in an appendix. For simplicity we have in the model assumed that spray and wind come from the same direction during the whole icing period, normal to the cylinder axis (as in Fig. 1). In order to reduce computer time, model calculations are made for the first
ACCEPTED MANUSCRIPT 30
spray periods and ice thickness and accreted mass are then extrapolated to the time ) by multiplying the results with the factor
T
.
IP
(50)
SC R
SEDCO 708, January 4-6, 1983
This icing event was documented by Minsk (1984) when the rig was operating in the Gulf of
NU
Alaska at 56.3 N and 161.8 W. Ice thickness profiles on a windward side of a column and diagonal truss were estimated by David Minsk from photos (the report contains some of these
MA
photos). The photos were taken between 1310 LST (Local Standard Time) and 1450 LST on January 6 and model input data were hence terminated at 1500 LST that day. Totally 60
D
hourly input data set were used starting at 0400 LST on January 4. The following
TE
meteorological and oceanographic registrations were made: air, dew-point and sea surface
CE P
temperature, sea surface salinity, air pressure (mb), significant wave height and period, wind direction, visibility and type and intensity of precipitation. The last four registrations were not used in model simulations. Calculation for a column and a truss were made from Z=12 m
AC
down to the nearest integer of
.
1. Large column D=9.15 m Time evolution of horizontal thickness profiles at Z=5.0 m for observation 36, 39…57, 60 are shown in Fig. 6. We observe that maximum ice thickness for the last observation is about 54 from the stagnation line. Fig. 7 shows the time evolution of vertical ice thickness distribution along this angle assuming the air entrapment is 5%. The model ice density for each separate 1h icing period for the 25 observations from obs.36 to obs.60 varied from 905
to 923
. This may be compared with the measured ice densities of accreted ice on the bulkhead and the deck of a vessel reported by Ryerson and Gow (2000). They found that the
ACCEPTED MANUSCRIPT 31
density was in the range 690-920 kg/m3. Assuming no air is entrapped in the ice and brine drainage has not taken place the density will normally be above 950 kg/m3. Due to the
T
uncertainty of the air volume it would be better test the model against accreted mass instead
SC R
IP
of ice thickness. But this is of course difficult for a single structure on a rig.
SEDCO 708, VERTICAL COLUMN, D=9.15 m, ICEMOD2.1
Ice thickness (cm)
8
MA
6
NU
10
-72
-54
-36
-18
2
0 0
CE P
-90
TE
D
4
18
36
54
72
90
Angle from stagnation line (deg)
AC
Fig. 6. Model ice thickness cumulative horizontal profiles for obs. 36,39,42…..57,60 (3h between each observation) at the height Z=5.0 m.
ACCEPTED MANUSCRIPT 32
SEDCO 708, VERTICAL COLUMN, D=9.15 m, ICEMOD2.1 12
T IP
10 9
SC R
8 7 6 5 4 1
2
3
4
5
6
7
8
9
10
MA
0
NU
Height above msl (m)
11
Ice thickness 54o from stagnation line (cm)
TE
D
Fig. 7. Model ice thickness cumulative vertical profiles -54⁰ from the stagnation line for obs.
CE P
36,39,42…..57,60 (3h between each observation).
In Fig.8 model results just after observation 60 are compared with the estimated results from
AC
Minsk (1984) for the large column. We observe that between 5 and 8 m above msl model ice thickness is 4-6 cm larger than reported values. But this is not too bad considering the rather long icing period and the uncertainty of several of the input parameters, for example spray mean mass flux, spray frequency, duration of one spray, spray flight time and heat transfer coefficient. The discrepancy between model and observed results documented in Fig. 8 may partly be explained by a wrong height distribution of the spray mass flux. We must however expect that model ice thickness at a certain location is larger than observed ice thickness due to the assumption of constant wind direction. Ideally the icing program should also be run for the whole icing period of 60 hours and not only for 10 spray periods (typically 1-1.5 min) for
ACCEPTED MANUSCRIPT 33
each observation. This is however too time-consuming for a PC. The field measurements are also not very exact and some ice may have fallen off during the icing period.
IP
12
ICEMOD2.1, obs 60 Estimated
SC R
11 10 9 8
NU
7 6 5
MA
Height above msl (m)
T
SEDCO 708, VERTICAL COLUMN, D=9.15 m
4 3 2 2
3
4
D
1
5
6
Ice thickness (cm)
7
8
9
10
TE
0
CE P
Fig. 8. Comparison between maximum model ice thickness at each level and estimated (from photos) maximum thickness just after observation no 60.
AC
2. Diagonal (70 ) truss
Ice thickness profiles at Z=4.2 m in the AB-plane (Fig.1) for obs. 33,36,…57,60 are shown in Fig.9. We now observe that maximum ice thickness for the last observation is 36-54 from the stagnation line. Fig.10 shows the time evolution of vertical ice thickness distribution along the angle 36 for the same set of observations. In Fig.11 model results for
just after
observation 60 are compared with the estimated results (also from Minsk, 1984) for a diagonal truss. Since the type of diagonal was not specified model results for both the 70 (D = 1.22 m) and the 50 (D = 1.52 m) truss are included in Fig. 11. The differences between the model results are very small. We observe that above 4.5 m the model thickness is larger than the estimated thickness (0 - 6 cm).
ACCEPTED MANUSCRIPT 34
SEDCO 708, 70o DIAGONAL TRUSS , D=1.22 m, ICEMOD2.1 14
T
Ice thickness (cm)
12
IP
10
SC R
8 6
NU
4 2
-72
-54
-36
MA
0
-90
-18
0
18
36
54
72
90
D
Angle from stagnation line (deg)
TE
Fig. 9. Ice thickness cumulative profiles in AB-plane for obs.36,39,….57,60 at Z=4.2 m.
CE P
SEDCO 708, 70o DIAGONAL TRUSS , D=1.22 m, ICEMOD2.1
10
AC
Height above msl (m)
12
8
6
4
2 0
2
4
6
8
10
12
14
Ice thickness (cm)
Fig. 10. Model ice thickness cumulative vertical profiles for obs. 36,39,42…..57,60 at φ= -36
ACCEPTED MANUSCRIPT 35
SEDCO 708, TRUSSES D=1.22 m AND D=1.52 m
D=1.22 m
D=1.52 m
12
Estimated
11
T IP
9 8
SC R
7 6 5 4
3 2 3
6
9
12
15
18
21
MA
0
NU
Height above msl (m)
10
Ice thickness (cm)
D
Fig. 11. Vertical profile of model maximum ice thickness on D=1.22 m (υ=70 ) and D=1.52
CE P
TE
m (υ=50 ) trusses and estimated (from photos) maximum ice thickness just after obs. 60.
AC
SEDCO 709, January 24-26, 1986 This icing event was documented by the drilling platform ice and weather observer (MacLaren Planseach,1986) when the rig was operating on the Scotian Shelf at 43.7 N and 59.9 W. Ice thickness profiles on a windward column and diagonal truss were estimated from photos. A number of photos were taken at 1800 GMT on January 25 when maximum ice build-up was reported. Model calculations are based on 25 input data sets from 1400 GMT on January 24 to 1600 GMT on January 25. Ice accretion before and after this time interval is negligible. Calculations for a column and a truss were made from the nearest integer to to Z=8.0 m. 1. Large column, D = 9.15 m
up
ACCEPTED MANUSCRIPT 36
Time evolution of horizontal thickness profiles at Z=2.5 m for observation 3, 5…23,25 are shown in Fig.12 We observe that maximum ice thickness is now about 36 from the
T
stagnation line. Fig.13 shows the time evolution of vertical ice thickness distribution along
IP
this angle (obs. 1 has here been added to show the effect of wave washing). In Fig.14 model
SC R
results just after observation 25 are compared with the estimated results from photos. We now observe that model ice thickness is larger than estimated values for
(0-6.5 m).
NU
SEDCO 709, LARGE COLUMN, D=9.15 m, Z=2.5 m, ICEMOD 2.1 12
Ice thickness (cm)
MA
10
-72
-54
AC
-90
CE P
TE
D
8 6
4 2 0
-36
-18
0
18
36
54
72
Angle from stagnation line (deg)
Fig. 12. Model ice thickness cumulative horizontal profiles for obs. 3,5,7…..23,25 (2h between each observation).
90
ACCEPTED MANUSCRIPT 37
SEDCO 709, LARGE COLUMN, D=9.15 m, ICEMOD2.1 5
T
4.5
IP SC R
3.5
3 2.5 2
NU
Height above msl (m)
4
1 0
2
4
MA
1.5
6
8
10
12
D
Ice thickness (cm)
TE
Fig. 13. Model ice thickness cumulative vertical profiles -36⁰ from the stagnation line for obs.
CE P
1,3,5,7…..23,25 (2h between each observation).
SEDCO 709, LARGE COLUMN, D=9.15 m
7
Height above msl
6
ICEMOD2.1, obs 25 Estimated
AC
8
5 4 3 2 1 0 0
1
2
3
4
5
6
7
Ice thickness (cm)
8
9
10
11
12
ACCEPTED MANUSCRIPT 38
Fig. 14. Comparison between maximum model ice thickness and estimated (from photos)
IP
T
maximum ice thickness just after observation 25.
SC R
2. Diagonal (70 ) truss
Ice thickness profiles at Z=2.0 m in the AB-plane (Fig.1) for obs. 3,5,…23,25 are shown in Fig. 15. We now observe that maximum ice thickness is about 36 from the stagnation line. for obs.
NU
Fig.16 shows the time evolution of vertical ice thickness distribution for
MA
3,5….23,25. In this case only the maximum thickness was documented from photos: 15.0 cm at the height 2.0 m above msl. Maximum height of ice accretion was reported to be 7.0 m.
TE
D
SEDCO 709, 70⁰ DIAGONAL TRUSS, D=1.22 m, Z=2.0 m, ICEMOD2.1 16
12 10 8
AC
Ice thickness (cm)
CE P
14
6 4 2 0
-90
-72
-54
-36
-18
0
18
36
54
72
Angle from stagnation line (deg)
Fig. 15. Model ice thickness cumulative profiles in the AB-plane for obs. 3,5,7,…..23,25
90
ACCEPTED MANUSCRIPT 39
SEDCO 709, 70⁰ DIAGONAL TRUSS, D=1.22 m,
T IP
4
SC R
Height above msl (m)
5
3
NU
2
0
2
MA
1 4
6
8
10
12
14
16
D
Ice thickness (cm)
CE P
TE
Fig. 16. Model ice thickness cumulative vertical profiles for obs. 1,3,5,7…23,25 at φ= -36 .
Icing event
AC
A summary of the two icing events presented above is shown in Table 1. Icing object
Corner column,
Max. thickness
Height of max.
Max. height of
(cm)
thickness (m)
icing (m)
Model
Model
Model
Obs.
Obs.
Obs.
9.7
7.5
5.0
4.0
~12.0
9.0
13.1
18.0*)
4.0-4.4
3.5*)
~12.0
9.0
13.0
18.0*)
3.6-3.8
3.5*)
~12.0
9.0
SEDCO D=9.15 m 708
70⁰ Truss, D=1.22 m 50⁰ Truss, D=1.52 m
ACCEPTED MANUSCRIPT SEDCO
Corner column,
709
D=9.15 m
10.0
2.5
3.5
~5.5
7.0
15.3
15.0
2.0
2.0
~5.4
7.0
15.2
-
2.0
~5.4
-
T
70⁰ Truss,
11.5
IP
40
D=1.22 m
-
SC R
50⁰ Truss, D=1.52 m
Table 1. Comparison of observed and predicted ice accretion thickness data. *) Type of
MA
NU
truss for ice thickness observations was not specified.
Although this cannot be compared with direct observations a summary of model ice mass on a
TE
D
column and a truss and corresponding “overall freezing fraction” (ratio of total accreted ice mass and total mass of spray hitting the cylinder during all the observations, excluding for
Icing event
CE P
each observation icing and water catch rate below the level Icing object
Corner column,
Overall freezing fraction
28451
0.28
70⁰ Truss, D=1.22 m
7832
0.28
50⁰ Truss, D=1.52 m
11170
0.27
16618
0.66
70⁰ Truss, D=1.22 m
3595
0.65
50⁰ Truss, D=1.52 m
5089
0.65
AC
SEDCO 708
Ice load (kg)
) is given in Table 2.
D=9.15 m
SEDCO 709 Corner column, D=9.15 m
Table 2. Predicted final ice mass and overall freezing fraction.
ACCEPTED MANUSCRIPT 41
4.2 Sensitivity analysis Sensitivity analyses will now be presented for three of the icing model parameters: air
T
temperature, wind speed at 10 m height and spray frequency. The constant and the variable
Constant parameters Variable parameters ) sensitivity
NU
1: Air temperature (
) sensitivity
(cont. spray)
D
3: Spray frequency (N) sensitivity
MA
2: Wind speed (
SC R
Case number
IP
parameters for each test are given in Table 3.
TE
Table 3. Parameter values to be used in the sensitivity tests. Bold types in the right
CE P
column is standard condition.
AC
We will consider icing on a column of diameter 5.0 m extending between the levels Z = 6.0 m and Z = 7.5 m. The following parameters are kept constant: sea surface salinity
,
relative humidity f = 0.8, atmospheric pressure p = 1000 mb and significant wave height Hs = 5.0 m. The significant will be calculated from a formula given in Horjen (1990): (51) This formula is based on 240 wave data (
collected on the rigs Bow drill 1 and 2
and SEDCO 709 operating on the Grand Banks area off the coast of Newfoundland. With the prescribed value of the significant wave height the spray period becomes 7.4 s. The angle and height step is 4.5 and 0.25 m respectively and the time step is 0.0005s. Model calculations
ACCEPTED MANUSCRIPT 42
are again made for the first
spray periods and ice thickness and accreted mass are
then extrapolated to 24 hours (
) by multiplying the results by the factor E
T
(Eq. 50). Wave washing of accreted ice is now neglected. The overall freezing fraction is then
IP
given by:
SC R
(52)
) and the icing area of the
NU
The integration is over the period of icing calculations ( cylinder.
MA
Vertical profiles at the stagnation line for Case no 1 (air temperature sensitivity) is shown in Fig.17. As expected more ice accretes when the air temperature decreases. This is mainly due
D
to the increase in convective and evaporative heat loss and a decrease in arriving spray
to
CE P
cylinder varies from
TE
temperature (the bulk impact-generated spray temperature upon impact at the top of the for Case no 1). An inspection of Fig. 17 gives an
approximate maximum sensitivity of 6.58 cm per 4 the range of -12
to -8
, i.e. 1.6 cm/ , for air temperature in
. Fig. 18 (vertical profiles) and Fig. 19 (horizontal profiles) show
AC
the effect of varying the wind speed. In our case ice accretion increases with increasing wind speed (but this is not generally valid). At the cylinder base the ice thickness for this particular test increases approximately 0.5 cm when the wind speed increases with 1 m/s in the range 17-25m/s. Although maximum Nusselt number is located at
(Eq. 39) maximum
ice accretion at the cylinder base is from Fig.19 about 30⁰ from the stagnation line. This is due to several mechanisms, for example a relative high heat input from the spray cloud (to be discussed in more detail in the next section) and the interaction between icing intensity and brine film salinity. Note that there are some small bumps around the locations of maximum heat transfer coefficient (but this may just be a coincidence or due to numerical noise).
ACCEPTED MANUSCRIPT 43
-8 deg
CASE 1: AIR TEMPERATURE SENSITIVITY
-10 deg 7.5
-12 deg
T IP
7
SC R
6.75 6.5 6.25 6 0
4
8
NU
Height above msl (m)
7.25
12
20
MA
Ice thickness (cm)
16
TE
D
Fig. 17. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 1.
17 m/s
CASE 2: WIND SPEED SENSITIVITY
CE P
7.5 7.25 7
25 m/s
AC
Height above msl (m)
21 m/s
6.75 6.5 6.25
6 0
3
6
9
12
15
18
Ice thickness (cm)
Fig. 18. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 2.
ACCEPTED MANUSCRIPT 44
17 m/s
CASE 2: WIND SPEED SENSITIVITY, Z=6.0 m
21 m/s 24
T
25 m/s
IP
Ice thickness (cm)
18
SC R
12
0 -90
-60
-30
0
NU
6
30
60
90
MA
Angle from stagnation line (⁰)
TE
D
Fig. 19. Horizontal ice thickness profile at the stagnation line extrapolated to 24h, Case no 2.
CE P
From the discussion in section 2.4 the observations of spray frequency on an offshore rig is rather deficient. A sensitivity analysis varying this parameter is hence included assuming
AC
standard condition. Somewhat surprisingly at the first sight it turns out that for the three frequencies studied,
( continuous spray, i.e. single
sprays that follow immediately after each other) the total 24 h accreted mass of ice are more or less equal: 1346 kg, 1286 kg and 1463 kg respectively. Overall freezing fraction was 76%, 71% and 83% respectively. The stagnation line ice thickness profiles (Fig. 20) are also almost the same for the three cases except at the upper level of the cylinder; here the ice thickness for continuous spray is much larger than for the other two cases. An attempt to explain this result is made below: At the top
of the cylinder one of the boundary conditions of the model is X=0. It
ACCEPTED MANUSCRIPT 45
then follows that the ice thickness at some point on the cylinder top is given by (Horjen,2013):
IP
) is the ice density close to the icing surface. From this equation it follows that
SC R
where
T
(53)
for periodic spray the ice will grow only during spray impingement (normally
). For
the Case 3 test the spray period was 8.30 s and 16.60 s for spray frequencies
. This
NU
respectively. The single spray duration was from Eq. 46
and
MA
means that for the major part of the icing period there is no growth in ice thickness at . During continuous spray, however, ice will accumulate during the whole icing (during spray
D
period with an icing intensity at the stagnation point of
and
TE
impingement the icing intensity for the high and low spray frequency is
CE P
respectively).
Contineous spray N=1/Ts
CASE 3: SPRAY FREQUENCY SENSITIVITY
N=0.5/Ts
AC
7.5
Height above msl (m)
7.25
7
6.75 6.5 6.25 6 0
4
8
12
16
Ice thickness (cm)
. Fig. 20. Vertical ice thickness profile at the stagnation line extrapolated to 24h, Case no 3. .
ACCEPTED MANUSCRIPT 46
From Fig. 20 it follows that the time integral of R extrapolated to the time
is almost
independent of spray frequency. An expression of the extrapolated integral at the stagnation ) of spray
T
line may be obtained from Eq.1 for periodic spray and an integral number (
IP
periods (formally using Fourier analysis for the local water catch rate):
is the heat loss by evaporation,
NU
where
SC R
(54)
is the extrapolated time,
does not change with spray
MA
defined by Eq.50. In the model the mean spray mass flux
and E is
frequency (since A in Eq. 41 is put equal to a constant). The first term on the right hand side
D
of Eq. 54 is hence independent of spray frequency. If the local brine film thickness after
TE
sprays is rather thin (let’s say smaller than 0.1 mm) this term will be the dominant one. For
CE P
Case 3 the brine film thickness after 10 spray periods at the mean height Z=6.75 m is 38μm for spray generated by both every significant wave and every second significant wave. For the
AC
continuous spray case the brine film thickness is literally zero. From Eq.54 it is likely that the accreted ice mass depends more on spray frequency for a thick brine film. To test this hypothesis we will run Case no 3 once more, but now with a mean spray mass flux 25 times larger than the value given by Eq. 41. Stagnation line ice thickness profiles are shown in Fig. 21 (for the continuous spray case there is no ice deposit here). Ice thickness now decreases about 50% when the spray frequency decreases to half the original value. Total accreted ice mass decreases 17%. At Z=6.75 m the brine film thickness after 10 spray periods is about 360 μm for the high frequency case ( the low frequency case (
).
) and about 460 μm for
ACCEPTED MANUSCRIPT 47
N=1/Ts
CASE 3: SPRAY FREQUENCY SENSITIVITY, 25*NSF
N=0.5/Ts
7.5
T IP
7
SC R
6.75 6.5 6.25 6 0
2
4
6
NU
Height above msl (m)
7.25
8
10
12
14
MA
Ice thickness (cm)
D
Fig. 21. Vertical ice thickness profile at the stagnation line extrapolated to 24h when the
CE P
TE
spray mass flux is 25 times the “normal” value (NSF) given by Eq.41 .
As mentioned before a more realistic spray model would be the use of the spray mass flux
AC
or the LWC of a single spray as input to the icing model. The time-average impact-generated spray mass flux is then generally given by (compare Eq. 41): (55) where
is a function depending on the sea state, height above msl and type of
marine structure. Assuming that the spray frequency is proportional to the meeting frequency of significant waves (actually this is most probable for a fixed structure) it then follows that the first term on the right hand side of Eq. 54 is independent of spray frequency only if where generally
. Expressing for example the significant
wave height and period by Eq. 44 and Eq. 45 this condition gives a LWC near the wave crest
ACCEPTED MANUSCRIPT 48
proportional to
. This is a rather special condition which we cannot expect to be valid.
Nevertheless the assumption made in the model that the mean spray mass flux does not
T
change with spray frequency has proven to give reasonable approximations of ice growth rate
IP
on some marine structures.
SC R
Details of the time evolution of icing intensity at the mean height and stagnation line of the cylinder are shown in Fig.22. We observe that the variations of icing intensity are quite
NU
different for the three spray frequency cases, but the time-average value of the icing intensity is more or less equal. Note that the left hand side of Eq. 54 may be written
, which
MA
again shows that total ice load per unit area for a specified icing period do not vary much with
D
spray frequency.
N=1/Ts
TE
ICING INTENSITY VARIATIONS FOR CASE 3, Z=6.75 m
CE P
Continuous spray
1.9
AC
Icing intensity (g/m 2 s)
2.4
N=0.5/Ts
1.4
0.9
0.4 0
10
20
30
40
50
Time (s)
Fig. 22. Detail of icing intensity time variations at the stagnation line and mean height of the cylinder. Six and three spray cycles are considered for the two periodic spray cases. The spray-on periods are the horizontal or monotonically increasing parts of the curves.
ACCEPTED MANUSCRIPT 49
For the continuous spray case at the level Z = 6.75 m the stagnation point spray mass flux is and the icing intensity is
(see Fig.22). The , which means that the assumption is justified
IP
that droplet salinity upon impact is equal to the sea surface salinity
T
median droplet diameter is from Eq. 47
SC R
(Horjen,1990). Neglecting the convective transport of the brine film and mass loss by evaporation the theoretical brine film salinity is simply given by (see section 4.3):
NU
(56)
close to the ICEMOD2.1 result of
, which is
MA
Putting now the parameters given above in this equation we find .
D
There are also other model parameters which are based on very limited or none experimental
TE
documentations, for example the droplet flight time. Tests have been done making the flight
CE P
time equal to half the wave length and the eighth part of the wave length. Total ice load for the standard condition then becomes 1351 kg and 1343 kg respectively, i.e. the variations with flight time seems to be small for normal sea spray flux (NSP) given by Eq. 41. The bulk
2
AC
impact-generated spray temperatures upon impact at the top of the cylinder were
and
for the two flight times considered.
Finally some diagrams showing the time evolution of icing intensity, brine film salinity, brine film thickness and cumulative accreted ice mass per unit area at the stagnation line and mean height Z=6.75 m of the cylinder are presented in Fig. 23-26 for Case no 1. The diagrams show the variation of the four quantities during the first 10 or 20 spray cycles for wind speed 21 m/s and air temperatures -8.0
and -12.0
. The spray period is 8.30 s (from Eq.45) and, as
mentioned above, the duration of a single spray is 2.56 s. Overall freezing fraction for -8.0 and -12.0
is 67% and 83% respectively.
ACCEPTED MANUSCRIPT 50
From Fig. 24 we see that after some few spray cycles the brine film salinity more or less fluctuates around 54
before any drainage has taken place. This result may be compared with the
T
about 18
when the air temperature is -8 . Corresponding ice salinities will be
IP
vessel icing investigations of Ryerson and Gow (2000). They found that the salinity of ice
7.5
SC R
accreted on the bulkhead varied between 11.5‰ and 14.7‰. The air temperature was from to -8.5 . Time evolution of the brine film thickness is shown in Fig.25. This is the only
example with 20 spray periods. The larger number of spray periods is chosen in order to show
NU
that the brine film thickness normally stabilize after a sufficient time period. In this case we
43μm for
MA
see that after some time the brine film thickness seems to fluctuate around 63μm and and
respectively. at the
D
Fig. 26 show the time evolution of the cumulative ice mass accretion per unit area
TE
stagnation point at Z=6.75 m (integration of the R-function in Fig.23). Due to the harmonic
CE P
wave nature of the icing intensity there is a periodic-linear increase with time of the accumulated ice mass. In this example
is practically proportional with time for the whole
time period considered (10 spray periods) which is actually a condition for the extrapolation
AC
method to give approximate correct answers.
ACCEPTED MANUSCRIPT 51
-8 deg
ICING INTENSITY VARIATION FOR CASE 1, Z=6.75 m
-12 deg
T IP
2
SC R
1.5 1 0.5
0 0
10
20
30
40
50
60
70
80
90
MA
Time (s)
NU
Icing intensity (g/m 2s)
2.5
D
Fig. 23. Time history of icing intensity at the stagnation line and mean height for the first 10
CE P
TE
spray periods for a =-8⁰C and a =-12⁰C.
-8 deg
BRINE FILM SALINITY VARIATION FOR CASE 1, Z=6.75 m 120
AC
Brine film salinity (‰)
105
-12 deg
90 75 60 45 30
0
10
20
30
40
50
60
70
80
90
Time (s)
Fig. 24. Time history of brine film salinity at the stagnation line and mean height for the first 10 spray periods for a =-8⁰C and a =-12⁰C.
ACCEPTED MANUSCRIPT 52
-8 deg
BRINE FILM THICKNESS VARIATION FOR CASE 1, Z=6.75 m
-12 deg
T
70
IP
60 50 40
SC R
Brine film thickness (μm)
80
30 20
0 0
20
40
60
80
NU
10 100
140
160
180
MA
Time (s)
120
D
Fig. 25. Time history of brine film thickness at the stagnation line and mean height for the
TE
first 20 spray periods for a =-8⁰C and a =-12⁰C
CE P
-8 deg
CUMULATIVE ICE MASS ACCRETION FOR CASE 1, Z=6.75 m 160
AC
Accreted ice mass (g/m2 )
140
-12 deg
120 100
80 60 40 20 0 0
15
30
45
60
75
90
Time (s)
Fig. 26. Time history of total accreted ice mass per unit area at the stagnation line and mean height for the first 10 spray periods for a =-8⁰C and a =-12⁰C.
ACCEPTED MANUSCRIPT 53
We will now study in more detail the mutual variation of icing intensity and brine film salinity. This is shown in Fig. 27 for the three last spray cycles when the wind speed is 21 m/s vary periodically synchronized with the
T
and the air temperature is -12⁰C. Both R and
IP
model’s spray on/off. The time period from A to B corresponds to the spray-on condition and
SC R
the time period from C to D corresponds to the model’s spray-off condition. Icing intensity is lowest during spray impingement, between for example A and B. This may be explained by the relatively large heat and mass input from the spray cloud. The temperature of the
. With the assumptions made in the model he bulk spray temperature upon impact then
MA
1.34
droplets just before impact at the upper level Z = 7.5 m of the cylinder is
NU
descending
becomes
. The overall freezing fraction for the whole cylinder
D
(Eq. 52) for the 10 spray periods is about 83% , but during spray impingement from for
TE
example A to B the local freezing fraction is only 18-20%. Brine film salinity has the lowest value at the end of a single spray impact, which is due to both mixing with sea water of lower
CE P
salinity and low icing intensity.
Just after the spray has arrived icing intensity is gradually increasing from the lowest value at
AC
A to a slightly higher value at B. During the same period the brine film salinity is decreasing. Icing intensity is gradually decreasing during the no-spray period from the highest value at C to a lower value at D. Now the brine film salinity is increasing. This reflects the interesting respond brine film salinity changes have upon icing intensity and vice versa. The reason for the abrupt changes in icing intensity immediately before A and immediately after B is the discontinuous nature of the spray mass flux.
ACCEPTED MANUSCRIPT 54
R
ICING INTENSITY AND BRINE FILM SALINITY, 3 SPRAY PERIODS
Sb 2.2
90
C
1.8
70 60
D
1.6
50 40
B
1.4
30 20
A
NU
1.2
55
60
65
70
MA
1
75
Brine film salinity (‰)
IP
T
80
SC R
Icing intensity (g/m2 s)
2
10 0 80
85
D
Time (s)
TE
Fig. 27. Time history of icing intensity and brine film salinity at the stagnation line and
CE P
Z=6.75 m. Detail of Fig. 24 and Fig.25 for the three last spray periods for a =-12⁰C. Spray
AC
on: A to B. Spray off: C to D.
For the sensitivity tests presented above the angle and lateral space steps were
and
giving a mesh length ratio of 1.27. A question arise: how sensitive is the model
for angle/space steps changes? To examine this question we shall run the standard condition test again, but now with four other sets of angle/space steps. The results are summarized in Table 4. We observe that the differences in results using the various sets of angle/space steps are small; compared with test no 1 (sensitivity test condition) maximum change in total accreted ice load is only 5.1 % (test no 4). A disadvantage of using very small angle/space steps is the increase in computer time, running for example test 4 takes for example about 5 times as long time as running test 3.
ACCEPTED MANUSCRIPT
Space
Total ice
“Mean” stag.
“Mean” front
Max. ice
Overall
number
steps
load
line thickness
thickness
thickness
freezing fraction
(cm)
4.5,
0.25
1346
14.4
2
9.0 , 0.25
1349
14.5
3
18.0, 0.25
1330
14.6
4
4,5, 0.125
1278
14.3
5
4.5,
1375
0.30
MA
1
(cm)
14.3
(cm)
(%)
10.8
19.5
76.2
10.7
19.4
76.5
10.1
18.6
75.9
10.8
19.7
77.9
10.7
19.4
75.7
SC R
(kg)
NU
(⁰,m)
IP
Test
T
55
D
Table 4. Comparison of standard condition results, extrapolated to 24h, for four different sets
TE
of space steps. “Mean” thickness is defined to be the mean thickness for the upper, middle is the overall freezing fraction defined by Eq. 52.
CE P
and lower section of the cylinder and
AC
4.3 How important is the inclusion of brine film movement in icing modeling? In this section we will take a closer look on the effect of brine film movement and also find out if the convective acceleration terms have any significance at all. We will consider icing on two vertical cylinders, the first one with the same diameter as the corner columns and the second one with the same diameter as the 70⁰ trusses of the SEDCO platforms. The following parameters are common for the two tests:
Air temperature:
Sea surface temperature:
Sea surface salinity:
ACCEPTED MANUSCRIPT 56
Relative humidity:
Significant wave height:
Air pressure:
Air porosity:
IP
T
. Icing will be calculated for
SC R
The angle step chosen is
. An
extended angle range is used to account for the possibility of brine flowing over to the lee side
NU
of the cylinder at high wind speeds. On the lee side only wind stress and gravity is acting on the brine film. Note that for a laminar boundary flow we can use the Blasius formula (Eq.14)
MA
up to 110⁰ from the stagnation line before instabilities occurs. The Nusselt number equation (Eq. 39) is also valid on some part of the lee side, up to 130 from the stagnation line. For
D
icing calculations at the cylinder shoulders and beyond it was necessary to choose a relatively
TE
small angle step. The time step chosen is the same as for the sensitivity tests (5
s).
CE P
Wave period is again calculated from Eq. 51. Like the sensitivity tests wave washing of accreted ice is neglected. Model calculations are also now made for the first 10 spray periods and ice thickness and accreted mass are then extrapolated to an icing period of 24 hour by
AC
multiplying this result by the factor E (Eq. 50). CASE no 1: D = 9.15 m, Z = 5.0-7.4 m,
, Wind speed:
At the level Z=5.00 m the instantaneous local water catch rate at the stagnation point during spray impingement is
approximately. There is no contribution from wind-
generated spray (overall collection efficiency equal to zero). The bulk spray temperature just before impact at the top of the cylinder is about
. It turns out that both ICEMOD2 and
ICEMOD2.1 give an overall freezing fraction of 37.1%. There is also no difference in total ice load and “mean” (defined under Table 4) windward side ice thickness: 2743 kg and 7.8 cm respectively. The observed wind speeds for the SEDCO708 and SEDCO709 icing cases were
ACCEPTED MANUSCRIPT 57
all less than 25 m/s. It is therefore reason to believe that the results given in section 4.1 would be the same using ICEMOD2.
T
We will now investigate a case with much larger mean water catch rate by multiplying the
IP
normal sea spray flux Eq.41 (NSF) by a factor of 50. This is not unrealistic since other data
SC R
sources of spray measurements indicates spray mass fluxes 10-1000 as large as the formula values of this work (Forest et al, 2005). The instantaneous local water catch rate at the . But even now
NU
stagnation point at Z=5.0 m during spray impingement is now
the two models give almost the same total ice load: 3755 kg and 3649 kg using ICEMOD2.1
MA
and ICEMOD2 respectively. There is however one parameter which is quite different for the two icing models for large spray mass fluxes: the horizontal brine film velocity
. This is
the ICEMOD2 horizontal velocity is about 28% higher than the ICEMOD2.1 result.
TE
υ
D
shown in Fig.28 (horizontal profiles at the cylinder base at the end of spray no 10). For
CE P
The vertical brine film velocity is however almost equal for the two models (Fig.29). Due to the increased brine film thickness at the cylinder shoulders (Fig.30) maximum absolute vertical velocity occurs at these locations (compare Eq. 8). We observe that this value is not
AC
much different for the two components. It has been common practice in icing modeling to only consider the gravity driven component of the brine mass flux (Kulyakhtin and Tsarau, 2014). Although laminar lubrication theory is an idealization this work shows that for a relatively high wind speed it is important to retain both velocity components.
ACCEPTED MANUSCRIPT 58
ICEMOD 2.1
MODEL COMPARISON CASE 1
ICEMOD 2
2400
800 0 -800
-108 -90 -72 -54 -36 -18
0
18
36
-1600
-2400
72
90 108
NU
-3200
54
IP
T
1600
SC R
Brine film velocity ub (mm/s)
3200
MA
Angle from stagnation line (⁰)
Fig. 28. Horizontal profiles of brine film mean lateral velocity across the brine film (defined
D
by Eq.5) at the cylinder base for ICEMOD2.1 and ICEMOD2. Model results refer to the end
CE P
TE
of spray no 10. Spray mass flux is 50 times as large as the normal value (NSF, Eq. 41).
ICEMOD 2.1
-108 -90
-72
AC
Brine film velocity (mm/s)
MODEL COMPARISON CASE 1 VERTICAL BRINE FILM VELOCITY, Z=5.0 m, 50*NSF
-54
-36
ICEMOD 2
0 -18
0
18
36
54
72
90
108
-400
-800
-1200
-1600
Angle from stagnation line (⁰)
Fig. 29. Horizontal profiles of brine film mean vertical velocity at the cylinder base for ICEMOD2.1 and ICEMOD2 at the end of spray no 10. Spray mass flux is 50 times as large as the normal value.
ACCEPTED MANUSCRIPT 59
MODEL COMPARISON CASE 1 BRINE FILM THICKNESS, Z=5.0 m, 50*NSF
ICEMOD 2.1 ICEMOD 2
Brine film thickness ( μm)
1000
IP
T
800
SC R
600 400
0 -108 -90
-72
-54
-36
-18
0
NU
200
18
36
54
72
90
108
MA
Angle from stagnation line (⁰)
D
Fig. 30. Brine film thickness profile at the base of the test cylinder for ICEMOD2.1 and
TE
ICEMOD2 at the end of spray no 10. Spray mass flux is 50 times as large as the normal
CE P
value.
Looking at Fig.28 and Fig.29 one may suspect that the brine film velocity has no effect at all on icing calculations. Perhaps we may as well put
? The icing equations
(57)
AC
then reduce to two simple ordinary differential equations:
Periodic spray,
A further simplification is the assumption of continuous spray. The brine film salinity is then constant in time and we just have to solve two algebraic equations (see Horjen,2013):
(58)
Continuous spray,
ACCEPTED MANUSCRIPT 60
We must here assume that the conductive heat flux to the accreted ice is independent of time, but as mentioned before we have generally put
in this work. Several other icing
T
models resemble this model, for example the old model of Stallabrass (1980) and RIGICE
IP
(Brown and Horjen, 1989). It turns out that for normal sea spray mass flux ICEMOD2.1 and
SC R
the periodic spray model with no brine film movement (called Model0l) results in no difference in the distribution of accreted ice. The continuous spray model with no brine film movement (called Model02) does however depict a somewhat different ice thickness
NU
distribution (Fig. 31 and Fig. 32). We observe that Model02 gives more ice deposit at the
MA
upper part (which again may be explained by Eq. 53) and less ice deposit at the lower part of the cylinder than the other two models. This results in a total ice load for the three models
D
differing not significantly. Maximum ice thickness for all the models is at
ICEMOD 2.1+ MODEL 01 MODEL 02
12 10 8 6
AC
CE P
TE
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, Z=6.20 m, NSF
Ice thickness (cm)
.
4 2 0
-108 -90
-72
-54
-36
-18
0
18
36
54
72
90
108
Angle from stagnation line (⁰)
Fig. 31. Ice thickness horizontal profiles at the mean height for three model versions: ICEMOD2.1, Model0l and Model02. Normal spray mass flux is assumed (Eq. 41).
ACCEPTED MANUSCRIPT 61
ICEMOD 2.1+ MODEL 01
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, ϕ=0⁰, NSF
MODEL 02
T
7.4
IP SC R
6.6 6.2 5.8 5.4
5 0
2
4
NU
Height above msl (m)
7
6
8
MA
Ice thickness (cm)
D
Fig. 32. Ice thickness vertical profiles for Case no 1 at the stagnation line for three model
CE P
TE
versions. Normal spray mass flux is assumed (Eq. 41).
Although the simplified models Model01 and Model02 in many cases approximates icing on a
AC
cylinder quite well these models bring about some logical problems about continuity principles: what happens with the unfrozen brine? Since there is no motion of the brine film we must conclude that unless brine is continuously blown off by wind action the brine film thickness will increase without limit when the freezing fraction is less than 1. The growth of brine film thickness in time for the three models ICEMOD2.1, Model01 and Model02 is shown in Fig. 33. We observe that Model02 gives a linear growth in time while Model01 gives a stepwise growth following the same trend as Model02. ICEMOD2.1 on the other hand gives a brine film thickness which after some spray cycles seems to stabilize and more or less fluctuate around a constant value. Fig. 34 shows the time history of brine film salinity for the same three models. As expected Model02 results in a constant salinity (40.3‰
freezing
ACCEPTED MANUSCRIPT 62
temperature -2.14 ) all the time while the other two models give fluctuating salinity values. Note that for periodic spray the salinity decreases during periods of spray-on and increases
SC R
IP
T
during periods of spray-off.
MODEL COMPARISON CASE 1 BRINE FILM THICKNESS, ϕ=0⁰, NSF
240
MODEL 01
NU
MODEL 02
MA
160 120 80
D
40 0 15
CE P
0
TE
Brine film thickness (μ m)
200
ICEMOD 2.1
30
45
60
75
90
Time (s)
AC
Fig. 33. Model comparison of time history of brine film thickness at the stagnation line and mean height for the first 10 spray periods. Normal sea spray mass flux.
ACCEPTED MANUSCRIPT 63
MODEL COMPARISON CASE 1. BRINE FILM SALINITY, ϕ=0O, NSF
ICEMOD 2.1 MODEL 01
44
MODEL 02
T IP
40
SC R
Salinity (‰)
42
38
34 0
15
30
45
60
75
90
MA
Time (s)
NU
36
D
Fig. 34. Model comparison of time history of brine film salinity at the stagnation line and
CE P
TE
mean cylinder height for the first 10 spray periods. Normal sea spray mass flux.
We will now make the same comparison as above for a very large spray mass flux case, 50
AC
times NSF. The result is shown in Fig.35. ICEMOD2.1 and Model01 give the same profile for . A small amount of brine is pushed over to the lee side which of course cannot be calculated by Model0l. Model02 gives quite a different result with no icing for
.
The constant mean spray mass flux is in this case too large on the major part of the cylinder upwind side to give any icing. For ICEMOD2.1 and Model0l almost all icing takes place during the spray-off part of the spray period. A summary of the results for both normal and very large spray mass flux for Case no 1 is shown in Table 5.
ACCEPTED MANUSCRIPT 64
ICEMOD 2.1
MODEL COMPARISON CASE 1 ICE THICKNESS PROFILES, Z=6.2 m, 50*NSF
MODEL 01
18
MODEL 02
T
16
IP
Ice thickness (cm)
14 12
SC R
10 8 6 4
NU
2 0 -72
-54
-36
-18
0
18
36
54
72
90
108
MA
-108 -90
Angle from the stagnation line (⁰)
TE
D
Fig.35. Ice thickness horizontal profiles at the mean height for three model versions. A very
CASE no 1
CE P
large spray mass flux is assumed (50*NSF).
ICEMOD2.1,
AC
Test number (
Model0l
, Periodic spray) (
M
Model02
, Periodic spray) (
M
, Cont. spray)
M
(kg)
(cm)
(%)
(kg)
(cm)
(%)
(kg)
(cm)
(%)
2743
7.78
37.1
2744
7.79
37.1
2903
7.74
39.5
2 ( 50*NSF) 3755
11.4
1.02
3447
10.6
0.93
295.7
0.753
0.08
1 (NSF)
Table 5. Summary of Case no 1 icing results for ICEMOD2.1 and two simplified models. =”mean” front thickness defined under Table 4,
CASE no 2: D =1.20 m, Z =6.0-6.48 m,
defined by Eq. 52.
, Wind speed:
ACCEPTED MANUSCRIPT 65
At the level Z=6.00 m the instantaneous local water catch rate at the stagnation point during spray impingement is approximately
.The water catch rate has now a small
T
contribution from wind-generated spray (overall collection efficiency 4%). The bulk spray . Both
IP
temperature just before impact at the top of the cylinder is also now about
SC R
ICEMOD2 and ICEMOD2.1 give an overall freezing fraction of 51.4% and there is no difference in total ice load and “mean” windward side ice thickness: 93.0 kg and 10.3 cm.
NU
We will now study the effect of a very large spray mass flux, but this time we multiply the NSF by a factor of 25. The instantaneous local water catch rate at the stagnation point at . Any significant difference
MA
Z=6.0 m during spray impingement then becomes
in ice load calculations for the two two-dimensional models is still not obtained, 129 kg
D
(ICEMOD2.1) and 131 kg (ICEMOD2). The single main difference between the two models
TE
is again the horizontal brine film velocity (Fig.36, some numerical noise occurs). Now the is about 62% lower than the ICEMOD2.1 result.
CE P
ICEMOD2 horizontal velocity at υ
The vertical brine film velocity was however almost equal at all angels for ICEMOD2.1 and ICEMOD2 (not shown here). Comparisons of ice thickness horizontal profiles for the three
AC
models ICEMOD2.1, Model0l and Model02 are shown in Fig.37. ICEMOD2.1 and Model01 give the same profile for
. More brine is now transported to the cylinder lee side
resulting in greater differences in ice loads between the three models compared with Case no 1 (Table 6). The strange Model02 ice profiles (only icing near the cylinder shoulders) is obviously again due to the constant large heat input from the spray cloud. For such cases this model must be discarded.
ACCEPTED MANUSCRIPT 66
MODEL COMPARISON CASE 2 BRINE FILM VELOCITY, Z=6.00 m, 25*NSF
ICEMOD 2.1 ICEMOD 2
400
T
300
IP
200 100 -100 -108 -90 -72 -54 -36 -18
0
SC R
0 18
36
-200 -300
-400
54
72
90 108
NU
Brine film velocity ub (mm/s)
500
-500
MA
Angle from stagnation line (⁰)
Fig. 36. Horizontal profiles of brine film mean lateral velocity at the cylinder base for
TE
D
ICEMOD2.1 and ICEMOD2. A very large spray mass flux is assumed (25*NSF).
ICEMOD 2.1
AC
Ice thickness (cm)
CE P
MODEL COMPARISON CASE 2 ICE THICKNESS PROFILES, Z=6.24 m, 25*NSF
MODEL 01 MODEL 02
25 20 15 10
5 0 -108 -90 -72 -54 -36 -18
0
18
36
54
72
90
108
Angle from stagnation line (⁰)
Fig.37. Ice thickness horizontal profiles for Case no 2 at Z=6.24 m for three model versions. A very large spray mass flux is assumed (25*NSF).
ACCEPTED MANUSCRIPT 67
ICEMOD2.1,
Test number (
Model0l
M
, Periodic spray) (
, Cont. spray)
T
, Periodic spray) (
Model02
M
M
IP
CASE no 2
(cm)
(%)
(kg)
(cm)
(%)
(kg)
(cm)
(%)
92.95
10.3
51.4
78.09
8.74
43.2
88.09
8.80
48.9
2 (25*NSF) 128.6
12.7
2.85
98.69
11.5
2.19
21.47
2.18
0.48
NU
1 ( NSF)
SC R
(kg)
MA
Table 6. Summary of Case no 2 icing results for ICEMOD2.1 and two simplified models.
D
5. CONCLUSIONS AND RECOMMENDATIONS
TE
In the old ice model ICEMOD only one space coordinate was considered. This implied that if
CE P
you wanted to calculate ice accretion along the whole circumference of for example a vertical rig column you had to specify in advance the shape the ice profile based on icing calculation at the stagnation line only. This is no longer necessary in the new two-dimensional models
AC
ICEMOD2 and ICEMOD2.1. For relative wet icing (freezing fraction much lower than 1) the shape of the ice profile is governed by mainly three factors: spray mass flux, heat balance of the brine film layer and movement of the brine film. Typical model ice profile is then like a saddle as shown in Fig. 6 and Fig. 9 for SEDCO 708 and Fig. 12 and Fig. 15 for SEDCO 709. For the vertical columns the model horizontal profiles are as expected symmetric about the stagnation line. For a truss we should expect that for very wet icing most ice accretes at the lower side of the column. This is due to the affect of gravity on the brine film movement. For the SEDCO 708 model result (Fig. 9) we observe that ice thickness is slightly higher on the lower side of the 70⁰ truss. For the SEDCO 709 case (Fig. 15) this tendency is not observed.
ACCEPTED MANUSCRIPT 68
We must however remember that the saddle form of the ice profile is here solely a theoretical result; it is impossible from the photos included in the report of Minsk, 1984, to make such a
T
conclusion. Unfortunately no photos were presented in the first theoretical analysis of the
IP
SEDCO icing cases 27 years ago (Brown and Horjen, 1988). At that time we certainly had
SC R
access to all photo material available. If these photos still exist it would of course be most interesting to make a comparison, if possible, with theoretical results. But most likely, due for example to changes in wind direction during the icing periods, the field ice profiles had not
NU
the nice symmetric form depicted in Fig. 6, 9, 12 and 15.
MA
We observe from Table 1 that model maximum ice thicknesses is up to 29% higher or lower than estimated values while the heights of maximum ice thickness are not very much different
D
from the estimated values from photos. Further we see from Fig. 8 and Fig. 11 that above
TE
Z=4.5-5.0 m maximum model ice thickness is larger than estimated thickness for the column and trusses of SEDCO 708, largest difference being about 6.0 cm for the column and about
CE P
6.4 cm for the truss. From Fig. 14 we see that model maximum ice thickness for the column of SEDCO 709 is larger than estimated thickness when
, the largest difference being
AC
about 6.5 cm. Comparing model and estimated ice thickness at specified heights show relative large differences. Assuming however that model thickness height distribution is displaced 1-2 m downwards many of the model result and the estimated results will more or less coincide; implying that there is some chance that total model and field ice load on a single column or truss will be comparable. This has not been analyzed in this work. For the sensitivity tests a constant significant wave height of 5.0 m is assumed while the wave period is calculated from an empirical relation based on field data from the Grand Bank region (not far from the icing location of SEDCO 709). From Fig. 17 we may conclude that for the base of the cylinder (Z=6.0 m) and an icing period of 24 h the ice thickness increase due to air temperature decrease is about 1.6 cm/
. When the wind speed increases by 1m/s
ACCEPTED MANUSCRIPT 69
the ice thickness increases about 0.5 cm at the same level (Fig. 18). Note however that these conclusions are only valid for the specified test conditions of the sensitivity tests. Finally the
T
effect of varying the spray frequency was investigated. Although not very realistic we have
IP
assumed that the time-average spray mass flux formula used in the icing model is independent
SC R
of spray frequencies (spray frequency was not registered simultaneously during the collection of spray mass flux data on Treasure Scout). The single spray mass flux is then a derived quantity obtained by dividing the time-average value by the product
. Reducing now the
NU
spray frequency to half the value used for rig icing modeling of this work does only have a
MA
minor effect upon ice thickness distribution (Fig.20) and total accreted ice mass (reduced by 4.5 %) for a specified icing period. For continuous spray (implying a stationary model for the
D
icing intensity) the changes are also rather small; total ice load then increases 8.7%. This
TE
conclusion is however not valid for sea spray mass flux 25 times larger than the normal value given by Eq.41 (Fig. 21). Now the total ice load and stagnation line ice thickness
CE P
increases with increasing spray frequency. Based on a theoretical consideration given in the text this is probably due to a large increase in brine film thickness (the example used in this
AC
work gives a factor of about 10). Details of the icing intensity, brine film salinity and brine film thickness variations at and mean cylinder height (Z=7.2 m) are shown in Fig. 23-26 for the initial 10 or 20 spray periods when the air temperature is -8.0
and -12.0 . All these quantities vary periodically
in time. It is interesting to note from Fig. 23 that maximum icing intensity occurs during the “no spray” portion of the on-off spraying cycle. After the first 2-3 spray cycles the brine film salinity seems to stabilize and fluctuate around a certain value. The brine film thickness does also stabilize and fluctuate around some value, but not before a longer time period has elapsed; in our examples 14-15 spray cycles (Fig. 25). As shown in Fig. 26 there is a periodiclinear increase with time of the total accumulated ice mass per unit area.
ACCEPTED MANUSCRIPT 70
The main goal of this article was the presentation of an extended icing model, ICEMOD2.1, based on four differential equations. In the beginning of this work it was expected that the
T
new model would give noticeable different and hopefully better result than the first two-
and Case 2: D = 1.20 m,
) we have
SC R
(Case 1: D = 9.15 m,
IP
dimensional model ICEMOD2. Based on the results of section 4.3 for two rig icing cases
however arrived at the following more or less unexpected conclusions:
NU
1. Concerning total ice load and thickness distribution ICEMOD2 gives the same result as ICEMOD2.1, either for the normal sea spray mass flux (Eq. 41) or a very large sea
or 25 (applied for Case 2).
MA
spray mass flux obtained by multiplying Eq. 41 by a factor of 50 (applied for Case 1)
D
2. Only one calculated parameter turns out to be noticeably different for the two models:
TE
the brine film lateral velocity for a very large sea spray mass flux (see above). Depending on icing conditions ICEMOD2.1 gives either larger (Case 2) or lower
CE P
(Case 1) maximum lateral velocity
than ICEMOD2.
3. For a very large spray mass flux the longitudinal velocity
is largest near the
AC
cylinder shoulders where we have maximum brine film thickness. 4. For NSF an icing model with no brine film movement (called Model0l) gives the same ice load and thickness distribution results as ICEMOD2.1. 5. For NSF a stationary model with no brine film movement (called Model02) gives a satisfactory total ice load approximation to ICEMOD2.1 results. There is however now a small difference in ice thickness distribution; Model02 stagnation line ice thickness being larger than ICEMOD2.1 and Model01 results at the upper section of the cylinder (a simple theoretical proof is given) and smaller at the lower section. 6. For a very large spray mass flux (in this analysis it was of order some brine may move to the lee side of the cylinder due to the action
ACCEPTED MANUSCRIPT 71
of wind stress; especially if the cylinder is slim such as for example a truss. On the lee side icing can be calculated only by using ICEMOD2 or ICEMOD2.1. For large
T
structures (rig columns) Model0l will still approximate total ice load quite well. The
(Case 2).
SC R
(Case 1) and
IP
horizontal thickness distribution of Model0l and ICEMOD2.1 are almost equal for
7. For a very large spray mass flux Model02 is useless since a large continuous mass flux will prevent icing on a major part of the upwind side of the cylinder. Using this model
NU
ice accretion is restricted to a small area near the cylinder shoulders, from about for Case 2,
MA
from the stagnation line for Case 1 and from about and this is obviously not correct.
D
With the latest extension of ICEMOD2 the main goal of making a relatively complete,
TE
theoretical structured model with a significant complement of supporting empirical parameterizations and empirical sub-models is obtained. The use of ICEMOD2.1 instead of
CE P
ICEMOD2 is however mainly of academic interest (brine film lateral velocity). In many cases a sufficient approximation of ice accretion on cylinders can be obtained by just solving two
AC
simple algebraic equations (Eq. 58). But it is after all satisfying that the complete set of icing equations in fact is solvable and gives realistic results. Improving some of the input parameters (for example based on turbulence theory) and considering other geometries may perhaps show that ICEMOD2.1 nevertheless is favorable compared with other more simple models. Further work on icing model should concentrate on studying in more detail some of the input parameters such as local heat transfer coefficient, local single spray mass flux (or alternatively liquid water content), spray frequency, single spray duration, droplet flight time and the interfacial distribution coefficient (is a constant value always a good approximation?). In most icing models heat conduction in the accreted ice close to the brine film/ice interface is
ACCEPTED MANUSCRIPT 72
neglected. A common argument for doing so has been the high liquid porosity of accreted ice near the icing front. But this is not necessarily true and a study of heat transport in the ice and
T
structure is a challenge for future research on spray icing. The problem is complicated for a
IP
time-dependent model due to the fluctuation of temperature at the icing surface.
SC R
A suggestion for future work is also an investigation of the influence of precipitation on sea spray icing. Snow showers were for example reported for several of the observations on the
NU
two SEDCO platforms during the spray icing periods. A combined spray+snow icing model was in fact introduced in the old one-dimensional version of ICEMOD (Horjen, 1990). Based
MA
on ICEMOD2 this model may be extended to a two-dimensional version. A mathematical model can however under no circumstances reproduce the nature completely.
D
A more complicated icing model will not necessarily give more reliable results for field
TE
conditions. Spray impingement on a rig structure is for example more or less a stochastic
CE P
process. The only way to really test the validity of an icing model is to compare model tests
AC
with controlled laboratory experiments.
References
Achenbach, E.G. (1977): “The effect of surface roughness on the cross flow of air”, J. of Heat and Mass Transfer, Vol. 20, pp 359-369 ANSYS ( 2009): ANSYS FLUENT 12.0 Theory guide, ANSYS Inc., Cannonsburg, 816 pp. Assur,A. (1958): “Composition of sea ice and its tensile strength”, in Arctic Sea Ice, U.S. National Academic Sciences – National Research Council Pub. 598, pp. 106-138
ACCEPTED MANUSCRIPT 73
Best, A.C. (1950): “The size distribution of raindrops”, Quart. J. Roy. Met. Soc.,Vol.77, pp.16-36
T
Best, A.C. (1951): “Drop-size distribution of clouds and fog”, Quart. J. Roy. Met. Soc., Vol.
IP
77, pp. 418-426
SC R
Brown, R. D. and Horjen, I. (1989): “Evaluation of state of the art drilling platform icing models”, Canadian Climate Centre, Report No. 89-10, Atmospheric Environment Service,
NU
Ontario, 80 pp.
MA
Brown, R. D., Horjen, I. , Jørgensen, T. and Roebber, P.(1988): “Evaluation of state- of- theart drilling platform ice accretion models”, In 4th Int. Workshop on Atmospheric icing of
D
Structures, Paris, pp. 208-213
TE
Cheney, W. and Kincaid, D. (1985): “Numerical mathematics and computing”, Brooks/Cole
CE P
Publishing Company, Monterey, California. Cox, G.F.N. and Weeks, W.F. (1983): “Equations for determining the gas and brine volumes
AC
of sea-ice samples”, Journal of Glaciology, 29(102), pp. 306-316 Cox, G.F.N. and Weeks, W.F. (1986): “Changes in the salinity and porosity of sea-ice samples during shipping and storage”, Journal of Glaciology, 32(112), pp. 371-375 Forest,T., Lozowski, E. and Gagnon, R. (2005): “Estimating marine icing on offshore structures using RIGICE04”, International Workshop of Atmospheric Icing on Structures (IWAIS), Montreal. Horjen, I. and Vefsnmo, S. (1985): “A kinematic and thermodynamic analysis of sea spray (in Norwegian), Offshore Icing – Phase II, Norwegian Hydrotechnical Laboratory (NHL). Report STF60 F85014
ACCEPTED MANUSCRIPT 74
Horjen, I. and Vefsnmo, S. (1986): “Calibration of ICEMOD - Extension to a time-dependent model”, Offshore Icing – Phase I/II, NHL- report STF60 F86040, 73 pp.
T
Horjen, I. and Vefsnmo, S. (1987): “Time-dependent sea spray icing on ships and drilling rigs
IP
– a theoretical analysis”, Offshore Icing – Phase IV, NHL- report STF60 F87130, 87 pp.
SC R
Horjen, I. (1990): “Numerical modeling of time-dependent marine icing, anti-icing and deicing”, doctor theses, Trondheim, Norway
NU
Horjen, I. (2013): “Numerical modeling of two-dimensional sea spray icing on vessel-
MA
mounted cylinders”, Cold Regions Science and Technology, 93, p. 20-35 Incropera, F. and De Witt, D.P. (1985): “Fundamentals of heat and mass transfer”, second ed.,
D
Wiley, New York
TE
Jones, K.F. and Andreas, E.L. (2012): “Sea spray concentration and the icing of fixed
144
CE P
offshore structures”, Quarterly Journal of the Royal Meteorological Society, 138(662):131-
AC
Jørgensen, T.S. (1984): “Sea spray measurements on the drilling rig “Treasure Scout” .Results from tests April 1983-February 1984” (in Norwegian), Offshore Icing-Phase II, OTTERreport STF88 F84038
Jørgensen, T.S. (1985): “Sea spray characteristics on a semi-submersible drilling rig”, Offshore Icing-Phase II, NHL-report STF60 F85015 Kays W.M. and Crawford M.E.(1980): “Convective heat and mass transfer”, McGraw-Hill, New York, 420 pp.
ACCEPTED MANUSCRIPT 75
Kulyakhtin, A., Kulyakhtin, S. and Løset, S. (2013):”Measurements of thermodynamic properties of ice created by frozen spray”, Proc. 23. International Offshore and Polar
T
Engineering, Anchorage, Alaska, p.1104-1111
IP
Kulyakhtin, A. and Tsarau, A. (2014):”A time-dependent model of marine icing with
SC R
application of computational fluid dynamics Measurements”, Cold Regions Science and Technology, 104-105, p. 33-44
NU
Kulyakhtin, A. (2014):”Numerical Modelling and Experiments on Sea Spray Icing”, doctor
MA
theses, Trondheim, Norway
Lai, R. J. and Shemdin, O.H. (1974): “Laboratory study of the generation of spray over
D
water”, J. Geophysical Res., 79(21), pp. 3055-3063
TE
Lax P.D. and Wendroff, B.(1964): “Difference schemes for hyperbolic equations with high
CE P
order of accuracy”, Comm. Pure Appl. Math., Vol. 17, pp. 381-398 Lozowski, E.P., Stallabrass, J.P. and Hearty, P.F. (1979): “The icing of an unheated nonrotating cylinder in liquid water droplet-ice crystal clouds”, Mech. Eng. Report LTR-LT-96,
AC
National Research Council, Canada, 109 pp. Makkonen, L. (1985): “Heat transfer and icing of a rough cylinder”, Cold Regions Science and Technology, vol. 10, pp.105-116 MacLaren Plansearch (1986): “Final report of the offshore icing measurement project, Scotian Shelf. Contractors final report prepared for Atmospheric Environment Service”, Bedford, Nova Scotia.
ACCEPTED MANUSCRIPT 76
Myers, T.G. and Charpin, J.P. (2004): A mathematical model for atmospheric ice accretion and water flow on a cold surface”, International Journal of Heat and Mass Transfer, 47(25):
T
5483-5500
IP
Minsk, L.D. (1984): “Ice observation program on the semisubmersible drilling vessel SEDCO
SC R
708”, CRREL Special Report 84-2, 14pp
Mitten, P. (1994): “Measurement and modeling of spray icing on offshore structures”, Final
NU
report, Atmospheric Environmental Service of Canada, Contract no. 07SE.KM169-8-7439
MA
Monaham, E.C. (1968): “Sea spray as a function of low elevation wind speed”, J. Geophysical Res., 73(4), pp. 1127-1137
D
Preobrazhenskii, L.Y. (1973): “Estimate of the content of spray drops in the near-water layer
TE
of the atmosphere”, Fluid Mechanics – Soviet Research, 2 (2): pp. 95-100
CE P
Richtmyer, R.D. (1963): “A survey of difference methods for non-steady fluid dynamics, N.C.A.R. Tech. notes 63-2
AC
Ryerson, C.R. and Gow, A.J. (2000):”Ship superstructure icing: Crystalline and physical properties”, ERDC/CRREL (Cold Regions Research and Engineering Laboratory) report TR00-11
Schlichting, H. (1955): ”Boundary layer theory”, McGraw-Hill, First English edition, 535 pp. Sobey, I. J. (2000): “Introduction to interactive boundary layer theory”, Oxford University Press, ISBN 0-19-850675-9. Stallabrass, R. (1980): “Trawler icing. A compilation of work done at N.R.C.”, Mech. Eng. Report MD-56, No 19327, Ottawa, Canada
ACCEPTED MANUSCRIPT 77
Weeks, W. and Lofgren, G. (1967): “The effective solute distribution coefficient during the freezing of NACL solutions”, In: Physics of snow and Ice, Institute of Low Temperature
T
Science, Hokkaido University, Sapporo, Japan, pp.579-597.
IP
Weeks, W. (2010): “On sea ice”, University of Alaska Press, Fairbanks, AK, 664 pp.
SC R
Wu, J., Murray, J.J. and Lai, R.J. (1984): “Production and distribution of sea spray”, Journal of Geophysical Res., 89 (C5), p. 8163-8169
NU
Zakrzewski, P. (1987): “Splashing a ship with collision-generated spray”, Cold Regions
MA
Science and Technology, 14, p. 65-83
TE
D
APPENDIX:
icing events
CE P
Weather parameters used as input to the SEDCO 708 and SEDCO 709
AC
The data in each line below is from left to right: Date (year, month, day, time), air temperature ( ), dew-point temperature ( ), air pressure (mb), sea surface temperature ( ) , wind speed (knot), significant wave height (m), significant wave period (s) and sea surface salinity (‰). For SEDCO 708 the time is LST while for SEDCO 709 the time is GMT. SEDCO 708 ICING EVENT (60 observations) 83010404, -1.7,
-3.3,
999.9,
3.8,
18.8,
1.5,
5.7, 32.5
83010405, -2.2,
-3.3,
999.9,
3.8,
19.1,
1.5,
5.2, 32.5
83010406, -2.2,
-3.3, 1000.5,
3.8,
23.4,
1.7,
5.2, 32.5
ACCEPTED MANUSCRIPT 78
-3.3, 1000.7,
3.8,
20.6,
1.6,
5.6, 32.5
83010408, -2.8,
-3.9, 1000.9,
3.8,
22.4,
1.5,
5.6, 32.5
83010409, -2.8,
-3.9, 1001.1,
3.8,
21.7,
1.6,
5.4, 32.5
83010410, -2.8,
-3.9, 1001.3,
3.8,
20.5,
83010411, -2.8,
-3.9, 1001.6,
3.8,
20.5,
83010412, -2.8,
-3.9, 1001.5,
3.8,
83010413, -3.9,
-5.0, 1001.1,
83010414, -3.3,
-4.4, 1000.2,
83010415, -2.8,
-4.4, 1000.4,
83010416, -1.7,
-3.3, 1000.4,
83010417, -2.2,
T
83010407, -2.2,
5.8, 32.5
1.6,
5.9, 32.5
21.4,
1.7,
6.1, 32.5
3.8,
17.3,
1.7,
6.1, 32.5
3.8,
17.7,
1.6,
5.8, 32.5
3.8,
14.9,
1.5,
5.9, 32.5
3.8,
19.3,
1.5,
6.1, 32.5
-3.3, 1000.1,
3.8,
20.6,
1.5,
5.7, 32.5
83010418, -2.3,
-3.4, 1000.4,
3.8,
18.5,
1.5,
5.5, 32.5
83010419, -2.3,
-3.4, 1000.4,
3.8,
16.2,
1.3,
5.4, 32.5
83010420, -2.4,
-3.5, 1000.4,
3.8,
15.9,
1.3,
5.3, 32.5
IP
SC R
NU
MA
D
TE
83010421, -2.4,
-3.5, 1000.1,
3.7,
14.7,
1.4,
5.5, 32.5
83010422, -2.5,
-3.6, 1000.2,
3.7,
16.9,
1.3,
5.3, 32.5
AC
CE P
1.6,
83010423, -2.5,
-3.6,
999.9,
3.7,
17.5,
1.3,
5.4, 32.5
83010500, -2.6,
-3.7,
999.6,
3.7,
20.3,
1.3,
5.2, 32.5
83010501, -2.6,
-3.7,
998.5,
3.8,
25.6,
1.5,
5.1, 32.5
83010502, -2.7,
-3.8,
998.0,
3.8,
30.0,
1.8,
5.5, 32.5
83010503, -2.7,
-3.8,
997.4,
3.8,
32.2,
2.2,
5.9, 32.5
83010504, -2.8,
-3.9,
997.6,
3.8,
31.6,
2.5,
6.1, 32.5
83010505, -2.8,
-3.9,
997.3,
3.8,
32.0,
2.3,
6.0, 32.5
83010506, -2.8,
-3.9,
996.7,
3.8,
35.2,
2.6,
6.4, 32.5
83010507, -2.2,
-3.3,
996.3,
3.7,
37.1,
2.8,
6.4, 32.5
83010508, -2.2,
-3.3,
996.3,
3.7,
38.9,
3.0,
6.6, 32.5
ACCEPTED MANUSCRIPT 79
-3.3,
996.4,
3.7,
39.6,
3.1,
7.2, 32.5
83010510, -2.2,
-3.3,
996.7,
3.7,
38.2,
3.5,
7.5, 32.5
83010511, -2.7,
-3.3,
996.9,
3.7,
38.8,
3.9,
8.0, 32.5
83010512, -2.8,
-3.4,
996.6,
3.7,
38.3,
83010513, -2.8,
-3.5,
996.0,
3.7,
37.9,
83010514, -2.8,
-3.6,
995.6,
3.7,
83010515, -2.8,
-3.7,
995.6,
83010516, -2.8,
-3.8,
995.8,
83010517, -2.8,
-3.9,
996.1,
83010518, -2.9,
-4.0,
996.5,
83010519, -2.9,
-4.1,
83010520, -3.0,
-4.2,
83010521, -3.1,
-4.3,
83010522, -3.2,
-4.4,
T
83010509, -2.2,
7.9, 32.5
4.2,
8.4, 32.5
39.3,
3.8,
7.9, 32.5
3.7,
35.2,
3.8,
8.1, 32.5
3.7,
33.4,
4.4,
8.4, 32.5
3.7,
34.8,
4.0,
8.0, 32.5
3.7,
33.7,
4.2,
7.9, 32.5
996.8,
3.7,
32.4,
4.4,
8.2, 32.5
997.3,
3.7,
31.7,
3.8,
8.0, 32.5
997.7,
3.7,
31.9,
3.8,
7.7, 32.5
998.0,
3.6,
30.3,
4.3,
8.0, 32.5
IP
SC R
NU
MA
D
TE
83010523, -3.3,
-4.4,
998.2,
3.6,
27.5,
4.1,
8.3, 32.5
83010600, -3.3,
-4.4,
998.6,
3.6,
29.1,
4.4,
8.6, 32.5
AC
CE P
3.5,
83010601, -3.3,
-4.4,
998.8,
3.6,
29.8,
4.1,
8.5, 32.5
83010602, -4.4,
-5.5,
998.9,
3.6,
30.4,
3.6,
8.2, 32.3
83010603, -5.6,
-6.7,
999.2,
3.7,
36.3,
3.8,
8.3, 32.3
83010604, -6.7,
-7.8,
999.2,
3.7,
38.1,
4.8,
9.0, 32.5
83010605, -6.7,
-7.8,
999.1,
3.7,
39.3,
5.2,
9.0, 32.5
83010606, -7.2,
-8.3,
999.1,
3.7,
38.8,
5.4,
9.3, 32.5
83010607, -8.3,
-9.4,
999.2,
3.7,
40.0,
5.3,
9.0, 32.5
83010608, -8.3,
-9.4,
999.2,
3.7,
39.4,
5.0,
8.8, 32.5
83010609, -8.3,
-9.4,
999.2,
3.7,
38.7,
4.4,
8.5, 32.5
83010610, -7.8,
-8.9,
999.0,
3.7,
38.6,
4.8,
8.9, 32.5
ACCEPTED MANUSCRIPT 80
-9.4,
998.8,
3.7,
38.7,
4.9,
8.7, 32.5
83010612, -8.3,
-9.4,
998.4,
3.7,
40.3,
5.2,
8.8, 32.5
83010613, -8.3,
-9.4,
997.9,
3.7,
39.1,
4.9,
9.2, 32.5
83010614, -8.3,
-9.4,
997.0,
3.7,
39.2,
83010615, -8.3,
-9.4,
996.3,
3.7,
39.1,
86012415, -5.0,
-8.0, 1028.0,
86012416, -5.0,
-8.0, 1027.5,
4.9,
9.0, 32.5
IP
SC R
-8.0, 1027.4,
1.9,
21.0,
1.4,
4.9, 34.0
1.8,
24.0,
1.4,
4.5, 34.0
1.8,
26.0,
1.5,
4.4, 34.0
86012417, -7.0, -10.0, 1027.4,
1.7,
27.0,
1.7,
4.4, 34.0
86012418, -7.0, -11.0, 1027.8,
1.7,
29.0,
1.7,
5.9, 34.0
86012419, -7.0, -11.0, 1028.2,
2.1,
28.0,
1.6,
4.4, 34.0
86012420, -7.0, -11.0, 1029.5,
1.9,
32.0,
1.8,
4.5, 34.0
86012421, -8.0, -11.0, 1030.9,
2.0,
35.0,
1.8,
4.4, 34.0
86012422, -8.0, -11.0, 1032.0,
2.0,
37.0,
2.1,
4.7, 34.0
86012423, -9.0, -11.0, 1033.6,
1.8,
34.0,
2.0,
4.5, 34.0
86012500, -9.0, -12.0, 1034.1,
1.6,
36.0,
2.0,
4.5, 34.0
86012501, -9.0, -12.0, 1034.9,
1.7,
32.0,
2.0,
4.5, 34.0
86012502,-10.0, -15.0, 1035.8,
1.8,
30.0,
2.1,
4.7, 34.0
86012503,-10.0, -15.0, 1035.8,
1.8,
30.0,
1.8,
4.5, 34.0
86012504,-10.0, -15.0, 1036.9,
2.1,
29.0,
1.8,
4.6, 34.0
86012505,-10.0, -14.0, 1037.3,
1.4,
27.0,
1.8,
4.8, 34.0
86012506,-10.0, -14.0, 1039.3,
0.8,
25.0,
1.7,
4.7, 34.0
86012507,-10.0, -14.0, 1040.4,
1.0,
21.0,
1.7,
4.7, 34.0
AC
CE P
TE
D
MA
86012414, -5.0,
8.9, 32.5
4.6,
NU
SEDCO 709 ICING EVENT (25 observations)
T
83010611, -8.3,
ACCEPTED MANUSCRIPT 81
1.4,
23.0,
1.8,
4.8, 34.0
86012509,-10.0, -16.0, 1041.8,
0.9,
21.0,
1.8,
4.9, 34.0
86012510,-10.0, -16.0, 1042.2,
0.5,
20.0,
1.6,
4.7, 34.0
86012511, -9.0, -14.0, 1043.1,
0.2,
24.0,
86012512, -9.0, -14.0, 1043.8,
1.1,
24.0,
86012513, -9.0, -14.0, 1044.9,
1.0,
86012516, -9.0, -14.0, 1044.8,
1.2,
T
86012508,-10.0, -14.0, 1041.2,
4.5, 34.0
1.4,
4.5, 34.0
22.0,
1.4,
4.5, 34.0
19.0,
1.1,
4.1, 34.0
IP
SC R
NU
1.6,
AC
CE P
TE
D
MA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ACCEPTED MANUSCRIPT 82
CE P
TE
D
MA
NU
SC R
IP
T
An extension of the two-dim. sea spray icing model ICEMOD2 is presented. Model results are compared with field results from two semi-submersible oil rigs. Sensitivity analysis of the effect of varying three input parameters is presented. The effect of including brine film motion in icing modeling has been discussed. The main difference between the new and old ICEMOD2 is the brine film velocity.
AC