NAnL INST OF STAND a TECH

ill 11 II 1 A111D7 3^0773

NBSIR 84-2882

NBS: Materials Measurements

U S. DEPARTMENT OF COMMERCE

National Bureau of Standards

Washington, DC 20234

June 1984

Annual Report

For the Period 1 April 1983 - 31 March 1984

NASA Government Order H-27954B

Properties of Electronic Materials

100

• U56

0^-2882 . 1984

. V

U S. DEPARTMENT OF COMMERCE

NATIONAL BUREAU OF STANDARDS

KATIOMAL r-UKEA.U OF t .A.-.IDARDS LiEPJvny ,

NBSIR 84-2882

'

NBS: MATERIALS MEASUREMENTS

)

J. R. Manning

U S. DEPARTMENT OF COMMERCE

National Bureau of Standards

Washington, DC 20234

June 1984

Annual Report

For the Period 1 April 1983 - 31 March 1984

NASA Government Order H-27954B

Properties of Electronic Materials

U.S. DEPARTMENT OF COMMERCE, Malcolm Baldrige, Secretary

NATIONAL BUREAU OF STANDARDS, EmMt Ambtor. Dirtctor

O

o

TABLE OF CONTENTS

Page

Summary 1

Task 1 - SURFACE TENSIONS AND THEIR VARIATIONS WITH TEMPERATURE

AND IMPURITIES

S. C. Hardy and J. Fine 5

Task 2 - CONVECTION DURING UNIDIRECTIONAL SOLIDIFICATION

S. R. Coriell, R. J. Schaefer, and G. B. McFadden .... 19

APPENDIX A: Asymmetric Instabilities in Buoyancy

Driven Flow in a Tall Vertical Annulus 42

APPENDIX B: Morphological Stability in the Presence

of Fluid Flow in the Melt 53

APPENDIX C: Effect of a Forced Couette Flow on Coupled

Convective and Morphological Instabilities During

Unidirectional Solidification 81

APPENDIX D: Convection-Induced Distortion of a Solid-Liquid Interface 105

Task 3 - MEASUREMENT OF HIGH TEMPERATURE THERMODYNAMIC

PROPERTIES

D. W. Bonnell 135

1
war"

.... %

Ml

'Mi m

Ei , 'vjj ^1

P^:p

’ •,\ ’*'.

’-f -y.- *' •'*- -A

Wmrn.

_: -r*i| ";^^>^' j® >' ^3?wfn ^ ^ , . ^ - a s: -v(^ ®

is:.,-

« i- t * • */-' ^^ *'*_ • * ®;'

SfStfjJi^VJK- v'"

"* '- r*'* bS^^ . P?-.

I:' is ^ *a^, -j^' ^r~^: „v ^ " m3Hl«ui.A3k- m ztumom '

!BT" # ,f'., ip^-, ' '•.^w;. *\'4 y/,-;jR' -: 5 .- > Tlan^<!>8 ,W *0

i&i'’ -.y-y

'•< ^wl

-><!»

^rr:- ?

, lf»i

-*!';' "'n

..ifK-:

se^

'

J-' 'i«i-

%ibk

\r

:3
NATIONAL BUREAU OF STANDARDS

MATERIALS MEASUREMENTS

Summary

This report describes NBS work for NASA in support of NASA's Microgravity

Science and Applications Program under NASA Government Order H-27954B

(Properties of Electronic Materials) covering the period April 1, 1983 to

March 31, 1984. This work is being carried out in three independent tasks:

Task 1. Surface Tensions and Their Variations with Temperature and

Impurities

Task 2. Convection During Unidirectional Solidification

Task 3. Measurement of High Tetrperature Thermodynami c Properties

The results for each task are given separately in the body of the report.

With the advent of the Space Shuttle and the accompanying availability

of space flight opportunities, it may become feasible to exploit the unique

microgravity environment of space flight to produce improved materials and

improved measurements of important materials properties. Reduction

of convection and enhanced capabilities for containerless processing are

two advantages frequently cited for materials processing in space. Both of

these topics are considered in the current work. Tasks 1 and 2 are directed

toward determining how the reduced gravity obtained in space flight can

affect convection and solidification processes. Emphasis in Task 3 is on

development of levitation and containerless processing techniques which can

be applied in space flight to provide thermodynamic measurements of reactive

materi al s

.

In materials processing on earth, gravity frequently produces densitydriven convection, thereby causing liquids to be stirred as they solidify.

This stirring disturbs the quiescent boundary layer at the solidifying
interface and can be very undesirable if nearly perfect crystals are required.

For example, it creates interface instabilities, introduces segregation of

components, and produces crystal defects in the resulting solid material.

These defects and inhomogeneities, which are particularly troublesome in

electronic technology and other advanced technical applications, might be

avoided in materials produced under microgravity conditions. Different

aspects of these convection processes are being considered in the first two

tasks. Surface tension gradients which may be the main source of convection

in space flight experiments are the measurement focus in Task 1; whereas

interactions between solidification processes, convection effects, and

interface stability, as influenced by gravity and the lack thereof, are

under investigation in Task 2.

In Task 1, measurements are being made of the surface tension of liquid

silicon. Silicon is one of the most important materials used in electronic

technology. Because of the need for a fully defect-free silicon, not

readily obtainable under earth-bound conditions where gravity driven

convection can be important, silicon is a major candidate material for

processing in space. In planning such space processing, it is important to

know the dependence of silicon surface tension on temperature and impurities,

particularly oxygen. Surface tension gradients produce Marangoni convection

and may be the main source of fluid flows which affect solidification

processes in space. Despite the technological importance of silicon, the

variations of its surface tension with temperatures and impurities were

not known prior to initiation of the current work. Since silicon is very

reactive, previous measurements had been for the most part non-systematic

and unreliable. The current sessile drop experiments have established

2
that the temperature (T) dependence of the surface tension (y) is given

by dy/dT = -0.28 mJ/m^K for liquid silicon in an atmosphere of purified

argon in the 200 degree temperature range above the silicon melting point.

Experiments are now underway to determine what appears to be an appreciable

dependence of y and dy/dT on oxygen impurities.

In Task 2, convective phenomena which arise during directional

solidification are being investigated by both theoretical calculations and

experiments. Direct measurements of convection and associated solidification

effects are being made during directional solidification of transparent

materials, particularly succi nonit ri le containing small amounts of ethanol;

and calculations are being made of convection, interface stability and

impurity segregation that should be expected during directional solidification.

The experimental observations in ground-based measurements at NBS of unidirectional upwards solidification of succi noni t ri le containing ethanol

have demonstrated that convective flow plays an important role in determining

the site at which interfacial stabilities develop. It was also determined

that the presence of a vapor bubble at the interface resulted in vigorous

convective flow, induced by surface tension gradients. Directional solidification techniques are the leading methods proposed for obtaining more

nearly defect-free materials in space. It is important in planning such

processing techniques to consider the role convection effects, including

previously unexpected effects found in the present work, could play under

reduced gravity conditions. The theoretical calculations thus specifically

include determinations of the effects a reduction in gravitational force

will have on convection, as will be important, for example, in space flight

experiments in the MEPHISTO project. Calculations have been made of the

3
solute, temperature, and flow fields at constant solidification velocity

taking into account solute rejection at the solidification interface. This

solute rejection can strongly influence convection effects. Linear stability

calculations for a flow parallel to a crystal -melt interface also have been

performed. These calculations demonstrated a coupling between hydrodynamic

and morphological instabilities and are in general agreement with experiments

carried out by Glicksman and colleagues.

In Task 3, assistance is being provided to a joint project involving

investigators from Rice University (Professor J. Margrave) and General

Electric Co. (Dr. R. T. Frost) in which a General Electric electromagnetic

levitation device is being used to develop levitation/drop calorimetry

techniques having space flight applications. The critical advantage of

such containerless processing experiments is the avoidance of specimen

contamination that would occur upon contact of a reactive material with

container walls during heating and melting. Measurements of the specific

heat of tungsten in earth-bound experiments have yielded some promising

results and also have allowed evaluation of certain limitations imposed

by gravity in such measurements. A number of automated techniques being

developed in this work should also be useful for space flight experiments.

4
Task 1 Surface Tensions and Their Variations with Tenperature and Impurities

S, C. Hardy

Metallurgy Division

Center for Materials Science

and

J. Fine

Surface Science Division

Center for Chemical Physics

SUMMARY

The surface tension of liquid silicon has been measured as a function

of temperature in purified argon atmospheres using the sessile drop technique.

The measurements show the surface tension is sensitive to low levels of an

impurity which is probably oxygen. The highest surface tension values

obtained under conditions which minimized the oxygen levels in the apparatus

are in good agreement with an isolated previous measurement in pure hydrogen.

The surface tension decreases linearly with increasing temperature and has

a temperature coefficient of -0.28 mJ/m^K.

5
Introducti on

The surface tension of silicon and its variation with temperature are

not known accurately. These data are of interest because it is now thought

that surface tension driving forces contribute significantly to fluid flow

in the Czochralski and float zone crystal growth of silicon. Fluid flows

arising from surface tensions are referred to as Marangoni or thermocapillary

flows. They are important in normal gravity and are the dominant type of

fluid flow in low gravity space experiments with silicon for which buoyancy

convection is negligibly small. Experiments involving free liquid silicon

surfaces are presently being developed for future space shuttle flights.

These surface tension measurements are essential for proper experimental

design and eventual quantitative i nterpretation of the low gravity observations.

The driving force for Marangoni flow is the surface tension gradient

rather than the surface tension itself. Thus the flow depends sensitively

on the value of the surface tension temperature coefficient K = dy/dT.

For silicon estimates of K range from approximately -0.1 mJ/m^K to -0.4 mJ/m^K.

Recent modelling of the Czochralski process for silicon found that this

variation in K gave significantly different flow patterns^. Thus a quantitative understanding of silicon processing is to some extent dependent on

accurate measurements of the temperature variation of the surface tension.

Although its high melting point, vapor pressure and reactivity complicate

the experiments, the uncertainty in K for silicon may be attributed largely

to the lack of systematic measurements of surface tension over a temperature

range. In general it is not possible to calculate K accurately from isolated

measurements of surface tension by different techniques or even by the same

technique in different laboratories. The small systematic errors in the

6
measurements are usually large enough to obscure the variation of surface

tension due to temperature. The only previous measurement of the silicon

surface tension over a temperature range is that of Lukin et al^. They

found K to be near -0.1 mJ/m^K. Although this is a reasonable value, the

density variation with tenperature used in their work is significantly

different from the generally accepted dependence. Their results are discussed

in more detail below.

Experimental

The surface tensions in this work were determined using the sessile

drop technique. This method is based on a comparison of the profile of a

liquid drop with the profile calculated by solving the Young-Laplace equation.

The comparison can be made in several ways; the traditional Bashforth-Adams^

procedure has been used here in conjunction with recently calculated drop

shape tables which virtually eliminate interpolation errors.

The high temperatures required to melt silicon are generated in an

inductively heated furnace. The silicon is contained in a shallow cup

which stabilizes the drop and imposes a circular cross section on the

drop base. The outside diameter of the cup provides the length standard

required in the data analysis. The cup and drop are at the center of a

cylindrical susceptor which lies along the axis of a quartz tube. This

tube connects through an 0-ring coupling to a vacuum system consisting of

a water cooled baffle and diffusion pump, trap, forepump, gauges, and

valves for admitting gases. Quartz windows at either end of the tube allow

simultaneous photography and temperature measurement of the drop with an

automatic optical pyrometer. After window corrections, the temperature

read by the optical pyrometer on slowly melting the sample was within a few

degrees of the silicon melting point.

7
The susceptor is constructed of three concentric cylinders of molybdenum

with an inner cylinder of tantalum. The tantalum liner has a much lower

vapor pressure than molybdenum and has been inserted to reduce the evaporation

of metal onto the sample during preliminary vacuum degassing. Several layers

of zirconia insulation surround but do not touch the susceptor. There is

no indication of any decomposition of the insulation at the temperatures

used in this work when the insulation does not touch the susceptor.

The induction heater has a power rating of 2.5 kW and operates at

450 kHz. The skin depth at this frequency is about 0.06 cm and calculation

indicates that about 10% of the field penetrates the four layer susceptor.

However, no difference was found in measured surface tension values when

the four layer susceptor was reduced to three layers. This suggests that

field penetration is not affecting the measurements significantly.

The choice of a cup material is critical because of the high reactivity

of liquid silicon and the sensitivity of the surface tension to even low

bulk concentrations of surface active impurities. Studies of the compatibility

of liquid silicon with various substances indicated that boron nitride was

a suitable container because of its insolubility in liquid silicon^.

Although the silicon may acquire some free boron when in contact with boron

nitride, this should not introduce serious errors in surface tension measurements because boron has a higher surface tension than silicon and does not

surface segregate. The development of a bulk boron concentration in the

silicon can measurably increase the surface tension, but only at high

concentrations^. We have found no change in surface tension with prolonged

high temperature contact with the BN cups. We conclude, therefore, that

boron contamination is not appreciably affecting these measurements. The

8
solid silicon does not bond to the cups, possibly because of the high

contact angle («130°) when liquid. Consequently the cup is not destroyed

by differences in thermal contraction on cooling and may be used for many

experiments.

The measurements are made in flowing argon purified by passage through

a titanium gettering furnace. The purity of the argon has not been measured

but the manufacturer claims that oxygen is reduced to the lO"^ ppm level

by the titanium getter. Additional purification of the argon is provided

by the molybdenum susceptor which is an efficient oxygen getter^. The

silicon used in these measurements is of six nines purity. The samples are

cut from a large boule using a carbide wheel and then ground to fit the

cup. The maximum diameter of the drops formed on melting is approximately

1.5 cm. Large drops are desirable to improve the precision of the measurement.^

The sessile drop technique is capable of highly accurate surface

tension measurements^. The potential of the technique can be realized,

however, only if the dimensions of the drop profile are measured to an

accuracy of a few parts in 10^. This requires a lens of large numerical

aperture to obtain sufficient resolution. We use an f/5.6 enlarging lens

with a 190 mm focal length in conjunction with a 4" x 5" camera body. The

image of the drop is formed on a ground glass plate which is observed with

a low power telescope for precise focussing. The drop is photographed

using its own radiation on a high contrast, fine grained film designed for

graphic arts work. The negatives are examined at lOOx on a measuring

microscope for analysis. The stage of the microscope can be rotated for

the precise orientation of the negative required in the analysis®.

In the normal experimental procedure the system is evacuated and

degassed at just below the silicon melting point. After cooling and starting

9
the argon flow, the temperature is slowly raised to the maximum (»1600°C)

and the sessile drop photographed. The power is then reduced and the drop

photographed again when temperature stability is achieved. This procedure

is repeated until the drop freezes.

After cooling, the silicon drop is removed from the apparatus and

subjected to Auger electron spectroscopic analysis. The only measurable

impurities are carbon and oxygen which is typical of samples exposed to

atmosphere. This layer is easily removed by argon ion sputtering. Although

it is possible that the solidification process redistributes impurities to

areas not accessible to the Auger analysis, we find no indication of any

specimen contamination by the cup or susceptor.

Results and Discussion

Although previous study has found little difference in measurements

with pure and oxygen doped silicon^, there is other evidence suggesting

that oxygen in dilute concentrations severely depresses the surface tension

of silicon. Kingery and HumeniklO observed the surface tension of liquid

silicon to increase dramatically when the atmosphere was changed from argon

to hydrogen. This increase in surface tension was attributed to a reduction

of oxygen activity. Our measurements also show a strong sensitivity to

what we believe to be a variation in oxygen pressure. In Figure 1 we show

the surface tension calculated from a number of sessile drops in three

different series of experiments. Although the atmosphere in all the experiments

was nominally ultra pure argon, somewhat different conditions in each group

of measurements produced different concentrations of impurities. Curve 1,

the lowest data, reflects the effects of some atmospheric contamination of

the flowing argon introduced through a small leak in the gas feed line

10
between the gettering furnace and the tube containing the silicon. In

these experiments the development of a yellow deposit on the quartz tube

near the molybdenum susceptor and downstream from the argon flow indicated

oxygen in the system. This deposit is most likely M0O 3 formed by reaction

with the susceptor and then sublimated^. The intermediate surface tension

data (curve II) was obtained in early experiments in which the zirconia

insulation was in contact with the molybdenum and may have been reduced

somewhat. Furthermore, the system was not heated during the preliminary

outgassing procedure. A much lighter deposit than in the previously described

experiments was formed on the tube under these conditions.

The highest data shown in Figure 1 we believe to be free of any significant depression due to oxygen in the system. The measurements were preceded

by a rigorous vacuum degassing at nearly 1400°C. There was no detectable

deposit formation during the measurement and the solidified silicon samples

had bright, smooth mi rror like surfaces. A linear regression analysis of

these data yields a K value of -0.279 mJ/m^K. The probable error in K is

.004 mJ/m^K. The K values found for the other groups of data are shown

on Fi gure 1

.

Previous measurements of the surface tension of liquid silicon are

listed in Table 1. These results scatter over the range of surface tensions

observed in our work, probably because of variations in oxygen levels. Our

high data lie slightly above that of Kingery and Humenik^O obtained in

pure hydrogen which supports our belief that the concentration of oxygen is

very low in our measurements. We suggest, therefore, the following expression

for the surface tension of clean silicon as a function of temperature:

Y

= 885 - 0.28(T - 1410)

11
As mentioned earlier, there has been one previous measurement of the

surface tension of silicon over a temperature ranged, Lukin et al found

the surface tension could be represented as

Y

= 727 - .1044(1 - 1410).

Their measurements were made using the sessile drop technique in purified

helium atmospheres. Some of the difference between their K value and ours

can be attributed to the use of different density data in the calculation

of the surface tension. Lukin et al determined the density from the volume

of the sessile drops. These measurements found a much smaller decrease in

density with increasing temperature than the generally accepted dependence

of Lucas^l which we have used. If their surface tension values are adjusted

to this latter density scale the surface tension temperature coefficient is

decreased to -0.15 mJ/m^K. Thus they find both a surface tension and |K|

significantly lower than the results of this work even after the density

correcti on.

The lowering of the surface tension by adsorption of a surface active

component is often accompanied by a reduction in the magnitude of K. This

