はじめに #
Rust製の最適制御向け最適化ライブラリOpEnに入門するためチュートリアルの非線形計画問題を解いたので、備忘録を兼ねてまとめた。
OpEnは非凸最適化問題を解くためのライブラリである。 PANOC (Proximal-Averaged Newton-type method for Optimal Control) という最適化手法が実装されている。 最適化問題のソースコードをRustで出力・コンパイルして解くため、計算の実行速度は速い。 また、OpEnはPythonやMATLABのインターフェイスを備えているため、Rustの知識がなくても利用できる。
環境 #
本記事ではOpEnのインターフェイスはPython+CasADiとした。 環境は以下の通り。また、OSはWindows 10である。
- Rust 1.55.0
- Python 3.8.8
- opengen 0.6.9
- OpEn 0.7.1
- CasADi 3.5.5
2022/01/30追記 opengen v0.6.8以前ではPythonからOpEnを実行するとエラーが発生します。
最適制御向け最適化ライブラリOpEnのRust build of TCP interface failedエラーについて – Helve Tech Blog
OpEnのインストール方法は以下の記事を参照。
Rust製最適化ライブラリOpEnのインストール – Helve Tech Blog
また、CasADiは以下でインストールできる。
pip install casadi
最適化問題の定義 #
ここからはチュートリアルの問題に沿って進める。 まず、CasADiとopengenをインポートする。
import casadi.casadi as cs
import opengen as og
次に、変数と目的関数を定義する。目的関数はRosenbrock関数である。
$$ \phi(u) = \sum_{i=0}^{4} \left( p_1 (u_{i+1}-u_i^2)^2 + (p_0 - u_i)^2 \right) $$u = cs.SX.sym("u", 5) # decision variable
p = cs.SX.sym("p", 2) # parameter
phi = og.functions.rosenbrock(u, p) # cost function
uは変数、pはパラメータである。パラメータは最適化計算の実行時に指定する。
なお、OpEnは最適化制御を解くために開発されたため、制御工学で入力を表す記号uを未知変数としていると思われる。
phiは目的関数である。ここでは用意された関数og.functions.rosenbrockを使っているが、自分で定義することも可能である。
例：
phi = sum([p[1]*(u[i+1]-u[i]**2)**2 + (p[0]-u[i])**2 for i in range(4)])
さらに、制約条件を定義する。
$$ u_1^2 + u_2^2 < 1.5^2, \\ -1 < u_3 < 0, -2 < u_4 < 10, -3 < u_5 < -1 $$ball = og.constraints.Ball2(None, 1.5) # ball centered at origin
rect = og.constraints.Rectangle(xmin=[-1,-2,-3], xmax=[0, 10, -1])
segment_ids = [1, 4]
bounds = og.constraints.CartesianProduct(segment_ids, [ball, rect])
problem = og.builder.Problem(u, p, phi).with_constraints(bounds)
Ball2(None, 1.5)は半径1.5の超球内に変数を収める制約である。最初の引数でNoneの代わりに数値を入れると、球の中心を指定できる。
Rectangle(xmin, xmax)は最小値、最大値を指定する制約である。
また、segment_idsで制約式を付加する変数のインデックスを指定している。
ビルドとソルバの設定 #
生成されるオプティマイザのメタ情報を設定する（省略可能）。
meta = og.config.OptimizerMeta() \
.with_version("0.0.0") \
.with_authors(["P. Sopasakis", "E. Fresk"]) \
.with_licence("CC4.0-By") \
.with_optimizer_name("the_optimizer")
生成されたソースコードをRustでビルドするときの設定を指定する。
with_build_modeは"debug"または"release"を指定できる。"release"を指定するとビルドは遅くなるが、最適化の実行は速くなる。
また、.with_tcp_interface_configはソルバをPythonまたはMATLABから呼び出すときに必要。
build_config = og.config.BuildConfiguration() \
.with_build_directory("python_build") \
.with_build_mode("debug") \
.with_tcp_interface_config()
ソルバの設定を行う。
solver_config = og.config.SolverConfiguration() \
.with_lbfgs_memory(15) \
.with_tolerance(1e-5) \
.with_max_inner_iterations(155)
最適化問題、メタ情報、ビルドとソルバの設定をbuilderにまとめる。
builder = og.builder.OpEnOptimizerBuilder(problem,
metadata=meta,
build_configuration=build_config,
solver_configuration=solver_config)
builder.build()
最適化計算の実行 #
PythonとオプティマイザはTCPで通信するため、TCPサーバを起動する。
mng = og.tcp.OptimizerTcpManager('python_build/the_optimizer')
mng.start()
pong = mng.ping() # check if the server is alive
print(pong)
OptimizerTcpManagerには、metaとbuild_configで指定したオプティマイザのパスを指定する。
TCPサーバとの通信が確立できれば、print(pong)で{'Pong': 1}と表示される。
次に、最適化計算を実行する。callにはパラメータの値を指定する。
response = mng.call([1.0, 50.0]) # call the solver over TCP
最適化の結果を確認する。
if response.is_ok():
# Solver returned a solution
solution_data = response.get()
u_star = solution_data.solution
exit_status = solution_data.exit_status
solver_time = solution_data.solve_time_ms
else:
# Invocation failed - an error report is returned
solver_error = response.get()
error_code = solver_error.code
error_msg = solver_error.message
response.get()で最適化の結果が得られる。主なプロパティは以下の通り。
| プロパティ | 説明 |
|---|---|
solution |
最適解 |
exit_status |
状態（Convergedなど） |
solve_time_ms |
実行時間（ミリ秒） |
cost |
評価関数の値 |
num_outer_iterations |
反復回数 |
最後にTCPサーバを停止する。
mng.kill()
最適化結果の確認 #
最後に、最適化結果を確認する。
>>> print(u_star)
[0.4940235276165833, 0.23381718884060976, 0.0, 0.006622525726116431, -1.0]
>>> print(exit_status)
Converged
>>> print(solver_time)
3.0322
u_starは最適解である。また、計算時間は3.0322ミリ秒であった。
