如何使用匿名函数在MATLAB中优化约束整数表达式? [英] How do I optimize constrained integral expressions in MATLAB using anonymous functions?

问题描述

我有一个综合错误表达式 E = int [abs(x- p)^ 2] dx ，范围为 x | 0 到 x | L .变量p是形式为 2 *(a * sin(x)+ b(a)* sin(2 * x)+ c(a)* sin(3 * x))的多项式.换句话说，系数 b 和 c 都是a的已知表达式.附加公式为 dE/da = 0 .如果定义了上限L，则方程组将关闭，我可以求解a，给出三个系数.

I have an integrated error expression E = int[ abs(x-p)^2 ]dx with limits x|0 to x|L. The variable p is a polynomial of the form 2*(a*sin(x)+b(a)*sin(2*x)+c(a)*sin(3*x)). In other words, both coefficients b and c are known expressions of a. An additional equation is given as dE/da = 0. If the upper limit L is defined, the system of equations is closed and I can solve for a, giving the three coefficients.

我设法获得了一个纯粹基于最大化L的优化例程来解决a.通过在下面的代码中设置optimize=0可以确认这一点.它提供的解决方案与我以解析方式解决问题的方式相同.因此，我知道求解系数a的方程式是正确的.

I managed to get an optimization routine to solve for a purely based on maximizing L. This is confirmed by setting optimize=0 in the code below. It gives the same solution as if I solved the problem analytically. Therefore, I know the equations to solve for the coefficent a are correct.

我知道我提出的示例可以用铅笔和纸解决，但我正在尝试构建针对此类问题的通用化优化功能(我有很多要评估的功能).理想情况下，polynomial作为函数的输入参数给出，然后输出xsol.显然，在担心泛化之前，我需要对本文介绍的polynomial进行优化.

I know the example I presented can be solved with pencil and paper, but I'm trying to build an optimization function that is generalized for this type of problem (I have a lot to evaluate). Ideally, polynomial is given as an input argument to a function which then outputs xsol. Obviously, I need to get the optimization to work for the polynomial I presented here before I can worry about generalizations.

无论如何，我现在需要在一些约束下进一步优化该问题.首先，选择L.这使我可以计算a.知道a后，polynomial仅是x的已知功能，即p(x).然后，我需要根据满足以下约束的0->x确定最大的INTERVAL:|dp(x)/dx - 1| < tol.这为我提供了系数为a的polynomial性能的度量.间隔就是我所说的带宽".我想强调两件事:1)带宽"与L不同. 2)带宽"内所有x的值都必须满足约束.函数dp(x)/dx确实会在公差标准内波动，因此无法对x的单个值进行测试.它必须经过一段时间的测试.违反的第一个实例定义了带宽.我需要最大化此带宽"/间隔.对于输出，我还需要知道哪个L会导致这种优化，因此我知道要为给定的约束选择正确的a.那是正式的问题陈述. (我希望这次我做对了)

Anyway, I now need to further optimize the problem with some constraints. To start, L is chosen. This allows me to calculate a. Once a is know, the polynomial is a known function of x only i.e p(x). I need to then determine the largest INTERVAL from 0->x over which the following constraint is satisfied: |dp(x)/dx - 1| < tol. This gives me a measure of the performance of the polynomial with the coefficient a. The interval is what I call the "bandwidth". I would like to emphasis two things: 1) The "bandwidth" is NOT the same as L. 2) All values of x within the "bandwidth" must meet the constraint. The function dp(x)/dx does oscillate in and out of the tolerance criteria, so testing the criteria for a single value of x does not work. It must be tested over an interval. The first instance of violation defines the bandwidth. I need to maximize this "bandwidth"/interval. For output, I also need to know which L lead to such an optimization, hence I know the correct a to choose for the given constraints. That is the formal problem statement. (I hope I got it right this time)

现在，我的问题是使用MATLAB的优化工具来设置整个过程.我尝试遵循以下文章的想法:

Now my problem is setting this whole thing up with MATLAB's optimization tools. I tried to follow ideas from the following articles:

• Optimization Toolbox™的教程

将if statement设置为optimize=1可以使用约束优化.我想过如何进行嵌套优化，但是我什么也做不了.我从IMSL优化库中提供了已知问题的解决方案以进行比较/检查.它们写在优化例程下面.无论如何，这是到目前为止我编写的代码:

Setting optimize=1 for the if statement will work with the constrained optimization. I thought some how nested optimization is involved, but I couldn't get anything to work. I provided known solutions to the problem from the IMSL optimization library to compare/check with. They are written below the optimization routine. Anyway, here is the code I've put together so far:

function [history] = testing()

% History
history.fval = [];
history.x = [];
history.a = []; 

%----------------
% Equations
polynomial = @(x,a) 2*sin(x)*a + 2*sin(2*x)*(9/20 -(4*a)/5) + 2*sin(3*x)*(a/5 - 2/15);
dpdx = @(x,a) 2*cos(x)*a + 4*cos(2*x)*(9/20 -(4*a)/5) + 6*cos(3*x)*(a/5 - 2/15);

% Upper limit of integration
IC = 0.8;       % initial
LB = 0;         % lower
UB = pi/2;      % upper

% Optimization
tol = 0.003;