is the case with silicon over the range of impurity concentrations encountered

in these measurements. In Figure 2 we plot |K| as a function of y at

the melting point for the data of Figure 1 and the adjusted result of

Lukin. The four points indicate an approximately linear variation of K

with Y for the conditions of these experiments. Although the adjustment

of the Lukin data into accord with our results may be fortuitous, their low

surface tensions are very likely due to surface active impurities.

In conclusion, we have measured the surface tension of liquid silicon

in purified argon atmospheres as a functicn of temperature. We find a

12
temperature coefficient near -0.28 mJ/m^K. Our experiments show a high

sensitivity of the surface tension to what we believe are low concentrations

of oxygen. Thus we cannot rule out some effect of low levels of oxygen in

our results. However, the highest surface tension values obtained in

conditions which minimized the residual oxygen pressure are in good agreement

with a previous measurement in pure hydrogen. We therefore believe that

depression of the surface tension by oxygen is insignificant in these

measurements

.

13
References

1. W. E. Langlois, J. Crystal Growth 56 (1982) 15.

2. S. V. Lukin, V. I. Zhuchkov and N. A. Vatolin, J. of the Less-Common

Metals 67 (1979) 399.

3. F. Bashforth and J. C. Adams, Capillary Action , University Press,

Cambridge, 1883.

4. T. O'Donnell, M. Leipold and M. Hagan, Compatibility Studies of

Various Refractory Materials in Contact with Molten Silicon,

DOE/JPL-1012-771 6, 1978.

5. S. Dushman, Scientific Foundations of Vacuum Technology (John Wiley

& Sons, Inc., New York, 1962) 646.

6. J. V. Najdich and V. N. Eremenko, Fiz. Metal, i Metalloved II (1961)

883.

7. D. W. G. White, Trans. ASM 55 (1962) 757.

8. C, A. Smolders and E. M. Duyvis, Recueil 80 (1961) 635.

9. Yu. M. Shashkov and T. P. Kolesnikova, Russ, J. of Phys. Chem. 37

(1963) 747.

10. W. D. Kingery and M. Humenik, J. Phys. Chem. 57 (1953) 359.

11. L. D. Lucas, Mem. Sci . Rev. Met., 61 (1964) 1.

12. P. H. Keck and W. Van Horn, Phys. Rev. 91 (1953) 51 2.

13. V. P. Eljutin, V. I. Kostikov and V. Ya. Levin, Izv. Vys. Uch. Sav.,

Tsvetn Met. 2 (1970) 131 .

14. N. K. Dshemilev, S. I. Popel and B. V. Zarevski, Fiz. Met. i Met. 18

(1964) 83.

15. P. V. Geld and M. S. Petrushevski , Izv. A. N., OTN, 3 (1961 ) 160.

14
Table 1

O

—1 O (mJ/m2)

atm. Ref.

1410 712 Vac. 9

1410 720 Argon 12

1500 825 Argon 13

1550 750 Vac. 14

1550 720 Vac. 15

1450 740 He. 10

1450 860 H2 10

1450 874 Argon Thi s

Table 1. Summary of previous silicon surface tension measurements.

15
900

700*—

1400

l_

1500 1600

Temperature (®C)

Figure 1. The surface tension of silicon as a function of temperature.

The surface tension and K variations in the three groups of

data are attributed to different concentrations of impurities.

See text.

16
700 800 900

Surface Tension (mJ/m^)

Figure 2. The absolute value of K as a function of the surface tension

at the melting temperature.

17
io she. ^ edT ,s
Task 2

Convection During Unidirectional Solidification

S. R. Coriell and R. J. Schaefer

Metallurgy Division

Center for Materials Science

and

G. B. McFadden

Mathematical Analysis Division

Center for Applied Mathematics

SUMMARY

During solidification of a binary alloy at constant velocity vertically

upwards, thermosolutal convection can occur if the solute rejection or

preferential i ncorporation at the crystal-melt interface decreases the density

of the melt. This convection can occur even in a microgravity environment,

though, since it is basically a gravity driven flow it will be much more

pronounced in earth-bound experiments than in space flight. A number of

numerical calculations of the solute, temperature, and flow fields have

been carried out using previously developed algorithms for Schmidt numbers

Sc = 1 , 10, 81. The time periodic flow, which occurs over a narrow range

of parameters, has been further investigated.

Experimental observations of the interaction of flow fields with the

crystal -melt interface during the unidirectional upward solidification of

succinonitri le containing ethanol have shown how the interface shape correlates

to different fluid flow patterns. Asymmetries in the thermal environment

are found to play an important role, which is dependent upon the crystal

growth velocity. In addition, the presence of a vapor bubble at the solidliquid interface is sometimes found to induce vigorous convection.

The general problem of the interaction of flow fields with crystal

melt interfaces has been investigated by linear stability calculations for

19
the case of a flow parallel to a crystal-melt interface. This research, in

collaboration with R. F. Sekerka of Carnegi e-Mel Ion University and M. E.

Glicksman and colleagues at Rensselaer Polytechnic Institute, has demonstrated

a coupling between hydrodynamic and morphological instabilities. The

calculations are in general agreement with experiments carried out by

Glicksman and colleagues on succi nonitrile.

20
Introduction

The general aim of this task is the theoretical and experimental

study of the fluid flow and solute segregation which occur during

solidification. The nature of the fluid flow, its effect on the shape of

the crystal-melt interface and on the resulting distribution of solutes

is examined. In particular, the role of solutal convection (gravity

driven flow due to solute gradients) during the vertical directional

solidification of a binary alloy is considered both theoretically and

experimental ly

.

We present a discussion of the numerical calculations first and then

a description of the experiments: each section can be read independently.

The results of a number of new calculations of thermosolutal convection

are presented. In particular it appears that the time periodic behavior

previously found is related to the second eigenmode in linear theory. A

few results for the lead-tin system are also given.

A number of manuscripts have been submitted for publication during the

year, viz., (1) Thermosolutal Convection During Directional Solidification,

6. B. McFadden, R. G. Rehm, S. R. Coriell, W. Chuck, and K. A. Morrish,

Met. Trans.; (2) Asymmetric Instabilities in Buoyancy-Driven Flow in a

Tall Vertical Annulus, G. B. McFadden, S. R. Coriell, R. F. Boisvert, and

M. E. Glicksman, Phys. Fluids; (3) Morphological Stability in the Presence

of Fluid Flow in the Melt, G. B. McFadden, S. R. Coriell, R. F. Boisvert,

M. E. Glicksman, and Q. T. Fang, Met. Trans.; (4) Effect of a Forced

Couette Flow on Coupled Convective and Morphological Instabilities During

Unidirectional Solidification, S. R. Coriell, G. B. McFadden, R. F.

Boisvert, and R. F. Sekerka, J, Crystal Growth; (5) Convection and

21
Interface Instability During Alloy Solidification, S. R. Coriell, G. B.

McFadden, R. G. Rehm, and M. E. Glicksman, Proceedings of Workshop on

Solidification and Fluid Dynamics in the Earth's and Space Laboratory,

Aachen; (6) Convection-Induced Distortion of a Solid-Liquid Interface,

R. J. Schaefer and S. R. Coriell, Met. Trans. The results in (1) and (5)

are based on previous NBS reports; the results in (2), (3), (4), and (6)

are new and these manuscripts are attached in Appendices (A-D).

The manuscript Coupled Convective Instabilities at Crystal -Melt

Interfaces, S. R. Coriell, G. B. McFadden, R. F. Boisvert, M. E. Glicksman,

and Q. T. Fang (Appendix in a previous report^) has been accepted for

publication in the Journal of Crystal Growth.

Numerical Results

The basic equations and methods for calculating the fluid velocity,

solute, and temperature fields during directional solidification of a

binary alloy have been previously reported^ *2. The notation used is

defined in Tables I and II.

Since it appears that the oscillations in time previously^ reported

are correlated with the second eigenmode of linear theory, we have carried

out linear stability calculations for Schmidt numbers Sc = 1 , 10, and 81.

Typically semiconductors have Sc = 10 and metals Sc = 100. The results

are shown in Figs. 1-3. For Sc = 1 , the first three eigenmodes are shown;

the dashed curve indicates how the wavelength of the fastest growing mode

varies with the degree of instability. Figures 2-3 are similar although

only the first two eigenmodes are shown. Previous nonlinear calculations

for Sc = 10 showed time periodic behavior for a narrow range of wavelengths

and values of Rs/Rs*. Comparison of Fig. 2 with Fig. 4 of reference 2

indicates that the non-linear time periodic behavior occurs in the vicinity

22
of the second eigenmode of linear theory.

We have carried out a number of calculations for the lead-tin system

for a growth velocity of 2(10"^) cm/s and a temperature gradient in the

liquid of 200 K/cm. The dimensionless parameters are then kL/kg = k/<s = 0.535,

k = 0.3, Sc = 81, Pr = 0.0225, Ra = 43.48, Rs* = 82.79. We report the

results of four computations:

(Pb 1) Rs/Rs* =1.51, U = 0.0741 , H/x = 1 .303, MESH (12, 54),

(Pb 2) Rs/Rs* = 2.13, U = 0.0741 , H/X = 1 .600, MESH(12, 54),

(Pb 3) Rs/Rs* = 3.86, U = 0.123, H/x = 3.223, MESH(12, 90),

(Pb 4) Rs/Rs* = 4.28, U = 0.123, H/x = 1 .790, MESH(12, 90).

The values given as MESH(N, M) indicate the number of intervals (N, M) in

the X and y direction, respectively; the calculation is carried out over

half a wavelength in the x-direction. In the first three calculations,

the relationship between wavelength x and degree of instability Rs/Rs*

was chosen to approximate the fastest growing mode of instability from

linear theory (see Fig. 3), and steady state behavior was found. For

calculations Pb 1, Pb 2, and Pb 3, the relative solute segregation at the

crystal melt -interface (the ratio of maximum minus minimum to average

concentration) is 0.378, 0.319, and 0.264 while the ratio of the maximum

flow velocity and the crystal growth velocity is 3.3, 4.7, and 8.5,

respecti vely

.

The wavelength and degree of instability Rs/Rs* for calculation Pb

4 were chosen in the region where time periodic flow might occur based on

our previous results for Sc = 10 and the location of the second eigenmode.

Oscillatory behavior did occur and the stream function and concentration

field at four different times during one period of oscillation is shown

23
in Figs. 4-5, respectively. The period of oscillation corresponds to about

880 s; in this time the crystal has grown 0.18 cm or 1.2 (D/V). Because

of the large Schmidt number, calculations for the lead-tin system require

very large computation times, e.g., about 3 hours of CPU time (Uni vac

1100/82) per dimensionless time unit for the oscillatory behavior shown

in Figs. 4-5, where the dimensionless time exceeds 35.

We have also initiated some calculations for Sc = 1 , where the

behavior appears qualitatively similar to the larger Schmidt number cases

(compare Figs. 1-3). The shorter computer times allow a more detailed

study of the variation of the flow fields with wavelength and solutal

Rayleigh number. One of the most interesting results is the dependence

of the steady state on the initial state (non-uniqueness of the steady

state). Figs. 6 and 7 show two different steady states with identical

values of the dimensionless variables, but with different initial

conditions. The relationship between non-uniqueness and oscillatory

behavior will be further investigated.

Several calculations have been carried out for reduced but constant

gravitational accelerations. We plan to modify the algorithm to permit

calculations for variation of the direction and magnitude of the gravitational acceleration with time. In some circumstances it appears that

growth vertically downward would be preferable. Linear and non-linear

computations for this case will be undertaken.

Experimental Results

Previously reported experiments^ demonstrated that convective flow

patterns present during the unidirectional growth of succinonitri 1

e

containing ethanol could result in the formation of a macroscopic pit in

24
the solid-liquid interface, with a cellular mi crostructure then developing

within the pit. The convective flow in these experiments was attributed

to the effect of radial teirperature gradients in the vicinity of the

solid-liquid interface. These gradients were due to heat loss to the

environment and to the generation of latent heat of fusion at the solidliquid interface when solidification occurred. It was demonstrated that

the activation of an auxiliary heater surrounding the solid-liquid

interface could be effective in supressing the formation of the interface

pit.

Additional measurements of the fluid flow above the solid-liquid

interface have indicated the importance of the asymmetry of the temperature

field. When heat losses are symmetrical, the fluid flow shows a toroidal

pattern; downward along the walls of the sample, radially inward above

the interface, and upward along the sample axis (Fig. 8a). When the

sample is held stationary in this condition, marker particles tend to

accumulate at the interface just above the point where the radial inward

flow converges, and to flow upward in a concentrated stream along the

axis of the torus (Fig. 9). Asymmetry of the thermal field causes the

interface to become tilted and the axis of the torus to be shifted to the

downhill side (Fig. 8b). Additional asymmetry causes the torus to pinch

off on one side and the flow then forms a single transverse roll (Fig. 8c).

The asymmetric part of the heat losses becomes relatively less

important when crystal growth is initiated, as a result of the large

symmetric generation of heat of fusion. In response to this effect, the

axis of the toroidal flow moves toward the center of the tube. Activation

of the auxiliary heater surrounding the interface generally resulted in

25
a more asymmetric flow field either because the auxiliary heater itself

was not symmetric or because the auxiliary heater cancelled the symmetric

part of the heat losses but left the asymmetric part. Thus the elimination

of the central interface depression when the auxiliary heater was activated

can be attributed in part to the decrease in radial gradients but also in

part to the change in the form of the flow field. Complete elimination

of horizontal gradients (and subsequent elimination of thermal convection)

by the use of auxiliary heaters is unlikely to be successful because of

the difficulty of producing perfectly uniform transfer of heat to the

specimen tube.

A very strong convective flow has been observed in the vicinity of a

vapor bubble which was trapped at the solid-liquid interface. This bubble

was formed when the sample was melted back into a region where the solid

contained a shrinkage void. The flow pattern above this bubble was

sometimes similar to that in samples containing no bubble but at other

times it showed a localized counterrotating toroidal pattern directly

above the bubble (Fig. 10). The velocities in this counterrotating flow

were at least an order of magnitude faster than those normally present,

and were fastest directly adjacent to the bubble surface. The counterrotating

flow appeared to be most vigorous shortly after the bubble shifted its

position. These observations indicate that the flow is driven by surface

energy gradients along the liquid-vapor interface; bubble motions would

cause sudden transient melting or solidification, thus increasing local

temperature and concentration gradients and inducing more rapid flow.

The effect would be especially significant under low gravity conditions,

where the location of vapor bubbles is difficult to control. More

26
observations of this effect are planned, to clarify the relative roles of

temperature and solute gradients.

27
References

R. J. Schaefer, S. R. Coriell, R. 6. Rehm, and G. B,

NBSIR 82-2560 (ed. by J. R. Manning) p. 25 (1982).

R. J. Schaefer, S. R. Coriell, G. B. McFadden, and F NBSIR 83-2772 (ed. by J. R. Manning) p. 19 (1983).

McFadden,

. G. Rehm,

28
Table I

Definition of Symbols Used for Physical (Xiantities

V (kinematic viscosity)

< (thermal diffusivity of melt)

<s (thermal diffusivity of crystal)

D (diffusion coefficient of melt)

g (acceleration of gravity)

a (thermal coefficient of expansion of melt)

etc (solutal coefficient of expansion of melt)

k|_ (thermal conductivity of melt)

kg (thermal conductivity of crystal)

k (distribution coefficient)

Coo (bulk concentration of solute)

Coo* (bulk concentration corresponding to onset of convection)

G|_ (temperature gradient in melt)

H (melt height)

V (crystal growth velocity)

X (horizontal wavelength)

t (time)

29
Table II

Dimensionless Quantities Used in Calculations

Ratio of Thermal Conductivities kL/kg

Ratio of Thermal Diffusivities k/ks

Distribution Coefficient k = 0.3

Schmidt number Sc = v/D

Prandtl number Pr = v/<; Pr$ = v/<s

Thermal Rayleigh number Ra = (goGL/vic) (D/V)^

Critical solutal Rayleigh number Rs* = (ga^C^vD) (D/V)^

Dimensionless crystal growth velocity U = VH/v

Solutal Rayleigh number Rs = (ga^Coo/vD) (D/V)^

Ratio of melt height to horizontal wavelength H/x

30
c/c«,

0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0

J I I I I I I I I

a>

Fig. 1. The normalized concentration at the onset of instability during

directional solidification as a function of the dimensionless

wavelength of a sinusoidal perturbation. Three eigenmodes and

the maximum growth rate are shown. For these curves (kL/k^) = {k/<s) = 1, k = 0.3, Sc = 1, Pr = 0.01, Ra = 1268, Rs* = 961.

31

UNEAK

THEORY

FOR

DA

=15

cm
XV/D

Fig. 2. The normalized concentration at the onset of instability during

directional Solidification as a function of the dimensionless

wavelength of a sinusoidal perturbation. Two eigenmodes and the

maximum growth rate are shown. For these curves (ki_/ks) = (</<s) = 1»

k = 0.3, Sc = 10, Pr = 0.01, Ra = 1268, Rs* = 963.

32

LINEAR

THEORY

FOR

DA

=

15

cm

Sc

=

10.

Pr=

.01
XV/D

C/Coo

0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0

I

I I I

I \ I

I I

Fig. 3. The normalized concentration at the onset of instability during

directional solidification of lead containing tin as a function

of the dimensionless wavelength of a sinusoidal perturbation.

Two eigenmodes and the maximum growth rate are shown. For these

curves (kL/ks) = (</<s) = 0.535, k = 0.3, Sc = 81 , Pr = 0.0225,

Ra = 43.48, Rs* = 82.79.

33

THEORY

FOR

LEAD-TIN

SYSTEM
RS= .35390+003 Rfl= .43480+002 SC= 81.0000 P= .0225 LV/D= 5.5866

STREAM FUNCTI0N C0NT0UR MflP AT T = 3S. 1342 STREAM FUNCTI0N CaNTOJR NAP AT T = 3S.3707

STREAM FUNCTION CONTOUR MAP AT T = 35.6664 STREAM FUNCTION CONTOUR NAP AT T = 35.9029

Fig. 4. The stream function as a function of position at four dimensionless

times for the directional solidification of lead containing 0.0015 wt%

tin at a growth velocity of 2(10"^) cm/s and a temperature gradient

in the liquid of 200 K/cm. The streamfuncton is periodic in time; one dimensionless time unit corresponds to 926 s.

34
RS= .35390+003 Rfl= .43480+002 SC= 81.0000 P= .0225 LV/D= 5.5866

C0NCENTRflT10N C0NT0UR MflP fiT T = 35.1342

C0NCENTRflT10N CflNT0UR Nff* RT T = 35.9029

Fig. 5. The concentration field as a function of position for the same conditions as Fig. 4.

35
RS= .33438+004 Rfl= .12679+004 SC= 1.0000 P= .0100 LV/D= 2.4644

ONTIUR

Fig

STREAM FUNCTION CONTOUR MAP AT T = 4.4741

.M999-001

6. The steady state stream function as a function of position; see

Fig. 1 caption for values of the dimensionless parameters not

given above (here L = x).

36
RS= .33438-*-004 Rfl 12679+004 SC= 1.0000 P 0100 LV/D= 2.4644

QWTIUR FIOH

STREAM FUNCTION CONTOUR MAP AT T = 5.6549

.37«S-001

Fig. 7. The steady state stream function as a function of position for

the same parameters as Fig. 6. There are two different steady

state stream functions (non-uniqueness) depending on the initial

conditions.

37
ig. 8. Convective flow patterns above the solid-liquid interface with

(a) symmetric, (b) slightly asymmetric, and (c) highly asymmetric

heat losses.

38
Fig. 9. Stream of particles along the axis of the toroidal flow field.

Fig. 10. Vapor bubble at the solid-liquid interface, with associated

rapid counterrotating convective flow.

39

APPENDIX A

ASYMMETRIC INSTABILITIES IN BUOYANCY-DRIVEN FLOW

IN A TALL VERTICAL ANNULUS

G.B. McFadden, S.R. Coriell, and R.F. Boisvert

National Bureau of Standards, Washington, D.C. 20234

and

M.E. Glicksman

Materials Engineering Department, Rensselaer Polytechnic

Institute, Troy, New York 12181

ABSTRACT

Linear stability of the one-dimensional flow between infinite vertical

coaxial cylinders induced by heating the inner cylinder is considered for

various ratios k of the inner radius to the outer radius, and for Prandtl

numbers P appropriate to air and water. For air with P=0.71 the least stable

disturbance is non- axi symmetric for k<0.44 and is axi symmetric for <>0.44,

and in either case the instability is due to the action of the shear forces.

For P=3.5, the situation is similar, except that the asymmetric shear mode

is superceded by an axi symmetric instability driven by buoyancy forces for

0.03<<<0.16. Wave speeds, v/ave numbers, and critical Grashof numbers for

these cases and for the case of zero Prandtl number are gi ven.

41
INTRODUCTION

A basic problem in the study of heat transfer concerns the stability of the

flov/ of a liquid placed between two vertical surfaces that are held at

different temperatures^. Due to the variation of the fluid density with

temperature, buoyancy forces typically cause the fluid to rise near the

hotter surface and sink near the colder surface. The geometry considered in

this note is a pair of infinite vertical coaxial cylinders with the inner

cylinder maintained at a higher temperature than that of the outer cylinder.

If the temperature difference is sufficiently small, the resulting flow is

purely vertical and the lateral heat transport is due solely to conduction.

If the temperature difference is large enough, however, the parallel flow

becomes unstable and a more complex flow develops; heat is then convected

1 aterally.

2 3 Several authors » have considered the linear stability of the basic parallel

flow to disturbances that are axially symmetric. The purpose of this note is

to extend these results to include asymmetric perturbations. The results

indicate that for small ratios ic of the inner radius to the outer radius, and

for small Prandtl numbers, the parallel flow first becomes unstable to

asymmetric disturbances.

Asymmetric (helical) instabilities of a different type have also been

observed experimentally and predicted by linear theory*^ for the situation in

which the outer cylinder is crystalline succinonitri le and the fluid is

molten succi nonitri le with a Prandtl number of 23. In that case the feature

of interest is the interaction between the non-rigid crystal-melt interface

and the buoyancy-dri ven flow.

42
FORMULATION AND SOLUTION

The governing equations are taken to be the Navier-Stokes equations in the

Oberbeck-Boussi nesq approximation^. The equations are made dimensionless by

scaling length by the radial gap L (outer radius minus inner radius), time

by L^/v, where v is the kinematic viscosity, temperature by the temperature

difference AT between the cylinders, fluid velocity by U = g a L aT/v, where

g is the magnitude of the gravitational acceleration which acts in the

negative z direction and a is the coefficient of thermal expansion, and

pressure by pU , where pis the density. The unperturbed temperature field

T^Cr) and the vertical base flow velocity Wjj(r) depend only on the radial

coordinate and are given by Choi and Korpela^. The perturbed quantities are

assumed to have the form f(r) exp(iwz + in<|» + at), and linear instability is

predicted when Re(a)>0 for any wavenumbers w and n. If the r-dependent

part of the perturbed radial, azimuthal, and axial velocity components, the

pressure, and the temperature are denoted by the complex quantities u,v,w,p,

and T, respectively, the linearized governing equations reduce to the

following set of ordinary differential equations:

Du + u/r + (inv)/r + iww = 0,

(a + ia)GWQ)u + GDp = D^u + Du/r - + (l+n^)/r^)u - 2inv/r^,

