Matrix Exponential to Approximate the Value of Matrizant

[Alexander122745]

Hello,

Consider the system of linear homogeneous differential equations of first order

dy/dx = A(x) y

where x denotes the independent variable, A(x) is a square matrix, and y is an unknown vector-function to be calculated.

Many-many years ago, when I was reading a good book on ordinary differential equations, I found a promising numerical method to approximately calculate the fundamental matrix (the matrizant or the Cauchy matrix) Ф(x) for the above linear system. The method is based on the replacement of the variable matrix A(x) by the constant matrix A(c) where c is the reasonably chosen numerical value for the variable x. Then the approximate value of the matrizant Ф(x) can be calculated by means of the formula

Ф(x) = Exp ( (x – a) A(c) )

where a is the initial value of x and Exp is so called matrix exponential that can be calculated by means of the existing very fast and precise computer algorithms.

Regrettably, I forgot the title of the book in which the idea of the method was briefly explained. It’s a petty because the method seems to be promising for stiff systems of differential equations.

My questions:

• Have you ever heard about the described method to approximate the matrizant?
• In what book or article was the method originally published?
• Is there any development or improvement of the original method to make it more accurate?
• Do you know any implementation of the method in any computing environment like Mathematica, Maple, or MATLAB?

My questions are caused by the necessity to accurately solve linear problems of hydrodynamic stability with a computer.

I would appreciate it if you could provide any additional information about the described numerical method. I mean names, titles, references, formulae, or Web links.

Thank you for your attention to this topic.

Respectfully,

Alexander

 

[Ray Vickson]

Hello,

Consider the system of linear homogeneous differential equations of first order

dy/dx = A(x) y

where x denotes the independent variable, A(x) is a square matrix, and y is an unknown vector-function to be calculated.

Many-many years ago, when I was reading a good book on ordinary differential equations, I found a promising numerical method to approximately calculate the fundamental matrix (the matrizant or the Cauchy matrix) Ф(x) for the above linear system. The method is based on the replacement of the variable matrix A(x) by the constant matrix A(c) where c is the reasonably chosen numerical value for the variable x. Then the approximate value of the matrizant Ф(x) can be calculated by means of the formula

Ф(x) = Exp ( (x – a) A(c) )

where a is the initial value of x and Exp is so called matrix exponential that can be calculated by means of the existing very fast and precise computer algorithms.

Regrettably, I forgot the title of the book in which the idea of the method was briefly explained. It’s a petty because the method seems to be promising for stiff systems of differential equations.

My questions:

• Have you ever heard about the described method to approximate the matrizant?
• In what book or article was the method originally published?
• Is there any development or improvement of the original method to make it more accurate?
• Do you know any implementation of the method in any computing environment like Mathematica, Maple, or MATLAB?

My questions are caused by the necessity to accurately solve linear problems of hydrodynamic stability with a computer.

I would appreciate it if you could provide any additional information about the described numerical method. I mean names, titles, references, formulae, or Web links.

Thank you for your attention to this topic.

Respectfully,

Alexander

To start: no, I don't know the answers to your specific questions, but I have spent some time, years ago, thinking about such issues.

Knowing that intuition can go badly wrong in such situations, I nevertheless assert that---intuitively---the method you outline might be "reasonable" for a slowly-varying matrix function $A(x)$. However, for better accuracy (again, intuitively) I would try to get a solution on interval $[a,b]$ by partitioning the interval into $[a,x_1], [x_1, x_2], \ldots, [x_n, b]$ and then using the outlined method on each sub-interval separately, stitching together the different pieces using continuity at the "join-points" $x_i.$ For "rapidly varying" $A(x)$ this would likely yield much better results.

As a Maple user, I usually do not bother with such fancy methods; as long as the system is not "too large" I would just go ahead and let Maple solve the system numerically if necessary. (However, just in case it works, I would first ask for a "symbolic" soltution, knowing full well that it is probably a waste of time.)

 

[Alexander122745]

To start: no, I don't know the answers to your specific questions, but I have spent some time, years ago, thinking about such issues.

Knowing that intuition can go badly wrong in such situations, I nevertheless assert that---intuitively---the method you outline might be "reasonable" for a slowly-varying matrix function $A(x)$. However, for better accuracy (again, intuitively) I would try to get a solution on interval $[a,b]$ by partitioning the interval into $[a,x_1], [x_1, x_2], \ldots, [x_n, b]$ and then using the outlined method on each sub-interval separately, stitching together the different pieces using continuity at the "join-points" $x_i.$ For "rapidly varying" $A(x)$ this would likely yield much better results.

As a Maple user, I usually do not bother with such fancy methods; as long as the system is not "too large" I would just go ahead and let Maple solve the system numerically if necessary. (However, just in case it works, I would first ask for a "symbolic" soltution, knowing full well that it is probably a waste of time.)

Dear Ray Vickson,
Thank you very much for reacting to my recent post at Physics Forum.
Alexander122745