Analysis of a model for the formation of fold-type oscillation marks in the continuous casting of steel

This paper investigates the different possible behaviours of a recent asymptotic model for oscillation-mark formation in the continuous casting of steel, with particular focus on how the results obtained vary when the heat transfer coefficient (m), the thermal resistance (Rmf ) and the dependence of the viscosity of the flux powder as a function of temperature, μf (T) , are changed. It turns out that three different outcomes are possible: (I) the flux remains in molten state and no solid flux ever forms; (II) both molten and solid flux are present, and the profile of the oscillation mark is continuous with respect to the space variable in the casting direction; (III) both molten and solid flux are present, and the profile of the oscillation mark is discontinuous with respect to the space variable in the casting direction. Although (I) gave good agreement with experimental data, it suffered the drawback that solid flux is typically observed during actual continuous casting; this has been rectified in this work via alternative (II). On the other hand, alternative (III) can occur as a result of hysteresis-type phenomenon that is encountered in other flows that involve temperature-dependent viscosity; in the present case, this manifests itself via the possibility of multiple states for the oscillation-mark profile at the instants in time when solid flux begins to form and when it ceases to form.


Introduction
Oscillation marks are more or less evenly spaced indentations along the surface of steel slabs produced via continuous casting. In this process, the solidification of the molten metal is initiated by a cooling mould that oscillates in the axis of the casting direction; furthermore, to prevent the solid shell from sticking to the mould, mould powder, often also termed flux or slag, is added at the top of the mould, melts and forms a lubricating layer in the gap between the steel and the mould walls, as well as a meniscus with the molten steel. A schematic of the situation is shown in Fig. 1(a).  The problem is then divided into a lower zone, where z > z 0 (t) , and an upper zone where z < z 0 (t) . In Section 2.1.1 below, we give the governing equations for the lower zone; after that, in Section 2.1.2, we give the equations for the upper zone.
2.1.1 Lower zone. For 0 < x < s (z, t) , heat transfer is governed by where T is the temperature and c (s) f , k (s) f and ρ (s) f are, respectively, the specific heat capacity, the thermal conductivity and the density of the solid flux; typically, the mould oscillates according to V (t) = V 0 cos ωt, (2.2) and the solid flux is thus assumed to oscillate with it. In writing (2.1), we have already assumed that the geometry we consider is slender, and we have therefore neglected heat conduction in the z-direction; we will do likewise for the remaining regions. For s (z, t) < x < s (z, t) + h (z, t) , conservation of mass, momentum in the z-direction and heat are expressed, respectively, by f and μ f are, respectively, the specific heat capacity, the thermal conductivity, the density and viscosity of the molten flux, u x and u z are the molten flux velocities in the x and z-directions respectively, p is the pressure and g is the gravitational acceleration. For (2.4), we have employed the usual lubrication approximation, whereby the inertia terms are neglected, which leads to ∂p/∂z being a function of z and t; note that the use of the lubrication approximation is consistent with the simplification made in writing (2.1) and (2.5) above and (2.6) below. We remark also that we will take μ f to be temperature-dependent, in contrast to Vynnycky et al. (2017), who assumed it to be constant; the exact nature of this temperature-dependence is discussed in detail in Section 2.3. For , the heat equation is written as where c (s) s , k (s) s and ρ (s) s are, respectively, the specific heat capacity, the thermal conductivity and the density of the solid steel.
Lastly, we point out that we do not solve any equation for the molten steel region but will simply include its effect in terms of a heat flux prescribed at x = s (z, t) + h (z, t) + f (z, t) .

