how can I calculate A^n when n is a symbolic positive integer?

Question

Songbai Jin 2023-1-23

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1898660-how-can-i-calculate-a-n-when-n-is-a-symbolic-positive-integer

评论： Songbai Jin 2023-1-25

A = [1/2,1/2,0,0;1/2,0,1/2,0;1/2,0,0,1/2;0,0,0,1]
syms n positive integer
A^n

This would not work because it got stuck and could never stop.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Walter Roberson 2023-1-24

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1898660-how-can-i-calculate-a-n-when-n-is-a-symbolic-positive-integer#answer_1154900

编辑：Walter Roberson 2023-1-24

在 MATLAB Online 中打开

syms N positive integer 
[V, D] = eig(sym(A)) 
result = V*diag(diag(D).^N)/V

Note those are the matrix operations * and / not element by element operations

8 个评论
显示 6更早的评论隐藏 6更早的评论

Walter Roberson 2023-1-25

在 MATLAB Online 中打开

The method does require that the matrix is diagonizable.

A = {[1/2,1/2,0,0;1/2,0,1/2,0;1/2,0,0,1/2;0,0,0,1], [1 3; 0 1]}
A = 1×2 cell array
    {4×4 double}    {2×2 double}
format long g
syms N positive integer 
for K = 1 : numel(A)
    [V, D] = eig(sym(A{K})); 
    if rank(V) ~= size(A{K},1)
        fprintf('diagonalization failed for matrix #%d\n', K);
    else
        result = V*diag(diag(D).^N)/V;
        result5 = subs(result,N,5);
        direct_version = sym(A{K})^5;
        fprintf('matrix #%d\n', K)
        double(result5 - direct_version)
    end
end
matrix #1
ans = 
       5.94278939548003e-71 -  6.91910479616603e-71i      -3.22608567183202e-71 +  1.17087100589517e-71i      -8.34820415079338e-72 -  6.46552480659025e-72i        7.3577392515467e-72 +  6.80939948034781e-75i
       2.16204338006988e-70 -  8.95310060594234e-71i       1.98092979849334e-71 +  1.79094167932363e-70i      -7.69732721700271e-71 -  2.37004100891168e-72i                          0 +  5.48748524668485e-75i
      -5.34851045593203e-71 -  2.10120053625901e-71i      -2.23562077258534e-71 +  4.81259819794677e-71i       5.23531446744669e-72 +  1.01119582283165e-72i       1.58474383879467e-71 -  5.10704961026086e-75i
                          0 +                     0i                          0 +                     0i                          0 +                     0i                          0 +                     0i
diagonalization failed for matrix #2

At the moment, I am not aware of any general formula for cases that are not diagonalizable.

Walter Roberson 2023-1-25

在 MATLAB Online 中打开

In cases where the matrix is numeric, you can approximate by using numeric eigenvalues, which removes problems with working with matrices 5 x 5 or larger (where it might not be possible to algebraically find the exact eigenvalues)

A = {[1/2,1/2,0,0;1/2,0,1/2,0;1/2,0,0,1/2;0,0,0,1], [1 3; 0 1]}
A = 1×2 cell array
    {4×4 double}    {2×2 double}
format long g
syms N positive integer 
for K = 1 : numel(A)
    [V, D] = eig((A{K})); 
    if rank(V) ~= size(A{K},1)
        fprintf('diagonalization failed for matrix #%d\n', K);
    else
        sV = sym(V);
        result = sV*diag(diag(D).^N)/sV;
        result5 = subs(result,N,5);
        direct_version = (A{K})^5;
        fprintf('matrix #%d\n', K)
        double(result5 - direct_version)
    end
end
matrix #1
ans = 
      -1.63526849241225e-16 - 2.96770310333167e-167i       -1.4643637431556e-17 + 2.47308591944306e-167i      -1.81451191767119e-16 - 3.09769808195142e-168i       3.59621678439901e-16 +                     0i
      -1.83427542145719e-16 + 5.72773437856958e-168i      -4.98395007893118e-17 + 2.38476142232009e-167i      -1.29593823909905e-16 -  1.0598939654756e-168i       3.62860866844936e-16 +                     0i
       -1.0923608999026e-16 - 1.69583034476095e-167i      -4.03961279131375e-17 - 1.48385155166583e-167i      -9.40591452819007e-17 +                     0i       2.43691363185299e-16 + 1.41319195396746e-168i
                          0 +                     0i                          0 +                     0i                          0 +                     0i                          0 +                     0i