(a + ioiGW^jv + inGp/r = D^v + Dv/r - + (l+n^)/r^)v + 2inu/r^,

(a + io)GVfQ)w + G(DWq)u + iwGp = T + D^w + Dv//r - + n^/r^)w,

P(a + i(oGW )T + PG(DT )u = D^j + DT/r - + n2/r2^T

n n

/

' ;

I .

43
Here the Grashof number G=g a AT L^/v^; the Prandtl number P=v/k; D=d/dr;

and k is the thermal diffusivity. The perturbed velocity components and

temperature field all vanish at the inner cylinder r=K/(l-i<) and at the outer

cyl i nder r=l/(l-K).

The equations can be recast as a set of 16 real first order ordinary

differential equations. They are solved using the linear boundary value

problem software SUPORT^ in combination with the nonlinear equation solver

SNSQE » . Following a method described by Keller , the complex boundary

condition u(l/(l-ic))=0 is replaced by Du(1/(1-k) )=1 in order to obtain a

nonsingular boundary value problem. The nonlinear solver then varies G and

Im(a) for fixed values of Re(a), w, n, k, and P until the computed value of

u(1/(1>k:)) is driven to zero. Values of G corresponding to marginal

stability are obtained by setting Re(a)=0 for fixed u, n, k, and P. For given

values of P and k, the critical Grashof number corresponds to the smallest

marginal value of G taken over all values of w and n. For n=0 this method

reproduces to within 1% the critical Grashof numbers obtained by Choi and

Korpela using a different method.

RESULTS

As discussed by Hart^°, energy can be transferred from the main flow to the

disturbance by the action of either the buoyancy forces or the shear forces

(Reynolds stress). Eigenmodes corresponding to each mechanism can be

identified. The Grashof number for marginal stability of the buoyant mode

increases markedly as the Prandtl number is decreased, whereas the shear mode

is relatively insensitive to the Prandtl number. The mechanism for

instability in effect at the critical Grashof number depends on the specific

value of P as well as k.

44
Figure 1 shows marginal Grashof numbers vs ic for air with Prandtl number

P=0.71. Two modes are shown: the dashed curve corresponds to azimuthal

wavenumber n=0 and the solid curve corresponds to n=±l. Each point on the

curves represents the minimum value of G with respect to w with the other

parameters fixed. The parallel flow is stable for Grashof numbers below both

curves, and the critical Grashof number for any < is given by the smaller

value of the two curves.

2

The dashed curve is the symmetric shear mode discussed by Choi and Korpela

and Shaaban and Ozisik . The solid curve is also a shear mode, and the

geometry requires that the two modes become indistinguishable in the limit as

< approaches 1, which represents a planar vertical slot. The asymmetric mode

lies slightly above the symmetric mode for 0.44<ic<l and the critical Grashof

number in this range corresponds to a symmetric disturbance, as was observed

experimentally for <=0.68 by Choi and Korpela^, For <<0.44 the critical

Grashof number corresponds to an asymmetric distortion of the flow. The

symmetric buoyant mode lies well above these curves at this Prandtl number

and is not shown.

Selected values of the critical Grashof number, axial wave speed, and axial

wavenumber are given in Table I. Positive wave speeds correspond to

upward-moving waves. The disturbance is stationary for <=1. Pxesults for P=0

are also included in Table I. In this limit the temperature field has no

effect on the flow and the problem reduces to the classical Orr-Sommerfeld

equati ons.

Figure 2 shows results for P=3.5. (The Prandtl number of water has a value

of 3.5 at about 50C but the thermophysical properties vary considerably with

45
temperature; the validity of the Boussinesq approximation for water therefore

restricts the magnitude of the temperature difference AT). For this value

the symmetric buoyant mode (dotted curve) has dropped sufficiently to become

dominant. The two shear modes bear roughly the same relation to each other

as before, but now the buoyant mode is the least stable mode over the range

0.03<<<0.16. Results for this case are also included in Table I. The buoyant

mode is increasingly destabilized as the Prandtl number increases, so that

for much higher values of P the critical Grashof numbers correspond to

symmetric disturbances.

In summary, for low Prandtl numbers and for low ratios of the inner cylinder

radius to the outer cylinder radius, the parallel flow first becomes unstable

to asymmetric perturbations. This effect should be observable in experiments

with air for radius ratios <<0.4.

46
ACKNOWLEDGEMENT

This work was conducted under the support of the Microgravity Sciences

Program, National Aeronautics and Space Administration.

47
REFERENCES

^G.Z. Gershuni, and E.M. Zhukhovitski i , Convective Stability of

Incompressible Fluids (Keter, Jerusalem, 1976),

^I.G. Choi, and S.A. Korpela, J. Fluid Mech. 99, 725 (1980).

^A.H, Shaaban, and M.N. Ozisik, 7th Inti. Conf. of Heat Transfer, 1982,

Munich, Germany.

S.R. Coriell, R.F. Boisvert, J.I. Mickalonis, and M.E. Glicksman, Advances

in Space Research (Proc. of 24th COSPAR Meeting) 3,95 (1983).

5

S. Chandrasekhar, Hydrodynamic and Hydromagnetic Stability (Oxford

University Press, London, 1961).

^M.R. Scott and H.A. Watts, SIAM J. Numerical Anal. 14, 40 (1977).

7

' SLATEC Common Math Library, National Energy Software Center, Argonne

National Lab., Argonne, IL, written by K. L. Hiebert, based on Powell®.

0

M.J.D. Powell, in: Numerical Methods for Nonlinear Algebraic Equations ,

ed. P. Rabinowitz, Gordon and Breach, NY (1970).

g H.B. Keller, Numerical Solution of Two Point Boundary Value Problems ,

Regional Conference Series in Applied Mathematics 24, SIAM, Philadelphia

(1976).

Hart, J. Fluid Mech. 35, 547 (1971).

48
TABLE I

k=.99 <=.60 o

CM .

II <=.10 <=.01

P=0.0 n=0 (Sh) n=0 (Sh) n=l (Sh) n=l (Sh) n=l (Sh)

G 7930.2 8319.4 10258. 10856. 14384.

0) 2.69 2.65 2.25 2.13 1.98

c .24088 12.568 29.529 30.712 33.113

P=.71 n=0 (Sh) n=0 (Sh) n=l (Sh) n=l (Sh) n=l (Sh)

G 8038.0 8512.4 10435. 10438. 10026.

0) 2.80 2.75 2.24 2.11 2.09

C .25432 13.400 32.316 32.320 24.770

P=3.5 n=0 (Sh) n=0 (Sh) n=l (Sh) n=0 (B) n=l (Sh)

G 7857.1 8347.5 10582. 9712.6 14767.

b) 2.75 2.75 2.25 1.21 1.90

C .24677 12.977 32.006 87.554 36.433

Table I. The critical Grashof number G, dimensionless axial wavenumber u, and

dimensionless axial wave speed c=|Im(a)|/aj (in units of v/L) as functions of

the radius ratio k. '(Sh)' denotes a shear instablity and '(B)' denotes a

buoyant instability.

49
FIGURE CAPTIONS

Figure 1. Marginal Grashof numbers vs radius ratio < for Prandtl number

P=.71. The dashed curve is the symmetric shear mode and the solid curve is

the asymmetric shear mode. The symmetric mode lies slightly below the

asymmetric mode for <>0,44.

Figure 2, Marginal Grashof numbers vs radius ratio < (logarithmic scale) for

Prandtl number P=3.5. The dashed curve is the symmetric shear mode and the

solid curve is the asymmetric shear mode. The dotted curve is the symmetric

buoyant mode. The symmetric shear mode lies slightly below the asymmetric

mode for large <.

50
P=.71

Figure 1

51
P=3.5

K

Figure 2

52
APPENDIX B

MORPHOLOGICAL STABILITY IN THE PRESENCE OF

FLUID FLOW IN THE MELT

G.B. McFadden, S.R. Coriell, and R.F. Boisvert

National Bureau of Standards, Washington, DC 20234

and

M.E. Glicksman and Q.T. Fang

Materials Engineering Department

Rensselaer Polytechnic Institute

Troy, New York 12181

ABSTRACT

Recent experiments have shown that the presence of a vertical buoyancy-driven

flow adjacent to an initially cylindrical crystal-melt interface may produce

a time-dependent helical deformation of the interface, with a rotation period

ranging from several minutes to many hours, depending upon the v/idth of the

melt. The temperature distribution is such that the interface is expected to

be morphologi cal ly stable in the absence of fluid flow. A linear stability

analysis reveals that the instability is due to a coupling betv/een a basic

hydrodynamic instability in the buoyant flow and the deformable boundary

separating the two phases. The crystal-melt interface lowers the critical

Grasliof number of an analogous rigid-v/alled system by an order of magnitude

for succi noni tri le v/ith a Prandtl number P=22.8; furthermore, the hydrodynamic mode that is actually destabilized by the interface is not the least

stable mode in the rigid-walled system for P=22.8. The results show that the

instability may be regarded either as a rather large alteration of a basic

hydrodynarni c i nstabi li ty by the crystal-melt interface, or as a significant

modification of the morphologi cal stability of the interface by the presence

of the buoyant flow.

53
INTRODUCTION

The occurrence of fluid flow in the melt during solidification has received

considerable attention, see e.g. the reviews by Hurle^, Carruthers^, and

Pimputkar and Ostrach^. It is well known that convection can have a significant effect on the structure and properties of the resulting solid. In

particular, the morphological stability^”^ of the crystal-melt interface may

be altered considerably by the presence of fluid flow. For example. Delves

has shown^ that a forced flow parallel to the interface helps stabilize the

interface against a Mul li ns-Sekerka instability^, V'/hile Coriell, et al.^

show that thermosol utal instabilities in the melt may interact with the

morphological instability in a complicated fashion and may in fact result in

segregation at solute levels much lower than those predicted on the basis of

% morphological stability theory alone.

A surprising effect of convection on an otherwise stable sol id- liquid

interface has recently been investigated experimentally by Fang, Glicksman

and Mi ckaloni They studied a pure sample of succi noni tri le (SCN)

with a vertical, stationary cylindrical crystal-melt interface in a radial

symmetric temperature field that increases into the liquid. Although

according to morphological stability theory the cylindrical interface should

be stable in the absence of convection, they found that the buoyancy-driven

flow induced by the temperature gradient in the melt v/ill under certain

conditions produce a slowly rotating helical interface, accompanied by a more

complex time-dependent flow. The critical Grashof number for the onset of

this instability^! is an order of magnitude lower than the critical Grashof

numbers associated with hydrodynamic instabilities in ri gi d-v/al led

contai ners!-»!-^.

54
Reference 14 contains a detailed description of the experiment and a comparison of experimental results vnth theoretical predictions of a linear stability analysis. A more detailed theoretical study for the related case of a

planar geometry has also been performed^^, shov\fi ng the dependence of the

critical Grashof number on the Prandtl number of the melt. In this paper the

experimental results are reviewed briefly, and some new theoretical results

concerning the nature of the crystal-melt instability are presented. It is

shown that the instability may be viewed either as a large alteration of a

hydrodynamic instability by the solidification process, or as a significant

modification to morphological stability theory due to the buoyancy-induced

f 1 ow

.

EXPERIMENT

% The experiment consists of observing the flow state and internal

configuration of a long glass cylinder containing solid succi noni tri le with

a platinum heating wire running along the axis of the cylinder. By passing a

current through the wire while cooling the outside of the container in a

isothermal bath of the proper temperature, a steady state cylindrical

crystal-melt interface may be produced which separates the outer annulus of

crystal from the core of melt near the wire, c.f. Figure la. The hotter

fluid near the wire is less dense than the fluid near the interface, and a

steady vertical flow is established such that the fluid flows upward near the

wire and downward near the interface. If the wire temperature is slowly

increased, the magnitude of the flow velocity increases until a critical melt

gap and a critical temperature difference AT between the wire and the

crystal-melt interface are reached for which the cylindrical interface is no

longer stable. The crystal-melt interface is then observed to assume a

55
helical shape, c.f. Figure lb.

The experimental results are summarized as follows:

1) A non-axi symmetri c mode was observed when the Grashof number exceeds

a critical value. The critical values range from 140 to 180,

depending on the average radial gap width L of the melt at the onset

of instability. (The width L at the onset may be varied by

adjusting the bath temperature and the wire temperature properly.)

2) The longitudinal wave speed associated with the crystal-melt

interface instability is about two orders of magnitude smaller than

the maximum base flow velocity, and was found to scale as L"'^. The

' wave length, on the other hand, was found to be proportional to L.

3) The helical interfacial waves were observed to rotate with periods

ranging from a few minutes to about ten hours, depending on the melt

gap. These rotation periods were found to scale as

4) The helical wave associated with the interface instability can

rotate with either a right-or a left-handed sense. It was found

that there is a definite preference for left-handed helices for

smaller value of L, whereas right-handed helices are preferred for

large values of L. In addition, the helical waves were always found

to travel upwards against gravity independent of their sense of

rotation. Experiments also shov/ed that the sense of rotation may be

influenced by rotating the specimen tube. For further details of

