Errors in using null command due to truncation error

17 Aug 2023
1 Answer
Updated 18 Aug 2023
2 Views (30 days)
Ran in:

0 votes

Ran in:
Hi,
I'm trying to find eigenvectors of a 9-by-9 square matrix, corresponding to its eigenvalues. The matrix consists of components with complex numbers and one symbolic 'a', so I found nine eigenvalues ('a'), via solving the determinant of the matrix. For some eigenvalues, I used the 'vpa' command since, without 'vpa', they are obtained as a form of 'root(eqn, z, integer)'. Here the issue seems to arise. Due to the truncation error, the 'null' shows empty eigenvectors corresponding to the eigenvalues. FYI, I don't know how to assign a specific variable in 'eig' and it takes forever to run. Is there a breakthrough other than Gauss elimination method with suppressing close-to-zero values?
clc;
close all;
clear;
syms a
R = 0.5234;
r = 0.0054;
s = 0.0084;
for p = 1:3
M = [a*R*r 1i*r 0 0 0 0 a*R*8.78+1i*(8.78-p*R) 0 0;
0 0 a*R*s 1i*s 0 0 0 a*R*79.88+1i*(79.88-p*R) 0;
0 0 0 0 a*R*0.4542 1i*0.4542 -(a*R*291.06+1i*(291.06+p*R)) -(a*R*291.06+1i*(291.06+p*R)) 0;
a*R*8.78+1i*(8.78-p*R) 0 0 0 0 0 0 0 a*R;
0 a*R*8.78+1i*(8.78-p*R) 0 0 0 0 0 0 1i;
0 0 a*R*79.88+1i*(79.88-p*R) 0 0 0 0 0 a*R;
0 0 0 a*R*79.88+1i*(79.88-p*R) 0 0 0 0 1i;
0 0 0 0 a*R*291.06+1i*(291.06+p*R) 0 0 0 a*R;
0 0 0 0 0 a*R*291.06+1i*(291.06+p*R) 0 0 1i];
Det = det(M);
DetEqn = Det == 0;
EigenVal1 = solve(DetEqn,a);
EigVal = vpa(EigenVal1);
for j=1:rank(M)
M_temp = subs(M,a,EigVal(j));
EigVec(:,j) = null(M_temp)
end
end
Unable to perform assignment because the indices on the left side are not compatible with the size of the right side.

Error in sym/privsubsasgn (line 1168)
L_tilde2 = builtin('subsasgn',L_tilde,struct('type','()','subs',{varargin}),R_tilde);

Error in indexing (line 999)
C = privsubsasgn(L,R,inds{:});

Sign in to answer this question.

Answers (1)

Walter Roberson
Walter Roberson on 17 Aug 2023

2 votes

Your code assumes that the null space is the same size each time, but most of the time the null space is empty. You cannot store an empty vector into a definite vector location.
You need to decide what you want to do when the null space is empty.

6 Comments
Show 4 older comments Hide 4 older comments

