matlab - Matlab中非线性方程组的一组解

Question

我在 Matlab 中有一个非线性方程组，我想找到它的解集。该系统如下所示。

想象一下，Matlab 函数mvncdf可以与solve. 我的方程组看起来像这样

clear
%SOME VALUES USED BELOW
U1true=1;
U2true=-1;
U3true=-2;
U4true=1;
rhotrue=-0.5;
mutrue=[0 0];
sigmatrue=[2*(1-rhotrue) 1-rhotrue; 1-rhotrue 2*(1-rhotrue)];

%SET VARIABLES
syms U1 U2 U3 U4 rho 

%SOLVE SYSTEM WITH 6 EQUATIONS
S=solve(mvncdf([ U1  -U2+U1], [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([ U1true  -U2true+U1true], mutrue, sigmatrue) ,...
        mvncdf([ U2   U2-U1], [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([ U2true   U2true-U1true], mutrue, sigmatrue),...
        mvncdf([-U1  -U2],    [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([-U1true  -U2true],        mutrue, sigmatrue),...
        mvncdf([ U3  -U4+U3], [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([ U3true  -U4true+U3true], mutrue, sigmatrue) ,...
        mvncdf([ U4   U4-U3], [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([ U4true   U4true-U3true], mutrue, sigmatrue),...
        mvncdf([-U3  -U4],    [0 0], [2*(1-rho) 1-rho; 1-rho 2*(1-rho)])==mvncdf([-U3true  -U4true],        mutrue, sigmatrue));

%DISPLAY SOLUTIONS    
M=[S.U1,S.U2,S.U3,S.U4, S.rho]; 

由于许多原因，这种方法显然不起作用，最重要的是因为mvncdf不能与 一起使用solve，如这里讨论的例子。

作为解决我的系统的替代方法，我考虑使用slicesample如下

这段代码实现了我的想法。

%MAIN
clear
rng default 

Utrue=[1;-1;-2;1];
rhotrue=-0.5;

number = 100000; 
initials=[Utrue;rhotrue]; %[U1;U2;U3;U4;rho]
rnd = slicesample(initials,number,'pdf',@dens);

%ANCILLARY FUNCTIONS TO SAVE
function density=dens(param)
         T=10^(-6); 

         if param(end)>=1 || param(end)<=-1 %rho outside (-1,1)
            density=0; %force density to be equal to zero
         else

            Utrue=[1;-1;-2;1];
            rhotrue=-0.5;

            mutrue=[0 0];
            sigmatrue=[2*(1-rhotrue) 1-rhotrue; 1-rhotrue 2*(1-rhotrue)];

            truth=mvncdf([ Utrue(1)  -Utrue(2)+Utrue(1);...
                        Utrue(2)   Utrue(2)-Utrue(1);...
                       -Utrue(1)  -Utrue(2)     ;...
                        Utrue(3)  -Utrue(4)+Utrue(3);...
                        Utrue(4)   Utrue(4)-Utrue(3);...
                       -Utrue(3)  -Utrue(4)], mutrue, sigmatrue);

           density=exp(-criterion(param,truth)/T); 
         end
end

function Q=criterion(param, truth)
         U=param(1:end-1);
         rho=param(end);
         mu=[0,0];
         sigma=[2*(1-rho) 1-rho; 1-rho 2*(1-rho)];
         candidates=mvncdf([ U(1)  -U(2)+U(1);...
                             U(2)   U(2)-U(1);...
                            -U(1)  -U(2)     ;...
                             U(3)  -U(4)+U(3);...
                             U(4)   U(4)-U(3);...
                            -U(3)  -U(4)], mu, sigma);

         Q=(candidates(1)-truth(1))^2+(candidates(2)-truth(2))^2+(candidates(3)-truth(3))^2+...
           (candidates(4)-truth(4))^2+(candidates(5)-truth(5))^2+(candidates(6)-truth(6))^2; %this is zero at a solution of the system

end

问题：使用我提出的方法，随着变量的增加，number我得到了越来越多的不同值的解rho，这与上面的评论一致，即系统对于每个值都有一个且只有一个关于 U 的解rho。但是，在我看来，要跨越整个(-1,1)区间rho，我需要将变量设置得number非常高。您是否有（i）可以帮助我改进上述代码以解决此问题的建议，或者（ii）更好的方法来提供我的方程组的解集？

标签： matlab