罚函数matlab代码

罚函数法是一种常用的非线性规划求解方法，其核心在于将原问题转化为一系列带约束的无约束问题。以下是一个简单的罚函数法的matlab代码实现，供参考：

function [x, fval, exitflag] = penalty(f, x0, constraints, options)
% f: 目标函数
% x0: 初始解
% constraints: 约束条件，形如 {c1, c2, ..., cm}，其中 ci(x) <= 0
% options: 选项，可以设置罚函数的初始值、放大因子、最大迭代次数等参数

% 默认选项
if nargin < 4
    options = struct();
end
if ~isfield(options, 'penalty')
    options.penalty = 1;
end
if ~isfield(options, 'factor')
    options.factor = 10;
end
if ~isfield(options, 'maxiter')
    options.maxiter = 100;
end
if ~isfield(options, 'tol')
    options.tol = 1e-6;
end

% 罚函数法求解
k = options.penalty; % 初始罚函数的值
x = x0;
fval = f(x);
for i = 1:options.maxiter
    % 构造罚函数
    penalty_f = @(x) f(x) + k * sum(max(constraints(x), 0).^2);
    % 无约束优化
    [x, fval, exitflag] = fminunc(penalty_f, x, options);
    % 判断是否满足约束条件
    if max(constraints(x)) <= 0
        break;
    end
    % 调整罚函数的值
    k = options.factor * k;
end
end

使用时，可以定义目标函数和约束条件，并调用 `penalty` 函数进行求解。例如：

% 目标函数
f = @(x) x(1)^2 + x(2)^2;
% 约束条件
c1 = @(x) x(1) + x(2) - 1;
c2 = @(x) x(1)^2 + x(2)^2 - 1;
constraints = {c1, c2};
% 初始解
x0 = [1; 0];
% 求解
options = struct('maxiter', 1000);
[x, fval, exitflag] = penalty(f, x0, constraints, options);
disp(['x = ', num2str(x')]);
disp(['fval = ', num2str(fval)]);

其中，约束条件以函数句柄的形式定义，并以一个 cell 数组的形式传入 `penalty` 函数。在选项中可以设置罚函数的初始值、放大因子、最大迭代次数等参数。