Métodos: discretización, muestreo y convergencia

Este capítulo expone, al estilo de una sección de métodos, los procedimientos computacionales mediante los cuales se extrae un número de la Schubart Master Formula. Procede del balance maestro continuo a su evaluación práctica en cinco pasos numerados: la integral exacta, su reducción a una rodaja representativa, su discretización como suma de cuadratura, su estimación por Monte Carlo y la propagación de la incertidumbre de entrada. Cada relación clave se numera para que el texto circundante y los capítulos posteriores puedan referirse a ella sin ambigüedad como «Ec. (n)». La derivación analítica completa que motiva estos pasos se da en la derivación completa; aquí el énfasis recae en cómo se calcula realmente la integral y en los errores que acarrea cada método.

Ecuación (1) — el balance expandido

El punto de partida es la representación física expandida de la Schubart Master Formula — la identidad contable de la potencia convertida instantánea escrita en su forma espectral–transporte,

$$P(t) = \iint \Phi(E,\theta,t)\,\sigma(E)\,\Delta p\,v_{\text{ph}}\,\eta_{\text{conv}}\;dE\,d\theta \tag{1}$$

donde la integración se extiende sobre la energía incidente $E$ y la coordenada angular $\theta$, y los cinco factores son la densidad de flujo, la sección eficaz de dispersión, la transferencia de momento efectiva, la velocidad de los fonones y la eficiencia de conversión. La Ec. (1) es una integral doble sobre un dominio bidimensional; todo método subsiguiente en este capítulo es una manera de evaluarla. La contabilidad dimensional que muestra que el lado derecho se reduce a vatios se lleva a cabo en el análisis dimensional, y se supone aquí.

Ecuación (2) — rodaja representativa

Cuando el integrando varía lentamente sobre una ventana bien elegida $\Delta E \times \Delta\theta$, o cuando se desea una primera estimación, la Ec. (1) se reduce a una única rodaja representativa reemplazando la integral por el integrando evaluado en valores representativos por el tamaño de la ventana,

$$P \approx \Phi\,\sigma\,\Delta p\,v_{\text{ph}}\,\eta_{\text{conv}}\;\Delta E\,\Delta\theta . \tag{2}$$

La Ec. (2) es la forma de valor medio de la Ec. (1): es exacta si los factores son constantes a lo largo de la ventana y una aproximación controlada en caso contrario. Su validez es exactamente la condición de promediado representativo del Supuesto C en supuestos y limitaciones — la ventana debe ser lo bastante ancha para ser representativa, pero lo bastante estrecha para que el integrando no oscile a lo largo de ella. La Ec. (2) es inestimable para estimaciones de orden de magnitud y para comprobar que un cálculo más completo no se ha desviado, pero no debe confiarse en ella donde el integrando esté agudamente concentrado, por ejemplo a través de una resonancia de la sección eficaz.

Ecuación (3) — discretización por cuadratura

Para evaluar la Ec. (1) sin suponer una única rodaja representativa, el dominio se particiona en una rejilla de celdas indexadas por $j$ (en energía) y $k$ (en ángulo), y la integral se aproxima por una suma de Riemann o de cuadratura sobre esas celdas,

$$P \approx \sum_{j,k} \Phi_{jk}\,\sigma_j\,\Delta p\,v_{\text{ph}}\,\eta_{\text{conv}}\;\Delta E_j\,\Delta\theta_k , \tag{3}$$

donde $\Phi_{jk}$ y $\sigma_j$ son los valores de los factores muestreados en los centros de las celdas y $\Delta E_j$, $\Delta\theta_k$ son las anchuras de las celdas. La Ec. (3) es simplemente la Ec. (2) sumada sobre muchas rodajas; a medida que las celdas se encogen, converge a la Ec. (1). Una regla uniforme del punto medio, tal como está escrita, tiene un error de discretización por celda de orden $(\Delta E)^3$ y $(\Delta\theta)^3$ en la segunda derivada del integrando, de modo que el error global de la regla del punto medio o del trapecio escala como $O(h^2)$ en el tamaño de celda $h$; la regla de Simpson lo mejora a $O(h^4)$ en integrandos suaves. Las rejillas no uniformes que colocan celdas más densamente donde el integrando varía más rápido — en torno a un pico espectral, por ejemplo — recuperan precisión de forma mucho más barata que el refinamiento uniforme. Esta evaluación basada en rejilla es la columna vertebral de la simulación numérica determinista.

Ecuación (4) — el estimador de Monte Carlo

La cuadratura en rejilla se vuelve costosa cuando el integrando es de alta dimensión, está agudamente concentrado o tiene una forma incómoda, porque el número de celdas crece como una potencia de la dimensión. La alternativa de Monte Carlo reescribe la integral como una esperanza y la estima por muestreo aleatorio. Extrayendo $N$ puntos $x_n = (E_n,\theta_n)$ del dominio y escribiendo $f(x)$ para el integrando completo de la Ec. (1) por el volumen del dominio, el estimador es la media muestral

$$\hat P = \frac{1}{N}\sum_{n=1}^{N} f(x_n), \qquad \mathrm{SE}\big[\hat P\big] \propto \frac{1}{\sqrt{N}} . \tag{4}$$

