2020 年 61 巻 10 号 p. 1974-1980
In this paper, a heat transfer model between the mold and the ingot with convex and concave surfaces for continuous casting (CC) process of copper and copper alloys is proposed and discussed. Conventionally, during CC process, vibration of the mold leads to ingot with convex and concave surfaces known as oscillation marks. These marks may cause heat resistance between the mold and the ingot. In the model, three areas where heat resistances occurs were considered: (ΔR1) non-contacting area, (ΔR2) concave area derived from decrease of thermal conductivity, and (ΔR3) area where non-effective heat flow exists in solid phase. The heat resistance values were obtained either analytically or by numerical methods for conditions typically observed in the CC process. Quantitative analyses and comparison of heat resistance values indicated that ΔR3 was the most significant factor and that ΔR1 and ΔR2 was negligible. Furthermore, it was found that slight changes in contact condition results in a large change in heat resistance.
This Paper was Originally Published in Japanese in J. Japan Inst. Copper 58 (2019) 109–115.

Fig. 5 Schematic drawing of three kinds of heat resistance.
In the continuous casting (CC) process of copper and copper alloys, how the heat transfer coefficient between the mold and ingot (the reciprocal of the thermal resistance) is determined has been an important research theme for a long time. This heat transfer coefficient is thought to fluctuate due to various factors. For example, it depends on the alloy type, the size of the ingot and the casting conditions. In some cases, even under the same casting conditions, it can also fluctuate over time due to some reasons. Fluctuations in the heat transfer coefficient have a large influence on the quality of the ingot, for instance, they can cause internal defects. However, the mechanism of determination and fluctuation of the value of the heat transfer coefficient has not been clarified.
In general, on the surface of an ingot produced by the CC process of copper or copper alloys, wave-like periodic convex and concave areas (hereinafter referred to as uneven areas or unevenness) are formed as the result of oscillation, that is, the vibration of the mold.1) The wavelength and amplitude of the unevenness is in the order of 10−1 to 101 mm and the formation of the unevenness may affect the heat transfer coefficient between the mold and ingot. Similar unevenness has been observed in continuous casting of steel, which has been reported to cause local delay of solidification and cracking of the ingot surface.2)
To date, several reports have been made on heat transfer models between the mold and ingot with an uneven surface, mainly in the steel industry. Ya Meng et al.3) assumed that the flux would enter the concave areas, and calculated the thermal resistance using the thermal conductivity of the flux. Shibata et al.4) assumed that gas is filled in the concave areas and proposed a model where the thermal resistance is calculated by using both thermal conductivity of air and radiation. However, there are not many studies that have examined the thermal resistance in detail by paying attention to the shape of unevenness.
The purpose of this study is to quantify the increase in thermal resistance of the ingot with an uneven surface according to the shape of the unevenness with respect to the thermal resistance of the ingot with a flat surface. In particular, we expected that some areas in solid phase close to concave areas become so-called thermal dead spaces that do not contribute to heat extraction, and the thermal resistance is increased by this. We named this thermal dead space “stagnation of heat” and carried out analysis in detail. Including this “stagnation of heat”, we decomposed the increased thermal resistance into three parts, and formulas that determine these have been established, using a total of eight parameters that mainly characterize the shape of the unevenness. Using these formulas, we clarified the dominant factors of the increase in thermal resistance, and discussed the causes of the thermal resistance fluctuations observed in the actual CC process.
The aim of this heat transfer model is to calculate the average thermal resistance in the direction perpendicular to the surface of the ingot (hereafter such a direction is referred to as heat extraction direction or y-direction) when two-dimensional, periodical and uniform unevenness is formed on the surface of the ingot in the casting direction (x-direction) during the CC process, as shown in Fig. 1(a). Here, two-dimensional means that, for example, in the case of a cylindrical ingot, the unevenness is formed in rotational symmetry with respect to the central axis of the cylinder. Similarly, in the case of a prismatic ingot, it was assumed that the shape of the unevenness does not change in the circumferential direction.

