Porting continues.
authorilb@NIH.GOV <ilb@NIH.GOV@ba61647d-9d00-f842-95cd-605cb4296b96>
Wed, 24 Jan 2018 23:22:45 +0000 (23:22 +0000)
committerilb@NIH.GOV <ilb@NIH.GOV@ba61647d-9d00-f842-95cd-605cb4296b96>
Wed, 24 Jan 2018 23:22:45 +0000 (23:22 +0000)
git-svn-id: https://citdcbmipav.cit.nih.gov/repos-pub/mipav/trunk@15350 ba61647d-9d00-f842-95cd-605cb4296b96

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

index 124750e..99ef653 100644 (file)
@@ -282,6 +282,10 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
     private double ZQUCA[] = new double[2]; // At first location of ZQUCA\r
     private double WPPQG[][] = new double[MQUCA][2];\r
     private double ZPPQG[][] = new double[MQUCA][2];\r
+    // COMMON /DSDTDA/PT,MD,HL\r
+    private int PT;\r
+    private double MD;\r
+    private double HL;\r
        public SymmsIntegralMapping() {\r
                \r
        }\r
@@ -12666,13 +12670,25 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
             \r
        //     LOCAL VARAIBLES\r
        \r
-       final int NIXINT = 200;\r
+       final int NXINT = 200;\r
        final int MAXSA = 100;\r
-       int I,IMNLA,J,L,MNSUA,NASYM,NCRVS,\r
-            NINFD,NPRVS,NXINT,NZERD,TNSUA;\r
+       int NASYM[] = new int[1];\r
+       int NCRVS[] = new int[1];\r
+       int NINFD[] = new int[1];\r
+       int NPRVS[] = new int[1];\r
+       int NZERD[] = new int[1];\r
+       int IMNLA = 0;\r
+       int I,J,L,MNSUA,TNSUA;\r
        final double BIG = 4.4;\r
-       double ANGSP,CCAPH,COCAP,COPHC,CPHCA,CR,EXCAP,EXPHC,LA,\r
-            OFLOW,PI,R1MACH,TOTLN,MCHEP;\r
+       double CCAPH[] = new double[1];\r
+       double COCAP[] = new double[1];\r
+       double COPHC[] = new double[1];\r
+       double CPHCA[] = new double[1];\r
+       double EXCAP[] = new double[1];\r
+       double EXPHC[] = new double[1];\r
+       double TOTLN[] = new double[1];\r
+       double ANGSP,CR,LA,\r
+            OFLOW,PI,R1MACH,MCHEP;\r
        String CHPC, CHCP;\r
        //CHARACTER OFLC*6,OFP0*6,OFP1*6,JBNM*4,CHPC*2,CHCP*2\r
        \r
@@ -12735,7 +12751,7 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
                LODP[I] = QUPTS[NARCS*NQPTS + I];\r
                LODW[I] = QUWTS[NARCS*NQPTS + I];\r
        }\r
