How do I do integral of an exponential matrix?

I am trying to find the integral of an equation that involves exponential matrices and I keep getting an error about "inner matrices must agree". Can any one see how to fix my problem?

Here is my function code: function [ y ] = Qkfunction(x)

F=[0 1 0;0 0 0;0 0 0];

L=[0 0;1 0;0 1];

Qc=[.2 0;0 .1];

dt=.01;

y=(expm(F*(dt-x)))*L*Qc*L'*(expm(F*(dt-x))');

end

and here is my script to take the integral:

dt=.01

Qk=integral(@Qkfunction,0,dt);

disp(Qk);

The error seems to occur with F*(dt-x) in the function but I do not know how to fix it since it should just be a matrix times a constant (dt-x) and not be an issue of matrix multiplication....

Answer 1

[Walter Roberson]

http://www.mathworks.com/help/matlab/ref/integral.html

and click on the "fun -- Integrand" link. It will show you

Integrand, specified as a function handle, defines the function to be integrated from xmin to xmax. For scalar-valued problems, the function y = fun(x) must accept a vector argument, x, and return a vector result, y. If you set the 'ArrayValued' option to true, fun must accept a scalar and return an array of fixed size.

Notice the "accept a vector argument".

[Savannah]

That worked! Thank you so much!!!

[Christoph Müller]

I'm relatively new to Matlab and have the same problem. However, despite the answer I cannot solve it. Could anyone please be more specific what the problem with the code is?

[Walter Roberson]

integral() passes a vector of x values of arbitrary length to the function and the function must calculate all of the results.

The code in the original question used

y=(expm(F*(dt-x)))*L*Qc*L'*(expm(F*(dt-x))');

and we have to read that with x being a vector of arbitrary length. With x being a vector, dt-x would be a vector of the same length. The code then has F*(dt-x) where F is a 3 by 3 matrix. The "*" operator is algebraic matrix multiplication, so that is a 3 x 3 matrix multiplied by a 1 x arbitrary matrix. Algebraic matrix multiplication requires that the "inner dimensions" must be the same, so if you use "*" to multiply an (m by n) matrix by (p by q) matrix, n must be the same as q, and the output size would be (m by q). But (3 x 3) * (1 x something) does not have "3" equal "1" so the "inner dimensions" do not agree and the matrix multiplication fails.

When you use integral() or integral2(), often you can simply vectorize the code. For example if you were integrating x*x - 5*x then that itself would fail because of the inner dimensions not agreeing for the "*" operation, but changing the expression to x.*x - 5*x would be all that would be needed to make it work. The .* operator is not algebraic matrix multiplication: it is multiplication of corresponding elements of two arrays of the same size (which does not have to be square).

But sometimes that is not enough. In the above example, there is no way to simply vectorize the F*(dt-x) expression. If you have a formula that calculates one scalar output value given one scalar input value but the calculation is too complicated to vectorize, then you will need to loop over the input x values calculating each in turn. The utility function arrayfun is often useful for that. For example instead of

y=(expm(F*(dt-x)))*L*Qc*L'*(expm(F*(dt-x))');

the original poster could have written

y = arrayfun(@(X) (expm(F*(dt-X)))*L*Qc*L'*(expm(F*(dt-X))'), x);

which calculates (expm(F*(dt-X)))*L*Qc*L'*(expm(F*(dt-X))') for each value in x individually.