Schematic drawing of the ingot and convex/concave.
In this study, the ingot size (radius or thickness of the ingot cross-section) was assumed to be more than one order larger than the wavelength and amplitude of the unevenness. Therefore, as the surface can be regarded as substantially planar macroscopically even in the case of a cylindrical ingot, it was assumed that the unevenness is formed on a two-dimensional flat plate that extends infinitely in z-direction of the Cartesian coordinate system as shown in Fig. 1(b). The shape of the unevenness was expressed by a two-dimensional cross section, and the thermal resistance and heat transfer coefficient were evaluated by the value per unit thickness unless otherwise specified.
Furthermore, the following assumptions were used for simplicity.
(1) The shape of the convex portion in uneven areas constituted by the solid phase is regarded as trapezoidal or rectangular, as shown in Fig. 2.

Schematic drawing of the model of convex/concave.
(2) The uneven areas have already been formed at the beginning of solidification, and the shape does not change thereafter.
(3) The interface between the solid phase and the melt is flat. Therefore, the solid-phase is divided into two parts, those are, the convex part and the flat plate part, as shown in Fig. 2.
(4) The concave areas outside the ingot are considered to be adiabatic, and radiation and thermal conduction are ignored for simplicity, since the purpose of this study is to analyze the effect of the shape of unevenness. Generally, in the CC process, concave areas are filled with lubricant or gas, so radiation and/or thermal conduction are not always negligible. However, the melting point of copper is lower than that of steel and the solidification of copper is a phenomenon that occurs in a temperature range where the effect of radiation is relatively small. And the thermal conductivity of the copper solid phase is more than one order higher than the flux used for lubricants and the thermal conductivity of the gas. Therefore, this assumption is not considered to be significantly different from the actual situation.
(5) The mold and solid phase are uniformly contacted with the upper side of the trapezoid. As shown in Fig. 3, the local heat transfer coefficient $\text{h}_{\ell }$[W/m2K] between the mold and ingot is set to h0 in the contacting positions (hereafter h0 is referred to as the heat transfer coefficient of the contacting positions), and is set to 0 in non-contacting positions, because the concave areas are adiabatic. In case of real CC process, the heat transfer coefficient might not be constant even within a contacting position but it might change according to the position as shown by the arrow and dotted line at the right end of Fig. 3. On the other hand, the depth of the unevenness studied in this model is in the order of mm, resulting in a large difference in the thermal resistance between the contacting and non-contacting positions. For this reason, we concluded that it is unnecessary to model the slight changes in heat transfer inside the contacting position in detail, and used the simple model. Moreover, the h0 value may vary with the progress of the solidification of the ingot, because it may change depending on, for example, the air gap formed between the mold and ingot by deformation of solidification shell. However, since the main purpose of this model is to analyze the effects of the shape of unevenness, fluctuations in the h0 value were not considered. Consequently, the value of h0 was assumed to be constant regardless of the shape of the unevenness and the time elapsed from the start of solidification.

The model of heat transfer coefficient.
(6) There is no thermal resistance between the melt and solid phase.
(7) The thermal conductivity k[W/mK] of the solid phase is constant regardless of temperature and position.
(8) The shape of the solid phase is represented by the following six parameters as shown in Fig. 4.

Six parameters describing the shape of convex/concave.
(9) Since the unevenness is repeated periodically in the casting direction (x-direction), it is sufficient for the modeling area to have only the length L in x-direction, as shown in Fig. 4. As the value of $\ell $ is equal to the solidification thickness, it increases with the progress of solidification, but it was considered to be uniform within the area of the model.
In this model, the increase in the resistance was expressed using eight parameters in total, consisting of six parameters related to the shape (L, D, $\ell $, A, B and C), a physical property k and an engineering parameter h0.
2.2 Assumptions about increase in resistance and formula derivation methodIn this model, as shown in Fig. 5, it was assumed that the following three types of thermal resistances (unit is m2K/W) are generated when unevenness is formed, compared to the case of the ingot with a flat surface.