% Coefficient
% --------------------------------------------------------------------------------------------
dpda = @(x,a) 2*sin(x) + 2*sin(2*x)*(-4/5) + 2*sin(3*x)*1/5;
dEda = @(L,a) -2*integral(@(x) (x-polynomial(x,a)).*dpda(x,a),0,L);
a_of_L = @(L) fzero(@(a)dEda(L,a),0);                                   % Calculate the value of "a" for a given "L"
EXITFLAG = @(L) get_outputs(@()a_of_L(L),3);                            % Be sure a zero is actually calculated


% NL Constraints
% --------------------------------------------------------------------------------------------
% Equality constraint (No inequality constraints for parent optimization)
ceq = @(L) EXITFLAG(L) - 1; % Just make sure fzero finds unique solution 
confun = @(L) deal([],ceq(L));

% Objective function
% --------------------------------------------------------------------------------------------
% (Set optimize=0 to test coefficent equations and proper maximization of L )
optimize = 1;

if optimize

%%%%  Plug in solution below

else 
    % Optimization options
    options = optimset('Algorithm','interior-point','Display','iter','MaxIter',500,'OutputFcn',@outfun);

    % Optimize objective
    objective = @(L) -L;
    xsol = fmincon(objective,IC,[],[],[],[],LB,UB,confun,options);

    % Known optimized solution from IMSL library
    % a = 0.799266;
    % lim = pi/2;
    disp(['IMSL coeff (a): 0.799266     Upper bound (L): ',num2str(pi/2)])
    disp(['code coeff (a): ',num2str(history.a(end)),'   Upper bound: ',num2str(xsol)])

end




    % http://stackoverflow.com/questions/7921133/anonymous-functions-calling-functions-with-multiple-output-forms
    function varargout = get_outputs(fn, ixsOutputs)
        output_cell = cell(1,max(ixsOutputs));
        [output_cell{:}] = (fn());
        varargout = output_cell(ixsOutputs);
    end

    function stop = outfun(x,optimValues,state)
        stop = false;

        switch state
            case 'init'
            case 'iter'
                % Concatenate current point and objective function
                % value with history. x must be a row vector.
                history.fval = [history.fval; optimValues.fval];
                history.x = [history.x; x(1)];
                history.a = [history.a; a_of_L(x(1))];

            case 'done'
            otherwise
        end
    end
end

我真的可以使用一些帮助来设置约束优化.我不仅是优化的新手，而且从未使用过MATLAB.我还应该注意，我上面的内容不起作用，并且对于约束优化来说是不正确的.

I could really use some help setting up the constrained optimization. I'm not only new to optimizations, I've never used MATLAB to do so. I should also note that what I have above does not work and is incorrect for the constrained optimization.

更新:我在if optimize部分中添加了一个for循环，以显示我要通过优化实现的目标.显然，我可以使用它，但是它看起来效率很低，尤其是当我提高range的分辨率并不得不多次运行此优化时.如果取消注释图，它将显示bandwidth的行为.通过遍历整个范围，基本上可以测试每个L，但是肯定有一种更有效的方法可以做到这一点?

UPDATE: I added a for loop in the section if optimizeto show what I'm trying to achieve with the optimization. Obviously, I could just use this, but it seems very inefficient, especially if I increase the resolution of range and have to run this optimization many times. If you uncomment the plots, it will show how the bandwidth behaves. By looping over the full range, I'm basically testing every L but surely there's got to be a more efficient way to do this??

更新:已解决

推荐答案

因此，看来fmincon并不是从事这项工作的唯一工具.实际上，我什至无法运行它.在下面，fmincon被卡住"在IC上，并且拒绝执行任何操作……为什么……这是另一篇文章！使用相同的布局和公式，fminbnd找到正确的解决方案.据我所知，唯一的区别是前者使用的是有条件的.但是我的条件是没有花哨，而且真的是不需要的.因此，它必须与算法有关.我想这就是使用黑匣子"时得到的.无论如何，经过长期的，痛苦的，痛苦的学习经历，这是一个解决方案:

So it seems fmincon is not the only tool for this job. In fact I couldn't even get it to work. Below, fmincon gets "stuck" on the IC and refuses to do anything...why...that's for a different post! Using the same layout and formulation, fminbnd finds the correct solution. The only difference, as far as I know, is that the former was using a conditional. But my conditional is nothing fancy, and really unneeded. So it's got to have something to do with the algorithm. I guess that's what you get when using a "black box". Anyway, after a long, drawn out, painful, learning experience, here is a solution:

options = optimset('Display','iter','MaxIter',500,'OutputFcn',@outfun);

% Conditional
index = @(L) min(find(abs([dpdx(range(range<=L),a_of_L(L)),inf] - 1) - tol > 0,1,'first'),length(range));

% Optimize
%xsol = fmincon(@(L) -range(index(L)),IC,[],[],[],[],LB,UB,confun,options);
xsol = fminbnd(@(L) -range(index(L)),LB,UB,options);

我要特别感谢@AndrasDeak的所有支持.如果没有帮助，我将无法解决！

I would like to especially thank @AndrasDeak for all their support. I wouldn't have figured it out without the assistance!

这篇关于如何使用匿名函数在MATLAB中优化约束整数表达式?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！