6175
J. Phys. Chem. 1993,97, 6175-6183
Population Inversion in a Multilevel System: A Model Study Bjarne Amstrup,td And& Liirincz,s and Stuart A. Rice'*? Department of Chemistry and The James Franck Institute, The University of Chicago, Chicago, Illinois 60637, Chemistry Department B, The Technical University of Denmark, DK-2800 Lyngby. Denmark, and Institute of Isotopes, The Hungarian Academy of Sciences, H- 1525 Budapest, Hungary Received: September 18, 1992; In Final Form: February 8, 1993
For most isolated molecules, the generation of a r-pulse which will create complete inversion of the population of two electronic states is impossible because the rotational contribution to the transition dipole moment is an explicit function of the projection of the total angular momentum along the electric field axis. It is, however, possible to obtain nearly complete inversion by use of phase- and amplitude-modulated pulses; we show examples of such inversion with scaled and swept hyperbolic secant pulses. We report the results of pulse shape calculations for model systems with one and eight transition dipole moments. The pulse shaping is based on optimal control theory, where both the time-dependent SchriMinger equation and the pulse energy are used as constraints and where the deviation of the pulse shape from a swept hyperbolic secant form is treated as a cost to be minimized. This use of a time-dependent penalty function is intended to bias the calculated optimal pulse shape in the direction of experimentally producible pulse shapes. The method is applied to a model system with two displaced potential energy surfaces which support harmonic vibrational motion.
I. Introduction It is-well known that one can generate complete population inversion in a two-level system (TLS) by applicationof a resonant electric field that satisfies the area condition'
where &o(t) is the envelope function of the radiation field and p is the matrix element of the transition dipole moment operator connecting the two levels. There are, in fact, an infinite number of so called *-phases that satisfy (1. l), since only the area of the pulse between the times to and tr is constrained. We note that a particular pulse will generate complete population inversion only for one value of the dipole matrix element; for any other value of the dipole matrix element, that particular pulse will generate less than 100% population inversion (see later). It is sometimes possible to find two "quasi-isolated" levels in a multilaevel system. Then, if the spectral width of the exciting pulse is less than the spacing between these levels and all others, it is possible to treat the dynamicsof population transfer between thequasi-isolated levels as if they formed a true two-level system. In general, however, the spectrum of the system of interest is such that the spectral width of an appropriate pulse is greater than the spacing between the pair of levels of interest and neighboring levels, whereupon application of the pulse to the system must result in excitation of many levels. This is, for example, the case in the Tannor-Rice*-3 scheme for guiding the evolution of a quantum system, since the pump and dump pulses could be, and normally should be, on the order of or shorter than a typical vibrationalperiod. Another instance where the isolation of a two-level subsystem from a many-level system breaks down occurs when the electric field has sufficient power to generate multiphoton processes. An example of this case, namely, the dynamics of population transfer to the upper state of a system with two displaced harmonic potential energy surfaces, has been studied by Soml6i and LiSrin~z.~ 7 The University of Chicago. t The Technical University of
Denmark.
8 The Hungarian Academy of Sciences.
We have already noted that a particular pulse shape will completely invert the population in a two-level system only for one value of the dipole matrix element. However, if the system of interest is a rotating molecule, there are many transitions between twovibronic states, each one correspondingto a particular allowed rovibronic transition. For example, a diatomic molecule composed of nuclei with zero spin has rovibrational states built on each electronic state which are characterized by the standard quantum numbers (v, J, mJ);if the nuclei have nonzero spin, F = J + I and mF must be added to the set. In the absence of an external electric or magnetic field, the (2J 1) levels labeled by mJ are degenerate but are connected to another state by allowed transitions withdifferentvaluesofii&.5 For R-branch transitions with linearly polarized light $ 8 0 0: ( ( J + 1)2 - m3)'I2,whereas for P-branch transitions we have G.20 a (P- m 3 ) l / 2 . Clearly, in such a case, because of the strong variation of the dipole matrix element with rotational substate, it is not obvious that one can find a pulse that will generate complete population inversion. Our approach to the design of pulse shapes that generate nearly complete population inversion in a system with many different dipole matrix elements connecting the states of interest is based on the use of optimal control theory. This approach was first introduced by Rabitz and co-workers6p7and exploited for the descriptionof processes involving two electronic potential energy surfaces by Kosloff et a1.8 Jin et aL9 have used optimal control theory to study electronic population inversion in a model system with two displaced harmonic potential energy surfaces, but they did not examine the effects associated with the existence of multiple transition dipole matrix elements, and they used somewhat different constraints on the electric field than weemploy in this paper. In this paper, we report the results of a study of population inversion in a model system with two displaced harmonic potential energy surfaces, with allowance for variation of the dipole matrix elements. Moreover, we introduce the deviation of the pulse shape from a prescribed form as a cost in the optimizationprocedure. The use of a time-dependent penalty function is intended to bias the optimization in the direction of experimentallyproduciblepulse shapes. We compare our results to those obtained using hyperbolic secant pulses, which are known to generate population inversion in an isolated two-level system in the high electric field limit.10 The notion that the hyperbolic
0022-365419312097-6175%04.00/0 0 1993 American Chemical Society
+
6176 The Journal of Physical Chemistry, Vol. 97, No. 23, 1993
secant is a good choice for the initial pulse shape comes from NMR spectroscopy. Warren et al.lO-l4 and Silver et ~ 1 . ~have 5 studied population inversion in NMR spectroscopy and laser spectroscopy for a variety of cases. They have shown that a hyperbolic secant pulse with a frequency sweep can induce a highly selective population inversion which, above a critical threshold, is independent of pulse amplitude for the two-level system; above that threshold, such a swept pulse will generate population inversion for a range of values of the dipole matrix element. We show that this is also the case for a system with two displaced harmonic potential energy surfaces, although the frequency sweep must be more rapid than for the case of the simple two-level system. We also show that it is possible to find optimized pulseshapesthat generate better than 99.8%population inversion in a model system with as many as eight degenerate levels with transition dipole matrix elements varying as if J = 8 or J = 7 in the P- or R-branches of the rotational spectrum.
Amstrup et al. Lagrange multiplier and a penalty factor, namely @ and A. In analogy with !P, @ is a two-component function defined by
and the objective functional with constraints takes the form
7' = ( *Afiqf)+
isff[ a, to
(@ih
H*)
- c.c.]
dt
+
Xt[J'flE(t)12 10 dt - El2 (2.8) It is sometimes desirable, say for the purpose of designing a field as close to an easily generated form as possible, to introduce constraints on the shape of the electric field. That can be done by introducing a weight function W(t)in the last term of (2.8). If the energy constraint is incorporated in the weight function, we find
II. General Considerations We now briefly outline the optimal control theory approach to the optical pulse shaping we have used. We consider a model system with two potential energysurfaces. Theequationofmotion in our system is (a, = a/&)
ih 13,*(3,t) = Hq(3,t)
(2.1)
Assuming that the Born-Oppenheimer approximation is valid and that the matter-field interaction can be treated semiclassically
X[J"1E(t)l2W(t) to dt - 112 (2.9)
The simplest weight function is a constant. In that special case, W(t)is the reciprocal of the target energy and (2.9) reduces to (2.8). However, in the general case W(t)also incorporates a constraint on the shape of the electric field. Application of the calculationsof variation in the usual manner (see, for example, refs 8,9,16,17) to (2.9) leads to the following set of coupled equations:
ih a,*(x,t) = H*(x,t)
(2.10)
i h d,@(x,t) = H@(x,t)
(2.1 1)
and (2.12) The subscript g refers to the "ground" state potential energy surface and u, correspondingly, to the "upper" state potential energy surface. Of course, Hi = Ti + K(3) where Ti = P / 2 m is the kinetic energy operator and V,is the corresponding potential energyoperator. pis the negative of the transitiondipolemoment operator,and & ( t ) represents the amplitudeof theclassical electric field. For convenience, in our calculations we allow the field to be complex and we use the rotating wave approximation. This procedure simplifies the numerical integration by removing the carrier frequency; the real part of the optimal electric field is readily extracted from the final results. As we are interested in the upper-state population at some final time tf, we define the objective functional to be optimized as a projection of the final state of the system (2.4) I = (#,(tf)l#,(tf)) = (*dPl*f) where P is a projection operator for the selected product and \kf is the wave function at the time tr. The projection operator in our case is
l@(tf))= Pl*(tf))
Equations 2.10-2.13 are further linked by the equation for the optimal field & ( t ) 4X[JfflE(f)12W(f)dt - 11W(t) & ( t ) - 2 q t ) = 0 4
In addition, it is always necessary that the system satisfy the equation of motion in (2.1). These two constraints introduce a
(2.14)
where
o(t)= -hi~ ~ 4 u ~ ~ ~ l P l-r (rc,,(t)IPl4&t))) c,*~o~
(2.15)
If the system has N,, different nonzero transition dipole moments, we must modify the objective functional from the form just dtscribed. We then have N,, different wave functions {eq, and each of these must satisfy the SchrMinger equation. As before, this constraint introduces a set of Lagrange multiplers {q. The objective functional equivalent to (2.9) is then
c.c.] dt We choose as one of the constraints on the system that the "energy" of the electric field is to be restricted in the optimization process:
(2.13)
+ A [ ~ t f l E ( t ) 1W(t) 2 dt - 112 (2.16) 10
This functional leads to a set of equations similar to the ones in (2.10)-( 2.15): i h a,'ki(~,t)= H,@(x,t); j = l...Nfi
(2.17)
i h at@ = HI@; j = l...N,,
(2.18)
I@(t0)) =
[yq;
j = I.."
(2.19)
The Journal of Physical Chemistry, Vol. 97, No.23, 1993 6177
Population Inversion in a Multilevel System
TABLE I: Potential Energy Surface Parameters for the Three Systems Consideredc I
I1 7Xleau
5XlVau 0 . 4 ~ 2.53
01
boa
bo/@
0 . 4 ~ 2.99
111
0.8 I
'
O
I
5x10-4au 0.6e 3.79
4*9h2/m&is the bohr radius. b u = ( h / M ~ ~ ) )C~The / 2 mass .
a
M = 80 OOOm, and the upper-state frequency o. = 5 X 10-4 au = 2.07 X 1013 s-l (-109.7 cm-l) are the same for all three systems.
0.3 j
I
I.
I
0.2
1
0.0
0.0
0.5
1.5
1.0
2.0
CLl(ea0) Figure 3. Population inversion of the TLS as a function of transition dipole moment for the 'sech-r-pulse" of Figure 1.
0.1
0.0
1.o j
0
8000
16000 time [au]
Figure 1. The "sech-*-pulse" for the TLS (24 OOO au
j
p
i L
24000
-
580.5 fs).
o
100 9(t)- - -
0.0 0.0
0.5
1.o
1.5
2.0
CLLl(ea0) Figure 4. Population inversion of the TLS as a function of transition dipole moment for the 'sech-lor-pulse" of Figure 2.
-0.2
-0.3
t
1 0
8000
16000 time [au]
24000
Figure 2. Real (R(t))and imaginary (3(t))parts of the frequencyswept 'sech-lh-pulse" (4.4) for the T U .
Id(tf)) = pl@(tf));j
(2.20) = I...IV~ which are further linked by the equation for the optimal field C 4X[S"(&(t)12 W(t)d t - 11W(t) 10
e(?)- 2 0 & )
=0
(2.21)
where
results as the Lanczos propagation scheme until 99% inversion was obtained. The second-order differencing method is faster than the Lanczos scheme, though it is generally not as precise as the latter. In fact, it could not be used for the final iterations in our study. In our numerical procedure, we find the electric field at a discrete set of times by the conjugate gradient method. Let us denote the field at time tk by &(tk) = 6, = a
k
+ iJk
gk
= 4X[Sff(a2(t) + J 2 ( t ) )W(t)dt - 11 x to
Wk(Rk
III. NumericalMethods We have used an iterative algorithm for the self-consistent determination of the optimal shape C(t), as given by (2.14) or (2.21), subject to satisfying the equations of motion for 9 and cli. The numerical solution of the time-dependent SchriSdinger equation is carried out using a second-order differencing method (see, for example, ref 18) or a short iterative Lanczospropagation scheme19J3together with the fast Fourier transform method for evaluating the kinetic energy. The second-order differencing method was used in the preliminary stages of this work. A test of the phase and amplitude of the wave function on both surfaces showed that the second-order differencingmethod gave the same
(3.1)
The gradient for the kth time step is then
+ iJk) - 20,
(3.2)
where 0, and w k are the overlap function ((2.15) or (2.22)) and the weight function, respectively, at time tk. Let i stand for the number of iterations. The direction conjugate to the gradient can then be defined for the kth time step in the ith iteration as
dk= -dk+ @-'d;1 with (Polak-Ribihe algorithm)20-22
(3.3)
6178 The Journal of Physical Chemistry, Vol. 97, No. 23, 1993
Amstrup et al.
TABLE Ik Amplitudes of Electric Fields Used' I
I
I n !
I
I
I
I
1
I
-.__----
I
111
I1
aQ ? = 2.5 X 10-4hartrees/(eao) = 1.286 X lo8 V m-'. They are found from the weak field limit using the dipole moment = 1 e@.
1
.
d
45 35
0.8 -
1.0
I
I
J
1
(4 -
2 0.611
-B
t
li
0.4
o.2ij
t
0.2
1
0.0 0.0
I
I
0.5
1.o
I
(111)
1.5
1 2.0
CLl(ea0) Figure 7. Population inversion of the displaced system given electric fields of the form in (4.8). The amplitudes used for the electric fields
are those from the first line of Table 11. The numbers in the figure are the values of the sweep parameter a used. 1-111 denotes the system.
For the initial step
Di = -gi
(3.5)
The new electric field is then defined by
+
6;' = 6 ; aiDi
(3.6) where a' is a coefficient chosen to maximize T(a? (the sum of the first and last term on the right-hand side of (2.9)) in the direction of 0:. As an initial guess for a' (see, for example, ref 20), we use the functional second derivatives of 1 with respect to the electric field:
where A is
A = X[~'r(s2z(t) + J2(t))W(t)dt - 11 '0
(3.8)
We have used a line search to find a'. This procedure involves finding I(aak) for a number, typically 4-10, of values of a. The value of a that gives the maximum value of ?(a&, together with the two closest a values and the target function values at these points, defines a second-order polynomial in a: T(aab) 'Y f,(a) = A d
+ Ba + C = A(a + ' / p / A ) 2
+ c- '/4Bz/A (3.9)
If A in (3.9) is negative, aa; = - (1/2)a$/A is used as the ai in (3.6) and C- (1/4)B2/A is the new estimated target function value. Otherwise, we use the value of aa; that gives the maximum value of 1as the a' in (3.6). In almost all cases, we
Population Inversion in a Multilevel System
The Journal of Physical Chemistry, Vol. 97, No. 23, 1993 6179
0.0005
0.0005
O.ooo4
O.OOO4
0.0003
0.0003
0.0002
0.0002
0.0001
0.0001
-0.0003
A
O.oo04
0.0004
0.0003
0.0003
0.0002
0.0002
0.0001
0.0001
‘ 1 . -0.OOO4
O.oo04
O.OOO4
0.0003
0.0003
0.0002 0.0001
a 0.0002
0.0001
-0.0005
O.OOO4
0.0004
0.0003
0.0003
0.0002
0.0002 I
0.0001 0
\
I
\
8000
time [au]
16000
0.0001
2z
0
0
Figure 8. Optimized electric fields (absolute valuw) and the initial field (top right) for the systems studied. The figure at the top left is for an energy-constraincdelectric field, whereas the other figures are for the shape-constrained electric field. The value of X is indicated for each shapeconstrained optimized field.
found that A was negative. This kind of line search is needed because, in the population inversion problem, both thedenominator and the numerator of (3.7) decrease toward zero as the inversion improves, and that behavior leads to numerical uncertainties. As is evident, a;is derived from a quadratic approximation to I with respect to the electric field. This approximation, however, can be very bad; in some cases, the range of ai used in the line search is less than 10-6 ah, whereas for other cases it has the same order of magnitude as a;. The line search speeds up the convergence compared to just taking the af directly from (3.7) and also compared to taking the value that gives the largest value of I (aa;) in the calculated points.
IV. Results of the Calculations The model used in the calculationsinvolves two electronicstates. The potential function for each electronic state is taken to be a
one-dimensional harmonic surface:
We assume that at t = tothe wave function is the eigenfunction of the ground vibrational state of the ground-state potential energy surface. When h o = 0 and w8 = o,,,the system reduces to the classical TLS because there exists only one allowed transition, from the initial state n, = 0 to nu = 0. Table I displays the parameters of the harmonic surfaces used in this study. The 0 4 energy is set to zero, as we have used the rotating wave approximation. The parameter set labeled I in Table I is also used in the TLS calculations, except that then A x 0 = 0. We consider, for comparison with the calculations reported later, a two-level system, which we drive with a T pulse of the
6180 The Journal of Physical Chemistry, Vol. 97, No. 23, 1993
Amstrup et al.
TABLE III: Results of Sbape-Restricted and Shape-Unrestricted (Bottom Row) Optimization Calculations Using the Dipole Moment C~ = 1 e& and System Ib ($u(rdhh(rd) x i J&14t)lzdr (%iW -0.0050 -0.0025 -0.0010 -0.0005 -0.0004 -0.0003 -0.3854
0.9887 0.9940 0.9941 0.9982 0.9987 0.9995 1.m
-
0.9856 0.9907 0.9893 0.9948 0.9959 0.9976 1.oooO
1.44X 1.51 X 1.22 X 1.28 X 1.29 X 1.48 X 0.88 X
1.000 1.Ooo 1.Ooo 1.Ooo 1.Ooo 1.000 1.000
lo-." lk3 lW3
lk3 10-3
Units of Au. 1 X 10-3 au e,,cJ&le(t)12 dr = 1.7 mJ/cm2. ('
[email protected]) is given to show conservation of the norm of the wave function.
TABLE I V Results of Shape-Restdcted and Shape-Unrestricted (Bottom Row) 0 tihation Calculations Using the Dipole Moments Given in &.IO) and System IC
x
CY
-0.0025 -0.0005 -0.0025 -0.0025 -0.3854
10 10 3 0 0
Ej=I(4(rf)l4(rf)) 7.9967 7.9975 7.9932b 7.9980 7.9981
-
I
S&14t)l2dt
7.9961 7.9972 7.9930 7.9972 7.9981
1.64 X O 1.72 X 1 t 2 1.27 X 1 k z 0.69 X 1p2 0.88 X
0 Units of au. 1 X 10-3 au cgJ&+(r)f dr = 1.7 mJ/cm2. Not fully converged. CY denota the swcep parameter used for the initial guess.
-0.004
-0.002
0 0.002 frequency [au]
0 004
0.006
Figure 9. Absolute values of the Fourier transforms of the initial field
(bottom), optimized field with energy constraint only (second from bottom), and optimized fields with shape-restricted constraint on the electric field (value of X indicated). The central frequency of the initial field is indicated by the dashed line (0.0015 au = 329.2 cm-l). following form:
Cl!/(t) = 6 ;: where 6 ;:
sech((t - t O ) / t p )
is determined from the condition
(4.3) with pc being the relevant matrix element of the transition dipole moment. Note the factor I/z in (4.3) which arises from the use of a complex electric field in the rotating wave approximation. The resulting electric field, using = 1 eao, is shown in Figure 1. The value of 1, is 2000 au, corresponding to a fwhm (intensity value) of 3525 au (85 fs). The population inversion as a function of the transition dipole moment for this electric field is shown in Figure 3. The real and imaginary parts of the scaled and frequency-swept sech-pulse lOC;?(sech((r - fO)/tp))l+Si (4.4) are shown in Figure 2, and the correspondingpopulation inversion as a function of transition dipole moment is shown in Figure 4. These results, which are all well-known (see, for example, ref 1 l), are given for comparison with the results described below. Note that we get the same curves with, e.g., a prefactor of 5 instead of 10, except for a change in the abscissa interval from [0,21 to t0,41. For the displaced oscillator systems labeled I, 11, and I11 in Table I, we consider an electric field with the envelope function: Cn,, hpX(t)
= C y sech((t - f O ) / t p X)
exp(i(nw, + l/Z(% - w&f) (4.5) This field is resonant with the nth vibrational state on the upper surface for system X (I, 11, or 111). The transition probability in theweak fieldlimit isa functionoftheoverlapsof thevibrational wave functions of the ground and excited electronic states as well as the center frequency of the electric field. To be more specific,
for very small values of the constant d, perturbation theory gives (+u,I+u(tf> ) = sin2(ap&)
(~Nc)~
(4.6)
Then, a is determined by calculating (+U(tf)l+u(t~))for some sufficiently small value of 6, and is found from
62'
ppx
=*/2w (4.7) These pulses will also be called sech-*-pulses. We have displayed the relevant values determined in this way, again for p, = 1 eao, in Table 11. For system I, the n = 3 vibrational state on the upper surface has the largest Franck-Condon overlap with the initial vibrational state on the lower surface. Note that the value for displ 6nS g3,, / O,r , amely 1.87,is smaller than the value taken from the reciprocal of the largest Franck-Condon factor, namely 2.12. The fact that this ratio is smaller than the reciprocal of thelargest Franck-Condon factor is a consequence of the use of short pulse excitations for the systems studied. These pulses have spectral bandwidths that are larger than the vibronic level spacing and hence will excite a coherent superpositionof severaleigenfunctions. As the excitation pulse approaches a delta function, the cited ratio approaches the two-level system value. For a square pulse, the ratio could not reasonably be smaller than 2.12.9 The population inversions in this displaced system as a function of transition dipole moment for the sech-a-pulse and the scaled frequency-swept sech-r-pulse are shown in Figure 5 and Figure 6, respectively. Figure 4 showsthe very high efficiencyof population inversion in the TLS for the frequency-swept hyperbolic secant pulse. It is evident from Figure 6 that such a pulse does not have the same efficiency of population inversion in the displaced surfaces case. Results similar to the ones displayed in Figure 6 are also found for fields resonant with the n = 1, n = 2, n = 4, hand n = 5 vibrational levels. It is, however, possible to achieve population inversion as efficient as that for the TLS by increasing the sweep parameter considerably from the value 5 used in (4.4). We have calculated the population inversion as a function of the transition dipole moment for a number of different sweep parameters a, where a is defined in the following equation for the electric field: n,,
Population Inversion in a Multilevel System
The Journal of Physical Chemistry, Vol. 97, No. 23, 1993 6181 0.0020
0.0020 A
0.0016
0.0016
0:0012
0.0012
0.0008
0.0008
0.0004
O.OOO4
0.0016
0.0016
0.0012
0.0012
0.0008
0.0008
0.0004
0.0004
0.0016 0.0012 0.0008 0.0004 0
A
C
Time [au]
0.0016 0.0012 0.0008 0.0004 0
i Time [au]
Figure 10. Absolute values of the shape-restricted optimized electric fields with eight transition dipole moments (B, C, E, and F). D shows that the target electric field for these calculations. A is the shape-unrestricted optimized electric field. Corresponding values of h and a (A, a)are A (-0,3854, 0), B (-0.0025,lo), C (-0.0005, lo), E (-0.0025,3),and F (-0.0025,O). &(t)
= 10&::pX(sech((t
a weight function which is the squared inverse of the field that was used as an initial guess divided by the time interval:
- to)/tp))l+ia X
exp(i(nw, + l/&”
- w,))t) (4.8)
The result is shown in Figure 7. Unfortunately, the frequency spectra for these strongly swept electric fields are unreasonably broad (see Figure 11). Can we find other electric fields that give the same strong population inversion of the TLS for a range of dipole moments? To answer this question we first test how much influence the shape-restricted optimization has on the resultant electric fields. We have found two electric fields that give essentially 100% inversion for one value of the dipole moment (pc = 1 eao) for system I. In both cases, we have used for the initial guess the hyperbolic secant field resonant with then = 3 state and amplitude as in (4.5). The absolute value of the initial field is shown in the top right side of Figure 8, and the corresponding Fourier transformis shown in thebottomof Figure9. The weight function for the first optimized field calculation is a constant; that is, we have a penalty on the energy of the electric field as in (2.8) and none on the shape of the field except that it is restricted to a time interval of 24 000 au (580.5 fs). The chosen A’ (as in (2.8)) was -500 000 au, and the “energy” E = 8.7795 X 10-4 au. Thus, the corresponding value of X (as in (2.9)) is X’EZ = 0.3854. The population inversion for theoptimized field was found to be 100% (within 1 X lo-’), and the final “energy” was 8.7795 X 10-4 au. The resulting field is shown in the top left side of Figure 8, and the corresponding Fourier transform is shown in Figure 9 (second from bottom). For the second optimization calculation, we used
W(t) = [e? sech((t - t o ) / t p ) ] - 2 ( t r -to)-’ (4.9) where the final X (as in (2.9)) is -0.0003. The resulting field is shown in Figure 8, and the corresponding Fourier transform is shown in Figure 9. The price one has to pay for restricting the shape of the pulse is an increase in the energy of the pulse. Thus, the energy of the shape-constrained optimized pulse is 69% larger than the energy of the energy-constrained optimized pulse. We have also carried out shape-restricted pulse optimization for other values of the parameter A. The results of these calculationsare summarizedin Table 111. (\kd\kf) is also displayed inTableIII toshow conservationofthenormofthewavefunction. The Fourier transforms of some of the final electric fields are shown in Figure 9, and the fields are shown in Figure 8. We note that the overall energy of the pulse does not change considerably as X decreases, but the pulse broadens considerably. That observation illustrates the importance of shape-restricted optimization. The somewhatlower energy of the energy-constrainedonly optimized field is due to the much higher value of A. Also, the frequency spectra of the optimized fields are clearly broader than the corresponding spectrum of the initial field. These fields do not have tails like the energy-constrained-onlyoptimized field, and they have somewhat higher energy. We have examined a simple model which illustrates the effect of varying a transition dipole matrix element on population inversion in a system with displayed potential energy surfaces. In our example, we consider that the electromagnetic field is
Amstrup et al.
6182 The Journal of Physical Chemistry, Vol. 97,No. 23, 1993
i
0.6
0.0I'
0.0
t
/
1
,
I
0.5
1.0
,
1.5
I
2.0
d(ea0)
Figure 12. Population inversion in the displaced system I as a functon of transition dipolemoment for theoptimized field with X = -0.0005. The
optimized points are denoted
Figure 11. Fourier transforms (absolute value) of the optimized fields from Figure 10 and of initial fields (dashed line).
linearly polarized, and we associate eight transitions with the excitation of the model system from the ground to the upper electronicenergy surface. These transitionsare assignedtransition dipole matrix elements as if they were associated with an initial level with J = 8 (P-branch) or an initial level with J = 7 (Rbranch). On the other hand, the model does not have any other sublevels, so it does not correspond to a real rotating diatomic. We have chosen to carry out calculations for the case that the maximum dipole matrix element has the value 2 eao, whereupon the eight transitions have dipole matrix elements equal to p,
= 2/8(82-
0'-
1)2
j = 1, ..., 8
(4.10)
The target function was chosen to be W(t) =
sech((t - to)/fp)]-2(tf-
to)-'
(4.1 1)
and we have set A equal to either 4.0025 or 4.0005. In all, four optimization calculations were carried out corresponding to an initial field with one shape but with different values for the sweep parameter, namely, a = 10, 10, 3, and 0, respectively. The calculations carried out with the same value of a had different values of A. We have also carried out a shape-unrestricted optimization of the field required to invert the population in the model system, using the same value of A as for the corresponding calculationwith only one transitiondipolematrixelement. Values of A and a for the several optimizations are displayed in Table IV. The calcualted optimal pulse shapes and their Fourier transforms are shown in Figures 10 and 11, respectively.
X.
The results displayed in Figures 10and 11have some interesting features. We note that the pulse shapes in parts B and C of Figure 10 are similar, even though the values of A are different. The pulse shapes in parts B, C, E, and F of Figure 10 show vanishing field for t = 0 and t = 24 000 au, whereas the pulse shape in Figure 10A (shape-unrestricted optimization) does not, clearly showing how the restriction of the pulse shape affects the wings of the calculated optimum pulse shape. The population inversion as a function of transition dipole moment for the model system with multipledipolematrix elements is shown in Figure 12. We note that nearly complete inversion is achieved for a considerable range of values of the dipole matrix elements.
V. Concluding Remarks The modification, hence also control, of the evolution of a quantum system has been demonstrated both formally and for model systems. However, the formal demonstration does not provide information about the extent of control possible, and the simple model systems considered greatly simplify the coupling between the system and the external field and also neglect important aspects of the molecular dynamics. This paper has addressed one of those omissions. Specifically, since all real molecules will have nonzero rotational excitation and the dipole transitions between rovibronic states have dipolematrix elements of different magnitude, it is not obvious that one can find a pulse shape that will invert all or nearly all of the population between two potential energy surfaces. We have shown how the use of optimal control theory permits construction of a pulse shape that does provide for population inversion in a model system. However, since we have not allowed the model system to have the full manifold of rotational levels, we have not allowed for depopulation via coupling between other pairs of levels as the frequency of the field is swept. Although our results provide an example of overcoming one complication associated with the properties of a real molecule, and we are confident that other complicationscan alsobeovercome,it remains tobeshown that,e.g.,onecangenerate nearly complete inversion in a real molecule, with all mechanical and electrical anharmonicity terms in the Hamiltonian accounted for, using pulse shapes that can be made in the laboratory. Acknowledgment. Bjarne Amstrup wishes to thank the Danish Natural Science Research Counsil for a grant (1 1-8929). This
Population Inversion in a Multilevel System research is part of a program supported by the National Science Foundation under CHE 8916046.
References and Notes (1) Allen, L.; Eberly, J. H. Opfical Resonance and Two-leuel Atoms; Dover, 1987 (J. Wiley & Sons: New York, 1975). (2) Tannor, D. J.; Rice, S.A. J . Chem. Phys. 1985,83, 5013-5018. (3) Rice, S. A.; Tannor, D. J.; Kosloff, R. J. Chem.Soc.,Faraday Trans. 2 1986,82, 2423-2444. (4) Sombi, J.; Ldrincz, A. Phys. Rev. 1991, A43, 2397-2401. ( 5 ) Townes, C. H.; Schawlow, A. L. Microwave Spectroscopy; Dover, 1975 (McGraw-Hill: New York, 1955). (6) Peirce, A. P.;Dahleh, M. A.; Rabitz, H. Phys. Rev. 1988,A37,49504964. (7) Shi, S.;Woody, A.; Rabitz, H. J. Chem. Phys. 1988,88,6870-6883. (8) Kosloff, R.; Rice, S.A,; Gaspard, P.; Tersigni, S.;Tannor, D. J. Chem. Phys. 1989, 139, 201-220. (9) Jin, Y .;Tannor, D. J.; Soml6i, J.; Mrincz, A., private communication. (10) Warren, W. S.;Bates, J. L.; McCoy, M. A.; Navratil, M.;Mueller, L. J. Opt. Soc. Am. B Opt. Phys. 1986, 3,488492.
The Journal of Physical Chemistry, Vol. 97, No. 23, 1993 6183 (11) Warren, W. S.;Silver, M. S.Adu. M a p . Reson. 1988,12,247-384. (12) Warren, W. S.J. Chem. Phys. lW, 81, 5437-5448. (13) Lin, C. P.; Bates, J.; Mayer, J. T.; Warren, W . S.J . Chem. Phys. 1987,86, 375G3751. Erratum in J. Chem. Phys. 1987,87,4241. (14) Spano, F. C.; Warren, W. S.Phys. Reu. 1988, A37, 1013-1016. (15) Silver, M. S.;Joseph, R. I.; Hoult, D. I. Phys. Reu. 1985, A N , 27532755. (16) Tersigni, S. H.; Gaapard, P.; Rice, S.A. J. Chem. Phys. 1990, 93, 1670-1680. (17) Amstrup, B.; Carlson, R.; Matro, A,; Rice, S . A. J. Phys. Chem. 1991,95, 8019-8027. (18) Kosloff, R. J. Phys. Chem. 1988,92, 2087-2100. (19) Leforestier, C.; et al. J. Comput. Phys. 1991, 91, 59-80. (20) Luenberger, D. G.Linear and Nonlinear Programming, 2nd 4.; Addison-Wesley: Reading, MA, 1984. (21) Glowinski, R. Numerical Methods for Nonlinear VariationalProblems; Springer-Verlag, 1984. (22) Reas, W. H.; Flannery, B. P.;Teukolsky, S.A.; Vetterlig, W. T. Numerical Recipes. (Fortran Version);Cambridge University Press, 1989. (23) Tal-Ezsr, H.; Kosloff, R.; Cerjan, C. J . Comput. Phys. 1992, 100, 179-187.