-       /*DIAGN4(CCAPH,COCAP,COPHC,CPHCA,EXCAP,EXPHC,ICRVS,IER,\r
+       DIAGN4(CCAPH,COCAP,COPHC,CPHCA,EXCAP,EXPHC,ICRVS,IER,\r
             IPRVS,NASYM,NCRVS,NINFD,NPRVS,NZERD,ARCLN,ASYMP,BCDMN,CORXX,\r
             TOTLN,VTARG,MAP11,DGPOL,JATYP,LOSUB,\r
             NARCS,NQPTS,NXINT,PARNT,TNSUA,AICOF,\r
@@ -12777,8 +12793,8 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
        \r
        OFLOW=Double.MAX_VALUE;\r
        MCHEP=EPS;\r
-       UPHYC[0]=MCHEP*CPHCA;\r
-       UCANP[0]=MCHEP*CCAPH;\r
+       UPHYC[0]=MCHEP*CPHCA[0];\r
+       UCANP[0]=MCHEP*CCAPH[0];\r
        System.out.println();\r
        Preferences.debug("\n", Preferences.DEBUG_ALGORITHM);\r
        System.out.println("PHYSICAL ROUNDOFF MAGNIFIES TO: " + UPHYC);\r
@@ -12796,33 +12812,33 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
        Preferences.debug("\n", Preferences.DEBUG_ALGORITHM);\r
        Preferences.debug("   MAP           ESTIMATED EVALUATION   ESTIMATED MAXIMUM \n", Preferences.DEBUG_ALGORITHM);\r
        \r
-       if (NINFD > 0) {\r
+       if (NINFD[0] > 0) {\r
            CHPC="**";\r
        }\r
        else {\r
             CHPC="  ";\r
        }\r
-       if (NZERD > 0) {\r
+       if (NZERD[0] > 0) {\r
            CHCP="**";\r
        }\r
        else {\r
            CHCP="  ";\r
        }\r
        \r
-       Preferences.debug("PHY --> CAN        " + CPHCA + CHPC +"          " + UPHYC + "\n", Preferences.DEBUG_ALGORITHM);\r
-       Preferences.debug("CAN --> PHY        " + CCAPH + CHCP + "          " + UCANP + "\n", Preferences.DEBUG_ALGORITHM);\r
+       Preferences.debug("PHY --> CAN        " + CPHCA[0] + CHPC +"          " + UPHYC + "\n", Preferences.DEBUG_ALGORITHM);\r
+       Preferences.debug("CAN --> PHY        " + CCAPH[0] + CHCP + "          " + UCANP + "\n", Preferences.DEBUG_ALGORITHM);\r
        \r
              \r
        Preferences.debug("* BASED ON UNIT ROUNDOFF IN DOMAIN OF MAP\n",Preferences.DEBUG_ALGORITHM);\r
-       if (NINFD > 0 || NZERD > 0) {\r
+       if (NINFD[0] > 0 || NZERD[0] > 0) {\r
            Preferences.debug("** CONDITION NUMBER DEPENDS ON UNIT ROUNDOFF,U" + "\n", Preferences.DEBUG_ALGORITHM);\r
-           if (NINFD > 0) {\r
-               Preferences.debug("   PHY --> CAN : CONDTN NO = " + COPHC + "*U**" + EXPHC + "\n",Preferences.DEBUG_ALGORITHM );\r
+           if (NINFD[0] > 0) {\r
+               Preferences.debug("   PHY --> CAN : CONDTN NO = " + COPHC[0] + "*U**" + EXPHC[0] + "\n",Preferences.DEBUG_ALGORITHM );\r
            }\r
-           if (NZERD > 0) {\r
-               Preferences.debug("   CAN --> PHY : CONDTN NO = " + COCAP + "*U**" + EXCAP + "\n", Preferences.DEBUG_ALGORITHM);\r
-           } // if (NZERD > 0)\r
-       } // if (NINFD > 0 || NZERD > 0)\r
+           if (NZERD[0] > 0) {\r
+               Preferences.debug("   CAN --> PHY : CONDTN NO = " + COCAP[0] + "*U**" + EXCAP[0] + "\n", Preferences.DEBUG_ALGORITHM);\r
+           } // if (NZERD[0] > 0)\r
+       } // if (NINFD[0] > 0 || NZERD[0] > 0)\r
        \r
        \r
        PI=4E+0*Math.PI;\r
@@ -12834,7 +12850,7 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
         Preferences.debug("SUBARC   % PHYSICAL        % CIRCLE\n",Preferences.DEBUG_ALGORITHM);\r
         for (I=1; I <= TNSUA; I++) {\r
             ANGSP=VTARG[I]-VTARG[I-1];\r
-            Preferences.debug(I + "     " + (ARCLN[I-1]/TOTLN) + "    " + (ANGSP/2.0/PI) + "\n", Preferences.DEBUG_ALGORITHM);\r
+            Preferences.debug(I + "     " + (ARCLN[I-1]/TOTLN[0]) + "    " + (ANGSP/2.0/PI) + "\n", Preferences.DEBUG_ALGORITHM);\r
         }\r
    \r
        \r
@@ -12848,7 +12864,7 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
                 LA=0.0;\r
             }\r
             else {\r
-                CR=2.0*PI*ARCLN[I-1]/Math.abs(ANGSP)/TOTLN;\r
+                CR=2.0*PI*ARCLN[I-1]/Math.abs(ANGSP)/TOTLN[0];\r
                 if (ERARC[I-1] == 0.0) {\r
                     LA=OFLOW;\r
                 }\r
@@ -12863,79 +12879,75 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
             Preferences.debug(I + " " + LA + " " + CR + "\n", Preferences.DEBUG_ALGORITHM);\r
         } // for (I=1; I <= TNSUA; I++)\r
        \r
-       Preferences.debug("MINIMUM SUBARC RESOLUTION IS " + RESMN + " ON SUBARC " + IMNLA + "\n",Preferences.DEBUG_ALGORITHM);\r
+       Preferences.debug("MINIMUM SUBARC RESOLUTION IS " + RESMN[0] + " ON SUBARC " + IMNLA + "\n",Preferences.DEBUG_ALGORITHM);\r
        System.out.println();\r
-       System.out.println("MINIMUM SUBARC RESOLUTION: " + RESMN);\r
+       System.out.println("MINIMUM SUBARC RESOLUTION: " + RESMN[0]);\r
        \r
        Preferences.debug("\n", Preferences.DEBUG_ALGORITHM);\r
-       if (!MAP11 || RESMN[0] < CRRES) {\r
+       if (!MAP11[0] || RESMN[0] < CRRES) {\r
 \r
            // MESSAGE TO STANDARD OUTPUT\r
        \r
            System.out.println("                   *** W A R N I N G  ***");\r
            Preferences.debug("                   *** W A R N I N G  ***\n", Preferences.DEBUG_ALGORITHM);\r
-               IF (RESMN.LT.CRRES) THEN\r
-                 WRITE(*,5) 'THE ABOVE RESOLUTION IS TOO SMALL:'\r
-               ENDIF\r
-               IF (.NOT.MAP11) THEN\r
-                 WRITE(*,5) 'BCF DERIVATIVE CHANGES SIGN:'\r
-               ENDIF\r
-       C\r
-               WRITE(CH0,*) '        ***  W A R N I N G ***'\r
-               WRITE(CH0,*)\r
-               I=0\r
-       C\r
-               IF (RESMN.LT.1) THEN\r
-                 I=I+1\r
-                 WRITE(CH0,190) I,'.  THE ABOVE SUBARC RESOLUTION MEANS THAT IT\r
-            + WILL BE PRACTICALLY'\r
-                 WRITE(CH0,*) '    IMPOSSIBLE FOR THE INVERSE MAP TO DISCRIMINA\r
-            +TE CORRECTLY'\r
-                 WRITE(CH0,200) '    BETWEEN NEIGHBOURING POINTS NEAR SUB ARC '\r
-            +,IMNLA \r
-               ELSE IF (RESMN.LT.CRRES) THEN\r
-                 I=I+1\r
-                 WRITE(CH0,190) I,'. THE ABOVE SUBARC RESOLUTION MEANS THAT THE\r
-            + INVERSE MAP MAY NOT'\r
-                 WRITE(CH0,*) '    BE ABLE TO RELIABLY DISCRIMINATE CORRECTLY B\r
-            +ETWEEN'\r
-                 WRITE(CH0,200) '    NEIGHBOURING POINTS NEAR ARC ',IMNLA \r
-               ENDIF\r
-       190   FORMAT(/,I1,A)\r
-       200   FORMAT(A,I3)\r
-       C\r
-               IF (NCRVS .GT. 0) THEN\r
-                 I=I+1\r
-                 WRITE(CH0,190) I,'.  THERE IS A COMPLETE REVERSAL OF DIRECTION\r
-            + ON THE FOLLOWING SUB ARCS:'\r
-                 WRITE(CH0,'(T10,I3)') (ICRVS(J),J=1,NCRVS)\r
-               ENDIF\r
-       C\r
-               IF (NPRVS .GT. 0) THEN\r
-                 I=I+1\r
-                 WRITE(CH0,190) I,'.  THERE IS A REVERSAL OF DIRECTION WITHIN T\r
-            +HE FOLLOWING SUB ARCS:'\r
-                 WRITE(CH0,'(T10,I3)') (IPRVS(J),J=1,NPRVS)\r
-                 WRITE(CH0,*) '    THE CORRESPONDING MINIMUM VALUES OF THE BOUN\r
-            +DARY CORRESPONDENCE'\r
-                 WRITE(CH0,*) '    DERIVATIVE ARE:'\r
-                 WRITE(CH0,'(T10,E9.2)') (BCDMN(J),J=1,NPRVS)\r
-               ENDIF\r
-       } // if (!MAP11 || RESMN < CRRES)\r
-             CLOSE(CH0)\r
-       999   CONTINUE\r
-       C\r
-       C**** WRITE CLOSING MESSAGE TO STANDARD OUTPUT CHANNEL AND LISTING FILE\r
-       C\r
-             CALL WRTAIL(5,0,IER)\r
-             CALL WRTAIL(5,CH0,IER)\r
-       C */\r
+           if (RESMN[0] < CRRES) {\r
+               System.out.println("THE ABOVE RESOLUTION IS TOO SMALL:");\r
+               Preferences.debug("THE ABOVE RESOLUTION IS TOO SMALL:\n", Preferences.DEBUG_ALGORITHM);\r
+           }\r
+           if (!MAP11[0]) {\r
+               System.out.println("BCF DERIVATIVE CHANGES SIGN:");\r
+               Preferences.debug("BCF DERIVATIVE CHANGES SIGN\n",Preferences.DEBUG_ALGORITHM);\r
+           }\r
+       \r
+           Preferences.debug("***  W A R N I N G ***\n", Preferences.DEBUG_ALGORITHM);\r
+           Preferences.debug("\n", Preferences.DEBUG_ALGORITHM);\r
+           I=0;\r
+       \r
+           if (RESMN[0] < 1) {\r
+               I=I+1;\r
+               Preferences.debug(I + ".  THE ABOVE SUBARC RESOLUTION MEANS THAT IT WILL BE PRACTICALLY\n", Preferences.DEBUG_ALGORITHM);\r
+               Preferences.debug("    IMPOSSIBLE FOR THE INVERSE MAP TO DISCRIMINATE CORRECTLY\n", Preferences.DEBUG_ALGORITHM);\r
+               Preferences.debug("    BETWEEN NEIGHBOURING POINTS NEAR SUB ARC " + IMNLA + "\n", Preferences.DEBUG_ALGORITHM);\r
+           } // if (RESMN[0] < 1)\r
+           else if (RESMN[0] < CRRES) {\r
+               I=I+1;\r
+               Preferences.debug(I + ". THE ABOVE SUBARC RESOLUTION MEANS THAT THE INVERSE MAP MAY NOT\n",Preferences.DEBUG_ALGORITHM);\r
+               Preferences.debug("    BE ABLE TO RELIABLY DISCRIMINATE CORRECTLY BETWEEN\n", Preferences.DEBUG_ALGORITHM);\r
+               Preferences.debug("    NEIGHBOURING POINTS NEAR ARC " + IMNLA + "\n", Preferences.DEBUG_ALGORITHM); \r
+           } // else if (RESMN[0] < CRRES)\r
+       \r
+           if (NCRVS[0] > 0) {\r
+               I=I+1;\r
+               Preferences.debug(I + ".  THERE IS A COMPLETE REVERSAL OF DIRECTION ON THE FOLLOWING SUB ARCS:\n", Preferences.DEBUG_ALGORITHM);\r
+               for (J = 0; J <  NCRVS[0]; J++) {\r
+                       Preferences.debug("         " + ICRVS[J] + "\n", Preferences.DEBUG_ALGORITHM);\r
+               }\r
+           } // if (NCARVS[0] > 0)\r
+       \r
+           if (NPRVS[0] > 0) {\r
+               I=I+1;\r
+               Preferences.debug(I + ".  THERE IS A REVERSAL OF DIRECTION WITHIN THE FOLLOWING SUB ARCS:\n", Preferences.DEBUG_ALGORITHM);\r
+               for (J = 0; J < NPRVS[0]; J++) {\r
+                       Preferences.debug("         " + IPRVS[J] + "\n", Preferences.DEBUG_ALGORITHM);\r
+               }\r
+               Preferences.debug("    THE CORRESPONDING MINIMUM VALUES OF THE BOUN+DARY CORRESPONDENCE\n",Preferences.DEBUG_ALGORITHM);\r
+               Preferences.debug("    DERIVATIVE ARE:\n", Preferences.DEBUG_ALGORITHM);\r
+               for (J = 0; J < NPRVS[0]; J++) {\r
+                   Preferences.debug("         " + BCDMN[J] + "\n", Preferences.DEBUG_ALGORITHM);\r
+               }\r
+           } // if (NPRVS[0] > 0)\r
+       } // if (!MAP11[0] || RESMN[0] < CRRES)\r
+       \r
+       // WRITE CLOSING MESSAGE TO STANDARD OUTPUT CHANNEL AND LISTING FILE\r
+       \r
+        WRTAIL(5,0,IER[0], null);\r
+            \r
     } // private void CNDPLT\r
 \r
