Tamaño de paso de descenso de gradiente adaptable cuando no puede hacer una búsqueda de línea

Tengo una función objetivo depende de un valor , donde es la solución para un PDE. Estoy optimizando por gradiente de descenso en la condición inicial de la PDE: . Es decir, actualizo y luego tengo que integrar el PDE para calcular mi residual. Eso significa que, si tuviera que hacer una búsqueda de línea para el tamaño del paso de descenso del gradiente (llámelo ), para cada valor potencial de tendría que integrar el PDE nuevamente.ϕ ( x , t = 1.0 ) ϕ ( x , t ) $E$$\phi(x, t = 1.0)$$\phi(x, t)$$E$ϕ ( x , t = 0.0 ) α $\phi(x, t = 0.0)$$\phi(x, t = 0.0)$$\alpha$$\alpha$

En mi caso eso sería prohibitivamente costoso. ¿Hay otra opción para el tamaño de paso de descenso de gradiente adaptativo?

No solo estoy buscando esquemas de principios matemáticos aquí (aunque, por supuesto, eso es mejor si algo existe), sino que estaría contento con cualquier cosa que generalmente sea mejor que un tamaño de paso estático.

¡Gracias!

Comenzaré con una observación general: la información de primer orden (es decir, usando solo gradientes, que codifican la pendiente) solo puede brindarle información direccional: puede decirle que el valor de la función disminuye en la dirección de búsqueda, pero no por cuánto tiempo . Para decidir qué tan lejos ir a lo largo de la dirección de búsqueda, necesita información adicional (el descenso del gradiente con longitudes de paso constantes puede fallar incluso en problemas cuadráticos convexos) Para esto, básicamente tienes dos opciones:

1. Use información de segundo orden (que codifica la curvatura), por ejemplo, usando el método de Newton en lugar del descenso de gradiente (para el cual siempre puede usar la longitud de paso suficientemente cerca del minimizador).$1$
2. Prueba y error (por lo que, por supuesto, me refiero a usar una búsqueda de línea adecuada como Armijo).

Si, mientras escribe, no tiene acceso a segundas derivadas, y evaluar la función de obediencia es muy costoso, su única esperanza es comprometerse: use suficiente información aproximada de segundo orden para obtener una buena longitud de paso de candidato, de modo que una línea la búsqueda solo necesita evaluaciones (es decir, como máximo un múltiplo constante (pequeño) del esfuerzo que necesita para evaluar su gradiente).$\mathcal{O}(1)$

Una posibilidad es utilizar las longitudes de paso Barzilai - Borwein (véase, por ejemplo, Fletcher: sobre el método Barzilai-Borwein . Optimización y control con aplicaciones, 235–256, Appl. Optim., 96, Springer, Nueva York, 2005 ). La idea es utilizar una aproximación de diferencia finita de la curvatura a lo largo de la dirección de búsqueda para obtener una estimación del tamaño del paso. Específicamente, elija arbitrario, establezca y luego para :g 0 : = ∇ f ( x 0 ) k = 0 , . . $\alpha_0>0$$g^0:=\nabla f(x^0)$$k=0,...$

1. Establezca yx k + 1 = x k + s $s^k = -\alpha_k^{-1} g^k$$x^{k+1}=x^k+s^k$
2. Evalúe y establezca$g^{k+1}=\nabla f(x^{k+1})$$y^k = g^{k+1}-g^{k}$
3. Establecer$\alpha_{k+1} = \frac{(y^k)^Ty^k}{(y^k)^Ts^k}$

Se puede demostrar que esta opción converge (en la práctica muy rápidamente) para funciones cuadráticas, pero la convergencia no es monótona (es decir, el valor de la función puede ser mayor que , pero solo de vez en cuando; vea el diagrama en la página 10 en el artículo de Fletcher). Para las funciones no cuadráticas, debe combinar esto con una búsqueda de línea, que debe modificarse para tratar la no monotonicidad. Una posibilidad es elegir (por ejemplo, retrocediendo) de modo que donde es el parámetro típico de Armijo y$f(x^{k+1})$$f(x^k)$$\sigma_k \in (0,\alpha_k^{-1})$

