finding non zero solution for homogeneous system of equations

6 次查看（过去 30 天）
MADHVI
MADHVI 2024-6-24
编辑： John D'Errico 2024-6-24
The following code is giving constant values of Nusselt number. The solutions are approximately close to zero.
D = 0.1;
L = 0.1;
B = 0.1;
xi = 0.1;
sigma1 = B^3 .* D^2 .* L^3;
sigma2 = D^2 .* L^3;
Ra_values = (10:100:100000).';
solutions = zeros(length(Ra_values), 7);
%Nu = zeros(length(Ra_values),1);
for i = 1:length(Ra_values)
Ra = Ra_values(i);
R = Ra * xi;
f1 = @(A) -((5.*A(2).*A(6).*pi^5.*B^2.*D^2.*L^4)./8 + (A(2).*A(6).*pi^5.*D^4.*L^4)./4)./sigma1 - (B .* ((A(3).*A(7).*D^4.*pi^5)/4 - (A(1).*L^4.*pi^4)./2 - (32.*A(5).*L^4.*Ra)./9 + (A(1).*L^4.*R.*pi^2)./2 + (5.*A(3).*A(7).*D^2.*L^2.*pi^5)./8))./sigma2;
f2 = @(A) -(B^2 .* (A(2).*A(5).*D^2.*L^4.*pi^5 - (A(2).*D^2.*L^4.*pi^4)./2 - (A(1).*A(6).*D^2.*L^4.*pi^5)./8 - (4.*A(6).*D^2.*L^4.*Ra)./9 + (A(4).*A(7).*D^2.*L^4.*pi^5)./2 + (3.*A(4).*A(7).*D^4.*L^2.*pi^5)./16 + (A(2).*D^2.*L^4.*R.*pi^2)./4) - (A(2).*D^4.*L^4.*pi^4)./4)./sigma1 - (B .* ((A(4).*A(7).*D^4.*pi^5)./8 - (A(2).*L^4.*pi^4)./4 - (16.*A(6).*L^4.*Ra)./9 + (A(2).*L^4.*R.*pi^2)./4 + (5.*A(4).*A(7).*D^2.*L^2.*pi^5)./16))./sigma2;
f3 = @(A) (B .* ((16.*A(7).*L^4.*Ra)./9 + (A(3).*D^4.*pi^4)./4 + (A(3).*L^4.*pi^4)./4 - (A(3).*L^4.*R.*pi^2)./4 + (A(3).*D^2.*L^2.*pi^4)./2 + (A(1).*A(7).*D^2.*L^2.*pi^5)./8 - A(3).*A(5).*D^2.*L^2.*pi^5 - (A(4).*A(6).*D^2.*L^2.*pi^5)./2))./sigma2 - (B^2 .* ((3.*A(4).*A(6).*pi^5.*D^4.*L^2)./16 + (5.*A(4).*A(6).*pi^5.*D^2.*L^4)./16) + (A(4).*A(6).*D^4.*L^4.*pi^5)./8)./sigma1;
f4 = @(A) (B .* pi^2 .* (2.*A(4).*D^4.*pi^2 - 2.*A(4).*L^4.*R + 2.*A(4).*L^4.*pi^2 + 4.*A(4).*D^2.*L^2.*pi^2 + A(2).*A(7).*D^2.*L^2.*pi^3 - 8.*A(3).*A(6).*D^2.*L^2.*pi^3 - 8.*A(4).*A(5).*D^2.*L^2.*pi^3))./(16.*D^2.*L^3) - ((B^2 .* pi^2 .* (2.*A(4).*D^2.*L^4.*R - 4.*A(4).*D^2.*L^4.*pi^2 - 4.*A(4).*D^4.*L^2.*pi^2 + 8.*A(2).*A(7).*D^2.*L^4.*pi^3 + A(2).*A(7).*D^4.*L^2.*pi^3 - A(3).*A(6).*D^2.*L^4.*pi^3 + A(3).*A(6).*D^4.*L^2.*pi^3 + 8.*A(4).*A(5).*D^2.*L^4.*pi^3))./16 - (A(4).*D^4.*L^4.*pi^4)./8)./sigma1;
f5 = @(A) (B^2 .* ((pi^5.*A(2)^2.*D^2.*L^4)./4 + (pi^5.*A(4)^2.*D^4.*L^2)./4 + (pi^5.*A(4)^2.*D^2.*L^4)./8) + (A(2)^2.*D^4.*L^4.*pi^5)./4 + (A(4)^2.*D^4.*L^4.*pi^5)./8)./sigma1 + (B .* ((A(3)^2.*D^4.*pi^5)/4 + (A(4)^2.*D^4.*pi^5)./8 + (8.*A(1).*L^4.*Ra)./9 + 8.*A(5).*L^4.*pi^4 - 2.*A(5).*L^4.*R.*pi^2 + (A(3)^2.*D^2.*L^2.*pi^5)./4 + (A(4)^2.*D^2.*L^2.*pi^5)./8))./sigma2;
f6 = @(A) (B^2 .* ((4.*A(2).*D^2.*L^4.*Ra)./9 + 2.*A(6).*D^2.*L^4.*pi^4 + (A(1).*A(2).*D^2.*L^4.*pi^5)./8 + (A(3).*A(4).*D^2.*L^4.*pi^5)./16 + (3.*A(3).*A(4).*D^4.*L^2.*pi^5)./16 - (A(6).*D^2.*L^4.*R.*pi^2)./4) + (A(6).*D^4.*L^4.*pi^4)./4)./sigma1 + (B .* (4.*A(2).*L^4.*Ra)./9 + 4.*A(6).*L^4.*pi^4 + (A(3).*A(4).*D^4.*pi^5)./4 - A(6).*L^4.*R.*pi^2 + (A(3).*A(4).*D^2.*L^2.*pi^5)./4)./sigma2;
f7 = @(A) (B .* ((4.*A(3).*L^4.*Ra)./9 + (A(7).*D^4.*pi^4)./4 + 4.*A(7).*L^4.*pi^4 - A(7).*L^4.*R.*pi^2 + 2.*A(7).*D^2.*L^2.*pi^4 + (A(1).*A(3).*D^2.*L^2.*pi^5)./8 + (A(2).*A(4).*D^2.*L^2.*pi^5)./16))./sigma2 + (B^2 .* ((3.*A(2).*A(4).*pi^5.*D^4.*L^2)./16 + (A(2).*A(4).*pi^5.*D^2.*L^4)./4) + (A(2).*A(4).*D^4.*L^4.*pi^5)./4)./sigma1;
F = @(A) [f1(A); f2(A); f3(A); f4(A); f5(A); f6(A); f7(A)];
A0 = [0.01 0 0 0 0.01 0 0];
A = fsolve(F, A0);
solutions(i, :) = A;
F(A)
end
A1_1_1 = solutions(:, 1);
A2_1_1 = solutions(:, 5);
Nu = xi - 1./2 - (pi^3 ./ Ra_values) .* A1_1_1 - (8 .* pi^3 ./ Ra_values) .* A2_1_1;
%Nu(i) = xi - 1./2 - (pi^3 ./ Ra_values(i)) .* A(1) - (8 * pi^3 ./ Ra_values(i)) .* A(5);
%end
plot(Ra_values, Nu);
xlabel('Ra');
ylabel('Nu');

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

