J. Phys. Chem. 1986, 90, 2264-2272
2264
TABLE VI: Measured Slope Enthalpies' of Zinc Sulfate DecomDositions
Zn0.2ZnS04 ZnSO, cell
ZnO-ZZnSO,
+ 4% Pt
AH',
TmIK
AH'r
TwIK
AH'r
T,,IK
7 3 2 5 1
254.7 257.2 257.4 261.7 263.4
817 800 796 795 799
278.7 253.0 258.8 261.8 262.3
868 860 848 831 844
236.4 245.2 245.6 249.0
812 796 787 789
equilibrium
232Sb
800
220.9c
850
213.0d
800
"In kJ-mol-'. bReaction 2. CReaction4. dReaction 5 been discussed in another paper,2 with reference to a similar mechanism for N 2 0 decomposition. It is worth noting here that Skeaff and Espelund'O found that additions of a few mole percent of Fez03 or Mn304were required in order to obtain steady and reproducible emf data on Z n S 0 4 and Zn0.2ZnS04, presumably because of induced catalytic effects. Although Fez03 and other transition-metal oxides were found to be ineffective in our studies a t 800-900 K, the emf measurements extended into the range 1000-1200 K. One can estimate the so-called vaporization coefficients associated with the several decomposition processes from the slope of the Whitman-M~tzfeldt'~-'~ plots. When the reciprocal of observed pressure is plotted against effective orifice area, =a, the intercept is the reciprocal of the equilibrium pressure, P,, and the slope is (aAPJ' where a is the vaporization coefficient and A is the sample surface area according to this approximate analysis. It is difficult to specify the true value of A , but if the cell cross sectional area is taken as a lower limit, then upper limit for reaction values of a at 800 K are calculated to be 5.3 X for reaction 4, and 6.0 X for the Pt-catalyzed 2, 1.5 X results on Zn0.2ZnS04. The kinetic barrier to the vaporization processes is such that the free surface vaporization fluxes of these
materials would be lower than the maximum kinetic theory values by approximately the derived values of a. These a values found in the present work are very similar in magnitude to those found for MgS0, decomposition.' In Table VI, the slope enthalpies derived from least-squares fitting of the vapor pressure data (where coefficient B = AH/ (2.303R)) for the three different samples are compared with each other and with the equilibrium values for the various reactions calculated from calorimetric and other thermodynamic data. Standard deviations for slope enthalpies can be derived from those for the B coefficients reported in Tables 1-111; the values range from 0.75 to 5.3 kJ-mol-I with a mean of 1.65 kJ.mol-'. With the exception of the smallest orifice (cell 7) data for Zn0.2ZnS04, the results generally show a decreasing AHo with decreasing hole size, in line with expectations as equilibrium is approached. Note, however, that all values of N o T are significantly larger than the equilibrium values, indicating the likelihood of a slow or ratelimiting step preceding the desorption of gaseous products in the decomposition process. This seems reasonable in view of the structural rearrangements associated with emission of SO3or SOz and O2 from the sulfate ion decomposition site. It may be significant that this difference between equilibrium and observed slope enthalpies was much smaller for MgSO, at 1000 K and that these processes with significant activation energies will manifest themselves increasingly at the lower temperatures. It should be noted that the true equilibrium enthalpies for the several decomposition reactions are probably larger than the values in Table VI by 10 to 12 kl-mol-' due to the expected small particle sizes as mentioned earlier, but this correction accounts for only a small part of the differences shown there. Acknowledgment. This research was supported by the Division of Chemical Science, Office of Basic Energy Sciences, U S . Department of Energy. Registry No. ZnSO,, 7733-02-0; ZnO-2ZnS04, 12037-14-8; Pt, 7440-06-4; CUO, 1317-38-0.
A Lie Approach to Global Senshivity Analysis of Systems Described by Ordinary Differential Equations Carl Wulfmant and Herschel Rabitz* Department of Chemistry, Princeton University, Princeton, New Jersey 08544 (Received: May 13, 1984)
A new form of global sensitivity analysis, based on Lie algebraic and group methods appropriate to many problems of chemistry, physics, and engineering, is outlined. An algorithm is provided that enables one to systematically determine, exactly or approximately, operators that convert solutions of ordinary differential equations containing parameters into families of solutions that arise when the parameters in the equation are changed in value. These techniques differ from conventional sensitivity analysis in that the new solutions generated by the present methods are guaranteed to satisfy the transformed differential equation. For example, if one knows a solution valid for some nonzero but restricted range of a parameter in the differential equation, then the operators can be used to obtain solutions valid outside of this range when such exist. If one knows a single solution valid for all values of the parameters, the operators can be used to obtain a many-parameter family of solutions valid for all values of the parameters. The operator determining equations of the present method are shown to reduce to those of conventional sensitivity analysis under appropriate restrictive conditions.
I. Introduction Most problems in chemistry, physics, and engineering may be expressed in terms of differential equations. The solutions to these equations, or their functionals, are quantities of physical interest. Typically, the underlying differential equations contain a number of parameters used to establish the physical model. In Some cases +Onleave from Department of Physics, University of the Pacific, Stockton, CA 95211.
0022-3654/86/2090-2264$01 .50/0
the model may be well-defined with parameters having essentially no uncertainty. In more typical circumstances, the model Parameters have significant uncertainties. Under both these conditions one would frequently like to know as much as possible about the solutions Of the differential equations as functions Of the Parameters. In the first case, the issue is not parameter uncertainty, but rather the desire to understand which portions of the model are responsible for Particular aspects of the observables. The second case involves the determination of the probability 0 1986 American Chemical Society
Lie Approach to Global Sensitivity Analysis distribution of observables from the prior known uncertainty distribution in the parameters. Problems of the type discussed above belong to the general realm of sensitivity analysis, which is a developing subject with a variety of practical and fundamental applications.’ This rather straightforward input/output sensitivity picture can become more involved and thereby take on a form outside of the most elementary applications of sensitivity analysis. For example, consider a coupled set of differential equations describing a flow reactor with N chemical species. Normally, this problem would be mathematically described in simplest terms as an initial value problem with N conditions specified at the input port of the reactor. However, a different situation arises if one desires to fix or optimize some of the product species leaving the reactor region. In particular, suppose that s of the products are desired to have fixed concentrations at the terminus of the reaction region. This constraint requires s of the input species now be allowed to float as unknowns; the problem then corresponds to the determination of the region in this s-dimensional space which gives rise to the desired fixed chemical outputs. Mathematically this problem may be viewed as a multidimensional boundary value problem which would be quite difficult to solve. From another perspective this situation may also be understood in terms of the desire to map out the s-dimensional initial concentration space and thereby identify the region of physical interest. An exactly analogous physical situation arises in semiclassical S-matrix scattering theory2 where one specifies the initial and final momenta on classical trajectories and desires to vary the initial coordinate phases in order to maintain the final momenta. Problems of this type are notoriously difficult to treat, but again they can be viewed as generalized sensitivity questions with the parameters being a subset of initial conditions. A problem related to that discussed in the last paragraph concerns how various parameters may be correlated in such a fashion as to produce exactly the same subset of observable values. For example, a solution to Schroedinger’s equation for inelastic or reactive scattering often produces a large number of distinct cross section^.^ Although presumably all of the cross sections together would imply a unique underlying Hamiltonian, this is likely not the case when only a subset are considered as available. Therefore, an important question involves seeking the family of Hamiltonians (or potential surfaces) that will give rise to exactly the same subset of observable cross sections. This question of uniqueness has a significant bearing on experimental design and the inversion of laboratory data back to the underlying Hamiltonian. Inversion questions of exactly the same nature arise in a variety of areas and in general the issue may be understood in terms of the goal of mapping out a region of parameter space while holding some constraints (observables) fixed. Sensitivity analysis may be conveniently divided into the broad categories of local and global question^.^ The issues raised in the first paragraph, especially concerning the understanding of which parameters control particular observables, may often be most prudently explored by the rapid computational techniques of local sensitivity analysis. Most prevalent amongst these techniques is the method of calculating gradients of the output with respect to the input variables and examining their magnitudes and trend^.^ Problems involving constrained or unconstrained mappings could in principle be treated through an extension of these techniques by considering higher order gradients in a Taylor series or its approximate analytic continued summation. This approach is generally a very tedious task when higher order terms in the series are needed. The traditional alternative of statistical ( I ) Tomovic, R.; Vukobratovic, M. “General Sensitivity Theory”; Elsevier: New York, 1972. (2) Miller, W. H. In “Advances in Chemical Physics”; Prigogine, I., Rice, S., Eds.; Wiley: New York, 1974; Vol. XXV. (3) Jochain, C. “Quantum Collision Theory”; Elsevier: New York, 1975. (4),Tilden, J.; Costanza, V.; McRae, G.; Seinfeld, J. In ‘Modelling of Chemical Reaction Systems”; Ebert, K., Deuflhard, P., Jager, W., Eds.; Springer-Verlag: West Berlin, 198 1. ( 5 ) Rabitz, H.; Kramer, M.; Dacol, D. Annu. Rev. Phys. Chem. 1983, 34, 419.
The Journal of Physical Chemistry, Vol. 90, No. 10, 1986
2265
Figure 1. The original and transformed solutions, 0, and 02, respectively, corresponding to eq 1.1 and 1.2 with parameter values specified in the text. The arrows connecting 0, and O2 show how some sample points are mapped from one solution to another.
or Monte Carlo-like sampling becomes prohibitively expensive for systems having more than a few parameters. Various importance-sampling techniques might be employed but again the expense can be very high. There seems to be no question that seeking global information will always be an arduous task, but presently available techniques may not be the optimal approach. This present paper considers a new view of these global mapping problems which shows promise for practical applications. It is too early at this stage to precisely estimate the labor involved in implementing these new techniques; however they present a perspective which deserves careful attention and exploration. The difference between previous approaches to global sensitivity analysis5 and that presented in this paper can perhaps be most easily understood with the aid of a simple example. Consider the equation dx/dt = -cx2
(1.1)
in which c is a parameter. In the usual approach to global sensitivity analysis one directly asks questions such as, “What is the effect on the solutions of (1.1) if one changes the parameter c to E?” In the approach set forth here we ask a question which is more restricted-in two senses. We first specify the manner in E = ac, with a being which c is to change, e.g., linearly as c a parameter. We then ask a question of the form, “If one changes c to ? = ac, in what manner must one change x to 1,and t to 7 so as to ensure that
-
dn/d7 =
?x2
(1.2)
when (1.1) is valid?” Insisting on specifying the manner in which c varies means that one is insisting on specifying a transformation. The insistence that (1.2) must follow from (1.1) is just the requirement that (1.1) is invariant under the transformation from c,x,t to ?,n,S. To appreciate the foremost consequences of asking our question in restricted form, let us determine the transformation of x,t required to convert (1 . l ) to (1.2) when we let ? = ac. Let us see if letting x = bx, t -*i = ht ( b and h are parameters), will work. (We note that there may be other such satisfactory transformations, and seeking out the largest number of “independent” transformations is a central part of the approach in this paper.) Making the substitutions -+
c=?/a,
x=.f/b,
t =i/h
(1.3a)
in (1.1) we obtain d(n/b) /d(i/ h) = -(?/a)(.f/ b)2
(1.3b)
i.e. (h/b)(dn/dt) = -(l/ab2)?X2
(1.3~)
2266 The Journal of Physical Chemistry, Vol. 90, No. 10, 1986
Wulfman and Rabitz In addition to its triviality, the example is atypical in a further way. We used some foreknowledge to guess a special form for the transformation of x and t required to compensate for the transformation of c and so leave the differential equation invariant. In the paragraphs below we will develop the logical procedure that determines the transformation of x and t required for each choice of transformation of c. As will become evident in later sections of the paper, choosing different transformations of c generaly transforms a solution of the initial equation into different solutions of the equation with the parameters altered. The use of invariance transformations of differential equations is, then, the first organizing principle of this paper. In local sensitivity analysis one considers derivatives of variables with respect to parameters. This, in the example above, is = c dc and Taylor-expanding the equivalent to letting c equation. In this paper we introduce more general types of infinitesimal changes. In the example of the previous section the infinitesimal change corresponding to c t = ac is c ? = (1 6a)c. As will appear in subsequent paragraphs, specification of this infinitesimal transformation of the parameters in the differential equation is, together with the invariance requirement, 3 = (1 + sufficient to ensure that there is a transformation x 6b)x, t i = (1 6h)t with
-
-
+
-
+
‘L
X Figure 2. Members of the one-parameter families of solutions to eq 1.1 with transformation c = ac, x = bc, and t = ht having the constraint abh = 1. The transformation parameter is a and new members of the family of solutions would correspond to different choices of b and h such that abh = 1. Some particular solutions are indicated by bold contours on the surface for selected values of a. We see that if h / b = l / a b 2 , Le., if abh = 1, eq 1.1 will indeed imply (1.2). Thus the transformations of c, x, and t are correlated in order to satisfy the criteria stated above. Now the solutions of (1.1) are x = l / ( c t + d) (1.4a) where d is a constant of integration. The solution curve for d = 1, c = 1, is plotted in Figure 1 and labeled as 8 , . The solutions of (1.2) may be written 2 = l / ( Z + d’) (1.4b) By way of illustration let us suppose our parameters have the values a = I / * , b = 2, h = 1. In Figure 1 we have plotted a series of points p obtained by moving points p on 8, with coordinates (x,t,c) to (Xi,?) with X = 2x, t = c/2, fX = t. The locus of the points p produces the curve which is obtained by substituting x/2 for x , 2? for c in the equation for 8, (1.4a). We thus have x / 2 = l(1
+ 2?i),
? =
i/*
(1.5a)
i.e. (1.5b) Comparing this with (1.4b) we see that the transformation has changed the constant of integration from d = 1 to d’ = I/*. We have converted a solution of (1. l ) , with parameter c = 1, into a solution of (1.2)-the same equation, but with a different value of the parameter (and, in this case, a different initial value). Note that on each curve the variables x,t and n,f, are generic, running, or dummy variables and the distinction between x and f and between t and 5 is relevant only to establishing the relation between the two curves. If we consider the parameter a as variable, then for each value of a and a choice for 6 and h such that abh = 1 , say the specific choice b = a-’ and h = 1, we obtain a transformed solution curve. In this manner, from a single solution curve of the original differential equation, one constructs a one-parameter family of solution curves of the transformed equation. This is indicated in Figure 2. Though the example is a trivial one, the result is general in the following sense: By requiring that the parameter be transformed in a specific way, and by requiring that the remaining variables be transformed so as to leave the differential equation invariant, one constructs a one-parameter family of solutions of the transformed equation from a single solution of the original differential equation.
-
+
(1
+ 6a)(l + 6b)(l + ah) = 1,
-
to O(6)
(1.6)
which converts solutions of (1.1) into solutions of (1.2). An emphasis on the primary role of infinitesimal transformations is the second organizing principle of this paper. Over a century ago Sophus Lie recognized the fundamental significance of the fact that the repetition of an infinitesimal transformation an arbitrarily large number of times gives rise to a one-parameter group of finite transformations.6 Lie established that if an infinitesimal transformation leaves a differential equation invariant at all points, then any finite transformation of the one-parameter group leaves the differential equation invariant, and vice versa. He developed a method for finding all the infinitesimal transformations that leave a given differential equation invariant, and he developed methods for finding the one-parameter group of transformations corresponding to any infinitesimal transformation. The global mapping procedure presented in this paper is thus based on the infinitesimal transformations of continuous groups which operate in the hyperspace defined by the system coordinates, state variables, and parameters. Considering this generalized space of all system variables the transformations change solutions of the original differential equation into solutions of the differential equation with altered parameters. Constraints on the physical system may be imposed either before relevant transformations are determined, or afterwards during their utilization. In many applications, examination of the infinitesimal transformations themselves gives valuable insight into the properties of the system, and the role of the parameters. In contrast to blind Monte Carlo parameter variation methods, the Lie approach will produce specialized transformation groups tailored to each particular physical or mathematical system. As one can choose to transform the parameters in an infinite variety of ways, there are an infinite number of one-parameter groups that carry out the corresponding invariance transformations. An important goal of this research is to determine appropriate independent sets from which all needed mappings may be created. Finally, we call attention to the fact that in Lie’s theory the generator of a one-parameter group plays a most useful role. The generator is a first-order differential operator (“Lie derivative”) of the form
u = p(z)(a/az’)
(1.7)
(Here and throughout the paper, we use the summation convention for repeated indices.) The quantum mechanical operators of linear ( 6 ) Lie’s collected works, Gesammelte Abhandlungen, Val. 1-7, are available from Johnson Reprint Corp., N.Y. His Differentialgleichungen, Vorlesungen Uber Continuerliche Grupen, and the three volumes of his Transformationsgruppen are all available from Chelsea Publishing Co., N.Y.
Lie Approach to Global Sensitivity Analysis and angular momentum are examples.' The finite transformations of a one-parameter group with parameter a and generator U may be carried out by the operator exp(aU). The role played by these exponential operators is very analogous to the role played by propagators in many quantum mechanical problem^.^ Indeed in the present case the finite transformation operators may be thought of as propagators in the hyperdimensional space of all system variables. With this analogy it is to be expected that many of the familiar approximations employed in quantum system propagation may also be applicable in the present circumstances. The approach presented here is a generalization of Lie techniques that have previously been found useful for investigating relations between the solutions of a given differential equation.*s9 However, the use of Lie groups and Lie algebras as an aid in analyzing the dependence of solutions of differential equations upon parameters in the differential equations is not in itself new. It has long been known that the dimensional analysis of differential equations provides a particularly simple example of Lie's analysis, and dimensional analysis often involves parameters in the equations.* Our analysis of eq 1.1 is a case in point. In general, dimensional consistency requires that every differential equation governing a physical system will admit transformations of the type considered in our example-combined dilatations of parameters and dynamical variables. Some years ago Bernal investigated Lie groups admitted by harmonic oscillators with force constants treated as variables (and arbitrary functions) on the same footing as the dynamical variables.I0 More recently Lie groups which vary parameters were introduced into the study of classical and quantum mechanical Hamiltonian systems, where they have special stability." Kitagawara utilized this analysis in a study of the dynamical origins of the regularities codified in the periodic chart of the elements.I2 An unpublished study of the problem of approximating the group generators admitted by a non-Hamiltonian system with limit cycle solutions showed that the approximation algorithm was stable to small variations in parameters in the equation and suggested the current in~estigati0n.l~We also call attention to a study by Ovsjannik~w'~ on a group theoretic analysis of nonlinear diffusion equations, treating an arbitrarily variable function on the same footing as other variables in the system. Finally, a recent paper by Huang et al.I5 considered a control theory approach to quantum mechanical wave packets with Lie techniques. The key innovations introduced in the present paper are two: 1. the introduction of the geometric viewpoints and methods of Lie into the foundations of sensitivity analysis; 2. the development and use of the theorem of section IV which enables one to systematize studies of parameter dependence and provides tests for the accuracy of approximations. The application of the theory to systems of ordinary linear differential equations is considered in a companion paper. 11. The Geometric View 1. The Spaces S of System Variables. Given an ordinary differential equation, or set of such equations, of arbitrary order we may rewrite them as a set of first-order equations in real varaibles dxJ/ds =f(x,s); j = 1, 2, ..., n' (2.1.1) (7) See, for example: Schiff, L. "Quantum Mechanics", 3rd ed.; McGraw-Hill: New Yuork, 1968; Chapter 7. (8) Bluman, G. W.; Cole, J. D. "Similarity Methods for Differential Equations"; Springer-Verlag: New York, 1974. (9) Ovsjanikov, L. V. "Group Properties of Differential Equations"; translated by G. W. Bluman, Department of Mathematics, University of British Columbia, Vancouver, 1962. ( I O ) Bernal, J. R.; M S c . Thesis, University of the Pacific, Stockton, CA, 1976. (1,l) Wulfman, C. E. In 'Recent Advances in Group Theory and Their Application to Spectroscopy"; Donini, J. C., Ed.; Plenum: New York, 1979. (12) Kitagawara, Y. MSc. Thesis, University of the Pacific, Stockton, CA, 1977. Cf. aiio ref 11. (13) Wulfman, C. E. "Machine Computed Approximate Generators Admitted by a Flow with a Limit Cycle", unpublished, 1980. (14) See ref 9. (15) Huang, G. M.; Tarn, T. J.; Clark, J. W. J . Math. Phys. 1983, 24, 2608.
The Journal of Physical Chemistry, Vol. 90, No. 10, 1986 2267
J
Figure 3. A schematic illustration of the projection of the vectors d x = f ( x , c ) dr onto the phase space defined with parameter value c,,. The vector projections at values c, and c2 produce a vector bundle in the phase space.
Temporal evolution problems in chemical kinetics and classical mechanics are already of this form, and others may be rewritten this way. Although partial differential equations will not be treated in this paper, their practical solution often reduces to solving equations of the above type. This is due to the fact that such problems are often spatially discretized to produce sets of ordinary differential equations. If the equations in (2.1.1) are not autonomous, or if the variable s is a periodic variable, we let ds = dt and write them as dxj/dt =f(x,s), ds/dt = 1 (2.1.2) We require t to run through all the reals. Referring to (2.1.1) and setting s = x" produces the standard autonomous set dxj/dt = j ( x ) ; j = 1, 2, ..., n; n = n'or n'+ 1 (2.1.3a) which is abbreviated as
x = f(x)
(2.1.3b)
We here require that the functions f are infinitely differentiable. A variable c may be treated as a parameter if it can be assigned an arbitrary value which does not change with time, dc/dt = & / a t + xj dc/dxJ = 0; j = 1, 2, ..., n (2.1.4) We then suppose that thef(x) are functions of a set of parameters ck, k = 1, 2, ...,N a n d that with respect to these parameters f is also infinitely differentiable. At each point p(x,c) in the space of parameters and state variables, we may consider the vector field dx defined by (2.1.5) dx = f(x,c) dt which may be projected onto the point p'(x) in the (x) subspace, the phase space, as schematically indicated in Figure 3. The usual phase space is generally defined with some fixed parameter value co. Through each point p'(x) there passes a phase curve of the flow whose tangent vector is defined by (2.1.5). As one varies c, the tangent vector p'(x) may be expected to vary and the phase curves through the point may be expected to change. They may in fact change from closed to open curves. No one-to-one continuous mapping can interconvert such topologically distinct curves. Consider however the variables XJ,Cand invertible one-to-one smooth transformations of XJ,Cthat carry contiguous points into contiguous points (diffeomorphisms). Through each point p(x,t) there may pass a solution curve with tangent vector defined by (dt, dx = f dt). As c is varied the tangent vectors and solution curves through p are, in general, varied. Now, however, because t ranges over all the reals there is a one-to-one mapping from every solution curve to the real line. There is thus no a priori reason for supposing that for those values of c for which solution curves exist they cannot be interconverted by diffeomorphisms of x,?, with the c's being treated as parameters in the sense previously defined.
2268 The Journal nf Physical Chemistry, Vol. 90, No. 10, 1986 We propose to study the effects of changes in the parameters c by investigating the mappings that convert solution curves into solution curves as one varies the c’s. There are no general algorithms for finding mappings that will interconvert solutions of differential equations. However Lie provided algorithms for finding all one-parameter groups of transformations that convert solutions of differential equations into solutions. The desire to have an algorithmic theory thus leads us to restrict our choice of mappings to those that may be carried out by a one-parameter continuous group of transformations or by a succession of such transformations. We will pay particular attention to situations in which a succession of one-parameter groups is a many-parameter group. The reason for this is that when a set of oneparameter groups comprise a many-parameter group, there exists a host of useful relations between different transformations of the system. However we must take cognizance of the fact that any differential equation admits an infinite number of one-parameter groups and only a small subset of these contain one-parameter groups that are interrelated i n such a way as to comprise many-parameter groups. We shall take the view that the effect of transformations upon the parameters c is at our command and time-independent while their effect upon the variables x and t is to be determined by making use of the requirement that the transformations convert solutions into solutions. We will show that once we have fixed our choice of the transformation of the parameters, the transformation generator U is essentially fixed and our theory becomes algorithmic. This perspective of the parameters as “controlling” the system is a natural one in many problems. There are cases in which this restricted view of the meaning of “parameter” would best be relaxed; however we will adhere to the restricted view in this paper. Because we wish to investigate the effect of any change of the parameters c, we insist that the succession of one-parameter groups be such as to enable the conversion of any c to any other c--that is, act transitively on (c),the space of all c’s. Finally, because the differential equations ,t = f ( x , c ) express the conditions that their solutions must satisfy, the list of variables is extended to ,t,x,c,t and we make the following definition: The spaces S are those defined by the real variables k,x,t,c and transformations of these sets into themselves by a succession of one-parameter continuous groups that act transitively in their c subspaces. The symbol Sois used to denote the subspace (x,t,c) of S, and S, is used to denote the subspace of (c). 2. One-Parameter Groups of’ Transformations. In this section we briefly review some salient properties of one-parameter transformation groups relevant to the following analysis. Extensive treatments of the material in this and in section III.1A can be found in the works of Lie,‘ in the introductory text of Cohen,I6 the text of Bluman and Cole,Rthe text of Ov~jannikov,’~ and in the book of Campbell.” To simplify notation we let (x,c,t)= ( z } ,with xl...x” = z’ ...z”; c‘ ...P = z~+’...z“++”;t = z”+”+I. We use a to denote the Lie group parameter-not contained in (21-and let zi Z’(z;a) represent a transformation with the following properties:
-
Z‘(z:0) = z’ Z’(Z(z;a);z’) = Z’(z;a’’),
(2.2.1) a” = a
+ a’
This is a one-parameter group of transformations of So in the standard form in which the identity transformation occurs when the parameter a is zero and the composition law is just the addition of parameters. We now define (2.2.2) (16) Cohen, A. “An Introduction to the Lie Theory of One-Parameter Groups”; Stechert: New York, 1931. (17) Campbell, J. E. “Introductory Treatise on Lie’s Theory of Finite Continuous Groups”; Oxford University Press: Oxford, England, 1903; reprinted, Chelsea. Bronx, N.Y., 1966
Letting a
-
Wulfman and Rabitz
-
6a gives the infinitesimal transformation zi
2’
=
2’
+ GaS(z)
(2.2.3)
The generator of the group is
u = S(Z) a/&’
a/axj
=
+
a/ack
+ T a/&
(2.2.4a)
where
r‘
=
UXJ,
T
=
ut, p = UCk
(2.2.4 b)
+
The operator 1 6aU is the operator of infinitesimal transformations: When it acts on any differentiable function g o n e has to first order in 6a (1
+ 6aU)g(z) = g((1 + 6aU)z) = g(2)
(2.2.5a)
where
z = (1
+ 6aU)z = z + 6a{(z)
(2.2.5 b)
Infinitesimal transformations themselves comprise a one-parameter group. We have, through first order in 6a
+ 6aU)(l + 6a‘U)z = (1 + 8a”J)t = i = 6a + sa’, 2 = z + 6af{(z), i = 2 + Ga{(r)
(1 where 6a”
(2.2.6a)
(2.2.6b)
Note that the last equation is obtained by treating Uas q(z) and therefore because of the group property the variables in U may be considered generic, running, or dummy variables. Thus, for example exp(-aU) exp(+aU)z = exp(-aU)? = z
(2.2.6~)
Also, the transformation inverse to (2.2.513) is z =P
-
ba{(t) = f - Gu{(z)
+ O(8a’)
(2.2.7)
If we have a function (2.2.8a)
g(z)
then as z = Z(z;-a), we must have, for an infinitesimal transformation g(Z(Z;-ba)) = g(z - Ga{(z))
(2.2.8b)
A function g(z) is invariant under the transformation z
g ( 4 = g(z)
-
2 if
(2.2.9)
for all z . If the function is a smooth one and the transformation is that of a one-parameter group with generator U, then g is invariant under the transformation if
Ug(z) = 0
(2.2.10)
for all z . The equation g(z) = 0 is invariant under the transformation if (2.2.10) is satisfied for all z for which g(z) = 0. Now consider the effect of one-parameter transformation groups upon infinitesimal displacements. Let dz be an infinitesimal displacement vector at p ( z ) . Let the transformation z i= Z(z;a) carry dz into d2 at ~ ( 2 ) The . corresponding infinitesimal transformation is ---t
dz’
-+
dZ’ = dz‘
+ 6 ~ - azj
dzj
(2.2.1 1)
This equation enables the determination of the effect of the transformation upon the ratio dxj:dt, which as Lie established, transforms in exactly the same manner as the derivative dxJ/dt, considered as a limit.’* (In other words, one need not consider xJ to be a function of t to derive the transformation law.) We henceforth treat xj interchangeably as the ratio dxj:dt and the derivative dxjldt. With these remarks in mind, the reader is (1 8) Lie, S. “Transformationsgruppen”; V.2, Abteilung I., Leipzig, 1890. Cf. also ref 9, p 35.
The Journal of Physical Chemistry, Vol. 90, No. 10, 1986 2269
Lie Approach to Global Sensitivity Analysis
G t
c
Figure 4. The subspace of x, x , and c of the differential equation x = Hypersurfaces such as this, but in a higher dimensional space, are defined by eq 2.1.3.
-cxz.
invited to show that as a becomes arbitrarily small, (2.2.1 1) implies that XJ = dxJ/dt is carried into XJ = dnj/di, where
and
D Dt
- = a/at
+
xi
a/ax‘
(2.2.12b)
We will also make use of the first “extension” or “prolongation” of U which is the Lie generator D=u+$ia/aXj (2.2.13) Later in this paper we will seek the generators 0 tha: leave the differential eq 2.1.3 invariant. Given such generators Uthe finite transformations they generate can be recovered in two ways. One can directly evaluate i = Z(z;a) = exp(aU)z (2.2.14) or obtain the finite transformation by integrating the set of ordinary differential equations (2.2.15)
-
Figure 5. A schematic illustration of an invariance transformation which takes a point p ( x , t ) into a point p = p(n,i) by virtue of the change c ?(a). Simultaneously, this may be thought of in terms of the tangent vector (dx,dt) to the solution at p being transformed into the tangent vector (d.t,di) at p’. Both tangent vectors can also be projected onto the surfaces t = 0, as illustrated.
define the original solution curves, then the inverse transformation converts them into the one-parameter family of curves defined by hi(.f,f,?:a)= g,’ (3.1.4a) where
hi(x,f,?,:a)= g‘(x(%,i,P;-a),t(R,t,?;-a),c(x,i,,T;-a)) (3.1.4b) It is worth noting that if eq 3.1.3 defines the general solution of the differential equations for all values of x,t,c, then it must be possible to rewrite eq 3.1.4b as
g(x,t,c) = C;(a),
CJ(0) = g;
(3.1.4~)
When the equations are written in this form the sole effect of the transformations is to change the value of the constants of integration 6,giving a one-parameter family of solutions with differing constants. For an infinitesimal transformation the equations of the transformed curves, (3.1.4b), become
g‘(x,t,c) - GaUg‘(x,t,c) = go‘, i = 1, 2, ..., N (3.1.5) As a’runs from 0 to a during the integration z’runs from z to 2. The same type of ordinary differential equation solvers capable of handling eq 2.12 would also be applicable here. Analogous relations hold for U in S, and all the relations of this section also hold for transformations of S for which U is the generator. 111. Invariance Transformations of the Differential Equations 1. Formulation and Interpretation. The set of differential eq 2.1.3 define a subspace of S, as suggested in Figure 4. The subspace is left invariant by a one-parameter Lie group with generator U if and only if D(X -f(x,c))l = 0 (3.1.1) Here the vertical line I is used to indicate that the equation is to hold identically for all values of the variables for which the original equations 2 = f ( x , c ) are satisfied. Thus one must have -f(x,c)) = -f(x,c))q (3.1.2)
o(x
where q is a function that has no factor (x - ~ ( x , c )with )~ b C 0. Note that U(x -f(x,c)) = 0 for all x,x is a special case of this. The transformations in Sogenerated by the U whose extension is Ucarry points p = p(x,t,c) into points p = p(%,t,T). At the same time the transformations, by virtue of (3.1.1) and (2.2.7), carry tangent vectors (dx,dt) at p’(x,t) defining solution curves through p’(x,t), into tangent vectors (dx,d?) at p’ = p’(%,f)defining solution curves through p’. These observations are illustrated in Figure 5. If the goi’s are a set of constants, and the equations
g‘(x,t,c) = g,’,
i = 1, 2, ..., N
(3.1.3)
The point p(x,2) obtained by infinitesimal transformation from p ( x , t ) will lie on a solution cruve through p ( x , t ) if ( E 6a, 7 ba), the tangent vector of the transformation, and Cfdt,dt), the tangent vector to the solution curve, are parallel at p(x,t). In the x,r subspace the variables x,t parameterize the original solution curve, while the variables %,tparameterize the transformed solution curve. Clearly, on any solution curve we may consider the variables x,t and n,i as generic ones so that the effect of the transformation is to convert a solution curve defined by (3.1.3) into a one-parameter family of solution curves defined by (3.1.4a). Phase curves ~ ( x c=) 0 are obtained by eliminating t from the solution curves defined by eq 3.1.3. Each phase curve is transformed into a one-parameter family of transformed phase curves $(%,?;a) = 0.
So far our remarks have been quite general. However, in this paper we shall henceforth consider each ck to be independent of x x, and t, and we shall only consider transformations f o r which the lClk depend solely upon c. This ensures that the transformed variables ck are not only independent of x and t , but that they are also independent of n,i. The subspace of (c) is therefore an invariant subspace in all that follows, and the space S is of a the structure particularly simple type. With this choice of the 1+5~ of the generators in (2.2.4a) will be shown to be dictated by the parameter excursions. 2. Derivation of the Determining Equations. We ?ow turn to the problem of determining the extended generators U and the base generators U of transformations that leave invariant the differential equations x = f(x,c). By the action of U on eq 2.1.3,
2270 The Journal of Physical Chemistry, Vol. 90, No. 10, 1986
Wulfman and Rabitz
one obtains for each component an equation
(p - U f ) l = 0
(3.2.1)
In view of eq 2.2.12 and 2.2.4a this gives
One also has the composition law = P”’(x)(“’’lPmlP”)
P“(x).P”(x)
To impose the requirement that this need vanish only for x =Ax) we eliminate factors of x from the equation by substitutingf(x) for iwherever it appears, and obtain the determining equations
(3.3.4)
These relations can now be used along with the orthonormality of the members of the basis to obtain an equation for each value of j
where j = 1, 2, ..., N (3.2.3) Despite the fact that eq 2.1.3 can be nonlinear, we see that the determining equations are a set of linear partial differential equations. Since no boundary or initial conditions are imposed the solutions depend upon arbitrary constants, or arbitrary functions, of integration. As a consequence, the integration of the determining equations yields functions l , T, $ connected via arbitrary constants (or functions) of integration in such a way that the admitted generator U is the form
u = cxru,
at
ax
”’ ax
- $k- aP = 0, j = 1, 2, ..., N
One may easily show that eq 3.3.la may be subsumed in the single equation (3.3.1 b)
where F is the evolution operator of (2.1.3), viz., F = ~ ( x , ca/axi )
There are a number of ways to solve or obtain approximate solutions to these equations. One can, for example, expandfand E in an appropriate basis, power series, Fourier series, etc. As a particular case consider the expansion in the orthonormal basis set (P‘(x))
P ( x ) = fiP’(x)
&‘(x,t)= t / ( t ) P ’ ( x )
ax + at
(3.3.2)
Kx = Ucf
(3.3.6a)
q k a/aCk
(3.3.6b)
where
uc
and the matrix K has elements M KJ,/J’,/‘
= I’U
N
- fi’ (Ilu/jlr’)1
(~Ill f~(’l‘’/ilu)’J’
The general solution of (3.3.6) is
x
= exp[K(t - to)] X ( t o )
+ s ‘ d s exp[K(t - s)]Uff
(3.3.7)
‘0
Solutions of the homogeneous equations determine initial condition transformations on 2.1.3. This topic will be discussed further in the companion paper. Here we shall focus on the particular solution to the inhomogeneous equation driven by Uff. It is important to also consider time-independent solutions of eq 3.3.6a
K x = Uff
(3.3.8)
In (3.3.8), and in the general case, the U have a particular form U, that depends upon the (free) choice $ : of the $k. Let us write
u, = u,x + u,c
(3.3.9a)
t; a/axJ
(3.3.loa)
ti‘ = Xi”P‘(X)
(3.3.1 Ob)
where u,x
(3.3.la)
aCk
au/at + [F,V = o
xJ9‘
(3.2.4)
where the A‘ are arbitrary constants or arbitrary functions and the U, are each independently admitted by the differential equations. Lie established that sets of ODE’Sadmit an infinite number of linearly independent U,, so that in practice it is important to solve the determining equations in such a way as to obtain a finite subset of solutions giving rise to the generators of the transformations of interest. We will make use of the freedom to choose initial conditions so as to bring this about. In the problems of interest here we want to obtain a set of one-parameter groups that act transitively in IC). It is preferable that this set of one-parameter groups is a many-parameter group, and for flexibility this many-parameter group should be as large as practically possible. Again, the freedom to impose conditions on the solutions of the determining equations will be of aid in bringing these desires to fruition. 3. Solution of the Determining Equations. The determining equations have the particular solution 7 = 1, l = 0 = $. That is, the system admits the generator U = a/&, the generator of time translation. (That an autonomous system is invariant under time translation is of course no surprise.) Thus the invariance operator exp(a 8/82) enables the arbitrary shift o f t to t + a. We therefore concentrate our attention on generators for which T = 0; nevertheless, the inclusion of transformations o f t can be important for some purposes, not considered here. The determining equations then simplify to
ap + --$ ap - p i
These equations are a set of inhomogeneous linear equations for the functions (f xJ*‘.They could have been written as the matrix representation of (3.3.la) on the chosen basis, P‘. However, we as choose to view them somewhat differently. Treating the the components of a vector x we have
=
with again using the expansion in (3.3.2). Here xi’ is thejl component of the x obtained when (3.3.7) is solved with the IC/k chosen to be .$ : Different choices for the expansion in (3.3.2) will produce results parallel to those obtained here. A power series expansion is considered in the accompanying paper. At this juncture it is appropriate to make a connection with conventional sensitivity analysis.’ The usual approach here consists of merely differentiating eq 2.1.3 with respect to ck. However, in the present paper it is more appropriate to consider differentiation with respect to the operator Uc = IC/k a/ack, Le., Lie differentiation, to give
We suppose that each P’(x) is a fixed function of x and that each
[/ is a function of t , and implicitly c, whose form is to be determined. Assuming the basis set is complete we have
where the x’s are considered to be functions of the e’s. Assuming
The Journal of Physical Chemistry, Vol. 90, No. 10, 1986 2271
Lie Approach to Global Sensitivity Analysis that eq 2.1.3 has been solved, the inhomogeneity in the present equation is known, and it is seen that these equations are linear in the unknown sensitivity coefficients UcxJ= qk (ax)/&) k . Equation 3.3.1 1 has a structural relation to the determining equations in (3.3.1), and the differences are largely associated with interpretation and procedural issues. First, the equations are mathematically equivalent when it is realized that the first two terms in eq 3.3.la may be rewritten to be the total time derivative occurring in the sensitivity eq 3.3.1 1. Furthermore, the sensitivity equations along with the original system in eq 2.1.3 effectively correspond to the method of characteristics for solving the determining equations. The central issue is how the determining equations and sensitivity equations are interpreted. In the latter case of sensitivity analysis, one is probing solution behavior around a nominal prescribed trajectory obtained from solving eq 2.1.3 with particular initial conditions. The inherently local nature of the analysis is apparent from examination of the vicinity of a given trajectory. In contrast, although the determining equations have incorporated eq 2.1.3 (cf., the transformation going from (3.1.2) to (3.2.3)), there is no longer any reference to a particular trajectory. Indeed, the perspective is one of exploring the state space with each of the variables xi treated as independent. This may also be understood considering that any point in the x space is accessible if all possible initial conditions are allowed. Therefore, although the determining equations are mathematically equivalent to those of sensitivity analysis, operationally they are treated in an entirely different fashion in order to obtain global information. This latter enhanced goal naturally leads to the more complex partial differential nature of the determining equations. The solution to these equations by repeated individual sensitivity trajectories does not appear to be an attractive approach given the goals of Lie theory. Other avenues and approximations are likely to be more productive. It is also anticipated that further insight into the similarities and differences between traditional sensitivity analysis and the present Lie approach will be most useful for fully understanding both techniques.
IV. Theorem The determining equations (3.3.1) are autonomous and so are themselves invariant under "time" translations. As a consequence t and the initial time to enter into their solution only through the difference At = t - to. We therefore are at liberty to consider the to be functions of the single variable At (and of course, x, and, implicitly, c). This observation lends utility to the following theorem which also has significant practical consequences. Theorem. Let the U,C be elements of a Lie algebra with commutation relations
e
[U,',U,'] = ep:U,C, p , q, r = 1, 2,
..., R
(4.1.1)
Here the e,; are the structure constants of a Lie algebra-an algebra which we are free to choose. This freedom arises because we may choose the coefficients $) as convenience dictates. Let the solutions in (3.3.7) establish a unique correspondence
U,C-U,,
r = 1 , 2 , ..., R
(4.1.2)
Then the full generators satisfy exactly the same commutation relations: (4.1.3) Proof. Since the Uc's are not functions of x, for each of them [a/axJ,Uc]vanishes for all j . Consequently since U, = U;X + U,C (4.1.4)
where and Up4only acts on the variable xi as a differential operator. Lie established that if two generators are admitted by a set of differential equations then so also is their c o m m ~ t a t o r . ' ~ From the
uniqueness of the correspondence between the U,C and the U, we must thus have UNX= ewrU?so that (4.1.3) follows from (4.1.1). (Lie also established that the commutation relations obeyed by extended generators are the same as those obeyed by the generators themselves.20) The theorem implies that the local Lie groups generated by the U,C, which act only in S,, and the U,, which act in So, are isomorphic. As a result of the liberty in choosing the transformation coefficients $, we may choose them to be those of some finite dimensional Lie group which acts transitively in S, and whose composition law is known. We may, for simplicity, choose the translation group in m variables, or for generality, choose the projective group in m variables. For the first, the order R of the group is m, while for the second it is ' 1 2 ( m+ l)(m + 2). The composition law for the translation group in m variables is in each dimension just that of the additive group of the real line, R(a?R(a) = R ( a a ? . Composition laws for the projective group in m variables have been discussed by Beckers et aL2' Knowledge of the composition law for the group of interest has the important feature that it enables the derivation of general relationships without knowing the exact form of the group generator. It enables one to choose a variety of paths through the space S each of which must lead from a given initial point p to a given final point p'. This has the consequence that if the generators or finite transformations are only known approximately, then comparing the coordinates of approximate points obtained in this way gives an estimate of the accuracy of the approximations. An estimate can also be obtained by making a corresponding comparison of infinitesimal transformations. The general operator of an R parameter group may also be written in the form (4.1.5) exp(a,Ul + a2U2 ... aRUR)
+
+ +
By choosing different values of the ratio a i / a Rone chooses operators that carry a point along different paths in the parameter space S, and the space So. Thus one has available an R - 1 parameter family of such paths. Each path explores the consequences of a differently coordinated change of parameters.
V. Invariance Transformations with Additional Requirements In this section we consider the following types of questions: 1. What changes in parameters will leave invariant a given property of the solutions of x = f ( x , c ) ? 2. What changes in parameters are most effective in changing some property of the solutions of x = f ( x c ) ? With the comments of the introduction in mind, these questions are seen to have a broad scope of applicability. A full treatment of them is not possible here, but some relevant aspects of the issues will be treated. Questions such as those above arise in a broad class of problems in which one has a specific goal while moving through the parameter space S,. A simple example might be the determination of the level curves of the functional J(x,c,t dt' = 0
or the level curves of a function. These questions introduce a possible new element into our previous considerations, for answering them can sometimes require relinquishing control over the transformation functions $k. In particular, a given question might only have an answer if $k is allowed to become a function of x,t. In this paper we will however only deal with cases in which $ can remain independent of x,t. There is another subtlety inherent in the questions of this section. It can happen that a transformation not admitted by a given set of equations is admitted if another equation or side condition is added to the set. Thus the function ( x ) ~ (y)2is not left invariant by the Lie transformation group with generator U = (1 - ( x ) ~ ) a / d x - ( x y ) a/dy. However if ( x ) ~+ Cy)2 is restricted to have the value 1, it is left invariant by such a transformation because
+
(20) For a modern proof, cf.: Ovsjannikov, L. V. ref 9, pp 38-40.
(21) Beckers, J.; Harnad, J.; Perroud, M.; Winternitz, P. J . Math. Phys. (19) Lie, S. 'Differentialgleichungen"; Kapitel 10, Leipzig, 1891.
1978, 19, 2126.
2272 The Journal of Physical Chemistry, Vol. 90, No. 10, 1986
+
+
U ( x ) * (y)2- 1)) = Z x ( ( x ) , b)2- 1) = 0. This paper will not consider cases in which the differential equations 1 = f ( x ) become restricted by auxiliary conditions in such a way as to admit an enlarged group necessary to the solution of the problem. Rather we will turn our attention to questions which can be answered with the aid of generators admitted by the equations x = f ( x ) alone. 1. Transformations That Leave Invariant a Function. Given the linear manifold of invariance generators U,available for a system, one can ask if there is any particular generator
U = Arc/,
(5.1.1)
of a one-parameter group that leaves invariant some given function F(x,t,c) of interest. If so,
0 = Ut: = ArUJ
X’F,
(5.1.2)
This equation may be solved for the coefficients X, with the aid of a projection of the functions F, onto a linearly independent set of functions G,. Letting Fr = (G,IFr)Gs
(5.1.3)
will yield a set of linear homogeneous equations for the A’:
(G,IF,)Ar = 0
(5.1.4)
One can also ask the “inverse” question, What functions are left invariant by a particular group with generator U? Let I be such a function so that we require UI = 0. This linear partial differential equation has as its characteristic equations
If these can be completely solved one obtains N + m - 1 functionally independent invariants I, of which m-1 depend only on C.
In a similar vein, we may also seek solutions of x = f ( x , c ) that are themselves invariant under the action of a Lie group with generator U. These are termed invariant solution^.^ Finally note that if the first order set x = f ( x , c ) has been obtained from a set of higher order equations then the invariant functions discussed in this section will, in general, be functions of derivatives in the original set, Le., differential invariants. 2. Transformations That Maximize the Change of Value of a Function. A transformation that leaves a function F(x,c) invariant moves a point in x,c space along a contour line of the function. At each point p ( x , c ) the vector of coefficients ( ( , T , $ ) in the generator U of the transformation is tangent to the contour. Given a set of one-parameter groups that act transitively in the space of x,t one can always find a linear combination of their generators that gave a generator U’whose coefficient vector points in any given direction in the space. In an N-dimensional space S , it is possible to find an N - I-dimensional subset of vectors each of which satisfy U’F = 0. If the (x,c,t)are coordinates of the point p in an orthogonal system of coordinates, each such C” will generate a transformation that moves points along paths orthogonal to the contour through p. Some linear combination of the U’must thus generate a transformation that moves p in a direction that will maximize the change in F. If the ( x , c ) are coordinates in an nonorthogonal system, then the argument must be given a geometrically evident modification. However, it should be noted that having a set of one-parameter groups that acts transitively in {c]by no means guarantees that one has a set that acts transitively in { x , t ]or { x ) . The usual situation will be that one will be able to use the group action to move points in some, but not all, directions “orthogonal” to the contours of constant value of F. The motion will necessarily be constrained by the differential equations and may be constrained by a lack of knowledge of a sufficiently large set of group generators. In this
Wulfman and Rabitz case the group action will only be able to be chosen to maximize the change in F subject to some constraint other than that required by the differential equations. One can always construct first-order differential operators that generate infinitesimal transformations orthogonal to the level curves of the function left invariant by a one-parameter group of transformations with generator U. VI. Conclusions The present paper should be viewed as an introduction of Lie’s geometric viewpoints and techniques into the study of sensitivity-related issues. Confining our attention to initial value problems described by eq 2.2.1 we seek group operators eauwith Lie generators U of the form given in eq 2.2.4a, which act on all the variables in the differential eq 2.2.1. The guiding philosophy throughout the paper is that of choosing, for the purpose at hand, the coefficients $k entering into the generators U. These are the coefficients that govern the transformation of the parameters in the differential equations. We demand that the operators e‘” transform solutions of the differential equations into solutions of the parametrically altered equations. We do this by demanding that the transformations carried out by the operators are such as to leave the equations invariant for all values of the parameters. We thereby obtain a set of linear equations, the Lie determining eq 3.3.1, that fix the form of the generators U. The general solution of these equations is given in eq 3.3.7. Again, taking the perspective that the coefficients fik are at our disposal, we prove in section IV that one is able to find invariance generators with a wide variety of commutation relations. This ensures that one can obtain a wide variety of many-parameter Lie groups, each of which enables one to change parameters in a host of differently coordinated ways. These then are the principal results and main working equations of the paper. An accompanying paper22illustrates their application to a system of coupled linear equations--the simplest general system of interest and one for which the determining equations may be solved analytically and all generators determined exactly. The task ahead is that of implementing the ideas for more complex systems-systems for which the governing differential equations are nonlinear. The Lie determining eq 3.3.1 remain linear and their general solution is given formally by (3.3.7). However, obtaining their general solution-but not some particular solutions-entails approximations. One is led to consider approximations to their solutions, (Le., approximations to group generators), and approximations to finite transformation operators. This we know involves a number of complexities. However, much is known about the approximate solution of linear equations. Furthermore, much of the research that has been done on approximating quantum mechanical evolution operators elHis directly applicable to the operators eau. The process of setting up the determining equations for any system is algorithmic and has been programmed by using symbol manipulating systems such as MACSYMA.’~ A computer program has been developed that both sets up and aprpoximates the solution to determining equations such as (3.3.1).24 In concluding, we note that the problems we have addressed in this paper are ubiquitous in chemistry, physics, and engineering, and that it appears the present approach holds considerable promise for advancing our understanding of them. Acknowledgment. We acknowledge many useful conversations and contributions from Drs. L. Hubbard and S. Shi. This work was supported by National Science Foundation Grant CHE8014165 and by the Office of Naval Research. _____.____
(22) Hubbard, L.; Wulfman, C.; Rabitz, H., following paper in this issue. (23) Wulfman, C. In “Group Theoretical Methods in Physics”; Serdaroglu, M., Inonu, E., Eds.; Springer-Verlag: New York, 1983; Lect. Notes Phys., Vol. 180, pp 187-91. (24) Nagao, G.M S c . Thesis, University of the Pacific, Stockton, CA, 1980.