-    private void DIAGN4(double CCAPH, double COCAP, double COPHC, double CPHCA, double EXCAP,\r
-        double EXPHC, int ICRVS[], int IER[], int IPRVS[], int NASYM, int NCRVS, int NINFD,\r
-        int NPRVS, int NZERD, double ARCLN[], double ASYMP[], double BCDMN[], double CORXX[],\r
-       double TOTLN, double VTARG[], boolean MAP11[], int DGPOL[], int JATYP[], int LOSUB[],\r
+    private void DIAGN4(double CCAPH[], double COCAP[], double COPHC[], double CPHCA[], double EXCAP[],\r
+        double EXPHC[], int ICRVS[], int IER[], int IPRVS[], int NASYM[], int NCRVS[], int NINFD[],\r
+        int NPRVS[], int NZERD[], double ARCLN[], double ASYMP[], double BCDMN[], double CORXX[],\r
+       double TOTLN[], double VTARG[], boolean MAP11[], int DGPOL[], int JATYP[], int LOSUB[],\r
        int NARCS, int NQPTS, int NXINT, int PARNT[], int TNSUA, double A1COF[], double ACOEF[],\r
        double B1COF[], double BCFSN[], double BCOEF[], double H0VAL[], double H1VAL[],\r
        double HALEN[], double JACIN[], double MIDPT[], double SOLUN[], double QUPTS[],\r
@@ -12948,349 +12960,463 @@ public class SymmsIntegralMapping extends AlgorithmBase  {
        // +BCDMN(*),BCFSN(*),BCOEF(*),CORXX(*),JACIN(*),MIDPT(*),H0VAL(*),\r
        // +H1VAL(*),HALEN(*),SOLUN(*),VTARG(*),QUPTS(*),QUWTS(*)\r
        // LOGICAL MAP11\r
-       /*C\r
-       C     IER=0  - NORMAL EXIT\r
-       C     IER=50 - LOCAL PARAMETER MXCOF MUST BE >= NQPTS.\r
-       C     IER=51 - NON-ANALYTIC ARC DETECTED\r
-       C\r
-       C**** LOCAL VARIABLES\r
-       C\r
-             INTEGER AJT,DG,I,I1,IA,JT,K,LOD,LOM,MININ,MXCOF,NINTS,PT,QP\r
-             REAL AL,ATOL,BT,CC,COF,D,DSDT,H0,HH,HL,JACSUM,MD,MCHEP,MPT,\r
-            +PHI,R1MACH,RTOL,SEND,SINC,SJT,SS,SUM,TERM,TINC,TT,TUPI,X,XX,YMAX,\r
-            +YMIN,YY\r
-             COMPLEX PARFUN,T1,T2\r
-             COMMON /DSDTDA/PT,MD,HL\r
-             PARAMETER (MININ=20,MXCOF=32,QP=4)\r
-             REAL JACOF(MXCOF)\r
-             EXTERNAL DSDT,JACSUM,PARFUN,R1MACH\r
-       C\r
-       C     INITIALISE SOME CONSTANTS\r
-       C\r
-             TUPI=8E+0*ATAN(1E+0)\r
-             MCHEP=R1MACH(4)\r
-             RTOL=1E+1*MCHEP\r
-             ATOL=1E+2*MCHEP\r
-             NCRVS=0\r
-             NPRVS=0\r
-             CCAPH=0E+0\r
-             CPHCA=0E+0\r
-             MAP11=.TRUE.\r
-             YMAX=R1MACH(2)\r
-             NASYM=0\r
-       C\r
-       C     START TO COMPUTE THE ARC LENGTHS OF EACH SUBARC (ARCLN) AND THE  \r
-       C     TOTAL LENGTH (TOTLN) OF THE BOUNDARY\r
-       C\r
-             TOTLN=0E+0\r
-             DO 10 IA=1,TNSUA\r
-               PT=PARNT(IA)\r
-               MD=MIDPT(IA)\r
-               HL=HALEN(IA)\r
-               T1=CMPLX(MD+HL)\r
-               T2=CMPLX(MD-HL)\r
-       C\r
-       C****   COMPOSITE QP-PANEL GAUSS-LEGENDRE ESTIMATE FOR ARCLN(IA)\r
-       C\r
-               HH=1E+0/QP\r
-               SUM=0E+0\r
-               DO 6 K=1,QP\r
-                 MPT=-1E+0+(2E+0*K-1E+0)*HH\r
-                 DO 3 I=1,NQPTS\r
-                   X=MPT+HH*QUPTS(I)\r
-                   SUM=SUM+QUWTS(I)*DSDT(X)\r
-       3         CONTINUE\r
-       6       CONTINUE\r
-               ARCLN(IA)=HH*SUM\r
-               TOTLN=TOTLN+ARCLN(IA)\r
-       10    CONTINUE\r
-       C\r
-       C     TEST FOR COMPLETE REVERSAL OF DIRECTION OF A BOUNDARY SUBARC ON \r
-       C     THE UNIT DISC.\r
-       C\r
-             DO 20 IA=1,TNSUA\r
-               IF (VTARG(IA+1) .LT. VTARG(IA)) THEN\r
-                 NCRVS=NCRVS+1\r
-                 ICRVS(NCRVS)=IA\r
-                 MAP11=.FALSE.\r
-               ENDIF\r
-       20    CONTINUE\r
-       C\r
-       C     COMPUTE THE NUMBERS *NINFD* (*NZERD*) OF POINTS WHERE THE \r
-       C     DERIVATIVE OF THE MAP PHYSICAL --> CANONICAL IS RESPECTIVELY\r
-       C     INFINITE (ZERO).     \r
-       C\r
-             NINFD=0\r
-             NZERD=0\r
-             DO 25 I=1,NARCS\r
-               IF (JACIN(I) .LT. 0E+0) THEN\r
-                 NINFD=NINFD+1\r
-               ELSE IF (JACIN(I) .GT. 0E+0) THEN\r
-                 NZERD=NZERD+1\r
-               ENDIF\r
-       25    CONTINUE\r
-       C\r
-       C     NOW START TO EVALUATE THE DIMENSIONLESS BOUNDARY CORRESPONDENCE\r
-       C     DERIVATIVE AT SELECTED VALUES OF DIMENSIONLESS ARC LENGTH;\r
-       C     OUTPUT RESULTS FOR SUBSEQUENT GRAPH PLOTTING IF REQUIRED AND\r
-       C     TEST FOR SIGN CHANGES IN THIS DERIVATIVE.\r
-       C\r
-             SS=0E+0\r
-             SEND=0E+0\r
-             DO 60 IA=1,TNSUA\r
-               NINTS=MAX(MININ,NINT(ARCLN(IA)*NXINT/TOTLN))\r
-               TINC=2E+0/NINTS\r
-               DG=DGPOL(IA)\r
-               IF (DG+1 .GT. MXCOF) THEN\r
-                 IER=50\r
-                 RETURN\r
-               ENDIF\r
-               JT=JATYP(IA)\r
-               AJT=ABS(JT)\r
-               H0=H0VAL(AJT)\r
-               BT=JACIN(AJT)\r
-               AL=1E+0/(1E+0+BT)\r
-               PT=PARNT(IA)\r
-               MD=MIDPT(IA)\r
-               HL=HALEN(IA)\r
-               LOM=LOSUB(IA)\r
-               LOD=(AJT-1)*NQPTS+1 \r
-               IF (JT.GT.0) THEN\r
-                 CC=VTARG(IA)-VTARG(1)\r
-               ELSE\r
-                 CC=VTARG(IA+1)-VTARG(1)\r
-               ENDIF\r
-               DO 30 I=1,DG+1\r
-                 I1=I+LOM-1\r
-                 JACOF(I)=SOLUN(I1)\r
-       30      CONTINUE\r
-               SJT=SIGN(1E+0,REAL(JT))\r
-               DO 40 I=2,DG+1,2\r
-                 JACOF(I)=SJT*JACOF(I)\r
-       40      CONTINUE\r
-               TT=-1E+0\r
-               D=DSDT(TT)\r
-               YMIN=YMAX\r
-               IF (IA .EQ. 1) THEN\r
-                 XX=0E+0  \r
-                 IF (BT .LT. 0E+0) THEN\r
-                   YY=YMAX\r
-                   NASYM=NASYM+1\r
-                   ASYMP(NASYM)=XX\r
-                 ELSE IF (BT .GT. 0E+0) THEN\r
-                   YY=0E+0\r
-                 ELSE\r
-                   PHI=JACSUM(TT,DG,ACOEF(LOD),BCOEF(LOD),H0,JACOF)\r
-                   IF (D .EQ. 0E+0) THEN\r
-                     IER=51\r
-                     RETURN \r
-                   ENDIF \r
-                   YY=TOTLN*PHI/D\r
-                 ENDIF\r
-                 IF (NINFD .EQ. 0E+0) THEN\r
-                   CPHCA=TUPI*ABS(YY)/TOTLN\r
-                 ENDIF\r
-                 IF (NZERD .EQ. 0E+0) THEN\r
-                   IF (YY .EQ. 0E+0) THEN\r
-                     CCAPH=YMAX\r
-                   ELSE\r
-                     CCAPH=TOTLN/TUPI/ABS(YY)\r
-                   ENDIF\r
-                 ENDIF\r
-                 WRITE(OUCH1,902) XX,YY\r
-                 YY=0E+0\r
-                 WRITE(OUCH0,902) XX,YY\r
-                 CORXX(1)=0E+0\r
-               ENDIF\r
-       C\r
-       C       ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS FOR INFINITE\r
-       C       DERIVATIVE CASES.\r
-       C\r
-               IF (BT .LT. 0E+0) THEN\r
-                 PHI=JACSUM(-1E+0,DG-1,A1COF(LOD),B1COF(LOD),H1VAL(AJT),\r
-            +               BCFSN(LOM+1))\r
-                 PHI=BCFSN(LOM)-2E+0*PHI\r
-                 COF=ABS(PHI)/D**(BT+1E+0)\r
-                 TERM=MCHEP**BT*COF\r
-                 IF (TERM .GT. CPHCA) THEN\r
-                   CPHCA=TERM\r
-                   COPHC=COF\r
-                   EXPHC=BT\r
-                 ENDIF\r
-               ENDIF\r
-               IF (BT .GT. 0E+0) THEN\r
-                 PHI=JACSUM(-1E+0,DG-1,A1COF(LOD),B1COF(LOD),H1VAL(AJT),\r
-            +               BCFSN(LOM+1))\r
-                 PHI=BCFSN(LOM)-2E+0*PHI\r
-                 IF (ABS(PHI) .EQ. 0E+0) THEN\r
-                   CCAPH=YMAX\r
-                   COCAP=YMAX\r
-                   EXCAP=AL-1E+0\r
-                 ELSE\r
-                   COF=D/ABS(PHI)**AL\r
-                   TERM=MCHEP**(AL-1E+0)*COF\r
-                   IF (TERM .GT. CCAPH) THEN\r
-                     CCAPH=TERM\r
-                     COCAP=COF\r
-                     EXCAP=AL-1E+0\r
-                   ENDIF\r
-                 ENDIF\r
-               ENDIF\r
-       C\r
-       C       "DO 50" LOOP FOR POINTS INTERIOR TO ARC NUMBER IA\r
-       C\r
-               DO 50 I=1,NINTS-1\r
-                 TT=TT+TINC\r
-       C\r
-       C****     ARC LENGTH INCREASE BY GAUSS-LEGENDRE\r
-       C\r
-                 SUM=0E+0\r
-                 DO 45 K=1,NQPTS\r
-                   X=TT+5E-1*TINC*(QUPTS(K)-1E+0)\r
-                   SUM=SUM+QUWTS(K)*DSDT(X)\r
-       45        CONTINUE\r
-                 SINC=5E-1*TINC*SUM\r
-                 SS=SS+SINC\r
-                 XX=SS/TOTLN\r
-       C\r
-       C         EVALUATE DIMENSIONLESS BCF DERIVATIVE *YY*\r
-       C\r
-                 PHI=JACSUM(SJT*TT,DG,ACOEF(LOD),BCOEF(LOD),H0,JACOF)\r
-                 D=DSDT(TT)\r
-                 IF (D .EQ. 0E+0) THEN\r
-                   IER=51\r
-                   RETURN \r
-                 ENDIF \r
-                 YY=TOTLN*(1E+0+SJT*TT)**BT*PHI/D\r
-                 WRITE(OUCH1,902) XX,YY\r
-                 YMIN=MIN(YY,YMIN)\r
-       C\r
-       C         ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS FOR FINITE\r
-       C         DERIVATIVE CASES.\r
-       C\r
-                 IF (NINFD .EQ. 0E+0) THEN\r
-                   CPHCA=MAX(CPHCA,TUPI*ABS(YY)/TOTLN)\r
-                 ENDIF\r
-                 IF (NZERD .EQ. 0E+0) THEN\r
-                   IF (YY .EQ. 0E+0) THEN\r
-                     CCAPH=YMAX\r
-                   ELSE\r
-                     CCAPH=MAX(CCAPH,TOTLN/TUPI/ABS(YY))\r
-                   ENDIF\r
-                 ENDIF\r
-       C\r
-       C         EVALUATE DIMENSIONLESS BCF *YY*\r
-       C\r
-                 PHI=JACSUM(SJT*TT,DG-1,A1COF(LOD),B1COF(LOD),H1VAL(AJT),\r
-            +               BCFSN(LOM+1))\r
-                 PHI=BCFSN(LOM)-(1E+0-SJT*TT)*PHI\r
-                 YY=(CC+SJT*(1E+0+SJT*TT)**(1E+0+BT)*PHI)/TUPI\r
-                 WRITE(OUCH0,902) XX,YY\r
-       50      CONTINUE\r
-       C\r
-       C       NEXT TAKE END POINT OF ARC NUMBER IA\r
-       C\r
-               TT=1E+0\r
-               D=DSDT(TT)\r
-               SEND=SEND+ARCLN(IA)\r
-               SS=SEND\r
-               XX=SS/TOTLN\r
-       C\r
-       C       EVALUATE DIMENSIONLESS BCF DERIVATIVE *YY*\r
-       C\r
-               IF (JT .LT. 0E+0) THEN\r
-                   IF (BT .LT. 0E+0) THEN\r
-                     YY=YMAX\r
-                     NASYM=NASYM+1\r
-                     ASYMP(NASYM)=XX\r
-                   ELSE IF (BT .GT. 0E+0) THEN\r
-                     YY=0E+0\r
-                   ELSE\r
-                     PHI=JACSUM(SJT*TT,DG,ACOEF(LOD),BCOEF(LOD),H0,JACOF)\r
-                     IF (D .EQ. 0E+0) THEN\r
-                       IER=51\r
-                       RETURN \r
-                     ENDIF \r
-                     YY=TOTLN*PHI/D\r
-                   ENDIF\r
-               ELSE\r
-                   PHI=JACSUM(TT,DG,ACOEF(LOD),BCOEF(LOD),H0,JACOF)\r
-                   IF (D .EQ. 0E+0) THEN\r
-                     IER=51\r
-                     RETURN \r
-                   ENDIF \r
-                   YY=TOTLN*2E+0**BT*PHI/D\r
-               ENDIF\r
-               WRITE(OUCH1,902) XX,YY\r
-               YMIN=MIN(YY,YMIN)\r
-               IF (YMIN.LT.0E+0 .AND. (VTARG(IA+1) .GE. VTARG(IA))) THEN\r
-                 NPRVS=NPRVS+1\r
-                 IPRVS(NPRVS)=IA\r
-                 BCDMN(NPRVS)=YMIN\r
-                 MAP11=.FALSE.\r
-               ENDIF\r
-       C\r
-       C       ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS\r
-       C\r
-               IF (NINFD .EQ. 0E+0) THEN\r
-                 CPHCA=MAX(CPHCA,TUPI*ABS(YY)/TOTLN)\r
-               ENDIF\r
-               IF (NZERD .EQ. 0E+0) THEN\r
-                 IF (YY .EQ. 0E+0) THEN\r
-                   CCAPH=YMAX\r
-                 ELSE\r
-                   CCAPH=MAX(CCAPH,TOTLN/TUPI/ABS(YY))\r
-                 ENDIF\r
-               ENDIF\r
-               IF (BT .LT. 0E+0) THEN\r
-                 PHI=JACSUM(-1E+0,DG-1,A1COF(LOD),B1COF(LOD),H1VAL(AJT),\r
-            +               BCFSN(LOM+1))\r
-                 PHI=BCFSN(LOM)-2E+0*PHI\r
-                 COF=ABS(PHI)/D**(BT+1E+0)\r
-                 TERM=MCHEP**BT*COF\r
-                 IF (TERM .GT. CPHCA) THEN\r
-                   CPHCA=TERM\r
-                   COPHC=COF\r
-                   EXPHC=BT\r
-                 ENDIF\r
-               ENDIF\r
-               IF (BT .GT. 0E+0) THEN\r
-                 PHI=JACSUM(-1E+0,DG-1,A1COF(LOD),B1COF(LOD),H1VAL(AJT),\r
-            +               BCFSN(LOM+1))\r
-                 PHI=BCFSN(LOM)-2E+0*PHI\r
-                 IF (ABS(PHI) .EQ. 0E+0) THEN\r
-                   CCAPH=YMAX\r
-                   COCAP=YMAX\r
-                   EXCAP=AL-1E+0\r
-                 ELSE\r
-                   COF=D/ABS(PHI)**AL\r
-                   TERM=MCHEP**(AL-1E+0)*COF\r
-                   IF (TERM .GT. CCAPH) THEN\r
-                     CCAPH=TERM\r
-                     COCAP=COF\r
-                     EXCAP=AL-1E+0\r
-                   ENDIF\r
-                 ENDIF\r
-               ENDIF\r
-       C\r
-       C       EVALUATE DIMENSIONLESS BCF *YY*\r
-       C\r
-               YY=(VTARG(IA+1)-VTARG(1))/TUPI\r
-               WRITE(OUCH0,902) XX,YY\r
-               IF (JT.LT.0) THEN\r
-                 CORXX(PT+1)=XX\r
-               ENDIF\r
-       C\r
-       60    CONTINUE\r
-       C          \r
-       901   FORMAT(2E16.8,1X,A3)\r
-       902   FORMAT(2E16.8)\r
-       C\r
-       C     NORMAL EXIT\r
-       C\r
-             IER=0\r
-       C */\r
+       \r
+       //     IER=0  - NORMAL EXIT\r
+       //     IER=50 - LOCAL PARAMETER MXCOF MUST BE >= NQPTS.\r
+       //     IER=51 - NON-ANALYTIC ARC DETECTED\r
+       \r
+        // LOCAL VARIABLES\r
+       \r
+       final int MININ = 20;\r
+       final int MXCOF = 32;\r
+       final int QP = 4;\r
+       int AJT,DG,I,I1,IA,JT,K,LOD,LOM,NINTS;\r
+       double AL,ATOL,BT,CC,COF,D,DSDT,H0,HH,JACSUM,MCHEP,MPT,\r
+            PHI,R1MACH,RTOL,SEND,SINC,SJT,SS,SUM,TERM,TINC,TT,TUPI,X,XX,YMAX,\r
+            YMIN,YY;\r
+       double T1[] = new double[2];\r
+       double T2[] = new double[2];\r
+       //COMPLEX PARFUN,T1,T2\r
+       //COMMON /DSDTDA/PT,MD,HL\r
+       double JACOF[] = new double[MXCOF];\r
+       //EXTERNAL DSDT,JACSUM,PARFUN,R1MACH\r
+       double A[];\r
+       double B[];\r
+       int N;\r
+       double CO[];\r
+       \r
+       //     INITIALISE SOME CONSTANTS\r
+       \r
+       TUPI=2.0*Math.PI;\r
+       MCHEP=EPS;\r
+       RTOL=10.0*MCHEP;\r
+       ATOL=100.0*MCHEP;\r
+       NCRVS[0]=0;\r
+       NPRVS[0]=0;\r
+       CCAPH[0]=0.0;\r
+       CPHCA[0]=0.0;\r
+       MAP11[0]=true;\r
+       YMAX=Double.MAX_VALUE;\r
+       NASYM[0]=0;\r
+       \r
+       //     START TO COMPUTE THE ARC LENGTHS OF EACH SUBARC (ARCLN) AND THE  \r
+       //     TOTAL LENGTH (TOTLN) OF THE BOUNDARY\r
+       \r
+       TOTLN[0]=0.0;\r
+       /*for (IA=1; IA <= TNSUA; IA++) {\r
+           PT=PARNT[IA-1];\r
+           MD=MIDPT[IA-1];\r
+           HL=HALEN[IA-1];\r
+           T1[0]=MD+HL;\r
+           T1[1]= 0.0;\r
+           T2[0]=MD-HL;\r
+           T2[1] = 0.0;\r
+       \r
+           // COMPOSITE QP-PANEL GAUSS-LEGENDRE ESTIMATE FOR ARCLN(IA)\r
+       \r
+           HH=1.0/QP;\r
+           SUM=0.0;\r
+           for (K=1; K <= QP; K++) {\r
+               MPT=-1.0+(2.0*K-1.0)*HH;\r
+               for (I=1; I <= NQPTS; I++) {\r
+                   X=MPT+HH*QUPTS[I-1];\r
+                   SUM=SUM+QUWTS[I-1]*DSDT(X);\r
+               } // for (I=1; I <= NQPTS; I++)\r
+           } // for (K=1; K <= QP; K++) \r
+           ARCLN[IA-1]=HH*SUM;\r
+           TOTLN[0]=TOTLN[0]+ARCLN[IA-1];\r
+       } // for (IA=1; IA <= TNSUA; IA++)\r
+       \r
+       // TEST FOR COMPLETE REVERSAL OF DIRECTION OF A BOUNDARY SUBARC ON \r
+       //     THE UNIT DISC.\r
+       \r
+       for (IA=1; IA <= TNSUA; IA++) {\r
+           if (VTARG[IA] < VTARG[IA-1]) {\r
+                 NCRVS[0]=NCRVS[0]+1;\r
+                 ICRVS[NCRVS[0]-1]=IA;\r
+                 MAP11[0]=false;\r
+           } // if (VTARG[IA] < VTARG[IA-1])\r
+       } // for (IA=1; IA <= TNSUA; IA++)\r
+       \r
+       //     COMPUTE THE NUMBERS *NINFD* (*NZERD*) OF POINTS WHERE THE \r
+       //     DERIVATIVE OF THE MAP PHYSICAL --> CANONICAL IS RESPECTIVELY\r
+       //     INFINITE (ZERO).     \r
+       \r
+       NINFD[0]=0;\r
+       NZERD[0]=0;\r
+       for (I=1; I <= NARCS; I++) {\r
+           if (JACIN[I-1] < 0.0) {\r
+               NINFD[0]=NINFD[0]+1;\r
+           }\r
+           else if (JACIN[I-1] > 0.0) {\r
+                NZERD[0]=NZERD[0]+1;\r
+           }\r
+       } // for (I=1; I <= NARCS; I++)\r
+       \r
+       //     NOW START TO EVALUATE THE DIMENSIONLESS BOUNDARY CORRESPONDENCE\r
+       //     DERIVATIVE AT SELECTED VALUES OF DIMENSIONLESS ARC LENGTH;\r
+       //     OUTPUT RESULTS FOR SUBSEQUENT GRAPH PLOTTING IF REQUIRED AND\r
+       //     TEST FOR SIGN CHANGES IN THIS DERIVATIVE.\r
+       \r
+       SS=0.0;\r
+       SEND=0.0;\r
+       for (IA=1; IA <= TNSUA; IA++) {\r
+       NINTS=Math.max(MININ,(int)Math.round(ARCLN[IA-1]*NXINT/TOTLN[0]));\r
+       TINC=2.0/NINTS;\r
+       DG=DGPOL[IA-1];\r
+       if (DG+1 > MXCOF) {\r
+           IER[0]=50;\r
+           return;\r
+       }\r
+       JT=JATYP[IA-1];\r
+       AJT=Math.abs(JT);\r
+       H0=H0VAL[AJT-1];\r
+       BT=JACIN[AJT-1];\r
+       AL=1.0/(1.0+BT);\r
+       PT=PARNT[IA-1];\r
+       MD=MIDPT[IA-1];\r
+       HL=HALEN[IA-1];\r
+       LOM=LOSUB[IA-1];\r
+       LOD=(AJT-1)*NQPTS+1; \r
+       if (JT > 0) {\r
+           CC=VTARG[IA-1]-VTARG[0];\r
+       }\r
+       else {\r
+           CC=VTARG[IA]-VTARG[0];\r
+       }\r
+       for (I=1; I <= DG+1; I++) {\r
+           I1=I+LOM-1;\r
+           JACOF[I-1]=SOLUN[I1-1];\r
+       } // for (I=1; I <= DG+1; I++)\r
+       if (JT >= 0) {\r
+               SJT = 1.0;\r
+       }\r
+       else {\r
+               SJT = -1.0;\r
+       }\r
+       for (I=2; I <+ DG+1; I +=2) {\r
+           JACOF[I-1]=SJT*JACOF[I-1];\r
+       }\r
+       TT=-1.0;\r
+       D=DSDT(TT);\r
+       YMIN=YMAX;\r
+       if (IA == 1) {\r
+           XX=0.0;  \r
+           if (BT < 0.0) {\r
+               YY=YMAX;\r
+               NASYM[0]=NASYM[0]+1;\r
+               ASYMP[NASYM[0]-1]=XX;\r
+           } // if (BT < 0.0)\r
+           else if (BT > 0.0) {\r
+               YY=0.0;\r
+           }\r
+           else {\r
+               A = new double[DG];\r
+               B = new double[DG];\r
+               for (N = 0; N < DG; N++) {\r
+                       A[N] = ACOEF[LOD+N-1];\r
+                       B[N] = BCOEF[LOD+N-1];\r
+               }\r
+               PHI=JACSUM(TT,DG,A,B,H0,JACOF);\r
+               if (D == 0.0) {\r
+                   IER[0]=51;\r
+                   return; \r
+               }\r
+               YY=TOTLN[0]*PHI/D;\r
+           } // else\r
+           if (NINFD[0] == 0.0) {\r
+               CPHCA[0]=TUPI*Math.abs(YY)/TOTLN[0];\r
+           }\r
+           if (NZERD[0] == 0.0) {\r
+               if (YY == 0.0) {\r
+                   CCAPH[0]=YMAX;\r
+               }\r
+               else {\r
+                   CCAPH[0]=TOTLN[0]/TUPI/Math.abs(YY);\r
+               }\r
+           } // if (NZERD[0] == 0.0)\r
+           Preferences.debug("XX = " + XX + " YY = " + YY + "\n", Preferences.DEBUG_ALGORITHM);\r
+           // WRITE(OUCH1,902) XX,YY\r
+           YY=0.0;\r
+           // WRITE(OUCH0,902) XX,YY\r
+           CORXX[0]=0.0;\r
+       } // if (IA == 1)\r
+       \r
+       //       ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS FOR INFINITE\r
+       //       DERIVATIVE CASES.\r
+       \r
+               if (BT < 0.0) {\r
+                       A = new double[DG-1];\r
+                       B = new double[DG-1];\r
+                       for (N = 0; N < DG-1; N++) {\r
+                               A[N] = A1COF[LOD+N-1];\r
+                               B[N] = B1COF[LOD+N-1];\r
+                       }\r
+                       CO = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               CO[N] = BCFSN[LOM+N];\r
+                       }\r
+                   PHI=JACSUM(-1E+0,DG-1,A,B,H1VAL[AJT-1],CO);\r
+                   PHI=BCFSN[LOM-1]-2.0*PHI;\r
+                   COF=Math.abs(PHI)/Math.pow(D,(BT+1.0));\r
+                   TERM=Math.pow(MCHEP,BT)*COF;\r
+                   if (TERM > CPHCA[0]) {\r
+                       CPHCA[0]=TERM;\r
+                       COPHC[0]=COF;\r
+                       EXPHC[0]=BT;\r
+                   }\r
+               } // if (BT < 0.0)\r
+               if (BT > 0.0) {\r
+                       A = new double[DG-1];\r
+                       B = new double[DG-1];\r
+                       for (N = 0; N < DG-1; N++) {\r
+                               A[N] = A1COF[LOD+N-1];\r
+                               B[N] = B1COF[LOD+N-1];\r
+                       }\r
+                       CO = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               CO[N] = BCFSN[LOM+N];\r
+                       }\r
+                   PHI=JACSUM(-1.0,DG-1,A,B,H1VAL[AJT-1],CO);\r
+                   PHI=BCFSN[LOM-1]-2.0*PHI;\r
+                   if (Math.abs(PHI) == 0.0) {\r
+                       CCAPH[0]=YMAX;\r
+                       COCAP[0]=YMAX;\r
+                       EXCAP[0]=AL-1.0;\r
+                   }\r
+                   else {\r
+                           COF=D/Math.pow(Math.abs(PHI),AL);\r
+                           TERM=Math.pow(MCHEP,(AL-1.0))*COF;\r
+                           if (TERM > CCAPH[0]) {\r
+                               CCAPH[0]=TERM;\r
+                               COCAP[0]=COF;\r
+                               EXCAP[0]=AL-1.0;\r
+                           } // if (TERM > CCAPH[0])\r
+                   } // else\r
+               // if (BT > 0.0)\r
+       \r
+                // "DO 50" LOOP FOR POINTS INTERIOR TO ARC NUMBER IA\r
+       \r
+               for (I=1; I <= NINTS-1; I++) {\r
+                   TT=TT+TINC;\r
+       \r
+                   // ARC LENGTH INCREASE BY GAUSS-LEGENDRE\r
+       \r
+                   SUM=0.0;\r
+                   for (K=1; K <= NQPTS; K++) {\r
+                       X=TT+0.5*TINC*(QUPTS[K-1]-1.0);\r
+                       SUM=SUM+QUWTS[K-1]*DSDT(X);\r
+                   } //  for (K=1; K <= NQPTS; K++)\r
+                   SINC=0.5*TINC*SUM;\r
+                   SS=SS+SINC;\r
+                   XX=SS/TOTLN[0];\r
+       \r
+                   // EVALUATE DIMENSIONLESS BCF DERIVATIVE *YY*\r
+                   A = new double[DG];\r
+                       B = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               A[N] = ACOEF[LOD+N-1];\r
+                               B[N] = BCOEF[LOD+N-1];\r
+                       }\r
+                   PHI=JACSUM(SJT*TT,DG,A,B,H0,JACOF);\r
+                   D=DSDT(TT);\r
+                   if (D == 0.0) {\r
+                       IER[0]=51;\r
+                       return;\r
+                   } // if (D == 0.0)\r
+                   YY=TOTLN[0]*Math.pow((1.0+SJT*TT),BT)*PHI/D;\r
+                   Preferences.debug("XX == " + XX + " YY = " + YY + "\n", Preferences.DEBUG_ALGORITHM);\r
+                   YMIN=Math.min(YY,YMIN);\r
+       \r
+                   // ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS FOR FINITE\r
+                   // DERIVATIVE CASES.\r
+       \r
+                   if (NINFD[0] == 0.0) {\r
+                       CPHCA[0]=Math.max(CPHCA[0],TUPI*Math.abs(YY)/TOTLN[0]);\r
+                   }\r
+                   if (NZERD[0] == 0.0) {\r
+                       if (YY == 0.0) {\r
+                           CCAPH[0]=YMAX;\r
+                       }\r
+                       else {\r
+                           CCAPH[0]=Math.max(CCAPH[0],TOTLN[0]/TUPI/Math.abs(YY));\r
+                       }\r
+                   } // if (NZERD[0] == 0.0)\r
+       \r
+                   // EVALUATE DIMENSIONLESS BCF *YY*\r
+                   A = new double[DG-1];\r
+                       B = new double[DG-1];\r
+                       for (N = 0; N < DG-1; N++) {\r
+                               A[N] = A1COF[LOD+N-1];\r
+                               B[N] = B1COF[LOD+N-1];\r
+                       }\r
+                       CO = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               CO[N] = BCFSN[LOM+N];\r
+                       }\r
+                   PHI=JACSUM(SJT*TT,DG-1,A,B,H1VAL[AJT-1],CO);\r
+                   PHI=BCFSN[LOM-1]-(1.0-SJT*TT)*PHI;\r
+                   YY=(CC+SJT*Math.pow((1.0+SJT*TT),(1.0+BT))*PHI)/TUPI;\r
+                   Preferences.debug("XX == " + XX + " YY = " + YY + "\n", Preferences.DEBUG_ALGORITHM);\r
+               } // for (I=1; I <= NINTS-1; I++)\r
+       \r
+               // NEXT TAKE END POINT OF ARC NUMBER IA\r
+       \r
+               TT=1.0;\r
+               D=DSDT(TT);\r
+               SEND=SEND+ARCLN[IA-1];\r
+               SS=SEND;\r
+               XX=SS/TOTLN[0];\r
+       \r
+               // EVALUATE DIMENSIONLESS BCF DERIVATIVE *YY*\r
+       \r
+               if (JT < 0.0) {\r
+                   if (BT < 0.0) {\r
+                       YY=YMAX;\r
+                       NASYM[0]=NASYM[0]+1;\r
+                       ASYMP[NASYM[0]-1]=XX;\r
+                   } // if (BT < 0.)\r
+                   else if (BT > 0.0) {\r
+                       YY=0.0;\r
+                   }\r
+                   else {\r
+                       A = new double[DG];\r
+                       B = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                                A[N] = ACOEF[LOD+N-1];\r
+                                B[N] = BCOEF[LOD+N-1];\r
+                       }\r
+                       PHI=JACSUM(SJT*TT,DG,A,B,H0,JACOF);\r
+                       if (D == 0.0) {\r
+                           IER[0]=51;\r
+                           return;\r
+                       }\r
+                       YY=TOTLN[0]*PHI/D;\r
+                   } // else\r
+               } // if (JT < 0.0)\r
+               else {\r
+                       A = new double[DG];\r
+                       B = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                                A[N] = ACOEF[LOD+N-1];\r
+                                B[N] = BCOEF[LOD+N-1];\r
+                       }\r
+                   PHI=JACSUM(TT,DG,A,B,H0,JACOF);\r
+                   if (D == 0.0) {\r
+                       IER[0]=51;\r
+                       return;\r
+                   }\r
+                   YY=TOTLN[0]*Math.pow(2.0,BT)*PHI/D;\r
+               } // else\r
+               Preferences.debug("XX == " + XX + " YY = " + YY + "\n", Preferences.DEBUG_ALGORITHM);\r
+               YMIN=Math.min(YY,YMIN);\r
+               if (YMIN < 0.0 && (VTARG[IA] >= VTARG[IA-1])) {\r
+                   NPRVS[0]=NPRVS[0]+1;\r
+                   IPRVS[NPRVS[0]-1]=IA;\r
+                   BCDMN[NPRVS[0]-1]=YMIN;\r
+                   MAP11[0]= false;\r
+               } // if (YMIN < 0.0 && (VTARG[IA] >= VTARG[IA-1]))\r
+       \r
+               // ESTIMATE FUNCTION EVALUATION CONDITION NUMBERS\r
+       \r
+               if (NINFD[0] == 0.0) {\r
+                   CPHCA[0]=Math.max(CPHCA[0],TUPI*Math.abs(YY)/TOTLN[0]);\r
+               }\r
+               if (NZERD[0] == 0.0) {\r
+                   if (YY == 0.0) {\r
+                       CCAPH[0]=YMAX;\r
+                   }\r
+                   else {\r
+                       CCAPH[0]=Math.max(CCAPH[0],TOTLN[0]/TUPI/Math.abs(YY));\r
+                   }\r
+               } // if (NZERD[0] == 0.0)\r
+               if (BT < 0.0) {\r
+                       A = new double[DG-1];\r
+                       B = new double[DG-1];\r
+                       for (N = 0; N < DG-1; N++) {\r
+                               A[N] = A1COF[LOD+N-1];\r
+                               B[N] = B1COF[LOD+N-1];\r
+                       }\r
+                       CO = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               CO[N] = BCFSN[LOM+N];\r
+                       }\r
+                   PHI=JACSUM(-1.0,DG-1,A,B,H1VAL[AJT-1],CO);\r
+                   PHI=BCFSN[LOM-1]-2.0*PHI;\r
+                   COF=Math.abs(PHI)/Math.pow(D,(BT+1.0));\r
+                   TERM=Math.pow(MCHEP,BT)*COF;\r
+                   if (TERM > CPHCA[0]) {\r
+                       CPHCA[0]=TERM;\r
+                       COPHC[0]=COF;\r
+                       EXPHC[0]=BT;\r
+                   } // if (TERM > CPHCA[0])\r
+               } // if (BT < 0.0)\r
+               if (BT > 0.0) {\r
+                       A = new double[DG-1];\r
+                       B = new double[DG-1];\r
+                       for (N = 0; N < DG-1; N++) {\r
+                               A[N] = A1COF[LOD+N-1];\r
+                               B[N] = B1COF[LOD+N-1];\r
+                       }\r
+                       CO = new double[DG];\r
+                       for (N = 0; N < DG; N++) {\r
+                               CO[N] = BCFSN[LOM+N];\r
+                       }\r
+                   PHI=JACSUM(-1.0,DG-1,A,B,H1VAL[AJT-1],CO);\r
+                   PHI=BCFSN[LOM-1]-2.0*PHI;\r
+                   if (Math.abs(PHI) == 0.0) {\r
+                       CCAPH[0]=YMAX;\r
+                       COCAP[0]=YMAX;\r
+                       EXCAP[0]=AL-1.0;\r
+                   }\r
+                   else {\r
+                       COF=D/Math.pow(Math.abs(PHI),AL);\r
+                       TERM=Math.pow(MCHEP,(AL-1.0))*COF;\r
+                       if (TERM > CCAPH[0]) {\r
+                           CCAPH[0]=TERM;\r
+                           COCAP[0]=COF;\r
+                           EXCAP[0]=AL-1.0;\r
+                       } // if (TERM > CCAPH[0])\r
+                   } // else\r
+               } // if (BT > 0.0)\r
+       \r
+               // EVALUATE DIMENSIONLESS BCF *YY*\r
+       \r
+               YY=(VTARG[IA]-VTARG[0])/TUPI;\r
+               Preferences.debug("XX == " + XX + " YY = " + YY + "\n", Preferences.DEBUG_ALGORITHM);\r
+               if (JT < 0) {\r
+                   CORXX[PT]=XX;\r
+               }\r
+       \r
+       } // for (IA=1; IA <= TNSUA; IA++)\r
+       \r
+       //     NORMAL EXIT\r
+       \r
+       IER[0]=0;*/\r
+        \r
     } // private void DIAGN4\r
+    \r
+    private double DSDT(double X) {\r
+\r
+        // TO COMPUTE THE PARAMETRIC DERIVATIVE OF THE ARC LENGTH FOR THE\r
+        // ARC SPECIFIED IN DSDTDA\r
+\r
+       // LOCAL VARIABLES\r
+\r
+       // INTEGER PT\r
+       // REAL MD,HL\r
+       // COMPLEX DPARFN\r
+       // COMMON /DSDTDA/PT,MD,HL\r
+       // EXTERNAL DPARFN\r
+        double DIN[] = new double[2];\r
+        DIN[0] = MD + HL *X;\r
+        DIN[1] = 0.0;\r
+        double DOUT[] = DPARFN(PT,DIN);\r
+        double result = zabs(DOUT[0],DOUT[1])*HL;\r
+        return result;\r
+\r
+    } // private double DSDT\r
+\r
 \r
       /**\r
        * zabs computes the absolute value or magnitude of a double precision complex variable zr + j*zi.\r