I have this system of linear equations and i need to solve it for unknowns. how i can solve it in MATLAB? we have 8 unknowns and 8 equations.

6 views (last 30 days)
javeria
javeria on 25 Mar 2025
Answered: javeria on 16 Apr 2025
$$
\begin{aligned}
& b_{0}^{l}=a_{0}^{l} A_{00}^{l}+\sum_{n=a_{n}}^{\infty} a_{n}^{l} A_{0} n^{\ell}+Z_{0}^{l} \\
& b_{m}^{l}=a_{0}^{l} A_{m 0}^{l}+\sum_{n=1}^{\infty} a_{n}^{l} A_{m n}^{l}+Z_{m}^{l} \\
& a_{0}^{l}=b_{0}^{l} B_{00}^{l}+\sum_{m=1}^{\infty} b_{m}{ }^{l} B_{0 m}^{l}+c_{0}^{l} C_{00}^{l}+\sum_{m=1}^{\infty} c_{m}^{l} C_{0 m}^{l}+Y_{0}^{l} \\
& a_{n}^{l}=b_{0}^{l} B_{n 0}^{l}+\sum_{m=1}^{\infty} b_{m}^{l} B_{n m}^{l}+c_{0}^{l} C_{n 0}^{l}+\sum_{m=1}^{\infty} c_{m}^{l} C_{n m}^{l}+Y_{n}^{l}
\end{aligned}
$$
and:
$$
\begin{aligned}
& c_{0}^{l}=d_{0}^{l} D_{00}^{\ell}+\sum_{n=1}^{\infty} d_{n}^{l} D_{0 n}^{l}+X_{0}^{l} \\
& c_{m}^{l}=d_{0}^{l} D_{m 0}^{\ell}+\sum_{n=1}^{\infty} d_{n}^{l} D_{m n}^{\ell}+X_{m}^{\ell} \\
& d_{0}^{\ell} = b_{0}^{l} E_{00}^{\ell}+\sum_{m=1}^{\infty} b_{m}^{\ell} E_{0 m}^{l}+c_{0}^{l} F_{00}^{l}+\sum_{m=1}^{\infty} c_{m}^{l} F_{0 m}^{l}+W_{0}^{l} \\
& d_{n}^{l}=b_{0}^{l} E_{n 0}^{l}+\sum_{m=1}^{\infty} b_{m}^{l} E_{n m}^{l}+c_{0}^{l} F_{n 0}^{l}+\sum_{m=1}^{\infty} c_{m}^{l} F_{n m}^{l}+W_{n}^{l}
\end{aligned}
$$
Unknowns are:-
$$
\begin{gathered}
\left(a_{0}^{l}, a_{n}^{l}\right) ;\left(b_{0}^{l}, b_{n}^{l}\right),\left(c_{0}^{l}, c_{n}^{l}\right) and \\
\left(d_{0}^{l}, d_{n}^{l}\right) .
\end{gathered}
$$
The above system of equations which is linear one is the latex code for these equations and i need to solve it for these unknowns. I am confused that how i can solved this system. how i can write these equations in matrix form to solved it for the unknowns? Help me.
  7 Comments
Show 5 older comments
javeria
javeria on 25 Mar 2025
@Sam Chak we can defined the number for its truncation by taking a big number i.e
M=20;
N=50;
My question is that i need to solve this system of linear equations. so i need to write 1st the system in matrix i.e AX=B then i can solve by using X=A^{-1}B. But I am confused about that how to write this system in AX=B form.

Sign in to comment.

Sign in to answer this question.

Accepted Answer