Walter Roberson
Walter Roberson on 17 Aug 2023
Ran in:
Ran in:
syms a
R = 0.5234;
r = 0.0054;
s = 0.0084;
for p = 1:3
M = [a*R*r 1i*r 0 0 0 0 a*R*8.78+1i*(8.78-p*R) 0 0;
0 0 a*R*s 1i*s 0 0 0 a*R*79.88+1i*(79.88-p*R) 0;
0 0 0 0 a*R*0.4542 1i*0.4542 -(a*R*291.06+1i*(291.06+p*R)) -(a*R*291.06+1i*(291.06+p*R)) 0;
a*R*8.78+1i*(8.78-p*R) 0 0 0 0 0 0 0 a*R;
0 a*R*8.78+1i*(8.78-p*R) 0 0 0 0 0 0 1i;
0 0 a*R*79.88+1i*(79.88-p*R) 0 0 0 0 0 a*R;
0 0 0 a*R*79.88+1i*(79.88-p*R) 0 0 0 0 1i;
0 0 0 0 a*R*291.06+1i*(291.06+p*R) 0 0 0 a*R;
0 0 0 0 0 a*R*291.06+1i*(291.06+p*R) 0 0 1i];
Det = det(M);
DetEqn = Det == 0;
EigenVal1 = solve(DetEqn,a);
EigVal = (EigenVal1);
for j=1:rank(M)
M_temp = subs(M,a,EigVal(j));
EV = null(M_temp);
if isempty(EV)
EigVec(:,j,p) = sym(NaN(size(EV,1),1));
else
EigVec(:,j,p) = EV;
end
end
end
format long g
EigVec = double(EigVec)
EigVec =
EigVec(:,:,1) = Columns 1 through 3 -0.0604434080666857 - 0.0568402099138265i -0.0604434080666857 + 0.0568402099138265i 0 + 0i -0.0568402099138265 + 0.0604434080666857i -0.0568402099138265 - 0.0604434080666857i 0 + 0i -0.00630053702136478 - 0.00625925383312013i -0.00630053702136478 + 0.00625925383312013i 0 + 0i -0.00625925383312013 + 0.00630053702136478i -0.00625925383312013 - 0.00630053702136478i 0 + 0i -0.00171477249055503 - 0.00171785608816912i -0.00171477249055503 + 0.00171785608816912i 0.998204973259795 + 0i -0.00171785608816912 + 0.00171477249055503i -0.00171785608816912 - 0.00171477249055503i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.06339171087373 + 0i NaN + 0i 0 + 0i 1 + 0i NaN + 0i 1.00659554466799 + 0i 0 + 0i NaN + 0i 1 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i Columns 7 through 9 NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i EigVec(:,:,2) = Columns 1 through 3 -0.0641387232677706 - 0.0564917511132487i -0.0641387232677706 + 0.0564917511132487i 0 + 0i -0.0564917511132487 + 0.0641387232677706i -0.0564917511132487 - 0.0641387232677706i 0 + 0i -0.00634195365255762 - 0.00625884452532305i -0.00634195365255762 + 0.00625884452532305i 0 + 0i -0.00625884452532305 + 0.00634195365255762i -0.00625884452532305 - 0.00634195365255762i 0 + 0i -0.00171169167540565 - 0.00171784779045346i -0.00171169167540565 + 0.00171784779045346i 0.996416379214725 + 0i -0.00171784779045346 + 0.00171169167540565i -0.00171784779045346 - 0.00171169167540565i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.13536440283453 + 0i NaN + 0i 0 + 0i 1 + 0i NaN + 0i 1.0132786693931 + 0i 0 + 0i NaN + 0i 1 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i Columns 7 through 9 NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i EigVec(:,:,3) = Columns 1 through 3 -0.0680253050205014 - 0.0558597772365388i -0.0680253050205014 + 0.0558597772365388i 0 + 0i -0.0558597772365388 + 0.0680253050205014i -0.0558597772365388 - 0.0680253050205014i 0 + 0i -0.00638363887049742 - 0.00625815577392186i -0.00638363887049742 + 0.00625815577392186i 0 + 0i -0.00625815577392186 + 0.00638363887049742i -0.00625815577392186 - 0.00638363887049742i 0 + 0i -0.0017086164151074 - 0.00171783399737567i -0.0017086164151074 + 0.00171783399737567i 0.99463418334813 + 0i -0.00171783399737567 + 0.0017086164151074i -0.00171783399737567 - 0.0017086164151074i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.21778690116231 + 0i NaN + 0i 0 + 0i 1 + 0i NaN + 0i 1.02005113025445 + 0i 0 + 0i NaN + 0i 1 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i 0 + 0i 0 + 0i NaN + 0i Columns 7 through 9 NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i NaN + 0i
Walter Roberson
Walter Roberson on 17 Aug 2023
Note that in this above code, the null() calculation is working on the symbolic solutions, so there is no truncation error going on.
Seung Hyeop Hyun
Seung Hyeop Hyun on 18 Aug 2023
Walter,
Thank you for leaving your comments. The issue here is that there must be eigenvectors corresponding to the eigenvalues, which can be obtained from the 'null', and also the dimension of the eigenvectors is definitely 9x1. The 'null' cannot find the eigenvectors and provides 0X1 empty vectors. I believe that it's because of the approximation due to using 'vpa'. The 'a' is already replaced by one of eigenvalues, and therefore, the null is not working on the symbolic solutions. Alternatively, I can do 'eig' for M after one of eigenvalues is entered.
Matt J
Matt J on 18 Aug 2023
Edited: Matt J on 18 Aug 2023
There doesn't appear to be any reason for you to be using symbolic operations for this problem. Eigenvalue problems generally do not have closed form solutions. Also, all of your problem input parameters (e.g. R,r,s) have a priori known values.
Torsten
Torsten on 18 Aug 2023
@Seung Hyeop Hyun
I don't know why you talk about "eigenvalues", but I agree that if "a" gives det(M(a)) = 0, null(M(a)) should be at least 1-dimensional and not empty.
Walter Roberson
Walter Roberson on 18 Aug 2023
Ran in:
Ran in:
@Seung Hyeop Hyun
You did not take into account that you use floating point constants and that some of the calculations take place in floating point instead of as symbolic numbers.
When you use symbolic numbers consistently then the problem does not show up.
Q = @(v) sym(v);
syms a
R = Q(5234)/Q(10)^4;
r = Q(54)/Q(10)^4;
s = Q(84)/Q(10)^4;
n8_78 = Q(878)/Q(10)^2;
n79_88 = Q(7988)/Q(10)^2;
n_4542 = Q(4542)/Q(10)^4;
n291_06 = Q(29106)/Q(10)^2;
for p = 1:3
M = [a*R*r 1i*r 0 0 0 0 a*R*n8_78+1i*(n8_78-p*R) 0 0;
0 0 a*R*s 1i*s 0 0 0 a*R*n79_88+1i*(n79_88-p*R) 0;
0 0 0 0 a*R*n_4542 1i*n_4542 -(a*R*n291_06+1i*(n291_06+p*R)) -(a*R*n291_06+1i*(n291_06+p*R)) 0;
a*R*n8_78+1i*(n8_78-p*R) 0 0 0 0 0 0 0 a*R;
0 a*R*n8_78+1i*(n8_78-p*R) 0 0 0 0 0 0 1i;
0 0 a*R*n79_88+1i*(n79_88-p*R) 0 0 0 0 0 a*R;
0 0 0 a*R*n79_88+1i*(n79_88-p*R) 0 0 0 0 1i;
0 0 0 0 a*R*n291_06+1i*(n291_06+p*R) 0 0 0 a*R;
0 0 0 0 0 a*R*n291_06+1i*(n291_06+p*R) 0 0 1i];
Det = det(M);
DetEqn = Det == 0;
EigenVal1 = solve(DetEqn,a);
EigVal = EigenVal1;
for j=1:rank(M)
M_temp = subs(M,a,EigVal(j));
EV = null(M_temp);
if isempty(EV)
EigVec(:,j,p) = sym(NaN(size(EV,1),1));
else
EigVec(:,j,p) = EV;
end
end
end
EigVec
EigVec(:,:,1) = 
EigVec(:,:,2) = 
EigVec(:,:,3) = 
format long g
EigVec = double(EigVec)
EigVec =
EigVec(:,:,1) = Columns 1 through 3 -0.0604434080666857 - 0.0568402099138265i -0.0604434080666857 + 0.0568402099138265i 0 + 0i -0.0568402099138265 + 0.0604434080666857i -0.0568402099138265 - 0.0604434080666857i 0 + 0i -0.00630053702136478 - 0.00625925383312013i -0.00630053702136478 + 0.00625925383312013i 0 + 0i -0.00625925383312013 + 0.00630053702136478i -0.00625925383312013 - 0.00630053702136478i 0 + 0i -0.00171477249055503 - 0.00171785608816912i -0.00171477249055503 + 0.00171785608816912i 0.998204973259795 + 0i -0.00171785608816912 + 0.00171477249055503i -0.00171785608816912 - 0.00171477249055503i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.06339171087373 + 0i -2.07917517265388 - 0.12335357734933i 0 + 0i 1 + 0i 2.08986241502569 + 0.131173171659899i 1.00659554466799 + 0i 0 + 0i -1.17431482491886 - 3.23249177899391i 1 + 0i 0 + 0i 1.16945872447306 + 3.25381182291119i 0 + 0i 0 + 0i 0.389692214298109 - 0.187973926292369i 0 + 0i 0 + 0i -0.392422256869245 + 0.187636508068213i 0 + 0i 0 + 0i 0.0467535877540691 + 0.00573056182938062i 0 + 0i 0 + 0i -0.15363302189684 + 0.127699515578846i 0 + 0i 0 + 0i 1 + 0i Columns 7 through 9 -2.07917517265388 + 0.12335357734933i -1.91125259969198 - 0.520193365405037i -1.91125259969198 + 0.520193365405037i 2.08986241502569 - 0.131173171659899i 1.91129494286941 + 0.553169312823223i 1.91129494286941 - 0.553169312823223i -1.17431482491886 + 3.23249177899391i -0.105661537118353 - 0.770245234313409i -0.105661537118353 + 0.770245234313409i 1.16945872447306 - 3.25381182291119i 0.0937570861782642 + 0.775325421161631i 0.0937570861782642 - 0.775325421161631i 0.389692214298109 + 0.187973926292369i 0.0716575680173222 - 0.185960430138347i 0.0716575680173222 + 0.185960430138347i -0.392422256869245 - 0.187636508068213i -0.0749584912828433 + 0.185626626193629i -0.0749584912828433 - 0.185626626193629i 0.0467535877540691 - 0.00573056182938062i 0.0363384221349967 + 0.0221560988269227i 0.0363384221349967 - 0.0221560988269227i -0.15363302189684 - 0.127699515578846i -0.00986540081559587 + 0.00258850280209401i -0.00986540081559587 - 0.00258850280209401i 1 + 0i 1 + 0i 1 + 0i EigVec(:,:,2) = Columns 1 through 3 -0.0641387232677706 - 0.0564917511132487i -0.0641387232677706 + 0.0564917511132487i 0 + 0i -0.0564917511132487 + 0.0641387232677706i -0.0564917511132487 - 0.0641387232677706i 0 + 0i -0.00634195365255762 - 0.00625884452532305i -0.00634195365255762 + 0.00625884452532305i 0 + 0i -0.00625884452532305 + 0.00634195365255762i -0.00625884452532305 - 0.00634195365255762i 0 + 0i -0.00171169167540565 - 0.00171784779045346i -0.00171169167540565 + 0.00171784779045346i 0.996416379214726 + 0i -0.00171784779045346 + 0.00171169167540565i -0.00171784779045346 - 0.00171169167540565i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.13536440283453 + 0i -1.03424396514104 - 0.0577669915193808i 0 + 0i 1 + 0i 1.04493120751284 + 0.0655865858299493i 1.0132786693931 + 0i 0 + 0i -0.589585462682334 - 1.60558586753832i 1 + 0i 0 + 0i 0.58472936223653 + 1.62690591145559i 0 + 0i 0 + 0i 0.193481085863487 - 0.0941556722582628i 0 + 0i 0 + 0i -0.196211128434622 + 0.0938182540341066i 0 + 0i 0 + 0i 0.0116310771845003 + 0.00138540923306503i 0 + 0i 0 + 0i -0.0380957587743605 + 0.0318852196758995i 0 + 0i 0 + 0i 1 + 0i Columns 7 through 9 -1.03424396514104 + 0.0577669915193808i -0.955605128257274 - 0.243608708993425i -0.955605128257274 + 0.243608708993425i 1.04493120751284 - 0.0655865858299493i 0.955647471434705 + 0.276584656411611i 0.955647471434705 - 0.276584656411611i -0.589585462682334 + 1.60558586753832i -0.058782994029221 - 0.382582523732593i -0.058782994029221 + 0.382582523732593i 0.58472936223653 - 1.62690591145559i 0.0468785430891321 + 0.387662710580816i 0.0468785430891321 - 0.387662710580816i 0.193481085863487 + 0.0941556722582628i 0.0341783223759006 - 0.0931471170415329i 0.0341783223759006 + 0.0931471170415329i -0.196211128434622 - 0.0938182540341066i -0.0374792456414217 + 0.0928133130968145i -0.0374792456414217 - 0.0928133130968145i 0.0116310771845003 - 0.00138540923306503i 0.00912923446797742 + 0.00536880051870816i 0.00912923446797742 - 0.00536880051870816i -0.0380957587743605 - 0.0318852196758995i -0.00244438930219064 + 0.000683128476076533i -0.00244438930219064 - 0.000683128476076533i 1 + 0i 1 + 0i 1 + 0i EigVec(:,:,3) = Columns 1 through 3 -0.0680253050205014 - 0.0558597772365388i -0.0680253050205014 + 0.0558597772365388i 0 + 0i -0.0558597772365388 + 0.0680253050205014i -0.0558597772365388 - 0.0680253050205014i 0 + 0i -0.00638363887049742 - 0.00625815577392186i -0.00638363887049742 + 0.00625815577392186i 0 + 0i -0.00625815577392186 + 0.00638363887049742i -0.00625815577392186 - 0.00638363887049742i 0 + 0i -0.0017086164151074 - 0.00171783399737567i -0.0017086164151074 + 0.00171783399737567i 0.99463418334813 + 0i -0.00171783399737567 + 0.0017086164151074i -0.00171783399737567 - 0.0017086164151074i 1 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 0 + 0i 1 + 0i 1 + 0i 0 + 0i Columns 4 through 6 0 + 0i 1.21778690116231 + 0i -0.685933562636758 - 0.0359047962427311i 0 + 0i 1 + 0i 0.696620805008562 + 0.0437243905532995i 1.02005113025445 + 0i 0 + 0i -0.394675675270157 - 1.06328389705312i 1 + 0i 0 + 0i 0.389819574824353 + 1.0846039409704i 0 + 0i 0 + 0i 0.12807737638528 - 0.0628829209135606i 0 + 0i 0 + 0i -0.130807418956415 + 0.0625455026894044i 0 + 0i 0 + 0i 0.00514395587663915 + 0.000594946349441334i 0 + 0i 0 + 0i -0.0167933653854523 + 0.0141531959048798i 0 + 0i 0 + 0i 1 + 0i Columns 7 through 9 -0.685933562636758 + 0.0359047962427311i -0.637055971112372 - 0.151413823522888i -0.637055971112372 + 0.151413823522888i 0.696620805008562 - 0.0437243905532995i 0.637098314289803 + 0.184389770941074i 0.637098314289803 - 0.184389770941074i -0.394675675270157 + 1.06328389705312i -0.0431568129995103 - 0.253361620205655i -0.0431568129995103 + 0.253361620205655i 0.389819574824353 - 1.0846039409704i 0.0312523620594214 + 0.258441807053877i 0.0312523620594214 - 0.258441807053877i 0.12807737638528 + 0.0628829209135606i 0.0216852404954267 - 0.0622093460092614i 0.0216852404954267 + 0.0622093460092614i -0.130807418956415 - 0.0625455026894044i -0.0249861637609478 + 0.061875542064543i -0.0249861637609478 - 0.061875542064543i 0.00514395587663915 - 0.000594946349441334i 0.0040759677296243 + 0.00231048172032896i 0.0040759677296243 - 0.00231048172032896i -0.0167933653854523 - 0.0141531959048798i -0.0010764184837051 + 0.000319388109428762i -0.0010764184837051 - 0.000319388109428762i 1 + 0i 1 + 0i 1 + 0i

Sign in to comment.

Categories

Find more on Linear Algebra in Help Center and File Exchange

Products

Release

R2017b

Tags

Asked:

Seung Hyeop Hyun
Seung Hyeop Hyun
on 17 Aug 2023

Commented:

Walter Roberson
Walter Roberson
on 18 Aug 2023

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!