function [root,count,err] = NewtonRoot(xin, p, q, tol, maxits)
%
% This function calculates the pth root of q ,
%
%   q^(1/p)
%
% using Newton's method. For simplicity, we assume that p is a 
% positive integer and q is a positive real number
%
% Input:
%   xin : initial guess
%   p : the positive integer degree of the root
%   q : the positive number whose pth root is to estimated
%   tol : stopping tolerance
%   maxits : the maximal number of iterations
%
% Output:
%   root : approximation of the root
%   count : number of newton iterations required to compute the 
%           root
%   err: = 0, if the algorithm proceeded to completion
%        = 1, if an error was encountered
%
  root = NaN;
  count = 0;
  err = 0;
  if int32(p) ~= p || p < 0
    disp('Error: p must be a positive integer'); 
    err = 1;
  return
  end
  if q < 0
    disp('Error: q must be positive');
    err = 1;
  return
  end
  diff = 1.0;
  x = xin;
  while diff > tol && count < maxits
    xo = x;
    x = xo - fn(xo,p,q)/dfn(xo,p,q);
    diff = abs(x-xo);
    count = count+1;
  end
  if count >= maxits
    err = 1;
  end
  root = x;
end

function y = fn(x,p,q)
  y = x^p-q;
end

function y = dfn(x,p,q)
  y = p*x^(p-1);
end

warning: using the gnuplot graphics toolkit is discouraged

The gnuplot graphics toolkit is not actively maintained and has a number
of limitations that are unlikely to be fixed.  Communication with gnuplot
uses a one-directional pipe and limited information is passed back to the
Octave interpreter so most changes made interactively in the plot window
will not be reflected in the graphics properties managed by Octave.  For
example, if the plot window is closed with a mouse click, Octave will not
be notified and will not update its internal list of open figure windows.
The qt toolkit is recommended instead.

xin = 3
p = 2
q = 81
tol = 1E-6
maxits = 100
[root,count,err] = NewtonRoot(xin, p, q, tol, maxits)

xin = 3

p = 2

q = 81

tol = 1.0000e-06

maxits = 100

root = 9
count = 6
err = 0