diagonalization failed for matrix #2

Walter Roberson 2023-1-25

在 MATLAB Online 中打开

The result is too complex and always contains imaginary unit i, which should not appear in the result.

Let us have a look at determining the roots of a real-valued polynomial of degree 3 that has real roots

C = [1 2 -3 -4]

C = 1×4

1 2 -3 -4

syms a3 a2 a1 a0 x

p = a3*x^3 + a2*x^2 + a1*x + a0

p =

sols = solve(p, 'maxdegree', 3)

sols =

The above is the general solution, including the case where the coefficients might be complex valued and where the roots might be complex valued. Let us narrow it down by substituting in real coefficients:

subsols = subs(sols, [a3 a2 a1 a0], C)

subsols =

Okay, so you would not necessarily expect the imaginary components to vanish immediately upon substituting in real coeficients. Let us simplify() the expressions; that should the trick:

ssubsols = simplify(subsols)

ssubsols =

No, not good enough.

Ummm, maybe we need to expand the expressions and simplify the results of that?

sss = simplify(expand(ssubsols))

sss =

Nope, there are still imaginary components in that. We have a cube root of a complex-valued expression. Do we have an reason to expect to be able to rewrite that (19+6*sqrt(51)*i)^(1/3) in a separable form A+B*i so that we would have a hope of being able to rationalize the expressions to eventually symbolically get rid of all of the *i components? I do not think so.

What are the roots?

double(sss)
ans = 3×1
    1.5616
   -1.0000
   -2.5616

Or in the more direct form,

solve(poly2sym(C))

ans =

double(ans)

ans = 3×1

-1.0000 -2.5616 1.5616

So the roots are real-valued... but can we really expect to create a generalized formula that does not involve i anywhere that would return the real roots just given the coefficients and the promise that the roots are all real?

Songbai Jin 2023-1-25

在 MATLAB Online 中打开

Thank you for your help! It really helps me a lot!

Although I do not know how to remove i from A^N result, I find a way to get real root result.

C = [1 2 -3 -4];

syms a3 a2 a1 a0 x

p = a3*x^3 + a2*x^2 + a1*x + a0;

sols = solve(p, 'maxdegree', 3);

subsols = subs(sols, [a3 a2 a1 a0], C)

subsols =

result = simplify(subsols,'Steps',100)

result =

This method does not work in A^N result, probably because it needs more steps, or MATLAB just can not simplify the A^N answer.

Songbai Jin 2023-1-25

在 MATLAB Online 中打开

What's more, I find the mpower function implemented by MATLAB needs improving.

A = [1,3;0,1];

syms N positive integer

A^N

ans =

A = [1,1;0,0];

syms N positive integer

A^N

Error using ^
Singularity.

The function throws error if A is singular matrix, but it can calculate A^N if A is not diagonalizable.

请先登录，再进行评论。

Answer 2

KSSV 2023-1-23

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1898660-how-can-i-calculate-a-n-when-n-is-a-symbolic-positive-integer#answer_1154075

在 MATLAB Online 中打开

A = [1/2,1/2,0,0;1/2,0,1/2,0;1/2,0,0,1/2;0,0,0,1]

A = 4×4

0.5000 0.5000 0 0 0.5000 0 0.5000 0 0.5000 0 0 0.5000 0 0 0 1.0000

syms n positive integer

A.^n

ans =

3 个评论
显示 1更早的评论隐藏 1更早的评论

Songbai Jin 2023-1-24

I hope to calculate matrix multiplication A^n not elementwise product A.^n.

Thank you for your help.

KSSV 2023-1-24

在 MATLAB Online 中打开

A = [1/2,1/2,0,0;1/2,0,1/2,0;1/2,0,0,1/2;0,0,0,1]
n = 10 ; 
B = eye(size(A)) ; 
for i = 1:n
    B = B*A ; 
end

请先登录，再进行评论。

how can I calculate A^n when n is a symbolic positive integer?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

8 个评论
显示 6更早的评论隐藏 6更早的评论

更多回答（1 个）

3 个评论
显示 1更早的评论隐藏 1更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

how can I calculate A^n when n is a symbolic positive integer?

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

8 个评论 显示 6更早的评论隐藏 6更早的评论

更多回答（1 个）

3 个评论 显示 1更早的评论隐藏 1更早的评论

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

WeChat

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

8 个评论
显示 6更早的评论隐藏 6更早的评论

3 个评论
显示 1更早的评论隐藏 1更早的评论