清空筛选条件
清空筛选条件

In solving two non-linear equations using fsolve, I am getting an error 'Failure in Initial objective function evaluation. FSOLVE cannot continue'. How to resolve?

3 次查看（过去 30 天）
Siva Rama Krishna Ayaluru
评论： Siva Rama Krishna Ayaluru 2023-5-11
%The code below solves two equations. The initial part of the code generates 10000 values of random data U,V,x_0,y_0 which are seed values for fsolve
clc
clear all
close all
n=10000;
%Generating 10000 values of Gaussian distributed random delay in the range
%5*10^-3 and 6*10^-3
a=5*10^-3;
b=6*10^-3;
X = (b-a).*randn(1,n) + a;
X1 = X;
d=10*10^-3;
%Generating T_1_i_A(:,1)
a=0.005;
b=5;
x = (b-a).*rand(1,n)+a;
T_1_i_A=zeros(n,30);
T_1_i_A(:,1)=transpose(x);
%Generating T_1_i_A for N=2:30 where N is no. of rounds of message exchange
for j=2:30
T_1_i_A(:,j)=T_1_i_A(:,j-1)+X(j);
end
%Generating T_2_i_B and T_2_i_P from T_1_i_A,d,X
T_2_i_B=zeros(n,30);
for j=1:30
T_2_i_B(:,j)=T_1_i_A(:,j)+20*10^-6+1.01*(T_1_i_A(:,j)-T_1_i_A(:,1))+15*10^-3+transpose(X);
end
T_2_i_P=zeros(n,30);
for j=1:30
T_2_i_P(:,j)=T_1_i_A(:,j)+30*10^-6+1.105*(T_1_i_A(:,j)-T_1_i_A(:,1))+5*10^-3+transpose(X1);
end
%finding U and V
U=T_2_i_P-T_2_i_B;
V=T_1_i_A-T_1_i_A(:,1);
fun = @root2d;
%Generating 10000 random values of pre-defined clock offset in the range
%5*10^-6 and 30*10^-6
a=5*10^-6;
b=30*10^-6;
n=10000;
x_0 = (b-a).*rand(1,n) + a;
%Generating 10000 random values of pre-defined clock skew in the range
%1.01 and 1.20
a1=1.01;
b1=1.20;
y_0=(b1-a1).*rand(1,n)+a1;
x0=zeros(2,n);
x=zeros(2,n);
for i=1:n
x0(:,i)=[x_0(i); y_0(i)];
end
for N=1:30
n=10000;
for i=1:n
x(:,i)= fsolve(@root2d,x0(:,i),[],U(i,N),V(i,N),N);
end
end
function f = root2d(x)
f(1)=(1.5./(U-x(1)-x(2)*V+10^-4)-0.05./(U-x(1)-x(2)*V+10^-4).^2)*N+0.1*sum(1:N);
f(2)=((1.5*V)./(U-x(1)-x(2)*+10^-4)-0.05*V./(U-x(1)-x(2)*V+10^-4).^2+0.1*V)*N;
return
end

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

回答（1 个）

