Implementation of Newton's Method with the function representation that depends on iteration number

Hello everyone in Matlab community ! I have been fighting with my problem for last year and now I think it is impossible to implement what I need it to do in Matlab. My problem is that I need to implement Newton-Raphson algorithm with help of two function Funk and Funkder (the functions are presented below ).I have no problem with that, however the biggest issue is that those two functions depend on number of iteration in Newton algorithm. I wonder if there is any way to submit number of iteration as a part of variable in the mentioned function . I would appreciate any hint , any advice, suggestion ! Thanks in advance, with all respect, Deemitry .

This is Newton method implementation :

      M=10;
      w=M;
      n=ones(M,1);
      S0 = 40; K = 45; r = 0.01; T =6/12;
      sigma = 0.3;
      delta=0.1;
      eps=0.000000001;
      counter=1;
      n(M)=(K*min(1,r/delta));
      while (abs(Funk(n))>eps)
           if (w>1)
          n(w-1)=n(w)-(Funk(n)/Funkder(n))
               w=w-1;
          else
          break
          end
      end

And here comes two functions:

function [o]=Funk(x1)
M=10;
S0 = 40; K = 45; r = 0.01; T =6/12;
sigma = 0.3;
delta=0.1;
delT=T/M;
x(M)=(K*min(1,r/delta));
i=8;
if (i+1)<(M-1)
   for j=(i+1):1:(M-1)
        fsum(j)=delT.*f1(x1(i),x1(j),(j-i)*delT);

     end
     summa1=sum(fsum);
          else
              for j= (M-1):1:(i+1)
                   fsum(j)=delT.*f1(x1(i),x1(j),(j-i)*delT);
              end
                summa1=sum(fsum(j));
     end
%end   
    o=x1(i)-K+p(x1(i),K,(M-i)*delT)+0.5*delT*(f1(x1(i),x(M),(M-i)*delT)+0.5*(r*K-delta*x1(i)))+summa1;
    if~nargout; 
    o;
    end
end

function [l]=Funkder(x3)
S0 = 40; K = 45; r = 0.01; T =6/12;
sigma = 0.3;
delta=0.1;
M=10;
delT=T/M;
i=8;
if (i+1)<(M-1)
for j=(i+1):1:(M-1)
       Fsum1(j)=exp(-delta*(j-i)*delT)*(c((log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)/(sigma*(((j-i)*0.5)^0.5)))*(delta-(r*K)/x3(j))-delta*normcdf(-(log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)));

end
summa=sum(Fsum1)
else
for j= (M-1):1:(i+1)
                 Fsum1(j)=exp(-delta*(j-i)*delT)*(c((log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)/(sigma*(((j-i)*0.5)^0.5)))*(delta-(r*K)/x3(j))-delta*normcdf(-(log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)));

end
summa=sum(Fsum1);
end
   l=1-exp(-delta*(M-i)*delT)*normcdf(-(log(x3(i)/K)+(r-delta+(sigma^2)/2)*(j-i)*delT))+delT*summa;
   if~nargout; 
    l;
   end
end

Answer 1

[Walter Roberson]

Your "while" statement evalutes Funk(n) . Does that count as an "iteration" distinct from when you invoke Funk(n) within the body of the while?

I would say that the Netwon-Raphson method is incompatible with having function values that depend upon the iteration. If you change according to iteration, then the implication is that if you were to change the initial conditions, that the point found would be different -- because the points "move" with each iteration. I could only see it working if the change in function values converged more quickly than N-R converges so that by the time you get down to a sufficiently narrow interval, the points have "frozen" to where they would have been if no movement had taken place at all. Such a thing might work, but I cannot see a point to it.

[Dimitry]

