Problem III a

function p=Q3a(p0, tol, N)

    i=1;

   while i<=N

       n=ceil(sqrt(((p0^4))/(6*sqrt(pi)*tol))^1/2)
       h=p0/n;
       u=(1:n-1)*h;
       fexp=exp((-u.^2)/2);
       s=sum(fexp);
       S=(h/2)*(1+2*s+exp((-p0^2)/2));

       f=(2/sqrt(pi))*S-1;
       df=(2/sqrt(pi))*exp((-p0^2)/2);

       p=p0-(f/df)
       if abs(p-p0)<=tol
           disp(['The error is less than the given tolerance after ' num2str(i) ' iterations'])
           return
       end

       i=i+1;
       p0=p;