10강 Optimization
사진과 같이 장애물이 존재하는 환경에서 로봇이 A point ⇒ B point로 이동하고자 할 때, 시간을 최소화시키는 최적 경로는 어떻게 구할 수 있을까요? 보행 로봇에게 있어 베터리 효율이 매우 중요합니다. 보행 로봇이 최소한의 에너지를 사용하도록 하는 경로는 어떻게 구할 수 있을까요?
⇒ 이번 시간에는 이러한 질문에 대한 답을 구할 수 있는 “최적화”에 대해 살펴보고자 합니다.
- 최적화를 위해 다음과 같은 방법들을 생각해 볼 수 있습니다.
 
| Method | Description | Drawback | 
|---|---|---|
| 함수의 그래프를 그리기 | 좌표평면 상에 해당 함수의 그래프를 그려 최솟값을 만족하는 해를 구해냅니다. | 차원이 3차원을 넘어가면 직관적으로 이해하기 어렵습니다. | 
| 추측해보기 | 직관적인 추론 | 추론이 들어맞는 특정 상황에서만 사용할 수 있습니다. | 
| 수치적 해 구하기 | 임의의 초기값과 반복적 검증을 통해 최적값을 도출합니다. | local minima에 빠질 위험이 있습니다. | 
⇒ 각각의 방법들이 장단점을 갖지만, local minima에 빠지지 않는다면 수치적 해를 구하는 것이 가장 범용적으로 사용될 수 있습니다.
여기서의 "local minima"라는 것에 대해 간단히 설명을 추가해보겠습니다.
아래 함수의 그래프가 곧 경로 실행에 필요한 비용이라고 생각해 보겠습니다. ( 시간이 될 수도 있고, 에너지가 될 수도 있으며, 혹은 이들을 모두 고려한 무언가일 수도 있지요.) 이 비용의 최솟값을 구하기 위해 그래프를 그렸을 때, 가장 아래에 위치하는 지점이 최소 지점이 될 것입니다.
우리는 파란색 지점이 모든 그래프를 통틀어 최고의 최솟값임을 알고 있지만, 수치적 해 기법은 아래와 같이 본인 근처에서의 최솟값을 발견한 뒤, 추가적인 탐색을 진행하지 않습니다. 이러한 최소 지점을 “local minima”라고 부릅니다.
최적 해를 구하는 문제는 변수의 개수에 따라 단일변수 / 다변수로, 제약조건의 유무에 따라 constrained / unconstrained로, 제약조건의 선형 여부에 따라 linear / non-linear로 나뉘게 됩니다.
Only One Variable?
- one variable
 - multi variable
 
Has Constraint?
- constrained
 - unconstrained
 
Is Linear?
- linear
 - non-linear
 
이번 시간, python의 scipy.optimize를 사용하여 각종 상황에 대한 최적해를 구해보고 약간의 이론적인 내용을 함께 다루어보겠습니다.
scipy.optimize.minimize
scipy.optimize.minimize기본 사용법
| Parameter | Description | 
|---|---|
| fun | 최적화시킬 함수식을 구현하여 전달합니다. | 
| x0 | 최적화를 시작할 초기값을 전달합니다. | 
| args | 함수식에 필요한 추가 매개변수가 있을 시 전달합니다. | 
| method | 사용되는 최적화 알고리즘을 선택합니다. | 
- 간단한 예시를 통해 감을 익혀봅시다. - method를 바꿔가며 실행해보세요!
 
- opt.minimize에 사용될 수 있는 method는 다음과 같은 option들을 갖습니다. 이들에 대해 모두 설명하는 것은 강의의 범위를 벗어나므로 훌륭한 보조자료들을 추천드리고 넘어가겠습니다.
 
- 다크 프로그래머 - 최적화 기법의 직관적 이해
 - 혁펜하임 - BFGS & DFP | 콰지-뉴턴법 (quasi-Newton) 끝판왕!
 - 모두의 연구소 - 모두를 위한 컨벡스 최적화
 
Constrained Optimization
경로 탐색에서의 장애물, 모터의 최대 토크, 속도 등 실제 로봇 개발 시 다양한 제약조건을 고려해야 합니다. 이러한 “제약조건”을 만족하는 최적 해를 구하는 상황을 상상하면서 학습을 이어나갑시다.
- 아래와 같은 함수 $f(x)$ 와 제약조건들이 있다고 가정해보겠습니다.
 
- 각각의 constraint들을 분류해보면 다음과 같습니다.
 
| Description | |
|---|---|
| f(x) | cost | 
| x1 x2 x3 x4 x5 | optimization variable | 
| x1 + x2 + x3 = 5 | linear equality constraint | 
| x3^2 + x4 = 5 | non-linear equality constraint (sin, cos 등 모두 비선형) | 
| x4^2 + x5^2 ≤ 5 | non-linear inequality constraint | 
| x1 ≥ 0.3 | Bound Conditions | 
⇒ 일반적으로 non-linear constraint 보다 linear constraint가 계산이 용이합니다. ( 미분을 계산함에 있어 용이하기 때문이지요!)
Method Adaptation
scipy에서 제공하는 constrained minimization method 중 3종류를 실습해보고자 하며, 각각에 대한 구현과 동작을 통해 최적화에 대한 감을 익혀보도록 하겠습니다.
| Method | Description | 
|---|---|
| COBYLA | 선형화 추론을 사용하며, equality constraint를 처리하지 못합니다. | 
| SLSQP | 구현 시 equal, inequal case를 dictionary로 구분합니다. | 
| trust-constr | 구현 시 linear, non-linear를 구분합니다. | 
cobyla.py
cobyla method는 equality constraint를 처리할 수 없기 때문에 inequality condition 2개를 통해 구현합니다.
💡 ineq condition은 기본적으로 > 0 이 기준입니다.
💡 Bounds를 통해 각 variable들의 최소, 최대값을 정의합니다.
- slsqp.py
 
constraint, bounds를 정의하는 방식은 COBYLA와 동일합니다. 다만 options에서 약간의 차이를 갖습니다.
- trust_constr.py
 
trust_constr는 linear, non-linear constraint에 대해 각기 다른 접근법을 취합니다. 따라서 opt.LinearConstraint, opt.NonlinearConstraint를 통해 constraint를 정의하고 전달하는 모습을 확인할 수 있습니다.
Torque Optimization Example
로봇의 하드웨어는, 지금 우리가 구현하고 있는 시뮬레이션과는 달리 최대 속도, 힘이 제한되어 있습니다. 따라서 로봇의 제어에 있어 이 제한을 넘지 않도록 Optimization하는 것은 매우 중요합니다. 이번 시간의 마지막 예시로, 간단한 pendulum의 예시를 통해 Torque Optimization을 맛보는 시간을 갖고자 합니다.
- 우선, 준비한 예시 프로그램을 실행시켜봅니다. PD 제어를 사용하여 pendulum을 수직으로 일으켜 세우는 제어 시뮬레이션 입니다.
 
cd ./lec10_optimization/torque_optimization python3 onelink_pd.py
- 아래 그래프에는 시간에 따라 변화하는 pendulum의 위치, 속도, 힘이 시각화되어 있습니다. 그런데, Torque 부분을 보면, 200을 넘어 약 250정도까지 도달하는 것을 볼 수 있습니다.
 
- 실제 pendulum과 같은 로봇 시스템에 많이 사용되는 모터의 스펙을 가져와 보았는데요, 최대 torque는 120Nm의 제한을 갖습니다. 따라서 시뮬레이션 상의 동작을 따라가지 못하게 되겠지요
 
image from : t-motor-store
따라서 우리는, 최대 torque의 제한을 고려하는, 그러면서도 pendulum을 원하는 위치에 도달시키는 새로운 제어기가 필요합니다.
- 이 예시를 구현하기 위해선 추가적으로 배워야 할 개념들이 조금 더 있습니다. 따라서 제가 작성한 실습 코드를 우선 실행시킨 뒤, Optimization의 중요성에 대해서 체감하는 시간으로 가볍게 생각하시면 되겠습니다 😊
 
python3 onelink_opt.py
💡 위 그래프에서 보시다시피, 최대 Torque가 훨씬 안정화 된 모습을 확인할 수 있습니다. 비록 목표 지점까지의 도달 시간을 0.5 ⇒ 1.2초로 늘어났지만 말이죠.
두번째 예시인 onelink_opt.py 중에서, optimization과 관련된 코드만 함께 살펴보겠습니다.- 초기 상태는 각도와 각속도 (-pi/2, 0)이며, 목표 상태는 (pi/2, 0)이 됩니다.
 
self.z_end = [np.pi/2, 0] ... z0 = [-np.pi/2, 0]
- 이번 예시에서는 시간과, 제어값인 torque에 대해 최대/최소 제약을 걸어주었습니다. 이 값은 여러분들께서 직접 바꿔보면서 실험해보세요!
 
    time_min, time_max = 1, 4
    u_min, u_max = -20, 120
    # object state (theta1, theta1dot)
    z_end = parms.z_end
    # temporal control inputs
    u_opt = (u_min + (u_max-u_min) * np.random.rand(1, N+1)).flatten()
    # prepare upper/lower bounds
    u_lb = (u_min * np.ones((1, N+1))).flatten()
    u_ub = (u_max * np.ones((1, N+1))).flatten()
    # state x (t, u)
    x0 = [1, *u_opt]
    x_min = [time_min, *u_lb]
    x_max = [time_max, *u_ub]
- 최대, 최소값 조건들과 더불어, pendulum_constraint이라는 제약 조건을 추가하고, 최적화 시 사용되는 Cost Function, cost를 사용하여 제어기를 작성해 보았습니다.
 
    limits = opt.Bounds(x_min, x_max)
    constraint = {
        'type': 'eq',
        'fun': pendulum_constraint
    }
    result = opt.minimize(
        cost, x0, args=(parms), method='SLSQP',
        constraints=[constraint],
        options={'ftol': 1e-6, 'disp': True, 'maxiter': 500},
        bounds=limits
    )
    opt_state = result.x
- Cost Function에 해당하는 함수, cost는 control의 제곱과 시간을 곱하고, 이것을 누적하고 있습니다.
 
def cost(x, args):
    N = OptParams().N
    time = x[0]
    dt = time / N
    u_opt = x[1:]
    tau_sum = sum([x*y for x, y in zip(u_opt, u_opt)]) * dt
    return tau_sum + time
- 수식적으로 나타내면 아래와 같은데요. 여기서 왜 제곱을 하고, 시간에 대해 적분을 해준 것일까요?
 
$$ \min\limits_{u(\cdot)} \int_{t_0}^{t_f} T^2 \,dt \\ $$
- 첫째로, 최적화의 계산 효율을 위해 우리는 시간을 10등분하여 이산적으로 만들어줄 것입니다. 때문에 전체 Torque를 계산하기 위해 적분의 개념이 도입되고, 이에 따라 dt의 곱셈이 추가된 것입니다.
 
  class OptParams:
    def __init__(self):
        self.N = 10
- 다음으로 제곱을 취한 이유는, 부호의 제거를 위함이며, 현재는 제어하는 모터가 1개 뿐이지만, 다수의 모터를 제어하는 경우 아래와 같이 벡터 형태로 표현될 것입니다.
 
$$ \min\limits_{u(\cdot)} \int_{t_0}^{t_f} (T^T T) \,dt \\ $$
💡 더불어 지금 사용한 Cost Function은 Torque에 대한 고려만 추가했지만, 시간, 충격량, 에너지 등 다양한 요소들을 추가하여 Cost Function을 제작할 수 있답니다.
- Time Bound와 Torque Bound, 그리고 Cost Function까지 살펴보았는데요, 그럼 목표 State에 도달시키는 것은 어떻게 구현할까요? - 이것 또한 하나의 Constraint가 됩니다.
 
def pendulum_constraint(x):
    parms = Parameters()
    z_end = parms.z_end
    z_aft, _, _, _ = simulator(x)
    theta1_diff = z_aft[0] - z_end[0]
    theta1dot_diff = z_aft[1] - z_end[1]
    return [theta1_diff, theta1dot_diff]
⇒ Simulation을 거친 결과가 우리의 목표 지점과 일치하도록 하는 Constraint가 되며, 여기서 사용할 수 있는 기법인 collocation, shooting method들은 다음 시간 계속해서 배워보도록 하겠습니다.
Exercise
- pendulum 시스템의 cost function을 다양하게 바꿔보고, 움직임과 그래프를 분석해 봅시다.
 - 강의 중 제시드린 다양한 레퍼런스 사이트들을 살펴보세요 😊
 