Boundary and initial conditions
A heat balance at x = 0 gives , T m,f < T < T m,s and s = 0, h = 0 , (2.8) where m is the heat transfer coefficient linking the temperature at x = 0 and T w , R mf is the interface thermal contact resistance between the solid flux and the mould and R ms is the interface thermal contact  Vynnycky et al. (2017); here, however, we do not make this assumption. Observe also that the second alternative in (2.8) suggests that molten flux would not give rise to a thermal contact resistance, which may in itself be an oversimplified description, in view of the rheology of the flux at these temperatures. We nevertheless keep this second alternative as it stands in (2.8), and as was presented in Vynnycky et al. (2017), before re-evaluating it in Section 5.1. Note that, after T has been determined, the temperature at the outer surface of the mould, T mould , and the heat flux at this surface, q, can be found via (2.11) where ( ) ± denotes the value of a function in the limit as x tends to s (z, t) from above and below, respectively, and ΔH f is the latent heat of fusion for the flux. Physically, equations (2.12)-(2.15) represent: respectively, the continuity of the z-component of velocity, so that the molten flux moves with the speed of the solid flux, which is in turn assumed to move with the speed of the mould wall, i.e. no slip; the continuity of the x-component of the velocity; the temperature being equal to the melting temperature of the flux; conservation of heat, relating the differences in heat flux to the latent heat released due to phase change. Since the effect of latent heat release at this interface is believed to be small (Hill et al., 1999), it is seldom, if ever, included in others' models for oscillation-mark formation; indeed, we were not even able to find a numerical value for ΔH f in the oscillation-mark literature. On the other hand, we have been able to find values in the context of electroslag refining (Kelkar et al., 2016;Kharicha et al., 2008Kharicha et al., , 2016Yanke et al., 2016) and we will thence proceed by retaining the right-hand side of equation ( The physical interpretations of ( where ΔH s is the latent heat of fusion for steel; note that using (2.16) and (2.17) reduces (2.22) to (2.23) The physical meaning of (2.21) and (2.22) can be deduced from the previous discussion and is therefore not repeated. However, we note that there is no governing equation to determine (∂T/∂x) + in equation (2.23), and we will therefore set, for the time being, and consider the effect of this term later. For the pressure, we have (2.26) The first of these is based on the idea that since molten flux and steel are in contact at z = z 0 (t) , their pressures are equal (Bland, 1984); the right-hand side of (2.25) is the molten steel metallostatic pressure, with p a as the atmospheric pressure. On the other hand, equation (2.26) represents the fact that at some distance L down the caster, the pressure will once again be atmospheric, either because an air gap will form between the solidified flux and the mould, or at the bottom of the mould itself, where the system is exposed to the ambient atmosphere. For simplicity, we will take L to be the length of the mould, since we can then prescribe it; otherwise, its value would need to be determined, although it would still be expected to be of the order of magnitude of the length of the mould (Bland, 1984). While (2.8)-(2.26) can be considered as boundary conditions for the lower zone, boundary conditions are also required for the upper zone and in particular for equation (2.7). For these, we have p = p a at z = 0, (2.27) (2.29) Respectively, these express the following: the meniscus is at atmospheric pressure; the continuity of the z-component of velocity, as in equation (2.12); zero shear stress at the interface of molten flux and molten steel, which arises from assuming continuity of shear stress at this interface, combined with the fact that the flux is much more viscous than the steel, giving the proposed simplification. Note also that, in the original development (Hill et al., 1999;Vynnycky et al., 2017), the momentum equations are solved for both upper and lower zones, whereas the heat equations are solved only in the lower zone; consequently, a boundary condition would be required at z = z 0 (t) for the temperature and, in Hill et al. (1999), this is calculated to be the profile obtained by assuming a conductive temperature profile in the flux and steel layers.
The problem also formally requires initial conditions, but since we will be considering the periodic behaviour of the system that is eventually established, and since our treatment here is analytical rather than numerical, there is no need to specify these explicitly here.

