Re: Re: Questions regarding MatrixExp, and its usage
- To: mathgroup at smc.vnet.net
- Subject: [mg63382] Re: [mg63355] Re: [mg63335] Questions regarding MatrixExp, and its usage
- From: Pratik Desai <pdesai1 at umbc.edu>
- Date: Wed, 28 Dec 2005 03:55:38 -0500 (EST)
- References: <BAY105-F241B75C5EDB643DA783A549A370@phx.gbl>
- Sender: owner-wri-mathgroup at wolfram.com
Michael Chang wrote:
> Hi,
>
>> From: Andrzej Kozlowski <akoz at mimuw.edu.pl>
To: mathgroup at smc.vnet.net
>> <michael_chang86 at hotmail.com>
>> Subject: [mg63382] Re: [mg63355] Re: [mg63335] Questions regarding MatrixExp,
>> and its usage
>> Date: Tue, 27 Dec 2005 15:52:53 +0900
>>
>> *This message was transferred with a trial version of CommuniGate(tm)
>> Pro*
>>
>> On 27 Dec 2005, at 09:58, Andrzej Kozlowski wrote:
>>
>>>
>>> On 27 Dec 2005, at 09:42, Andrzej Kozlowski wrote:
>>>
>>>>
>>>> On 27 Dec 2005, at 08:19, Andrzej Kozlowski wrote:
>>>>
>>>>> Now you should know that in general this is not going to hold
>>>>> in all of complex plane (but will hold in most).
>>>>
>>>>
>>>>
>>>> I wrote the above quite thoughtlessly: obviously there no sense in
>>>> which the equality "holds in most of the complex plane". Clearly
>>>> in every sense it holds just as often as it does not. Sorry about
>>>> that; I replied too quickly.
>>>> Andrzej Kozlowski
>>>>
>>>>
>>>
>>>
>>> One more correction is needed: in a certain obvious sense the
>>> equation does not hold "more often" than it holds since Exp is a
>>> surjective mapping of the compelx plane to itself which covers it
>>> infinitely many times (the fibre is Z - the integers).
>>>
>>> Because of the holidays I am now constanlty in a rush and can't
>>> find enough free time even to write a proper reply!
>>>
>>> Andrzej Kozlowski
>>
>>
>>
>> Actually, even the above is not strictly correct: Exp is a
>> surjective mapping form the complex plane to the complex plane minus
>> the point 0. Now that I have a little bit of time I can try to
>> analyse the entire problem more carefully. (In fact, I have not
>> taught complex analysis for over 15 years and I have become a little
>> bit rusty. So when I first saw this post I thought the problem lied
>> in the branch discontinuity of Log, which is why I wrote the
>> relations was true in "most of the complex plane". Of course I was
>> completely wrong in this respect).
>>
>> Let's again define the function
>>
>>
>> f[x_, y_] := E^(x*y) - E^(y*Log[E^x])
>>
>> We want to investigate where in the complex plane this is 0. This is
>> by definition the same as
>>
>>
>> f[x, y]
>>
>>
>> E^(x*y) - (E^x)^y
>>
>>
>> First, this is going to be zero for any real x and and an arbitrary y:
>>
>>
>> ComplexExpand[f[x,y],{y}]
>>
>> 0
>>
>> Secondly, suppose we have any pair of complex numbers a,b where
>> f[a,b] ==0. That is:
>>
>> a /: f[a, b] = 0;
>>
>> Then we have
>>
>>
>>
>> ExpandAll[FullSimplify[
>> f[a + 2*Pi*I, b]]]
>>
>>
>> E^(a*b + 2*I*b*Pi) - (E^a)^b
>>
>>
>> This will be zero if an only if b is an integer:
>>
>>
>> Simplify[%,Element[b,Integers]]
>>
>> 0
>>
>> So for every pair (a,b) for which the identity holds and b is not an
>> integer we can generate uncountably many pairs for which it does not
>> hold by simply adding 2*Pi*I to a. For example:
>>
>>
>> f[2,3/4]
>>
>> 0
>>
>>
>> FullSimplify[f[2 + 2*Pi*I,
>> 3/4]]
>>
>>
>> (-1 - I)*E^(3/2)
>>
>> On the other hand, we can get pairs of complex numbers for which the
>> identity holds provided the imaginary part of the first complex
>> number is not large:
>>
>>
>> f[1,2+3I]
>>
>>
>> 0
>>
>> However, for complex numbers with large imaginary part:
>>
>>
>> Simplify[f[1 + 12*I, 2 + 3*I]]
>>
>> (-E^(-34 + 27*I))*
>> (-1 + E^(12*Pi))
>>
>> it is easy in this way to give a complete description of the pairs
>> (a,b) for which f is 0, but I will skip it and turn to matrices.
>>
>> In this case, while I am not 100% sure, I tend to believe the
>> situation to be quite analogous. We are interested in the equation
>>
>> MatrixExp[B*p]==MatrixPower[MatrixExp[B],p]
>
>
> Many thanks to Pratik, Daniel, and Andrzej for their very insightful
> and expert feedback! :)
>
>> I believe this will hold for real matrices B and (probably) all
>> complex p but will not hold in general. In fact I believe most what
>> I wrote above can be generalised to this case, although the
>> statements and proofs would be more complicated.
>
>
> Hmm ... actually, from the sample example listed below, I don't
> believe that it will hold *in general* for real B *and* real p:
>
> In[1]: params={theta->Pi^Pi,p->Sqrt[2]};
> In[2]: B=theta {{Cot[theta],Csc[theta]},{-Csc[theta],-Cot[theta]}};
> In[3]: test1=Simplify[MatrixExp[B p]/.params];
> In[4]: test2=Simplify[MatrixPower[MatrixExp[B],p]/.params];
> In[5]: Simplify[test1 == test2]
> Out[5]: False
>
> Daniel has suggested that for (square matrix) B and (scalar) p both
> being real-valued, this only will hold if B is positive definite
> (although I suspect that this also may hold with B being positive
> semi-definite too). Using the above example:
>
> In[6]: BLim = Limit[B,theta->0];
> In[7]: Eigenvalues[BLim]
> Out[7]: {0, 0}
> In[8]: MatrixPower[MatrixExp[BLim],p]==MatrixExp[BLim,p]
> Out[8]: True
This is the special case where your matrix is Nilpotent refer to
http://mathworld.wolfram.com/MatrixExponential.html
The other thing that might be of interest is the diagonalizability of
the matrices
One can easily try this
Clear[aa, theta, p]
theta = Pi^Pi
p = Rationalize[1.414]
aa = theta {{Cot[theta], Csc[theta]}, {-Csc[theta], -Cot[theta]}};
{mata, matb} = JordanDecomposition[aa]
test1 = mata.MatrixExp[matb, p].Inverse[mata] // N // Chop
test3 = MatrixExp[aa, p]
test1 - test3 // Chop // MatrixForm
>>{{0, 0}, {0, 0}}
The advantage of using the Jordan Decomposition is that matb will be
diagonal for full Rank case. Hence matrix exponential is essentially the
exponential of each element of the diagonal
In fact this may be a way to find Powers of Matrix with non integer n
with out controversies and they hold
test4 = MatrixPower[aa, Sqrt[2]] // N
test5 = mata.MatrixPower[matb, Sqrt[2]].Inverse[mata] // N
test4 - test5 // Chop
>>{{0, 0}, {0, 0}}
>
> (By the way ... does anyone know *exactly* what the second argument
> for MatrixExp does? I've emailed Wolfram, since they only document
> MatrixExp with one argument, but I've seen their own documentation
> examples using *two* arguments; empirically, thus far, it seems that:
>
> MatrixExp[B,p]==MatrixExp[B p]
>
> with p being a scalar. But I digress ...)
>
> Anyways, for B and p real, I can 'sorta' see this point from what I
> (trivially) understand of the Spectral Mapping Theorem, since, as
> Andrzej has pointed out,
>
> Exp[a]^b !=Exp[a b]
>
> in general, with a complex, and b real; hence, any (strictly
> complex-valued) eigenvalue of
>
> MatrixPower[MatrixExp[B],p]
>
> will in general *not* be equal to
>
> MatrixExp[B p]
>
> Does this seem reasonable?
>
> Overall, too, I guess that I'm still kinda perplexed by what
> MatrixPower[B,pi] *means*? Somehow, I can feel 'comfortable' with
>
> MatrixExp[B pi]
>
> but not with
>
> MatrixPower[MatrixExp[B],pi]
>
> since I tend to think of MatrixExp[B pi] as Limit[MatrixExp[B
> t],t->Pi] (and can even revert back to an infinite power series matrix
> sum for an additional 'ease of understanding'), but, unlike general
> 'x^y' for scalars, can't quite grasp what the MatrixPower[*,*]
> equivalent signifies ...) :(
>
>> Let's just illustrate this in the case of a 2 by 2 random matrix.
>>
>>
>> B=Array[Random[Integer,{1,6}]&,{2,2}]
>>
>>
>> {{6,1},{5,1}}
>>
>> Let's take some complex p, e.g. 1+I
>>
>> In[65]:=
>> N[MatrixExp[B*(1 + I)]]==N[MatrixPower[MatrixExp[B],1+I]]
>>
>> Out[65]=
>> True
>>
>> To produce a case where the relationship does not hold just imitate
>> the procedure for complex numbers given above. First we add to 2Pi
>> * times the identity matrix to B:
>>
>> Z = B + 2 Pi*I IdentityMatrix[2];
>>
>> For p take any non-integer number, real or complex:
>>
>>
>> N[MatrixExp[Z*(1/2)]]==N[MatrixPower[MatrixExp[Z],1/2]]
>>
>>
>> False
>>
>> Andrzej Kozlowski
>
>
> Do any of my comments above make sense (or does anyone have a better
> explanation of what exactly how MatrixPower[B,p] can be interpreted
> with p not an integer)? My musings simply are from a 'layman's'
> perspective, and probably not very mathematically 'strict' ... :(
>
> Regards,
>
> Michael
>
>