Integral en el espacio log-log

Estoy trabajando con funciones que, en general, son mucho más suaves y se comportan mejor en el espacio log-log --- así que ahí es donde realizo la interpolación / extrapolación, etc., y eso funciona muy bien. ¿Hay alguna manera de integrar estas funciones numéricas en el espacio log-log?

es decir, espero usar algún tipo de regla trapezoidal simple para realizar una integral acumulativa (por ejemplo, en python, use scipy.integrate.cumtrapz), para encontrar alguna st$F(r)$

$$F(r) = \int_0^r y(x) \, dx$$

Pero espero usar los valores y , en lugar de y (cuando sea posible).l o g ( x ) y $log(y)$$log(x)$$y$$x$

Solo puedes cambiar las variables. Configuración de , . La integral se convierteb ( a ) = l o g ( y ( x ) $a=log(x)$$b(a)=log(y(x))$

$F(r)=\int^{log(r)}_{-\infty} exp(a+b) da$

Debe tener un poco de cuidado porque se está integrando desde . Lo que debe hacer exactamente dependerá de cómo se vea .y ( x $-\infty$$y(x)$

No uso Python, pero si entiendo correctamente, entonces por estás pensando en algo como donde es un vector que muestrea la integral sobre una cuadrícula .F = i n t e g r un t e ( y , x ) F = [ F 1 , . . . , F n ]

$$F(r)=\int_0^ry(x)dx$$ $$\mathbf{F}=\mathrm{integrate}(\mathbf{y},\mathbf{x})$$ $\mathbf{F}=[F_1,...,F_n]$$\mathbf{x}$

Sin embargo, no tiene muestras de e , sino que tiene muestras de y .y x = log ( x ) y = log ( Y $x$$y$$\hat{x}=\log(x)$$\hat{y}=\log(y)$

Por supuesto, el enfoque más simple sería pero esto sería propenso a errores, porque no es uniforme, aunque es.y ( x ) y ( x

$$\mathbf{F}=\mathrm{integrate}(\exp(\mathbf{\hat{y}}),\exp(\mathbf{\hat{x}})) ,$$ $y(x)$$\hat{y}(\hat{x})$

Ahora la regla trapezoidal esencialmente supone que su entrada es lineal por partes. Entonces, la generalización simple sería asumir que es lineal por partes.y ( x$y(x)$$\hat{y}(\hat{x})$

En este caso, al definir , tiene Δ F k = ∫ x k + 1 x k y ( x ) d x = ∫ x k + 1 x k e y ( x ) e x d x = ∫ x k + 1 x$\Delta F_k=F_{k+1}-F_k$

$$\Delta F_k=\int_{x_k}^{x_{k+1}}y(x)dx =\int_{\hat{x}_k}^{\hat{x}_{k+1}}e^{\hat{y}(\hat{x})}e^\hat{x}d\hat{x} =\int_{\hat{x}_k}^{\hat{x}_{k+1}}\tilde{y}(\hat{x})d\hat{x}$$

Luego, definiendo , tiene y , con y .y k + t ≈ y k + t Δ y k ~ y ( t ) ≈ un e b t un = e y k + x k b = Δ y k + Δ x$t=(\hat{x}-\hat{x}_k)/\Delta \hat{x}_k$

$$\hat{y}_{k+t} \approx \hat{y}_k + t \Delta \hat{y}_k$$ $\tilde{y}(t) \approx a e^{bt}$$a=e^{\hat{y}_k+\hat{x}_k}$$b=\Delta \hat{y}_k+\Delta \hat{x}_k$

Entonces la integral se convierte en

$$\Delta F_k \approx a \Delta \hat{x} \int_0^1e^{bt}dt = a \Delta \hat{x} \frac{e^b-1}{b}$$

En Matlab esto se vería algo así

dlogx=diff(logx); dlogy=diff(logy); k=1:length(logx)-1;  
b=dlogx+dlogy; a=exp(logx+logy);  
dF=a(k).*dlogx.*(exp(b)-1)./b;  
F=cumsum([0,dF]);

¡Espero que esto ayude!

(Editar: Mi respuesta es esencialmente la misma que la respuesta mucho más concisa que Damascus Steel dio mientras escribía. La única diferencia es que intenté dar una solución particular para el caso donde la " es un poco -linear discretizado sobre una malla discreta , con )$y(x)$$\hat{y}(\hat{x})$$\mathbf{\hat{x}}$$F(\hat{x}_1)=0$

Si la función se ve fluida en un gráfico log-log, puede interpolar usando una ley de potencia en cada intervalo (las leyes de potencia son, por supuesto, lineales en log-log). Por lo tanto, entre los puntos y bajo el supuesto de que dentro del intervalo , se obtiene y . La contribución a la integral desde el intervalo es entonces $(x_i, y_i)$$(x_{i+1}, y_{i+1})$$y = C_i x^{n_i}$$i$$n_i = \ln(y_i/y_{i+1})/\ln(x_i/x_{i+1})$$C_i = \ln(y_i) - n_i\ln(x_i)$$i$

$$\Delta F_i = \int_{x_i}^{x_{i+1}} C_i x^{n_i}dx = \left\{ \begin{array}{l@{\quad}l} \frac{C_i}{n_i + 1} (x_{i+1}^{n_i+1} - x_i^{n_i+1}), & n_i \neq -1 \\ C_i(\ln x_{i+1} - \ln x_i), & n_i=-1, \end{array} \right.$$ donde obviamente necesita cierta tolerancia para identificar el caso especial en su implementación.$n_i=-1$

Creo que hay un poco de confusión con el cambio de variables en algunas de las respuestas anteriores, así como algunos errores. La integral de una función de registro no es el registro de la integral. Creo que en general es difícil escribir la integral de una función conociendo la integral de su registro. Si alguien sabe cómo hacerlo, estaría interesado.

Mientras tanto, la solución anterior de @ Stefan es la forma de moverse integrando una función en el espacio log-log. El punto de partida es que la función que está tratando es lineal en el espacio log-log para segmentos lo suficientemente pequeños.

Luego se puede escribir la ecuación de la línea en los puntos finales del segmento:

$$\log(y_1)=m_1 \log(x_1) + n_1$$ $$\log(y_2)=m_1 \log(x_1) + n_1$$

donde es la pendiente de la línea y es su intersección en y.$m_1$$n_1$

Restando los dos, uno puede encontrar:

$$m_1=\frac{\log(y_1) -log(y_2)}{\log(x_1)-log(x_2)}$$

y de la sustitución:

$$n_1=\log(y_1)-m_1 \log(x_1)$$

Si en el espacio log-log la ecuación de un segmento está cerca de una línea, en el espacio normal (lineal) la ecuación del segmento está cerca de un exponencial:

$$y(x)\simeq x^{m} e^{n}$$

Si tenemos una formulación analítica para este segmento, es fácil de integrar:

$$\int_{x_1}^{x_2}y(x)dx = \frac{e^{n_1}}{m_1+1}(x_2^{m_1+1}-x_1^{m_1+1}), \,\, \text{for } \, m\neq -1$$

y

$$\int_{x_1}^{x_2}y(x)dx = e^{n_1} \log \frac{x_2}{x_1},\, \,\text{for } \, m = -1$$

Esto se siente un poco como hacer trampa, pero esto es un muestreo en el espacio log-log de modo que podamos aproximar la función en el espacio lineal a un exponencial con parámetros derivados del espacio log-log.

La solución que uso es básicamente una implementación de la regla del trapecio y hace uso de la scipy.misc.logsumexpfunción para mantener la precisión. Si tiene alguna función lnyque devuelve el logaritmo de yentonces puede hacer, por ejemplo:

desde scipy.misc import logsumexp
importar numpy como np

xmin = 1e-15
xmax = 1e-5

# obtener valores de x espaciados logarítmicamente
xvs = np.logspace (np.log10 (xmin), np.log10 (xmax), 10000)

# evalúa tu función en xvs
lys = lny (xvs)

# realizar la integración de la regla de trapecio
deltas = np.log (np.diff (xvs))
logI = -np.log (2.) + logsumexp ([logsumexp (lys [: - 1] + deltas), logsumexp (lys [1:] + deltas)])

El valor logIes el registro de la integral que desea.

Obviamente, esto no funcionará si necesita configurar xmin = 0. Pero, si tiene algún límite inferior positivo distinto de cero para la integral, puede jugar con el número de puntos xvspara encontrar un número donde converja la integral.