Schematic drawing of three kinds of heat resistance.
Then, the increase in total thermal resistance ΔR is expressed as ΔR = ΔR1 + ΔR2 + ΔR3.
Here, “stagnation of heat” means the solid phase portions that become dead spaces from the viewpoint of the heat flow from the molten metal to the mold and do not sufficiently contribute to the heat extraction. It is thought to be generated mainly in the place near the concave area shown by the circle in Fig. 5(c).
The equations that determine ΔR1, ΔR2 and ΔR3 were established using the 8 parameters defined in the previous section. The equations for ΔR1 and ΔR2 were expressed using the basic formulas for heat transfer. Regarding ΔR3, it was difficult to express it using an analytical solution, so an approximate expression was obtained by performing a numerical analysis using the two-dimensional finite difference method.
When the surface of the ingot is uneven, the average heat transfer coefficient between the mold and the solid phase by contact is represented by hav = A · h0. Therefore, the thermal resistance R1 per unit cross section is expressed as R1 = 1/(A · h0). On the other hand, since the thermal resistance R1p in case of the ingot with a flat surface is 1/h0, increased thermal resistance ΔR1 can be expressed as follows.
| \begin{equation} \Delta R_{1} = R_{1} - R_{1p} = \frac{(1 - A)}{A \cdot h_{0}} \end{equation} | (1) |
When the shape of the convex portion is rectangular instead of trapezoidal, the solid phase cross-sectional area in the uneven areas effective for thermal conduction in the heat extraction direction is thought to decrease from 1 to B because concave areas are assumed to be adiabatic. Therefore, the apparent thermal conductivity k2 of the uneven areas can be expressed as follows.
| \begin{equation} k_{2} = B \cdot k\quad (B < 1) \end{equation} | (2) |
| \begin{equation} R_{2} = \frac{D}{B \cdot k} \end{equation} | (3) |
The thermal resistance R2p in the absence of unevenness is expressed as R2p = D/k. Consequently, the increase in the resistance ΔR2 can be expressed as follows.
| \begin{equation} \Delta R_{2} = R_{2} - R_{2p} = \frac{D(1 - B)}{B \cdot k} \end{equation} | (4) |
In the case of a trapezoid, since the width changes in the height direction (y-direction), the thermal resistance also changes continuously in the height direction in proportion to the reciprocal of the width. The thermal resistance between the upper and lower sides of the trapezoid can be approximately obtained by integrating the reciprocal of the width in the height direction, assuming that continuously changing thermal resistors are connected in series. Here, the ratio of the trapezoidal thermal resistance value to that of the rectangle with the same area and height, r is determined. As shown in Fig. 6(a), the shape of the trapezoid is defined as base 1 + d, top 1 − d and height 1. The width of this trapezoid w at an arbitrary height position v is expressed by the following equation.
| \begin{equation*} w = -2d \cdot v + (1 + d) \end{equation*} |