Walter Roberson
Walter Roberson 2023-5-10
You really need to recheck your root2d code. You are passing scalar U and V entries but your for i loop involves U(:,j) which requires non-scalar U. But if U is non-scalar then your accumulated k is going to be non-scalar and then storing it into f(1) and f(2) would fail.
Your for j loop accumulating into k uses expressions that are independent of j except for the final 0.1*j -- so you might as well only calculate the expression once and then add on the cumulative total sum(1:N) * 0.1 .
Your for i loop is independent of i so it is going to give a result that is N times the value of the expression once.
Note that your use of U(:,j) inside for i is not indexing by i and is instead using whatever value of j happens to be sitting around in the workspace -- in particular it would be using U(:,N) since for j=1:N is going to leave j set to N after the loop.
%The code below solves two equations. The initial part of the code generates 10000 values of random data U,V,x_0,y_0 which are seed values for fsolve
n=10000;
%Generating 10000 values of Gaussian distributed random delay in the range
%5*10^-3 and 6*10^-3
a=5*10^-3;
b=6*10^-3;
X = (b-a).*randn(1,n) + a;
X1 = X;
d=10*10^-3;
%Generating T_1_i_A(:,1)
a=0.005;
b=5;
x = (b-a).*rand(1,n)+a;
T_1_i_A=zeros(n,30);
T_1_i_A(:,1)=transpose(x);
%Generating T_1_i_A for N=2:30 where N is no. of rounds of message exchange
for j=2:30
T_1_i_A(:,j)=T_1_i_A(:,j-1)+X(j);
end
%Generating T_2_i_B and T_2_i_P from T_1_i_A,d,X
T_2_i_B=zeros(n,30);
for j=1:30
T_2_i_B(:,j)=T_1_i_A(:,j)+20*10^-6+1.01*(T_1_i_A(:,j)-T_1_i_A(:,1))+15*10^-3+transpose(X);
end
T_2_i_P=zeros(n,30);
for j=1:30
T_2_i_P(:,j)=T_1_i_A(:,j)+30*10^-6+1.105*(T_1_i_A(:,j)-T_1_i_A(:,1))+5*10^-3+transpose(X1);
end
%finding U and V
U=T_2_i_P-T_2_i_B;
V=T_1_i_A-T_1_i_A(:,1);
fun = @root2d;
%Generating 10000 random values of pre-defined clock offset in the range
%5*10^-6 and 30*10^-6
a=5*10^-6;
b=30*10^-6;
n=10000;
x_0 = (b-a).*rand(1,n) + a;
%Generating 10000 random values of pre-defined clock skew in the range
%1.01 and 1.20
a1=1.01;
b1=1.20;
y_0=(b1-a1).*rand(1,n)+a1;
x0=zeros(2,n);
x=zeros(2,n);
for i=1:n
x0(:,i)=[x_0(i); y_0(i)];
end
options = optimoptions(@fsolve, 'display', 'none');
for N=1:30
n=10000;
for i=1:n
x(:,i,N)= fsolve(@(x)root2d(x, U(i,N),V(i,N),N), x0(:,i), options);
end
end
function f = root2d(x, U, V, N)
k=0;
for j=1:N
k=k+1.5./(U-x(1)-x(2)*V+10^-4)-0.05./(U-x(1)-x(2)*V+10^-4).^2+0.1*j;
end
f(1)=k;
k=0;
for i=1:N
k=k+(1.5*V)./(U-x(1)-x(2)*+10^-4)-0.05*V./(U-x(1)-x(2)*V+10^-4).^2+0.1*V;
end
f(2)=k;
return
end
  4 个评论
