J. Phys. Chem. B 1999, 103, 1393-1401
1393
ARTICLES Qualitative and Quantitative Analysis of Solid State Free Induction Decay (1H NMR) Curves Using a Combination of the Methods of Gardner and Prony: Isotactic Polypropylene as a Case Study S. Schreurs and J.-P. Franc¸ ois* Limburgs UniVersitair Centrum, Institute of Materials Science, Department SBG, Research Group Analytical Chemistry, UniVersitaire Campus, B-3590 Diepenbeek, Belgium
P. Adriaensens and J. Gelan Limburgs UniVersitair Centrum, Institute of Materials Science, Department SBG, Research Group Organic and Polymer Chemistry, UniVersitaire Campus, B-3590 Diepenbeek, Belgium ReceiVed: July 10, 1998; In Final Form: NoVember 3, 1998
A combination of the methods of Gardner and Prony is used for the qualitative and quantitative analysis of solid-state free induction decay (1H NMR) curves of isotactic polypropylene, recorded at 40, 60, 80, and 100 °C at constant time intervals. These multicomponent decay curves contain information about the physicochemical structure of the polymer sample: each component is characterized by a specific model function containing a spin-spin relaxation time T2 and an intensity as parameters depending on the properties and amount of the corresponding phase. The number of phases is still an unanswered question. Each FID curve is analyzed by using an improved method of Gardner et al. in order to reveal the number of phases, to predict the best fitted model function, and to estimate the T2 value of each component. The quantification of the parameters (T2 and fraction) with Prony’s improved method is discussed. The presence of an intermediate phase is unambiguously demonstrated. The influence of temperature on the physicochemical structure of the polymer is examined. The precision of the computed parameters is discussed for the analyses of three consecutively measured FID curves of isotactic polypropylene at 60 °C.
Introduction Solid state proton wide-line (broad-band) NMR provides a sensitive probe for the molecular state of a nuclear environment through the short range nature of magnetic dipolar interactions. For heterogeneous polymer systems (physical or chemical), the molecular mobility of each component and the fraction of the component can be directly estimated from the measurement of the spin-spin relaxation time T2.1 A solid echo technique provides a way to collect the complete response of the system to a 90° pulse near the resonance frequency without influence of the system recovery time. The recorded free induction decay (FID) curve is a multicomponent decay curve having for each component a characteristic line shape, a spin-spin relaxation time T2, and an intensity depending on the characteristics of the corresponding phase. Unfortunately, the analysis of the recorded FID curve is a difficult task. Up till now, there is still uncertainty about the number of components and on their mathematical model functions in the decay curve. A lot of techniques have been proposed yet to analyze multicomponent decay curves based on, e.g., graphical extraction,2 linear3,4 and nonlinear5 (NLLS) least-squares methods, singular value decomposition,6-8 Gardner’s method,9-14 or * To whom correspondence should be sent E-mail:
[email protected].
Prony’s method.15-18 The success of these techniques is often hampered by the ill conditioning property of the exponentialsum fitting problem. Furthermore, the final results seem to depend strongly upon the input of the correct number of exponential components and upon the initial estimates of the preexponential factors N°i and decay constants λi. The original method of Gardner et al.9 has been developed to decompose a multicomponent decay curve, which is a sum of independent exponential functions. The technique is based on nonlinear transformations of the variables and on Fourier transformations. From the theoretical point of view, the method results in a spectrum of delta functions whose number is equal to the number of components and their positions are equal to the decay constants. From the practical point of view, this method did not gain much success because at that time the numerical evaluation of the Fourier transformations was very cumbersome and tedious. After development of the fast Fourier transform (FFT) technique19 the method was revised by Schlesinger10 and, independently, by Smith and Cohn-Sfetcu.11-14 Additionally, they combined the technique with an adequate low-pass filter in order to decrease the effects of computational noise and of the high frequency noise normally present in experimental data. They also showed that the improved method of Gardner et al. was able to analyze all kinds of multicomponent decay curves as long as the model functions were known. The analysis
10.1021/jp9829694 CCC: $18.00 © 1999 American Chemical Society Published on Web 02/17/1999
1394 J. Phys. Chem. B, Vol. 103, No. 9, 1999 resulted in a continuous spectrum in which the peaks are characterized by a finite width and intensity. The number of real peaks is equal to the number of components, their positions are good estimates for the decay constants, and their intensities are proportional to the preexponential factors, each divided by the corresponding decay constant. Prony’s method15-18 is a technique that has also been adjusted and modified over the years. It was published for the first time in 1795. The basic idea in this method is a separation of the linear parameters, N°i, and the nonlinear ones, λi. The original version has been evaluated extensively by Lanczos.20 He concluded that the method is numerically unstable, is not able to analyze a triexponential function unambiguously, and fails in the presence of noise. Two hundred years after Prony, Osborne and co-workers17,18 published two versions of an improved method of Prony: the recurrence and difference versions. Osborne et al. published an extensive mathematical study of these new versions and proved the asymptotic stability of Prony’s improved method. In a small simulation study they compared the technique with the Marquardt-Levenberg algorithm21 and concluded that the improved methods of Prony converge rapidly and are remarkably tolerant of poor starting values. In this contribution, we present the analysis of free induction decay curves of an isotactic polypropylene (iPP) film obtained at different temperatures. The improved method of Gardner et al. has been used in order to obtain information about the number of components (phases), to predict the best fitted model function (line shape), and to estimate the decay constants (T2 values). In addition, the improved method of Prony and the Marquardt-Levenberg algorithm were applied to quantify the contributions (fractions) of each component in the complete FID curve. The influence of the temperature on the fractions of the components is examined and related to the physicochemical structure of the polymer. The precision of the computed parameters is discussed for the analyses of three consecutively measured FID curves of iPP at 60 °C.
Schreurs et al. fit (Marquardt-Levenberg algorithm) or the MatLab 4.2c1 programs for the improved methods of Gardner et al. and Prony. Computational Methods The model function
NMR. Proton wide-line measurements of Solid State isotactic polypropylene (iPP) films were performed on a Varian Inova 400 spectrometer in a wide-line probe equipped with a 10 mm coil. Measurements as a function of temperature were done on the same sample, starting after a 30 min thermal stabilization time for each temperature. To determine the precision (stability) of the proposed method, three consecutive FID curves were recorded under identical spectrometer settings on a fresh sample. On-resonance FID’s were acquired by means of the “solid echo” technique (90°x-τ-90°y) developed by Powles and Strange22 in an effort to overcome the effects of the deadtime of the receiver, with an echo delay of 8.5 µs. According to these authors, the NMR signal, τ µs after the second pulse, is a good approximation for a standard free induction decay after a 90° pulse. The 90° pulses were of 2.3 µs duration. The dwell time during acquistion (sampling rate) is 0.5 µs, allowing a very accurate determination of the data point with maximum intensity. This maximum of the “echo” FID curve has been calibrated to time zero. All FID curves were accumulated with 350 scans and a preparation delay of 5 s. After recording, the number of experimental points (sampling rate) has been reduced to 25% (only one point out of four points has been used). These data-reduced FID curves were converted into an ASCII file and read into the Mac Kaleidagraphe 3.0.2 program for the NLLS
1
Mi exp ∑ a i)1
i
t
ai
T2i
(1)
can be used for fitting the NMR signal, where Mi, T2i, ai are the magnetization, spin-spin relaxation time, and shape parameter for the ith component, respectively. If all ai are equal to 1, the model function M(t) simplifies to a sum of exponential functions, which will be noted hereafter in the general form of a multicomponent exponential decay curve f(t) n
f(t) )
N°i exp(-λit) ∑ i)1
(2)
where N°i and λi are the preexponential factor and decay constant of component i, respectively. Comparison of eqs 1 and 2 shows that Mi ) N°i and T2i ) 1/λi for each component i, if ai ) 1. The multicomponent decay curve analysis consists of two parts. First, the number of components n is determined and the decay constants λi are estimated using the improved method of Gardner et al. (further denoted as the improved GGM method). This method can also validate the exponential (ai ) 1) or Gaussian (ai ) 2) character of each component i. Second, Prony’s improved method is used to quantify the decay constants λi and contributions N°i of each exponential component i. Both methods are explained in this section, starting from eq 2 as a general model function. 1. Improved Method of Gardner et al. The original method of Gardner et al.9 considered a multicomponent exponential decay curve f(t) (eq 2) as a special form of the Laplace integral equation
Experimental Section 1H
( ( )( ) )
n
M(t) )
f(t) )
∫0∞exp(-λt) g(λ) dλ,
t>0
(3)
provided n
g(λ) )
N°i δ(λ - λi) ∑ i)1
(4)
where δ(λ) is the Dirac delta function. Gardner et al. used a nonlinear transformation of the variables, x ) ln t and λ ) exp(-y), and Fourier transformations in order to obtain the solution of eq 3 as
g(e-x) ) F
{
-1
F [exf(ex)] x
F [exe(-e )]
}
(5)
where F is the Fourier transform operator. The resulting spectrum g(e-x) vs x is equivalent to the spectrum of g(λ)/λ vs λ. The number of peaks in the continuous spectrum equals the number of components n, the positions of the peaks correspond to the decay constants λi, and the intensities of the peaks are proportional to N°i/λi. Schlesinger10 showed that the cumbersome numerical evaluations of the Fourier integrals were no longer needed after the development of the fast Fourier algorithm.19 The GGM technique was further improved by using a low-pass filter10-14 in
Analysis of Solid State FID Curves
J. Phys. Chem. B, Vol. 103, No. 9, 1999 1395
order to reduce the effect of experimental and computational noise. The best results could be obtained with a Gaussian lowpass filter
( )
H(µ) ) exp -
µ2 µd2
(6)
ff(x) ) gg(x) X kk(x)
∫-∞+∞gg(z) kk(x - z) dz
(7)
with ff(x) ) exf(ex), gg(x) ) g(e-x) and the model function kk(x) ) ex exp(-ex). Equation 7 can be deconvolved for gg(x) by applying Fourier transform techniques:
gg(x) ) F
-1
{
}
F [ff(x)] F [kk(x)]
(12)
where δk are the recurrence Prony parameters and the characteristic polynomial Pδ(z)
(9)
k)1
n+1
∑ξkzk-1
(10)
k)1
Expression 9 can then be transformed as
(13)
has n roots Fj ) e-λj/N. The problem is now to find the optimal set of Prony parameters δk in the recurrence equations (12), followed by the calculation of the roots of the corresponding polynomial Pδ(z) (eq 13). The recurrence equations (12) can also be stated in matrix notation. Let fi ) f(ti), i ) 1, ..., N, f ) (f1, ..., fN)T, and let Xδ be the N × (N - n) matrix
[ ] δ1 · · ·
·
·
·
δ1
·
·
·
(11)
(14)
· · ·
δn+1
where δk are the recurrence parameters, then f satisfies
XTδ f ) 0
(15)
The method that Osborne et al.17,18 developed for finding the optimal set of Prony parameters is based on the treatment of the exponential fitting problem as a separable regression. This constitutes the second formalism. By definition, a least squares problem is called separable if the fitting function can be written as a linear combination of functions involving further parameters in a nonlinear manner. The model function f(t) in eq 2 can also be written as
f ) A(λ)N°
(16)
where f ) (f1, ..., fN)T, A is a (N × n) matrix function of λ with elements Aij ) e-λjti, and N° ) (N°1, ..., N°n)T. Osborne16 showed that such a separable problem can be transformed to a minimization problem with one equality condition involving the nonlinear parameters only. Let y ) (y1, ..., yN)T be the vector of observations that can be represented best by the model function f in eq 16. For any fixed value of λ the sum of squares
Φ(N°,λ) ) (y - A(λ)N°)T(y - A(λ)N°)
In eq 9 D is the differential operator and ξk are constant coefficients. The characteristic polynomial (10) has n roots (-λj). In practice, f(t) is often known in a finite number N of discrete points t1, ..., tN in a finite interval. Assuming that the interval is ]0,1], then ti ) i/N for i ) 1, ..., N. Osborne et al. represented eq 9 in two equivalent discrete formulations: the difference and recurrence versions. Only the recurrence version, which we have applied in the FID analyses, is explained further. The recurrence version uses the forward shift operator Π:
Π f(t) ≡ f(t + 1/N)
δkzk-1 ∑ k)1
Xδ ) δn+1
n+1
∑ξkDk-1f(t) ) 0
Pδ(z) )
(8)
Cohn-Sfetcu et al.12 proved that this method can also be used to analyze multicomponent Gaussian, Lorentzian, and (sin x)/x type functions by using another model function kk(x). They13 applied the improved GGM method for the analysis of pulsed NMR measurements. In the present work, a MatLab 4.2c1 program is used for the deconvolution of multicomponent exponential and Gaussian functions with the improved GGM method combined with a Gaussian low-pass filter. As input, the program requires a data set of f(t) at equal intervals of ln(t) and a µd value for the Gaussian low-pass filter. As output, it generates the g(λ)/λ vs λ spectrum. The number of peaks in this spectrum is equal to the number of components n in the decay curve. The shapes of the peaks are symmetric and the peak positions are good estimates for the decay constants λi when the correct model function is used. 2. Prony’s Improved Method. The improved method of Prony16-18 is based on two formalisms. The first states that the coefficients of the differential equation (eq 9) with characteristic polynomial Pξ(z) (eq 10) comprise a more stable parametrization than those of eq 2.
Pξ(z) )
∑δkΠk-1f(t) ) 0
k)1
n+1
characterized by its parameter µd. Smith et al.14 recognized that the solution offered by Gardner et al. implies a mapping of the integral equation (eq 3) into the convolution integral equation (eq 7)
)
n+1
(17)
is minimized by
N ˆ °(λ) ) (ATA)-1ATy
(18)
Substituting this back into Φ gives the reduced sum of squares
Ψ(λ) ) Φ(N ˆ °(λ),λ) ) yT(I - PA)y
with
PA ) A(ATA)-1AT
(19)
where PA is the orthogonal projection onto the column space of A. The least squares problem may be solved by minimizing Ψ with respect to λ, and recovering the least squares estimate
1396 J. Phys. Chem. B, Vol. 103, No. 9, 1999
Schreurs et al.
of N° from eq 18. Note that XTA ) 0, so that the columns ofX and A span orthogonal spaces. Therefore, eq 19 can be rewritten in terms of the Prony parameters
Ψ ) yTPXy
PX ) X(XTX)-1XT
with
(20)
where PX is the orthogonal projection onto the column space of Xδ. Hence, Ψ is a function of the recurrence Prony parameters δk. The derivative of Ψ with respect to δ can be written as
Ψ ˙ δ ) 2Bδ(δ)δ
(21)
with Bδ a symmetric (n + 1) × (n + 1) matrix function of δ with elements T yBδij ) yTXδi(XTδ Xδ)-1Xδj T Xδj(XTδ Xδ)-1XTδ y (22) yTXδ(XTδ Xδ)-1Xδi
and with Xδj ) ∂Xδ/∂δj. The function Ψ is now minimized with respect to the Prony parameters δ with an additional constraint δTδ ) 1 so that the elements of δ are normalized. The necessary condition for a stationary value leads to a nonlinear eigenvalue problem:
(B(δ) - νI)δ ) 0 δTδ ) 1
in
δ ) δ*
(23)
Let W(δ,ν) ) Ψ(δ) + ν(1 - δTδ) be the Lagrange auxiliary function with ν the Lagrange multiplier for the constraint and Wδ ) 2B(δ)δ - 2νδ ) 0. Because Ψ does not depend on ||δ||, Ψ 4 (δ) must be in a direction orthogonal to δ, so that
˙ (δ) ) 2δTB(δ)δ ) 0 δTΨ
(24)
Premultiplying W 4 δ ) 0 by δT and using eq 24 leads to the result ν ) 0. Thus, the nonlinear eigenvalue problem (23) can be solved iteratively as a linear eigenvalue problem. Given δk, solve the linear eigenvalue problem
{
(B(δk) - νk+1I)δk+1 ) 0 T
δk+1 δk+1 ) 1
(25)
for δk+1 and νk+1, with νk+1 the nearest to zero of such solutions. Convergence can be accepted when
νk+1 , ||B|| This iterative scheme is the improved Prony algorithm for the recurrence equations. In the present work, we used a MatLab 4.2c1 program for the recurrence version of Prony’s improved method. As input, the program requires a data set of the function f(t) at equal intervals of t, the number of exponential components n, and optionally, starting values for the decay constants λi. As output, the program generates the optimized set of (N°i, λi) parameters for each component i. Results and Discussion Before the data analysis, the FID curves of the iPP sample at 40, 60, 80, and 100 °C were reduced by a factor of 4 (see Experimental Section). Figure 1 shows the adapted FID curve
Figure 1. FID curve of an iPP sample at 60 °C on both linear and semilogarithm (natural) scales.
obtained at 60 °C on both linear and semilogarithmic (natural) scales. These multicomponent decay curves contain information about the physicochemical structure of the polymer sample: each component is characterized by a specific model function containing a spin-spin relaxation time T2 and an intensity as parameters depending on the properties of the corresponding state. In general, in T2 measurements the FID for polymer segments below their glass-transition temperature, Tg (the rigid glassy or crystalline state), can be approximated by a Gaussian function, while the decay for polymer segments well above their glass-transition temperature can be approximated by an exponential function. Only for very rigid systems or low temperatures can the fastest decaying part be fitted by an Abragamian function23 rather than by a Gaussian one. Hence, for systems having rigid and mobile components, a linear combination of these functions is used in order to obtain an NLLS fit of the free induction decay. The form of the fitted equation for an NMR signal M(t), given by eq 1, consists of several decay function types depending on the shape parameter ai of the ith component (ai ) 1, exponential; ai ) 2, Gaussian; 1 < ai < 2, Weibullian). Still a lot of discussion is going on concerning the determination of the number of components n used to fit the FID curve. For iPP, regions with high structural order (crystalline phase (C)) and low molecular mobility behave differently from other regions with low structural order (elastomeric phase (E)) and high molecular mobility. As a transition between these regions, an intermediate phase (I) is assumed, having a lower molecular mobility compared to that of the E phase but a lower ordering compared to that of the C phase. According to Tanaka24,25 the FID curve of iPP is a superposition of three distinct curves, corresponding to the C, E, and I regions. The existence of a third intermediate component was proposed by Fujimoto et al.26 for low-density polyethylene (LDPE), in agreement with the two spin-lattice relaxation times T1 found in the amorphous phase by Peterlin et al.27 Furthermore, it is generally accepted that in the 40-100 °C thermal range the component having the shortest T2 shows a Gaussian type decay and may be ascribed to the crystalline phase while the component having the largest T2 decays exponentially and can be attributed to the elastomeric phase (glass phase well above Tg). The complete FID curve
Analysis of Solid State FID Curves
J. Phys. Chem. B, Vol. 103, No. 9, 1999 1397
can then be represented by either eq 26 or 27:
f(t) ) N°Ce-0.5(t/T2C) + N°Ie-0.5(t/T2I) + N°Ee-t/T2E
(26)
f(t) ) N°Ce-0.5(t/T2C) + N°Ie-t/T2I + N°Ee-t/T2E
(27)
2
2
2
The analysis of these multicomponent decay curves consists of two separate parts. In the first part, the improved GGM method is used as a qualitative technique for the determination of the number of components, the best model function of the intermediate phase and for the estimation of the decay constants. In the second part, the recurrence version of Prony’s improved method is used to quantify the parameters T2i ) 1/λi and the fractions N°i/N° (N° ) N°C + N°I + N°E) for each exponential component i. 1. Qualitative Analysis: The Improved Method of Gardner et al. Combined with a Gaussian Filter. In the improved GGM method, the choice of the model function kk(x) depends on the functional relation between the variable t and the response f(t) of the decay curve. Equations 26 and 27 suggest that a mixture of Gaussian and exponential components can be expected. Therefore, first, a bicomponent decay curve (eq 28) is simulated and analyzed by using the improved GGM method with different model functions kk(x).
f(t) ) 1000e-(t/10) + 100e-t/100 2
(28)
Parts a-c of Figure 2 show the resulting g(λ)/λ vs λ spectra. Figure 2a is obtained by using the model function kk(x) ) exp(-exp(x)). Two peaks are present: a symmetric one at λ ) 0.01, corresponding to the decay constant of the exponential component, and an asymmetric one with a large negative lob at the right side belonging to the Gaussian component. At its right side a small positive maximum can be found. A closer look at the spectrum (Figure 2b) shows also very small error ripples at both sides of the exponential component. These additional signals are computational artifacts due to extrapolation and rounding off errors. They can be distinguished from real signals because the relative intensity between a real and an artifactual peak changes when other µd values for the Gaussian filter are used. Figure 2c is the output of the improved GGM method using the model function kk(x) ) exp(-exp(x2)). Again two peaks can be seen: a symmetric one at λ ) 0.1 belonging to the Gaussian component, and a broader and somewhat asymmetric peak belonging to the exponential function. Hence, if the model function kk(x) is chosen correctly, a symmetric peak can be found at the corresponding value of λ in the g(λ)/λ vs λ spectrum. The experimental FID curves are collected at constant time intervals. Since the input data set of the function f(t) in the improved GGM method combined with a Gaussian filter should be gathered at constant intervals of x ) ln t ) 0.1, a new data set is calculated using a cubic spline interpolation method on the original data set. The model function kk(x) in eq 7 is set to exp(-exp(x)) in order to evaluate the exponential components and to exp(-exp(x2)) in order to locate the Gaussian components. In all analyses, the implementation of a Gaussian lowpass filter was necessary for reducing most of the experimental and computational noises. The g(λ)/λ vs λ spectra obtained from a first analysis of the FID curves, to determine the exponential components are given in Figure 3a-d. Besides the expected error ripples (vide supra), which are now additionally amplified by noise present in the experimental data, three significant peaks can be observed by considering the four parts (a-d) of Figure 3. A first exponential
Figure 2. g(λ)/λ vs λ spectra obtained in the analysis of f(t) ) 1000e-(t/10)2 + 100e-t/100 by using the improved GGM method with a model function (a, b) kk(x) ) exp(-exp(x)) and (c) kk(x) ) exp(exp(x2)). Figure 2b is an expanded plot of (a) in the region λ ) 0.0010.1.
component is situated at λ ≈ 0.0095-0.015 µs-1 (or a T2 ≈ 67-105 µs) and can therefore be assigned to the elastomeric phase (E). The central peak at λ ≈ 0.03-0.05 µs-1 can for the same reason be attributed to the intermediate phase (I). The third peak is asymmetric and has a large negative lobe at the right side. Analogous to the case illustrated in Figure 2a, the third signal should be related to a Gaussian component. To confirm
1398 J. Phys. Chem. B, Vol. 103, No. 9, 1999
Schreurs et al.
Figure 3. Qualitative analysis of the FID curves obtained at (a) 40, (b) 60, (c) 80, and (d) 100 °C by using the improved GGM method with kk(x) ) exp(-exp(x)) in order to identify the exponential components.
that conclusion, a second analysis of the FID curves is carried out using the improved GGM method, in which the model function kk(x) ) exp(-exp(x2)) is used in conjunction with the same Gaussian filter. The resulting g(λ)/λ vs λ spectra are plotted in Figure 4a-d. A narrow symmetric peak at λ ≈ 0.08 µs-1 or T2 ) (λx2)-1 ≈ 8.8 µs is superimposed on a broad composite feature. This symmetric peak is related to the Gaussian component in the FID curve and can be assigned to the crystalline phase (C). The results obtained with the improved GGM method indicate that eq 27 is the best model function for the studied FID curves of the iPP film over the temperature range 40-100 °C. The multicomponent decay curves are superpositions of a single Gaussian component and two exponential components. The estimated T2 values for each component at the specific temperatures are summarized in Table 1. 2. Quantitative Analysis: Prony’s Improved Method. It should be noted that Prony’s improved method can only unravel multicomponent exponential decay curves. Therefore, the Gaussian function belonging to the crystalline phase (C), was first quantified using the Marquardt Levenberg algorithm with model function (27) and subtracted from the original signal. Then, the intermediate (I) and elastomeric (E) components in the remaining signal have been quantified using the recurrence version of Prony’s improved method. 2.1. Crystalline Phase. The improved GGM method shows that the crystalline phase is represented by the Gaussian
component in eq 27. Furthermore, starting values for the T2 values could be determined (Table 1). However, the MarquardtLevenberg algorithm requires also starting values for the parameters N°i of each component i. These are empirically set to 50, 30, and 20% of the initial intensity of the FID curve at a specific temperature for the crystalline, intermediate, and elastomeric phases, respectively. The optimized Gaussian functions obtained in this fitting problem are collected in Table 2. A T2 value of ca. 8 µs remaining almost constant over the 40-100 °C temperature range agrees nicely with that reported by Tanaka et al.25 The Gaussian components are subtracted from the initial responses. The complete removal of the Gaussian component was checked by applying the improved GGM method with kk(x) ) exp(exp(x2)). The g(λ)/λ vs λ spectrum at 60 °C is shown in Figure 5. Comparing Figures 4b and 5 shows that only a very small fraction of the Gaussian component remains. The curves, obtained after subtracting the Gaussian components in Table 2, are therefore good approximations for superpositions of two pure exponential functions. Hence, they meet the requirements for Prony’s improved method. 2.2. Intermediate and Elastomeric Phases. Prony’s improved method requires the number of exponential components, which is for all the analyzed FID curves equal to 2 as known from the results obtained with the improved GGM method. This latter method provides also the starting values of the decay constants listed in Table 1. The decay constants are optimized
Analysis of Solid State FID Curves
J. Phys. Chem. B, Vol. 103, No. 9, 1999 1399
Figure 4. Qualitative analysis of the FID curves obtained at (a) 40, (b) 60, (c) 80, and (d) 100 °C by using the improved GGM method with kk(x) ) exp(-exp(x2)) in order to identify the Gaussian components.
TABLE 1: Estimates of the T2 Values Obtained with the Improved GGM Method T2 value (µs) phase
40 °C
60 °C
80 °C
100 °C
C I E
8.8 20.0 66.7
8.8 31.3 100
8.8 33.3 100
8.3 33.3 111
TABLE 2: Optimized Gaussian Functions, Obtained with the Marquardt-Levenberg Algorithm, Describing the Crystalline Phase in a FID Curve Recorded at a Number of Selected Temperatures T (°C)
fC(t)
40 60 80 100
2.26 × 105 exp(-0.5(t/8.03)2) 2.29 × 105 exp(-0.5(t/7.93)2) 2.05 × 105 exp(-0.5(t/7.90)2) 1.92 × 105 exp(-0.5(t/7.79)2)
with Prony’s method and the preexponential factors are evaluated. The results obtained from these analyses are shown in Table 3. The T2 values and fractions of the intermediate (I) and elastomeric (E) phases are calculated. They are found to be in good agreement with the results obtained from the nonlinear regression analyses presented in Table 4. Moreover, it is worthwhile to mention that also for the NLLS fit, the errors and especially the χ2 are smaller when eq 27 is used instead of eq 26. With the combined application of the improved GGM technique and Prony’s method, the existence of an intermediate
Figure 5. g(λ)/λ vs λ spectrum at 60 °C after subtraction of the Gaussian component.
phase has been proved now and quantified. The results show that its participation is not negligible in a FID signal recorded in a temperature range of 40-100 °C. Furthermore, the combination of both methods in the analysis of FID curves performs equally well as the common used nonlinear regression
1400 J. Phys. Chem. B, Vol. 103, No. 9, 1999
Schreurs et al.
TABLE 3: Results Obtained Using Prony’s Improved Method T (°C) 40 60 80 100
phase
λ (µs-1)
105 N°
T2 (µs)
fraction (%)
I E I E I E I E
0.0452 0.0166 0.0372 0.0119 0.0379 0.0106 0.0431 0.00980
1.90 0.55 1.92 1.16 1.41 1.64 1.05 1.81
22.1 60.2 26.9 84.0 26.4 94.8 23.2 102.0
40.4 11.6 35.8 21.6 27.7 32.2 21.9 37.8
TABLE 4: Results Obtained with the MarquardtLevenberg Algorithm T (°C) 40 60 80 100
phase
T2 (µs)
fraction (%)
C I E C I E C I E C I E
8.04 22.2 60.6 7.93 26.9 84.1 7.90 27.5 96.3 7.79 23.3 102.1
48.0 40.5 11.5 42.6 35.8 21.6 40.2 28.6 31.2 40.2 22.0 37.8
comparable rigidity, and therefore, the fractions found are rather totals of the fractions of the C phase and parts of the I phase with very low mobility. Only at higher temperatures (from 80 °C for the C phase, constant fraction of 40% crystallinity) does the fraction represent the real C phase content. Concerning the molecular mobility (determining the line shape) below 80 °C, part of the E phase is clearly not far enough above Tg to give a pure Lorentzian line shape, resulting in an overestimation of the C and I phases. From 80 °C on, the Gaussian line shape is better defined, allowing a subtraction of the response of the crystalline species. At higher temperatures, a decrease of the I phase contribution and an increase of the fraction of the E phase can be seen toward their two limiting values, yielding the intrinsic contribution of the E and I phases. 2.4. Precision of the Parameters. When several FID curves are recorded under the same spectrometer conditions and are analyzed using the same methods, one can get an idea about the precision of the parameters found, which accounts for both instrumental and computational errors. Therefore, a series of three FID curves of iPP is recorded at 60 °C. The T2 values and the fractions of the corresponding phases are quantified using the procedure described above. The mean values and the 95% confidence intervals of the parameters obtained in these three analyses are (a) for the crystalline phase, T2 ) 7.12 ( 0.01 µs and a fraction of 40.5 ( 0.2%, (b) for the intermediate phase, T2 ) 22.8 ( 0.2 µs and a fraction of 39.55 ( 0.05%, and (c) for the elastomeric phase, T2 ) 80.3 ( 0.7 µs and a fraction of 19.9 ( 0.2%. Thus, for all the parameters a high precision could be reached. However, the values obtained in the thermal study described previously (see Tables 3 and 4) for an iPP sample at 60 °C are quite different, e.g., in Table 3, T2 ) 26.86 µs for the intermediate phase (I) at 60 °C, which is far outside the 95% confidence interval that has been calculated here. This could be due to experimental conditions such as, e.g., heterogeneities in the iPP film, folding of the film into the NMR tube, thermal stabilization, etc. Therefore, the confidence level that can be obtained in individual measurements is much lower (deviations of more than 1% are realistic). Thus, it can be concluded that the procedure used to analyze the FID curves is numerically stable since it is able to determine the parameters within the entire experimental precision of the recording of the multicomponent decay curves. Conclusions
Figure 6. Influence of temperature on the fraction of a phase.
methods. However, in our technique no a priori information is needed concerning the model functions describing the different phases, their number, characteristic decay constants and contributing fractions. 2.3. Influence of the Temperature on the Fractions of the Phases. Table 4 shows nearly constant values for T2C with increasing temperature, indicating that molecular motion is almost unchanged in the crystalline regions over the temperature range covered in this experiment. On the other hand, T2E increases with rising temperature, indicating that molecular mobility in the nonconstrained (decoupled from crystalline regions) amorphous regions is enhanced as the temperature is raised. Figure 6 shows the influence of temperature on the fractions of each phase. The contribution of the elastomeric (E) phase increases rapidly when the temperature increases. An opposite tendency can be observed for the intermediate (I) phase. The contribution of the crystalline (C) phase asymptotically tends to the real crystallinity. At lower temperatures the C and I phases are difficult to distinguish from each other due to their
The combination of two methods, programmed in MatLab 4.2c1, is proposed as a new technique for analyzing free induction decay curves. The improved GGM method is able to identify the three components, as well as their model functions, and to estimate the decay constants related to the T2 values. The recurrence version of Prony’s improved method is used to optimize the decay constants and to calculate the preexponential factors. Free induction decay curves of iPP films, recorded at temperatures between 40 and 100 °C, are analyzed using these methods. The results are in good agreement with those obtained with the Marquardt-Levenberg algorithm; however, less a priori information is needed in our proposed method. The presence of an intermediate phase (I) could be proved unambiguously and quantified. The influence of temperature on the fractions of the phases is discussed. The precision of the quantification of the parameters is shown to be much smaller than the one that can be reached with the recording of the FID curves. Acknowledgment. The authors wish to thank Profs. G. K. Smyth and M. R. Osborne for providing us the MatLab program
Analysis of Solid State FID Curves of the recurrence version of Prony’s modified method. This work constitutes part of research results of project IUAP-IV-P4/11 sponsored by the Belgian Federal Office for Scientific, Technical, and Cultural Affairs (OSTC). References and Notes (1) Fukumori, K.; Sato, N.; Kurauchi,T. Rubber Chem. Technol. 1990, 64, 522. (2) Mancini, P.; Pilo, A. Compt. Biomed. Res. 1970, 3, 1. (3) Hunt, G. E.; Grant, I. P. J. Atmospheric Sci. 1969, 26, 963. (4) Paatero, P. Annales Academiae Scientiarum Fennicae series A: VI Physica 1969, 319, 1. (5) Varah, J. M. SIAM J. Sci. Comput. 1985, 6 (1), 30. (6) Kung, S. Y.; Arun, K. S.; Bhaskar, D. V. J. Opt. Am. Soc. 1983, 73, 1799. (7) Barkhuijsen, H.; De Beer, R.; van Ormondt, D. J. Magn. Reson. 1987, 73, 553. (8) Lupu, M.; Todor, D. Chem. Int. Lab. Sys. 1995, 29, 11. (9) Gardner, D. G.; Gardner, J. C.; Meinke, W. W. J. Chem. Phys. 1959, 31 (4), 978. (10) Schlesinger, J. Nucl. Instrum. Methods 1973, 106, 503. (11) Smith, M. R.; Cohn-Sfetcu, S. Nucl. Instrum. Methods 1974, 114, 171.
J. Phys. Chem. B, Vol. 103, No. 9, 1999 1401 (12) Cohn-Sfetcu, S.; Smith, M. R.; Nichols, S. T. Proc. IEEE 1975, 63, 326. (13) Cohn-Sfetcu, S.; Smith, M. R.; Nichols, S. T.; Henry, D. L. Proc. IEEE 1975, 63, 1460. (14) Smith, M. R.; Cohn-Sfetcu, S.; Buckmaster, H. A. Technometrics 1976, 18 (4), 467. (15) Prony, R. J. Ecole Polytech. 1795, 1, 24. (16) Osborne, M. R. SIAM J. Numer. Anal. 1975, 12 (4), 571. (17) Osborne, M. R.; Smyth, G. K. SIAM J. Sci. Stat. Comput. 1991, 12 (2), 362. (18) Osborne, M. R.; Smyth, G. K. SIAM J. Sci. Comput. 1995, 16 (1), 119. (19) Cooley, J. W.; Tukey, J. W. Math. Comput. 1965, 19, 297. (20) Lanczos, C. Applied analysis; Prentice-Hall: Englewood Cliffs, NJ, 1956; p 272. (21) Levenberg, K. Quart. Appl. Math. 1944, 2, 164. (22) Powles, J. G.; Strange J. H. Proc. Phys. Soc. London 1963, 82, 6. (23) Fedotov, V. D.; Schneider, H. Structure and Dynamics of Bulk Polymers by NMR-Methods; Springer-Verlag: Berlin, Heidelberg, 1989; p 58. (24) Tanaka, H. J. Appl. Polym. Sci. 1982, 27, 2197. (25) Tanaka, H.; Inoue, Y. Polym. Int. 1993, 31, 9. (26) Fujimoto, K.; Nishi, T. Polym. J. 1972, 3 (4), 448. (27) Crist, B.; Peterlin, A. J. Polym. Sci., Part A2 1969, 7, 1165.