Heat resistance of trapezoid.
When the thermal conductivity is 1, the relationship between d and r can be expressed by eq. (5) and Fig. 6(b).
| \begin{align*} r &= \int\limits_{0}^{1}\frac{dv}{w} = \int\limits_{0}^{1}\frac{dv}{-2d \cdot v + (1 + d)} \\ &= \left[-\frac{1}{2d}\ln(-2d \cdot v + 1 + d)\right]_{0}^{1} \end{align*} |
| \begin{equation} r = \frac{\mathit{ln}\{(1 + d)/(1 - d)\}}{2d} \end{equation} | (5) |
As can be seen from eq. (5) and Fig. 6(b), even for a very sharp trapezoid such as d = 0.9 (the ratio of the upper to the lower sides is 1:19), the value of r is at most 1.6, and the difference in the thermal resistance between the trapezoid and rectangle is not very large. In other words, as shown in Fig. 5(b′), the difference between the calculated value and the true value is not so large even if calculations are performed using rectangles instead of trapezoids. For this reason, ΔR2 was evaluated using eq. (4) thereafter. In addition, the value of B can be changed formally independently of A and C so that this formula can be applied to shapes other than trapezoids.
3.3 Increase in thermal resistance due to stagnation of heat: ΔR3 3.3.1 Numerical analysis modelAs shown in Fig. 5(c), the area of numerical analysis for obtaining ΔR3 was limited to the rectangular area in the flat plate part of solid phase. For the boundary conditions, it was assumed that heat flows in from the entire lower side of the rectangular area which is contacted to the molten metal, and heat flows out from a part of the upper side that is connected to the convex part. The temperature of the lower side of the rectangular was set to 1000°C, the part of the upper side, which is connected to the convex part was set to 0°C, and the part of the upper side which is not connected to the convex part was set to adiabatic. Since the temperature distribution in x-direction repeats periodically, the vertical sides of the rectangular were set as the adiabatic boundary condition. Under these boundary conditions, two-dimensional thermal analysis was performed to obtain the internal temperature distribution. In addition, using these results, the heat flux at each position was calculated and the thermal resistance between the top and bottom sides of the trapezoid were obtained. The values of four parameters related to this area were changed in the series of calculation, those were, solid phase thickness ($\ell $), solid phase thermal conductivity (k), one half of the wavelength (L) and connection rate (C).
This numerical analysis is not solidification analysis, but numerical analysis for calculating the thermal resistance of a certain shape.
3.3.2 Numerical analysis resultsFigure 7(a) shows the temperature distribution of the flat plate part when C = 0.5 and $\ell = 2\text{L}$, as an example of numerical analysis results. The width of 4L is shown repeatedly in x-direction in Fig. 7(a). For comparison, Fig. 7(b) shows the temperature distribution in the case without unevenness and the heat flux is completely uniform and one-directional. The isothermal line means that the wider the interval, the smaller the temperature gradient and the smaller the heat flux. In Fig. 7(a), the area shown by a black circle has a wider isotherm interval and smaller heat flux than the area shown in Fig. 7(b), and it indicates that there is a “stagnation of heat” there. On the contrary, the interval of the isothermal lines is narrow around the area shown by the black arrow, and due to the influence of this, the interval of the isothermal lines is wide and the heat flux is small around the area shown by the white arrow. Around the area shown by the white arrow, the isothermal line is perpendicular to the heat extraction direction and almost flat, so the heat flux around this area as equal to the average heat flux of this shape. This indicates that the total thermal resistance between the upper and lower sides is increased due to the unevenness.

Temperature distribution of the plate area.
As described above, by comparing Fig. 7(a) and Fig. 7(b) it is clearly shown that the stagnation of heat is formed in the flat plate part, and the mechanism by which this phenomenon occurs was clarified intuitively and visually.
Figure 8 shows the relationship between the thickness of the flat plate part $\ell $ and ΔR3 when the values of the thermal conductivity of the solid phase k and the wavelength L are changed. The value of C is fixed at 0.5. Here, as shown in the Fig. 8, the value of $\ell $ and ΔR3 were made dimensionless and are expressed as $\ell ^{ \sim } = \ell /\text{L}$ ($\ell ^{ \sim }$ is the dimensionless solid phase thickness) and ΔR3∼ = k/ΔR3/L (ΔR3∼ is the dimensionless increase in resistance), respectively. It can be seen that all values of ΔR3∼ is almost independent of L and k.

The relationship between ΔR3∼ and $\ell ^{ \sim }$.
From Fig. 8, it can be seen that the ΔR3∼ value increases with the increase of $\ell ^{ \sim }$ when $\ell ^{ \sim }$ is small, and reaches the maximum value (denoted as ΔR3∼*) at $\ell ^{ \sim } \fallingdotseq 1$ ($\ell \fallingdotseq \text{L}$), with no increase thereafter. This is because, as shown in Fig. 7(a), “stagnation of heat” exists only in the vicinity of uneven areas, and does not exist at a distance L or more from uneven areas.
As described above, since ΔR3∼* does not depend on L, k and $\ell $, ΔR3∼* is considered to be a function of only the connection rate C. Then, when the relationship between C and ΔR3∼* was obtained from the results of numerical analysis, as shown in Fig. 9, It turned out that ΔR3∼*/(1 − C) and −ln (C) have a nearly linear relationship. Consequently, it was found that ΔR3* can be expressed by the following approximate equation.
| \begin{equation} \varDelta R_{3}^{*} = -\frac{0.6071 \cdot L(1 - C) \cdot \mathit{ln}\,C}{k} \end{equation} | (6) |

