0

我想找到使函数最大化的向量p（即）的值。 p(1),p(2),p(3),...A(p)

我正在使用 MATLAB 来做到这一点，我发现fsolve我认为它可以帮助我。所以我做了功能A：

function A = myfun(p)

    R  = 0.1; 
    u1 = 500;
    u2 = 400;
    u3 = 300;

    A = ( (p(1)+p(2)+p(3)) * (1/u1+1/u2+1/u3)) * ...
        (1 + R*(p(1)^2+p(2)^2+p(3)^2) * (1/u1+1/u2+1/u3) );

然后我需要求解一个方程组，它将是：

diff(A,p(1))==0

diff(A,p(2))==0

diff(A,p(3))==0

结果p向量将解决我的问题。

如何fsolve求解这个方程组 ( p0=[1 1 1]) ？

4

1 回答 1

0

以下是如何使用符号工具箱执行此操作的示例：

%// (there's probably a way to generalize this as well)
dAdp = cellstr(char(...
    [char(diff('C * (p1+p2+p3) * (1 + R*C*(p1^2+p2^2+p3^2))', 'p1')) ';'],...
    [char(diff('C * (p1+p2+p3) * (1 + R*C*(p1^2+p2^2+p3^2))', 'p2')) ';'],...
    [char(diff('C * (p1+p2+p3) * (1 + R*C*(p1^2+p2^2+p3^2))', 'p3')) ';']))

%// convert to proper vector equation
dAdp = regexprep(dAdp, 'p([0-9])', 'p\($1\)');

%// convert to function handle
F = str2func( strcat('@(p) [', dAdp{:}, ']') );

但我不建议这样做（它完全不可读且容易出错）。

您可以使用符号工具箱编写一个适当的函数并dAdp在上面进行评估，但我也不建议这样做（它非常慢）。

我是一个数字人，因为我只相信计算机的用途：计算。我在纸上做的推导，除了过于简单但又冗长乏味的推导（你可以说是这样，我仍然喜欢在纸上做，因为我也喜欢锻炼我的大脑:)。

我建议您这样做，和/或使用符号工具箱不断检查自己。恕我直言，您应该将其用作辅助工具，而不是作为主要引擎。

所以，我们开始吧。这又是你的函数，这次以不同的形式，应用了更多的大脑：

function A = myfun(p)    
    R = 0.1; 
    u = [500; 400; 300];     
    C = sum(1./u);
    B = C * sum(p) * (1 + R*C*sum(p.^2));

所以，你想解决 ∇A( p ) = 0。最好的办法是多动脑筋。您可以验证向量导数是否等于：

function F = Aprime(p)
    R = 0.1;
    u = [500; 400; 300];
    C = sum(1./u);
    F = C*( C*R*sum(p.^2) + 2*C*R*p*sum(p) + 1 );

你可以用更多的大脑来解决：

C·( CR·Σp² + 2CR·p·Σp + 1 ) = 0
                Σp² + 2p·Σp = -1/(CR)

矢量 == 标量：这意味着 p 中的所有元素都是相等的。
代入 q = p1 = p2 = p3，则

 3q² + 2q·3q = -1/(CR)
         9q² = -1/(CR)

⇒ q = ⅓·√(-1/(CR))

（表示麻烦），或者fsolve像这样：

fsolve(@Aprime, [1 1 1])

这将立即导致

fsolve 已完成，因为在初始点处的函数值向量按函数容差的默认值测量接近于零，并且根据梯度测量，问题看起来很规则。

这确实表明麻烦。

由于您已表示您对最大值感兴趣，并且该函数似乎没有固定点，因此唯一的结论是该函数没有最大值也没有最小值，除非您对p施加界限。事实上，如果你将维度减少到 2，并绘制一个图：

R = 0.1;
u = [500; 400; 300];
C = sum(1./u);
B = @(p) C * sum(p) .* (1 + R*C*sum(p.^2));

[p1,p2] = meshgrid(-10:0.1:10);
surf(p1,p2,reshape(B([p1(:) p2(:)].'), size(p1)), 'edgecolor', 'none')

在此处输入图像描述

...没有极值。

于 2014-03-20T16:28:02.713 回答