the experiment, see reference 10 and 14.

56
The remainder of this paper consists of a collection of new theoretical

results obtained by linear stability analysis. A brief description of the

governing equations and the method of numerical solution is also provided.

NUMERICAL METHOD

The experiment is modelled theoretically by assuming that the aspect ratio of

the container is infinite. The governing equations are taken to be the

Boussinesq approximati on^^ to the Navi er- Stokes equations for the velocity

field, a convection-diffusion equation for the temperature in the melt, and a

diffusion equation for the temperature in the crystal. The unperturbed

problem then has a simple solution depending only on the radius, with the

temperature field and vertical velocity field as given by Choi and Korpela^^.

The equations are made dimensionless by choosing the length scale to be the

radial distance L between the unperturbed interface and the wire, the time

scale to be L^/v, where v is the kinematic viscosity of the melt, and the

temperature scale to be aT, the temperature difference betv/een the wire and

the unperturbed crystal-melt interface. The velocity scale is taken to be

v/L, v/hich differs from the scaling used in earlier papers^^»l^>^5. For the

purposes of a linear stability analysis the equations and boundary conditions

are linearized about the unperturbed velocity and temperature fields and

about the cylindrical interface. The z,4>,and t dependence in cylindrical

coordinates of the perturbed quantities is assumed to have the form exp(ot +

iojz + in({)), where w is the axial wavenumber, n is the azimuthal wavenumber,

and the complex time constant a determines the temporal stability of the

system; instability is predicted when ar=Re(a)>0. The resulting equations

are

57
Du + u/r + (inv)/r + iww = 0

(a + io)Wo)u + Dp = D^u + Du/r - + (l+n^)/r^)u - 2i nv/r^

(a + iwWo)v + inp/r = D^v + Dv/r - + (l+n^)/r2)v + 2i nu/r^

(o + i(jjVjQ)v/ + (DWo)u + iwp = D^w + Dw/r - + n^/r^)w + GT

P(a + i(oWo)T + P(DTq)u = D^T + DT/r - (o)2 + n^/r^)!

i n the melt, and

PsaTs = d2Ts+ DTs/r - (co^ + n2/r2)Ts

in the crystal, vvhere D=d/dr, Wo(r) is the unperturbed vertical velocity, and

To(r) is the unperturbed temperature in the liquid. The r-dependent part of

the perturbed radial, azimuthal, and a>tial velocity components, the pressure

field, and the liquid and solid temperature fields are denoted by u(r), v(r),

v/(r), p(r), T(r), and T 5 (r), respectively. The Prandtl number P=v/K, where K

is the thermal diffusivity in the liquid, Ps=v/Ks, where Ks is the thermal

diffusivity in the solid, and the Grashof number G=gaATL^/v2, where a is the

coefficient of volume expansion. With this choice of units, Wo(r) is proportional to G.

The ratio of the wire radius to the iinperturbed interface radius is denoted

by K. The boundary conditions at the wire radius r=</(l-<) are u=v=w=T=0.

The solidification boundary conditions at the interface r=l/(l-K) are^

u = 0

V = 0

58
w = -6 DWq

T = -6DTo

Ts = -6(kL/ks)DTo

-o6H = (kL/ks)DT - DTj

where kL and are the liquid and solid thermal conductivities, and

H=vLv/(ksAT) , where Ly is the latent heat released per unit volume upon

solidification. The perturbation of the interface has the form 5exp(at + iwz

+ in(})). It is also convenient to define A=GH. These boundary conditions

have been simplified by assuming that the dimensionless capillarity constant

Y=GTn^r/(ATL)=0, where T^^ is the melting point of SCN and r is the ratio of

the 5 olid-liquid surface tension to Ly> and the dimensionless density change

upon solidification e=(p 5 /p|_)-l= 0 , where ps and pL are the densities of the

solid and liquid, respectively. This is an excellent approximation for the

experiment with SCN. The crystal extends to r=Ls> where the perturbed

temperature vanishes.

The equations can be rewritten as a system of sixteen real first-order ordinary differential equations in the liquid. The solution procedure is similar

to that employed in previous linear stability analyses^. Briefly, the eigenvalue problem is solved using the linear boundary value problem software

SUPORT^^ in conjunction with the nonlinear equation solver SNSQE^^>^^ following a method outlined by Keller^O.

It is also of interest to consider the case of a planar geometrylS^ ^an

be described by the above equations with r formally taken to be infinitely

large and with the azimuthal velocity set to zero. The crystal is assumed to

59
extend to infinity. The unperturbed vertical flow field induced by buoyancy

has a simple cubic profile and the unperturbed temperature field in the

liquid is linear. Two methods were used for the planar case: the procedure

outlined above and another method based on discretizing the equations using

second-order accurate finite differences on a uniformly-spaced mesh to obtain

a matrix eigenvalue problem. The latter method is useful for obtaining

initial guesses for the nonlinear root finder in the first procedure as well

as for perform! ng quick scans in parameter space to identify various

eigenmodes. The eigenvalues were computed using the subroutine F02GJF from

the NAG library^l, which is based on the QZ algorithm^^^

NUMERICAL RESULTS

A coijiparison of theory and experiment^'*^ shows that the critical Grashof

number, wavelength, and phase velocity of the instability are predicted

rather well by linear stability theory. (However, the theory is unable to

account for the observed tendency for left-handed helices to be favored for

small liquid gaps, even when Coriolis force terms are included in the

differential equations!^). In this section the nature of the interfacial

instability and its relationship to other unstable modes in the system are

considered.

Marginal stability curves for a radius ratio of k=0.02 are shov/n in Figure 2.

The parameters used in computing these curves are P=Ps=22.8, A=10170.0,

kf/ks^l.O, £-0.0283, y=4. 53(10"^) , and Ls=1.0. The numerical results

indicate that setting e=Y=0 has a negligible effect on the marginal stability

curves, and the rest of the results will all correspond to the case e=y=0.

The cylindrical interface is stable for Grashof numbers less than 176, and is

60
first unstable to an asymmetric, helical disturbance (n=l). For sufficiently

large Grashof numbers, the system is actually unstable to several different

types of disturbance, and five of these modes are indicated in Figure 2. The

curves labelled "n=0" represent axi symmetric disturbances, and those labelled

"n=l" represent asymmetric disturbances. For each of the curves labelled

"n=l" there is also an identical curve with n=-l which corresponds to an

eigenmode that differs only in its symmetry; that is, n=l represents a

left-handed helical disturbance and n=-l represents a right-handed

disturbance. Linear combinations of n=l and n=-l modes which are not helical

are also possible.

The upper three curves correspond to instabilities that would also occur if

the crystal-melt interface were replaced by a rigid isothermal wal

% The top two curves correspond to hydrodynamic instabilities driven by energy

transferred from the main flow to the disturbance by the action of Reynolds

stress^^ ("shear modes"), which represents momentum flux due to convection.

The third curve (n=0) corresponds to instabilities due to energy transfer by

the buoyancy forces ("buoyant modes"). For large k there are both n=0 and

n=l buoyant modes, but for the n=l case as < is decreased the unstable region

inside the closed marginal stability curves ap=0 (cf. Figure 3) becomes

smaller and smaller. Eventually op becomes negative everywhere and the

marginal stability curves vanish altogether.

The lower two curves in Figure 2 represent the new modes that appear when a

crystal-melt interface that can change shape is present rather than a rigid

isothermal wall. The critical Grashof number G=176 occurs at a v/avenumber

w=1.35. The disturbance takes the form of a travelling wave that moves in

the upward direction, with oi=Im(a)=-0.977(10"2) . (Mote that the unperturbed

61
flow is directed downward near the interface). This value of a-j corresponds

to a phase velocity that is two orders of magnitude lower than the maximum

flow velocity. By contrast, the wavespeed for the buoyant mode is comparable

to the base flow velocity. Values of G, w, and a-j at the minimum of each of

the five curves are given in Table I.

The value of k= 0.02 used in Figure 2 is representative of the values studied

by Fang, et al.^^. A detailed numerical examination of the dependence of

the crystal-melt instability on < shows that for <>0.8 the least stable mode

is axi symmetri c; that is, the relative positions of the lowest n=l and n=0

modes shown in Figure 2 are reversed. A qualitatively similar dependence on

< is also obtained for the top tv/o shear modes in the case of rigid

boundari es^^.

The marginal stablity curves in Figure- 2 give values of G and w for which

op=0. More generally, a grov/th rate is determined by each value of G and w,

with different growth rates a^(G,w) associated with each distinct mode of

instability. Surpri si ngly , the crystal-melt instability and the shear

instability v/ith the same value of n share the same growth rate ar(G,o)); that

is, the growth rate varies continuously along a path connecting the

marginal stability curves for the crystal-melt mode and the shear mode. This

result implies that the crystal -melt mode is in fact a branch of the ordinary

hydrodynamic shear mode. The five marginal stability curves in Figure 2

arise from just three distinct relationships ap(G,tu) correspondi ng to the n=0

and n=l shear modes and the n=0 buoyant mode. The crystal-melt interface

apparently destabilizes the shear modes and allows instabilities to occur at

si gni fi cantly lower Grashof numbers.

62
This relationship can be observed more easily if the Prandtl number is

decreased, which apparently stabilizes the crystal-melt mode while leaving

the shear mode relatively unchanged^^. Figure 4 shows marginal stability

curves for P=0.148, <=0.02, and n=l for both a crystal-melt interface (solid

curve) and a rigid isothermal wall (dashed curve). The other parameters have

the same values as for Figure 2, except that e and y have been set to zero.

The curves are similar for u>l, while for io<l the development of the

crystal-melt mode is clear. Values of ai for both cases are also included in

the bottom half of the plot. Note the large variation in ai near w=l,

indicating the relatively slow wave speed of the crystal-melt mode. As the

Prandtl number increases from P=0.148 to P=22.8 the marginal stability curve

of the crystal -melt mode drops down and to the right while the interior

maximum near u=l rises sharply.

In this context the crystal-melt mode may be viewed as a modification of a

hydrodynamic instablity by the effect of a phase transition at one boundary.

On the other hand, an equally valid i nterpretati on is that the

gravity-induced flow in the melt alters the morphological stabilty of the

soli d-li quid interface. Vlith the temperature increasing from the interface

into the liquid, morphological stability theory predicts that the cylindrical

sol id- liquid interface is stable to small disturbances in the absence of

fluid flow. The theory also predicts that in other circumstances a

sol id- liquid interface may become unstable if the liquid is supercool ed^^ or

if enough solute is added to the liquid^ (constitutional supercooling). The

crystal-melt mode described here is yet another instance of interfacial

instability, in this case due to a buoyancy-driven flow.

63
The curve labelled G=0 in Figure 5 shows the growth rate or<0 of asymmetric

disturbances to the cylindrical crystal-melt interface in the absence of

fluid flow. For all wavenumbers w the disturbances are damped out in time.

This curve represents the n=l version of the mode of interfacial deformation

studied in classical morphological stability theory. The other curves show

how this mode is destabilized by the presence of the gravi ty-i nduced flow

parallel to the interface, with G=176 again correspondi ng to the critical

Grashof number for which the interface is first unstable. In this figure the

temperature difference between the interface and the outer cylinder is the

same for all curves, and the variation in Grashof number may be considered to

be caused by different values of the gravitational acceleration g. The point

to be emphasized here is that it is the same eigenmode as is examined in

morphological stability theory that is excited by the applied flow field.

For this calculation P=Ps=22.8,Ls=!<L/ks~l‘^> H=57.8.

The crystal-melt mode is thus due to the interaction of the non-rigid

interface with the buoyant flow parallel to the interface. The dependence of

this mode on the specific form of the flow field is considered next. For

simplicity a planar geometry will be considered. The planar case also

exhibits a crystal-melt mode, a buoyant mode, and a shear mode^^, with a

critical Grashof number G=251, wavenumber u=0.8, and time constant

ai=-3. 47(10-2) for P=Ps=22.8, and A=1.017(104)

.

A natural question from the metallurgist's point of view is whether the

iriterfacial instability depends only on the fluid flow in the immediate

vicinity of the interface or whether it depends instead on details of tlie

flow field farther out from the interface. For example, do similar instabilites occur when the cubic velocity profile of the buoyancy-dri ven flow is

64
replaced by simple Couette or Poiseuille flows having linear and quadratic

velocity profiles, respectively? For H=16.95, P=22.8, Ps=0, and o)=0.7, in

the absence of all fluid flow the growth rate of the morphological mode is