$$f(x^k - \sigma_k g^k) \leq \max_{\max(k-M,1)\leq j\leq k} f(x^j) - \gamma \sigma_k (g^k)^Tg^k,$$ $\gamma\in(0,1)$$M$controla el grado de monotonicidad (p. ej., ). También hay una variante que usa valores de gradiente en lugar de valores de función, pero en su caso el gradiente es aún más costoso de evaluar que la función, por lo que no tiene sentido aquí. (Nota: por supuesto, puede intentar aceptar ciegamente las longitudes de los pasos BB y confiar en su suerte, pero si necesita algún tipo de solidez, como escribió en sus comentarios, sería una muy mala idea).$M=10$

Un enfoque alternativo (y, en mi opinión, mucho mejor) sería utilizar esta aproximación de diferencia finita ya en el cálculo de la dirección de búsqueda; Esto se llama un método cuasi-Newton . La idea es construir una aproximación incremental de la de Hesse utilizando diferencias de gradientes. Por ejemplo, podría tomar (la matriz de identidad) y para resolver y establece con como arriba y . (Esto se llama actualización Broyden$\nabla^2 f(x^k)$$H_0=\mathrm{Id}$$k=0,\dots$

$$H_{k}s^{k} = -g^{k},\label{cc1}\tag{1}$$ ykxk+1=xk+s $$H_{k+1} = H_k + \frac{(y^k-H_ks^k)^T(s^k)^T}{(s^k)^Ts^k}$$ $y^k$$x^{k+1} = x^k +s^k$y rara vez se usa en la práctica; Una actualización mejor pero un poco más complicada es la actualización BFGS , para lo cual, y más información, me refiero al libro de Optimización numérica de Nocedal y Wright .) La desventaja es que a) esto requeriría resolver un sistema lineal en cada paso (pero solo del tamaño de lo desconocido, que en su caso es una condición inicial, por lo tanto, el esfuerzo debe dominarse resolviendo PDE para obtener el gradiente; también, existen reglas de actualización para aproximaciones de la arpillera inversa , que solo requieren calcular una matriz única -vector producto) yb) todavía necesita una búsqueda de línea para garantizar la convergencia ...

Afortunadamente, en este contexto existe un enfoque alternativo que hace uso de cada evaluación de la función. La idea es que para simétrico y positivo definido (que está garantizado para la actualización BFGS), resolver es equivalente a minimizar el modelo cuadrático En un método de región de confianza , lo haría con la restricción adicional que , donde es un radio de región de confianza elegido adecuadamente (que desempeña el papel de la longitud del paso ). La idea clave ahora es elegir este radio de forma adaptativa, en función del paso calculado. Específicamente, miras la relación (1) q k ( s ) = $H_k$$\eqref{cc1}$‖S‖≤Δ

$$q_k(s) = \frac12 s^T H_k s + s^T g^k.$$ $\|s\| \leq \Delta_k$$\Delta_k$$\sigma_k$ $$\rho_k := \frac{f(x^k)-f(x^k+s^k)}{f(x^k)-q_k(s^k)}$$ de la reducción real y prevista del valor de la función. Si es muy pequeño, su modelo era malo y descarta e intenta nuevamente con . Si está cerca de , su modelo es bueno y configura y aumenta . De lo contrario, simplemente establezca y deje solo . Para calcular el minimizador real de$\rho_k$$s^k$$\Delta_{k+1}<\Delta_k$$\rho_k$$1$$x^{k+1}=x^k+s^k$$\Delta_{k+1}>\Delta_k$$x^{k+1}=x^k+s^k$$\Delta_k$$s^k$$\min_{\|s\|\leq \Delta_k} q_k(s)$, existen varias estrategias para evitar tener que resolver el problema de optimización con restricciones completas; mi favorito es el método de CG truncado de Steihaug . Para más detalles, nuevamente me refiero a Nocedal y Wright.