Added public abstract f, fQ, g, and Jac for external calls. Use fTestMode, fQTestMod...
authorilb@NIH.GOV <ilb@NIH.GOV@ba61647d-9d00-f842-95cd-605cb4296b96>
Fri, 9 Mar 2018 17:19:04 +0000 (17:19 +0000)
committerilb@NIH.GOV <ilb@NIH.GOV@ba61647d-9d00-f842-95cd-605cb4296b96>
Fri, 9 Mar 2018 17:19:04 +0000 (17:19 +0000)
git-svn-id: https://citdcbmipav.cit.nih.gov/repos-pub/mipav/trunk@15409 ba61647d-9d00-f842-95cd-605cb4296b96

mipav/src/gov/nih/mipav/model/algorithms/CVODES.java

index e971c19..a42a086 100644 (file)
@@ -558,6 +558,7 @@ public abstract class CVODES {
 \r
     final int cvsRoberts_dns = 1;\r
     int problem = cvsRoberts_dns;\r
+    boolean testMode = true;\r
        \r
        \r
        public CVODES() {\r
@@ -935,10 +936,10 @@ public abstract class CVODES {
         * @param user_data\r
         * @return\r
         */\r
-       private int f(double t, NVector yv, NVector ydotv, CVodeMemRec user_data) {\r
+       private int fTestMode(double t, NVector yv, NVector ydotv, CVodeMemRec user_data) {\r
                double y[] = yv.getData();\r
                double ydot[] = ydotv.getData();\r
-               if (problem == 1) {\r
+               if (problem == cvsRoberts_dns) {\r
                    ydot[0] = -0.04*y[0] + 1.0E4*y[1]*y[2];\r
                    ydot[2] = 3.0E7*y[1]*y[1];\r
                    ydot[1] = -ydot[0] - ydot[2];\r
@@ -946,10 +947,14 @@ public abstract class CVODES {
                return 0;\r
        }\r
        \r
-       private int fQ(double t, NVector x, NVector y, CVodeMemRec user_data) {\r
+       public abstract int f(double t, NVector yv, NVector ydotv, CVodeMemRec user_data);\r
+       \r
+       private int fQTestMode(double t, NVector x, NVector y, CVodeMemRec user_data) {\r
                return 0;\r
        }\r
        \r
+       public abstract int fQ(double t, NVector x, NVector y, CVodeMemRec user_data);\r
+       \r
        /**\r
         * g routine.  Compute functions g_i(t,y) for i = 0,1.\r
         * @param t\r
@@ -958,15 +963,17 @@ public abstract class CVODES {
         * @param user_data\r
         * @return\r
         */\r
-       private int g(double t, NVector yv, double gout[], CVodeMemRec user_data) {\r
+       private int gTestMode(double t, NVector yv, double gout[], CVodeMemRec user_data) {\r
                double y[] = yv.getData();\r
-               if (problem == 1) {\r
+               if (problem == cvsRoberts_dns) {\r
                    gout[0] = y[0] - 0.0001;\r
                    gout[1] = y[2] - 0.01;\r
                }\r
                return 0;\r
        }\r
        \r
+       public abstract int g(double t, NVector yv, double gout[], CVodeMemRec user_data);\r
+       \r
        /**\r
         * Jacobian routine.  Compute J(t,y) = df/dy.\r
         * @param t\r
@@ -977,24 +984,26 @@ public abstract class CVODES {
         * @param tmp2\r
         * @return\r
         */\r
-       private int Jac(double t, NVector yv, NVector fy, double J[][], CVodeMemRec data, NVector tmp1, NVector tmp2, NVector tmp3) {\r
+       private int JacTestMode(double t, NVector yv, NVector fy, double J[][], CVodeMemRec data, NVector tmp1, NVector tmp2, NVector tmp3) {\r
                double y[] = yv.getData();\r
-           J[0][0] = -0.04;\r
-           J[0][1] = 1.0E4 * y[2];\r
-           J[0][2] = 1.0E4 * y[1];\r
-           \r
-           J[1][0] = 0.04;\r
-           J[1][1] = -1.0E4 * y[2] - 6.0E7*y[1];\r
-           J[1][2] = -1.0E4 * y[1];\r
-           \r
-           J[2][0] = ZERO;\r
-           J[2][1] = 6.0E7 * y[1];\r
-           J[2][2] = ZERO;\r
+               if (problem == cvsRoberts_dns) {\r
+                   J[0][0] = -0.04;\r
+                   J[0][1] = 1.0E4 * y[2];\r
+                   J[0][2] = 1.0E4 * y[1];\r
+                   \r
+                   J[1][0] = 0.04;\r
+                   J[1][1] = -1.0E4 * y[2] - 6.0E7*y[1];\r
+                   J[1][2] = -1.0E4 * y[1];\r
+                   \r
+                   J[2][0] = ZERO;\r
+                   J[2][1] = 6.0E7 * y[1];\r
+                   J[2][2] = ZERO;\r
+               }\r
            \r
            return 0;\r
        }\r
        \r
-       \r
+       public abstract int Jac(double t, NVector yv, NVector fy, double J[][], CVodeMemRec data, NVector tmp1, NVector tmp2, NVector tmp3);\r
        \r
        \r
        \r
@@ -2391,8 +2400,14 @@ public abstract class CVODES {
        * If computing quadr. sensi., call fQS at (t0,y0,yS0), set znQS[1][is] = yQS'(t0), is=1,...,Ns.\r
        */\r
 \r
-      retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+      if (testMode) {\r
+         retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+                  cv_mem.cv_zn[1], cv_mem.cv_user_data);   \r
+      }\r
+      else {\r
+          retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                             cv_mem.cv_zn[1], cv_mem.cv_user_data); \r
+      }\r
       cv_mem.cv_nfe++;\r
       if (retval < 0) {\r
         cvProcessError(cv_mem, CV_RHSFUNC_FAIL, "CVODES", "CVode",\r
@@ -2406,8 +2421,14 @@ public abstract class CVODES {
       }\r
 \r
       if (cv_mem.cv_quadr) {\r
-        retval = fQ(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+       if (testMode) {\r
+               retval = fQTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+                    cv_mem.cv_znQ[1], cv_mem.cv_user_data);    \r
+       }\r
+       else {\r
+            retval = fQ(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                                cv_mem.cv_znQ[1], cv_mem.cv_user_data);\r
+       }\r
         cv_mem.cv_nfQe++;\r
         if (retval < 0) {\r
           cvProcessError(cv_mem, CV_QRHSFUNC_FAIL, "CVODES", "CVode",\r
@@ -4228,15 +4249,27 @@ public abstract class CVODES {
    \r
    /* tempv <- f(t+h, h*y'(t)+y(t)) */\r
 \r
-   retval = f(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
+   if (testMode) {\r
+          retval = fTestMode(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
+               cv_mem.cv_tempv, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = f(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
                          cv_mem.cv_tempv, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nfe++;\r
    if (retval < 0) return(CV_RHSFUNC_FAIL);\r
    if (retval > 0) return(RHSFUNC_RECVR);\r
 \r
    if (cv_mem.cv_quadr && cv_mem.cv_errconQ) {\r
-     retval = fQ(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
+        if (testMode) {\r
+                retval = fQTestMode(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
+                 cv_mem.cv_tempvQ, cv_mem.cv_user_data);        \r
+        }\r
+        else {\r
+         retval = fQ(cv_mem.cv_tn+hg, cv_mem.cv_y,\r
                             cv_mem.cv_tempvQ, cv_mem.cv_user_data);\r
+        }\r
      cv_mem.cv_nfQe++;\r
      if (retval < 0) return(CV_QRHSFUNC_FAIL);\r
      if (retval > 0) return(QRHSFUNC_RECVR);\r
@@ -4513,8 +4546,14 @@ else                return(snrm);
      cv_mem.cv_uround*HUNDRED;\r
 \r
    /* Evaluate g at initial t and check for zero values. */\r
-   retval = g(cv_mem.cv_tlo, cv_mem.cv_zn[0],\r
+   if (testMode) {\r
+          retval = gTestMode(cv_mem.cv_tlo, cv_mem.cv_zn[0],\r
+               cv_mem.cv_glo, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = g(cv_mem.cv_tlo, cv_mem.cv_zn[0],\r
                             cv_mem.cv_glo, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nge = 1;\r
    if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -4532,8 +4571,14 @@ else                return(snrm);
    smallh = hratio*cv_mem.cv_h;\r
    tplus = cv_mem.cv_tlo + smallh;\r
    N_VLinearSum_Serial(ONE, cv_mem.cv_zn[0], hratio, cv_mem.cv_zn[1], cv_mem.cv_y);\r
-   retval = g(tplus, cv_mem.cv_y,\r
+   if (testMode) {\r
+          retval = gTestMode(tplus, cv_mem.cv_y,\r
+               cv_mem.cv_ghi, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = g(tplus, cv_mem.cv_y,\r
                             cv_mem.cv_ghi, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nge++;\r
    if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -4578,8 +4623,14 @@ else                return(snrm);
    if (cv_mem.cv_irfnd == 0) return(CV_SUCCESS);\r
 \r
    CVodeGetDky(cv_mem, cv_mem.cv_tlo, 0, cv_mem.cv_y);\r
-   retval = g(cv_mem.cv_tlo, cv_mem.cv_y,\r
+   if (testMode) {\r
+          retval = gTestMode(cv_mem.cv_tlo, cv_mem.cv_y,\r
+               cv_mem.cv_glo, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = g(cv_mem.cv_tlo, cv_mem.cv_y,\r
                             cv_mem.cv_glo, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nge++;\r
    if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -4606,8 +4657,14 @@ else                return(snrm);
    } else {\r
      CVodeGetDky(cv_mem, tplus, 0, cv_mem.cv_y);\r
    }\r
-   retval = g(tplus, cv_mem.cv_y,\r
+   if (testMode) {\r
+          retval = gTestMode(tplus, cv_mem.cv_y,\r
+               cv_mem.cv_ghi, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = g(tplus, cv_mem.cv_y,\r
                             cv_mem.cv_ghi, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nge++;\r
    if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -4731,8 +4788,14 @@ else                return(snrm);
    }\r
 \r
    /* Set ghi = g(thi) and call cvRootfind to search (tlo,thi) for roots. */\r
-   retval = g(cv_mem.cv_thi, cv_mem.cv_y,\r
+   if (testMode) {\r
+          retval = gTestMode(cv_mem.cv_thi, cv_mem.cv_y,\r
+               cv_mem.cv_ghi, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = g(cv_mem.cv_thi, cv_mem.cv_y,\r
                             cv_mem.cv_ghi, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nge++;\r
    if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -4925,8 +4988,14 @@ else                return(snrm);
      }\r
 \r
      CVodeGetDky(cv_mem, tmid, 0, cv_mem.cv_y);\r
-     retval = g(tmid, cv_mem.cv_y, cv_mem.cv_grout,\r
+     if (testMode) {\r
+        retval = gTestMode(tmid, cv_mem.cv_y, cv_mem.cv_grout,\r
+                 cv_mem.cv_user_data);  \r
+     }\r
+     else {\r
+         retval = g(tmid, cv_mem.cv_y, cv_mem.cv_grout,\r
                               cv_mem.cv_user_data);\r
+     }\r
      cv_mem.cv_nge++;\r
      if (retval != 0) return(CV_RTFUNC_FAIL);\r
 \r
@@ -5134,8 +5203,14 @@ else                return(snrm);
         * If f() fails recoverably, treat it as a convergence failure and \r
         * attempt the step again */\r
 \r
-      retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
+       if (testMode) {\r
+          retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_y,\r
+                   cv_mem.cv_ftemp, cv_mem.cv_user_data);   \r
+       }\r
+       else {\r
+           retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
                              cv_mem.cv_ftemp, cv_mem.cv_user_data);\r
+       }\r
        cv_mem.cv_nfe++;\r
        if (retval < 0) return(CV_RHSFUNC_FAIL);\r
        if (retval > 0) {\r
@@ -5931,8 +6006,14 @@ else                return(snrm);
    /* Initialize delS and Delp to avoid compiler warning message */\r
    delS = Delp = ZERO;\r
 \r
-   retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+   if (testMode) {\r
+          retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+               cv_mem.cv_tempv, cv_mem.cv_user_data);   \r
+   }\r
+   else {\r
+       retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                          cv_mem.cv_tempv, cv_mem.cv_user_data);\r
+   }\r
    cv_mem.cv_nfe++;\r
    if (retval < 0) return(CV_RHSFUNC_FAIL);\r
    if (retval > 0) return(RHSFUNC_RECVR);\r
@@ -6033,8 +6114,15 @@ else                return(snrm);
 \r
      Delp = Del;\r
 \r
-     retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
+     if (testMode) {\r
+        retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_y,\r
+                 cv_mem.cv_tempv, cv_mem.cv_user_data);\r
+        \r
+     }\r
+     else {\r
+         retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
                            cv_mem.cv_tempv, cv_mem.cv_user_data);\r
+     }\r
      cv_mem.cv_nfe++;\r
      if (retval < 0) return(CV_RHSFUNC_FAIL);\r
      if (retval > 0) return(RHSFUNC_RECVR);\r
@@ -6128,8 +6216,14 @@ else                return(snrm);
 \r
    for(;;) {\r
 \r
-     retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+     if (testMode) {\r
+        retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+                 cv_mem.cv_ftemp, cv_mem.cv_user_data);         \r
+     }\r
+     else {\r
+            retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                            cv_mem.cv_ftemp, cv_mem.cv_user_data);\r
+     }\r
      cv_mem.cv_nfe++; \r
      if (retval < 0) return(CV_RHSFUNC_FAIL);\r
      if (retval > 0) return(RHSFUNC_RECVR);\r
@@ -6262,8 +6356,14 @@ else                return(snrm);
          return(-1);\r
        }\r
 \r
-       retval = Jac(cv_mem.cv_tn, y, fy, cvdls_mem.A, \r
+       if (testMode) {\r
+          retval = JacTestMode(cv_mem.cv_tn, y, fy, cvdls_mem.A, \r
+                   cvdls_mem.J_data, tmp1, tmp2, tmp3);   \r
+       }\r
+       else {\r
+           retval = Jac(cv_mem.cv_tn, y, fy, cvdls_mem.A, \r
                                cvdls_mem.J_data, tmp1, tmp2, tmp3);\r
+       }\r
        if (retval < 0) {\r
          cvProcessError(cv_mem, CVDLS_JACFUNC_UNRECVR, "CVSDLS", \r
                          "cvDlsSetup",  MSGD_JACFUNC_FAILED);\r
@@ -6605,8 +6705,14 @@ else                return(snrm);
        \r
        /* Save norm of correction, evaluate f, and loop again */\r
        Delp = Del;\r
-       retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
+       if (testMode) {\r
+          retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_y,\r
+                   cv_mem.cv_ftemp, cv_mem.cv_user_data);   \r
+       }\r
+       else {\r
+           retval = f(cv_mem.cv_tn, cv_mem.cv_y,\r
                              cv_mem.cv_ftemp, cv_mem.cv_user_data);\r
+       }\r
        cv_mem.cv_nfe++;\r
        if (retval < 0) return(CV_RHSFUNC_FAIL);\r
        if (retval > 0) {\r
@@ -6993,8 +7099,14 @@ else                return(snrm);
           cv_mem.cv_qwait = LONG_WAIT;\r
           cv_mem.cv_nscon = 0;\r
 \r
-          retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+          if (testMode) {\r
+                  retval = fTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+                   cv_mem.cv_tempv, cv_mem.cv_user_data);   \r
+          }\r
+          else {\r
+              retval = f(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                                 cv_mem.cv_tempv, cv_mem.cv_user_data);\r
+          }\r
           cv_mem.cv_nfe++;\r
           if (retval < 0) return(CV_RHSFUNC_FAIL);\r
           if (retval > 0) return(CV_UNREC_RHSFUNC_ERR);\r
@@ -7003,8 +7115,14 @@ else                return(snrm);
 \r
           if (cv_mem.cv_quadr) {\r
 \r
-            retval = fQ(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+            if (testMode) {\r
+                retval = fQTestMode(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
+                     cv_mem.cv_tempvQ, cv_mem.cv_user_data);    \r
+            }\r
+            else {\r
+                    retval = fQ(cv_mem.cv_tn, cv_mem.cv_zn[0],\r
                                    cv_mem.cv_tempvQ, cv_mem.cv_user_data);\r
+            }\r
             cv_mem.cv_nfQe++;\r
             if (retval < 0) return(CV_QRHSFUNC_FAIL);\r
             if (retval > 0) return(CV_UNREC_QRHSFUNC_ERR);\r
@@ -7072,8 +7190,14 @@ else                return(snrm);
           int retval;\r
 \r
           /* Save quadrature correction in acorQ */\r
-          retval = fQ(cv_mem.cv_tn, cv_mem.cv_y,\r
+          if (testMode) {\r
+                  retval = fQTestMode(cv_mem.cv_tn, cv_mem.cv_y,\r
+                   cv_mem.cv_acorQ, cv_mem.cv_user_data);   \r
+          }\r
+          else {\r
+              retval = fQ(cv_mem.cv_tn, cv_mem.cv_y,\r
                                  cv_mem.cv_acorQ, cv_mem.cv_user_data);\r
+          }\r
           cv_mem.cv_nfQe++;\r
           if (retval < 0) return(CV_QRHSFUNC_FAIL);\r
           if (retval > 0) return(QRHSFUNC_RECVR);\r
@@ -7831,14 +7955,24 @@ else                return(snrm);
             N_VLinearSum_Serial(ONE,y,Delta,yS,ytemp);\r
             cv_mem.cv_p[which] = psave + Delta;\r
 \r
-            retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ySdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
             \r
             N_VLinearSum_Serial(ONE,y,-Delta,yS,ytemp);\r
             cv_mem.cv_p[which] = psave - Delta;\r
 \r
-            retval = f(t, ytemp, ftemp, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ftemp, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ftemp, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
@@ -7853,25 +7987,45 @@ else                return(snrm);
             \r
             N_VLinearSum_Serial(ONE,y,Deltay,yS,ytemp);\r
 \r
-            retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ySdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
             N_VLinearSum_Serial(ONE,y,-Deltay,yS,ytemp);\r
 \r
-            retval = f(t, ytemp, ftemp, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ftemp, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ftemp, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
             N_VLinearSum_Serial(r2Deltay, ySdot, -r2Deltay, ftemp, ySdot);\r
             \r
             cv_mem.cv_p[which] = psave + Deltap;\r
-            retval = f(t, y, ytemp, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, y, ytemp, cv_mem.cv_user_data);   \r
+            }\r
+            else {\r
+                retval = f(t, y, ytemp, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
             cv_mem.cv_p[which] = psave - Deltap;\r
-            retval = f(t, y, ftemp, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, y, ftemp, cv_mem.cv_user_data);   \r
+            }\r
+            else {\r
+                retval = f(t, y, ftemp, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
@@ -7889,7 +8043,12 @@ else                return(snrm);
             N_VLinearSum_Serial(ONE,y,Delta,yS,ytemp);\r
             cv_mem.cv_p[which] = psave + Delta;\r
 \r
-            retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ySdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
             \r
@@ -7901,14 +8060,24 @@ else                return(snrm);
             \r
             N_VLinearSum_Serial(ONE,y,Deltay,yS,ytemp);\r
 \r
-            retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, ytemp, ySdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = f(t, ytemp, ySdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
             N_VLinearSum_Serial(rDeltay, ySdot, -rDeltay, ydot, ySdot);\r
             \r
             cv_mem.cv_p[which] = psave + Deltap;\r
-            retval = f(t, y, ytemp, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fTestMode(t, y, ytemp, cv_mem.cv_user_data);   \r
+            }\r
+            else {\r
+                retval = f(t, y, ytemp, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
@@ -7997,14 +8166,24 @@ else                return(snrm);
             N_VLinearSum_Serial(ONE, y, Delta, yS, tmp);\r
             cv_mem.cv_p[which] = psave + Delta;\r
 \r
-            retval = fQ(t, tmp, yQSdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fQTestMode(t, tmp, yQSdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = fQ(t, tmp, yQSdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
             \r
             N_VLinearSum_Serial(ONE, y, -Delta, yS, tmp);\r
             cv_mem.cv_p[which] = psave - Delta;\r
 \r
-            retval = fQ(t, tmp, tmpQ, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fQTestMode(t, tmp, tmpQ, cv_mem.cv_user_data);         \r
+            }\r
+            else {\r
+                retval = fQ(t, tmp, tmpQ, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
 \r
@@ -8020,7 +8199,12 @@ else                return(snrm);
             N_VLinearSum_Serial(ONE, y, Delta, yS, tmp);\r
             cv_mem.cv_p[which] = psave + Delta;\r
 \r
-            retval = fQ(t, tmp, yQSdot, cv_mem.cv_user_data);\r
+            if (testMode) {\r
+                retval = fQTestMode(t, tmp, yQSdot, cv_mem.cv_user_data);       \r
+            }\r
+            else {\r
+                retval = fQ(t, tmp, yQSdot, cv_mem.cv_user_data);\r
+            }\r
             nfel++;\r
             if (retval != 0) return(retval);\r
             \r