or=-7. 44(10"^). With buoyant flow the crystal-melt instability occurs for a

critical Grashof number that produces a velocity shear dw/dr=20.17 at the

interface. Setting G=0 and adding either enough Couette or enough Poiseuille

flow to produce the same shear dw/dr=20.17 at the interface actually stabilizes this mode, with growth rates ap=-1.29(10"l) and a|-=-9,05(10"2) , respectively. The stabilizing effect is consistent with the results of Delves^

showing that a forced flow has a stabilizing influence on a constitutionally

supercooled, advancing interface. The matrix stability code was also

employed to search for instabilities over some other parameter ranges, but no

instances of interfacial instability were observed for either Couette or

Poiseuille flow.

A related question is whether the purely hydrodynamic instability that is

known to occur for planar Poiseuille flow^*^ at a Reynolds number of 5772 is

affected by the presence of a crystal -melt interface. This Reynolds number

is based on the maximum flow velocity and the half-width of the liquid and

corresponds to a much larger shear rate dv//dr=4.62(10'^) at the interface.

The effect seems to be negligible: with H=16.95 and P=Ps=22.8, the critical

Reynolds number is altered by only 0.003%. Figure 6 shows a plot of vs. w

for both the hydrodynamic mode (solid curve) and the morphological mode

(dashed curve) for a Reynolds number R=5800. In contrast to the buoyant flow

case (Figure 5) for which the shear mode coincides with the morphological

mode, these two modes remain distinct, with the morphological mode stable and

the hydrodynamic mode just unstable at this Reynolds number.

65
A more conservative approach to the question of the dependence of the

crystal-melt mode on the velocity profile is taken in Figure 7, where a

Couette flow is superposed on a buoyancy-driven flow v/ith a Grashof number of

G=600. The Couette flow is generated by moving the hot wall opposite the

interface with a dimensionless velocity V (again based on v/L). With no

additional Couette flow the planar interface is unstable for 0.28<cu<0.93 at

G=600, and this unstable range of wavenumbers corresponds to the interval

inside the closed solid curve with V=0. The solid curve is the marginal

stability curve with V and w considered as parameters, and the planar interface is unstable outside the curve and is unstable inside it. The figure

indicates that the interface becomes stable if enough Couette flow is

present. For G=600 the net flow is nowhere downward for V>50 and is nowhere

upward for V<-6.25, so that the figure indicates that the crystal-melt

instablity can occur even when the flow is unidirectional (-11.23<V<-6.25)

.

The stabilization for large |V| is consistent with the argument that the

existence of the crystal -melt instability hinges upon the presence of the

shear instability, since in the limit of pure Couette motion it is known that

the flow field is stable to all infinitesimal disturbances^^. The exact

details of the interaction are not straightforward, however, since, for example, as V decreases from zero the range of unstable values of to decreases,

then increases, and finally vanishes. Values for a-j are included as the

dashed curve in Figure 7. The Couette flow alters the wavespeed considerably;

note, however, that remains negative, meaning that the travelling wave

still moves upv/ard.

Tt is perhaps worth emphasizing that the existence of the crystal-melt mode

does not seem to be a simple consequence of the cubic velocity profile. The

66
instability exhibits a strong dependence on Prandtl number^^, whereas both

the unperturbed velocity field and the unperturbed temperature field are

independent of Prandtl number. For small Prandtl numbers, in which case the

perturbed temperature becomes negigible, the crystal-melt mode is found to

revert to the hydrodynamic shear instability that is present for rigid

boundaries. For large Prandtl numbers the critical Grashof number of the

crystal-melt mode is inversely proportional to the Prandtl number, so that

the Rayleigh number given by the product PG approaches a constant for large

P. Some simplification of the governing equations can be acheived for large

Prandtl number, but the analysis of the resulting system still appears to be

compl i cated.

SUMMARY

The observed helical instability results from a coupling between the

hydrodynamic shear instability associated with the buoyant flow and the

morphological mode of interface deformation that arises in ordinary

morphological stability theory. The morphological mode apparently does not

couple to the buoyant mode, which is actually less stable than is the shear

mode for SCN, nor does it couple to the hydrodynamic instability associated

with a Poiseuille flow. In one example, the addition of enough Couette flow

stabilizes the crystal-melt interface, although it is possible to retain the

helical instability while adding enough Couette flow to ensure that the flow

is everywhere unidirectional. The observed tendency for left-handed helices

to occur at small liquid gaps is not explained by the theory. A convincing

physical argument that explains the coupling between the shear mode and the

interface is also lacking.

67
ACKNOWLEDGEMENT

The authors gratefully acknowledge the support of the Microgravity Sciences

Program, National Aeronautics and Space Administration. We also thank Drs.

J. I. Mickalonis, S. Ostrach, R. 6. Rehm, R. J. Schaefer, R. F. Sekerka, and

J. A. Simmons for helpful discussions.

68
REFERENCES

^D. T. Hurle, in: Current Topics in Materials Science, 1977, vol 2, p. 549,

ed. by E. Kaldis and H. J. Scheel, North-Holland , Amsterdam.

2j. R. Carruthers, in: Preparation and Properties of Solid State Materials,

1977, vol. 3, p. 1, ed. by W. R. Wilcox and R. A. Lefever, Dekker, New

York.

^S. M. Pimputkar and S* Ostrach: J. Crystal Growth, 1981, vol. 55, p. 614.

%

^J. S. Langer: Rev. Mod. Phys., 1980, vol. 52, p. 1.

^R. T. Delves, in: Crystal Growth, 1974, vol. 1, p. 40., ed. by B. R.

Pamplin, Pergamon Press, Oxford.

^W. W. Mullins and R. F. Sekerka: J. Appl. Phys., 1964, vol. 35, p. 444.

^S. R. Coriell, M. R. Cordes, W. J. Boettinger, and R. F. Sekerka: J.

Crystal Growth, 1980, vol. 49, p. 13.

®J. I. Mickalonis: M. S. Thesis , Materials Engineering Dept., Rensselaer

Polytechnic Institute, 1982.

69
E. Glicksman and J. I. Mickalonls: Proc. of the 16th Southestern

Seminar on Thermal Sci., 1982, vol, 1, p. 505, Miami, FL.

T. Fang: M. S. Thesis , Materials Engineering Dept., Rensselaer

Polytechnic Institute, 1983.

R. Coriell, R. F. Boisvert, J. I. Mickalonis, and M. E. Glicksman:

Advances in Space Research (Proc. of 24th COSPAR Meeting), 1983, vol. 3,

no. 5, p. 95.

G. Choi and S. A. Korpella: J. Fluid Mech. , 1980, vol. 99, p. 725.

^^A. H. Shaaban and M. N. Ozisik; 7th Inti. Conf. of Heat Transfer, 1982,

Munich, Germany.

T. Fang, M. E. Glicksman, S. R. Coriell, G. B. McFadden, and R. F.

Boisvert: submitted to J. Fluid Mech.

^^S. R. Coriell, G. B. McFadden, R. F. Boisvert, M. E. Glicksman, and Q. T.

Fang: to appear in J. Crystal Growth.

^^S. Chandrasekhar: Hydrodynamic and Hydromagnetlc Stability , Oxford

University Press, London, 1961.

^^M. R. Scott and H. A. Watts: SIAM J. Numerical Anal., 1977, vol. 14,

p.40.

70
^®SLATEC Common Math Library, National Energy Software Center, Argonne

National Lab., Argonne, IL, written by K. L. Hiebert, based on Powell^®.

I

J. D. Powell, in: Numerical Methods for Nonlinear Algebraic Equations ,

ed. P. Rabinowitz, Gordon and Breach, IfSf, 1970.

B. Keller, Numerical Solution of Two Point Boundary Value Problems ,

Regional Conference Series in Applied Mathematics, 1976, vol. 24, SIAM,

Philadelphi, PA.

NAG Fortran Library Manual: Mark 9, vol. 4, Numerical Algorithms Group,

Ltd., Oxford, 1982.

%

22c. B. Moler and G. W. Stewart: SIAII J. Numer. Anal., 1973, vol. 10, p.

241.

2^g. B. McFadden, S. R. Coriell, R. F. Boisvert, and M. E. Glicksman,

submitted to Phys. Fluids.

2^P. G. Drazin and W. H. Reid: Hydrodynamic Stability , Cambridge University

Press, Cambridge, 1981.

2^W. W. Mullins and R. F. Sekerka: J. Appl. Phys., 1963, vol 34, p. 323.

71
TABLE I

G u n

175.72 1.35 -9.7735(10-3) i

464.14 1.45 -3.6457(10-2) 0

2151.8 2.30 -43.729 0

13216.3 2.00 -69.049 1

60308.3 2.10 -708.84 0

Table I. Approximate minimum Grashof numbers G, axial wavenumbers «, time

constant and azimuthal wave number n for the five neutral stability

curves shown in Figure 2.

72
FIGURE CAPTIONS

Figure 1. Side view of the stable cylindrical mode (a), and helical mode (b).

The black spot on the right hand side is the profile of a thermistor that

detects temperature variations in the crystal.

Figure 2, The Grashof number at the onset of instability as a function of the

axial wavenumber w of a sinusoidal perturbation. The lowest two modes occur

only for a crystal- melt interface, while the higher three modes occur for

both a crystal-melt interface and a rigid interface. The curves are labelled

with the azimuthal wavenumber n.

Figure 3. The Grashof number at the onset of instability as a function of the

axial wavenumber m of a sinusoidal perturbation for the n=l buoyant mode with

rigid boundaries. The solid curve corresponds to a radius ratio ic=0.325, and

the dashed curve to <=0.292.

Figure 4. The Grashof number G (top curves) and the time constant (bottom

curves) at the onset of instability as a function of the axial wavenumber w of a asymmetric perturbation for Prandtl number P=0.148, radius ratio

<=0.148. The solid curves correspond to crystal-melt boundary conditions,

and the dashed carves to rigid wall boundary conditions.

Figure 5. The growth rate Oj- as a function of axial wavenumber w for

asymmetric perturbations for G=0, 90, 176, and 350. The lowest curve with

G=0 corresponds to the nK>rphological mode in the Mullins-Sekerka theory (in

this case stable).

Figure 6. The growth rate Oj- as a function of the axial wavenumber to for

disturbances in a planar geometry with Reynolds number 5800. The unstable

hydrodynamic mode (solid curve) is distinct from the stable morphological

mode (dashed curve).

Figure 7. The Couette flow velocity V and the time constant -a^ as a function

of the vertical wavenumber to at the onset of instability in a planar

geometry. The Grashof number G is fixed at G=600. The region inside the

closed solid curve corresponds to unstable values, the region outside

corresponds to stable values.

73
cm
CYLINDRICRL GEOMETRY

CRY5TRL-MELT INTERFACE

P=22.8 /C=0.02
CY

5000 n

4000 -1

3000 -J

2000 -

1000 -

0 4

- 0

INDRICBL GEOMETRY

RIGID BOUNDRRY

R=22.8 n=l

/C='.325

/C=.292

0.5 1 1.5 2 2.5 3 3.5

00

Figure 3
CYLINDRICRL GEOMETRY

P=.H8 /C=.02 n=l

CRYSTRL-MELT INTERFRCE

RIGID BOUNDRRY

00

Figure 4

77
CYLINDRICRL GEOMETRY

CRYSTflL-nELT INTERFACE

P-22.8 /C-0.02 n-1

Figure 5

78
PLflNRR GEOMETRY

POISEUILLE FLOW R-5800

HYDRODYNAMIC MODE

MORPHOLOGICAL MODE

Figure 6

79
PLRNRR GEOMETRY

CRYSTRL-MELT INTERERCE

P=22.8 G=600

0.4

0.2

0.0

b

1

80
APPENDIX C

Effect of a Forced Couette Flow on Coupled Convective and

Morphological Instabilities During Unidirectional Solidification

S. R. Coriell, G. B. McFadden, R. F. Boisvert

National Bureau of Standards

Washington, D.C. 20234

and

R. F. Sekerka*

Carnegi e-Mel Ion University

Pittsburgh, PA 15213

*Consultant, National Bureau of Standards
ABSTRACT

The effect of a forced Couette flow, parallel to a horizontal crystalmelt interface during directional solidification of an alloy of lead

containing tin, on the onset of convective and morphological instabilities

is calculated numerically via a linear stability analysis. Such a flow

does not affect perturbations with wave vectors perpendicular to the

flow. The onset of morphological instability is somewhat suppressed and

thermosolutal convection is greatly suppressed. When instabilities occur,

they are oscillatory and correspond to travelling waves. For values of

the crystal growth velocity for which mixed morphological and convective

modes occur, the presence of a forced flow produces sufficient decoupling

to allow formerly degenerate branches to be identified.

82
1. INTRODUCTION

The onset of natural convection and morphological instability during

vertical directional solidification of a binary alloy has been calculated

previously by linear stability analysis [1-6]. For given growth conditions,

it is well known that the crystal -melt interface is morphologically

unstable in the presence of sufficiently large solute concentrations [7].

For growth vertically upwards toward a hotter fluid whose density decreases

with increasing temperature, and in the absence of horizontal gradients,

the fluid is convectively unstable (thermosolutal convection) for solute

concentrations above some critical value provided that a solute which is

rejected at the crystal -melt interface is less dense than the bulk melt

or a solute that is preferentially incorporated at the crystal-melt

interface is more dense than the bulk melt. It is well known that such

thermosolutal convection can occur even if the unperturbed fluid density

decreases with height [8] because of the different rates at which heat

and solute diffuse. For the solidification of lead containing tin with a

temperature gradient in the melt of 200 K/cm, it was found [1] that the

mode of instability at onset corresponds to thermosolutal convection for

growth velocities less than about 40 un/s and to morphological instability

above this velocity. In the vicinity of this transition velocity, there

are mixed modes that involve complex interaction between convective and

morphological instabilities.

Delves [9-11] has calculated the effect on morphological stability

of a forced flow parallel to the unperturbed crystal -melt interface. Such

a flow has no effect on interface shape perturbations with wave vector

perpendicular to the direction of the flow, but usually increases the

83
critical concentration for morphological instability, and gives rise to

traveling waves, for perturbations with wave vector parallel to the

direction of the flow. Similarly, it is known that a Couette shear flow

stabilizes perturbations with wave vector parallel to the direction of

flow for Ray lei gh-Benard convection, i.e., a pure liquid heated from

below T12-14] and for thermosolutal convection with linear solute

and temperature fields [15].

In this article we report numerical results for the effect of a

forced fluid velocity, parallel to the crystal-melt interface and varying

linearly with distance from the interface (Couette flow), on the onset of

convective and morphological instabilities during the directional solidification of lead containing tin. Since perturbations with wave vector

perpendicular to the flow direction are unaffected by the flow, the

results obtained here are for perturbations with wave vector parallel to

the direction of flow. Our results extend those of Delves to include

thermosolutal convection and extend our previous work on thermosolutal

convection to include a forced flow. When a forced flow is included, the

instability is oscillatory and corresponds to traveling waves. In general,

we find that the forced flow stabilizes perturbations with wave vector

parallel to the flow in the sense that the critical concentration for

their instability is larger than in the absence of such flow* However,

unless perturbations with wave vector perpendicular to the flow are

absent or can be otherwise suppressed, the onset of instability will

occur at the same concentration as in the absence of forced flow but will

correspond to two-dimensional rolls with the axes of the rolls parallel

to the flow direction.

84
2. THEORY

We consider unidirectional solidification at constant velocity V in

the vertical direction of a binary alloy with an initially planar interface.

We choose a moving coordinate system (x, y, z) attached to this planar

crystal-melt interface located at z = 0 with the z-axis directed into the

melt and antiparallel to the gravitational acceleration. The basic

equations and boundary conditions governing the convective and morphological

stability of this system have been derived previously [1,3] and will

simply be stated here. The fluid velocity U', the temperature T', the

concentration C, and the crystal-melt interface shape z$l be expressed

in the form

ii' = Uo(z) F(x, y, t),

T* = To(z) + T(z) F(x, y, t),

C = Co(z) + C(z) F(x, y, t),

z$L = <5 ^{x. y, t),

where the subscript o indicates the unperturbed fields which are functions

of z alone, 5 is the initial amplitude of the perturbed interface, U, T and C

express the dependence of the perturbed fluid velocity, temperature, and

concentration, on z and F(x, y, t) = exp [at + i(ojxx + ayy)]. The complex

number a governs the dependence on time t of the perturbed fields, and

ux and ay are spatial frequencies. The system is unstable if the

real part of a is positive for any values of ojx and ay. The

unperturbed temperature and concentration fields depend exponentially on

z. The unperturbed flow field is of the formJUo(z) =^q(z)

with e = (ps/pl) - where ps and pL are densities of crystal

and melt, respectively; we measure fluid velocities in the laboratory

85
frame. The function q represents a forced flow in the x-direction

parallel to the unperturbed crystal-melt interface; we shall later take

q to be linear in z in our calculations.

The linearized differential equations that govern stability are

{L[vL + V*(d/dz) - a] - iwx(qL - d^q/dz^) } W = guj2(aT + ct^C),