La Ec. (4) es insesgada — su esperanza es igual a la integral verdadera — y su error estándar cae como $1/\sqrt{N}$, independientemente de la dimensión del dominio. Esa independencia de la dimensión es la razón por la que Monte Carlo supera a la cuadratura en rejilla tensorial para integrandos complicados, y es también por lo que el método se extiende naturalmente a la propagación de incertidumbre: extrayendo los propios factores de entrada de sus distribuciones de incertidumbre en un bucle externo, cada pasada de Monte Carlo produce una realización de $P$, y el conjunto reproduce la distribución de salida. El procedimiento anidado completo se expone en la validación de Monte Carlo.

Ecuación (5) — propagación de incertidumbre

El método final no concierne al valor de la integral, sino a su incertidumbre. Dado que la Ec. (1) es un producto de factores, las varianzas relativas de los factores se suman. Para $P = \prod_i x_i$ con factores estadísticamente independientes, una propagación a primer orden da

$$\left(\frac{\sigma_P}{P}\right)^{2} = \sum_i \left(\frac{\sigma_{x_i}}{x_i}\right)^{2} . \tag{5}$$

La Ec. (5) se sigue del log-diferencial $\mathrm{d}\ln P = \sum_i \mathrm{d}\ln x_i$ y se deriva por completo, junto con los términos de covarianza de entradas correlacionadas, en el análisis de incertidumbre. Su papel entre los métodos es decirnos, antes de cualquier muestreo costoso, qué entrada domina el error de salida: dado que los términos se suman en cuadratura, la mayor incertidumbre relativa controla el resultado, y el esfuerzo de refinamiento debe dirigirse allí. En el límite de incertidumbre pequeña, la dispersión de Monte Carlo de la Ec. (4) reproduce la Ec. (5), de modo que ambas sirven como comprobaciones mutuas de consistencia.

Convergencia

Los tres métodos de evaluación — las Ecs. (2), (3) y (4) — acarrean un comportamiento de convergencia cualitativamente diferente, y elegir entre ellos es cuestión de ajustar el método al integrando.

La rodaja representativa, Ec. (2), no converge en ningún sentido de refinamiento; es una aproximación fija cuyo error lo fija enteramente cuán representativa sea la ventana elegida. Es una primera estimación, no un esquema convergente.

La cuadratura, Ec. (3), converge como una potencia del tamaño de celda: $O(h^2)$ para las reglas del punto medio y del trapecio, $O(h^4)$ para la regla de Simpson, en integrandos suaves. Halvar el tamaño de celda recorta, por tanto, el error del trapecio en un factor de cuatro. Esta convergencia rápida y de alto orden hace de la cuadratura el método de elección para integrandos suaves y de baja dimensión — precisamente el régimen de la Ec. (1) bidimensional cuando los factores se comportan bien. Su debilidad es que el recuento de celdas crece como una potencia de la dimensión, y que los picos agudos fuerzan rejillas muy finas.

Monte Carlo, Ec. (4), converge como $1/\sqrt{N}$: reducir el error diez veces exige un aumento de cien veces en las muestras. Esto es lento comparado con la cuadratura en bajas dimensiones, pero el exponente $-\tfrac12$ es independiente de la dimensión, de modo que Monte Carlo gana para integrandos de alta dimensión o agudamente estructurados donde la rejilla de la cuadratura explota. Su ritmo lento se mejora rutinariamente mediante reducción de varianza: el muestreo por importancia, que concentra las extracciones donde el integrando es mayor y reponderara por $f/g$; el muestreo estratificado, que garantiza la cobertura de subregiones importantes; y las variables de control, que restan una magnitud correlacionada de integral conocida. Cada una de ellas rebaja la varianza efectiva $\mathrm{Var}[f]$ en el numerador del error estándar sin cambiar el escalado $1/\sqrt{N}$, de modo que un $N$ modesto logra lo que de otro modo requeriría órdenes de magnitud más de muestras uniformes.

En la práctica, una evaluación robusta de la Ec. (1) usa las tres: la Ec. (2) para un orden de magnitud de comprobación, la Ec. (3) para el valor determinista cuando el integrando es suave, y la Ec. (4) con reducción de varianza para la propagación de incertidumbre y para cualquier régimen en el que la rejilla se vuelva impracticable. La ley de propagación Ec. (5) suministra entonces el presupuesto de error analítico que se espera que el muestreo reproduzca.

Lecturas adicionales

• W. H. Press, S. A. Teukolsky, W. T. Vetterling y B. P. Flannery, Numerical Recipes: The Art of Scientific Computing, 3.ª ed. (Cambridge University Press, 2007) — sobre las reglas de cuadratura, sus órdenes de error y la integración de Monte Carlo con reducción de varianza.
• N. Metropolis y S. Ulam, «The Monte Carlo Method», Journal of the American Statistical Association 44 (1949), 335–341 — el artículo fundacional del método que subyace a la Ec. (4).
• C. P. Robert y G. Casella, Monte Carlo Statistical Methods, 2.ª ed. (Springer, 2004) — un tratamiento moderno y exhaustivo de los estimadores, la convergencia y la reducción de varianza.
• J. R. Taylor, An Introduction to Error Analysis, 2.ª ed. (University Science Books, 1997) — sobre la ley de propagación de la Ec. (5).