javeria
javeria on 14 Apr 2025
clc;
close all;
clear all;
tic;
N = 20; % Number of terms in region E and I
M= 50; % Number of term in region B
RE = 1.0; % Radius of external cylinder
ratios = [1.031,1.5, 2.0, 4.0]; % Ratios RE/RI you want to test
b= 2.0*RE; % The distance from seabed to the bottom of cylinder (h-d)=b
h = 4.0*RE; % Water depth.
g = 9.81; % Gravity acceleration
% Preallocate wave numbers k
nu =0.01; % nu
k=linspace(0,0,N);
% % Use the algorithm of Chamberlin & Porter
a2=nu*h;
p=a2*(1-4*(1-(1+a2)*exp(-2*a2))/(2*a2+sinh(2*a2)))^(-0.25); % approximation for k_0 (Eq. 19)
k(1)=p;
for i2=1:N-1
b2=i2*pi-pi/2*tanh(2*a2/(i2*pi*pi));%aproximation for k_n (Eq. 26)
for j=1:2
u=b2/(1-2*(b2*tan(b2)+a2)*sin(2*b2)/(a2*(sin(2*b2)+2*b2)))^0.5; % From article (Eq. 23), an upper bound for k_n, iterated to refine the roots
b2=u;
end
k(i2+1)=b2;
end
k=k/h;
k(:,2:N)=k(:,2:N)*complex(0,1) % print all roots where first one is real and remaining ones are purely imaginary
% Compute the matrix A for different ratios
for r = 1:length(ratios)
RI = ratios(r) * RE; % Update RI based on the current ratio
%% finding the root \lemda_m =m*pi/b%%%%%%
for j=1:M
m(j)=(pi*(j-1))/b;
end
% Initialize matrix A
A = zeros(M, N);
% Set A_00 (first element)
A(1, 1) = cosh(k(1) * h) * sinh(k(1) * b) / (k(1) * b * (2 * k(1) * h + sinh(2 * k(1)* h))) * besselh(0, 1, k(1) * RE) / (k(1) * besselh(-1, 1, k(1) * RE));
% Set A_0n for j = 2:N
for j = 2:N
A(1, j) = cos(k(j) * h) * sin(k(j) * b) / (k(j) * b * (2 * k(j) * h + sin(2 * k(j) * h))) * besselk(0, k(j) * RE) / (-k(j) * besselk(-1, k(j) * RE));
end
% Set A_m0 for i = 2:M
for i = 2:M
A(i, 1) = 2 * cosh(k(1) * h) * (-1)^(i - 1) * k(1) * sinh(k(1) * b) / (b * (2 * k(1) * h + sinh(2 * k(1) * h)) * (k(1)^2 + m(i)^2)) * besselh(0, 1, k(1) * RE) / (k(1) * besselh(-1, 1, k(1)* RE));
end
% Set A_mn for i = 2:M, j = 2:N
for i = 2:M
for j = 2:N
A(i, j) = 2 * (-1)^(i - 1) * sin(k(j) * b) / (b * (2 * k(j) * h + sin(2 * k(j) * h)) * (k(j)^2 - m(i)^2)) * besselk(0, k(j) * RE) / (-k(j) * besselk(-1, k(j) * RE));
end
end
end
% this is for one value of nu and now i want it for nu = linspace(0.01, 10, 100);% Loop over the range of nu

More Answers (4)

Sam Chak
Sam Chak on 25 Mar 2025
I see, @javeria. So, you would like to learn how to solve a system of linear equations. You can refer to the example at the following link and begin writing the code accordingly.
https://www.mathworks.com/help/symbolic/solve-a-system-of-linear-equations.html
If you are able to create the equations in LaTeX, I believe you can certainly write those equations in MATLAB. You can also search for more examples on Google: "how to solve a system of linear equations with MATLAB".
  1 Comment
javeria
javeria on 25 Mar 2025
@Sam Chak for the simple equations i solved it, but here the system is very complicated because these coefficients of unknowns are also matrix not a single value. so i am little bit confuse that how i can define these equations in matlab.

Sign in to comment.

Torsten
Torsten on 25 Mar 2025
Edited: Torsten on 25 Mar 2025
@javeria
If you want to avoid writing your system as A*X = B, you can use a nonlinear solver like "fsolve". As the system is linear, it should converge in one step.
Otherwise - as @Sam Chak already said - you can use "equationsToMatrix" to let MATLAB compute A and B from your equations:
https://uk.mathworks.com/help/symbolic/sym.equationstomatrix.html
  28 Comments
Show 26 older comments
Torsten
Torsten on 3 Apr 2025
Did you copy my full code ?
I changed
b = sym('b', [1 M]); % b1, b2, b3
a = sym('a', [1 N]); % a1, a2
c = sym('c', [1 M]); % c1, c2, c3
d = sym('d', [1 N]); % d1, d2
to
a = sym('a', [N 1]); % a1, a2
b = sym('b', [M 1]); % b1, b2, b3
c = sym('c', [M 1]); % c1, c2, c3
d = sym('d', [N 1]); % d1, d2
Maybe this is the problem.

Sign in to comment.

