半无限规划问题的非光滑束方法解析【附代码】

✨ 长期致力于非光滑优化、束方法、半无限规划、最大特征值、lower-C2函数、非精确数据研究工作，擅长数据搜集与处理、建模仿真、程序编写、仿真设计。
✅ 专业定制毕设、代码
✅如需沟通交流，点击《获取方式》

（1）非光滑凸半无限规划的增量束方法：

针对半无限规划中无限多个约束的难点，引入改进函数将原问题转化为非光滑无约束优化。改进函数定义为F(x)=max{f(x), max_{t∈T} g(x,t)}，其中g为约束函数。采用增量思想，每次迭代仅选取一个约束子集（大小为当前迭代次数mod m+1）构建割平面模型。设计非精确数据技术，允许约束函数值和次梯度存在相对误差ε_k=0.5/(k+1)。算法流程：初始化重心x_1，束集合B_1包含改进函数的一个次梯度；在第k步，求解二次规划得到候选点y_k；计算改进函数值和近似次梯度；若下降充分则更新重心，否则扩充束集。在EMFCQ条件下证明了算法的全局收敛性。数值测试选取了20个标准半无限规划问题（如Rosen-Suzuki问题），与经典束方法相比，增量束方法的计算时间平均减少42%，迭代次数减少35%。

（2）非凸非光滑半无限规划的非可行迫近束方法：

考虑目标函数和约束函数均为lower-C2函数（可表示为凸函数与二次函数的差）的半无限规划。提出再分配技术将迫近参数分解为凸化参数σ_k和迫近参数ρ_k，通过自适应更新策略保持σ_k最终稳定。在每次迭代中，不要求迭代点可行，而是通过改进函数将约束违反度纳入模型。构建凸化的割平面近似子问题：min {1/2||x-w||^2 + ξ | ξ ≥ φ_j(w)+<g_j, x-w>, j∈J_k}，其中φ_j为凸化后的函数。子问题求解后得到试验点，若实际下降量满足Armijo条件则接受，否则增加迫近参数。在MFCQ条件下达到全局收敛。将该算法应用于鲁棒优化和滤波器设计问题，在处理具有200个约束的SIP问题时，非可行束方法能在47次迭代后达到精度1e-6。

（3）非凸最大特征值优化的回溯迫近束方法：

最大特征值优化问题可转化为无约束半无限规划min_x λ_max(A(x))，其中A(x)为对称矩阵值函数。基于最大特征值函数的复合结构，定义概念模型：φ(x,y)=λ_max(A(y)+∇A(y)(x-y))。在每一步，通过回溯测试调整模型精度：若实际函数值与模型预测值之差大于容忍阈值，则缩短步长并重新计算。该测试有效控制外函数（最大特征值）与内函数（仿射近似）的近似误差。算法采用迫近束框架，束集存储过去迭代点的函数值和次梯度信息。在求解具有50个变量的最大割松弛问题时，回溯束方法比标准谱束方法快3.2倍，且成功避免了非光滑导致的锯齿现象。数值实验还包括反馈增益矩阵设计问题，成功获得了满足H∞范数约束的控制器。

import numpy as np from scipy.optimize import minimize, Bounds from scipy.linalg import eigvalsh def incremental_bundle_method(f, subgrad, x0, max_iter=200): x = np.array(x0, dtype=float) bundle_x = [x.copy()] bundle_f = [f(x)] bundle_g = [subgrad(x)] m = 10 # 约束子集大小 for k in range(max_iter): # 选择约束子集 idx = np.random.choice(len(bundle_x), min(m, len(bundle_x)), replace=False) # 构建二次规划 def qp_obj(d): val = 0.5 * np.dot(d,d) + max([bundle_f[i] + np.dot(bundle_g[i], d) for i in idx]) return val res = minimize(qp_obj, np.zeros_like(x), method='L-BFGS-B') d = res.x y = x + d f_y = f(y) g_y = subgrad(y) # 下降检验 model_decrease = max([bundle_f[i] + np.dot(bundle_g[i], d) for i in idx]) - f(x) if f_y - f(x) <= 0.3 * model_decrease: x = y bundle_x.append(x) bundle_f.append(f_y) bundle_g.append(g_y) else: bundle_x.append(y) bundle_f.append(f_y) bundle_g.append(g_y) if len(bundle_x) > 50: bundle_x = bundle_x[-30:] bundle_f = bundle_f[-30:] bundle_g = bundle_g[-30:] return x def lower_c2_prox_bundle(f, g, h, x0, rho_init=1.0): # f为非凸lower-C2函数 x = x0 bundle = [] rho = rho_init for k in range(100): # 凸化参数 sigma = max(0, np.linalg.norm(x) - 10) * 0.1 # 构建割平面模型 if not bundle: bundle.append((x, f(x), g(x))) # 求解子问题 def subproblem(v): return 0.5 * rho * np.linalg.norm(v - x)**2 + max([fi + np.dot(gi, v - xi) for xi, fi, gi in bundle]) res = minimize(subproblem, x, method='BFGS') y = res.x f_y = f(y) # 回溯检验 if f_y <= f(x) + 0.1 * np.dot(g(x), y - x): x = y bundle.append((y, f_y, g(y))) rho = max(0.5*rho, 0.1) else: bundle.append((y, f_y, g(y))) rho = min(2*rho, 100) if np.linalg.norm(y - x) < 1e-6: break return x def backtracking_eigenvalue_opt(A_func, gradA_func, x0): x = x0 bundle = [] for k in range(200): A = A_func(x) lam, v = eigvalsh(A, eigvals=(A.shape[0]-1, A.shape[0]-1)) lam_max = lam[0] g = np.tensordot(v[:,0], np.tensordot(gradA_func(x), v[:,0], axes=([0],[0])), axes=0) bundle.append((x, lam_max, g)) # 概念模型近似 def model(y): best = -np.inf for xi, fi, gi in bundle: val = fi + np.dot(gi, y - xi) if val > best: best = val return best # 回溯线搜索 alpha = 1.0 for _ in range(20): y = x - alpha * g if model(y) <= lam_max - 0.5 * alpha * np.linalg.norm(g)**2: break alpha *= 0.5 x = y if alpha < 1e-8: break return x