11 votos

Solución numérica de una ecuación integral

Tengo problemas con una solución de una ecuación integral en MATLAB: todas las condiciones son de un doble control, pero la respuesta es incorrecta. Permítame expresar la ecuación: $$ x(s) = g(s)+\int\limits_0^1x(t)f(t-s)\,dt\quad(1) $$ donde $$ f(t) = \frac{\sqrt{2}}{\pi(1+(t+1)^4)} $$ y $$ g(s) = \int\limits_{1-s}^\infty f(t)\,dt. $$ Para la solución de $x$ I sabe que debe estar delimitada $0\leq x(s)\leq 1$ y monótonamente creciente. Para resolver esta ecuación que se utiliza fie caja de herramientas que también es apoyado por muchas publicaciones de K. Atkinson y proporciona a los estrictos límites en el error. La gráfica de la solución: fie

Nota un feo saltar cerca de $0$, así como la no-monotónica comportamiento de la función. Un error que dice ser de menos de $10^{-4}$.

Este resultado no me satisfizo y luego encontré $x$, con un enfoque ingenuo: se hizo una cuadrícula en $[0,1]$ con un paso de $10^{-3}$ y reemplazado $x(s)$ por $x_i = x(s_i)$, $g(s)$ por $g_i = g(s_i)$, $f(t-s)$ por $f_{ij} = f(s_j-s_i)$ y la integral de una suma. Como resultado obtuve $(\mathbf I - \delta \mathbf f)\mathbf w = \mathbf g$. La solución de este sistema es en esta imagen. Con este método un error es menor que $0.01$.

Naive solution

He utilizado fie caja de herramientas para el problema que es muy difícil y el resultado fue perfecto, pero para la actual ecuación de que algo va mal. Espero que sea por el hecho de que he codificado algo incorrecto.

He pasado ya tres días tratando de encontrar un error - pero no hay mucho código y me pregunto si usted me puede ayudar. Tal vez, será mejor si pongo el código aquí para que alguien me ayude a facturar? En el caso de que esta pregunta está fuera superior, lo siento - solo dime, yo voy a borrar.

P. S. Mathematica es capaz de dar una forma explícita para $g$, así que puedo poner aquí por si ayuda.

Editado: Añadido El Código. Yo use el siguiente procedimiento para llamar fie

[soln,errest,cond] = Fie(1,0,1,1,@kPdfNew,@rhsNew)

Primeros números son $\lambda =1$, $a = 0$, $b = 1$ y el comportamiento de $=1$ - suave. Para la función de $g$, uso

function output = rhsNew(x)
n = length(x);
output = x;
for i=1:n
        output(i) = quad(@Poi,1-x(i),100);
end

y para el kernel que puse

function output = kPdfNew(s,t)
z = t-s; 
n = length(z);
output = z;
for i = 1:n
    output(i) = Poi(z(i));
end

donde Pdi es una densidad de $f(t)$ dada por

function output = Poi(z)
n = length(z);
output = z;
for i=1:n
    output(i) = sqrt(2)/pi/(1+(z(i)+1)^4);
end

Puede comentar sobre la integración hasta la $100$ lugar $\infty$ -, pero, en primer lugar, la diferencia de $\int\limits_{100}^\infty f(t)\,dt$ es más pequeño que el de una máquina de precisión y en segundo lugar, antes de utilizar la forma explícita de $g$ y el resultado fue el mismo.

5voto

Grant Puntos 116

Sólo en el caso de alguien tendrá el mismo problema, aquí está la respuesta. El problema estaba en el trabajo con argumento de vectores, funciones de valor de vector en MATLAB, por lo que el kernel debería darse de una manera ligeramente diferente:

function output = kPdfNew(s,t)
z = t-s; 
[nx,ny] = size(s);
output = zeros(size(s));
    for ii = 1:nx
      for jj = 1:ny
        output(ii,jj) = sqrt(2)/pi/(1+(z(ii,jj)+1)^4);
      end
    end
end

5voto

thelsdj Puntos 3344

En caso de rendimiento y legibilidad son cuestiones, puede lograr lo que su respuesta se hizo con el siguiente código, así:

function output = kPdfNew(s,t)
output = sqrt(2)/pi./(1+(t-s+1).^4);
end

(El problema principal con el código original era que itera for i=1:length(z) , en vez de length se podría haber usado numel)

i-Ciencias.com

I-Ciencias es una comunidad de estudiantes y amantes de la ciencia en la que puedes resolver tus problemas y dudas.
Puedes consultar las preguntas de otros usuarios, hacer tus propias preguntas o resolver las de los demás.

Powered by:

X