formapolinomioinverso.m

function F = formapolinomioinverso(nombrefuncion, K, grado, X)
%formapolinomioinverso
%   F = formapolinomioinverso(funcion)
%
%   Construye una aproximaci�n en forma de recetas en tiempo discreto de
%   unas recetas en tiempo continuo, definidas en un script o con el
%   comando inline; un ejemplo de este ultimo caso seria  
%      F = inline('[-280+575/k,-12; -120,-8+20/k; -280+400/k, -12]')
%
%   F = formapolinomioinverso(funcion, K, grado)
%      funcion   nombre de la funci�n con las recetas en tiempo continuo
%      K         factores sobre los que se calcula la aproximaci�n
%      grado     n�mero de pasos temporales para el tiempo discreto
%      F         aproximaci�n con recetas en tiempo discreto
%
%   Si K es un escalar la aproximaci�n se hace en el intervalo [K/2, K*3/2]

if nargin < 2, K = []; end
if nargin < 3, grado = []; end
if nargin < 4, X = []; end

if isempty(K), K = 1; end                                %si K no est� definido es 1
if size(K(:),1) == 1, K = (K/2:K/5:K*3/2); end           %si K es un escalar la aproximaci�n se efectua en [K/2, K*3/2] con 6 puntos
puntos = size(K, 2);                                     %determina el n�mero de factores para los que se hace la aproximaci�n
if isempty(grado), grado = fix((puntos+3)/2); end        %define el grado del polinomio inverso por defecto
if grado > puntos, grado = puntos; end                   %el grado debe ser menor o igual al numero de puntos analizados
                                                         %
if isempty(X)                                            %si el problema est� escrito con una matriz de flujos actualizados
    F0 = feval(nombrefuncion, K(1));                         %calcula la magnitud de los flujos actualizados para K(1)
else                                                     %si el problema est� escrito con funciones de producci�n actualizadas
    F0 = feval(nombrefuncion, K(1), X);                      %calcula la magnitud de los flujos actualizados para K(1) y X
end                                                      %
[filas, columnas] = size(F0);                            %establece el tama�o del problema
Fk = zeros(filas, columnas, puntos);                     %crea la variable para almacenar los resultados
Fk(:, :, 1) = F0;                                        %almacena la magnitud de los flujos actualizados para K(1)
for c1 = 2:puntos                                        %para el resto de K analizadas
   if isempty(X)                                             %si el problema est� escrito con una matriz de flujos actualizados
       Fk(:,:,c1) = feval(nombrefuncion, K(c1));                         %calcula la magnitud de los flujos actualizados para K(1)
   else                                                      %si el problema est� escrito con funciones de producci�n actualizadas
       Fk(:,:,c1) = feval(nombrefuncion, K(c1), X);                      %calcula la magnitud de los flujos actualizados para K(1) y X
   end                                                       %
end                                                      %
x = 1 ./ K;                                              %define x como 1 / K analizadas
y = zeros(1, puntos);                                    %define y
for c1 = 1:filas                                         %para cada fila
    for c2 = 1:columnas                                      %para cada columna
        y(:) = Fk(c1, c2, :);                                     %extrae y
        F(c1, c2, :) = ajusta(x, y, puntos, grado);               %calcula los coeficientes del polinomio inverso	
    end                                                      %
end                                                      %






function A = ajusta(x, y, fil, col)                      %funci�n para el c�lculo de los coeficientes de un polinomio inverso
M = zeros(fil, col);                                     %crea la matriz del sistema lineal
for c1 = 1:fil                                           %para cada fila
    for c2 = 1:col                                           %para cada columna
        M(c1, col-c2+1) = x(c1) .^ (col-c2);                     %establece la matriz
    end                                                      %
end                                                      %
A = (M \ y')';                                           %resuelve el sistema lineal para establecer los coeficientes del polinomio