Commit 5a7e535e authored by Jakob Gabriel's avatar Jakob Gabriel
Browse files

misc.place: different to builtin place, arbitrary eigenvalues can be set, that...

misc.place: different to builtin place, arbitrary eigenvalues can be set, that includes eigenvalues of multiplicity greater than rank(B)
+ some comments.
parent 2ed23e7c
...@@ -31,7 +31,9 @@ for n=1:size(A,1) ...@@ -31,7 +31,9 @@ for n=1:size(A,1)
adjA(m,n)=cofactor(n,m,A); % By switching indices, the resultant matrix is transposed. adjA(m,n)=cofactor(n,m,A); % By switching indices, the resultant matrix is transposed.
end end
end end
end % adj()
% Calculate elements of the cofactor-matrix % Calculate elements of the cofactor-matrix
function [a] = cofactor(i,j,A) function [a] = cofactor(i,j,A)
a=(-1)^(i+j)*det(A([1:i-1, i+1:end],[1:j-1,j+1:end])); a=(-1)^(i+j)*det(A([1:i-1, i+1:end],[1:j-1,j+1:end]));
end % cofactor()
function K = place(A, B, p)
%MISC.PLACE calculates the state feedback K such that all eigenvalues of (A-B*K) will result in
%p, i.e. eig(A-B*K) == p. Diffrent to the built-in place, this method supports eigenvalues of
%multiplicity greater than rank(B);
%
% Example
% K = misc.place(magic(2), [0; 1], [-1, -1])
% eig(magic(2) - [0; 1] * K)
arguments
A (:, :) double;
B (:, :) double;
p (:, 1) double = zeros([size(A, 1), 1]);
end
try
% try matlab built-in place first, as it is more accurate
K = place(A, B, p);
catch
assert(size(B, 2) == 1, "only tested for SISO-case yet");
K = zeros(size(B, 2), size(A, 1));
Atmp = A;
for it = 1 : size(A, 1)
[~, D, W] = eig(Atmp, "nobalance", "vector");
idx = chooseIdx(D, p);
if ~isempty(idx)
Ktmp = (W(:, idx)' * B) \ (D(idx) - p(it)) * W(:, idx)';
K = K + Ktmp;
Atmp = A - B*K;
else
break;
end
end
K = real(K);
end % try...catch
end % misc.place
function idx = chooseIdx(D, p)
% helper function
toBeMoved = setdiff(round(D, 15), round(p, 15));
if isempty(toBeMoved)
idx = [];
else
idx = find(round(D, 15) == toBeMoved(1), 1);
end
end
\ No newline at end of file
...@@ -38,7 +38,7 @@ if size(B,1)~=size(A,1) ...@@ -38,7 +38,7 @@ if size(B,1)~=size(A,1)
error('size(B,1)~=size(A,1)!') error('size(B,1)~=size(A,1)!')
end end
if size(B,2) ==1 if size(B, 2) == 1
K = place(A,B,pD); K = place(A,B,pD);
p = pD; p = pD;
else % mimo case: else % mimo case:
......
...@@ -4,6 +4,14 @@ function [tests] = testMisc() ...@@ -4,6 +4,14 @@ function [tests] = testMisc()
tests = functiontests(localfunctions()); tests = functiontests(localfunctions());
end end
function testPlace(tc)
A = magic(2);
B = ones(2, 1);
K = misc.place(A, B, 0*B);
tc.verifyEqual(eig(A - B*K), 0*B, "AbsTol", 1e-6);
end % testPlace()
function testReporter(testCase) function testReporter(testCase)
m = misc.Message("hallo"); m = misc.Message("hallo");
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment