26.4 형상·설계 최적화와 역설계

시뮬레이션의 가장 큰 쓰임 중 하나는 "이 설계가 어떻게 동작하는가"가 아니라 "어떤 설계가 가장 좋은가"에 답하는 것이다. 날개 단면, 방열판 핀 배치, 구조물의 재료 분포를 바꿔 가며 성능을 높이는 설계 최적화(design optimization)와, 원하는 성능에서 거꾸로 설계를 찾는 역설계(inverse design)가 그 일이다. 최적화 알고리즘은 10장에서 다뤘다. 이 장은 그 도구를 물리 시뮬레이션에 붙일 때 생기는 판단, 즉 기울기를 어떻게 싸게 얻는가와 대리모델이 찾은 최적해를 어떻게 믿을 수 있게 하는가를 다룬다.

최적화는 모델의 오차를 찾아내는 데 능하다. 대리모델이 어딘가에서 성능을 과대 예측하면 최적화는 정확히 그곳으로 간다. 그래서 설계 최적화에서 대리모델의 오차는 평균적인 크기보다 최적해 근처에서의 부호가 중요하다. 이 장의 계산이 그 현상을 직접 보인다.

학습 목표

설계 최적화 문제를 설계 변수·목적·제약으로 정의한다
adjoint 기울기와 유한차분 기울기의 비용 차이를 설명하고 기울기 기반 설계에 쓴다
대리모델 기반 최적화의 위험을 고충실도 재검증으로 관리한다

설계 문제의 정의

왜 필요한가. "더 좋은 방열판"은 최적화 문제가 아니다. 무엇을 바꿀 수 있고(설계 변수), 무엇을 좋게 하며(목적), 무엇을 지켜야 하는지(제약)를 정해야 알고리즘이 일을 할 수 있다. 이 정의에서 빠진 것은 최적화가 반드시 이용한다. 질량 제약을 빼면 무거운 설계가, 응력 제약을 빼면 부러지는 설계가 "최적"이 된다.

정의. 설계 변수 $\mathbf{p}$ , 상태 $\mathbf{u}$ (시뮬레이션 해), 목적 $J$ , 제약 $g_i$ 에 대해

\min_{\mathbf{p}}\ J(\mathbf{u}, \mathbf{p}) \quad \text{subject to} \quad \mathbf{R}(\mathbf{u}, \mathbf{p}) = \mathbf{0},\quad g_i(\mathbf{u}, \mathbf{p}) \le 0

이다. $\mathbf{R} = \mathbf{0}$ 은 이산화된 지배 방정식(상태 방정식)으로, 설계가 바뀔 때마다 시뮬레이션을 다시 풀어야 한다는 뜻이다. 이것이 일반 최적화 문제와의 가장 큰 차이다. 목적과 제약을 한 번 계산하는 데 시뮬레이션 한 번이 든다.

파라미터화. 설계 변수를 어떻게 잡느냐가 탐색할 수 있는 설계의 범위와 문제의 차원을 정한다.

치수 최적화(sizing): 두께, 폭, 지름 같은 몇 개의 치수. 차원이 작고 해석이 쉽다.
형상 최적화(shape): 경계의 위치를 spline 제어점, 자유형 변형(free-form deformation) 등으로 표현한다. 메시를 형상에 따라 움직여야 한다.
위상 최적화(topology): 영역의 모든 요소에 재료 밀도를 하나씩 두어, 구멍의 개수와 연결 구조까지 바꾼다. 설계 변수가 메시 요소 수만큼(수만–수백만 개) 생긴다(Bendsøe, Sigmund, 2003).

다목적. 공학 설계의 목적은 대개 둘 이상이고 서로 충돌한다. 10.7장의 Pareto 개념을 강철 외팔보에 적용해 본다. 길이 1 m, 끝단 하중 1 kN, 단면은 폭 $b$ ·높이 $h$ 의 직사각형이다. 질량 $m = \rho b h L$ 과 끝단 처짐 $\delta = FL^3/(3EI)$ ( $I = bh^3/12$ )을 모두 줄이고 싶고, 뿌리의 굽힘 응력 $\sigma = 6FL/(bh^2)$ 이 허용 응력 160 MPa를 넘으면 안 된다. 폭 1–5 cm, 높이 2–10 cm의 격자 설계를 모두 계산한다.

import numpy as np

F_n, Lb_m, E_pa, rho_kgm3, sig_allow_pa = 1000.0, 1.0, 200e9, 7850.0, 160e6   # 강철 외팔보, 끝단 하중
b_m, h_m = [a.ravel() for a in np.meshgrid(np.linspace(0.01, 0.05, 41), np.linspace(0.02, 0.10, 81))]
mass_kg = rho_kgm3 * b_m * h_m * Lb_m                        # (3321,)
defl_mm = 4 * F_n * Lb_m**3 / (E_pa * b_m * h_m**3) * 1e3     # 끝단 처짐 δ = F L³ / (3 E I), I = b h³ / 12
ok = 6 * F_n * Lb_m / (b_m * h_m**2) <= sig_allow_pa          # 뿌리의 굽힘 응력 제약
Fobj = np.stack([mass_kg, defl_mm], 1)[ok]                    # (M, 2) 두 목적, 모두 최소화
front = np.array([not np.any(np.all(Fobj <= p, 1) & np.any(Fobj < p, 1)) for p in Fobj])
bf, hf = b_m[ok][front], h_m[ok][front]
print(f"응력 제약을 만족하는 설계 {ok.sum()}개 중 Pareto 최적 {front.sum()}개")
print(f"그중 폭이 최소(1 cm)이거나 높이가 상한 근처(9.5 cm 이상)인 설계: {np.sum(np.isclose(bf, 0.01) | (hf >= 0.095))}개")
for d_max in [8.0, 2.0, 1.0, 0.3]:                            # epsilon-제약: 처짐 <= d_max에서 질량 최소
    feas = Fobj[:, 1] <= d_max
    if not feas.any():
        print(f"처짐 ≤ {d_max:3.1f} mm: 설계 범위 안에 해 없음")
        continue
    i = np.flatnonzero(feas)[np.argmin(Fobj[feas, 0])]
    print(f"처짐 ≤ {d_max:3.1f} mm: 폭 {b_m[ok][i] * 100:.1f} cm, 높이 {h_m[ok][i] * 100:.1f} cm, 질량 {Fobj[i, 0]:5.2f} kg")

응력 제약을 만족하는 설계 2564개 중 Pareto 최적 110개
그중 폭이 최소(1 cm)이거나 높이가 상한 근처(9.5 cm 이상)인 설계: 110개
처짐 ≤ 8.0 mm: 폭 1.0 cm, 높이 6.3 cm, 질량  4.95 kg
처짐 ≤ 2.0 mm: 폭 1.0 cm, 높이 10.0 cm, 질량  7.85 kg
처짐 ≤ 1.0 mm: 폭 2.0 cm, 높이 10.0 cm, 질량 15.70 kg
처짐 ≤ 0.3 mm: 설계 범위 안에 해 없음

Pareto 최적 설계는 모두 두 경계 위에 있다. 처짐이 큰 쪽(가벼운 설계)은 폭을 최소로 두고 높이를 키우고, 높이가 상한에 닿은 뒤에야 폭을 키운다. 굽힘 강성이 $bh^3$ 에 비례하고 질량은 $bh$ 에 비례하므로, 같은 질량이면 높이를 키우는 것이 훨씬 효율적이기 때문이다. 가장 가벼운 설계(처짐 8 mm 이하)는 응력 제약이 정한다. 높이 6.3 cm에서 응력이 허용값 바로 아래다. 처짐 0.3 mm는 설계 범위 안에서 불가능하다.

이 결과에는 공학적으로 중요한 함정이 있다. 최적해가 "높이 상한"에 붙어 있다는 것은 상한을 넓히면 설계가 더 높고 얇아진다는 뜻이다. 그런데 높고 얇은 보는 이 모델에 없는 실패 모드, 즉 옆으로 비틀리며 무너지는 횡비틀림 좌굴(lateral-torsional buckling)에 약하다. 최적화는 모델링하지 않은 실패 모드 쪽으로 설계를 민다. 범위의 상한이 실제로는 그런 제약을 대신하고 있었다면, 그 제약을 명시적으로 넣어야 한다.

gradient 기반 설계

왜 필요한가. 설계 변수가 수십 개를 넘으면 기울기 없는 탐색은 비현실적이다. 기울기를 유한차분으로 구하면 변수 하나마다 시뮬레이션을 한 번 더 돌려야 하므로, 변수 $n$ 개에 시뮬레이션 $n+1$ 번이 든다. 위상 최적화처럼 변수가 수십만 개면 반복 한 번에 시뮬레이션 수십만 번이다.

adjoint 방법. 상태 방정식이 선형 $\mathbf{K}(\mathbf{p})\mathbf{u} = \mathbf{f}$ 이고 목적이 $J = \mathbf{c}^\top\mathbf{u}$ 라 하자. 상태 방정식을 $p_e$ 로 미분하면 $\mathbf{K}\,\partial\mathbf{u}/\partial p_e = -(\partial\mathbf{K}/\partial p_e)\mathbf{u}$ 이므로

\frac{dJ}{dp_e} = \mathbf{c}^\top\frac{\partial\mathbf{u}}{\partial p_e} = -\mathbf{c}^\top\mathbf{K}^{-1}\frac{\partial\mathbf{K}}{\partial p_e}\mathbf{u} = -\boldsymbol{\lambda}^\top\frac{\partial\mathbf{K}}{\partial p_e}\mathbf{u}, \qquad \mathbf{K}^\top\boldsymbol{\lambda} = \mathbf{c}

이다. $\boldsymbol{\lambda}$ 를 adjoint(수반) 변수라 한다. $\boldsymbol{\lambda}$ 는 설계 변수의 개수와 상관없이 선형계 하나로 구해지고, 그 뒤 각 $p_e$ 에 대한 기울기는 값싼 행렬 곱이다. 시뮬레이션 1회 + adjoint 해석 1회로 모든 설계 변수의 기울기가 나온다. 이것은 8.4장의 reverse mode 자동미분과 같은 계산이다(Giles, Pierce, 2000). 대신 목적·제약이 $m$ 개면 adjoint 해석이 $m$ 번 필요하므로, adjoint는 "설계 변수는 많고 출력은 적은" 문제에 유리하다.

직관. $\lambda_i$ 는 "절점 $i$ 에 열을 조금 더 넣으면 목적이 얼마나 바뀌는가"라는 감도다. 목적이 한 지점의 온도라면 $\boldsymbol{\lambda}$ 는 그 지점에서 본 Green 함수다. 그러니 $\boldsymbol{\lambda}$ 를 한 번 구해 두면 어느 요소의 전도도를 바꾸든 그 효과를 곧바로 계산할 수 있다.

길이 1 m인 막대를 요소 50개로 나누고, 요소마다의 열전도도 $k_e$ 를 설계 변수로 둔다. 막대 전체에서 균일하게 열이 나고(1 W/m), 왼쪽 끝은 온도 0의 방열부, 오른쪽 끝은 단열이다. 목적은 가장 뜨거운 오른쪽 끝의 온도 $J = T(L)$ 이다.

n, L_m, q_wm = 50, 1.0, 1.0                       # 요소 수, 막대 길이, 단위 길이당 발열 [W/m]
h = L_m / n                                       # 요소 길이 [m]

def stiffness(k):                                 # 1차 유한요소 전도 행렬, x=0의 T=0인 절점은 뺀다
    main = np.append(k[:-1] + k[1:], k[-1]) / h   # (n,)
    return np.diag(main) - np.diag(k[1:] / h, 1) - np.diag(k[1:] / h, -1)   # (n, n)

f = np.full(n, q_wm * h); f[-1] /= 2              # (n,) 절점 열부하, 오른쪽 끝은 단열
def peak_temp(k):                                  # 목적: 단열된 끝의 온도 J = T(L)
    return np.linalg.solve(stiffness(k), f)[-1]

def grad_adjoint(k):                               # 해석 1회 + adjoint 해석 1회
    K = stiffness(k)
    T = np.concatenate([[0.0], np.linalg.solve(K, f)])               # (n+1,)
    lam = np.concatenate([[0.0], np.linalg.solve(K.T, np.eye(n)[-1])])   # K^T λ = ∂J/∂T
    return -np.diff(lam) * np.diff(T) / h          # (n,) dJ/dk_e = −λ^T (∂K/∂k_e) T

k0 = np.ones(n)                                    # 균일 열전도도 [W/(m·K)]
g_adj = grad_adjoint(k0)
g_fd = np.array([(peak_temp(k0 + 1e-6 * np.eye(n)[e]) - peak_temp(k0)) / 1e-6 for e in range(n)])
print(f"J = T(L) = {peak_temp(k0):.4f} K (해석해 q L²/(2k) = 0.5)")
print("adjoint와 유한차분의 최대 상대 차이 < 1e-5:", bool(np.max(np.abs(g_adj - g_fd)) / np.max(np.abs(g_adj)) < 1e-5))
print(f"선형계 풀이 횟수: adjoint 2회, 전진 유한차분 {n + 1}회")

J = T(L) = 0.5000 K (해석해 q L²/(2k) = 0.5)
adjoint와 유한차분의 최대 상대 차이 < 1e-5: True
선형계 풀이 횟수: adjoint 2회, 전진 유한차분 51회

두 기울기는 전진 유한차분의 절단 오차 수준에서 일치한다(실행해 보면 상대 차이는 약 $10^{-6}$ 이다). 같은 기울기를 adjoint는 선형계 2회, 유한차분은 51회로 얻었다. 설계 변수가 10만 개인 위상 최적화라면 2회 대 100,001회다. 유한차분은 이런 비교에서 adjoint 구현을 검증하는 용도로만 쓴다. 몇 개 성분을 골라 비교하는 이 검사를 빼먹으면, 부호나 인덱스가 틀린 기울기로 최적화가 조용히 엉뚱한 곳에 수렴한다.

이제 이 기울기로 설계한다. 재료 예산 $\sum_e k_e h = 1$ 을 지키면서 $T(L)$ 을 줄인다. 위상 최적화에서 널리 쓰이는 최적성 기준(optimality criteria) 갱신을 쓴다. 재료 단위당 효과 $-\partial J/\partial k_e / h$ 가 큰 요소의 전도도를 키우고, 예산을 맞추는 Lagrange 승수는 이분법으로 찾는다(10.2장의 KKT 조건이 이 갱신의 근거다). 이 문제는 연속 극한에서 해석해가 있다. 위치 $x$ 를 지나는 열유속이 $q(L - x)$ 이므로 $T(L) = \int_0^L q(L-x)/k(x)\,dx$ 이고, Cauchy–Schwarz 부등식으로 최적 분포는 $k \propto \sqrt{L - x}$ , 최적값은 $4qL^3/(9V)$ 다.

V = 1.0                                            # 재료 예산 Σ k_e h = V
k = np.full(n, V / L_m)
for it in range(30):
    B = -grad_adjoint(k) / h                       # (n,) 재료 단위당 온도 감소 효과 (>= 0)
    lo, hi = 1e-9, 1e9
    while hi / lo > 1 + 1e-10:                     # 예산을 맞추는 Lagrange 승수를 이분법으로 찾는다
        mid = np.sqrt(lo * hi)
        k_new = np.clip(k * np.sqrt(B / mid), np.maximum(0.01, 0.5 * k), 2 * k)
        lo, hi = (mid, hi) if k_new.sum() * h > V else (lo, mid)
    k = k_new
    if it in (0, 1, 4, 29):
        print(f"반복 {it + 1:2d}  J = {peak_temp(k):.5f} K")
x = (np.arange(n) + 0.5) * h
k_star = np.sqrt(L_m - x); k_star *= V / (k_star.sum() * h)   # 연속 문제의 최적해 k ∝ sqrt(L − x)
print(f"연속 문제의 최적값 4qL³/(9V) = {4 * q_wm * L_m**3 / (9 * V):.5f} K")
print("찾은 설계가 k ∝ sqrt(L − x)와 일치(최대 차이 < 1e-6):", bool(np.abs(k - k_star).max() < 1e-6))

반복  1  J = 0.44874 K
반복  2  J = 0.44502 K
반복  5  J = 0.44466 K
반복 30  J = 0.44466 K
연속 문제의 최적값 4qL³/(9V) = 0.44444 K
찾은 설계가 k ∝ sqrt(L − x)와 일치(최대 차이 < 1e-6): True

다섯 번 만에 수렴했고, 찾은 설계가 해석해의 분포와 일치한다. 균일 설계(0.5 K)보다 최고 온도가 11% 낮다. 이산 최적값 0.44466 K가 연속 최적값 0.44444 K보다 약간 큰 것은 요소 50개의 이산화 오차다. 열이 많이 지나는 방열부 쪽에 재료를 몰고 끝으로 갈수록 줄이는 설계로, 방열판 핀이 뿌리 쪽에서 두껍고 끝으로 갈수록 가늘어지는 것과 같은 원리다.

미분 가능한 시뮬레이션. solver 전체를 자동미분 프레임워크로 짜면 adjoint를 손으로 유도하지 않아도 된다(24.3장). 그러나 비선형 solver의 반복, 적응적 시간 간격, 메시 재생성처럼 미분하기 까다로운 단계가 있고, 긴 시간 적분의 reverse mode는 메모리를 많이 쓴다. 손으로 쓴 adjoint와 자동미분 중 어느 쪽이든, 유한차분과의 비교 검증은 생략하지 않는다.

대리모델 기반 설계

왜 필요한가. adjoint를 쓸 수 없는 경우가 많다. 상용 solver라 코드를 고칠 수 없거나, 목적이 시뮬레이션 출력의 미분 불가능한 함수이거나, 시뮬레이션 1회가 너무 비싸 기울기가 있어도 반복 수십 번을 감당할 수 없는 경우다. 이때는 몇십 번의 시뮬레이션으로 대리모델(26.2장)을 만들고 대리모델 위에서 최적화한다.

위험: 가짜 최적해. 대리모델 위의 최적화는 대리모델이 성능을 과대 예측하는 곳을 찾아간다. 오차가 평균적으로 작아도, 최적화가 그 오차 중 가장 유리한 것을 골라내기 때문이다. 앞의 막대 문제를 형상 파라미터 하나로 줄여 본다. 열전도도를 $k(x) \propto (L - x + 0.02)^p$ 로 두고(예산으로 정규화), $p \in [0, 0.3]$ 의 네 점에서만 고충실도 해석을 한 뒤 3차 다항식 대리모델로 $p \in [0, 1.5]$ 를 탐색한다.

design = lambda p: (lambda kk: kk * V / (kk.sum() * h))((L_m - x + 0.02) ** p)   # 형상 파라미터 p 하나
P = [0.0, 0.1, 0.2, 0.3]                                       # 고충실도 해석을 한 설계
Js = [peak_temp(design(p)) for p in P]
grid = np.linspace(0.0, 1.5, 151)
coef = np.polyfit(P, Js, 3)                                    # 3차 다항식 대리모델
p_s = grid[np.argmin(np.polyval(coef, grid))]
print(f"대리모델 최적 p = {p_s:.2f}: 예측 J = {np.polyval(coef, p_s):.5f} K, 재검증 J = {peak_temp(design(p_s)):.5f} K")
for it in range(3):                                            # 재검증 결과를 데이터에 넣고 탐색 범위를 제한한다
    lo, hi = min(P), max(P)
    cand = np.linspace(max(0.0, lo - 0.5 * (hi - lo)), hi + 0.5 * (hi - lo), 301)
    coef = np.polyfit(P, Js, 3)
    p_new = cand[np.argmin(np.polyval(coef, cand))]
    P.append(p_new); Js.append(peak_temp(design(p_new)))
    print(f"반복 {it + 1}: p = {p_new:.3f}, 예측 {np.polyval(coef, p_new):.5f} K, 재검증 {Js[-1]:.5f} K")
true_J = [peak_temp(design(p)) for p in grid]
print(f"참고: 이 설계군의 참 최적 p = {grid[np.argmin(true_J)]:.2f}, J = {min(true_J):.5f} K")

대리모델 최적 p = 0.63: 예측 J = 0.44294 K, 재검증 J = 0.44599 K
반복 1: p = 0.450, 예측 0.44621 K, 재검증 0.44661 K
반복 2: p = 0.576, 예측 0.44448 K, 재검증 0.44506 K
반복 3: p = 0.556, 예측 0.44499 K, 재검증 0.44496 K
참고: 이 설계군의 참 최적 p = 0.55, J = 0.44496 K

대리모델은 $p = 0.63$ 에서 0.44294 K를 예측했다. 그러나 이 값은 앞 셀의 연속 문제 최적값 0.44444 K보다 낮다. 어떤 열전도도 분포로도 도달할 수 없는 값이므로, 재검증을 하기 전에도 가짜임을 알 수 있다. 실제로 재검증하면 0.44599 K다. 예측된 개선(11.4%)의 일부는 대리모델의 외삽 오차였다. 재검증 결과를 데이터에 넣고 탐색 범위를 이미 해석한 점들 근처로 제한하는 반복을 세 번 하자, 예측과 재검증이 일치하며 참 최적( $p = 0.55$ , 0.44496 K)에 도달했다. 이것이 신뢰 영역(trust region) 방식의 대리모델 최적화이고, 다음 표본을 예측값과 불확실성의 균형으로 고르는 베이지안 최적화(10.7장; Jones, Schonlau, Welch, 1998)도 같은 "예측 → 고충실도 확인 → 갱신" 루프다.

생성모델 기반 역설계

왜 필요한가. 최적화는 목적 하나에 설계 하나를 찾는다. 그런데 설계자는 "이 성능을 내는 설계들"의 여러 후보를 원할 때가 많다. 제조 가능성, 비용, 미적 요구처럼 목적함수에 넣기 어려운 기준으로 마지막에 고르기 위해서다. 원하는 성능을 조건으로 설계를 생성하는 조건부 생성모델(19장)이 이 역할을 한다.

구성. 시뮬레이션으로 만든 (설계, 성능) 쌍으로 $p(\text{설계} \mid \text{성능})$ 을 학습하고, 목표 성능을 조건으로 설계를 표본 추출한다. 역문제는 대개 해가 여러 개이므로(서로 다른 형상이 같은 양력을 낼 수 있다), 평균 하나를 내는 회귀보다 분포를 내는 생성모델이 맞는 구조다.

제약 만족. 생성모델의 표본은 학습 데이터와 비슷할 뿐, 제약(응력, 최소 두께, 연결성)을 만족한다는 보장이 없다. 그래서 (1) 생성 후 제약을 검사해 걸러 내거나, (2) 생성 결과를 시작점으로 기울기 기반 최적화를 짧게 돌려 제약을 맞추거나, (3) 제약을 만족하는 설계만 표현하도록 파라미터화를 설계한다. 어느 쪽이든 생성된 설계의 성능은 고충실도 시뮬레이션으로 다시 계산한 값으로 보고한다. 생성모델이 "조건으로 받은 성능"은 그 설계의 성능이 아니다.

검증

왜 필요한가. 최적화 결과는 그 결과를 만든 모델(대리모델, 거친 격자, 단순화한 물리) 위에서 최적일 뿐이다. 최적해는 정의상 학습 데이터나 기존 설계에서 가장 먼 곳에 있는 경우가 많아, 모델 오차가 가장 클 가능성이 높은 곳이기도 하다.

검증 절차.

고충실도 재계산: 최적 설계를 학습 데이터 생성에 쓴 것보다 정확한 설정(더 고운 격자, 더 완전한 물리)으로 다시 계산한다. 예측값과 재계산값을 함께 보고한다.
물리적 하한·상한 검사: 셀 4처럼 이론적으로 도달 불가능한 성능을 예측했다면 즉시 가짜다. 에너지 보존, 열역학 법칙, 알려진 이론 한계(Betz 한계 같은)를 검사 기준으로 둔다.
모델링하지 않은 실패 모드 점검: 셀 1의 높고 얇은 보처럼, 최적해가 범위 경계에 붙어 있으면 그 경계가 실제로 무엇을 대신하는지 묻는다.
강건성: 제조 공차, 하중·운전 조건의 변동에서 성능이 급격히 나빠지지 않는지 본다. 최적점이 뾰족한 골짜기에 있으면 실제 제품에서 그 성능을 내지 못한다.
기준 비교: 기존 설계, 단순한 규칙 기반 설계와 같은 고충실도 평가로 비교한다.

실습

대리모델로 찾은 최적 설계를 고충실도 해석으로 재검증해 대리모델 오차가 최적화에 미친 영향을 측정한다. 코드는 labs/physics-and-engineering/design-optimization/surrogate_design.py다. 막대의 열전도도 분포를 제어점 세 개의 로그값으로 파라미터화하고(설계 변수 3개), 요소 200개 해석을 고충실도로 둔다. 표본 8·16·32개의 GP 대리모델과 요소 4개짜리 저충실도 해석 각각으로 최적화한 뒤, 찾은 설계를 고충실도로 재계산한다. 기준은 고충실도 해석을 직접 최적화한 값이다. 집필 시점의 실행 결과는 다음과 같다.

고충실도 직접 최적화(기준): J = 0.44783  z = [ 1.04  0.77 -0.3 ]
균일 설계: J = 0.50000
대리모델                       예측 J     재검증 J     예측 오차     기준 대비 손실  예측 < 이론 하한 4/9
GP, 표본 8개               0.35558   0.51264  -0.15706     +0.06481  예
GP, 표본 16개              0.26201   0.55387  -0.29185     +0.10603  예
GP, 표본 32개              0.40965   0.45948  -0.04983     +0.01165  예
저충실도 해석(요소 4개)          0.45314   0.44800  +0.00514     +0.00016  아니오

결과를 읽는 법은 다음과 같다.

GP 대리모델이 찾은 "최적해"는 모두 가짜다. 예측값이 모두 이론 하한 4/9 ≈ 0.444보다 낮다. 표본 8개와 16개로 찾은 설계는 재검증하면 균일 설계(0.5)보다도 나쁘다. 최적화가 대리모델이 과소 예측한 영역을 정확히 찾아간 결과다.
표본을 늘린다고 단조롭게 나아지지 않는다. 표본 16개가 8개보다 나빴다. 3차원 설계 공간에서 수십 개 표본은 여전히 성기고, 어느 표본 배치가 최적해 근처를 덮느냐에 따라 결과가 크게 달라진다.
물리를 담은 저충실도 모델이 데이터 기반 대리모델보다 나았다. 요소 4개짜리 해석은 예측 오차가 0.005로 작고 재검증 성능이 기준과 거의 같다. 거칠어도 물리 방정식을 푸는 모델은 설계 공간 전체에서 행동이 합리적이어서, 최적화가 이용할 "구멍"이 적다. 26.2장의 다중 충실도가 설계 최적화에서 특히 유용한 이유다.

스크립트를 고쳐 다음을 확인한다. (1) GP 대신 셀 4의 재검증 루프(예측 → 고충실도 확인 → 데이터 추가 → 탐색 범위 제한)를 3차원에서 구현하면 표본 몇 개로 기준에 도달하는가. (2) GP의 예측 평균 대신 "예측 평균 + 2 × 예측 표준편차"를 최소화하면(비관적 탐색) 가짜 최적해가 줄어드는가. (3) 저충실도 모델의 요소 수를 2개로 줄이면 언제부터 최적해가 틀려지는가.

이어서 adjoint 기울기를 직접 구현하는 연습문제를 푼다.

실패 사례

대리모델의 오차를 이용한 '가짜 최적해'를 그대로 보고하는 것. 최적화는 대리모델이 성능을 과대 예측한 곳을 찾아간다. 대리모델 위에서의 예측값이 아니라 고충실도 재계산값을 보고하고, 둘의 차이를 함께 적는다. 예측값이 물리적 한계를 넘으면 그 자체로 실패의 증거다.
adjoint 기울기를 검증하지 않고 쓰는 것. 부호나 인덱스가 틀린 기울기도 최적화는 어딘가로 수렴시킨다. 무작위 방향 몇 개에서 유한차분과 비교하는 검사를 구현의 일부로 둔다.
목적·제약에서 빠진 것을 잊는 것. 질량을 빼면 무거운 설계가, 좌굴을 빼면 얇은 설계가 나온다. 최적해가 설계 범위의 경계에 붙어 있으면 그 경계가 대신하는 실제 제약을 찾는다.
설계 변수가 많은 문제에 유한차분 기울기를 쓰는 것. 변수 $n$ 개에 시뮬레이션 $n+1$ 회가 들어 반복 한 번의 비용이 감당할 수 없게 된다. adjoint나 자동미분을 쓰고, 유한차분은 검증에만 쓴다.

통과 기준

다루는 설계 문제의 설계 변수·목적·제약·파라미터화를 정의하고, 다목적이면 Pareto front와 그 위에서 고른 근거를 보고한다.
설계 변수 $n$ 개, 목적·제약 $m$ 개인 문제에서 adjoint와 유한차분의 시뮬레이션 횟수를 비교하고, adjoint 기울기를 유한차분으로 검증한 결과를 보고한다.
최적화 결과를 고충실도 재검증과 함께 보고한다. 예측값과 재검증값의 차이, 물리적 한계 검사 결과를 적는다.
연습문제 adjoint-gradient의 테스트를 통과한다.

참고 문헌

Giles, M. B., Pierce, N. A. (2000). An introduction to the adjoint approach to design. Flow, Turbulence and Combustion, 65(3–4), 393–415. doi:10.1023/A:1011430410075 — 설계 최적화를 위한 adjoint 방법의 입문
Jameson, A. (1988). Aerodynamic design via control theory. Journal of Scientific Computing, 3(3), 233–260. doi:10.1007/BF01061285 — 공기역학 형상 설계에 adjoint를 쓴 초기 연구
Bendsøe, M. P., Sigmund, O. (2003). Topology Optimization: Theory, Methods, and Applications. Springer. doi:10.1007/978-3-662-05086-6 — 위상 최적화와 최적성 기준 갱신
Jones, D. R., Schonlau, M., Welch, W. J. (1998). Efficient global optimization of expensive black-box functions. Journal of Global Optimization, 13(4), 455–492. doi:10.1023/A:1008306431147 — 비싼 함수의 대리모델 기반 전역 최적화(EGO)