Function to solve unknowns in Ka=b. Where K is a square matrix fully known, however, knowns and unknowns are scattered in a & b.

Hello,

I am trying to write a function to solve unknowns in Ka=b. The tricky part is unknowns are scattered in a & b and it has to work for different compatible matrix sizes. The inputs are K, a, b and I am trying to output the solved variables in a & b. Are there functions already similar to this? How would you go about this process wise?

Assume:

1. K is a square matrix fully known
2. a & b are column matrices with knowns and unknowns
3. dimensions can vary but assume they are always compatible.

[KSSV]

Showing one exmple will help.

[JT]

Right thats a good idea. Here is an example.

Ksys=[1,1,1,1,1,1;...

2,2,2,2,2,2;...

3,3,3,3,3,3;...

4,4,4,4,4,4;...

5,5,5,5,5,5;...

6,6,6,6,6,6];

F={1,2,'Ukn',4,'Ukn','Ukn'}';

Delta={'Ukn','Ukn',3,'Ukn',5,6}';

[F_sol,Delta_sol] = f_Ksys_solveF_Delta(Ksys,F,Delta);

Ideally, I am trying to make it be able solve symbolically where some of the "knowns" are symbolic.

I actually wrote the function out and it works for some scenarios but I am questioning its correctness in others. I have attached the code. If anyone has any tips to make it more efficient or if you see any situations where it wont work, your input is much appreciated. Thanks!

[Bruno Luong]

This is a careless example given: rank Ksys is 1, so not invertible.

[JT]

Apologies. Assume Ksys is an invertible matrix.

Answer 1

[Bruno Luong]

You might need to carry out some normalization of a, b so that they are unitless before doing this.

% Generate test data 
m = 12;
n = 10;
K = rand(m,n);
a = rand(n,1);
b = K*a;
known_a=rand(size(a))>0.4;
known_b=rand(size(b))>0.4;

% Solve a(~known_a) and b(~known_b)
% from K, a(known_a) and b(known_b)
tlsqflag = true; % recommend true
[m,n] = size(K);
M = [K, -eye(m)];
known_ab = [known_a(:); known_b(:)];
X = zeros(m+n,1);
X(known_ab,1) = [a(known_a); b(known_b)];
if ~tlsqflag
    % not recommended
    Y = -M*X;
    X(~known_ab) = M(:,~known_ab)\Y;
else
    % recommended method
    [~,~,V] = svd([M(:,~known_ab), M*X],0);
    Vend = V(:,end);
    X(~known_ab) = Vend(1:end-1)/Vend(end);
end
a = X(1:n)
b = X(n+1:end)

% Check residu
norm(K*a-b)/norm(b)