javeria
javeria on 3 Apr 2025
Ran in:
% Define dimensions
M = 4; % For m = 0, 1, 2, 3
N = 3; % For n = 0, 1, 2
% Initialize coefficient matrices
A = zeros(M, N); % 4x3 (Aij)
B = zeros(N, M); % 3x4 (Bij)
C = zeros(N, M); % 3x4 (Cij)
D = zeros(M, N); % 4x3 (Dij)
E = zeros(N, M); % 3x4 (Eij)
F = zeros(N, M); % 3x4 (Fij)
% Define A (Aij)
A(1,1) = 2;
for j = 2:N
A(1,j) = 3;
end
for i = 2:M
A(i,1) = 4;
end
for i = 2:M
for j = 2:N
A(i,j) = 5;
end
end
% Define B (Bij)
B(1,1) = 1;
for j = 2:M
B(1,j) = 2;
end
for i = 2:N
B(i,1) = 3;
end
for i = 2:N
for j = 2:M
B(i,j) = 1;
end
end
% Define C (Cij)
C(1,1) = 2;
for j = 2:M
C(1,j) = 1;
end
for i = 2:N
C(i,1) = 0;
end
for i = 2:N
for j = 2:M
C(i,j) = 2;
end
end
% Define D (Dij)
D(1,1) = 1;
for j = 2:N
D(1,j) = 2;
end
for i = 2:M
D(i,1) = 3;
end
for i = 2:M
for j = 2:N
D(i,j) = 1;
end
end
% Define E (Eij)
E(1,1) = 1;
for j = 2:M
E(1,j) = 2;
end
for i = 2:N
E(i,1) = 3;
end
for i = 2:N
for j = 2:M
E(i,j) = 1;
end
end
% Define F (Fij)
F(1,1) = 2;
for j = 2:M
F(1,j) = 1;
end
for i = 2:N
F(i,1) = 0;
end
for i = 2:N
for j = 2:M
F(i,j) = 2;
end
end
% Define right-hand side vectors
Z = zeros(M, 1); % Z0, Z1, Z2, Z3
Y = zeros(N, 1); % Y0, Y1, Y2
X = zeros(M, 1); % X0, X1, X2, X3
W = zeros(N, 1); % W0, W1, W2
% Define Z (start at 1 and increment by 1)
Z(1) = 1;
for i = 2:M
Z(i) = 2;
end
% Define Y (start at 5 and increment by 1)
Y(1) = 5;
for i = 2:N
Y(i) = 9;
end
% Define X (start at 8 and increment by 1)
X(1) = 8;
for i = 2:M
X(i) = 9;
end
% Define W (start at 12 and increment by 1)
W(1) = 12;
for i = 2:N
W(i) = 10;
end
% Define identity and zero submatrices
I1 = zeros(M, M); % 4x4 identity
I2 = zeros(N, N); % 3x3 identity
I3 = zeros(M, M); % 4x4 identity
I4 = zeros(N, N); % 3x3 identity
for i = 1:M
I1(i,i) = 1;
I3(i,i) = 1;
end
for i = 1:N
I2(i,i) = 1;
I4(i,i) = 1;
end
O1 = zeros(M, M); % 4x4 zero
O2 = zeros(M, N); % 4x3 zero
O3 = zeros(N, M); % 3x4 zero
O4 = zeros(N, N); % 3x3 zero
O5 = zeros(M, M); % 4x4 zero
O6 = zeros(M, N); % 4x3 zero
O7 = zeros(N, N); % 3x3 zero
% Assemble the full matrix G
G = [I1, -A, O1, O2;
-B, I2, -C, O4;
O5, O6, I3, -D;
-E, O7, -F, I4];
% Assemble the vector b
y = [Z; Y; X; W];
b = sym('b', [1 ,M]); % b0, b1, b2, b3
a = sym('a', [1, N]); % a0, a1, a2
c = sym('c', [1 ,M]); % c0, c1, c2, c3
d = sym('d', [1, N]); % d0, d1, d2
x = [b, a, c, d]'; % Vector of unknowns
% Solve the system
x_solution = inv(G) * y;
% Display all components
% Display results
disp('Matrix M (14x14):'); disp(M);
Matrix M (14x14): 4
disp('Vector b (14x1):'); disp(b);
Vector b (14x1):
disp('Solution vector x:'); disp(x_solution);
Solution vector x: -3.1481 -11.8201 -11.8201 -11.8201 -10.3598 2.7619 2.7619 19.6878 6.4444 6.4444 6.4444 -3.3598 3.7619 3.7619
% Extract bm, an, cm, dn from x_solution
bm = x_solution(1:M); % b0, b1, b2, b3
an = x_solution(M+1:M+N); % a0, a1, a2
cm = x_solution(M+N+1:M+N+M); % c0, c1, c2, c3
dn = x_solution(M+N+M+1:end); % d0, d1, d2
disp('bm (b0, b1, b2, b3):'); disp(bm);
bm (b0, b1, b2, b3): -3.1481 -11.8201 -11.8201 -11.8201
disp('an (a0, a1, a2):'); disp(an);
an (a0, a1, a2): -10.3598 2.7619 2.7619
disp('cm (c0, c1, c2, c3):'); disp(cm);
cm (c0, c1, c2, c3): 19.6878 6.4444 6.4444 6.4444
disp('dn (d0, d1, d2):'); disp(dn);
dn (d0, d1, d2): -3.3598 3.7619 3.7619
when i run my this code online it works but on my matlab it give me error i donot why
  20 Comments