回答（1 个）

Aquatris
Aquatris 2024-6-24
编辑：Aquatris 2024-6-24
Ran in:
Since zeros(7,1) is kind of a solution for you, you should select your initial value away from that point. Below you can find a simple adjustment to your code that shows some of the nonzero solutions. If you want to find nonzero solutions for all Ra_values, then you should implement a simple algorithm that changes the initial values for specific Ra until a nonzero solution is found
D = 0.1;
L = 0.1;
B = 0.1;
xi = 0.1;
sigma1 = B^3 .* D^2 .* L^3;
sigma2 = D^2 .* L^3;
Ra_values = (10:100:100000).';
solutions = zeros(length(Ra_values), 7);
%Nu = zeros(length(Ra_values),1);
options = optimoptions('fsolve','Display','none');
for i = 1:length(Ra_values)
Ra = Ra_values(i);
R = Ra * xi;
f1 = @(A) -((5.*A(2).*A(6).*pi^5.*B^2.*D^2.*L^4)./8 + (A(2).*A(6).*pi^5.*D^4.*L^4)./4)./sigma1 - (B .* ((A(3).*A(7).*D^4.*pi^5)/4 - (A(1).*L^4.*pi^4)./2 - (32.*A(5).*L^4.*Ra)./9 + (A(1).*L^4.*R.*pi^2)./2 + (5.*A(3).*A(7).*D^2.*L^2.*pi^5)./8))./sigma2;
f2 = @(A) -(B^2 .* (A(2).*A(5).*D^2.*L^4.*pi^5 - (A(2).*D^2.*L^4.*pi^4)./2 - (A(1).*A(6).*D^2.*L^4.*pi^5)./8 - (4.*A(6).*D^2.*L^4.*Ra)./9 + (A(4).*A(7).*D^2.*L^4.*pi^5)./2 + (3.*A(4).*A(7).*D^4.*L^2.*pi^5)./16 + (A(2).*D^2.*L^4.*R.*pi^2)./4) - (A(2).*D^4.*L^4.*pi^4)./4)./sigma1 - (B .* ((A(4).*A(7).*D^4.*pi^5)./8 - (A(2).*L^4.*pi^4)./4 - (16.*A(6).*L^4.*Ra)./9 + (A(2).*L^4.*R.*pi^2)./4 + (5.*A(4).*A(7).*D^2.*L^2.*pi^5)./16))./sigma2;
f3 = @(A) (B .* ((16.*A(7).*L^4.*Ra)./9 + (A(3).*D^4.*pi^4)./4 + (A(3).*L^4.*pi^4)./4 - (A(3).*L^4.*R.*pi^2)./4 + (A(3).*D^2.*L^2.*pi^4)./2 + (A(1).*A(7).*D^2.*L^2.*pi^5)./8 - A(3).*A(5).*D^2.*L^2.*pi^5 - (A(4).*A(6).*D^2.*L^2.*pi^5)./2))./sigma2 - (B^2 .* ((3.*A(4).*A(6).*pi^5.*D^4.*L^2)./16 + (5.*A(4).*A(6).*pi^5.*D^2.*L^4)./16) + (A(4).*A(6).*D^4.*L^4.*pi^5)./8)./sigma1;
f4 = @(A) (B .* pi^2 .* (2.*A(4).*D^4.*pi^2 - 2.*A(4).*L^4.*R + 2.*A(4).*L^4.*pi^2 + 4.*A(4).*D^2.*L^2.*pi^2 + A(2).*A(7).*D^2.*L^2.*pi^3 - 8.*A(3).*A(6).*D^2.*L^2.*pi^3 - 8.*A(4).*A(5).*D^2.*L^2.*pi^3))./(16.*D^2.*L^3) - ((B^2 .* pi^2 .* (2.*A(4).*D^2.*L^4.*R - 4.*A(4).*D^2.*L^4.*pi^2 - 4.*A(4).*D^4.*L^2.*pi^2 + 8.*A(2).*A(7).*D^2.*L^4.*pi^3 + A(2).*A(7).*D^4.*L^2.*pi^3 - A(3).*A(6).*D^2.*L^4.*pi^3 + A(3).*A(6).*D^4.*L^2.*pi^3 + 8.*A(4).*A(5).*D^2.*L^4.*pi^3))./16 - (A(4).*D^4.*L^4.*pi^4)./8)./sigma1;
f5 = @(A) (B^2 .* ((pi^5.*A(2)^2.*D^2.*L^4)./4 + (pi^5.*A(4)^2.*D^4.*L^2)./4 + (pi^5.*A(4)^2.*D^2.*L^4)./8) + (A(2)^2.*D^4.*L^4.*pi^5)./4 + (A(4)^2.*D^4.*L^4.*pi^5)./8)./sigma1 + (B .* ((A(3)^2.*D^4.*pi^5)/4 + (A(4)^2.*D^4.*pi^5)./8 + (8.*A(1).*L^4.*Ra)./9 + 8.*A(5).*L^4.*pi^4 - 2.*A(5).*L^4.*R.*pi^2 + (A(3)^2.*D^2.*L^2.*pi^5)./4 + (A(4)^2.*D^2.*L^2.*pi^5)./8))./sigma2;
f6 = @(A) (B^2 .* ((4.*A(2).*D^2.*L^4.*Ra)./9 + 2.*A(6).*D^2.*L^4.*pi^4 + (A(1).*A(2).*D^2.*L^4.*pi^5)./8 + (A(3).*A(4).*D^2.*L^4.*pi^5)./16 + (3.*A(3).*A(4).*D^4.*L^2.*pi^5)./16 - (A(6).*D^2.*L^4.*R.*pi^2)./4) + (A(6).*D^4.*L^4.*pi^4)./4)./sigma1 + (B .* (4.*A(2).*L^4.*Ra)./9 + 4.*A(6).*L^4.*pi^4 + (A(3).*A(4).*D^4.*pi^5)./4 - A(6).*L^4.*R.*pi^2 + (A(3).*A(4).*D^2.*L^2.*pi^5)./4)./sigma2;
f7 = @(A) (B .* ((4.*A(3).*L^4.*Ra)./9 + (A(7).*D^4.*pi^4)./4 + 4.*A(7).*L^4.*pi^4 - A(7).*L^4.*R.*pi^2 + 2.*A(7).*D^2.*L^2.*pi^4 + (A(1).*A(3).*D^2.*L^2.*pi^5)./8 + (A(2).*A(4).*D^2.*L^2.*pi^5)./16))./sigma2 + (B^2 .* ((3.*A(2).*A(4).*pi^5.*D^4.*L^2)./16 + (A(2).*A(4).*pi^5.*D^2.*L^4)./4) + (A(2).*A(4).*D^4.*L^4.*pi^5)./4)./sigma1;
F = @(A) [f1(A); f2(A); f3(A); f4(A); f5(A); f6(A); f7(A)];
A0 = randn(7,1)*1e3;
A = fsolve(F, A0,options);
solutions(i, :) = A;
end
A1_1_1 = solutions(:, 1);
A2_1_1 = solutions(:, 5);
Nu = xi - 1./2 - (pi^3 ./ Ra_values) .* A1_1_1 - (8 .* pi^3 ./ Ra_values) .* A2_1_1;
nonzeroSolutions = solutions(all(abs(solutions)>1e-8,2),:)
nonzeroSolutions = 557x7
-11.1439 -2.4841 0.5966 -0.3815 0.1024 -0.1500 0.0153 -1.1769 1.1402 0.9476 -0.1999 -0.1267 0.0423 -0.1175 -161.9110 -3.7509 -0.9152 15.7362 3.5740 -39.5678 39.4107 -5.9693 -4.2321 0.1918 0.0064 -0.0264 -0.8981 -0.0094 -6.0161 0.2271 -0.9818 -4.0342 -0.5380 1.5874 -1.7717 -7.1687 5.3215 2.1127 -0.2700 -0.7532 0.3789 -0.1887 -10.5027 -7.6824 -1.8260 -0.1935 -1.0856 -0.5672 0.1584 -11.6710 -8.5242 1.4140 0.1407 -1.1966 -0.6349 -0.1218 -11.7617 0.0084 10.3163 -0.0083 0.4097 -0.0118 2.9077 -12.3794 0.2234 10.8350 0.0211 0.2714 0.0482 2.8422
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
%Nu(i) = xi - 1./2 - (pi^3 ./ Ra_values(i)) .* A(1) - (8 * pi^3 ./ Ra_values(i)) .* A(5);
%end
plot(Ra_values, Nu);
xlabel('Ra');
ylabel('Nu');
  1 个评论
John D'Errico
John D'Errico 2024-6-24
编辑：John D'Errico 2024-6-24
+1 of course. I very much appreciate seeing posts like yours. Well reasoned, good explanations and in good depth. I hope we continue to see you posting for a long time on the forum.
I would only suggest an idea that when you start with a random point at each pass for subtly different values of the parameter Ra, especially since this is a nonlinear system, you may be finding totally different solutions. It may be better to use the previous solution as a start point for fsolve at the next value of Ra.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Loops and Conditional Statements 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by