www.gusucode.com > mbcdata 工具箱 matlab 源码程序 > mbcdata/mbcOSmodalopt.m
function out = mbcOSmodalopt(action, in) %MBCOSMODALOPT CAGE Optimization library function % MBCOSMODALOPT contains the configuration and evaluation routines for % optimal selection of operating modes in CAGE. This function is called % by CAGE using the syntaxes below. % % OPTOBJ = MBCOSMODALOPT('OPTIONS', OPTOBJ) is called by CAGE when an % Mixed optimization is created in CAGE. This call defines the % structure of the optimization problem, using a CGOPTIMOPTIONS object. % The completed configuration for the Optimization is returned to CAGE % via OPTOBJ. % % OPTIMSTORE = MBCOSMODALOPT('EVALUATE', OPTIMSTORE) is called by CAGE % when a user runs an FOPTCON optimization. Optimization parameters are % passed to the evaluation routine through OPTIMSTORE. The results of the % optimization are stored in the OPTIMSTORE object which is returned to % CAGE for further processing. % % See also MBCOSFMINCON % Copyright 2009-2015 The MathWorks, Inc. and Ford Global Technologies, Inc. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Deal with the two action cases %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% switch lower(action) case 'options' % Return the updated options object out = i_Options(in); case 'evaluate' % Return a handle to the main evaluation routine out = i_Evaluate(in); end end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % options Function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function opt = i_Options(opt) % Set the evaluation function interface to always be the latest one opt = setRunInterfaceVersion(opt, 0); % Add a name opt = setName(opt, 'Modal optimization'); % Add a description opt = setDescription(opt, 'Modal optimization using a single objective optimization subject to constraints'); % Set up the free variables %%%% Any number are allowed %%%% % Set up the objective functions %%%% Any objective is allowed %%%% opt = setObjectivesMode(opt,'any'); % Set up the constraints %%%% Any number of constraints are allowed %%%% % Allow objectives and constraints to be renamed opt = setRenameMode(opt, true); % allow equality constraints opt = setEqConMode(opt,'on'); % Set up the operating point sets %%%% No operating point sets are allowed %%%% opt = setOperatingPointsMode(opt, 'fixed'); % Set up the optimization parameters opt = addParameter(opt, 'Display', ... {'list', {'none', 'final', 'iter'}}, 'none'); opt = addParameter(opt, 'MaxIter', ... {'integer', 'positive'}, 500, ... 'Maximum iterations'); opt = addParameter(opt, 'MaxFunEvals', ... {'integer', 'positive'}, 2000, ... 'Maximum function evaluations'); opt = addParameter(opt, 'TolX', ... {'number', 'positive'}, 1e-6, ... 'Variable tolerance'); opt = addParameter(opt, 'TolFun',... {'number', 'positive'}, 1e-6, ... 'Function tolerance'); opt = addParameter(opt, 'TolCon', ... {'number', 'positive'}, 1e-6, ... 'Constraint tolerance'); opt = addParameter(opt, 'Algorithm', ... {'list', {'active-set','sqp', 'interior-point'}}, 'interior-point',... 'Constrained optimization algorithm'); opt = addParameter(opt, 'DiffMinChange', ... {'number', 'positive'}, 1e-8, ... 'Minimum change in variables for gradient'); opt = addParameter(opt, 'DiffMaxChange', ... {'number', 'positive'}, 0.1, ... 'Maximum change in variables for gradient'); opt = addParameter(opt, 'ModeVariable',... {'integer', 'positive'},1,... 'Index to the mode free variable'); opt = addParameter(opt, 'ModeObjective',... {'integer', 'positive'},1,... 'Index to the objective to determine best mode'); end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Evaluation Function %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function optimstore = i_Evaluate(optimstore) % Optimization routine must conform to the following API % optimstore = <OPTIMIZATION_FUNCTION>(optimstore) % % Inputs: optimstore: store of all the CAGE data needed % Outputs: optimstore: updated object containing all the CAGE data, including outputs % % Get bounds and linear constraints A = getA(optimstore); B = getB(optimstore); lb = getLB(optimstore); ub = getUB(optimstore); ModeVariable = getParam(optimstore, 'ModeVariable'); try ModeObjective = getParam(optimstore, 'ModeObjective'); catch % old versions don't have ModeObjective ModeObjective = 1; end % Get the user defined options and put them in an optimization toolbox % options structure for fmincon options = optimset(optimstore); doDisplay = ~strcmp(options.Display,'none') && ~getStopState(optimstore); x0 = getInitFreeVal(optimstore); if isScalarFreeVariables(optimstore) % modal optimization only works for scalar problems % Set some extra options in fmincon options = optimset(options, ... 'OutputFcn',@i_Output,... 'GradObj', 'off', ... 'GradCon', 'off'); if ModeVariable>numel(x0) error(message('mbc:mbcOSmodalopt:InvalidValue')) end NumModes = ub(ModeVariable); MinMode = lb(ModeVariable); if NumModes~=fix(NumModes) || NumModes<1 || ~isfinite(NumModes) || ... MinMode~=fix(MinMode) || MinMode~=1 || ~isfinite(MinMode) error(message('mbc:mbcOSmodalopt:InvalidValue1')) end IndexCont = (1:length(x0))~=ModeVariable; % Initialise return quantities bestx = repmat(x0,NumModes,1); bestfun = Inf(NumModes, 1); exitFlag = zeros(1,NumModes); % make output structure OUTPUT(NumModes).iterations = []; OUTPUT(NumModes).funcCount = []; OUTPUT(NumModes).stepsize = []; OUTPUT(NumModes).algorithm = []; OUTPUT(NumModes).firstorderopt = []; OUTPUT(NumModes).message = []; % remove Mode Variable from bounds and linear constraints LB = lb(IndexCont); UB = ub(IndexCont); if ~isempty(A) A = A(IndexCont,:); B = B(IndexCont); end Xfull = x0; runNo = getRunIndex(optimstore); if doDisplay fprintf('Run %d\n',runNo); end Objectives = getObjectives(optimstore); ObjectiveFuncTypes = getObjectiveType(optimstore); if ModeObjective>length(Objectives) error(message('mbc:mbcOSmodalopt:InvalidValue2')) end for Mode = 1:NumModes Xfull(ModeVariable) = Mode; Xmode0 = x0(:, IndexCont); if getStopState(optimstore) y0= NaN; else y0 = i_evalObj(Xmode0,Xfull,IndexCont,optimstore); end bestx(Mode, ModeVariable) = Mode; if isfinite(y0) % optimize for each mode [bestx(Mode, IndexCont), bestfun(Mode), exitFlag(Mode), outstruct] = fmincon(@i_evalObj, ... Xmode0, A, B, [],[], LB, UB, @i_evalCon, options,Xfull,IndexCont,optimstore); if ModeObjective>1 % evaluate a different objective as mode objective bestfun(Mode) = evaluateObjective(optimstore, bestx(Mode,:),Objectives(ModeObjective)); if strcmp(ObjectiveFuncTypes{ModeObjective},'max') bestfun(Mode) = -bestfun(Mode); end end if doDisplay if exitFlag(Mode)<0 fprintf(' Mode %d: Not optimal\n', Mode); elseif exitFlag(Mode)==0 fprintf(' Mode %d: Incomplete %g\n', Mode,bestfun(Mode)); else fprintf(' Mode %d: Optimal %g\n', Mode,bestfun(Mode)); end end OUTPUT(Mode).iterations = outstruct.iterations; OUTPUT(Mode).funcCount = outstruct.funcCount; OUTPUT(Mode).stepsize = outstruct.stepsize; OUTPUT(Mode).algorithm = outstruct.algorithm; OUTPUT(Mode).firstorderopt = outstruct.firstorderopt; OUTPUT(Mode).message = sprintf('Mode %d: %s',Mode,iGetShortMessage(outstruct.message)); else bestx(Mode, IndexCont) = NaN; bestfun(Mode) = Inf; if getStopState(optimstore) % user stopped optimization OUTPUT(Mode).message = sprintf('Mode %d: Modal optimization stopped by the output or plot function.',Mode); exitFlag(Mode) = -1; else % most likely the optimization didn't run due to an invalid % operating point for the switch model if doDisplay fprintf(' Mode %d: Invalid\n', Mode); end OUTPUT(Mode).message = sprintf('Mode %d: Invalid\n', Mode); exitFlag(Mode) = -7; end OUTPUT(Mode).iterations = 0; OUTPUT(Mode).funcCount = 0; OUTPUT(Mode).stepsize = 0; OUTPUT(Mode).algorithm = 'fmincon'; OUTPUT(Mode).firstorderopt = []; OUTPUT(Mode).cgiterations = []; end end % Choose the best run and get the results at that run if any(exitFlag>0) % only use successful optim modes IndOK = find(exitFlag>0); [~, bestInd] = min(bestfun(IndOK)); SelectedSolution = IndOK(bestInd); elseif any(exitFlag>-7) % use all modes IndOK = find(exitFlag>-7); [~, bestInd] = min(bestfun(IndOK)); SelectedSolution = IndOK(bestInd); else [~, SelectedSolution] = min(bestfun); end if doDisplay fprintf('Best Mode %d: %g \n\n', SelectedSolution,bestfun(SelectedSolution)); end else % Not required to run the optimization. Return the initial conditions % as the "solution". bestx = getInitFreeVal(optimstore); exitFlag = -8; OUTPUT.iterations = 0; OUTPUT.funcCount = 0; OUTPUT.stepsize = 0; OUTPUT.algorithm = 'fmincon'; OUTPUT.firstorderopt = []; OUTPUT.cgiterations = []; OUTPUT.message = 'Modal optimization cannot be used for nonscalar problems. Algorithm not run.'; SelectedSolution = []; end % Save the results optimstore = setFreeVariables(optimstore, bestx,SelectedSolution); optimstore = setExitStatus(optimstore, exitFlag, {OUTPUT.message}); end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Options plus objective and constraint evaluation % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function [y, yg] = i_evalObj(x,Xfull,IndexCont,optimstore) Xfull(IndexCont) = x(:); yg = []; Objectives = getObjectives(optimstore); y = evaluateObjective(optimstore, Xfull,Objectives(1)); ObjectiveFuncTypes = getObjectiveType(optimstore); switch ObjectiveFuncTypes{1} case 'max' % Function and gradient need to be negated y = -y; case 'min' otherwise error(message('mbc:mbcOSfmincon:InvalidState')); end end function [c, ceq, cg, cgeq] = i_evalCon(x,Xfull,IndexCont,optimstore) Xfull(IndexCont) = x(:); if nargout>2 % evaluate nonlinear inequality constraints [c, cg] = evaluateIneqCon(optimstore, Xfull); % evaluate nonlinear equality constraints [ceq, cgeq] = evaluateEqCon(optimstore, Xfull); else c = evaluateIneqCon(optimstore, Xfull); ceq = evaluateEqCon(optimstore, Xfull); end end function stop= i_Output(varargin) % get stop state so the optimization can stop in the middle of run stop= getStopState(varargin{end}); end function shortstr = iGetShortMessage(str) % Process output message from fmincon and extract the first short sentence. shortstr = regexp(str, '\n*(.*?)\n\n', 'once', 'tokens'); if ~isempty(shortstr) shortstr = shortstr{1}; else % Regexp did not match anything. Forward the string as-is. shortstr = str; end end