The relationship between ΔR3∼* and C.
As described above, the total thermal resistance increase due to the formation of the unevenness ΔR can be calculated using eqs. (1), (4) and (6).
In order to analyze fluctuations in the thermal resistance occurring in the CC process using this model, the mold temperature in the actual CC process was measured and the changes in the thermal resistance depending on the casting conditions were calculated.
4.1 The method and the results of temperature measurementThe measurement of the mold temperature was carried out on a low-alloy cylindrical ingot under conditions 1 and 2. Under conditions 1 and 2, the casting speed, molten metal temperature (i.e. total amount of heat input), and the mold used were not changed while other conditions were changed. The temperature was measured at three locations, those are, upper, middle and lower points of the mold by changing the mold height in the casting direction. Measurement was performed using a K type thermocouple at 1-second intervals for 20 minutes, and the average value was taken as the temperature at each position.
Figure 10 shows the measurement results. Under conditions 1 and 2, the mold temperature changed by several tens of K. This is considered to be due to the change in the heat transfer coefficient (thermal resistance) between the mold and ingot.

The result of mold temperature measurement.
In order to calculate the difference in thermal resistance between condition 1 and condition 2 from the measured temperature data, a two-dimensional axisymmetric thermal conductivity analysis of the CC process was performed using casting conditions as input conditions. First, in order to determine the overall heat transfer coefficient (heat transfer coefficient between the cooling water and ingot) under condition 1, numerical calculations were performed so that the mold temperature measured under condition 1 and the mold temperature calculated under condition 1 match.
As the overall heat transfer coefficient includes two unknown values; one is the heat transfer coefficient between the mold and ingot and the other is the heat transfer coefficient between the mold and cooling water (the water cooling transfer coefficient), in order to determine the both heat transfer coefficients respectively, more information was needed. Consequently, calculation was also performed so that the depth of the molten pool (distance from the molten metal surface position to the deepest point of the pool) measured under condition 1 and the depth of the molten pool calculated under condition 1 become almost same. Next, in order to obtain the heat transfer coefficient between the mold and ingot under the condition 2, calculation was performed so that the mold temperature measured under the condition 2 and the mold temperature calculated under the condition 2 were almost the same value, without changing the water cooling transfer coefficient that was obtained before.
As the result, the values of the overall heat transfer coefficient were as follows.
From these results, the difference in thermal resistance ΔRx between condition 1 and condition 2 was calculated as ΔRx = 0.90 × 10−4 [m2K/W]. According to the simulation results, the depth of the molten pool under condition 2 was estimated to be about 8.2% deeper than that under condition 1. Although detailed description is omitted, the ingot produced under the condition 1 had internal defects, and the ingot produced under the condition 2 was sound. As explained so far, in the actual continuous casting process, the heat transfer coefficient of about 10% easily fluctuates, which may cause a great change in the quality of the ingot.
Using eqs. (1), (4) and (6), the eight parameters were varied to calculate the values ΔR1, ΔR2 and ΔR3*, and the magnitude of the values was compared. For reference, the thermal resistance of the solid phase $\text{R}_{\text{s}} = \ell /\text{k}$ was also calculated and compared. The range of variation of the parameters was as follows, which is usually seen for copper and copper alloys.
Since a low-concentration copper alloy with a conductivity of 20% IACS or higher was assumed, k was set to 50 to 400 W/mK. L was set to 10 mm at maximum because it is almost the same as the wavelength of the unevenness. D was also set to 10 mm at maximum because it is the depth of the unevenness, and $\ell $ was set to 0 to 200 mm because it is the radius of the ingot or one half of the thickness of the ingot. In the range of the general shape of unevenness, the value of A, B and C varied from 0.1 to 0.9, and A ≦ B ≦ C, because the shape of the convex part is trapezoidal and it becomes thinner from the bottom to the top. Since the h0 value is equivalent to the average heat transfer coefficient between the mold and ingot without unevenness, it was estimated to be approximately 1000 to 2000 W/m2K from actual measurement experience. It is unknown whether the value is the same for the ingot with an uneven surface, but in this analysis, the maximum was set to 2000 W/m2K.
First, as a common example of copper, the combination of parameter values shown in set 1 of Table 1 is called the default condition. As a result of this combination, ΔR2 and ΔR3* were extremely small and negligible compared to ΔR1 and Rs.

