Nuclear Instruments and Methods in Physics Research A278 (1989) 425-430 North-Holland, Amsterdam
425
SUPERCONDUCTING DETECTOR FOR MINIMUM IONIZING PARTICLES A. GABUTTI, R.G . WAGNER, K.E . GRAY and R.T . KAMPWIRTH Argonne National Laboratory, Argonne, IL 60439, USA
R.H. ONO
National Institute of Standards and Technology, Boulder, CO 80303, USA
Received 19 December 1988
Although the detection of a-particles by thin superconducting films has been demonstrated previously, detection of the significantly smaller energy deposited by minimum ionizing particles (mips) requires independent verification . In this paper, experiments using - 6 keV X-rays, which deposit energies comparable to mips, are used to show that the switching probability is smaller than predicted by the simplest energy balance model. As a result, the simplifying assumptions of this model are in question, and it is demonstrated that practical detectors of mips will greatly benefit from, and most probably require, superconducting transition temperatures which are close to a low operating temperature (e .g ., 4.2 K) . In addition, an existing thermal propagation model is shown to adequately describe the behavior of the normal region after switching.
1 . Introduction Recently considerable attention has been devoted to superconductive detectors and their applications to high energy physics experiments [1]. For example, superconducting strips of NbN exhibit [2] very high radiation hardness and could be used for a microvertex detector [3] to provide "close-in" track reconstruction with point-to-point spatial resolution of a few microns. Although successful switching results [4] have been obtained on indium [5], tin [6] and aluminum [7] strips, the transition from the superconducting to normal state across the full width of a film was caused by ionization from a - 5 MeV a-particle . In a microvertex detector, the energy loss by a minimum ionizing particle (mip) traversing a thin film is approximately 200 times smaller than the energy deposited by a 5 MeV a-particle . In addition, a-particles deposit appreciable energy in the substrate which assists in switching the film . Therefore, the capability of superconducting strips to detect mips must be independently verified, and it is expected that films narrower than the 5-50 ~tm used in the a-particle experiments will be required. We report the beginning of an in-depth experimental study of the energy balance model [8]. For convenience in testing this model, the complication of the large energy deposited by a-particles in the substrate is avoided by irradiating films with - 6 keV X-rays, which deposit energies comparable to mips : Of course, there is 0168-9002/89/$03 .50 U Elsevier Science Publishers B.V . (North-Holland Physics Publishing Division)
a substantial difference between these absorption mechanisms : mips completely traverse the film depositing energy reasonably uniformly along their path, whereas 6 keV X-rays deposit their energy via the photoelectric effect . For Nb the absorption creates a - 3.4 keV photoelectron and an L-shell hole which creates, with very high probability, another 1 .93 keV Auger electron during de-excitation. The ranges [9] of the photoelectron and Auger electron are relatively small (- 0.07 and 0.03 wm, respectively), so that most or all of the energy, which is about five times the energy loss of a mip traversing a 1 [Lm thick film, is absorbed . Thus X-ray experiments can be used to test the model and evaluate the prospects for detecting mips. In this paper, experiments using - 6 keV X-rays are used to show that the switching probability is smaller than predicted by the simplest energy balance model [8]. This model is therefore in question . It is also shown that practical detectors of mips will greatly benefit from, and most probably require, superconducting transition temperatures, Tc, which are close to a low operating temperature (e .g., 4.2 K) . In addition, measurements of the evolution of the normal region after switching are shown to be in good agreement with the predictions of a thermal propagation model [10] . These results indicate that a self-recovering detector will be difficult to achieve. Finally, the sensitivity to soft X-rays is encouraging for the detector's use in high resolution X-ray imaging.
42 6
A . Gabutti et al. / Superconducting detector for minimum ionizing particles
2. Energy balance model
3. Experimental results
It is first necessary to review the switching principle and energy balance model [8] for mip interactions with a thin film. The ionization process produces an excess of high energy electrons as primary excitations in 3 X 10 -15 s: these come [111 into local equilibrium with the phonons, after a time of a few ps, at a temperature Te > Tbath. The system returns to Tbath by the much slower process of heat diffusion, eventually into the substrate. During this relaxation, a cylindrical hot spot, with T > T, will form around the mip's path . Starting from the initial primary excitations, the hot spot will grow to a maximum size determined by details of the thermal properties of the film and substrate, and subsequently decay. The hot spot destroys the superconductivity, thereby forcing the bias current to flow in a smaller cross section, where, if it can be made to exceed the critical value, the section of the film near the particle path is driven into its normal, nonsupefconducting state. The condition for the transition from the superconducting to normal state to occur for the full cross section of the film is : I> _ ],(w- 2ro )/w, where w is the film width and ro the maximum hot spot radius . Note that since the energy loss per unit path, dE/d x, is uniform throughout the film thickness, the normal region is cylindrical to a good approximation. Then ro can be evaluated using the energy balance relation : (H[TJ - H[ TbathI )1rro =dE/dx,
(2)
where H(T) is the temperature-dependent specific enthalpy of the superconducting material . In the absorption of low energy X-rays, the photoelectric effect predominates, and we assume that the primary excitations are created along the Auger and photoelectron paths. Although this energy is not deposited uniformly along their paths, the ranges are so short [9] that we will assume that the primary excitations are created in a small volume centered at the absorption point. Thus the hot spot will more closely resemble a sphere, rather than a cylinder, and a simple energy balance gives the radius from :
where i1 E is the Auger and photoelectron energy, which is - 5.33 keV in this case. In this spherical approximation, the switching threshold for the bias current becomes:
where t is the film thickness.
35 30 E
25
A
20
15
{H[TJ-H[Tbath]}41Tro/3=AE,
I >> IC (wt-ITrô )/wt,
For convenience in testing the simple energy balance model [8], the complication of the large energy deposited by a-particles in the substrate was avoided by irradiating films with X-rays . Our source of - 6 keV X-rays, from 55 Fe, deposits 5.33 keV, which is about five times the energy loss of a mil? traversing a 1 p m thick film . Thin-film detector strips of NbN and Nb were formed by a standard photolithographic process and patterned by either an isotropic plasma etch or a reactive ion etch (RIE) technique. They are 3.5 mm long, 0.2-1 .5 [Lm wide and 0.4-0 .8 p m thick, and have a number of electrical contacts . The NbN films were sputtered onto a sapphire substrate at a low substrate temperature of -- 200 ° C, resulting in T, = 11 .6 K and at 4.2 K, J~ -- 1 .8 x 10 5 A/cm2. Unfortunately the switching rate was barely measurable . This may have been because the cross-sectional area was too large (the widths and thicknesses of the NbN strips were - 1 .5 Wm and -- 0.8 pm, respectively) . However, in order to also eliminate the possibility that the random nature of the weak intergrain coupling of NbN resulted in anomalous weak spots (with very low J,) along the very narrow strips, pure Nb films, which exhibit no weak grain boundaries, were also studied. Note that except for radiation hardness, all the properties of Nb films should closely emulate NbN, including switching studies and testing of the energy balance model. Niobium films were sputtered onto an unheated substrate using the clean room environment at NIST to minimize contaminations which could also form weak spots along the strip. The Nb films had T, -- 9 K and at 4.2 K, Jc -- (2 .2-3 .3) X 10 7 A/cm2, if the geometrical cross sections found by SEM analysis are used . In fig. 1, Jc is plotted versus temperature for a representative sample . Unfortunately, the initial Nb strips, which were isotropically etched by a plasma in CF4 and OZ to widths narrower than 1 R m, showed excessive lateral etch and a triangular cross section . Starting from a 0.8 pm thick
(4)
Sample
10
AA
;1
5 0
4
Fig. 1. The typical temperature dependence of the critical current for Nb strips. Sample A, shown here, has a cross-sectional area of 0.1 win2 .
A . Gabutti et al. / Superconducting detector for minimum ionizing particles
Table 1 Summary of the properties of various Nb detectors and the effective normal zone radius, ro, determined from eq . (4) Strip A B C D
Etch
Shape
plasma plasma RIE RIE
triangular triangular rectangular rectangular
Area J~ [Wmz] [107A/cmz ] 0 .10 0 .18 0.08 0.24
3 .36 2 .2-2 .6 3 .8-5 .2 3 .5
ro [Itm] 0 .11 0 .07 0 .08
Nb film, the detector strip exhibited a maximum thickness in the range 0.3-0.4 R m. This was confirmed by both scanning electron micrographs and thickness profilometer analysis . This is a significant problem because it could lead to a rather strong variation in cross-sectional area along the length of the strip, thereby reducing the sensitivity of a majority of the film, except where the cross section is a minimum . To test this, electrical contacts were made at several locations along the strip, and Ic typically varied between 36 and 40 mA for four different sections of a strip that was - 0.8 p m wide at the substrate. A second set of Nb samples was patterned using an anisotropic reactive ion etch, for which the vertical etch rate was higher than the lateral one. Films of thickness - 0.4 wm exhibited rectangular cross sections with only a small spreading of the width very close to the substrate interface . Although these are expected to be more uniform, 1c varied between 30 and 42 mA for three different sections of a strip that was - 0.2 ~tm wide. Note that another strip patterned to a width of - 0.6 R.m by reactive ion etching had no measurable variation in I., indicating that the difficulties of maintaining uniform cross-sectional areas with our lithography techniques occur for ultrasmall dimensions. The properties of the detector strips used in these studies are summarized in table 1. All switching tests were done in liquid ° He at 4.2 K, and individual switching events were recorded on a multichannel counter, with the count rate converted to detector efficiency by dividing by the incident X-ray beam intensity. Fig. 2 shows the detector efficiency as a function of reduced bias current for three different Nb strips tested. Samples A and B had quasitriangular cross sections of area 0.1 and 0.18 Wmz, respectively, while sample C had a rectangular cross section of width 0.2 pm and thickness 0.4 ~tm for an area of 0.08 wmz. Measurements on the detectors shown in fig. 2 without the X-ray source indicated no switching events up to -99% of I, The distance from the X-ray source to detector is 3 ± 1 mm, leading to systematic uncertainties of up to 50% of the reported efficiencies for each sample . However, the highest efficiency for III, = 0.9 is only - 20% and thus five times smaller than an ideal detector .
42 7
30 25 20 v 15 5
0 .' ----» a a a, _"c',. ; 0 .8 0 .5 0 .6 0 .7 Fig. 2 . The switching probability, or efficiency, for three Nb strip detectors tested with - 6 keV X-rays. Specific data on these samples are indicated in table 1. The switching thresholds for X-rays in each sample shown in fig. 2 can be compared directly to the energy balance model [8]. Using eq . (4), the values of I/Ic of 0.58, 0.9 and 0.74 for samples A, B and C, respectively, lead to values of ro of 0.11, 0.07 and 0 .08 ~Lm. In eq. (3), ro is estimated to be 0.17 pm using values of H[Tbath] and H[T] for Nb taken from standard specific heat tables [12] . Therefore, although the experimental values of the effective ro determined from the switching threshold are reasonably the same, at - 0.1 lt m, for samples which had very different cross sections, they are consistently smaller than the predictions [8]. After switching, the hot spot can either collapse or propagate [10], depending on the power dissipated in the hot spot and the size of the hot spot relative to the
Fig. 3 . Digital oscilloscope traces of the voltage (upper trace) and current (lower trace) for switching with - 6 keV X-rays in sample A, which is biased at 33 mA . Note that the middle horizontal line is the zero for both I and V, and I is plotted as a negative quantity . The sensitivities, per large division, are: voltage, 0 .5 V ; current, 10 mA ; and time, 0 .2 Vs. The current was externally reduced at 0.4 lrs after the switching event.
42 8
A . Gabutti et al. / Superconducting detector for minimum ionizing particles
thermal healing length [13] . Although this property is independent of the switching efficiency, the possibility of a collapsing hot spot is of importance to the operation of a practical detector. To test this model [101, both the bias current, which was monitored by the voltage drop across a 50 0 series resistor, and the voltage across the strip were measured with a digital oscilloscope (having a 1 MS2 input impedance for the voltage measurement) . Typical results are shown in fig. 3 : a very fast risetime voltage pulse (-- 40 mV/ns) is observed across the strip, which is followed by a slower rise (- 1 [Ls) to the power supply voltage, Yo . The current decays during the first 0 .4 [Ls from its initial value of 33 mA to - 8 mA, after which the power supply reduces the current to a much lower value . Superconductivity recovers in - 1 ~ts, during which time the small residual current increases to - 2 mA and the detector is reset. For repetitive switching, the current is restored to 33 mA after 2 .5 ms . The thermal propagation model [10] makes very specific predictions about the switching, against which the decay of the current in fig . 3 can be checked without adjustable
parameters .
4 . Switching threshold The poor agreement of the calculated switching threshold with experiment no doubt results from various inadequacies of the simple model [8] . However, since several of these compensate each other, the net effect is, a priori, not simple to predict. For example, the relationships implied by eqs . (2) and (3) assume that the excess energy distributes itself uniformly within cylindrical or spherical volumes without extending beyond . This is clearly incorrect for thermal diffusion . However, in the real sample there will be a reduction of the critical current density, J., for any increase in temperature . Therefore, although the "normal" zone with T > T, will be smaller than ro , the effective radius will be larger than this "normal" zone due to the decrease in J. . A possible improvement in eqs . (2) and (3) is to replace Tc by Tc (I), obtained by inverting the relationship Ic (T). This gives the condition for the onset of additional heat generated by the transport current, but again suffers from assuming that the excess energy distributes itself uniformly . If H[Tc ] in eq . (3) is replaced by H[T.(I)], ro is found to be 0 .22 [.m for sample A which is larger than the above estimate and in poorer agreement with the model . Even larger values are found for the other samples . A calculation of the switching threshold, without such sweeping simplifications, requires numerical evaluation of the heat flow equation, self-consistently with the electrodynamics of the redistributed, nonuniform current flow . The latter is required to guarantee
that equipotentials are established across the film at locations far up- and downstream of the hot spot . Although such numerical evaluations are planned, for the present, the experiments must guide us, and eventually be checked against future calculations . A possible cause of the low switching efficiency is that some particle interactions result in self-recovering hot spots . These may result in a small, very fast voltage spike which does not trigger the counter nor reset the current in the present experiments . These could be studied experimentally if a greater voltage sensitivity for very short pulses can be realized . There are some additional considerations, not within the scope of the above numerical evaluations . These have to do with the possibility of nonequilibrium effects relating to the spatial uniformity of the sharing of energy between electrons and phonons [8,141 . Initially, because the energies are so high, the mean free paths and scattering times for each are short, and a good approximation to equilibrium is obtained at a high effective temperature. However, as the excess excitations cool down towards T, differences in the mean free paths and scattering times for the electrons and phonons may be important and the mean free paths become comparable with ro . Such effects can again be calculated in principle, but add an enormous burden, especially considering the different time scales involved . The next step beyond this paper is to compare experiments with the above, more modest numerical evaluation to see if the differences warrant further inclusion of nonequilibrium effects in the theoretical model . The shapes of the curves in fig. 2 as well as the efficiency are important factors in detector practicality . Both of these are no doubt influenced by several geometrical effects, but in ways which are not easily quantified in our preliminary experiments . The first is the nonuniformity of Ic along the length of the strip due to the inevitable variations of the cross-sectional area . Although these variations cannot be adequately measured in the present experiments, their effect, together with the potential for nonuniform material properties (e .g ., in granular NbN), should be less significant if larger areas (i .e ., wider films) are used. The second geometrical effect results from the fact that the effective diameter of the hot spot, 2ro , is close in magnitude to the sample dimensions, so that edge effects must be considered . For example, if the X-ray energy is deposited near a surface, the hot spot will be hemispherical, but will have the same volume as the spherical hot spot . As a result, the effective hemispherical radius will be larger, but the effective fraction of the film's crosssectional area will be smaller . For X-rays this corresponds to a reduction in area of 63% and for mips interacting near the film edges it is 71% . Obviously, this reduction will depend on the exact location of the interaction with respect to the film's boundaries, e.g .,
429
A . Gabutti et al. / Superconducting detector for minimum ionizing particles
corners are worse. However, this effect will also be significantly less important for large area strips (i .e ., wider films) . The third geometrical effect, occurring only for Xrays, is the escape from the film of either the photoelectron or Auger electron, which have ranges of 0.07 and 0.03 wm, respectively, in this experiment . Therefore part of the energy can occasionally escape for our small geometries, resulting in smaller values of AE which will lead to higher threshold currents for these interactions . However, films with larger cross-sectional areas will again minimize this effect. It is very important to realize that wider films can alleviate the geometrical effects which contribute to the low values of efficiency. Referring to eqs. (1)-(4), it is clear that a smaller value of enthalpy difference is required for wider films, for a given dE/dx or i1 E. As a result, practical detectors of mips will greatly benefit from, and most probably require, superconducting transition temperatures, Tc, which are close to a low operating temperature (e .g ., 4.2 K) .
Since the values of pNJzd are always much greater than 2WO = 310 W/cm2 , the thermal propagation velocity is given simply by [10] : v= [PNJ2dK/Cvd(Tc-
Tbath)1 1/2 .
(5 )
Note that v is proportional to J, so it will decrease as the hot spot grows and J is reduced. Values for C,, and K for these films can be taken from tables, using the Wiedemann-Franz law and the measured pN to determine the latter . This procedure gave excellent agreement for previous [10] propagation velocity measurements. In order to determine I(t), we first calculate the time dependence of the strip resistance, R(t) . Recall from eq . (5) that v(t) is proportional to J(t) and, if inductive effects are neglected, the simple circuit requires J(t)/J(0)=1/(1+R(t)/RO), where R o is the 50 0 series resistor used to measure the current. Then : dR(t)/dt=2vRN
(6)
1+R(t)/Ro
is trivially solved to give : I(t)/j(o) _ {1/(1 + at) ) 1/2,
5. Thermal propagation The thermal propagation model [10] makes very specific predictions about the switching, against which the decay of the current in fig. 3 can be checked without adjustable parameters. In the following it is assumed that the hot spot propagates in the film alone, and the substrate acts only as an additional thermal reservoir, at the bath temperature, Tbath, to the actual 4 He liquid covering the film. It is then possible to estimate the propagation velocity of the normal-superconductor boundary using: the total thermal boundary conductance to the substrate and liquid 4 He ; the specific heat of Nb at Tc, C ; the thermal conductivity of Nb at Tc, K; the film thickness, d; and the power dissipation per unit film area, pNJzd, where pN is the normal state resistivity at T= 10 K. The boundary conductances, Wo , reflect the power per unit area necessary for a temperature rise of Tc - Tbath in the film . Measurements show a value of 125 W/cm2 for NbN/sapphire substrate boundaries [10] and 1 W/cm2 for Nb/liquid 4 He boundaries [15] . Little change is expected for Nb/sapphire boundaries, but the conductance of 4 He boundaries is known to be time dependent because of the formation time for nucleate boiling at the surface [15] . As a result, a value of 30 W/cm2 is appropriate [15] for times less than -- 10 ws, and the total Wo is 155 W/cm2 . Using the experimental values of pN = 1.9 ltQ cm and d = 0.5 ~.m, one finds pNJz d= 8 x 10 4 W/cm2 immediately after switching in fig. 3 when J = 3.3 x 10 7 A/cm2 , and pNJzd = 0.47 x 10 4 W/cm2 after 0.4 [.s when J = 0.8 x 10 7 A/cm2 .
(7)
where a = 2 v [R N/l ]/R o ; R N is the normal resistance (at T = 10 K) of the detector of length 1; and the coefficient 2 reflects the propagation of the hot spot in both directions from the particle interaction. It is easily shown that V(t)/VO = 1 - I(t)/I(0). The results of this calculation, which uses no adjustable parameters, are shown in fig. 4 and obviously correctly reproduce the overall features of the experimental data shown in fig. 3, including the reduction of I by a factor of 4 after 0.4 Ws . Although the ringing, caused by inadequate terminations, makes the interpretation difficult, the voltage invariably seems to show a much 1 "-- .-.,0 .8
Sample A
I~ "
0 .6
UO~I`
"
V(t)1Vo
0 .4
I(t)/I(0)
O 0 .1
_- _ 0 .2
-.
0 .3
time(microsec)
0 .4
Fig. 4. Calculated time dependence for the reduced current and voltage using the thermal propagation model. This calculation uses only experimentally measured parameters corresponding to sample A shown in fig. 3, and no adjustable parameters. Vo =1 .65 V and I(0) = 33 mA .
43 0
A. Gabutti et al. / Superconducting detector for minimum ionizing particles
faster rise and a peak at about 20 ns . The peak is not explainable by the thermal propagation model and is thought to result from stray reactances and/or an inevitable ground loop in the measuring circuit. However, two effects could result in a faster risetime for the voltage. Although inductance in the circuit could delay the current decrease relative to the increase in R(t), calculated inductance values yield a negligible effect . A more plausible explanation is that the initial power density of pNJ zd = 8 x 10 4 W/cmz causes a temperature rise considerably in excess of T, so the resistance is larger than that determined above using the lowtemperature pN . The room-temperature/low-temperature resistance ratio was 6. Future studies will address this point. It should also be noted that the thermal recovery of the hot spot was also observed in fig. 3 when, after 400 its, the current was reduced, but not to zero . After superconductivity was recovered in - 1 ws, the current of 2 mA corresponded to p N J zd = 294 W/cmz, which is just less than 2WO = 310 W/cm z. Thus the thermal propagation model predicts negative velocities for I < 2 mA, and experimentally both I and V of the hot spot are seen to recover in fig. 3. Note that the maximum size of the hot spot shown in fig. 3 is - 1000 wm, assuming the low-temperature PN . Using the thermal propagation model, it is a simple matter to determine the conditions necessary for a self-recovering hot spot. The thermal healing length, rl, must be evaluated for Nb films on sapphire substrates, using the values of WO appropriate to the fast times expected, where A value of 1 .3 gm is found for 0.4 l.Lm thick Nb films. The critical hot spot length, L, below which the hot spot collapses, is given by [10] :
z L- =2 t1WolPNJ d, 9 which is only 5 nm! It will be difficult to increase Lnc sufficiently (to >_ 0.2 ~tm) to be practical for a selfrecovering hot spot . However, decreases in PN, which also increases in thermal conductivity K, would have to be only -- 12 . This is attainable in pure crystals of Nb, but not necessarily in such thin films, and certainly not in disordered materials like NbN which are radiation hardened . Reducing J, and therefore necessarily J., in order to reach the switching threshold, by a factor of 6 will be sufficient for Nb, but for radiation hardened materials there must be an additional reduction in J to overcome the larger PN and lower K.
6. Conclusions The main conclusion of this research is that a more detailed model calculation of switching thresholds will be necessary to adequately analyze the results of switching experiments for detectors of minimum ionizing radiation. It also seems clear that practical detectors of mips will require superconducting films with low Tc so that wider films, which are less susceptible to geometrical effects, can be used . The thermal propagation model seems to adequately describe the overall behavior after switching, but it implies that a self-recovering switch will be very difficult to achieve. In addition, the present results already show promise for the development of a very high spatial resolution soft X-ray detector . Acknowledgements The authors thank Rex Craig for help during sample preparation. This work was supported by the US Department of Energy, Divisions of High Energy Physics and Basic Energy Sciences-Materials Sciences under contract #W-31-109-ENG-38, and the National Institute of Standards and Technology . References [1] E.g., Superconductive Particle Detectors, ed . A. Barone (World Scientific, Singapore, 1988). [2] P. Gregshammer, H.W . Weber, R.T . Kampwirth and K.E. Gray, J. Appl. Phys . 64 (1988) 1301 . [3] R .G. Wagner and K.E . Gray, in ref. [1], p. 204. [4] N.K . Sherman, Phys . Rev. Lett . 8 (1962) 438. [5] D.E. Spiel, R.W . Boom and E.C . Crittendon, Jr ., Appl . Phys . Lett . 7 (1965) 292 ; E.C . Crittendon, Jr . and D.E . Spiel, J. Appl . Phys. 42 (1971) 3182.
N. Ishihara, R. Arai, T. Kohriki and N. Ujiie, Jpn. J. Appl . Phys. 23 (1984) 735. K.W. Shepard, W.Y . Lai and J.E. Mercereau, J. Appl . Phys . 46 (1975) 4664 . [8] K.E . Gray, in ref. [1], p. 1 . [91 K.-H. Weber, Nucl . Instr. and Meth . 25 (1964) 261. [10] K.E. Gray, R.T . Kampwirth, J.F . Zasadzinski and S.P. Ducharme, J. Phys . F: Met. Phys . 13 (1983) 405. [11] R.W. Schoenlein, W.Z . Lin, J.G . Fujimoto and G.L . Eesley, Phys. Rev. Lett. 58 (1987) 1680 . [12] Y.S . Touloukian and E.H .Buyco, Thermophysical Properties of Matter (Plenum, New York, 1970) p. 153. [13] W.J. Skocpol, M.R. Beasley and M. Tinkham, J. Appl. Phys. 45 (1974) 4054 . [14] N .E . Booth, in ref. [1], p. 18 . [15] W.G. Steward, Int. J. Heat Mass Transfer 21 (1978) 863. [6]