Model parameters
Most of the model parameters are reasonably well known and are given in Table 1, apart from m, R mf , R ms and μ f (T) . We give the background to these in turn.
For m, literature values vary between around 1×10 4 and 4.5×10 4 Wm −2 K −1 (Hill et al., 1999;Jonayat & Thomas, 2014) and it is a difficult quantity to specify for a number of reasons. It is a reflection of the heat transfer across the mould and a thermal boundary layer, possibly turbulent, at the mould surface in the water channels. Moreover, the channels are the space between fins that are aligned along the z-direction; the structure can be as in Fig. 7(b) of Meng & Thomas (2003). Consequently, m is unlikely to be constant, although we follow common practice and take it to be so.
The value of R mf is difficult to gauge, since solid mould slags can exist as glassy or crystalline phases or as mixtures of the two, i.e. slag films, and the properties for the various phases can vary considerably. Since the crystalline phase has a higher density than the glassy phase within the film, crystallization is accompanied by shrinkage, which produces rugosity and loss of contact with the mould, creating a thermal resistance, as noted by Holzhauser et al. (1999). However, values for R mf as high as 6×10 −4 m 2 K −1 W −1 have been suggested (Ramirez-Lopez et al., 2010).   Vynnycky et al. (2017). However, the values for ΔH f are taken from Kelkar et al. (2016); Kharicha et al. (2008Kharicha et al. ( , 2016 and Yanke et al. (2016) Symbol Value Unit A suitable value for R ms is similarly difficult to prescribe; Vynnycky et al. (2017) simply set it equal to R mf . In fact, as we show later, its value is insignificant as regards to determining the oscillation-mark profiles, and we therefore leave its value undetermined.
Numerous models have been developed to estimate molten slag viscosity based on its composition and temperature during cooling, based mainly on Arrhenius or Weymann relations (Iida et al., 2000;Mills & Sridhar, 1999;Nagano & Koyama, 1987;Riboud et al., 1981); the vast area of casting powders has been covered recently in Mills & Däcker (2017). A widely used model, by Riboud et al. (1981) based on 45 slags, gives the slag viscosity as where T is temperature in Kelvin and A,B > 0 are parameters that depend on the composition of the slag; this can typically be a vast array of oxides, although predominantly SiO 2 and CaO (Jonayat & Thomas, 2014). However, a limitation of the Riboud model is that it does not predict the abrupt increase in viscosity observed at some temperature during cooling, termed as the break temperature. An alternative Table 2 Model parameters for viscosity profiles is then to use a power-law relation of the form where n > 0 and T fsol are chosen empirically to fit measured data and μ 0 is the viscosity measured at a reference temperature T 0 (Jonayat & Thomas, 2014;Meng & Thomas, 2003). In what follows, we will use these two expressions, although in the case of the first one, we approximate it by the form which is commonly referred to as the Reynolds law (Bland, 1984;Hill et al., 1999). The replacement of exp (B/T) by exp(−B * T) is clearly suggestive of the Frank-Kamenetskii approximation (King, 2009;Vynnycky & O'Brien, 2012), and although (2.30) contains T as a pre-multiplicative factor, it turns out that (2.32) fits (2.30) very well in a least-squares sense for all 45 slags in Riboud et al. (1981); further details of this are given in appendix A. Moreover, the significance of this is that, although (2.31) and (2.32) are both nonlinear functions of T, they will still permit quasi-analytical progress with regard to determining the oscillation-mark profiles. To fix ideas still further, we will consider in particular (2.31) and the two most extreme cases from Riboud et al. (1981), termed as 'high' and 'low' and written in the form (2.32); the relevant constants are given in Table 2, and the three profiles for μ f (T) are plotted in Fig. 2.

Nondimensionalization
We nondimensionalize the model equations, (2.1)-(2.8) and (2.12 )-(2.29), with a view to identifying the key dimensionless parameters and being able to propose a reduced model later on. For this purpose, we set , and ΔT are, respectively, x-length, z-length, time and temperature difference scales that are taken to be also, [Q] is a heat flux scale that we assume to be known. The numerical values for all model parameters are given in Table 1; these consist of the data used in Vynnycky et al. (2017), Kharicha et al. (2008Kharicha et al. ( , 2016, Kelkar et al. (2016) and Yanke et al. (2016). For later use, we note that the data lead to these values serve for the purposes of later comparison with the dimensions of the experimentally measured and theoretically calculated oscillation marks. Apart from the variable viscosity, the nondimensionalized equations are, for the most part, as given in Section 4 in Vynnycky et al. (2017) and are therefore not repeated here; the only differences arise in distinguishing between the thermal contact due to solid flux and solid steel in equations (4.8) and (4.10), respectively, and in the retention of the latent heat term in (4.14) in Vynnycky et al. (2017), which will appear in equation (4.10) below. Instead, we consider the values of the model's nondimensional parameters.

Nondimensional parameters
The model equations contain 15 nondimensional parameters: , Their typical values, based on the values of the dimensional parameters given in Table 1, are presented in Table 3.

Reduced model governing equations
We now proceed by proposing a reduced model by making use of the sizes of the dimensionless parameters to simplify the governing equations, as appropriate. In view of the values given in Table 3, we observe that In addition, although values for Bi mf and Bi ms are not given in Table 3, as we will be varying m, R mf and R ms , we will take them as being O (1) quantities, as was the case in Vynnycky et al. (2017). Thus, we arrive at the following reduced model equations.

Boundary conditions
At X = 0, (2.8) becomes Note that, in (4.10), we have now included the latent heat term for completeness. However, Table 3 indicates St f > 1, and possibly even St f 1; for this reason, we will drop this term in what follows in Section 5, although its inclusion would certainly make for interesting future work. However, a number of physical and mathematical comments can be made here. Physically, this term represents a release of thermal energy, meaning that it can contribute to the re-melting of the solid flux. On the other hand, mathematically, even though it is possible that (4.14) (4.16)

Analysis
We proceed in Sections 5.1 and 5.2 by treating the heat and momentum equations, respectively. The analytical forms found for the temperature and velocity field are then used in Section 6 to constitute solutions for the actual oscillation-mark profile.

Heat
For heat transfer, the analysis is similar to that in Vynnycky et al. (2017), although with some deviation because R mf = R ms . The key results are as follows: 1. Because K fs 1, it is found that in the solid steel. This was based on the fact that St s 1 and that if the superheat of the molten steel, i.e. the difference between its input and melting temperatures, is sufficiently small, then the second term on the left-hand side of (4.16) will be small since K 1; thence, (4.5) and (4.14) will lead to (5.1). Although K was not quantified in Vynnycky et al. (2017), we do so here. A preliminary estimate for [Q] would be, on considering a conductive heat flux scale, where T in is the input temperature and l is a characteristic distance from the nozzle that delivers molten steel, which would be located on the right of Figs 1(a) and 1(b). With l ∼ 0.1 m and T in − T m,s ∼ 30 K, we obtain [Q] ∼ 1.4 × 10 4 Wm −2 , giving K ∼ 10 −4 . However, this estimate would certainly be too low, since the flow of the molten steel is turbulent, and the usual way to account for this is to use an effective value for the molten steel thermal conductivity that is several times that of the value of k (l) s ; e.g. Lait et al. (1974) and Hill et al. (1999) take the factor to be seven. Even so, we find that K 1, meaning that the contribution in (4.16) from the molten steel can be safely neglected.
2. To maintain analytical tractability, it was assumed that α f 1, so that the left-hand side of (4.4) can be neglected; the left-hand side of (4.1) was also neglected for reasons of tractability, as was the heat flux from the molten steel. These assumptions, which imply that heat transfer is dominated by conduction as a consequence of the slenderness of the flux layer, were the ones adopted in the earlier model by Hill et al. (1999), from which the current model derives. While equations (4.1) and (4.4) suggest that this might not be the case, it is certainly plausible that neglecting these terms does not change the qualitative behaviour of the solution, e.g. as would happen if α f 1. Moreover, it is well established that experience suggests that asymptotic solutions are useful numerically far beyond their nominal range of validity, as discussed by Crighton (1994) and Andrianov & Awrejcewicz (2006). Indeed, a recent example of this was seen in the context of the oscillation-mark problem itself in the two-dimensional computations for the flow of flux above the solidification point, wherein a time-derivative term multiplied by an O (1) constant was dropped from the Navier Stokes equations in order to expedite the numerical solution; as seen in Fig. 7 in Vynnycky & Zambrano (2018), the consequences were negligible. Lastly on this point, the fact that the resulting solutions are found to compare favourably against experimental data, as discussed in Section 6.3, gives a stronger basis of validity for neglecting the terms on the left-hand sides of (4.1) and (4.4) than if there had been no data.
3. Three situations can now arise: • solid and molten flux are present, in which case i.e. the flux temperature at X = 0 is greater than or equal to the flux melting temperature. Interestingly, however, requiring S ≥ 0 in equation (5.5) does not to lead to (5.7) unless Bi mf , then either the temperature at the flux phase-change interface is greater than the melting temperature or it is less than the melting temperature, i.e. there is a discontinuity in θ . If it is greater, we should then in fact not even have solid flux at all; if it is less, then we would get a jump in the oscillation-mark profile. In fact, we will find later on that there is also another mechanism in this reduced model for jumps in the oscillation-mark profile, and this issue will have to be considered in detail anyway. However, to simplify matters here, we will simply set Bi where F L 1 is a function to be determined. Noting that, in view of (5.3) and (5.6),μ f =μ f (X, Z, τ ) , (5.11) can be integrated with respect to X to give with Π := ∂P/∂Z and ; (5.13) From (2.12) and (2.16), we have that F L 1 and F L 2 satisfy and Π (Z, τ ) f 1 (S + H, Z, τ ) + f 2 (S + H, Z, τ ) F L 1 (Z, τ ) + F L 2 (Z, τ ) = 1, where (Z, τ ) is suppressed on S and H when they appear in the limits of the integrals or as arguments of f 1 and f 2 , whence In order to be able to connect the upper and lower zones, we will need to consider the mass flux of molten flux across Z = Z 0 (τ ) ; denoting the dimensionless counterpart of this for the lower zone by Q L l (Z, τ ), we will have , (5.14) On using the nondimensional form of (2.20) and integrating with respect to Z, we have where Q * R (τ ) can be determined from boundary conditions (2.25) and (2.26); here, Q * R (τ ) is the timedependent function of integration that arises from the integration with respect to Z of the dimensionless version of (2.20). We find Note that Z 0 (τ ) ≤ 1 and Z L 1, and that Γ ∼ Λ, indicating that equation (5.18) can be reduced to (5.19) i.e. the unknown Z 0 (τ ) is removed from this part of the model, as it has negligible impact on the integrals in (5.18), since Z L Z 0 (τ ). The implication is that, for this fold-type oscillation-mark model, there is no need to determine Z 0 (τ ) , as the interface between the molten flux and molten steel above the solidification point has been assumed to be planar, as shown in Fig. 1(b). On the other hand, in a model that takes into account a curved molten flux-molten steel interface, as shown in Fig. 1(a) would be necessary to determine Z 0 (τ ) ; this is discussed at length in Vynnycky & Zambrano (2018), and recapped shortly in section 5.2.2. Finally, we arrive at (5.20)

Upper zone.
In the upper zone, we have (5.29) where Q U l (Z, τ ) denotes the dimensionless mass flux of flux for the upper zone. Now, since we expect the mass flux of flux to be continuous on passing from the upper to the lower zone, we must have that Also, in view of (4.12), we must have Note that this condition is independent of Z 0 (τ ) ; as mentioned in Section 5.2.1, it is not needed for this fold-type oscillation-mark model. However, in a study that considers curved molten flux-molten steel interface, Vynnycky & Zambrano (2018) suggest that it would be determined by requiring the total heat flux to be continuous across Z = Z 0 (τ ) , by analogy with (5.30).
Recalling that we see that (5.32) can be formulated just in terms of H 0 . By analogy with Vynnycky et al. (2017), it becomes clear how (5.32) will simplify even when μ f is not constant. The key point is whether or notF L 3 (Z, τ ) ever vanishes. We see that when V (τ ) = 1, which clearly occurs when τ = τ * , where we will haveF L 1 = 0 andF L 2 = 1 from (5.14) and (5.15), respectively, leading to (5.34) where S * 0 = S 0 (τ * ) and H * 0 = H 0 (τ * ) . It turns out that the only possibility is that S * 0 = 0, H * 0 = 0; this is explained in detail in appendix B. In turn, this implies thatF L 3 (Z, τ ) vanishes, and which means that (5.32) can be simplified to just (5.36)

Results
To proceed further, we focus on whether m, R mf , R ms and μ f (T) can be chosen so that the model can reproduce the oscillation-mark profile that was obtained experimentally in Vynnycky et al. (2017); the experimental data and the method used to obtain it are documented in more detail in Section 2 of Vynnycky et al. (2017) and in the thesis by Saleem (2016). We recall that good agreement was achieved in Vynnycky et al. (2017), but with the drawback that the model did not predict the formation of solid flux. In fact, this can be done by considering m and R mf first, and then considering μ f (T) .

Effect of heat transfer coefficient, m, and interface thermal contact resistance, R mf
The starting point is that we must have where d exp max denotes the greatest depth of the experimentally obtained oscillation mark, and h crit = [x] H crit , with the latter then simply reducing to In effect, this gives an inequality that m and R mf must satisfy to ensure that any solid flux forms at all. This is shown in Fig. 3, which gives h crit as a function of m for different values of R mf , where we have taken d exp max = 0.33 mm from Vynnycky et al. (2017). Hence, for each value of R mf , there is clearly a critical value of m, which we shall call m cr , below which h crit > d exp max ; thus, Fig. 4 shows m cr as a function of R mf . Moreover, solutions are not possible for all values of R mf , but only for which, for the current parameter set, corresponds to 3.7 × 10 −4 m 2 K W −1 . A further subtlety is that if h crit d exp max , the theoretical profile is more likely to resemble the experimental profile. The reason for this is that it is inevitable that the solution produced by solving (5.36) will have a kink in the profile at the transition point where H 0 = H crit , as seen in Fig. 7 in Vynnycky et al. (2017). The original experimental result did not display any discernible kink, and the best way to reproduce this situation theoretically would be if the transition is close to the location where the oscillation mark begins and ends. This is achieved by having m and R mf so that d exp max − h crit is maximized. From Fig. 3, this appears to be when R mf = 0 and m is as great as possible. Note that this discussion is independent of the value of R ms . Thus, even though we have not prescribed it, its value does not affect the oscillation-mark profile; on the other hand, it would be necessary to prescribe its value in order to determine F, as seen from  (5.10), but determining F is not required for determining either S or H, because of the decoupling that occurs as a consequence of K fs 1. Observe also that the above is independent of the viscosity profile being used, and we now proceed to consider these in turn, with a view to reconstructing the experimental oscillation mark profile. the form (6. 2) The reason this happens is because μ f is assumed to be a function of the temperature alone, and we are applying constant temperature boundary conditions at the solid flux-molten flux and solid steel-molten flux interfaces, i.e. (4.9) and (5.1), respectively. It should be observed, however, that (5.1) arose because K fs 1, i.e. the thermal conductivity of the flux is much smaller than that of the steel, which then renders the steel to be at constant temperature at leading order. Indeed, the original conditions at the molten flux/solid steel interface were (2.18) and (2.19), and there was no sign, at that stage, that we would be able to reduce this to just a constant-temperature boundary condition.
In fact, since (6.1) is a depressed cubic, it is straightforward to determine the number of real roots and what they are. Writing (6.1) as there are three distinct real roots if and only one distinct real root and two nonreal complex conjugate roots if Φ < 0. Now, with Furthermore, the fact that C < 0 and D > 0 indicate that there can only be zero or two positive roots. If there are no positive roots, then the only real root will be negative, which clearly cannot satisfy H 0 > H crit ; thus the formulation can only be self-consistent if Φ > 0. (6.4) can now be narrowed down further: the formulation is self-consistent only if To find the three real roots for H 0 when this inequality is satisfied, we can use the well-known triple-angle identity, cos 3χ = 4 cos 3 χ − 3 cos χ , to find the solution for H 0 > H crit in closed form. Putting where n is an integer. The three roots can then be obtained by setting n = 0, 1, 2, although it is not clear which two of these will correspond to the positive roots, or whether one or both of these will be greater than H crit . Moreover, setting V = V 0 in (6.6), we obtain the theoretical maximum oscillation mark depth; requiring this to be equal to the experimental value, we arrive in turn at a relation between m and R, which both appear in Bi mf . This relation is given graphically in Fig. 5 for each of the viscosity profiles. Comparing these with the curve in Fig. 4 indicates that only the curve labelled 'Reynolds law (low)' lies above it, meaning that only in this case is the theoretical thickness of the solid flux layer less than the experimentally measured oscillation-mark depth; hence, of the three viscosity profiles, only the lower version of the Reynolds law would be capable of reproducing the oscillation-mark profile. However, 24 K. M. DEVINE ET AL. while we will only consider that case when we compare the model and experimental results later in Fig. 10, we will include the other two cases in the ongoing discussion.
6.2.2 H 0 < H crit . It remains to construct H 0 for H 0 < H crit , which will require the numerical solution to P H 0 , τ = 0, with P being given by (C.10) in appendix C for the Reynolds-law viscosity and (D.9) in appendix D for the power-law viscosity; the corresponding form for P when the viscosity is constant, which we shall also need in what follows, is Let τ 0 denote the value of τ at which V = 1 andV > 0, and hence H 0 = 0, so that and consider the behaviour of P for τ > τ 0 ; thus, τ = τ 0 is the start of what is known as the negative strip time, during which the mould moves downwards faster than the solidified steel. Moreover, we denote by τ crit the value of τ at which H 0 = H crit . The key point can now be illustrated by comparing the behaviour of P for the constant and power-law cases as τ varies; this is done in Figs 6 and 7, where we have arbitrarily set H crit = 0.5, although its exact value will in general depend on Bi mf , as will the exact form of P, rendering what happens in an actual case even more delicate. In Fig. 6, the behaviour is comparatively straightforward: as τ increases from τ 0 to τ crit , the curve moves down, but with P H crit , τ > 0 for τ < τ crit . It is evident that, for τ < τ crit , there is a unique root for H 0 , such that H 0 < H crit . Furthermore, for τ crit < τ < 1, there is a unique root for H 0 , such that H 0 > H crit . However, for Fig. 7 the situation is surprisingly different. For τ = τ crit , we see that P has two roots for H 0 ≤ H crit . Moreover, when τ > τ crit , it seems that there would be a further time interval, τ crit < τ < τ 2 say, during which there would be two possible roots for H 0 < H crit and one root for H 0 > H crit . In particular, this leads to the situation where H 0 cannot be a continuous function of τ , although the discontinuity can manifest itself in one of two ways: (I) at τ = τ crit , i.e. at the earliest possible opportunity, the solution jumps up to H crit from the lower strictly positive root of P H 0 , τ = 0; (II) even after τ = τ crit , one continues to use P as defined for H 0 < H crit to compute the root, meaning that H 0 would be continuous at τ = τ crit , but the solution would have to jump to the root for which H 0 > H crit some time later instead, i.e. at the latest possible opportunity.
How this would manifest itself in the prediction of an oscillation-mark profile is shown in Figs 8 and 9. These results have been generated using R mf = 0, m = 3994 Wm −2 K −1 , which give h crit = 0.21 mm, and equation (2.31), but with μ 0 = 0.02 kgm −1 s −1 instead of the value given in Table 2. A lower value of μ 0 is used in order to make the discontinuity more distinct on a plot; in fact, the discontinuity is present for all values of μ 0 but is much less visible for the value given in Table 2. Figure 8 shows the possible mark profiles, depending on whether mode I or II is adopted, with Fig. 9 showing blow-ups of the regions where a discontinuity can occur. While the discontinuity is scarcely discernible in the z-direction, it constitutes around 20% of the maximum mark depth in the x-direction. It is worth exploring what is causing this behaviour, how it might be avoided and whether it is in some way inherent in the power-law profile, since it can clearly cannot occur for the constant viscosity case. In fact, in the context of heat conduction problems in fluids having temperature-dependent viscosity, this behaviour is not surprising and is akin to the hysteresis phenomenon considered by Skul'skiy et al. (1999) and Costa & Macedonio (2002). This is explored further in appendix E, where it is shown that whether the behaviour occurs or not depends on the value of θ m,f and can even occur for a Reynolds-law viscosity profile. On another tack, however, we can observe that this also appears 26 K. M. DEVINE ET AL.  to be a nonobvious consequence of reducing the model to its current form, and in particular neglecting the terms on the left-hand sides of (4.1) and (4.4). Had these been retained, the model would have been analytically intractable requiring these equations to be solved numerically, but would not have introduced discontinuities.

Oscillation-mark profiles
Having explored how the different ways that the model can behave, we finally turn to whether it can predict an oscillation-mark profile that agrees with the experimentally obtained one focused on Downloaded from https://academic.oup.com/imamat/advance-article-abstract/doi/10.1093/imamat/hxaa010/5822756 by guest on 01 May 2020  in Vynnycky et al. (2017). Although we could arguably have tried to formulate this formally as an optimization problem over m, R mf and μ f (T) , there is little point in doing this. Instead, we proceed as indicated earlier by: • ensuring that solid flux is formed and that the maximum depths agree; • any kink that appears in the profile is restricted to the opening and closing of the mark, meaning that, for most of the negative strip time, solid flux is present.
The first of these means that (m, R mf ) lies on the curve labelled 'Reynolds law (low)' in Fig. 5, while the second implies that we should take R mf = 0. Figure 10 compares the profile generated this way with the experimentally obtained profile from Vynnycky et al. (2017), as well as the theoretical result from Vynnycky et al. (2017), for which the viscosity was assumed to be constant. As is evident, the agreement between the new model result and the experimental profile is reasonable, but with the advantage over the old model result that the current model predicts the formation of solid flux. In addition, Fig. 11 shows an enlargement of the two oscillation-mark profiles in Fig. 10; from this, it becomes clearer that the experimentally measured profiles are not identical, although it should be noted the discrepancies are in fractions of a millimetre.

Conclusions
This paper has investigated the different possible behaviours of a recent asymptotic model for oscillation-mark formation in the continuous casting of steel, with particular focus on how the results obtained vary when the heat transfer coefficient (m), the thermal resistance between the mould and the flux (R mf ) and the dependence of the viscosity of the flux powder as a function of temperature, μ f (T) , are changed. It turns out that three different outcomes are possible: (a) the flux remains in molten state and no solid flux ever forms; (b) both molten and solid flux are present, and the profile of the oscillation mark is continuous with respect to the space variable in the casting direction; (c) both molten and solid flux are present, and the profile of the oscillation mark is discontinuous with respect to the space variable in the casting direction.
Although alternative (a) gave good agreement with experimental data in Vynnycky et al. (2017), it suffered the drawback that solid flux is typically observed during continuous casting; this has been rectified in this work via alternative (b). On the other hand, alternative (c) is found to arise when the flux viscosity profile used, which has a power-law dependence on the temperature and contains a break temperature at which the viscosity becomes infinite, is one that is frequently advocated in the literature (Jonayat & Thomas, 2014); subsequently, it was found that it may arise for Reynoldslaw profiles also, although not if the viscosity is assumed constant. Moreover, behaviour (c) appears to be related to the hysteresis-type phenomenon encountered in other flows that involve temperaturedependent viscosity, although in the present case it appears to be a consequence of neglecting heat accumulation and convection terms in the energy equation for the flux; this warrants further investigation.
Lastly, we note that while it may always prove difficult to obtain values for m and R for a particular industrial casting process, it is evident that the quality of the model will improve if both μ f (T) and the final oscillation-mark profiles are available; unusually, here we had access to the latter, but unfortunately not to the former. and found that the only possibility was that H * 0 = 0. This was because if H * 0 > 0, then clearly (B.1) cannot be satisfied, since Γ > Λ and the second term must be positive, since H 0 ≥ 0 for all τ . However, we see that H 0 must also have vanished at τ = τ * − 1, τ * − 2, since V (τ * − 1) = 1 However, when V = 1, in which case we have the result when V = 1 is of significance, since it means that the integrals in (B.2), and hence (5.32), will be singular. From (B.2), it appears that the denominator of the last term on the left-hand side of (B.2) is more singular than the numerator, since the most singular terms in the integrands in the numerator and denominator behave as 1/H 2 0 and 1/H 3 0 , respectively; this suggests that the entire term involving the integrals can be neglected. Moreover, the fact that H 0 vanishes was also suggested by the computations in Hill et al. (1999), which were also for constant viscosity. Thus, the task now is to show that H 0 = 0 at τ = τ * whenμ f is not constant.