Next, we investigated whether the combinations of parameters may or may not exist such that ΔR2 and ΔR3* equal to or greater than ΔR1. As shown in Set2, when A < B < C with respect to the default condition, ΔR2, ΔR3* decreased significantly more than ΔR1. Moreover, as shown in Set3, when values A, B and C increased simultaneously while maintaining A = B = C with respect to the default condition, ΔR1 and ΔR2 decreased simultaneously, and ΔR3* decreased significantly. As shown in Set4 to Set7, the value of ΔR1 was always maximum even when L, D, k and h0 were appropriately varied from the default value. Furthermore, as shown in Set 8, ΔR1 was maximum even when all the parameters were selected so that the values of ΔR2 and ΔR3* became the largest relative to ΔR1.
From the above results, for copper and copper alloys, the increase in resistance when unevenness is formed on the surface of the ingot is approximately ΔR1 ≫ ΔR2 ≒ ΔR3*, indicating that ΔR1 was dominant. The effect of “stagnation of heat” focused on in this study was not so great.
5.2 Estimation of factors of thermal resistance change due to differences in casting conditionsThe reasons for the change in the thermal resistance of casting conditions 1 and 2 calculated in Chapter 4 were estimated using the results obtained so far. First, from the results of the discussion in Section 5.1, the change in ΔR is considered to be due to the fact that only ΔR1 changed. Consequently;
| \begin{equation} \Delta R_{x} = \Delta R_{1}(\textit{Cond.}\ 2) - \Delta R_{1}(\textit{cond.}\ 1) \end{equation} | (7) |
Here, substituting the values of the contact rate A under conditions 1 and 2 at A1 and A2, respectively, and the values of the respective heat transfer coefficients h0 to be h01 and h02, the following is established.
| \begin{equation} \Delta R_{x} = \frac{1}{A_{2} \cdot h_{02}} - \frac{1}{A_{1} \cdot h_{01}} = 0.90 \times 10^{-4} \end{equation} | (8) |
Since eq. (8) included only two kinds of parameters A and h0, the change in the change in thermal resistance from condition 1 to condition 2, ΔRx is thought to be the results of the change in A or h0. Therefore, the two cases in which either A or h0 has changed were reviewed as described below.
First, in Case 1, the change in A was considered. Here, when A1 = 0.5 and h01 = h02 = 2000 W/m2K, the value of A2 became 0.458, and the difference in the value of A between conditions 1 and 2 became about 10%. Even if the value of A1 was changed to a value other than 0.5, the change in A2 was almost the same. Although the shape of the unevenness of the ingot produced in this casting test was not been measured, such a change in the shape of unevenness may occur generally in the continuous casting of copper and copper alloys.
Therefore, it was clearly demonstrated that the cause of the fluctuation obtained by this actual CC process data can be explained by the change in the shape of the unevenness.
Next, in Case 2, only h0 was changed. Here, when A1 = A2 = 0.5, and h01 = 2000 W/m2K, the value of h02 became 1833 W/m2K, and the change was found to be about 10%. The value of h0 changes due to some cause other than the formation of unevenness, and is considered to change depending on, for example, the air gap due to the detachment of the ingot from the mold. The possibility that h0 changes by about 10% due to these factors cannot be denied at present. Therefore, the fluctuations in thermal resistance observed in the actual CC process this time could be due to the change in the contact ratio A, but it was not confirmed that it was the only cause.
The increase in thermal resistance when unevenness is formed on the surface of an ingot produced by the CC process was divided into three factors (ΔR1, ΔR2 and ΔR3), and each thermal resistance increase value was expressed using a simple equation consisting of eight parameters. Then, using these equations, the causes of the fluctuation in the thermal resistance observed in the actual CC process were discussed.
As a result, it was found that the most dominant factor of the fluctuation in the thermal resistance was the contact ratio, A.
In addition, it was found that if the cause of the fluctuation in thermal resistance observed in the actual CC process is limited only to the change of the shape of the unevenness, the cause can be explained by the change in A.