显示 2更早的评论
Walter Roberson
Walter Roberson 2023-5-10
Minor change to improve performance by pre-allocating the output array.
I ran the previous version without error.
%The code below solves two equations. The initial part of the code generates 10000 values of random data U,V,x_0,y_0 which are seed values for fsolve
n=10000;
%Generating 10000 values of Gaussian distributed random delay in the range
%5*10^-3 and 6*10^-3
a=5*10^-3;
b=6*10^-3;
X = (b-a).*randn(1,n) + a;
X1 = X;
d=10*10^-3;
%Generating T_1_i_A(:,1)
a=0.005;
b=5;
x = (b-a).*rand(1,n)+a;
T_1_i_A=zeros(n,30);
T_1_i_A(:,1)=transpose(x);
%Generating T_1_i_A for N=2:30 where N is no. of rounds of message exchange
for j=2:30
T_1_i_A(:,j)=T_1_i_A(:,j-1)+X(j);
end
%Generating T_2_i_B and T_2_i_P from T_1_i_A,d,X
T_2_i_B=zeros(n,30);
for j=1:30
T_2_i_B(:,j)=T_1_i_A(:,j)+20*10^-6+1.01*(T_1_i_A(:,j)-T_1_i_A(:,1))+15*10^-3+transpose(X);
end
T_2_i_P=zeros(n,30);
for j=1:30
T_2_i_P(:,j)=T_1_i_A(:,j)+30*10^-6+1.105*(T_1_i_A(:,j)-T_1_i_A(:,1))+5*10^-3+transpose(X1);
end
%finding U and V
U=T_2_i_P-T_2_i_B;
V=T_1_i_A-T_1_i_A(:,1);
fun = @root2d;
%Generating 10000 random values of pre-defined clock offset in the range
%5*10^-6 and 30*10^-6
a=5*10^-6;
b=30*10^-6;
n=10000;
x_0 = (b-a).*rand(1,n) + a;
%Generating 10000 random values of pre-defined clock skew in the range
%1.01 and 1.20
a1=1.01;
b1=1.20;
y_0=(b1-a1).*rand(1,n)+a1;
x0=zeros(2,n);
x=zeros(2,n);
for i=1:n
x0(:,i)=[x_0(i); y_0(i)];
end
options = optimoptions(@fsolve, 'display', 'none');
maxN = 30;
n=10000;
x = zeros(size(x0,1), n, maxN);
for N=1:maxN
for i=1:n
x(:,i,N)= fsolve(@(x)root2d(x, U(i,N),V(i,N),N), x0(:,i), options);
end
end
function f = root2d(x, U, V, N)
k=0;
for j=1:N
k=k+1.5./(U-x(1)-x(2)*V+10^-4)-0.05./(U-x(1)-x(2)*V+10^-4).^2+0.1*j;
end
f(1)=k;
k=0;
for i=1:N
k=k+(1.5*V)./(U-x(1)-x(2)*+10^-4)-0.05*V./(U-x(1)-x(2)*V+10^-4).^2+0.1*V;
end
f(2)=k;
return
end
Siva Rama Krishna Ayaluru
Siva Rama Krishna Ayaluru 2023-5-11
Thank you very much. I want to find the Mean square error of estimated clock offset and clock skew values across the pre-defined clock offset and skew values as shown in the code below. But the MSE is fluctuating, I mean not decreasing steadily. So I tried the same program with replacing fsolve by lsqnonlin at two places. But I am getting the error 'LSQNONLIN requires the following inputs to be of data type double: 'LB'.'
%The code below solves two equations. The initial part of the code generates 10000 values of random data U,V,x_0,y_0 which are seed values for fsolve
n=100;
%Generating 10000 values of Gaussian distributed random delay in the range
%5*10^-3 and 6*10^-3
a=5*10^-3;
b=6*10^-3;
X = (b-a).*randn(1,n) + a;
X1 = X;
d=10*10^-3;
%Generating T_1_i_A(:,1)
a=0.005;
b=5;
x = (b-a).*rand(1,n)+a;
T_1_i_A=zeros(n,30);
T_1_i_A(:,1)=transpose(x);
%Generating T_1_i_A for N=2:30 where N is no. of rounds of message exchange
for j=2:30
T_1_i_A(:,j)=T_1_i_A(:,j-1)+X(j);
end
%Generating T_2_i_B and T_2_i_P from T_1_i_A,d,X
T_2_i_B=zeros(n,30);
for j=1:30
T_2_i_B(:,j)=T_1_i_A(:,j)+20*10^-6+1.01*(T_1_i_A(:,j)-T_1_i_A(:,1))+15*10^-3+transpose(X);
end
T_2_i_P=zeros(n,30);
for j=1:30
T_2_i_P(:,j)=T_1_i_A(:,j)+30*10^-6+1.105*(T_1_i_A(:,j)-T_1_i_A(:,1))+5*10^-3+transpose(X1);
end
%finding U and V
U=T_2_i_P-T_2_i_B;
V=T_1_i_A-T_1_i_A(:,1);
fun = @root2d;
%Generating 10000 random values of pre-defined clock offset in the range
%5*10^-6 and 30*10^-6
a=5*10^-6;
b=30*10^-6;
n=100;
x_0 = (b-a).*rand(1,n) + a;
%Generating 10000 random values of pre-defined clock skew in the range
%1.01 and 1.20
a1=1.01;
b1=1.20;
y_0=(b1-a1).*rand(1,n)+a1;
x0=zeros(2,n);
x=zeros(2,n);
for i=1:n
x0(:,i)=[x_0(i); y_0(i)];
end
options = optimoptions(@fsolve, 'display', 'none');
maxN = 30;
n=100;
x = zeros(size(x0,1), n, maxN);
for N=1:maxN
for i=1:n
x(:,i,N)= fsolve(@(x)root2d(x, U(i,N),V(i,N),N), x0(:,i), options);
end
mse_theta_offset_hat_final(N)=mse(x(1,:,N),x_0);
mse_theta_skew_hat_final(N)=mse(x(2,:,N),y_0);
end
%Plots
k=1:30;
subplot(2,1,1),plot(k,mse_theta_offset_hat_final,'r+-')
title('clock offset')
xlabel('No. of Observations')
ylabel('MSE')
grid on
subplot(2,1,2),plot(k,mse_theta_skew_hat_final,'-*')
title('clock skew')
xlabel('No. of Observations')
ylabel('MSE')
grid on
function f = root2d(x, U, V, N)
k=0;
for j=1:N
k=k+1.5./(U-x(1)-x(2)*V+10^-4)-0.05./(U-x(1)-x(2)*V+10^-4).^2+0.1*j;
end
f(1)=k;
k=0;
for i=1:N
k=k+(1.5*V)./(U-x(1)-x(2)*+10^-4)-0.05*V./(U-x(1)-x(2)*V+10^-4).^2+0.1*V;
end
f(2)=k;
return
end

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Sources and Sinks 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by