Thank you very much for your answer. The idea that I am trying to implement is that the function value at time T is the starting approximation of Newton method. Then I need to calculate all values for function from T-1,T-2,..... 0. I am trying to implement it by creating a vector with size of unknown function's value and with assigning the approximation value to the last element of the vector. And then by assigning obtained during iterations values to W-1,W-2 .... elements of vector I plan to obtain all function values written in vector. The problem is I can not skip Newton's method because it is required by the task conditions and Funk and Funkder as you probably have seen contain partial sums which depend on which element of n vector ( variables for these function) is evaluated. I wonder if it possible to assign some variable that is responsible for calculating which element of vector is going into the function and in such way assign starting value for partial sums exist in Funk and Funkder.

Thanks again

[Walter Roberson]

Unfortunately I do not understand what you are trying to do.

Could you post the task statement / problem description ? The above is your planned _implementation_, but I do not grasp what _problem_ you are trying to solve.

[Dimitry]

I am very sorry that I did not respond earlier , I simply did not see this comment. Now my problem is solved and the solution was quite easy. Thank you very much for advices and effort!!!!

Answer 2

[Dimitry]

Hello again I changed the functions a little bit.But now when I try to run I x1 at Funkder consists only of one element while it should be 10-element vector. I wonder if someone can help with this one. Thanks in advance!

function [l]=Funkder(x3,y3)
    S0 = 40; K = 45; r = 0.01; T =6/12;
    sigma = 0.3;
    delta=0.1;
    M=10;
    delT=T/M;
    i=y3;
    if i==M
        l=1;
    else
    if (i+1)<(M-1)
    for j=(i+1):1:(M-1)
           Fsum1(j)=exp(-delta*(j-i)*delT)*(c((log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)/(sigma*(((j-i)*0.5)^0.5)))*(delta-(r*K)/x3(j))-delta*normcdf(-(log(x3(i)/x3(j))+(r-delta+(sigma^2)/2)*(j-i)*delT)));       
    end
    summa=sum(Fsum1);
     l=1-exp(-delta*(M-i)*delT)*normcdf(-(log(x3(i)/K)+(r-delta+(sigma^2)/2)*(M-i)*delT)/(sigma*((M-i)*delT)^0.5))+delT*summa;
    else    
               summa=0;
       l=1-exp(-delta*(M-i)*delT)*normcdf(-(log(x3(i)/K)+(r-delta+(sigma^2)/2)*(M-i)*delT)/(sigma*((M-i)*delT)^0.5))+delT*summa;
    end
       if~nargout; 
        l;
       end
    end
    end                

function [o]=Funk(x1,y1)
M=10;
S0 = 40; K = 45; r = 0.01; T =6/12;
sigma = 0.3;
delta=0.1;
delT=T/M;
x(M)=(K*min(1,r/delta));
i=y1;
if i==M
    summa1=0;
    o=x(M);
else
if (i+1)<(M-1) 
   for j=(i+1):1:(M-1)
        fsum(j)=delT.*f1(x1(i),x1(j),(j-i)*delT);

     end
     summa1=sum(fsum);
     o=x1(i)-K+p(x1(i),K,(M-i)*delT)+0.5*delT*(f1(x1(i),x(M),(M-i)*delT)+0.5*(r*K-delta*x1(i)))+summa1;
          else
                   summa1=0;
      o=x1(i)-K+p(x1(i),K,(M-i)*delT)+0.5*delT*(f1(x1(i),x(M),(M-i)*delT)+0.5*(r*K-delta*x1(i)))+summa1;
  end
  end
   if~nargout; 
      o;
      end
  end
  % The implementation of Newton method
  M=10;
  w=M;
  n=M;
  u=rand(M,1);
  S0 = 40; K = 45; r = 0.01; T =6/12;
  sigma = 0.3;
  delta=0.1;
  eps=0.000000001;
  u(w)=(K*min(1,r/delta));
  u
while (abs(Funk(u,n))>eps)
     if (w>1) && (n>1)
    u(w-1)=u(w)-(Funk(u(w),n)/Funkder(u(w),n))
         w=w-1;
         n=n-1

      else
      break
      end
  end