{kcL + V*(d/dz) - a - i(Uxq}T = (dTo/dz)W,

{DL + V*{d/dz) - a - i(oxq}C = (dCo/dz)W, (2)

where W is the z-conponent of U, v is the kinematic viscosity, V* = (1 + e)V,

L = (d^/dz^ - u?) ^ u? = (uj2 + (u^), g is the magnitude of the gravitational

acceleration, a and ac are coefficients used to express the

change in fluid density p from a reference density pq, viz., p = pqCI - aT' - ocC'],

< is the thermal diffusivity and D is the diffusion coefficient. The

unperturbed temperature and concentration gradients are

(dTo/dz) = Gi_ exp (-V*z/k)

and

(dCo/dz) = Gc exp (-V*z/D),

where Gl and G^- are the temperature and concentration gradients at the

unperturbed crystal -melt interface z = 0 .

The linearized boundary conditions at z = 0 are [1, 3]

W = -eo6,

(dW/dz) = [eoj^ + iojj^(dq/dz)] 6 ,

ksujsT - kL(dT/dz) = 6 {oLv - (V*k|_GL/K) + (VksGs/Kj) - ' ^s)

VC + (D/[(l + e)(l - k)])(dC/dz) = - 6 {[VGck/(k - 1)] +Cooa/k>.

mC - T = 6 {G|^ - mG^ + T^ruj^}, (3)

where kL and kg are the thermal conductivities of the melt and crystal,

respectively, <5 is the thermal diffusivity of the crystal.

86
u:^ = -(V/ 2 <

3

) + [(V/ 2 <

5 )^ + 4- a/<s]^'^^, Ly is the latent heat per unit

volume of the crystal, G5 = (VLy + ki_GL)/k 5

, k is the distribution coefficient,

Coo is the solute concentration far from the interface, Gc = V*Cco(k - l)/(Dk),

m is the liquidus slope, is the melting point of the pure crystal with

a planar interface, and r is a capillary constant, i.e., the ratio of the

crystal-melt surface tension and Ly. The initial interface amplitude 6

can be eliminated from eq. (3) to obtain four boundary conditions involving

only W, T and C.

The boundary conditions far from the interface at z = zl are

W = (dW/dz) = (dT/dz) = C = 0. (4)

As long as zl is sufficiently large, the numerical results with (dq/dz)

constant are independent of zl and independent of whether Dirichlet or

Neumann boundary conditions are specified for the perturbed temperature

and concentration fields.

It is apparent from the differential equations and boundary conditions

that q and its derivatives appear only when multiplied by Therefore,

defining e such that wx = oj cos e and = cj sin 0 , the

solution depends only on w and q cos 0 . Without any loss of

generality, we take 0 = 0 in our calculations since results for other

values of 0 can be obtained by replacing q by q cos 0 . Clearly,

if 0 = tt/ 2, i.e., ujx = 0, the imposed flow field has no

effect on the critical concentration for the onset of instability.

The numerical methods have been described previously [1]. The

differential equations are rewritten as a set of 16 real first order

equations. They are solved by using the software SUPORT [16]

developed for linear boundary value problems. We set op. the real

87
part of a, equal to zero and iterate on c„ and a], the imaginary part

of a, using the non-linear equation solver SNSQE [17-18], until we find

a solution of the differential equations that satisfies the boundary

conditions. In addition to the tests of the code described previously

[1], we have reproduced results of Gersting and Jankowski [19] for the

Orr-Sommerfeld equation and of Gallagher and Mercer [14] for RayleighBenard convection with shear.

3. RESULTS

Numerical calculations have been carried out for lead containing tin for

crystal growth velocities of 2, 40, and 200 un/s with Gl = 200K/cm, ay = 0,

and e = 0. The properties of the lead-tin system used in the calculations

are given in Table 1 of reference 1. We take q(z) = Sz, where the constant

S determines the amount of shear in the forced flow field.

We first give results for a relatively low growth velocity (V = 2 un/s)

for which the mode of instability at onset corresponds to thermosol utal

convection [1]. Then we turn to a relatively high velocity (V = 200 un/s)

for which the mode is morphological. Finally, we consider an intermediate

velocity at which mixed modes prevail.

For V = 2 un/s, the mode of instability at onset is thermosolutal

convection which occurs at a bulk concentration c„ = 3.615(10"^) wt.% for

S = 0, i.e., for the absence of forced convection. In the absence of all

convection, morphological instability would not occur until c« = 3.899 wt.%.

We characterize the amount of shear, S, of the forced flow field by the

effective velocity SD/V which is the value of q at a distance D/V into

the liquid. The effects of shear on the critical concentration, wavelength, 2tt/oj, and wave speed, -ai/uj, at the onset of instability are

88
shown in Figures 1-3, respectively. The critical concentration and wave

speed increase with increasing shear. The wavelength decreases, reaches

a minimum, and then slowly increases as the shear increases. For

SD/V = 0.087 cm/s, the critical concentration is about 100 times that for

no forced flow. The wave speed is comparable to SD/V and is about 350

times the crystal growth velocity for SD/V = 0.087 cm/s.

For V = 200 uT)/s, the mode of instability is morphological and

occurs at a bulk concentration of 0.0554 wt.% for S = 0. For this

case, an overstable convective mode [5] does not occur until c« = 3.7 wt.%.

The critical concentration, wavelength, and wave speed as a function of SD/V

are shown in Figures 4-6, respectively. As SD/V increases, the critical

concentration and wave speed increase while the wavelength decreases.

Extremely large shear is required to increase the critical concentration,

e.g., S = 10^ s"^ (which corresponds to SD/V = 1.5 cm/s) increases the

critical concentration by about 50%. For SD/V = 1.5 cm/s, the wave speed

of 6.5(1 0"3) cm/s is about one third the crystal growth velocity. The

magnitude of the shear required to increase the critical concentration

for morphological stability is in general agreement with the calculations

of Delves [11] for gallium-doped germanium.

At a crystal growth velocity of 40 u^i/s, the onset of convective

and morphological instability occur at about the same concentration, and

the mode of instability (see Fig. 1 of reference 1 and Figs. 2-4 of

reference 3) is mixed. The effect of shear on this instability is

shown in Figures 7-9 for S = 10‘2, 10"^, and 1 s“^, respectively. The

critical concentration and ai are plotted as functions of the spatial

frequency u. The onset of instability occurs at the minimum value of

89
Coo as a function of cj . We identify the minimum at larger values of

tu with morphological instability (dashed curves) and the minimum at

smaller values of w with convective instability (solid curves). For

S = 0.01 s"l, the curves are similar to the S = 0 curves (Fig. 3 of

reference 3). For finite S, ai is generally non-zero, while for

S = 0, oi vanishes identically for a range of w values. As S increases,

the convective and morphological modes appear to decouple. For S = 1 s”l

,

the concentration curves appear to cross each other with little interaction; however, ai of the morphological mode changes sign near the

crossing point. In Fig. 7, the minimum of the convective branch

occurs at a smaller value of c« than in Fig. 9; the minima in the

morphological branches occur at essentially the same value of c^.

DISCUSSION

Extension of our previous work on coupled convective and morphological

instabilities to include a forced flow parallel to the interface affects

perturbations with wave vectors along the forced flow in the following

ways: (1) the onset of morphological instability is somewhat suppressed

and the onset of thermosolutal convection is greatly suppressed; (2). the

resulting instabilities are oscillatory and correspond to travelling waves;

(3) for values of the growth velocity for which mixed morphological and

convective modes occur, the presence of a forced flow produces sufficient

decoupling to allow formerly degenerate branches to be identified.

Roughly speaking, a forced flow tends to smooth inhomogeneities in the

temperature, concentration, and flow fields, along the direction of the

flow, that would otherwise form spontaneously and lead to instabilities.

90
It should be emphasized that the above results were calculated for the

case of a very simple flow velocity q= Sz (Couette flow) which is known

to be hydrodynami cal ly stable to infinitesimal disturbances for all

values of S [20]. Other flows that can be hydrodynami cally unstable,

such as those recently investigated [21], can lead to strong destabilizing

interactions between convective and morphological modes, even though the

flow field near the interface itself could be approximated by Couette

flow. Flows with shear rates comparable to those found to be significant

in the present paper have been produced during eutectic solidification by

rotation of the crystal [22]. Other flows related to horizontal temperature

gradients have been calculated by Chang and Brown [23] and shear rates at

the interface can be calculated from such flows. These flows are, however,

much more complex than the simple Couette flow treated here and great

caution should be excercised in the application of our results to such

situations.

Acknowledgements

This work was partially sponsored by the Microgravity Sciences

Program of the National Aeronautics and Space Administration. One of us

(R. F. Sekerka) received partial support from the National Science

Foundation under grant DMR78-22462. We would like to thank R. A. Brown,

R. G. Rehm and R. J. Schaefer for helpful discussions.

91
REFERENCES

1. S. R. Coriel 1, M. R. Cordes, W. J. Boettinger, and R. F. Sekerka,

J. Crystal Growth 49 (1980) 13.

2. S. R. Coriel 1, M, R. Cordes, W. J. Boettinger, and R. F. Sekerka,

Adv. Space Res 1 (1981) 5.

3. S. R. Coriel 1 and R. F. Sekerka, PhysicoChetnical Hydrodynamics 2

(1981 ) 281 .

4. R. J. Schaefer and S. R. Coriell, in: Materials Processing in

the Reduced Gravity Environment of Space, Ed, G. E. Rindone

(North-Hol land, Amsterdam, 1982) p. 479.

5. D. T. J. Hurle, E. Jakeman, and A. A. Wheeler, J. Crystal Growth 58

(1982) 163.

6. 0. T. J. Hurle, E. Jakeman, and A. A. Wheeler, Phys. Fluids 26

(1983) 624.

7. W. W. Mullins and R. F. Sekerka, J. Appl. Phys 35 (1964) 444.

8. J. S. Turner, Buoyancy Effects in Fluids (Cambridge Univ. Press,

London, 1973).

9. R. T. Delves, J. Crystal Growth 3/4 (1968) 562.

10. R. T. Delves, J. Crystal Growth 8 (1971) 13.

11. R. T. Delves, in: Crystal Growth, Vol. 1, Ed. B. R. Pamplin

(Pergamon, Oxford, 1974) p. 40.

12. G. Z. Gershuni and E. M. Zhukhovitskii , Convective Stability

of Incompressible Fluids (Keter, Jerusalem, 1976).

13. J. W. Deardorff, Phys. Fluids 8 (1965) 1027.

14. A. P. Gallagher and A. McD. Mercer, Proc. Roy. Soc. A286 (1965) 117.

15. P. F. Linden, Geophys. Fluid Dyn. ^ (1974) 1 .

16. M. R. Scott and H. A. Watts, SIAM J. Numerical Anal. 14 (1977) 40.

17. SLATEC Common Math Library, National Energy Software Center,

Argonne National Lab., Argonne, IL. The program SNSQE was written

by K. L. Hiebert and is based on an algorithm of Powell.

92
1ft. M. J. D. Powell, in: Numerical Methods for Nonlinear Algebraic

Equations, Ed. P. Rabinowitz (Gordon and Breach, New York 1970) p. ft7.

19. J. M. Gersting, Jr. and D . R. Jankowski, Int. J. Num. Methods Eng.

4 (1975) 195.

?0. P. G. Drazin and W. H. Reid, Hydrodynamic Stability (Cambridge

Univ. Press, London, 1981).

21. S. R. Coriell, R. F. Boisvert, J. I. Mickalonis, and M. E. Glicksman,

Adv. Space Research 2_, No. 5 (1983) 95.

22. J. M. Quenisset, R. Sokolowski, and M. E. Glicksman, J. Crystal

Growth 63^ (1983) 389.

23. C. J. Chang and R. A. Brown, J. Crystal Growth (1983) 343.

93
FIGURE CAPTIONS

Fig. 1

.

Fig, 2.

Fig. 3.

Fig. 4.

Fig. 5.

Fig. 6.

Fig. 7.

The critical concentration c* of tin in lead for the onset

of thermosolutal convection as a function of SO/V for a growth

velocity V of 2 ijn/s and a temperature gradient G|_ in the

liquid of 200 K/cm.

The wavelength 2Tr/oj of sinusoidal perturbations at the

onset of thermosolutal convection as a function of SD/V for

conditions corresponding to those given in Fig. 1.

The wave speed (-ci/oj) of sinusoidal perturbations at the

onset of thermosolutal convection as a function of SD/V for

conditions corresponding to those given in Fig. 1.

The critical concentration c« of tin in lead for the onset

of morphological instability as a function of SD/V for a growth

velocity V of 200 jin/s and a temperature gradient Gl in the

liquid of 200 K/cm.

The wavelength 2 tt/oj of sinusoidal perturbations at the

onset of morphological instability as a function of SO/V for

conditions corresponding to those given in Fig. 4.

The wave speed (-ai/u) of sinusoidal perturbations at the

onset of morphological instability as a function of SD/V for

conditions corresponding to those given in Fig. 4.

The concentration c«o of tin in lead and ai at the onset

of instability during directional solidification at V = 40 j/n/s

as a function of spatial frequency w of a sinusoidal perturbation

for a shear S = 0.01 s“^ and a temperature gradient G^ in the

94
liquid of 200 K/cm. The dashed curves correspond to morphological

modes and the solid curves correspond to thermosolutal convection.

Fig. 8. The concentration c« of tin in lead and ai at the onset of

instability during directional solidification at V = 40 i/n/s as

a function of spatial frequency u of a sinusoidal perturbation

for a shear S = 0.1 s"^ and a temperature gradient in the

liquid of 200 K/cm. The dashed curves correspond to morphological

modes and the solid curves correspond to thermosolutal convection.

Fig. 9. The concentration c« of tin in lead and ai at the onset of

instability during directional solidification at V = 40 ^n/s as

a function of spatial frequency w of a sinusoidal perturbation

for a shear S = 1.0 s"^ and a temperature gradient in the

liquid of 200 K/cm. The dashed curves correspond to morphological

modes and the solid curves correspond to thermosolutal convection.

95
(wt

%)

Linear Stability for Lead—Tin

V = 2 fJLTCi/s Gl = 200 K/cm

96
Wavelength

(c

Linear Stability for Lead-Tin

V = 2 jjLxa./s Gl = 200 K/cm

97
Wave

Speed

(cm/s)

Linear Stability for Lead-Tin

V = 2 /im/s Gl = 200 K/cm

98
Linear Stability for Lead-Tin

V = 200 yum/s Gl = 200 K/cm

99
Linear Stability for Lead—Tin

V = 200 fxm./s Gl = 200 K/cm

100
Linear Stability for Lead-Tin

V = 200 /2m/s Gl = 200 K/cm

101
(wt

%)

Linear Stability for Lead-Tin

V = 40 /zm/s Gt = 200 K/cm

S = .01 s-'

102
(wt

%)

Linear Stability for Lead—Tin

V = 40 yLim/s Gt = 200 K/cm

S = 0.1 s-‘

103
(wt

%)

8

o

Linear Stability for Lead—Tin

V = 40 jjunx/s Gt = 200 K/cm

S = 1.0 s-‘

2.0-1

0*0 1 1—I

r I I I I

I

1 1—I—

f

-

| -

1 Tn

10* l(f lO"

CO (cm“‘)

104
APPENDIX D

CONVECTION-INDUCED DISTORTION OF A SOLID-LIQUID INTERFACE

R. J. Schaefer and S, R. Coriell

Metallurgy Division

National Bureau of Standards

Washington, D.C. 20234

ABSTRACT

Measurements of convective flow fields and solid-liquid interface

shapes during the solidification of a pure and a slightly alloyed

transparent material reveal that the convective transport of solute can

cause a macroscopic depression to develop in the solid-liquid interface.

This effect occurs under conditions close to those which are predicted to

produce morphological instability of a planar interface. A cellular or

dendritic mi crostructure later develops within the interface depression.

The convection is attributed to the effect of radial temperature gradients

in the crystal growth apparatus.

105
I. INTRODUCTION

During unidirectional growth of a crystal by a Bridgeman or similar

process, convection in the liquid phase can be driven be several forces.

In the presence of gravitational forces, the most important source of

convective flow is ordinarily the density differences due to temperature

and conposition differences within the sample. Additional flow can be

driven by the volume change which accompanies the phase change, and if

any free surfaces (liquid-vapor or liquid-liquid) are present, surface

energy gradient (Marangoni ) flows are possible. The latter two types of

flow are present even in the absence of gravitational forces.

In the presence of gravitational forces, density -induced convective

flow can be avoided during crystal growth if special conditions are met.

In a pure material it is sufficient that the temperature increase with

height (for a normal material which expands with increasing temperature)

and that no gradients be present in the horizontal directions.

If a solute is present, the solidification process will lead to

composition differences ahead of the solid-liquid interface. If the

rejected solute is lower in density than the solvent, the conposition

differences can lead to a thermo-solutal convective instability. The

conditions which lead to this instability are not simply that the density

in some region increase with height: the conditions are determined by an

analysis of the behavior of elements of the fluid which are displaced

vertically [1 , 2 , 3].

When convective flow occurs, it redistributes both heat and solute

in the liquid phase, and thus causes the tenperature and concentration

fields in this phase to differ from those which would be calculated on

106
the basis of diffusive transport alone. As a result, the composition of

the crystal as it grows will be less homogeneous than would be possible

in the absence of convection. A steady state flow pattern can result in

lateral segregation while a non-steady flow pattern will produce additional

longitudinal non-uni formi ty

.

Because convective flow can modify the temperature and solute fields

in the liquid adjacent to the solid-liquid interface, it can cause local

variations in the conditions which determine the morphological stability

[4] of the solid-liquid interface. Thus for exanple it can lead to the

growth of crystals which contain localized regions of cellular mi crostructure

within an otherwise cell-free matrix. When a highly uniform crystal is

needed, it may thus be essential to keep these convective effects under

control

.

In this paper we present observations of convective phenomena during

the solidification of a transparent "alloy" system, succinonitri le

containing ethanol. The conditions which lead to the onset of thermosolutal convective and interfacial (morphological) instabilities in this

system are calculated, and other sources of convective flow are considered.

II. THEORY

The stability of the unidirectional upward solidification of

succinonitri le containing ethanol, in the absence of horizontal temperature

gradients, can be analyzed by the same methods previously used for the

solidification of lead containing tin [1]. The material properties used

in the calculations are listed in Table I: the measurement of some of

these properties is described in the experimental section. The results

of the calculations are qualitatively similar to those previously obtained

107
for the lead-tin system in that two types of instability are observed.

Instabilities of short wavelength, which are associated with deviations

of the solid-liquid interface from planarity, are termed interfacial

instabilities. Instabilities of longer wavelength, which are associated

with motion of the liquid above the interface, are termed convective

instabilities. Figure 1 shows, as a function of the velocity of

solidification, the concentration of ethanol in succinonitri le required

to produce the two types of instability. The calculations apply to a

system growing with a temperature gradient in the liquid of 10 K/cm, a

typical value for this type of system.

The concentration of solute required to produce interfacial

instabilities decreases with increasing velocity, as is the familiar

situation described by the theory of morphological instability except at

extremely high growth velocities. The convective instabilities show the

opposite velocity dependence. Unidirectional upward solidification with

a planar interface is stable at concentration-velocity combinations lying

below both curves in Figure 1. At concentrations near 3 x lO"^ wt.% ethanol

there exists a range of velocities in which the solidification is stable,

above which interfacial instabilities are predicted and below which convective instabilities are predicted. At concentrations above 5.4 x 10"3 wt.%

ethanol, growth is not stable at any velocity with this temperature gradient.

The dotted portion of the convective instability curve in Figure 1

represents an oscillatory convective instability which occurs in the

region where the convective and interfacial instability curves cross.

The onset of instability by an oscillatory mode is the most striking

difference between the results for lead-tin and succinonitri le-ethanol

.

108
Figure 2 shows the wavelength at the onset of instability relative

to the diffusion boundary layer thickness D/V, for a tenperature gradient

in the liquid of 10 K/cm. Whereas the wavelength of the interfacial

instabilities is approximately equal to the thickness of the diffusion

layer, the wavelength of the convective instabilities is typically about

an order of magnitude larger and at velocities of less than 1 ym/s the

predicted wavelength of the convective instabilities are greater than 1 cm.

Therefore one cannot expect to observe these instabilities in small

contai ners

.

Increasing the tenperature gradient tends to stabilize the system

with respect to both interfacial and convective instabilities. Figure 3

shows interfacial and convective stability curves for several values of

the temperature gradient in the liquid. One of the samples used in this

study contained 2.6 x lO"^ wt.% ethanol, as indicated in Figure 3, and it

is seen that at this composition the gradient in the liquid has only a

small effect on the velocity at which convective instability occurs, but

has a strong effect on the velocity at which interfacial instability

occurs.

In a real crystal growth process, horizontal temperature gradients

are difficult to avoid. This is especially true for materials such as

succi nonitri le which have low thermal conductivities. In such cases the

thermal conductivity of the container walls is comparable to that of the

material under study, and the requirement of maintaining good visibility

of the material in the solid-liquid interface region precludes the use of

surrounding insulation.

Radial gradients resulting from lateral heat losses can be minimized

109
by careful design of the crystal growth system, but in the vicinity of

the solid-liquid interface there is a thermal interaction with the chamber

walls which cannot be completely eliminated. Conservation of heat at the

interface leads to the relation

ksGs - k^Gi = vL (1)

where kg and k] are the thermal conductivities of the solid and liquid,

Gg and G] are the temperature gradients in the solid and liquid, v is the

growth velocity and L is the latent heat of fusion. In general there is

therefore a veloci ty -dependent change of the tenperature gradient across

the interface. Because there is no correspondi ng change of gradient in

the container walls, the isotherms in the vicinity of the container walls

are distorted from the ideal horizontal planes which would be needed to

avoid a driving force for convection.

In succi nonit ri le the thermal conductivity of the solid and liquid

are almost identical, so that at zero growth velocity there is virtually

no change in tenperature gradient across the interface and it is possible

to achieve almost planar interfaces. However, Eq. 1 indicates that there

will be a gradient change of 2,1 K/cm at the interface for each 1 ym/s

of growth velocity. As a result, at growth velocities of 1 ym/s or more

there will be significant horizontal gradients near the container walls

and thermally driven convective flow must be expected.

III. EXPERIMENT

Many of the properties of succi nonitri le have been reported in the

literature, as indicated in Table I. We have determined the relevant end

of the succinonitrile-ethanol phase diagram (Figure 4) and find the

liquidus and solidus slopes at low concentrations to be 3.6 K/wt.% and 81

no
K/wt.% respectively. In addition, we have measured the density of

succi noni tri 1 e-ri ch liquids (Figure 4). Tlie only property used in the

calculations which has not actually been measured is the diffusion

coefficient for ethanol in succinonitri le, for which an estimated value

has been used.

Samples for crystal growth experiments were sealed under vacuum in

borosilicate glass tubes, approximately 45 cm long and with 17 mm inside

diameter. The composition of the samples was evaluated from measurements

of the solidus and liquidus tenperatures . A high purity sample was found

to have a melting range of 1.2 x 10~2 K, corresponding to the effect of

1.5 X 10"^ wt.% ethanol or equivalent impurity, while an ethanol-doped

sample had a melting range of 0.215 K, corresponding to 2.6 x 10"^ wt.%

ethanol

.

For crystal growth the samples were drawn downwards by a tracking

mechanism through a temperature gradient produced by a cooling water

jacket and a small electric furnace. In the region between the water

jacket and the furnace, the interface was surrounded by a transparent

block with flat walls to reduce the optical distortion of the solid-liquid

interface region (Figure 5). A flat auxiliary heater was embedded in

this block to provide some control over the temperature gradients in the

interface region.

Convection was measured from multiple-exposure photographs of latex

microspheres floating in the liquid. The spheres were 10-20 <yn in

diameter and were close in density to the liquid succinonitri le so that

they tended to settle out only at velocities much less than the convective

flow velocities which were being measured. By using unequal intervals

m
between the multiple exposures, the direction of motion of individual

particles could be determined. The particles were best seen in dark

field illumination whereas bright field images showed the solid-liquid

interface most clearly.

IV. RESULTS

When the sample is held stationary (zero growth velocity), any

convection present in the liquid is attributable to the effect of radial

heat losses. Stationary solid-liquid interfaces are observed to be

slightly concave upwards, indicating that the tenperature is higher in

the center than around the edges. As expected under these circumstances,

a convective flow pattern is present above the interface with flow downward

near the walls, radially inward directly above the interface, and upwards

in the center of the tube. The flow pattern above the interface can

thus be described as a toroidal roll. Any thermal asymmetries of the

system result in displacement of the node and axis of this torus from the

central axis of the sample tube. If the asymmetry is sufficiently great,

the torus pinches off to form a single transverse roll (Figure 6).

When the toroidal flow remains for several hours in a stationary

sample, one frequently observes an accumulation of marker particles in

the liquid just above the interface at the node where the radial inward

flow converges and turns upward. A dense trail of particles streams

upward with the liquid from this point and along the axis of the toroidal

roll (Figure 7). The velocity of this upward flow near the center of the

tube, at a distance 2.5 mm above the interface, is typically 10-15 vm/s.

Smaller velocities could be attained when the interface region was well

insulated but only at the cost of reduced visibility.

112
In a typical experiment the interface was equilibrated for several

hours at zero growth velocity and the downward motion of the sample tube

at constant velocity would then be started. The solidification velocity

does not of course immediately attain the velocity at which the tube is

withdrawn: instead, the interface starts to move downward with the tube

and as the gradients change solidification commences and the solidification

velocity accelerates toward the tube withdrawal velocity, which is reached

when the interface reaches a steady state position. When crystal growth

starts, the radial gradients increase both due to the effect of drawing

the hotter material downward through the viewing block and due to the

latent heat emitted from the solid-liquid interface. The convective flow

pattern is correspondingly altered. In cases where the axis of the

toroidal roll was initially displaced from the tube center (Figure 6b),

the axis moves toward the center of the tube. In radially symmetric

cases, the upward flow velocity near the center of the tube increases

when crystal growth starts. Tlie upward flow velocities at the tube

center, 2.5 mm above the interface, were typically 20-40 j/n/s for a

sample growing at 2 ijm/s.

Figure 8 shows examples of measured convective flow tracks in the

sample containing 2.6 x lO"^ wt.% ethanol. Figure 8a shows the flow

field before the start of crystal growth, and Figure 8b shows the flow

field 19 minutes after the start of sample motion at 2 vin/s. It is seen

that in Figure 8b the flow velocities have increased significantly compared

to those in Figure 8a. Note that in Figure 8b the time is still small

compared to the time (D/v^k) required to establish the steady state solute

distribution, which in this case is approximately 95 minutes.

113
the slowest velocities at which it was seen (about 2 un/s) it attained

its most dramatic appearance (Figure 10), becoming sharply pointed and

about 1 mm deep before breaking down into a cellular structure. At a

growth velocity of 3 vm/s, the pit was considerably wider and shallower

when the cellular structure appeared (Figure 11) and at 4 im/s the

pit was hardly discernable from the overall concavity of the interface

when the cellular structure appeared. The range of growth velocities

over which this effect is seen is therefore very narrow, but it is a

critical range in that it is near the limits of the concentration and

velocity for which the growth of microsegregation-free crystals is

possible.

The marker particles indicate low flow rates in the liquid within

the pit, and in all cases show that the pit lies under the point at which

the flow converges and turns upward.

The formation of such pits has not been observed in the pure saiiple,

nor has it been observed under conditions where the flow forms a single

roll (Figure 6c) instead of a toroidal roll (Figure 6a or 6b).

The auxiliary heater surrounding the sample could be used to alter

the radial temperature gradients in the vicinity of the solid-liquid

interface. Ideally it would act as a substitute source of latent heat,

establishing a tenperature gradient discontinuity within the transparent

viewing block to match that created within the succinonitrile by the

emission of latent heat (Eq. 1). In practice it was generally found that

as the radial gradients were decreased by the use of the auxiliary heater,

the toroidal flow field could be suppressed but there were always enough

asymmetries in the system to yield a single roll flow pattern. By

114
Together with the change in the fluid flow, there is a change in

the interface shape which develops when crystal growth starts. In the

pure material the interface assumes the shape of a smooth upwardly concave

bowl (Figure 9). This sample was sufficiently pure that no conditions

(other than bulk supercooling of the liquid) were found which resulted in

interfacial instabilities. Qualitatively, this is the interface shape

change which would be expected even in the absence of fluid flow. It is

not known how much additional interface shape change is produced by the

fluid flow which does occur,

A much more conplex shape change occurs during crystal growth in the

ethanol -doped sample at velocities close to that which leads to the

formation of interfacial instabilities. For growth at 1 ym/s, the

interface develops a gentle bowl shape similar to the pure material, and

remains featureless. For growth at 5 jiii/s or more the interface

rapidly develops an overall cellular structure. At intermediate velocities,

however, the interface starts to develop a smooth bowl shape but after an

incubation period which varies in length from a few minutes to as much

as two hours, an additional macroscopic depression develops in the

interface directly under the point where the radial inward flow above the

interface converges and turns upward. Once the depression appears, it

deepens rapidly until interfacial instabilities develop within it. These

instabilities then develop into cells or dendrites which grow rapidly to

fill the pit up to the level of the original bowl-shaped interface. This

interface then advances with a localized cellular region in the area

where the pit had been, surrounded by a featureless (cell-free) area.

The shape of the pit varied with the solidification velocity. At

115
activation of the auxiliary heater it was found that the formation of the

interface pit could be suppressed during growth at 2 vm/s, but the

effect may better be attributed to the change of the flow pattern rather

than to an actual suppression of the flow.

An additional convective flow was occasionally observed in cases

where a vapor bubble was trapped at the solid-liquid interface. This

cold occur when the sample was melted down into a region where an axial

shrinkage pipe had been formed by a previous freezing operation. When

the melt contacted the shrinkage pipe, the vapor bubble which formed

would sometimes be held at the interface by surface tension forces before

breaking away and rising to the top of the tube. The flow in the liquid

near the bubble was often similar to that which occurred in its absence,

but at other times a counter-rotating toroidal flow would be present in

the liquid immediately above the bubble (Figure 12). The flow rates in

the counter-rotating flow were at least an order or magnitude faster than

those in the normal flow pattern, and immediately adjacent to the bubble

surface they were so rapid that they were not recorded on the photographs.

The counter-rotating flow was most vigorous shortly after slight changes

in the position of the bubble.

V. DISCUSSION

The formation of a transient macroscopic pit in a solidifying sample

surface, in the presence of both convection and solute, appears to be an

important stage in the development of a localized cellular structure.

The convective flows observed in these experiments can be attributed

to the effects of radial temperature gradients, which in the case of

stationary samples are the result of heat losses to the environment and

116
in the case of growing samples have an additional contribution from the

interaction of the latent heat with the chamber walls.

In stationary samples, the conposition of the liquid is uniform and

therefore is unaffected by the convective flow. In a solidifying sample

with a solute concentration Cq far from the interface, a solute-enriched

layer is present in the liquid adjacent to the interface. At steady

state and in the absence of convection, the solute concentration Cl in

the liquid ahead of a planar interface would have a simple exponential

form

C|_ = Cq [1 + (2)

and the composition immediately adjacent to the interface would be Co/k

.

Because of the low value of k for the succinonitrile-ethanol system

(0.044), the concentration differences are large. Thus in the sample

containing Cq = 2.6 x 10"^ wt.% ethanol, the composition of the liquid

adjacent to the interface would be 5.9 x 10"2 wt.% ethanol and the

interface tenperature would be depressed 0.2 K below the temperature for

equilibrium of solid with liquid of composition Cq.

Convective flow will change the temperature and solute fields in the

liquid ahead of the interface. However, because the thermal diffusivity

is approximately 100 times as large as the solute diffusivity, the solute

field can be much more severely distorted than is the temperature field.

The distortion of the liquid temperature and composition fields will lead

to changes of the interface shape.

In a crystal growth system such as that shown in Figure 5, in the

absence of convection the solid-liquid interface of a growing crystal

will not be planar but will be macroscopical ly concave due to the radial

117
heat losses and the interaction with the chamber walls. When solute is

present, this concavity will cause the solute distribution in the liquid

adjacent to the interface to be somewhat non-uniform even in the absence

of interfacial instabilities [8],

If we add in the effect of convection, the interface shape will be

modified in the pure material by the convective transport of heat. In a

crystal growth apparatus such as ours, we can expect the convective flow

to reduce the deviation of the interface from planarity because the

downward flow on the cooler parts of the tube (near the walls) will bring

warmer material from above to the interface region. The actual magnitude

of the interface shape change due to the convective flow in the pure

material is not known, although some indication of it could be obtained

by comparing interface shapes in our experiments to shapes observed when

the solid lies above the liquid, in which case the convection would

increase the deviation from planarity. This type of geometry was

approximated during the zone refining of the samples, and the deviation

from planarity was much more pronounced.

During solidification of the sample containing ethanol, the convective

flow can change the distribution of solute ahead of the solid-liquid

interface. A simple uniform horizontal shear flow above the interface

would not alter the exponential solute distribution, but substantial

effects can be expected with the fluid flow patterns which prevail in a

finite sized system. For the toroidal flow as illustrated in Figure 6a,

the material flowing down near the walls will have composition close to

that of the bulk liquid and the concentration of ethanol ahead of the

interface near the walls will therefore be less than that predicted in

118
the absence of convection. As the liquid moves radially inward above the

interface, it will accumulate the rejected solute. Because of the radial

nature of the flow, solute will be relatively concentrated in the vicinity

of the node where the flow converges and turns upward, in a manner somewhat

analogous to the process by which the marker particles accumulate in this

vi ci ni ty

.

We feel that the interface pit which develops at certain growth

velocities is the result of this localized concentration of solute by the

thermally-induced convective flow field. The evidence supporting this

hypothesis is that

1) The pit is seen only in the sairple doped with ethanol; and

2) The pit occurs directly under the point where the flow field

converges, even when this point is displaced from the tube center by

asymmetries of the thermal field, and formation of the pit is not observed

when the asymmetries are strong enough to cause the flui.d. flow to change

from a toroidal (Figure 6a) to a single roll (Figure 6c) configuration.

Chang and Wilcox [8] have reported observations of solid-liquid

interface shapes as affected by radial heat flow and convection during

zone melting of napthalene. Their heat flow geometry, and hence their

convective flow fields, were different from that which occurred in our

experiments. They observed localized cellular break-down of the solidifying

interface, which in some cases occurred within an indentation of the

interface. They attribute the localization of interface break-down to

the redistri bution of heat and solute by convection currents, but do not

report the effect which we find in which the convective redistribution of

solute causes a macroscopic change of interface shape.

119
VI. CONCLUSIONS

A small concentration of ethanol in succi noni tri le is sufficient to

induce interfacial or thermosolutal convective instabilities during

unidirectional upward solidification. Radial temperature gradients

present during solidification in a cylindrical tube also induce convective

flow, and when solute is present it is redi stri buted by the convective

flow to produce changes in the shape and interfacial stability of the

solid-liquid interface.

Acknowledgement

The authors acknowledge the support of NASA's Microgravity Science and

Applications Division.

120
References

1. S. R. Coriel 1, M. R. Cordes. W. J. Boettinger, and R. F. Sekerka, J.

Crystal Grov/th 49, 13

2. S. R. Coriel 1, M. R. Cordes, W. J. Boettinger, and R. F. Sekerka,

Adv. Space Res. 1, 5 (1981).

3. S. R. Coriell and R. F. Sekerka, PhysicoChemi cal Hydrodynamics 2, 281

(1981 ).

4. W. W. Mullins and R. F. Sekerka, J. Appl. Phys. 35, 444 (1964).

5. K. J. Timmermans and Mme. Hennaut-Roland, Journal de Chimie Physique

34, 693 (1937).

6. C. A. Wulff and E. F. Westrum, J. Phys. Chem. 67, 2375 (1963).

7. R. J. Schaefer, M. E. Glicksman and J. D. Ayers, Phil. Mag. 32, 725

(1 975).

8. C. E. Chang and W. R. Wilcox, J. Crystal Growth 21, 182 (1974).

121
Figure Captions

Fig. 1

Fig. 2

Fig. 3

Fig. 4

Fig. 5

Fig. 6

Fig. 7

Fig. 8

Fig. 9

Fig. 1

Fig. 1

Fig. 1

. The concentration of ethanol in succi noni tri 1 e above which

interfacial or convective instability occurs.

. The ratio of the instability wavelength, "A , to the diffusion

boundary layer thickness, D/v, at the onset of instability.

. Effect of the temperature gradient in the liquid on the

concentration of ethanol in succi noni tri 1 e required to induce

instabi 1 ity

.

. Succinonitrile-ethanol phase diagram, with lines of constant

liquid density.

. Exploded view of crytal growth apparatus.

. Convective flow patterns above solid-liquid interface for (a)

symmetrical, (b) slightly asymmetrical and (c) highly asymmetrical

heat losses.

. Plume of particles marking the axis of a toroidal flow pattern.

. Measured flow velocities above the solid-liquid interface (a)

before the start of growth and (b) 19 minutes after the start of

growth.

. Concave interface resulting from growth in pure succinonitri le.

0. Macroscopic depression in the solid-liquid interface of a sample

containing 2.6 x 10’4 wt.% ethanol, growing at 2 un/s.

1. Macroscopic interface depression in sample growing at 3 un/s.

2. Vigorous counter-rotating toroidal flow above vapor bubble in

sample containing 2.6 x 10~^ wt.% ethanol.

122
Table. I. Properties of Succinonitrile

Property Symbol Value Reference

Density of Liquid 983 kg/m^** [5]

Thermal Expansivity of Solid

“s -5.6 X IO'Vk [6]

Thermal Expansivity of Liquid -8.1 X 10-Vk [5*1

Shear Vi scosity 0 2.6 X 10-3 Pa-s**

Kinematic Viscosity V 2.6 X 10-Sm2/s** [5]

Latent Heat of Fusion L 4.70 X 10^ J/kg [6]

Thermal Conductivity of Solid 0.225 J/mKs** [7]

Thermal Conductivity of Liquid k^ 0.223 J/mKs** [7]

Thermal Diffusivity of Solid ^s 1 .16 X lO-'^m^/s**

Thermal Diffusivity of Liquid 1 .12 X 10-"7m2/s** ***

Distribution Coefficient for

Ethanol in Succinonitrile

0.044 ****

Liquidus slope for Ethanol

in Succinonitrile

mi 3.6 K/wt.%

Solidus slope for Ethanol

in Succinonitrile

ms 81 K/wt.%

* Indicates value derived fron properties given in cited reference.

** Indicates value of property at melting point, 331 .23K. *** Indicates property derived from other properties in the table. **** Indicates properties measured in current work.

123
Ethanol

Concentration

(wt.,%)

124
100

a/ YA MlBueieABM pazMBuuoN

125

1.0

10.0

Solidification

Velocity

V

((im/s)
SuccinonitrileEthanol

(%'IM) °®Q uo!)ej}uaouoo jouem^

126

Solidification

Velocity

V

(|im/s)
Ethanol

Concentration

(wt.%)

Temperature (°C)

lO

CO __

cn

127
Furnace

128
129
Fig. 7

130
131
Fig. 9
Fig. n

/

Fig. 12

133

Task 3

Measurement of High Temperature Thermodynamic Properties

D. W. Bonne!

1

Inorganic Materials Division

Center for Materials Science

SUMMARY

The newly developed General Electric Space Systems Division (GE)

imaging optical pyrometer was completed, tested, and calibrated in inservice operating conditions. This pyrometer is designed to allow accurate

temperature measurement of a levitated high temperature sample whose specific

heat then can be determined by dropping it into a calorimeter, or in space

flight by moving the calorimeter to enclose the sample. The pyrometer has

the capability for two color operation, with the alternate colors sampled

at several hertz. The unit proved to be stable in operation, and appears

to be a reasonable prototype. Special techniques were developed to connect

the temperature scale of the pyrometer to that of the standard lamp.

For the electromagnetic induction levitation system, a new, tighter

coil design was tested with major improvement in holding ability and stable

levitation for high melting point metals. Successful coupling of this

levitation coil and the electron-beam heat source needed to melt the sample

proved more difficult than anticipated, with the result that successful

complete melting of tungsten has not yet been achieved.

An additional major effort of this reporting period is the preparation

of a chapter on Levitation Calorimetry, which represents a comprehensive

review of the subject of electromagnetic levitation and thermophysical

property measurement, the first such effort since the early 1960's. Levitation

activities at Rice University have resumed in full force, with complete

reactivation of their facilities. A number of new materials (silver.

135
gallium, etc.) are being studied, and regular discussions indicate that

publication quality results are now being obtained.

Introduction

The primary objective of this task has been to develop, in conjunction

with General Electric Space Systems Division (GE), and Rice University (RICE),

the techniques, methodology, and experimental apparatus for measuring

thermophysical properties (e.g., heat of fusion, enthalpy increment, heat

capacity, etc.) of materials at extremely high temperatures. A primary

technique has been the development of a coupled electron-beam heating/

electromagnetic levitation-calorimeter system. The e-beam heating/levitation apparatus is based on equipment developed originally by GE specifically

for melting tungsten. As the heat of fusion and enthalpy increment functions

for tungsten are known only through estimation and indirect measurements

[Bonnell, 1983], and tungsten occupies a unique position as the most refractory

element, adapting this apparatus to calorimetric measurements seemed a

particularly significant opportunity to expand the limits of ground-based

property measurements.

The technical difficulty of merging these two techniques has proved to

be a measure of the dividing line between experiments which can be performed

on the ground, and those which will require reduced gravity for accomplishment.

The significance of these experiments lies in the observation that the

temperatures require heater powers (10 to 100 kW) for reasonable sized samples,

far in excess of that expected to be available in space in the near future.

Thus, to obtain the scientific data, and to continue to develop methods for

use in space, the clearest path which will allow development of hardware

able to take advantage of the low-gravity environment of space has been to

pursue limiting ground-based cases.

136
GE Pyrometry System

The new GE imaging pyrometer was optically similar to the previous

system, and installed in the same manner. This arrangement and the techniques

developed for calibration and operation are discussed in detail elsewhere

[Bonnell, 1981]. The new GE imaging pyrometer system was designed around a

commercially available photodiode array camera system. Custom circuitry to

allow the gating of any individual pixel to two conpanion high quality

sample-and-hold operators was included. This allowed dc sampling and

display of the selected single pixel output for each of two colors passed

by two narrow hand filters. The filters were mounted in an arm which could

be alternated between the two filters by a 180° stepper motor. The location

of the arm in the two alternate positions was sensed by switch closures,

and the camera image collection (pixel gating) restricted to frame times

when the filters were stationary. The normal 30 hz framing rate was used,

and the filters could be changed with the loss of only one frame time. The

filter change frequency was limited by the motor to a few hertz.

Since the camera uses silicon diode based detection, and is thus

sensitive to wavelengths beyond 1000 nm, it was practical to select

filters straddling the sensitivity peak by choosing one of the band-pass

wavelengths in the infrared. The selected filters were 850 and 656 nm

center frequency, 10 nm Full -Width-at-Half -Maxi mum band pass, approximately

13 mm diameter dichroics. The extreme narrow band pass choice was based on

GE's expectation of radiance at the melting point, and the available maximum

aperture, approximately 10 mm. The choice was correct and at melting,

output from the pyrometer for both wavelengths was about mid-scale. However,

there was no direct way to obtain significant output from the system at

137
standard lamp temperatures which did not jeopardize the stability of the

standard lamp (about 1600°C for vcauum [ASTM, 1975]). This determination

is essential to calibrate the unit to the brightness temperature scale, to

obtain scale intercomparisons among the various aperture/fi Iter combinations,

and to perform running checks on outgas/vaporization coating of the vacuum

system windows and mirrors.

One way to deal with this problem would be to use wider band-width

filters, and insert neutral density filters for coarse scale control. This

alternative is still possible, but was not available at test time.

A modification of the camera control circuitry so that image collection

could be spread over multiples of the frame time allowed the camera to be

used as an integrator. Since calibration is a static process, this was a

reasonable alternative, although it made alignment of the camera somewhat

more awkward. Even though digital multiples of the frame time were used

(typically 10 frame times), the actual integration ratio had to be determined

empirically, as two frame times includes the retrace time, for example,

which is not part of one frame time. Also, there is the possibility that

the sensor array may not store the image without losses. The empirical

determination yielded a gain value for ten frame integration of 10.47 ± 0.3.

The error was estimated based on the digital granularity of the smallest

signal

.

A disk with a series of apertures was located behind the filter arm

to provide known aperture limiting, and wavelength independent gain. The

apertures were measured with a reticule comparator. Even using the integration technique, it was not possible to get enough signal to intercompare

the smaller aperatures with the standard lamp. If the current filters are

138
to remain as is, it would be wise to dismount the apertures and measure

their areas photometrically. Table 1 lists the aperture plate dimensions.

The lineal measurement pair represent the largest and the smallest cross

sections. The holes were not visibly oval so the aperture area ratios

were based on the average of the extrema measurements. For the temperatures

observed (maximum ~3800 K), it was never necessary to use other than the

A and B apertures.

Table 1

Pyrometer Apertures

Aperature Di ameter

i n.3

Area^

Sq. in.

Ratio

A 0.239

0.243 0.04562 > 2.13

B 0.162

0.168 0.02138

> 1.96

C 0.117

0.119 0.01094 > 2.00

D 0.082

0.085 0.00548

1

3 in. = 2.54 cm.

Calibration of the pyrometer was carried out in the manner outlined

before [Bonnell, 1981]. Briefly, the pyrometer output was compared to a

calibrated micro-optical pyrometer reading the temperature of a miniature

strip lamp inside the vacuum system. The levitaton coil system flange was

139
replaced for this purpose by a pyrex plate to allow direct viewing of the

filament. The pyrex plate absorbance was measured separately by measuring

the brightness temperature of an external lamp with and without the plate

in the optical path.

The current (measured by the voltage drop across a quality 0.0100 ohm

shunt) applied to the lamp during the process of calibrating the pyrometer

served to calibrate the standard lamp for later use in determining vacuum

window and mirror absorption changes. It is to be noted that the strip lamp

was used at brightness temperatures much higher than the long term limits

recommended for pyrometers [ASTM, 1975]. The stability of the lamp was checked

by reproducing temperatures below the safe limit after use at the higher values.

The observed output of the pyrometer (V, in volts) is related to the

brightness temperature by the expression.

derived from the familiar Wein approximation to the Planck radiation law.

This form indicates that the best fit to experimental T versus photometer

output is a one parameter fit to determine k. Two parameter fits, in which

the leading term (X/C 2 ) is a parameter, produced fits in which the leading

term was within statistically estimated error limits of the theoretical

slope, but such fits actually represent the data more poorly, as they do

not have correct extrapolation behavior in the original data space. The

one parameter fit for the 656 nm filter is

In V - k

^2

IOOOO/T3

= 3.759 ± 0.02 - 0.4561 In V

as shown in Figure 1. For the 850 nm filter, the result is

10000/Tg = 3.362 ± 0.02 - 0.5908 In V

140

I
Note that the functions are for V measured in the xlO integrate mode. Since

the true gain, g, due to integration was measured (g = 10.47), the relationship

for the single frame high temperature operation mode is obtained by substituting

V = Vg in the above. This gives the working expression for brightness

temperature as,

10000/Tg = 2.688 ± 0.03 - 0.4561 In V (656 nm filter)

and

1000/Tg = 1 .975 ± 0.03 - 0.5908 In V (850 nm filter)

Corrections for aperture are gain corrections, and made in the same

manner. Ln(E), the emissivity, and window obscuration corrections are also

treated as gain terms, and the corrections made accordingly.

In operation, the pyrometer system showed a slight but disconcerting

image shift with the changes of the two filters. The shift was never large

enough to significantly affect normal temperature measurements, but as the

standard filament image is only about twice the width of a measurement

pixel, calibrations were made with the automatic filter switching off.

The system was free of the interference and much of the noise of the

previous system. The window obscuration problem, due to vapor deposition

at high temperature, is still a major correction, and changes rapidly

during an experiment. This makes it virtually impossible to determine

temperature even approximately during an experiment and remains a significant

problem. A method for dealing with this problem in real time has been

suggested [Bonnell, 1981], but not yet prototyped.

Levitation Calorimetry

An extensive effort to obtain more data on liquid tungsten was made

during this reporting period. The changes suggested from the last effort

141
[Bonnell, 1983] were made, with mixed results. A newly designed coil was

installed (see Fig. 2). This coil provided improved sample stability, with

little tendency to lose the samples due to inadvertent contact with the

coil, or due to sample ejection during excessive oscillations, once the

upper to lower coil gap was decreased. The upper turn was deliberately

left larger than optimum to minimize effects of the RF field on electron

beam focus, and to insure that it could be protected by the molybdenum

shade ring. This coil was nearly optimum for the samples. Although

liquefying tungsten is a serious problem, and has become increasingly

difficult as the apparatus ages, in two cases where other system failures

prevented capture, samples which were apparently completely molten were

contained by the coil.

The primary problem with this coil system was that it intercepts more

radiant energy from the levitated sample because its turns are much closer

to the sample. Perhaps even more importantly, beam electrons, forward-scattered

from the sample, are scattered onto the lower turns of the coil more effectively.

The combination of these two effects caused a series of coil failures, and

subsequent severe contamination of the aluminum walled vacuum system with

water. Although the system has excellent safety interlocks, the rupture of

a coil turn under 50 to 60 psig water pressure results in ounces of water

inside the vacuum system. When the problem was finally correctly diagnosed,

higher pressure water pumps were not available. We were able to divert the

water supply for another GE experiment to the task. That water system was

designed to deliver about 50 gallons to a linear levitation system over the

period of minutes. It consisted of a large 125 psig rated pressure vessel

filled with water and pressurized to 120 psig with nitrogen gas. At this

142
pressure, the levitation coil used about 4 liters of water per minute, an

approximately 1.5 improvement in water flow. This allowed 35 minutes

operation time with margin for orderly shutdown. This rate of water flow

was adequate for the increased heat dissipation requirements. It would be

desirable in the future to have a pump to supply 125 psig water, although

no problem ensued from the pressurized system. Higher pressure would also

be recommended, up to the working limit of the reinforced hose, as it is

still not certain how much safety margin remains.

Two additional problems accompanied the coil failures. The first was

repeated failure of the automatic gate system to function. This trouble

was finally identified as being caused by insufficient clearance for the

radiation gate hub, brought about by repeated minor splatters of molten

tungsten, aluminum (from the upper gate housing and vacuum system walls)

and tantalum (from the gate system heat shields). Curiously, the system

would work on the bench, but jam erratically under vacuum. The problem was

rectified only by a complete overhaul of the entire system. The other

problem was a return of the sample contamination problem. The problem was

not immediately identifiable due to difficulty with electron beam melting.

When finally correctly diagnosed (as due again to soft solder vapor contamination)

it was found that a small run of soft solder from a feedthrough cooling tube

was in line of sight to the coil from below. The problem was rectified by

removing all excess soft solder. However, extreme care will be necessary in

the future as critical joints in the coil flange are held together with

soft solder, which cannot be removed.

Scanning Electron Micrographs (SEM) examination of a sectioned contaminated

specimen with x-ray fluorescence showed of the order of 1/2 percent zinc present.

Quantitation was difficult, as 30 KeV excitation had to be used, which produces

143
some spurious peaks from the steen vacuum chamber. It was not possible to

tell if there was significant concentration along grain boundaries. The

sample showed large grains, indicative of near melting. Except for the

oblate shape change without surface melting, it would be difficult to

differentiate this specimen from a good run. These results are tentative

at present, with more work planned.

The electron beam continued to be the limiting part of the system. A

molybdenum shade ring was reinstalled to protect the top turns of the coil.

This ring was in place during the coil failures, implying that the direct

electron beam melting of the coil was not a likely cause of coil failure.

Extensive tests of the interaction of the beam and coil were carried out at

all power levels of beam and levitation field. Tests involved samples

mounted on the previously developed "sting," [Bonnell, 1983] as well as

freely levitated specimens. No discernable defocusing of the beam by the

RF coil occurred under any operating conditions. The electron beam image

consisted of two similar oval shaped spots separated by approximately 1 cm.

It was not possible to merge the spots without defocusing the beam. Any

focus degradation made it impossible to melt the sample, even with 100 mA,

30,000 V beam currents. In addition, when placing one spot on the sample,

the other spot could not be seen, and presumably impinged the side of the

sample, or the coil. This uncertainty makes analysis of the problem more

difficult. Another difficulty appears to be related to the need for a

precisely shaped electron beam spot. If the spot is too small, drilling of

the specimen occurs [Bonnell, 1983]; if too large, the sample will not

melt. Considerable more effort is required here by GE.

The new coil provides good holding for the current 0.9 cm diameter

samples. A possible improvement in melting ability should be attainable

144
by using even smaller samples. Samples 0.65 cm in diameter are recommended.

The main requirement seems to be to improve the electron beam image

shape, and to eliminate the extra spot. This will probably require

dismounting the beam source and adjusting it in successive tries.

Tungsten Heat of Fusion

A selection of recovered samples were sectioned and examined metal lographically. All the sectioned samples showed interior voids, suggesting

that significant 1 iquefication had occurred. However, no sample showed

clear regions with no grain boundaries. It is apparent that during long

heating in the levitation coil, grain growth occurs. Figure 3 shows the

original grain structure, both of the original rod stock and a ground

sphere. Figure 4 shows micrographs of a typical run, labeled A [see Bonnell,

1983]. There appear to be two distinct regions. The smaller region with

finer grain structure is in contact with the central void, and the region

(upper right in figure) where the electron beam impinged. The difference

in grain structure between this region and the rest, which shows major

grain growth, may be evidence for separation of melted versus not melted

tungsten. However, a more reasonable explanation is the grain growth occurs

throughout the heating process, and much of the grain structure could be

due to multiple nucleation sites upon resolidification. The sectioned

samples also suggest that different parts of the sample could have been

molten at different times. This is consistent with the difficulty in

melting noted above.

The sectioned samples do support the previously assigned extent of

melting, but are not conclusive enough to improve those assignments.

Figure 5 shows a plot of the available melting tungsten data, fitted

145
without constraint. The agreement with JANAF [1971] at the solidus is

excellent. The disagreement at the liquidus simply represents the previous

poor state of data for liquid tungsten, as has been previously discussed

[Bonnell, 1983],

Although the present data are a major improvement over previously

available data, measurements of the completely molten material are still

needed, both for the heat of fusion and to obtain the liquid heat capacity.

Levitation of Tantalum

During the latest tungsten series, a series of tantalum samples were

levitated. Tantalum proved easy to melt with the electron beam, requiring

less than 60 percent (of 100 mA available) beam current. The liquid could

not be held in the coil as was anticipated. The automatic gate and capture

system functioned perfectly, and two molten runs were captured. Post run

analysis, assuming the observed thermal arrest was the melting temperature,

yielded a heat of fusion, relative to the JANAF evaluated solid, of less

than 1 kJ/mol . Even correcting to the reported melting temperature with

generous error estimates, the observed heat of fusion was less than 8 to

10 kJ/mol , in serious disagreement with 36.6 kJ/mol heat of fusion value

from the exploding wire measurements of Lebedev et al., [1971]. Our measurements yield an unrealistically low value for AS(fus). Other tantalum

samples in the same series showed the shape deformation characteristic of

vaporized soft solder contamination. A number of tungsten samples previous

to the tantalum runs were also suspected of being contaminated. However,

GE assured us that no soft solder had been used in the coil region, and

sample preparation was suspected. With the evidence of bad tantalum runs

as well, a very careful examination revealed a tiny smear of soft solder on

146
the underside of a coolant line in line of sight to the coil, a flow which

had occurred in installing a new coil.

All available tantalum samples had been exposed to the soft solder, so

no more runs were attempted. Tantalum does appear a good candidate for

further study with this system.

Acknowledgements

Mr. C. Brady of the NBS Metallurgy Division prepared the metal sections,

and Dr. D. Butrymowicz and Mr. C. Olsen of that Division helped with the

SEM examinations.

147
References

1. ASTM Subcommittee E 20.02 on Optical Pyrometers, G. D. Nutter, Chairman,

"Recommended Practice in Optical Pyrometry" (1975).

2. D. W. Bonnell, in "NBS: Materials Measurements," ed. J. R. Manning,

NBSIR 81-2295 (1981 ).

3. D. W, Bonnell, in "NBS: Materials Measurements," ed. J. R. Manning,

NBSIR 82-2560 (1982).

4. D. W. Bonnell, in "NBS: Materials Measurements," ed. J. R. Manning,

NBSIR 83-2772 (1983).

5. D. W. Bonnell, R. L. Montgomery, B. Stephenson, P. C. Sundareswarean,

and J. L. Margrave, "Levitation Calorimetry" in Specific Heat of Solids , ed. A. Cezairliyan, in press, (1984).

6. A. Cezairliyan, NBS Report 10326, July 1970, Chs. 7 and 8.

7. JANAF Thermochemical Tables, NSRDS-NBS 37 (1972) and updates through 1983.

8. S. V. Lebedev, A. I. Savvatimskii , and Yu. B. Smirnov, High Temperature

USSR 9, 635 (1971 ).

148
Temperature

(K,

corrected)

2000

2500

Fig. 1. Graph of 1 parameter fit to pyrometer calibration data with 656 nm

filter. V is observed voltage for integration over ten frame times.

149
Fig. 2. a. Side view of coil, with inch scale. Note that the shape is

nearly solenoidal, and that the top reverse turn is closely

spaced.

Fig. 2.b. Top view, showing throat diameter, and large upper reverse turn.

Scale tic marks are l/16th inch.

150
Fig. 3. a. Original grain structure of tungsten rod. Scale tic marks are

mi 1 1 i meters.

Fig. 3.b. Somewhat enlarged grain structure of large sphere after heating.

Maximum temperature less than 3000 K. The flattened area at

lower left is where this sphere melted into the aluminum vacuum

wall on impact.

151
Fig. 4. Grain pattern for typical drop (run "A"). The electron beam

impinged along the long axis of the sample from upper right.

Note void, (tear drop shaped region). This sample is estimated

to have been 75 percent molten at impact in the calorimeter.

Scale tics are mm.

152
MOLTEN

-

3695 K

ENTHALPY INCREMENT (kJ/M)

Fig. 5. Plot of fit to tungsten drops, showing the agreement with the

JANAF [1972] selection. In the liquid region, the JANAF curve

has been raised 18 kJ from its estimated value. H(fus) determined

here is 53.0 ± 2.3 kJ/mol.

153

3400

3600

3800

4000
HBS-114A (REV. 2 -aoi U.S. DEPT. OF COMM.

BIBLIOGRAPHIC DATA

SHEET (See in struct/on s)

1. PUBLICATION OR

REPORT NO.

2. Performint Ortan. Report No,

NBSIR 84-2882

3. Publication Date

June 1984

4. TITLE AND SUBTITLE

NBS: Materials Measurement

5. AUTHOR(S)

J. R. Manning

6. PERFORMING ORGANIZATION (If Joint or other than NBS. see Instructions)

NATIONAL BUREAU OF STANDARDS

DEPARTMENT OF COMMERCE

WASHINGTON, D.C. 20234

7. Contract/Grant No. NASA Govt. Order

H-27954B

t. Type of Report & Period Covered

Annual Report

Apr. 1, 1983-Mar. 31, 1984

9. SPONSORING ORGANIZATION NAME AND COMPLETE ADDRESS (Street. City, State, ZIP)

Microgravity Science and Applications Program

National Aeronautics and Space Administration

Washington, D.C. 20546

10. SUPPLEMENTARY NOTES

I—

I

Document describes a computer program; SF-185, FIPS Software Summary, Is attached.

11. ABSTRACT (A 200-word or less factual summary of most significant information. If document Includes a significant

bibliography or literature survey, mention it here)

This report describes NBS work for NASA in support of NASA's Microgravity

Science and Applications Program under NASA Government Order H-27954B (Properties

of Electronic Materials) covering the period April 1, 1983 to March 31, 1984.

This work is being carried out in three independent tasks:

Task 1. Surface Tensions and Their Variations with Temperature and

Impurities

Task 2. Convection During Unidirectional Solidification

Task 3. Measurement of High Temperature Thermophysical Properties.

The results for each task are given in detail in the body of the report.

Tasks 1 and 2 are directed toward determining how the reduced gravity obUined

in space flight can affect convection and solidification processes.^ Emphasis

in Task 3 is on development of levitation and containerless processing^

techniques which can be applied in space flight to provide thermodynamic

measurements of reactive materials.

12. KEY WORDS fSIx Co twelve entries; alphabetical order; capitalize only proper names; and separate key words by semicolons)

binary alloy; convection; directional solidification; silicon; surface tension,

thermodynamic properties; tungsten

13. AVAILABILITY

[ XI Unlimited

I—

I

For Official Distribution. Do Not Release to NTIS

r~] Order From Superintendent of Documents. U.S. Government Printinf Office, Washintton, D.C.

20402. gOrder From National Technical Information Service (NTIS), Sprintfleld, VA. 22161

14. NO. OF

PRINTED PAGES

155

15. Price

$16.00

USCOMM'OC •04S-P«0