498
WILHELM GODSCHALK
Mathematical Formulation of Rotor Deceleration Experiments in Ultracentrifugation
by Wilhelm Godschalk Department of Biochemistry and Department of Physics, Unhersitv of Virginia, Charbttesville, Virginia 22001 (Received June 6 , 1967)
Equilibrium centrifugation of large particles is only possible when the angular speed of the rotor is kept uniformly constant. The degree of accuracy required depends on the relaxation time of the system, which is defined in this paper, and ultimately on the diffusion coefficient of the macromolecule. A decrease in rotor speed will generally result in a reversal of the net flow in the cell; a t sufficiently high rates of deceleration, the flow becomes completely diffusion controlled. The flow equation for this case has been solved for an ideal two-component system, using a special form of a finite Hankel transform. The exact solution for the concentration as a function of radial distance and time is obtained as an infinite sum of composite Bessel functions. Convergence has been tested and numerical solutions have been evaluated on the digital computer for simulated runs on ribonuclease and turnip yellow mosaic virus. This method opens the possibility of measuring diffusion coefficients from the slope of a logarithmic plot of the integral transform of the concentration us. time. It is expected that this treatment can be extended and used to advantage for the study of multicomponent systems.
Introduction The original model of the magnetically suspended ultracentrifuge, developed by Beams and coworkers,'J has been in use in our laboratory for several years. A well-known property of this type of centrifuge is the slight decrease in rotor speed through the absence of an external drive mechanism. The rate of deceleration depends on the residual air pressure in the evacuated chamber and the occurrence of stray currents when ferromagnetic material is present in the vicinity of the centrifuge. Deceleration rates as low as O.l%/day have been found, but, under less ideal conditions, this figure can easily go up to several per cent per day. Quite recently, a magnetic drive system has been devised by Beams,3which is capable of keeping the rotor speed constant with suitable precision during the run, without making any mechanical contact with the freely suspended rotor. I n spite of this new and very useful development, a freely coasting rotor is still employed for some applications. The beneficial effect of deceleration is the time-saving factor in reaching equilibrium. This effect has been known to us as an experimental fact for some time but has been put on a theoretical basis only recently by Billick, et d 4 9 6 For molecules of small or medium size a state of quasi-equilibrium can be maintained experimentally, in spite of the continuous loss of angular speed. I n this case, the system can respond instantly to the changes in centrifugal field, as relatively small molecules have The Journal of Physical Chemistry
large diffusion coefficients. With large particles, like the plant viruses we are studying (molecular weight >lo8), complications arise when the diffusion coefficients are small and the speed loss is appreciable. The response of the system in this case will lag behind the decrease in centrifugal field; this means that the flow in the cell is controlled by diffusion only; Le., the system can be treated as if there were no sedimentation at all. Such an experiment can yield additional information about diffusion coefficients and might be very useful for the determination of the composition of paucidisperse systems. Accordingly, it became desirable to develop a theory for a deceleration experiment of this kind. I n this paper, is presented the simplest case of a monodisperse ideal system, where the diffusion coefficient of the macromolecule is independent of concentration and pressure. The partial differential equation is solved by means of a special integral transform. Though other approaches are possible, the use of this transform was preferred because of its potential for a uniform treatment of related problems.6 (1) J. W. Beams, R.D. Boyle, and P. E. Hexner, Rev. Sci. Instr., 3 2 , 646 (1961). (2) J. W. Beams, R. D. Boyle, and P. E. Hexner, J . Polymer Sci., 57, 181 (1962). (3) J. W. Beams, Rev. Sci. Instr., 37, 687 (1966). (4) I. H. Billick, M. Dishon, M. Schultz, G. H. Weias, and D. A. Yphantis, Proc. Natl. Acad. Sci. U.S., 56, 399 (1966). (6) I. H. Billiok, M. Sohults, and G. H. Weiss, J . Phya. Chem., 71, 2496 (1967).
ROTORDECELERATION EXPERIMENTS IN ULTRACENTRIFUGATION The system to be considered in the following derivations is thought to be at equilibrium at an angular speed w a t the start of the experiment. At zero time, the rotor starts to decelerate at such a rate that the response of the system is slow by comparison. This situation is similar to the one treated by Klenin, et UZ.,’ where the speed was thought to decrease instantly to a lower value. In the present case, this lower speed is assumed to be zero. The concentration distribution in the cell as a function of time will be evaluated through the use of exact solutions to the diffusior, equation, rather than the rectangular approximation.s
The Integral Transform A function, f , of variable r is assumed to possess an integral transform.&), defined by b
.J(P)
=
J f(dBo(pr)r dr
499
and p2 represent two different roots of eq 4,the integral equals zero. On the other hand, if pl and p , designate the same root p , the value of the integral can be found by applying 1’Hbpital’s rule to eq 6. First, eq 6 is differentiated with respect t o p 2 , then p 2 and p l are both put equal to p . The resulting expression is evaluated using the well-known recurrence formula for Bessel functions
pl
zJn’(s)- nJ,(x)
(7)
= -XJ,+1(X)
and the Wronskian relation
With simple, though tedious algebra, the value of the integral is found to be
(1)
in which
Bo(pr) ==Jo(pr)Yo’(pa) - Yo(pr)Jo’(pa)
(2)
The terms Jo and Yoare Bessel functions of zero order and of the first and second kind, respectively. We define further
(3) The parameter p represents the positive roots of the equation
Jo’(pb)Yo’(pa) - Yo’(pb>Jo’(~a) = 0
(4)
with Jo’(pb) and Yo’(@)defined as in eq 3, for r = b. Although this type of finite Hankel transform is known,O this particular specimen, with p as roots of eq 4,was not found in the literature. However, the inverse transform can be found from the theory of Fourier-Bessel series. For two values, p l and p2, of parameter p , the following relationship ho1ds”J for any cylinder function B, JrB,(p14B,(p~!4r dr =
When the function Bo is substituted for B , and the integration is performed between a and b, eq 5 is equivalent to (p2
- p l Z ) s b B , ( p ~ r ) D ~ ( p 2 rdr) r = b [piBo(ph)Bo’(pib)- ~zBo(pib)B0’(~2b) 1a [piBo(pza)Bo’(pia>- ~zBo(~ia)Bo’(pza) I (6)
As Bo’(pb) and Bo’(pa) are both zero, the right-hand member of thia equation vanishes. Accordingly, if
I f we now assume that a function can be expanded in a series
+ CamBo(pmr) m
f(r) = K
(10)
m=l
in which K is a constant and the pmparameters designate all the positive roots of eq 4,the coefficients a, can be found by multiplying eq 10 by rBo(p,r) and integrating between a and b l b f ( r ) B o ( p , r ) r dr = K
s,”.
o (p,r)r dr
+
)
LbBo(prnr)r(g p ~ o ( p m r ) dr
(11)
The left-hand member is nothing else butf(p), the transform of f(r), which was defined in eq 1. The first integral on the right equals zero, as can be easily verified by integration of the recurrence formula
XB,‘(X)
+ nB,(z) =
XB,-l(X)
(12)
and using eq 4. When the integration of the second term on the right is carried out, all mixed integrals involving two different values of p , vanish, so that only one integral remains, the value of which is given by eq 9. Accordingly, a typical coefficient in eq 10 is
and the inversion formula for the transform becomes (6) W.Godschalk, unpublished data. (7) 8.Klenin, H. Fujita, and R. A. Albright, J. Phys. Chem., 70, 946 (1966). (8) D. A. Yphantis and D. F. Waugh, {bid., 60, 623 (1956). (9) E.C.Titchmarsh, Proc. London Math. Soc., 2 2 , 13 (1923). (10) G. N.Watson, “A Treatise on the Theory of Bessel Functions,” 2nd ed, Cambridge University Press, London, 1945,p 134. Volume 76, Number 2 February 1968
WILHELM GODSCHALK
500
Now Bo(pr)is a linear combination of Bessel functions and is itself a solution of the Bessel equation
(14) The Diffusion Problem Consider a macromolecule with sedimentation coefficient s and diffusion coefficient D, in equilibrium at rotor speed w, in a sector-shaped cell with inner radius r = a and outer radius r = b. Both s and D are treated as independent of the concentration C. The condition of zero flow in the cell, plus the consideration that a t the hinge point” r h the concentration is equal to Co, the original uniform concentration of the solution, yieIds the initial concentration distribution in the cell Co = Co exp[yMw2(r2 -
rh’)]
(15)
S
1 - Bp
2DM
2RT
y=-=-
c = A exp[--Dpzt]
(23)
The initial condition, eq 15, can be transformed to yield the coefficient A
A = 0 = Coiexp[yMw2(r2- rha)]Bo(pr)rdr
where 8 is the partial specific volume of the solute and p is the density of the solution. The location of ?+h can easily be determined experimentally, by means of a white-light fringe, when a compensation cell is used in the interferometer system, as suggested by Beams.2J2 When other optical systems are used, the equilibrium distribution can be expressed in any other suitable way, without affecting the treatment given in this paper. Starting from equilibrium, the rotor is decelerated at such a rate that sedimentation becomes negligible. The partial differential equation which has to be solved reduces then to Fick’s second law for cylindrical diffusion.
While the initial condition is given by eq 15, the final condition, as t approaches infinity, is uniform concentration Co throughout the cell. The boundary conditions to eq 17 can be formulated as
because the flow is zero at the ends of the cell and is proportional to the concentration gradient in the absence of a sedimentation force. Application of the transformation of eq 1 to the partial differential equation eliminates the space variable r
lb
g B o ( p r ) r dr
+ s,” gBo(pr) dr
(19)
in which c is the transform of C. Using eq 4 and integrating the first term on the right twice by parts, we obtain
D I_ dt
with the solution
b
in which
1 dC D -dt =
Substituting eq 21 into eq 20, we obtain the simple differential equation
= pLbCBo’(pr)dr
The Journal of Physical Chemistry
+ p21bCBo”(pr)r dr
(20)
(24)
With the inversion formula (14) and substituting Co for the constant K , the concentration in the cell as a function of radial distance and time becomes
BO(P,~>C~ exp[--Dpm2t1 (25) in which the transform of the initial dishribution, Co,is given by eq 24. Equation 25 meets the requirements for a good solution, as it satisfies all the conditions that were imposed. At very large t, C becomes Co throughout the cell. On putting t = 0, the formula represents the inverse transform of the transform of the initial condition, i.e., the initial condition itself. Differentiation of eq 25 with respect to r shows that it satisfies the boundary conditions. Finally, multiplication of eq 25 by r and integration between a and b gives the expression of conservation of mass
s,”Cr dr = 1/zCo(b2- az)
(26)
Computer Analysis of the Results A program in extended ALGOL was written for the Burroughs B5500 digital computer, to evaluate the terms of the series in eq 25. The roots p , of eq 4 are only dependent on the length of the liquid column in the centrifuge cell. The Bessel functions were calculated on the computer for an actual cell, with a = 6.4228 cm and b = 6.7344 cm. A plot of the cylinder function Ba’(pb)us. a continuously variable parameter p was generated by means of an appropriately programmed Calcomp plotter. The function is single valued and has a half-period of about 10 for these values of a and b, as is shown in Figure 1. The approximate values of the zeros, read from this graph, (11) W. J. Archibald, J. Phys. Colloid Chem., 51, 1204 (1947). (12) J. W. Beams, Rev. Sci. Instr., 34, 139 (1963).
ROTOR DECELERATION EXPERIMENTS IN ULTRACENTRIFUGATION
501
few minutes of deceleration. However, experiments suggest that it does not take more than 10 min for the gradients a t the boundaries to drop to zero, as eq 18 requires, provided the speed loss is sufficiently great. Numerical solutions of eq 25 have been obtained for simulated runs on ribonuclease ( M = 1.36 X lo4, D = 1.19 X and turnip yellow mosaic virus (TYMV) ( M = 5.4 X lo6, D = 1.51 X lo-'). For intervals over 4 min, 20 terms were ample for the summation in the case of TYMV; for t exceeding 4 hr, the concentration distribution was almost entirely determined by the first two terms. With ribonuclease, even fewer terms sufficed. The results of the calculations, as direct reproductions of the Calcomp plotter output, are presented in Figures 2 and 3. From the curves it is obvious that the plot of (l/rC)(dC/dr) us. T , as suggested by Archibald," is the most sensitive to changes in rotor speed. All three ways of representing the data show that the effect of deceleration is most pronounced near the boundaries, which is the expected behavior. The larger diffusion coefficient of ribonuclease is reflected in a faster approach to a uniform concentration, as is apparent from Figure 2.
Relaxation Times
Figure 1. The cylinder function Bo'(pb) as a function of p, for b = 6.7344cm.
were used as starting values in a computer routine (BO020 of the University of Virginia Computer Science Center), which makes use of the iterative Muller procedurela for finding the roots of any function. The first 100 positive roots p , of eq 4 were obtained to a precision of eight figures. The required processing time was about 40 min; the results are universally applicable, as long as a and b remain unchanged. The values of pm were plugged into the main program,l* which evaluated the terms of eq 25 including the integrals and carried out the summation and the plotting. Once the values of p were known, it took an additional 40 min of processing time to calculate and add up 100 terms, at 61 T values in the cell, for 11 different time intervals. I n addition to the concentrations, values for the concentration gradient at each point were obtained by analytical diflerentiation of the terms of eq 25. As these derivatives showed a tendency to oscillate for small time intervals, numerical differentiation of the C 11s. T plot was also employed, with notably better results. For small t values, the series of eq 25 needs many terms to converge properly. This is no serious drawback, as it is at least questionable whether the boundary conditions are even approximately valid during the first
The deceleration experiment can be treated as a relaxation phenomenon. Assume a system at sedimentation equilibrium; the concentration distribution in the cell is given by eq 15. This initial distribution will be denoted by Co. At time zero the system starts to relax toward a final condition at w = 0, where the concentration is uniformly equal to Co, the original bulk concentration of the solution. A relaxation time r for the whole process can now be defined in the usual way
c - co
c*- co
=
(27)
From Archibald's famous paper," we know that for a monodisperse, ideal system the hinge point is constant during sedimentation. That it remains constant on deceleration can easily be shown in the following way. A combination of eq 15 and 27 yields
for any intermediary time 1; r is not necessarily constant, but may be a function of t. By equating C for an arbitrary value of t , with Co, at time t = 0, we can find the value of r at which the concentration is at all times equal to the initial equilibrium concentra(13) W. L. Frank, J . Assoc. Computing Machinery, 5 , 154 (1958). (14) This program, named DECEL, oan be obtained from the author. Volume 72,Number 2
February 1968
502
WILHELM GODSCHALK
@I d
si d
Figure 2A. Relaxation of a 0.4 g/100 ml solution of ribonuclease, starting from sedimentation equilibrium at 300 rps. concentration 08. radial distance.
Figure 2B. See caption for Figure 2A. (l/TC)(dC/dr) ua. r; this function is a horizontal straight line a t equilibrium ( t = 0). The Journal of PhysCal Chemistry
503
ROTOR DECELERATION EXPERIMENTS IN ULTRACENTRIFUQATION
El
'1 17MIN 4 S E C
;I : (P.
4 H R S 33MIN 4 S E C
8: 0
9.
l8HRS 12MIN I 6 S E C
z
4: I.
w.700
Figure 2c. See caption for Figure 2A.
w.200
Ln C us. r3; this plot is linear at equilibrium.
17MIN 4 S E C
I H R BMlN
I L 1 I I \ Y
T d 1.111,
1 4
TU&"
291 HRS 16MIN 16 S E C
!7. 5!
:
6.Y20
I
6.W
6.W
6.SYO
.
6.580
,
,
,.
6.620
RROIRL OISTRNCE
, ,
I
.
6.660
.
,
, .
6.700
,
. , . 6.N
7 6.
(CM)
Figure 3A. Relaxation of a 0.4 g/100 ml solution of TYMV, starting from sedimentation equilibrium at 14 rps. The initial concentration distribution is approximately equal to that of Figure 2. The time intervals are the same. C us. r. Volume 78, Number 2 February 1068
WILHELM GODSCHALK
504
17MIN 4SEC kHR 8 M I N I G S E C
d
4HRS 3 3 M I N 4 S E C
\ \I 72HRS 4 9 M l N 4 S E C
RRDIRL DISTRNCE
(CM)
Figure 3B. See caption for Figure 3A. (l/rC)(dC/dr) vs. r.
l7MlN 4 S 9
Figure 3C.
See caption for Figure 3A.
The Journal of Physical Chemiatra,
Ln C vs. rt.
ROTOR DECELERATION EXPERIMENTS IN ULTRACENTRIFUQATION tion at that point. we get
exp[
-$]I>
For t = ti and t = 0, respectively,
= Coexp[yMw2(r2-
rh2)]
(29)
7
It is obvious that this only holds at the point r = r h , as has been chosen different from zero. To get an idea of the order of magnitude of the relaxation time for a real system, we can write eq 27 in the following form, using eq 25 tl
2PBo(pmr)pee
l
DpmPt
2 PBo(pmr)p
=
e--t/7
(30)
lll-1
in which P represents the first factor under the summation sign in eq 25. From eq 30 it appears that r is some “mean value” of the variable 1/Dpm2. Statistically, the exponential represents the most probable decay frequency of the spectrum. Qualitatively, we can immediately infer that 7 must be of the order of lo4 sec for macromolecules such as ribonuclease or TYMV. For ribonuclease,’b with D = 1.19 X the largest possible value of r is 8.4 X loa, as the first positive root of eq 4 is approximately 10 for a normal centrithe fuge cell. For TYMV,l* with D = 1.51 X upper limit of 7 is 6.6 X lo4. As t approaches infinity, r tends to this limit, because for large values of t the sum in the numerator of eq 30 is almost entirely determined by the first term of the series.
Discussion The treatment given in this paper is only applicable at relatively large deceleration rates. If we define a relaxation time, r,, for the rotor deceleration
in which wo is the initial rotor speed, then the above treatment will break down whenever r and 7, become of comparable magnitude. In terms of the two examples cited, this means that, starting at 14 rps, a TYMV run has to be decelerated more than 2.7%/ hr, or 48%/day. Ribonuclease, with its larger diffusion coefficient, requires higher deceleration rates, better than 19.3%/hr, or 99.5%/day. These rates are considerably higher than the normal speed loss of a freely coasting rotor in the magnetic ultracentrifuge. However, with the recently developed magnetic drive system,8 it is possible to introduce artificial deceleration by operating the rotating magnetic field in reverse. With mechanical centrifuges, as the Spinco Model E, the brake system can be used, but the
505
experimenter should be continuously aware of the possibility of convection in the cell. The mathematical treatment given in this paper is applicable to the field relaxation method of Kegeles and Sia.17 The deceleration rates used by these authors were as high as 5%/min for chymotrypsin. Although their method did not involve sedimentation equilibrium, the present theory can be used without extensive modifications. In this paper, the equilibrium distribution has been chosen as a well-defined initial condition. However, in eq 25 the initial condition is only represented by its integral transform P , so the equation can be made valid for an arbitrary initial concentration distribution by making Coequal to the numerical value of the integral transform of this actual distribution. As this relaxation method is the reverse of approach to equilibrium, it can be employed for the determination of diffusion coefficients. For small time intervals, it is advantageous to use the expression for the transform of the concentration in the form of eq 23 = coe-D~m2t
(32)
Again, Co may represent any arbitrary initial distribution. Thus, for any eigenvalue p,, the slope of a plot of In C vs. time is -Dpm2, while the intercept yields In 0. The transform of C can be calculated by numerical integration; this is most conveniently done on the computer. The same program generates the appropriate Bessel functions which are needed for the calculation of the transform. If desired, the plotting and evaluation of D can be carried out in the same operation. For very large values of t, eq 25 can be used. As t increases, finally all terms higher than the first become negligible, and D can be obtained from the limiting slope of a plot of In C vs. t at any fixed point in the cell. In general, the time required to reach this state will be excessively long, but the value of D obtained in this case corresponds to conditions of very small concentration gradients and might be different from that resulting from the first method or from approach to sedimentation equilibrium.18 Preliminary measurements have been performed in our laboratory and will be reported elsewhere. It is expected that this treatment can be extended to multicomponent systems and will provide a useful tool for the study of viruses and their protein subunits. Acknowledgments. Thanks are due to Professor J. W. (15) J. T.Edsall, Proteins, 1, 549 (1953). (16) R. Markham, Discussions Faraday Soc., 11, 221 (1951). (17) G.Kegeles and C . L.Sia, Biochemistry, 2,906 (1963). (18) K. E. Van Holde and R. L. Baldwin, J . Phys. Chem., 62, 734 (1958). Volume 7.9, Number 8 February 1068
506
M. KRUMPELT, J. FISCHER, AND I. JOHNSON
Beams for his kind interest and for providing the instrumentation which made this work relevant. Drs. D. W. Kupke and G. K. Ackers and Mr. G. T. Mayer were very helpful as discussion partners. The author is especially grateful to Dr. G. H. Weiss of the National Institutes of Health for making two manuscripts
available prior to publication and for solving the problem under study in the rectangular approximation. This work could not have been done without the outstanding service of the University of Virginia Computer Science Center. This work was supported by U. S. Public Health Service Grant No. GM-12569.
The Reaction of Magnesium Metal with Magnesium Chloride by Michael Krumpelt, Jack Fischer, and Irving Johnson Argonne National Laboratory, Argonne, Illinois 60489
(Received June 6 , 1967)
The reaction of metallic magnesium with molten magnesium chloride has been studied by the determination of the magnesium(1) content of molten MgC12 in equilibrium with liquid Mg and Mg-Cu alloys at 800". An almost linear dependence of the Mg(1) content of the molten NlgCl2 on the magnesium activity in the alloy was found. This was taken to confirm the existence of Mg22+as the solute species. The small deviations from exact linearity have been discussed in terms of a concentration dependence of the activity coefficient of the dissolved magnesium. The activity coefficient was found to decrease with increase in concentration. The temperature dependence of the Mg(1) content of the MgC12 in equilibrium with Mg was determined between 750 and 850".
Introduction It has been suggested' that solutions of the alkaline earth metals in their dihalides consist of an equilibrium mixture of monovalent ions, divalent ions, and F-centerlike electrons. The results from studies of the electrical conductivity have been used to deduce the relative concentrations of monovalent ions and electrons in a few cases. However, the cryoscopic and phase-diagram data have not established unambiguously whether the monovalent ion exists as the monomer, M +, or the dimer, M22+, or a mixture of the two. I n the present study, the latter question has been investigated for the MgMgClz system by measuring the magnesium metal solubility in molten magnesium chloride which is in equilibrium with various liquid magnesium-copper alloys. The term magnesium metal solubility, while used by most workers in this field, is misleading, since the liquid salt solution very likely does not contain a significant concentration of metal atoms. I n the present report, the magnesium content of Mg-MgCl2 mixtures as determined by the measurement of the quantity of hydrogen released by reaction with aqueous HC1 will be referred to as the Mg(1) content. The over-all reaction for the equilibration of liquid magnesium or liquid magnesium-copper alloys with The Journal of Physical Chemistry
molten magnesium chloride may be described by one, or possibly both, of the equations Mg Mg
+ Mg2+
2Mg+
(1)
+ Mg2+I_ Mgz2+
(2)
If eq 1 represents the reaction, then the activity of Mg+ in the molten salt would be directly proportional to the square root of the magnesium activity in the alloy, whereas, if eq 2 represents the reaction then the activity of Mg22+would be directly proportional to the magnesium activity in the alloy. If it is assumed that the activity of Mg+ or Mg22+is directly proportional to the Mg(1) content of the magnesium chloride, then the applicability of eq 1 or 2 may be established. Three previous attempts were made to determine the parabolic or linear d e ~ e n d e n c e . ~ -The ~ first two groups of investi(1) A. 9. Dworkin, H. R. Bronstein, and M. A. Bredig, J . Phys. Chem., 70, 2384 (1966). (2) N. G. Bukun and E. A. Ukshe, Russ. J . Inorg. Chem., 6 , 466 (1961). (3) P.S. Rogers, J. W. Tomlinson, and F. D. Richardson, "Physical Chemistry of Process Metallurgy," Interscience Publishers, Inc., New York, N. Y.,1961,p 909. (4) J. D. Van Norman and J. J. Egan, J . Phys. Chern., 67, 2460 (1963).