Show 18 older comments
Torsten
Torsten on 15 Apr 2025
% % This is the code for finding roots of the dispersion relation
% % K = k * tanh(kh)
% % here K = nu^2/g, nu is angular frequency and h is water depth
% function [ k1 ] = dispersion( nu,h,N)
RE = 1.0; % Radius of external cylinder
ratios = [1.031, 1.5, 2.0, 4.0]; % Ratios RE/RI you want to test
b = 2.0 * RE; % Distance from seabed to the bottom of cylinder (h-d)=b
h = 4.0 * RE; % Water depth
N=3; % how many roots you need? Change your N to get roots
M = 4;% Change your M according to u
Nu = [0.1;0.2]; % angular frequency
h=4.0; % water depth
%% finding the root \lemda_m =m*pi/b%%%%%%
for j=1:M
m(j)=(pi*(j-1))/b;
end
for nn = 1:numel(Nu)
nu = Nu(nn);
k = zeros(1, N);
% % Use the algorithm of Chamberlin & Porter
a2=nu*h;
p=a2*(1-4*(1-(1+a2)*exp(-2*a2))/(2*a2+sinh(2*a2)))^(-0.25); % approximation for k_0 (Eq. 19)
k(1)=p;
for i2=1:N-1
b2=i2*pi-pi/2*tanh(2*a2/(i2*pi*pi));%aproximation for k_n (Eq. 26)
for j=1:2
u=b2/(1 - (2*(b2*tan(b2) + a2)*sin(2*b2))/(a2*(sin(2*b2)+2*b2)))^0.5; % From article (Eq. 23), an upper bound for k_n, iterated to refine the roots
b2=u;
end
k(i2+1)=b2;
end
k=k/h;
k(:,2:N)=k(:,2:N)*complex(0,1); % print all roots where first one is real and remaining ones are purely imaginary
for r = 1:length(ratios)
RI = RE / ratios(r); % Update RI based on the current ratio
% Example: Initialize matrix A as a zero matrix of size NxN
A = zeros(M, N);
% Set A_00 (first element)
A(1, 1) = cosh(k(1) * h) * sinh(k(1) * b) / (k(1) * b * (2*k(1)*h + sinh(2 * k(1)* h))) * besselh(0, 1, k(1) * RE) / (k(1) * besselh(-1, 1, k(1) * RE));
% Set A_0n for j = 2:N
for j = 2:N
A(1, j) = cos(k(j) * h) * sin(k(j) * b) / (k(j) * b * (2*k(j)*h + sin(2*k(j) * h))) * besselk(0, k(j) * RE) / (-k(j) * besselk(-1, k(j) * RE));
end
% Set A_m0 for i = 2:M
for i = 2:M
A(i, 1) = 2 * cosh(k(1) * h) * (-1)^(i - 1) * k(1) * sinh(k(1) * b) / (b*(2 * k(1) * h + sinh(2 * k(1) * h)) * (k(1)^2 + m(i)^2)) * besselh(0, 1, k(1) * RE) / (k(1) * besselh(-1, 1, k(1)* RE));
end
% Set A_mn for i = 2:M, j = 2:N
for i = 2:M
for j = 2:N
A(i, j) = 2 * (-1)^(i - 1) * sin(k(j) * b) / (b * (2 * k(j) * h + sin(2 * k(j) * h)) * (k(j)^2 - m(i)^2)) * besselk(0, k(j) * RE) / (-k(j) * besselk(-1, k(j) * RE));
end
end
% After calculating matrix A, you can process or print it as required
fprintf('Matrix A for ratio %.3f and RI %.3f:\n', ratios(r), RI);
disp(A);
end
end

Sign in to comment.

javeria
javeria on 16 Apr 2025
@Torsten @Sam Chak I need some help regarding my research I am a Phd student. I derive analytical equations for my system to find the unknown coefficient however my results are not corrected. If i provide my analytical expression can you do help me in coding part. Or if you have mathematics background can you also help me in analytical part.
Thanks

Categories

Find more on Applications in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!