diff -r 20017c4feae6 src/odepkg_auxiliary_functions.cc
--- a/src/odepkg_auxiliary_functions.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_auxiliary_functions.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -129,18 +129,18 @@
    octave_value vy, octave_value_list vextarg, octave_idx_type vdeci) {
 
   ColumnVector vresult  (vy.vector_value ());
-  ColumnVector vreduced (vy.length ());
+  ColumnVector vreduced (vy.numel ());
 
   // Check if the user has set the option "OutputSel" then create a
   // reduced vector that stores the desired values.
-  if (vsel.is_empty ()) {
+  if (vsel.isempty ()) {
     for (octave_idx_type vcnt = 0; vcnt < vresult.numel (); vcnt++)
       vreduced(vcnt) = vresult(vcnt);
   }
   else {
-    vreduced.resize (vsel.length ());
+    vreduced.resize (vsel.numel ());
     ColumnVector vselect (vsel.vector_value ());
-    for (octave_idx_type vcnt = 0; vcnt < vsel.length (); vcnt++)
+    for (octave_idx_type vcnt = 0; vcnt < vsel.numel (); vcnt++)
       vreduced(vcnt) = vresult(static_cast<int> (vselect(vcnt)-1));
   }
 
@@ -153,25 +153,25 @@
   else if (vdeci == 1) varin(2) = "";
   else if (vdeci == 2) varin(2) = "done";
 
-  for (octave_idx_type vcnt = 0; vcnt < vextarg.length (); vcnt++)
+  for (octave_idx_type vcnt = 0; vcnt < vextarg.numel (); vcnt++)
     varin(vcnt+3) = vextarg(vcnt);
 
   // Evaluate the output function and return the value of the output
   // function to the caller function
   if ((vdeci == 0) || (vdeci == 2)) {
     if (vplt.is_function_handle () || vplt.is_inline_function ())
-      feval (vplt.function_value (), varin, 0);
+      octave::feval (vplt.function_value (), varin, 0);
     else if (vplt.is_string ()) // String may be used from the caller
-      feval (vplt.string_value (), varin, 0);
+      octave::feval (vplt.string_value (), varin, 0);
     return (true);
   }
 
   else if (vdeci == 1) {
     octave_value_list vout;
     if (vplt.is_function_handle () || vplt.is_inline_function ())
-      vout = feval (vplt.function_value (), varin, 1);
+      vout = octave::feval (vplt.function_value (), varin, 1);
     else if (vplt.is_string ()) // String may be used if set automatically
-      vout = feval (vplt.string_value (), varin, 1);
+      vout = octave::feval (vplt.string_value (), varin, 1);
     return (vout(0).bool_value ());
   }
 
@@ -201,7 +201,7 @@
 
   // If vjac is a cell array then we expect that two matrices are
   // returned to the caller function, we can't check for this before
-  if (vjac.is_cell () && (vjac.length () == 2)) {
+  if (vjac.iscell () && (vjac.numel () == 2)) {
     varout(0) = vjac.cell_value ()(0);
     varout(1) = vjac.cell_value ()(1);
     if (!varout(0).is_matrix_type () || !varout(1).is_matrix_type ()) {
@@ -218,10 +218,10 @@
     varin(1) = vy;  // varin(1).print_with_name (octave_stdout, "vy", true);
     varin(2) = vdy; // varin(2).print_with_name (octave_stdout, "vdy", true);
     // Fill up RHS arguments with extra arguments that are given
-    for (octave_idx_type vcnt = 0; vcnt < vextarg.length (); vcnt++)
+    for (octave_idx_type vcnt = 0; vcnt < vextarg.numel (); vcnt++)
       varin(vcnt+3) = vextarg(vcnt);
     // Evaluate the Jacobian function and return results
-    varout = feval (vjac.function_value (), varin, 1);
+    varout = octave::feval (vjac.function_value (), varin, 1);
   }
 
   // In principle this is not possible because odepkg_structure_check
@@ -266,10 +266,10 @@
     varin(0) = vt;
     varin(1) = vy;
     // Fill up RHS arguments with extra arguments that are given
-    for (octave_idx_type vcnt = 0; vcnt < vextarg.length (); vcnt++)
+    for (octave_idx_type vcnt = 0; vcnt < vextarg.numel (); vcnt++)
       varin(vcnt+2) = vextarg(vcnt);
     // Evaluate the Jacobian function and return results
-    varout = feval (vjac.function_value (), varin, 1);
+    varout = octave::feval (vjac.function_value (), varin, 1);
     vret = varout(0);
   }
 
@@ -313,25 +313,25 @@
   else if (vmass.is_function_handle () || vmass.is_inline_function ()) {
     octave_value_list varin;
     octave_value_list varout;
-    if (vstate.is_empty () || !vstate.is_string ())
+    if (vstate.isempty () || !vstate.is_string ())
       error_with_id ("OdePkg:InvalidOption",
         "If \"Mass\" value is a handle then \"MStateDependence\" must be given");
  
     else if (vstate.string_value ().compare ("none") == 0) {
       varin(0) = vt;
-      for (octave_idx_type vcnt = 0; vcnt < vextarg.length (); vcnt++)
+      for (octave_idx_type vcnt = 0; vcnt < vextarg.numel (); vcnt++)
         varin(vcnt+1) = vextarg(vcnt);
     }
 
     else { // If "MStateDependence" is "weak" or "strong"
       varin(0) = vt; varin(1) = vy;
       // Fill up RHS arguments with extra arguments that are given
-      for (octave_idx_type vcnt = 0; vcnt < vextarg.length (); vcnt++)
+      for (octave_idx_type vcnt = 0; vcnt < vextarg.numel (); vcnt++)
         varin(vcnt+2) = vextarg(vcnt);
     }
 
     // Evaluate the Mass function and return results
-    varout = feval (vmass.function_value (), varin, 1);
+    varout = octave::feval (vmass.function_value (), varin, 1);
     vret = varout(0);
   }
 
@@ -368,13 +368,13 @@
 
   octave_scalar_map vretval;
 
-  if (vstats.length () < 5)
+  if (vstats.numel () < 5)
     error_with_id ("OdePkg:InvalidArgument",
       "C++ function odepkg_auxiliary_makestats error");
   else {
     vretval.assign ("nsteps",   vstats(0));
     vretval.assign ("nfailed",  vstats(1));
-    vretval.assign ("nfevals",  vstats(2));
+    vretval.assign ("noctave::fevals",  vstats(2));
     vretval.assign ("npds",     vstats(3));
     vretval.assign ("ndecomps", vstats(4));
     vretval.assign ("nlinsols", vstats(5));
@@ -411,7 +411,7 @@
   // in vsel if vdeci == (0 || 1).
   RowVector vrow;
   if (vdeci != 2) {
-    vrow.resize (vy.length ());
+    vrow.resize (vy.numel ());
     vrow = RowVector (vy.vector_value ());
   }
 
diff -r 20017c4feae6 src/odepkg_octsolver_ddaskr.cc
--- a/src/odepkg_octsolver_ddaskr.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_ddaskr.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -111,28 +111,28 @@
   // function that keeps the set of implicit differential equations
   octave_value_list varin;
   varin(0) = T; varin(1) = A; varin(2) = APRIME;
-  for (F77_INT vcnt = 0; vcnt < vddaskrextarg.length (); vcnt++)
+  for (octave_idx_type vcnt = 0; vcnt < vddaskrextarg.numel (); vcnt++)
     varin(vcnt+3) = vddaskrextarg(vcnt);
@@ -228,10 +222,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;                      // The OdePkg options structure
 
@@ -277,7 +272,7 @@
 
   // Check if the third and the fourth input argument (check for
   // vector already was successful before) have the same length
-  if (args(2).length () != args(3).length ()) {
+  if (args(2).numel () != args(3).numel ()) {
     error_with_id ("OdePkg:InvalidArgument",
       "Third and fourth input argument must have the same length");
     return (vretval);
@@ -289,8 +284,8 @@
   if (nargin >= 5) {
 
     // Fifth input argument != OdePkg option, need a default structure
-    if (!args(4).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(4).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();        // Create a default structure
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
         vddaskrextarg(vcnt-4) = args(vcnt); // Save arguments in vddaskrextarg
@@ -300,7 +295,7 @@
     else if (nargin > 5) { 
       octave_value_list varin;
       varin(0) = args(4); varin(1) = "odekdi";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();        // Create structure of args(4)
       for (octave_idx_type vcnt = 5; vcnt < nargin; vcnt++)
@@ -311,7 +306,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(4); varin(1) = "odekdi"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -319,7 +314,7 @@
   } // if (nargin >= 5)
 
   else { // if nargin == 4, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -331,13 +326,13 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"RelTol\" not set, new value 1e-6 is used");
   }
   else if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -347,13 +342,13 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"AbsTol\" not set, new value 1e-6 is used");
   }
   else if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -362,9 +357,9 @@
 
   // Setting the tolerance type that depends on the types (scalar or
   // vector) of the options RelTol and AbsTol
-  F77_INT vitol;
+  octave_idx_type vitol;
   if (vreltol.is_scalar_type () && vabstol.is_scalar_type ()) vitol = 0;
-  else if (vreltol.length () == vabstol.length ()) vitol = 1;
+  else if (vreltol.numel () == vabstol.numel ()) vitol = 1;
   else {
     error_with_id ("OdePkg:InvalidOption", 
       "Options \"RelTol\" and \"AbsTol\" must have same length");
@@ -374,7 +369,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"NormControl\" will be ignored by this solver");
@@ -382,14 +377,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   octave_value vplot = vodeopt.contents ("OutputFcn");
-  if (vplot.is_empty () && nargout == 0) vplot = "odeplot";
+  if (vplot.isempty () && nargout == 0) vplot = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -409,13 +404,13 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
-    if (!vinitstep.is_empty ())
+  if (args(1).numel () > 2) {
+    if (!vinitstep.isempty ())
       warning_with_id ("OdePkg:InvalidOption",
        "Option \"InitialStep\" will be ignored if fixed time stamps are given");
     vinitstep = args(1).vector_value ()(1);
   }
-  else if (vinitstep.is_empty ()) {
+  else if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value 1e-6 is used");
@@ -424,7 +419,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -439,19 +434,19 @@
 
   // The options 'Jacobian', 'JPattern' and 'Vectorized'
   octave_value vjac = vodeopt.contents ("Jacobian");
-  if (!vjac.is_empty ()) vddaskrjacfun = vjac;
+  if (!vjac.isempty ()) vddaskrjacfun = vjac;
 
   // The option Mass will be ignored by this solver. We can't handle
   // Mass-matrix options with IDE problems
   octave_value vmass = vodeopt.contents ("Mass");
-  if (!vmass.is_empty ())
+  if (!vmass.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"Mass\" will be ignored by this solver");
 
   // The option MStateDependence will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmst = vodeopt.contents ("MStateDependence");
-  if (!vmst.is_empty ())
+  if (!vmst.isempty ())
     if (vmst.string_value ().compare ("weak") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -459,14 +454,14 @@
   // The option MvPattern will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MvPattern will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"MassSingular\" will be ignored by this solver");
@@ -474,14 +469,14 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // Implementation of the option MaxOrder has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (vmaxder.is_empty ()) {
+  if (vmaxder.isempty ()) {
     vmaxder = 3;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxOrder\" not set, new value 3 is used");
@@ -494,7 +489,7 @@
 
   // The option BDF will be ignored because this is a BDF solver
   octave_value vbdf = vodeopt.contents ("BDF");
-  if (!vbdf.is_empty ())
+  if (!vbdf.isempty ())
     if (vbdf.string_value () != "on") {
       vbdf = "on"; warning_with_id ("OdePkg:InvalidOption", 
         "Option \"BDF\" set \"off\", new value \"on\" is used");
@@ -504,12 +499,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -523,39 +518,39 @@
   NDArray vRTOL   = vreltol.array_value ();
   NDArray vATOL   = vabstol.array_value ();
 
-  vddaskrneqn = TO_F77_INT (args(2).length ());
+  vddaskrneqn = args(2).numel ();
   double T       = vTIME(0);
   double TEND    = vTIME(vTIME.numel () - 1);
   double *Y      = vY.fortran_vec ();
   double *YPRIME = vYPRIME.fortran_vec ();

:
:
:
   
   INFO[0]   = 0;      // Define that it is the initial first call
   INFO[1]   = vitol;  // RelTol/AbsTol are scalars or vectors
   INFO[2]   = 1;      // An intermediate output is wanted
   INFO[3]   = 0;      // Integrate behind TEND
-  if (!vjac.is_empty ()) INFO[4]   = 1;
+  if (!vjac.isempty ()) INFO[4]   = 1;
   else INFO[4] = 0;   // Internally calculate a Jacobian? 0..yes
   INFO[5]   = 0;      // Have a full Jacobian matrix? 0..yes
   INFO[6]   = 1;      // Use the value for maximum step size
@@ -571,14 +566,14 @@
   // initialize these IO-functions for further use
   octave_value vtim (T); octave_value vsol (vY); octave_value vyds (vYPRIME);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vplot.is_empty ()) 
+  if (!vplot.isempty ()) 
     odepkg_auxiliary_evalplotfun (vplot, voutsel, args(1), args(2), vddaskrextarg, 0);
 
   octave_value_list veveideargs;
   veveideargs(0) = vsol; 
   veveideargs(1) = vyds;
   Cell veveidearg (veveideargs);
-  if (!vevents.is_empty ())
+  if (!vevents.isempty ())
     odepkg_auxiliary_evaleventfun (vevents, vtim, veveidearg, vddaskrextarg, 0);
 
   // We are calling the core solver here to intialize all variables
@@ -618,17 +613,17 @@
       }
       vsol = vcres; vyds = vydrs; vtim = T;
 
-      if (!vevents.is_empty ()) {
+      if (!vevents.isempty ()) {
         veveideargs(0) = vsol;
         veveideargs(1) = vyds;
         veveidearg = veveideargs;
         veveres = odepkg_auxiliary_evaleventfun (vevents, vtim, veveidearg, vddaskrextarg, 1);
-        if (!veveres(0).cell_value ()(0).is_empty ())
+        if (!veveres(0).cell_value ()(0).isempty ())
           if (veveres(0).cell_value ()(0).int_value () == 1) {
             ColumnVector vttmp = veveres(0).cell_value ()(2).column_vector_value ();
             Matrix vrtmp = veveres(0).cell_value ()(3).matrix_value ();
@@ -638,7 +633,7 @@
           }
       }
 
-      if (!vplot.is_empty ()) {
+      if (!vplot.isempty ()) {
         if (odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vddaskrextarg, 1)) {
           error ("Missing error message implementation");
           return (vretval);
@@ -655,16 +650,16 @@
 
   veveideargs(0) = vsol;
   veveideargs(1) = vyds;
   veveidearg = veveideargs;
-  if (!vevents.is_empty ())
+  if (!vevents.isempty ())
     odepkg_auxiliary_evaleventfun (vevents, vtim, vsol, vddaskrextarg, 2);
-  if (!vplot.is_empty ())
+  if (!vplot.isempty ())
     odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vddaskrextarg, 2);
 
   // Return the results that have been stored in the
@@ -697,7 +692,7 @@
     vretmap.assign ("solver", "odekdi");
     if (vstats.string_value ().compare ("on") == 0)
       vretmap.assign ("stats", vstatinfo);
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretmap.assign ("ie", veveres(0).cell_value ()(1));
       vretmap.assign ("xe", veveres(0).cell_value ()(2));
       vretmap.assign ("ye", veveres(0).cell_value ()(3));
@@ -715,7 +710,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretval(2) = veveres(0).cell_value ()(2);
       vretval(3) = veveres(0).cell_value ()(3);
       vretval(4) = veveres(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_mebdfdae.cc
--- a/src/odepkg_octsolver_mebdfdae.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_mebdfdae.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -274,10 +268,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;                      // The OdePkg options structure
 
@@ -319,8 +314,8 @@
   if (nargin >= 4) {
 
     // Fourth input argument != OdePkg option, need a default structure
-    if (!args(3).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(3).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();       // Create a default structure
       for (octave_idx_type vcnt = 3; vcnt < nargin; vcnt++)
         vmebdfdaeextarg(vcnt-3) = args(vcnt); // Save arguments in vmebdfdaeextarg
@@ -330,7 +325,7 @@
     else if (nargin > 4) {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "odebda";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();       // Create structure from args(4)
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
@@ -341,7 +336,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "odebda"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -349,7 +344,7 @@
   } // if (nargin >= 4)
 
   else { // if nargin == 3, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -361,14 +356,14 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"RelTol\" not set, new value %3.1e is used", 
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -378,14 +373,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"AbsTol\" not set, new value %3.1e is used", 
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -403,7 +398,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"NormControl\" will be ignored by this solver");
@@ -411,14 +406,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   octave_value vplot = vodeopt.contents ("OutputFcn");
-  if (vplot.is_empty () && nargout == 0) vplot = "odeplot";
+  if (vplot.isempty () && nargout == 0) vplot = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -438,13 +433,13 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
-    if (!vinitstep.is_empty ())
+  if (args(1).numel () > 2) {
+    if (!vinitstep.isempty ())
       warning_with_id ("OdePkg:InvalidOption",
        "Option \"InitialStep\" will be ignored if fixed time stamps are given");
     vinitstep = args(1).vector_value ()(1);
   }
-  if (vinitstep.is_empty ()) {
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -454,7 +449,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -469,16 +464,16 @@
 
   // The options 'Jacobian', 'JPattern' and 'Vectorized'
   octave_value vjac = vodeopt.contents ("Jacobian");
  octave_idx_type vmebdfdaejac = 22; // We need to set this if no Jac available
-  if (!vjac.is_empty ()) {
+  if (!vjac.isempty ()) {
     vmebdfdaejacfun = vjac; vmebdfdaejac = 21;
   }
 
   // Implementation of the option 'Mass' has been finished, these
   // options can be set by the user to another value than default
   vmebdfdaemass = vodeopt.contents ("Mass");
   octave_idx_type vmebdfdaemas = 0;
-  if (!vmebdfdaemass.is_empty ()) {
+  if (!vmebdfdaemass.isempty ()) {
     vmebdfdaemas = 1;
     if (vmebdfdaemass.is_function_handle () || vmebdfdaemass.is_inline_function ())
       warning_with_id ("OdePkg:InvalidOption",
@@ -488,7 +483,7 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   vmebdfdaemassstate = vodeopt.contents ("MStateDependence");
-  if (!vmebdfdaemassstate.is_empty ())
+  if (!vmebdfdaemassstate.isempty ())
     if (vmebdfdaemassstate.string_value ().compare ("weak") != 0) // 'weak' is default
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -496,14 +491,14 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MassSingular will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -511,14 +506,14 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // Implementation of the option MaxOrder has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (vmaxder.is_empty ()) {
+  if (vmaxder.isempty ()) {
     vmaxder = 3;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxOrder\" not set, new value %1d is used", 
@@ -542,12 +537,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -560,34 +555,34 @@
   NDArray vATOL   = vabstol.array_value ();
   NDArray vY0     = args(2).array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
+  octave_idx_type N = args(2).numel ();
   double T0 = vTIME(0);
   double HO = vinitstep.double_value ();
   double *Y0 = vY0.fortran_vec ();
   double TOUT = T0 + vinitstep.double_value ();
   double TEND = vTIME(vTIME.numel ()-1);

:
:
:   
   
+  if (!vmebdfdaemass.isempty ()) MASBND[0] = 1;

@@ -598,9 +593,9 @@
   // etc. and initialize the plot, events and the solstore functions
   octave_value vtim (T0); octave_value vsol (vY0);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vplot.is_empty ()) odepkg_auxiliary_evalplotfun 
+  if (!vplot.isempty ()) odepkg_auxiliary_evalplotfun 
     (vplot, voutsel, args(1), args(2), vmebdfdaeextarg, 0);
-  if (!vevents.is_empty ()) odepkg_auxiliary_evaleventfun 
+  if (!vevents.isempty ()) odepkg_auxiliary_evaleventfun 
     (vevents, vtim, args(2), vmebdfdaeextarg, 0);
 
   // We are calling the core solver here to intialize all variables
@@ -648,12 +643,12 @@
         vcres(vcnt) = Y0[vcnt];
       vsol = vcres; vtim = TOUT;
-      if (!vevents.is_empty ()) {
+      if (!vevents.isempty ()) {
         veveres = odepkg_auxiliary_evaleventfun (vevents, vtim, vsol, vmebdfdaeextarg, 1);
-        if (!veveres(0).cell_value ()(0).is_empty ())
+        if (!veveres(0).cell_value ()(0).isempty ())
           if (veveres(0).cell_value ()(0).int_value () == 1) {
             ColumnVector vttmp = veveres(0).cell_value ()(2).column_vector_value ();
             Matrix vrtmp = veveres(0).cell_value ()(3).matrix_value ();
@@ -662,7 +657,7 @@
             TOUT = TEND; // let's get out here, the Events function told us to finish
           }
       }
-      if (!vplot.is_empty ()) {
+      if (!vplot.isempty ()) {
         if (odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfdaeextarg, 1)) {
           error ("Missing error message implementation");
           return (vretval);
@@ -695,12 +690,12 @@
 
         vcres(vcnt) = Y0[vcnt];
       vsol = vcres; vtim = TOUT;
-      if (!vevents.is_empty ()) {
+      if (!vevents.isempty ()) {
         veveres = odepkg_auxiliary_evaleventfun (vevents, vtim, vsol, vmebdfdaeextarg, 1);
-        if (!veveres(0).cell_value ()(0).is_empty ())
+        if (!veveres(0).cell_value ()(0).isempty ())
           if (veveres(0).cell_value ()(0).int_value () == 1) {
             ColumnVector vttmp = veveres(0).cell_value ()(2).column_vector_value ();
             Matrix vrtmp = veveres(0).cell_value ()(3).matrix_value ();
@@ -709,7 +704,7 @@
             vtimecnt = vtimelen; // let's get out here, the Events function told us to finish
           }
       }
-      if (!vplot.is_empty ()) {
+      if (!vplot.isempty ()) {
         if (odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfdaeextarg, 1)) {
           error ("Missing error message implementation");
           return (vretval);
@@ -725,11 +720,11 @@
 
   vsol = vcres; vtim = TOUT;
-  if (!vevents.is_empty ())
+  if (!vevents.isempty ())
     odepkg_auxiliary_evaleventfun (vevents, vtim, vsol, vmebdfdaeextarg, 2);
-  if (!vplot.is_empty ())
+  if (!vplot.isempty ())
     odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfdaeextarg, 2);
 
   // Return the results that have been stored in the
@@ -763,7 +758,7 @@
     vretmap.assign ("solver", "odebda");
     if (vstats.string_value () == "on") 
       vretmap.assign ("stats", vstatinfo);
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretmap.assign ("ie", veveres(0).cell_value ()(1));
       vretmap.assign ("xe", veveres(0).cell_value ()(2));
       vretmap.assign ("ye", veveres(0).cell_value ()(3));
@@ -781,7 +776,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretval(2) = veveres(0).cell_value ()(2);
       vretval(3) = veveres(0).cell_value ()(3);
       vretval(4) = veveres(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_mebdfi.cc
--- a/src/odepkg_octsolver_mebdfi.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_mebdfi.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -319,23 +319,17 @@
@@ -350,10 +344,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;                      // The OdePkg options structure
 
@@ -399,7 +394,7 @@
 
   // Check if the third and the fourth input argument (check for
   // vector already was successful before) have the same length
-  if (args(2).length () != args(3).length ()) {
+  if (args(2).numel () != args(3).numel ()) {
     error_with_id ("OdePkg:InvalidArgument",
       "Third and fourth input argument must have the same length");
     return (vretval);
@@ -411,8 +406,8 @@
   if (nargin >= 5) {
 
     // Fifth input argument != OdePkg option, need a default structure
-    if (!args(4).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(4).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();        // Create a default structure
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
         vmebdfiextarg(vcnt-4) = args(vcnt); // Save arguments in vmebdfiextarg
@@ -422,7 +417,7 @@
     else if (nargin > 5) { 
       octave_value_list varin;
       varin(0) = args(4); varin(1) = "odebdi";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();        // Create structure of args(4)
       for (octave_idx_type vcnt = 5; vcnt < nargin; vcnt++)
@@ -433,7 +428,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(4); varin(1) = "odebdi"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -441,7 +436,7 @@
   } // if (nargin >= 5)
 
   else { // if nargin == 4, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -453,14 +448,14 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"RelTol\" not set, new value %3.1e is used", 
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -470,14 +465,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"AbsTol\" not set, new value %3.1e is used", 
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption", 
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -495,7 +490,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"NormControl\" will be ignored by this solver");
@@ -503,14 +498,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   octave_value vplot = vodeopt.contents ("OutputFcn");
-  if (vplot.is_empty () && nargout == 0) vplot = "odeplot";
+  if (vplot.isempty () && nargout == 0) vplot = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -530,13 +525,13 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
-    if (!vinitstep.is_empty ())
+  if (args(1).numel () > 2) {
+    if (!vinitstep.isempty ())
       warning_with_id ("OdePkg:InvalidOption",
        "Option \"InitialStep\" will be ignored if fixed time stamps are given");
     vinitstep = args(1).vector_value ()(1);
   }
-  if (vinitstep.is_empty ()) {
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -546,7 +541,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -561,22 +556,22 @@
 
   // The options 'Jacobian', 'JPattern' and 'Vectorized'
   octave_value vjac = vodeopt.contents ("Jacobian");
   octave_idx_type vmebdfijac = 22; // We need to set this if no Jac available
-  if (!vjac.is_empty ()) {
+  if (!vjac.isempty ()) {
     vmebdfijacfun = vjac; vmebdfijac = 21;
   }
 
   // The option Mass will be ignored by this solver. We can't handle
   // Mass-matrix options with IDE problems
   octave_value vmass = vodeopt.contents ("Mass");
-  if (!vmass.is_empty ())
+  if (!vmass.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"Mass\" will be ignored by this solver");
 
   // The option MStateDependence will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmst = vodeopt.contents ("MStateDependence");
-  if (!vmst.is_empty ())
+  if (!vmst.isempty ())
     if (vmst.string_value ().compare ("weak") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -584,14 +579,14 @@
   // The option MvPattern will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MvPattern will be ignored by this solver. We
   // can't handle Mass-matrix options with IDE problems
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -599,14 +594,14 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // Implementation of the option MaxOrder has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (vmaxder.is_empty ()) {
+  if (vmaxder.isempty ()) {
     vmaxder = 3;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxOrder\" not set, new value %1d is used", 
@@ -630,12 +625,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -649,7 +644,7 @@
   NDArray vY0     = args(2).array_value ();
   NDArray vYPRIME = args(3).array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
+  octave_idx_type N = args(2).numel ();
   double T0 = vTIME(0);
   double HO = vinitstep.double_value ();
   double *Y0 = vY0.fortran_vec ();
@@ -740,16 +735,16 @@
         vcres(vcnt) = Y0[vcnt]; vydrs(vcnt) = YPRIME[vcnt]; 
       }
       vsol = vcres; vyds = vydrs; vtim = TOUT;
-      if (!vevents.is_empty ()) {
+      if (!vevents.isempty ()) {
         veveideargs(0) = vsol;
         veveideargs(1) = vyds;
 	veveidearg = veveideargs;
         veveres = odepkg_auxiliary_evaleventfun (vevents, vtim, veveidearg, vmebdfiextarg, 1);
-        if (!veveres(0).cell_value ()(0).is_empty ())
+        if (!veveres(0).cell_value ()(0).isempty ())
           if (veveres(0).cell_value ()(0).int_value () == 1) {
             ColumnVector vttmp = veveres(0).cell_value ()(2).column_vector_value ();
             Matrix vrtmp = veveres(0).cell_value ()(3).matrix_value ();
@@ -758,7 +753,7 @@
             TOUT = TEND; // let's get out here, the Events function told us to finish
           }
       }
-      if (!vplot.is_empty ()) {
+      if (!vplot.isempty ()) {
         if (odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfiextarg, 1)) {
           error ("Missing error message implementation");
           return (vretval);
@@ -791,16 +786,16 @@
         vcres(vcnt) = Y0[vcnt]; vydrs(vcnt) = YPRIME[vcnt];
       }
       vsol = vcres; vyds = vydrs; vtim = TOUT;
-      if (!vevents.is_empty ()) {
+      if (!vevents.isempty ()) {
         veveideargs(0) = vsol;
         veveideargs(1) = vyds;
 	veveidearg = veveideargs;
         veveres = odepkg_auxiliary_evaleventfun (vevents, vtim, veveidearg, vmebdfiextarg, 1);
-        if (!veveres(0).cell_value ()(0).is_empty ())
+        if (!veveres(0).cell_value ()(0).isempty ())
           if (veveres(0).cell_value ()(0).int_value () == 1) {
             ColumnVector vttmp = veveres(0).cell_value ()(2).column_vector_value ();
             Matrix vrtmp = veveres(0).cell_value ()(3).matrix_value ();
@@ -809,7 +804,7 @@
             vtimecnt = vtimelen; // let's get out here, the Events function told us to finish
           }
       }
-      if (!vplot.is_empty ()) {
+      if (!vplot.isempty ()) {
         if (odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfiextarg, 1)) {
           error ("Missing error message implementation");
           return (vretval);
@@ -825,16 +820,16 @@
   vsol = vcres; vyds = vydrs; vtim = TOUT;
   veveideargs(0) = vsol;
   veveideargs(1) = vyds;
   veveidearg = veveideargs;
-  if (!vevents.is_empty ())
+  if (!vevents.isempty ())
     odepkg_auxiliary_evaleventfun (vevents, vtim, veveidearg, vmebdfiextarg, 2);
-  if (!vplot.is_empty ())
+  if (!vplot.isempty ())
     odepkg_auxiliary_evalplotfun (vplot, voutsel, vtim, vsol, vmebdfiextarg, 2);
 
   // Return the results that have been stored in the
@@ -868,7 +863,7 @@
     vretmap.assign ("solver", "odebdi");
     if (vstats.string_value () == "on") 
       vretmap.assign ("stats", vstatinfo);
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretmap.assign ("ie", veveres(0).cell_value ()(1));
       vretmap.assign ("xe", veveres(0).cell_value ()(2));
       vretmap.assign ("ye", veveres(0).cell_value ()(3));
@@ -886,7 +881,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vevents.is_empty ()) {
+    if (!vevents.isempty ()) {
       vretval(2) = veveres(0).cell_value ()(2);
       vretval(3) = veveres(0).cell_value ()(3);
       vretval(4) = veveres(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_radau.cc
--- a/src/odepkg_octsolver_radau.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_radau.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -197,10 +197,10 @@
   octave_value vy = octave_value (A);
 
   // Check if an 'Events' function has been set by the user
-  if (!vradauevefun.is_empty ()) {
+  if (!vradauevefun.isempty ()) {
     vradauevesol = odepkg_auxiliary_evaleventfun 
       (vradauevefun, vt, vy, vradauextarg, 1);
-    if (!vradauevesol(0).cell_value ()(0).is_empty ())
+    if (!vradauevesol(0).cell_value ()(0).isempty ())
       if (vradauevesol(0).cell_value ()(0).int_value () == 1) {
         ColumnVector vttmp = vradauevesol(0).cell_value ()(2).column_vector_value ();
         Matrix vrtmp = vradauevesol(0).cell_value ()(3).matrix_value ();
@@ -216,15 +216,15 @@
 
   // Check if an 'OutputFcn' has been set by the user (including the
   // values of the options for 'OutputSel' and 'Refine')
-  if (!vradaupltfun.is_empty ()) {
+  if (!vradaupltfun.isempty ()) {
     if (vradaurefine.int_value () > 0) {
       ColumnVector B(N); double vtb = 0.0;
@@ -322,8 +317,8 @@
   if (nargin >= 4) {
 
     // Fourth input argument != OdePkg option, need a default structure
-    if (!args(3).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(3).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();       // Create a default structure
       for (octave_idx_type vcnt = 3; vcnt < nargin; vcnt++)
         vradauextarg(vcnt-3) = args(vcnt); // Save arguments in vradauextarg
@@ -333,7 +328,7 @@
     else if (nargin > 4) {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "ode2r";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();       // Create structure from args(4)
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
@@ -344,7 +339,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "ode2r"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -352,7 +347,7 @@
   } // if (nargin >= 4)
 
   else { // if nargin == 3, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -364,14 +359,14 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"RelTol\" not set, new value %3.1e is used",
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -381,14 +376,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"AbsTol\" not set, new value %3.1e is used",
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -397,10 +392,10 @@
 
   // Setting the tolerance type that depends on the types (scalar or
   // vector) of the options RelTol and AbsTol
-  F77_INT vitol = 0;
-  if (vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  octave_idx_type vitol = 0;
+  if (vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 0;
-  else if (!vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  else if (!vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 1;
   else {
     error_with_id ("OdePkg:InvalidOption",
@@ -411,7 +406,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"NormControl\" will be ignored by this solver");
@@ -419,14 +414,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   vradaupltfun = vodeopt.contents ("OutputFcn");
-  if (vradaupltfun.is_empty () && nargout == 0) vradaupltfun = "odeplot";
+  if (vradaupltfun.isempty () && nargout == 0) vradaupltfun = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -443,11 +438,11 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
+  if (args(1).numel () > 2) {
     error_with_id ("OdePkg:InvalidOption",
       "Fixed time stamps are not supported by this solver");
   }
-  if (vinitstep.is_empty ()) {
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -457,7 +452,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -472,20 +467,20 @@
   // Implementation of the option 'Jacobian' has been finished, these
   // options can be set by the user to another value than default
   vradaujacfun = vodeopt.contents ("Jacobian");
  octave_idx_type vradaujac = 0; // We need to set this if no Jac available
-  if (!vradaujacfun.is_empty ()) vradaujac = 1;
+  if (!vradaujacfun.isempty ()) vradaujac = 1;
 
   // The option JPattern will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vradaujacpat = vodeopt.contents ("JPattern");
-  if (!vradaujacpat.is_empty ())
+  if (!vradaujacpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"JPattern\" will be ignored by this solver");
 
   // The option Vectorized will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vradauvectorize = vodeopt.contents ("Vectorized");
-  if (!vradauvectorize.is_empty ())
+  if (!vradauvectorize.isempty ())
     if (vradauvectorize.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"Vectorized\" will be ignored by this solver");
@@ -493,8 +488,8 @@
   // Implementation of the option 'Mass' has been finished, these
   // options can be set by the user to another value than default
   vradaumass = vodeopt.contents ("Mass");
-  F77_INT vradaumas = 0;
-  if (!vradaumass.is_empty ()) {
+  octave_idx_type vradaumas = 0;
+  if (!vradaumass.isempty ()) {
     vradaumas = 1;
     if (vradaumass.is_function_handle () || vradaumass.is_inline_function ())
       warning_with_id ("OdePkg:InvalidOption",
@@ -504,7 +499,7 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   vradaumassstate = vodeopt.contents ("MStateDependence");
-  if (!vradaumassstate.is_empty ())
+  if (!vradaumassstate.isempty ())
     if (vradaumassstate.string_value ().compare ("weak") != 0) // 'weak' is default
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -512,14 +507,14 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MassSingular will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -527,21 +522,21 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // The option MaxOrder will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (!vmaxder.is_empty ())
+  if (!vmaxder.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxOrder\" will be ignored by this solver");
 
   // The option BDF will be ignored by this solver, the core Fortran
   // solver doesn't support this option
   octave_value vbdf = vodeopt.contents ("BDF");
-  if (!vbdf.is_empty ())
+  if (!vbdf.isempty ())
     if (vbdf.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"BDF\" will be ignored by this solver");
@@ -550,12 +545,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -567,31 +562,31 @@
   NDArray vRTOL = vreltol.array_value ();
   NDArray vATOL = vabstol.array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
+  octave_idx_type N = args(2).numel ();
   double X = args(1).vector_value ()(0);
   double* Y = vY0.fortran_vec ();
   double XEND = args(1).vector_value ()(1);
@@ -602,9 +597,9 @@
   octave_value vtim = args(1).vector_value ()(0);
   octave_value vsol = args(2);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vradaupltfun.is_empty ()) odepkg_auxiliary_evalplotfun 
+  if (!vradaupltfun.isempty ()) odepkg_auxiliary_evalplotfun 
     (vradaupltfun, vradauoutsel, args(1), args(2), vradauextarg, 0);
-  if (!vradauevefun.is_empty ())
+  if (!vradauevefun.isempty ())
     odepkg_auxiliary_evaleventfun (vradauevefun, vtim, args(2), vradauextarg, 0);
 
   // We are calling the core solver and solve the set of ODEs or DAEs
@@ -641,9 +636,9 @@
   octave_value vted = octave_value (XEND);
   octave_value vfin = octave_value (vlastline);
 
-  if (!vradaupltfun.is_empty ()) odepkg_auxiliary_evalplotfun
+  if (!vradaupltfun.isempty ()) odepkg_auxiliary_evalplotfun
     (vradaupltfun, vradauoutsel, vted, vfin, vradauextarg, 2);
-  if (!vradauevefun.is_empty ()) odepkg_auxiliary_evaleventfun
+  if (!vradauevefun.isempty ()) odepkg_auxiliary_evaleventfun
     (vradauevefun, vted, vfin, vradauextarg, 2);
   
   // Get the stats information as an octave_scalar_map if the option 'Stats'
@@ -669,9 +664,9 @@
     vretmap.assign ("x", vtres);
     vretmap.assign ("y", vyres);
     vretmap.assign ("solver", "ode2r");
-    if (!vstatinfo.is_empty ()) // Event implementation
+    if (!vstatinfo.isempty ()) // Event implementation
       vretmap.assign ("stats", vstatinfo);
-    if (!vradauevefun.is_empty ()) {
+    if (!vradauevefun.isempty ()) {
       vretmap.assign ("ie", vradauevesol(0).cell_value ()(1));
       vretmap.assign ("xe", vradauevesol(0).cell_value ()(2));
       vretmap.assign ("ye", vradauevesol(0).cell_value ()(3));
@@ -689,7 +684,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vradauevefun.is_empty ()) {
+    if (!vradauevefun.isempty ()) {
       vretval(2) = vradauevesol(0).cell_value ()(2);
       vretval(3) = vradauevesol(0).cell_value ()(3);
       vretval(4) = vradauevesol(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_radau5.cc
--- a/src/odepkg_octsolver_radau5.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_radau5.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -114,29 +109,29 @@
   // function that keeps the set of differential algebraic equations
   octave_value_list varin;
   varin(0) = X; varin(1) = A;
-  for (F77_INT vcnt = 0; vcnt < vradau5extarg.length (); vcnt++)
+  for (octave_idx_type vcnt = 0; vcnt < vradau5extarg.numel (); vcnt++)
     varin(vcnt+2) = vradau5extarg(vcnt);
-  octave_value_list vout = feval (vradau5odefun.function_value (), varin, 1);
+  octave_value_list vout = octave::feval (vradau5odefun.function_value (), varin, 1);
 
   // Return the results from the function evaluation to the Fortran
   // solver, again copy them and don't just create a Fortran vector
@@ -202,10 +197,10 @@
   octave_value vy = octave_value (A);
 
   // Check if an 'Events' function has been set by the user
-  if (!vradau5evefun.is_empty ()) {
-    vradau5evesol = odepkg_auxiliary_evaleventfun
+  if (!vradau5evefun.isempty ()) {
+    vradau5evesol = odepkg_auxiliary_evaleventfun 
       (vradau5evefun, vt, vy, vradau5extarg, 1);
-    if (!vradau5evesol(0).cell_value ()(0).is_empty ())
+    if (!vradau5evesol(0).cell_value ()(0).isempty ())
       if (vradau5evesol(0).cell_value ()(0).int_value () == 1) {
         ColumnVector vttmp = vradau5evesol(0).cell_value ()(2).column_vector_value ();
         Matrix vrtmp = vradau5evesol(0).cell_value ()(3).matrix_value ();
@@ -221,15 +216,15 @@
 
   // Check if an 'OutputFcn' has been set by the user (including the
   // values of the options for 'OutputSel' and 'Refine')
-  if (!vradau5pltfun.is_empty ()) {
+  if (!vradau5pltfun.isempty ()) {
     if (vradau5refine.int_value () > 0) {
       ColumnVector B(N); double vtb = 0.0;
@@ -307,10 +271,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;               // The OdePkg options structure
 
@@ -352,8 +317,8 @@
   if (nargin >= 4) {
 
     // Fourth input argument != OdePkg option, need a default structure
-    if (!args(3).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(3).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();       // Create a default structure
       for (octave_idx_type vcnt = 3; vcnt < nargin; vcnt++)
         vradau5extarg(vcnt-3) = args(vcnt); // Save arguments in vradau5extarg
@@ -363,7 +328,7 @@
     else if (nargin > 4) {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "ode5r";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();       // Create structure from args(4)
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
@@ -374,7 +339,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "ode5r"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -382,7 +347,7 @@
   } // if (nargin >= 4)
 
   else { // if nargin == 3, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -394,16 +359,16 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"RelTol\" not set, new value %3.1e is used",
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
        "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
      }
   }
@@ -411,14 +376,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"AbsTol\" not set, new value %3.1e is used",
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -427,10 +392,10 @@
 
   // Setting the tolerance type that depends on the types (scalar or
   // vector) of the options RelTol and AbsTol
  octave_idx_type vitol = 0;
-  if (vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  if (vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 0;
-  else if (!vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  else if (!vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 1;
   else {
     error_with_id ("OdePkg:InvalidOption",
@@ -441,7 +406,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"NormControl\" will be ignored by this solver");
@@ -449,14 +414,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   vradau5pltfun = vodeopt.contents ("OutputFcn");
-  if (vradau5pltfun.is_empty () && nargout == 0) vradau5pltfun = "odeplot";
+  if (vradau5pltfun.isempty () && nargout == 0) vradau5pltfun = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -470,31 +435,14 @@
 +  if (args(1).numel () > 2) {
+    error_with_id ("OdePkg:InvalidOption",
+      "Fixed time stamps are not supported by this solver");
+  }
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -504,10 +452,10 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && tSpan.numel() == 2) {
-    vmaxstep = (tSpan(1) - tSpan(0)) / 12.5;
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
+    vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption",
-     "Option \"MaxStep\" not set, new value %3.1e is used",
       vmaxstep.double_value ());
   }
 
@@ -519,20 +467,20 @@
   // Implementation of the option 'Jacobian' has been finished, these
   // options can be set by the user to another value than default
   vradau5jacfun = vodeopt.contents ("Jacobian");
  octave_idx_type vradau5jac = 0; // We need to set this if no Jac available
-  if (!vradau5jacfun.is_empty ()) vradau5jac = 1;
+  if (!vradau5jacfun.isempty ()) vradau5jac = 1;
 
   // The option JPattern will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vradau5jacpat = vodeopt.contents ("JPattern");
-  if (!vradau5jacpat.is_empty ())
+  if (!vradau5jacpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"JPattern\" will be ignored by this solver");
 
   // The option Vectorized will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vradau5vectorize = vodeopt.contents ("Vectorized");
-  if (!vradau5vectorize.is_empty ())
+  if (!vradau5vectorize.isempty ())
     if (vradau5vectorize.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"Vectorized\" will be ignored by this solver");
@@ -540,8 +488,8 @@
   // Implementation of the option 'Mass' has been finished, these
   // options can be set by the user to another value than default
   vradau5mass = vodeopt.contents ("Mass");
  octave_idx_type vradau5mas = 0;
-  if (!vradau5mass.is_empty ()) {
+  if (!vradau5mass.isempty ()) {
     vradau5mas = 1;
     if (vradau5mass.is_function_handle () || vradau5mass.is_inline_function ())
       warning_with_id ("OdePkg:InvalidOption",
@@ -551,7 +499,7 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   vradau5massstate = vodeopt.contents ("MStateDependence");
-  if (!vradau5massstate.is_empty ())
+  if (!vradau5massstate.isempty ())
     if (vradau5massstate.string_value ().compare ("weak") != 0) // 'weak' is default
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -559,14 +507,14 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MassSingular will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -574,44 +522,44 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // The option MaxOrder will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (!vmaxder.is_empty ())
+  if (!vmaxder.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxOrder\" will be ignored by this solver");
 
   // The option BDF will be ignored by this solver, the core Fortran
   // solver doesn't support this option
   octave_value vbdf = vodeopt.contents ("BDF");
-  if (!vbdf.is_empty ())
+  if (!vbdf.isempty ())
     if (vbdf.string_value ().compare ("off") != 0)
      warning_with_id ("OdePkg:InvalidOption",
         "Option \"BDF\" will be ignored by this solver");
 
   // Implementation of the option NewtonTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vNTOL = vodeopt.contents ("NewtonTol");
-  if (vNTOL.is_empty ()) {
+  if (vNTOL.isempty ()) {
     vNTOL = 0;
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" not set, default value is used");
   }
 
   // Implementation of the option MaxNewtonIterations has been finished, this
   // option can be set by the user to another value than default value
  octave_value vmaxnit =
     vodeopt.contents ("MaxNewtonIterations");
-  if (vmaxnit.is_empty ()) {
-    vmaxnit = 7; warning_with_id ("OdePkg:InvalidOption",
+  if (vmaxnit.isempty ()) {
+    vmaxnit = 7; warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" not set, default value 7 is used");
   }
   else if (vmaxnit.int_value () < 1) {
-    vmaxnit = 7; warning_with_id ("OdePkg:InvalidOption",
+    vmaxnit = 7; warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" is zero, default value 7 is used");
   }
 
@@ -623,31 +571,31 @@
   NDArray vRTOL = vreltol.array_value ();
   NDArray vATOL = vabstol.array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
-  double X = tSpan(0);
+  octave_idx_type N = args(2).numel ();
+  double X = args(1).vector_value ()(0);
   double* Y = vY0.fortran_vec ();
-  double XEND = tSpan(1);
+  double XEND = args(1).vector_value ()(1);
   double H = vinitstep.double_value ();
   double *RTOL = vRTOL.fortran_vec ();
   double *ATOL = vATOL.fortran_vec ();
@@ -657,12 +605,12 @@
 
   octave_value vtim = args(1).vector_value ()(0);
   octave_value vsol = args(2);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vradau5pltfun.is_empty ()) odepkg_auxiliary_evalplotfun
-    (vradau5pltfun, vradau5outsel, octave_value(tSpan), args(2), vradau5extarg, 0);
-  if (!vradau5evefun.is_empty ())
+  if (!vradau5pltfun.isempty ()) odepkg_auxiliary_evalplotfun 
+    (vradau5pltfun, vradau5outsel, args(1), args(2), vradau5extarg, 0);
+  if (!vradau5evefun.isempty ())
     odepkg_auxiliary_evaleventfun (vradau5evefun, vtim, args(2), vradau5extarg, 0);
 
   // octave_stdout <<  "X VALUE=" << X << XEND << std::endl;
@@ -700,11 +648,11 @@
   octave_value vted = octave_value (XEND);
   octave_value vfin = octave_value (vlastline);
 
-  if (!vradau5pltfun.is_empty ()) odepkg_auxiliary_evalplotfun
+  if (!vradau5pltfun.isempty ()) odepkg_auxiliary_evalplotfun
     (vradau5pltfun, vradau5outsel, vted, vfin, vradau5extarg, 2);
-  if (!vradau5evefun.is_empty ()) odepkg_auxiliary_evaleventfun
+  if (!vradau5evefun.isempty ()) odepkg_auxiliary_evaleventfun
     (vradau5evefun, vted, vfin, vradau5extarg, 2);
-
+  
   // Get the stats information as an octave_scalar_map if the option 'Stats'
   // has been set with odeset
   // "nsteps", "nfailed", "nfevals", "npds", "ndecomps", "nlinsols"
@@ -725,18 +673,12 @@
     vretmap.assign ("solver", "ode5r");
-    if (!vstatinfo.is_empty ()) // Event implementation
+    if (!vstatinfo.isempty ()) // Event implementation
       vretmap.assign ("stats", vstatinfo);
-    if (!vradau5evefun.is_empty ()) {
+    if (!vradau5evefun.isempty ()) {
       vretmap.assign ("ie", vradau5evesol(0).cell_value ()(1));
       vretmap.assign ("xe", vradau5evesol(0).cell_value ()(2));
       vretmap.assign ("ye", vradau5evesol(0).cell_value ()(3));
@@ -744,38 +686,23 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vradau5evefun.is_empty ()) {
+    if (!vradau5evefun.isempty ()) {
       vretval(2) = vradau5evesol(0).cell_value ()(2);
       vretval(3) = vradau5evesol(0).cell_value ()(3);
       vretval(4) = vradau5evesol(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_rodas.cc
--- a/src/odepkg_octsolver_rodas.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_rodas.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -288,21 +288,21 @@
   // function that keeps the set of differential algebraic equations
   octave_value_list varin;
   varin(0) = X; varin(1) = A;
-  for (F77_INT vcnt = 0; vcnt < vrodasextarg.length (); vcnt++)
+  for (octave_idx_type vcnt = 0; vcnt < vrodasextarg.numel (); vcnt++)
     varin(vcnt+2) = vrodasextarg(vcnt);
-  octave_value_list vout = feval (vrodasodefun.function_value (), varin, 1);
+  octave_value_list vout = octave::feval (vrodasodefun.function_value (), varin, 1);
 
   // Return the results from the function evaluation to the Fortran
   // solver, again copy them and don't just create a Fortran vector
@@ -398,10 +398,10 @@
   octave_value vy = octave_value (A);
 
   // Check if an 'Events' function has been set by the user
-  if (!vrodasevefun.is_empty ()) {
+  if (!vrodasevefun.isempty ()) {
     vrodasevesol = odepkg_auxiliary_evaleventfun 
       (vrodasevefun, vt, vy, vrodasextarg, 1);
-    if (!vrodasevesol(0).cell_value ()(0).is_empty ())
+    if (!vrodasevesol(0).cell_value ()(0).isempty ())
       if (vrodasevesol(0).cell_value ()(0).int_value () == 1) {
         ColumnVector vttmp = vrodasevesol(0).cell_value ()(2).column_vector_value ();
         Matrix vrtmp = vrodasevesol(0).cell_value ()(3).matrix_value ();
@@ -417,15 +417,15 @@
 
   // Check if an 'OutputFcn' has been set by the user (including the
   // values of the options for 'OutputSel' and 'Refine')
-  if (!vrodaspltfun.is_empty ()) {
+  if (!vrodaspltfun.isempty ()) {
     if (vrodasrefine.int_value () > 0) {
       ColumnVector B(N); double vtb = 0.0;
      for (octave_idx_type vcnt = 1; vcnt < vrodasrefine.int_value (); vcnt++) {
@@ -491,10 +485,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;                      // The OdePkg options structure
 
@@ -536,8 +531,8 @@
   if (nargin >= 4) {
 
     // Fourth input argument != OdePkg option, need a default structure
-    if (!args(3).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(3).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();       // Create a default structure
       for (octave_idx_type vcnt = 3; vcnt < nargin; vcnt++)
         vrodasextarg(vcnt-3) = args(vcnt); // Save arguments in vrodasextarg
@@ -547,7 +542,7 @@
     else if (nargin > 4) {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "oders";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();       // Create structure from args(4)
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
@@ -558,7 +553,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "oders"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -566,7 +561,7 @@
   } // if (nargin >= 4)
 
   else { // if nargin == 3, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -578,14 +573,14 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"RelTol\" not set, new value %3.1e is used",
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -595,14 +590,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"AbsTol\" not set, new value %3.1e is used",
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -611,10 +606,10 @@
 
   // Setting the tolerance type that depends on the types (scalar or
   // vector) of the options RelTol and AbsTol
  octave_idx_type vitol = 0;
-  if (vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  if (vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 0;
-  else if (!vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  else if (!vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 1;
   else {
     error_with_id ("OdePkg:InvalidOption",
@@ -625,7 +620,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"NormControl\" will be ignored by this solver");
@@ -633,14 +628,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   vrodaspltfun = vodeopt.contents ("OutputFcn");
-  if (vrodaspltfun.is_empty () && nargout == 0) vrodaspltfun = "odeplot";
+  if (vrodaspltfun.isempty () && nargout == 0) vrodaspltfun = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -657,11 +652,11 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
+  if (args(1).numel () > 2) {
     error_with_id ("OdePkg:InvalidOption",
       "Fixed time stamps are not supported by this solver");
   }
-  if (vinitstep.is_empty ()) {
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -671,7 +666,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -686,20 +681,20 @@
   // Implementation of the option 'Jacobian' has been finished, these
   // options can be set by the user to another value than default
   vrodasjacfun = vodeopt.contents ("Jacobian");
  octave_idx_type vrodasjac = 0; // We need to set this if no Jac available
-  if (!vrodasjacfun.is_empty ()) vrodasjac = 1;
+  if (!vrodasjacfun.isempty ()) vrodasjac = 1;
 
   // The option JPattern will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vrodasjacpat = vodeopt.contents ("JPattern");
-  if (!vrodasjacpat.is_empty ())
+  if (!vrodasjacpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"JPattern\" will be ignored by this solver");
 
   // The option Vectorized will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vrodasvectorize = vodeopt.contents ("Vectorized");
-  if (!vrodasvectorize.is_empty ())
+  if (!vrodasvectorize.isempty ())
     if (vrodasvectorize.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"Vectorized\" will be ignored by this solver");
@@ -707,8 +702,8 @@
   // Implementation of the option 'Mass' has been finished, these
   // options can be set by the user to another value than default
   vrodasmass = vodeopt.contents ("Mass");
  octave_idx_type vrodasmas = 0;
-  if (!vrodasmass.is_empty ()) {
+  if (!vrodasmass.isempty ()) {
     vrodasmas = 1;
     if (vrodasmass.is_function_handle () || vrodasmass.is_inline_function ())
       warning_with_id ("OdePkg:InvalidOption",
@@ -718,7 +713,7 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   vrodasmassstate = vodeopt.contents ("MStateDependence");
-  if (!vrodasmassstate.is_empty ())
+  if (!vrodasmassstate.isempty ())
     if (vrodasmassstate.string_value ().compare ("weak") != 0) // 'weak' is default
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -726,14 +721,14 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MassSingular will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -741,21 +736,21 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // The option MaxOrder will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (!vmaxder.is_empty ())
+  if (!vmaxder.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxOrder\" will be ignored by this solver");
 
   // The option BDF will be ignored by this solver, the core Fortran
   // solver doesn't support this option
   octave_value vbdf = vodeopt.contents ("BDF");
-  if (!vbdf.is_empty ())
+  if (!vbdf.isempty ())
     if (vbdf.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"BDF\" will be ignored by this solver");
@@ -763,12 +758,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -780,33 +775,33 @@
   NDArray vRTOL = vreltol.array_value ();
   NDArray vATOL = vabstol.array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
+  octave_idx_type N = args(2).numel ();
  octave_idx_type IFCN = 1;
   double X = args(1).vector_value ()(0);
   double* Y = vY0.fortran_vec ();
@@ -815,9 +810,9 @@
   octave_value vtim = args(1).vector_value ()(0);
   octave_value vsol = args(2);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vrodaspltfun.is_empty ()) odepkg_auxiliary_evalplotfun 
+  if (!vrodaspltfun.isempty ()) odepkg_auxiliary_evalplotfun 
     (vrodaspltfun, vrodasoutsel, args(1), args(2), vrodasextarg, 0);
-  if (!vrodasevefun.is_empty ())
+  if (!vrodasevefun.isempty ())
     odepkg_auxiliary_evaleventfun (vrodasevefun, vtim, args(2), vrodasextarg, 0);
 
   // We are calling the core solver and solve the set of ODEs or DAEs
@@ -854,9 +849,9 @@
   octave_value vted = octave_value (XEND);
   octave_value vfin = octave_value (vlastline);
 
-  if (!vrodaspltfun.is_empty ()) odepkg_auxiliary_evalplotfun
+  if (!vrodaspltfun.isempty ()) odepkg_auxiliary_evalplotfun
     (vrodaspltfun, vrodasoutsel, vted, vfin, vrodasextarg, 2);
-  if (!vrodasevefun.is_empty ()) odepkg_auxiliary_evaleventfun
+  if (!vrodasevefun.isempty ()) odepkg_auxiliary_evaleventfun
     (vrodasevefun, vted, vfin, vrodasextarg, 2);
   
   // Get the stats information as an octave_scalar_map if the option 'Stats'
@@ -881,9 +876,9 @@
     vretmap.assign ("x", vtres);
     vretmap.assign ("y", vyres);
     vretmap.assign ("solver", "oders");
-    if (!vstatinfo.is_empty ()) // Event implementation
+    if (!vstatinfo.isempty ()) // Event implementation
       vretmap.assign ("stats", vstatinfo);
-    if (!vrodasevefun.is_empty ()) {
+    if (!vrodasevefun.isempty ()) {
       vretmap.assign ("ie", vrodasevesol(0).cell_value ()(1));
       vretmap.assign ("xe", vrodasevesol(0).cell_value ()(2));
       vretmap.assign ("ye", vrodasevesol(0).cell_value ()(3));
@@ -901,7 +896,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vrodasevefun.is_empty ()) {
+    if (!vrodasevefun.isempty ()) {
       vretval(2) = vrodasevesol(0).cell_value ()(2);
       vretval(3) = vrodasevesol(0).cell_value ()(3);
       vretval(4) = vrodasevesol(0).cell_value ()(1);
diff -r 20017c4feae6 src/odepkg_octsolver_seulex.cc
--- a/src/odepkg_octsolver_seulex.cc	Mon Jan 02 19:24:34 2017 -0500
+++ b/src/odepkg_octsolver_seulex.cc	Sat Dec 29 18:52:54 2018 +0100
@@ -112,29 +112,29 @@
   // function that keeps the set of differential algebraic equations
   octave_value_list varin;
   varin(0) = X; varin(1) = A;
-  for (F77_INT vcnt = 0; vcnt < vseulexextarg.length (); vcnt++)
+  for (octave_idx_type vcnt = 0; vcnt < vseulexextarg.numel (); vcnt++)
     varin(vcnt+2) = vseulexextarg(vcnt);
-  octave_value_list vout = feval (vseulexodefun.function_value (), varin, 1);
+  octave_value_list vout = octave::octave::feval (vseulexodefun.function_value (), varin, 1);
 
   // Return the results from the function evaluation to the Fortran
   // solver, again copy them and don't just create a Fortran vector
@@ -202,10 +202,10 @@
 
   vseulexevebrk = false;
   // Check if an 'Events' function has been set by the user
-  if (!vseulexevefun.is_empty ()) {
+  if (!vseulexevefun.isempty ()) {
     vseulexevesol = odepkg_auxiliary_evaleventfun 
       (vseulexevefun, vt, vy, vseulexextarg, 1);
-    if (!vseulexevesol(0).cell_value ()(0).is_empty ())
+    if (!vseulexevesol(0).cell_value ()(0).isempty ())
       if (vseulexevesol(0).cell_value ()(0).int_value () == 1) {
         ColumnVector vttmp = vseulexevesol(0).cell_value ()(2).column_vector_value ();
         Matrix vrtmp = vseulexevesol(0).cell_value ()(3).matrix_value ();
@@ -223,15 +223,15 @@
   // Check if an 'OutputFcn' has been set by the user (including the
   // values of the options for 'OutputSel' and 'Refine')
   vseulexpltbrk = false;
-  if (!vseulexpltfun.is_empty ()) {
+  if (!vseulexpltfun.isempty ()) {
     if (vseulexrefine.int_value () > 0) {
       ColumnVector B(N); double vtb = 0.0;
      for (octave_idx_type vcnt = 1; vcnt < vseulexrefine.int_value (); vcnt++) {
@@ -285,10 +279,11 @@
 
-  octave_idx_type nargin = args.length (); // The number of input arguments
+  octave_idx_type nargin = args.numel (); // The number of input arguments
   octave_value_list vretval;               // The cell array of return args
   octave_scalar_map vodeopt;                      // The OdePkg options structure
 
@@ -330,8 +325,8 @@
   if (nargin >= 4) {
 
     // Fourth input argument != OdePkg option, need a default structure
-    if (!args(3).is_map ()) {
-      octave_value_list tmp = feval ("odeset", tmp, 1);
+    if (!args(3).isstruct ()) {
+      octave_value_list tmp = octave::feval ("odeset", tmp, 1);
       vodeopt = tmp(0).scalar_map_value ();       // Create a default structure
       for (octave_idx_type vcnt = 3; vcnt < nargin; vcnt++)
         vseulexextarg(vcnt-3) = args(vcnt); // Save arguments in vseulexextarg
@@ -341,7 +336,7 @@
     else if (nargin > 4) {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "odesx";
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value ();       // Create structure from args(4)
       for (octave_idx_type vcnt = 4; vcnt < nargin; vcnt++)
@@ -352,7 +347,7 @@
     else {
       octave_value_list varin;
       varin(0) = args(3); varin(1) = "odesx"; // Check structure
-      octave_value_list tmp = feval ("odepkg_structure_check", varin, 1);
+      octave_value_list tmp = octave::feval ("odepkg_structure_check", varin, 1);
       if (error_state) return (vretval);
       vodeopt = tmp(0).scalar_map_value (); // Create a default structure
     }
@@ -360,7 +355,7 @@
   } // if (nargin >= 4)
 
   else { // if nargin == 3, everything else has been checked before
-    octave_value_list tmp = feval ("odeset", tmp, 1);
+    octave_value_list tmp = octave::feval ("odeset", tmp, 1);
     vodeopt = tmp(0).scalar_map_value (); // Create a default structure
   }
 
@@ -372,14 +367,14 @@
   // Implementation of the option RelTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vreltol = vodeopt.contents ("RelTol");
-  if (vreltol.is_empty ()) {
+  if (vreltol.isempty ()) {
     vreltol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"RelTol\" not set, new value %3.1e is used",
       vreltol.double_value ());
   } // vreltol.print (octave_stdout, true); return (vretval);
   if (!vreltol.is_scalar_type ()) {
-    if (vreltol.length () != args(2).length ()) {
+    if (vreltol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"RelTol\" must be the same as the number of equations");
       return (vretval);
@@ -389,14 +384,14 @@
   // Implementation of the option AbsTol has been finished, this
   // option can be set by the user to another value than default value
   octave_value vabstol = vodeopt.contents ("AbsTol");
-  if (vabstol.is_empty ()) {
+  if (vabstol.isempty ()) {
     vabstol = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"AbsTol\" not set, new value %3.1e is used",
         vabstol.double_value ());
   } // vabstol.print (octave_stdout, true); return (vretval);
   if (!vabstol.is_scalar_type ()) {
-    if (vabstol.length () != args(2).length ()) {
+    if (vabstol.numel () != args(2).numel ()) {
       error_with_id ("OdePkg:InvalidOption",
         "Length of option \"AbsTol\" must be the same as the number of equations");
       return (vretval);
@@ -405,10 +400,10 @@
 
   // Setting the tolerance type that depends on the types (scalar or
   // vector) of the options RelTol and AbsTol
  octave_idx_type vitol = 0;
-  if (vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  if (vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 0;
-  else if (!vreltol.is_scalar_type () && (vreltol.length () == vabstol.length ()))
+  else if (!vreltol.is_scalar_type () && (vreltol.numel () == vabstol.numel ()))
     vitol = 1;
   else {
     error_with_id ("OdePkg:InvalidOption",
@@ -419,7 +414,7 @@
   // The option NormControl will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnorm = vodeopt.contents ("NormControl");
-  if (!vnorm.is_empty ())
+  if (!vnorm.isempty ())
     if (vnorm.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"NormControl\" will be ignored by this solver");
@@ -427,14 +422,14 @@
   // The option NonNegative will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vnneg = vodeopt.contents ("NonNegative");
-  if (!vnneg.is_empty ())
+  if (!vnneg.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"NonNegative\" will be ignored by this solver");
 
   // Implementation of the option OutputFcn has been finished, this
   // option can be set by the user to another value than default value
   vseulexpltfun = vodeopt.contents ("OutputFcn");
-  if (vseulexpltfun.is_empty () && nargout == 0) vseulexpltfun = "odeplot";
+  if (vseulexpltfun.isempty () && nargout == 0) vseulexpltfun = "odeplot";
 
   // Implementation of the option OutputSel has been finished, this
   // option can be set by the user to another value than default value
@@ -451,11 +446,11 @@
   // Implementation of the option InitialStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vinitstep = vodeopt.contents ("InitialStep");
-  if (args(1).length () > 2) {
+  if (args(1).numel () > 2) {
     error_with_id ("OdePkg:InvalidOption",
       "Fixed time stamps are not supported by this solver");
   }
-  if (vinitstep.is_empty ()) {
+  if (vinitstep.isempty ()) {
     vinitstep = 1.0e-6;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialStep\" not set, new value %3.1e is used",
@@ -465,7 +460,7 @@
   // Implementation of the option MaxStep has been finished, this
   // option can be set by the user to another value than default value
   octave_value vmaxstep = vodeopt.contents ("MaxStep");
-  if (vmaxstep.is_empty () && args(1).length () == 2) {
+  if (vmaxstep.isempty () && args(1).numel () == 2) {
     vmaxstep = (args(1).vector_value ()(1) - args(1).vector_value ()(0)) / 12.5;
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxStep\" not set, new value %3.1e is used", 
@@ -480,20 +475,20 @@
   // Implementation of the option 'Jacobian' has been finished, these
   // options can be set by the user to another value than default
   vseulexjacfun = vodeopt.contents ("Jacobian");
   octave_idx_type vseulexjac = 0; // We need to set this if no Jac available
-  if (!vseulexjacfun.is_empty ()) vseulexjac = 1;
+  if (!vseulexjacfun.isempty ()) vseulexjac = 1;
 
   // The option JPattern will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vseulexjacpat = vodeopt.contents ("JPattern");
-  if (!vseulexjacpat.is_empty ())
+  if (!vseulexjacpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"JPattern\" will be ignored by this solver");
 
   // The option Vectorized will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vseulexvectorize = vodeopt.contents ("Vectorized");
-  if (!vseulexvectorize.is_empty ())
+  if (!vseulexvectorize.isempty ())
     if (vseulexvectorize.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"Vectorized\" will be ignored by this solver");
@@ -501,8 +496,8 @@
   // Implementation of the option 'Mass' has been finished, these
   // options can be set by the user to another value than default
   vseulexmass = vodeopt.contents ("Mass");
  octave_idx_type vseulexmas = 0;
-  if (!vseulexmass.is_empty ()) {
+  if (!vseulexmass.isempty ()) {
     vseulexmas = 1;
     if (vseulexmass.is_function_handle () || vseulexmass.is_inline_function ())
       warning_with_id ("OdePkg:InvalidOption",
@@ -512,7 +507,7 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   vseulexmassstate = vodeopt.contents ("MStateDependence");
-  if (!vseulexmassstate.is_empty ())
+  if (!vseulexmassstate.isempty ())
     if (vseulexmassstate.string_value ().compare ("weak") != 0) // 'weak' is default
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MStateDependence\" will be ignored by this solver");
@@ -520,14 +515,14 @@
   // The option MStateDependence will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmvpat = vodeopt.contents ("MvPattern");
-  if (!vmvpat.is_empty ())
+  if (!vmvpat.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MvPattern\" will be ignored by this solver");
 
   // The option MassSingular will be ignored by this solver, the
   // core Fortran solver doesn't support this option
   octave_value vmsing = vodeopt.contents ("MassSingular");
-  if (!vmsing.is_empty ())
+  if (!vmsing.isempty ())
     if (vmsing.string_value ().compare ("maybe") != 0)
       warning_with_id ("OdePkg:InvalidOption",
         "Option \"MassSingular\" will be ignored by this solver");
@@ -535,21 +530,21 @@
   // The option InitialSlope will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vinitslope = vodeopt.contents ("InitialSlope");
-  if (!vinitslope.is_empty ())
+  if (!vinitslope.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"InitialSlope\" will be ignored by this solver");
 
   // The option MaxOrder will be ignored by this solver, the core
   // Fortran solver doesn't support this option
   octave_value vmaxder = vodeopt.contents ("MaxOrder");
-  if (!vmaxder.is_empty ())
+  if (!vmaxder.isempty ())
     warning_with_id ("OdePkg:InvalidOption",
       "Option \"MaxOrder\" will be ignored by this solver");
 
   // The option BDF will be ignored by this solver, the core Fortran
   // solver doesn't support this option
   octave_value vbdf = vodeopt.contents ("BDF");
-  if (!vbdf.is_empty ())
+  if (!vbdf.isempty ())
     if (vbdf.string_value ().compare ("off") != 0)
       warning_with_id ("OdePkg:InvalidOption", 
         "Option \"BDF\" will be ignored by this solver");
@@ -557,12 +552,12 @@
   // this solver, IT NEEDS TO BE CHECKED IF THE FORTRAN CORE SOLVER
   // CAN HANDLE THESE OPTIONS
   octave_value vntol = vodeopt.contents ("NewtonTol");
-  if (!vntol.is_empty ())
+  if (!vntol.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"NewtonTol\" will be ignored by this solver");
   octave_value vmaxnewton = 
     vodeopt.contents ("MaxNewtonIterations");
-  if (!vmaxnewton.is_empty ())
+  if (!vmaxnewton.isempty ())
     warning_with_id ("OdePkg:InvalidOption", 
       "Option \"MaxNewtonIterations\" will be ignored by this solver");
 
@@ -574,32 +569,32 @@
   NDArray vRTOL = vreltol.array_value ();
   NDArray vATOL = vabstol.array_value ();
 
-  F77_INT N = TO_F77_INT (args(2).length ());
-  F77_INT IFCN = 1;
+  octave_idx_type N = args(2).numel ();
+  octave_idx_type IFCN = 1;
   double X = args(1).vector_value ()(0);
   double* Y = vY0.fortran_vec ();
   double XEND = args(1).vector_value ()(1);
@@ -610,9 +605,9 @@
   octave_value vtim = args(1).vector_value ()(0);
   octave_value vsol = args(2);
   odepkg_auxiliary_solstore (vtim, vsol, 0);
-  if (!vseulexpltfun.is_empty ()) odepkg_auxiliary_evalplotfun 
+  if (!vseulexpltfun.isempty ()) odepkg_auxiliary_evalplotfun 
     (vseulexpltfun, vseulexoutsel, args(1), args(2), vseulexextarg, 0);
-  if (!vseulexevefun.is_empty ())
+  if (!vseulexevefun.isempty ())
     odepkg_auxiliary_evaleventfun (vseulexevefun, vtim, args(2), vseulexextarg, 0);
 
   // We are calling the core solver and solve the set of ODEs or DAEs
@@ -653,9 +648,9 @@
   octave_value vted = octave_value (XEND);
   octave_value vfin = octave_value (vlastline);
 
-  if (!vseulexpltfun.is_empty ()) odepkg_auxiliary_evalplotfun
+  if (!vseulexpltfun.isempty ()) odepkg_auxiliary_evalplotfun
     (vseulexpltfun, vseulexoutsel, vted, vfin, vseulexextarg, 2);
-  if (!vseulexevefun.is_empty ()) odepkg_auxiliary_evaleventfun
+  if (!vseulexevefun.isempty ()) odepkg_auxiliary_evaleventfun
     (vseulexevefun, vted, vfin, vseulexextarg, 2);
   
   // Get the stats information as an octave_scalar_map if the option 'Stats'
@@ -680,9 +675,9 @@
     vretmap.assign ("x", vtres);
     vretmap.assign ("y", vyres);
     vretmap.assign ("solver", "odesx");
-    if (!vstatinfo.is_empty ()) // Event implementation
+    if (!vstatinfo.isempty ()) // Event implementation
       vretmap.assign ("stats", vstatinfo);
-    if (!vseulexevefun.is_empty ()) {
+    if (!vseulexevefun.isempty ()) {
       vretmap.assign ("ie", vseulexevesol(0).cell_value ()(1));
       vretmap.assign ("xe", vseulexevesol(0).cell_value ()(2));
       vretmap.assign ("ye", vseulexevesol(0).cell_value ()(3));
@@ -700,7 +695,7 @@
     vretval(2) = vempty;
     vretval(3) = vempty;
     vretval(4) = vempty;
-    if (!vseulexevefun.is_empty ()) {
+    if (!vseulexevefun.isempty ()) {
       vretval(2) = vseulexevesol(0).cell_value ()(2);
       vretval(3) = vseulexevesol(0).cell_value ()(3);
       vretval(4) = vseulexevesol(0).cell_value ()(1);
