text:
CCCCCCCCCCCCCCCCCCCCCCCC MIDACO FORTRAN HEADER CCCCCCCCCCCCCCCCCCCCCCCCC
C
C     _|      _|  _|_|_|  _|_|_|      _|_|      _|_|_|    _|_|
C     _|_|  _|_|    _|    _|    _|  _|    _|  _|        _|    _|
C     _|  _|  _|    _|    _|    _|  _|_|_|_|  _|        _|    _|
C     _|      _|    _|    _|    _|  _|    _|  _|        _|    _|
C     _|      _|  _|_|_|  _|_|_|    _|    _|    _|_|_|    _|_|
C
C                                          Version 4.0 (limited)
C
C    MIDACO - Mixed Integer Distributed Ant Colony Optimization
C    ----------------------------------------------------------
C
C    MIDACO solves the general Mixed Integer Non-Linear Program (MINLP):
C
C
C       Minimize     F(X)           where X(1,...N-NI)   is *CONTINUOUS*
C                                   and   X(N-NI+1,...N) is *DISCRETE*
C
C       Subject to:  G_j(X)  =  0   (j=1,...ME)     Equality Constraints
C                    G_j(X) >=  0   (j=ME+1,...M) Inequality Constraints
C
C       And bounds:  XL <= X <= XU
C
C
C    MIDACO is a (heuristic) global optimization solver that stochastically
C    approximates a solution 'X' to the above displayed MINLP problem. MIDACO
C    is based on an extended Ant Colony Optimization framework (see [1]) in
C    combination with the Oracle Penalty Method (see [2]) for constraints 'G(X)'.
C
C    In case of mixed integer problems, the continuous variables are stored
C    first in 'X(1,...N-NI)', while the discrete (also called integer or
C    combinatorial) variables are stored last in 'X(N-NI+1,...N)'.
C    As an example consider:
C
C       X = (  0.1234,  5.6789,  1.0000,  2.0000,  3.0000)
C
C       where 'N' = 5 and 'NI' = 3
C
C    Note that all 'X' is of type double precision. Equality and inequality
C    constraints are handled in a similar way. The vector 'G' stores at first
C    the 'ME' equality constraints and behind those, the remaining 'M-ME'
C    inequality constraints are stored.
C
C    MIDACO is a derivative free black box solver and does not require the
C    relaxation of integer variables (this means, integer variables are
C    treated as categorical variables). MIDACO does not require any user
C    specified parameter tuning as it can run completely on 'Autopilot'
C    (all parameter set equal to zero). Optionally, the user can adjust
C    the MIDACO performance by setting some parameters explained below.
C
C
C    List of MIDACO subroutine arguments:
C    ------------------------------------
C
C    P  :   (Parallelization Factor)
C            If no parallelization is desired, set P = 1.
C
C    N  :   Number of optimization variables in total (continuous and integer ones).
C           'N' is the dimension of the iterate 'X' with X = (X_1,...,X_N).
C
C    NI :   Number of integer optimization variables. 'NI' <= 'N'.
C           Integer (discrete) variables must be stored at the end of 'X'.
C
C    M  :   Number of constraints in total (equality and inequality ones).
C           'P*M' is the dimension of a constraint vector 'G' with G = (G_1,...,G_M).
C
C    ME :   Number of equality constraints. 'ME' <= 'M'.
C           Equality constraints are stored in the beginning of 'G'.
C           Inequality constraints are stored in the end of 'G'.
c
C    X(P*N) :   Array containing the iterate 'X'. For P=1 (no parallelization)
C               'X' stores only one iterate and has length 'N'. For P>1
C               'X' contains several iterates, which are stored one after
C               another.
C
C    F(P)   :   Array containing the objective function value 'F' corresponding
C               to the iterates 'X'. For P=1 (no parallelization), 'F' is a single
C               value. For P>1 F stores several values, corresponding to 'X' one
C               after another.
C
C    G(P*M) :   Array containing the constraint values 'G'.For P=1 (no parallelization)
C               'G' has length 'M'. For P>1 'G' has length 'P*M' and stores the vectors
C               of constraints, corresponding to 'X' one after another.
C
C     XL(N) :   Array containing the lower bounds for the iterates 'X'.
C               Note that for integer dimensions (i > N-NI) the bounds should also be
C               discrete, e.g. X(i) = 1.0000.
C
C     XU(N) :   Array containing the upper bounds for the iterates 'X'.
C               Note that for integer dimensions (i > N-NI) the bounds should also be
C               discrete, e.g. X(i) = 1.0000.
C
C    IFLAG :    Information flag used by MIDACO. Initially MIDACO must be called with IFLAG=0.
C               If MIDACO works correctly, IFLAG flags lower than 0 are used for internal
C               communication. If MIDACO stops (either by submitting ISTOP=1 or automatically
C               by the FSTOP or AUTOSTOP parameter), an IFLAG value between 1 and 9 is returned
C               as final message. If MIDACO detects at start-up some critical problem setup, a
C               ***WARNING*** message is returned by IFLAG as value between 10 and 99. If
C               MIDACO detects an ***MIDACO INPUT ERROR***, an IFLAG value between 100 and 999
C               is returned and MIDACO stops. The individual IFLAG flags are as follows:
C
C               Final Messages:
C               ---------------
C               IFLAG = 1 : Feasible solution,   MIDACO was stopped by MAXEVAL or MAXTIME
C               IFLAG = 2 : Infeasible solution, MIDACO was stopped by MAXEVAL or MAXTIME
C               IFLAG = 3 : Feasible solution,   MIDACO stopped automatically by AUTOSTOP
C               IFLAG = 4 : Infeasible solution, MIDACO stopped automatically by AUTOSTOP
C               IFLAG = 5 : Feasible solution,   MIDACO stopped automatically by FSTOP
C
C               WARNING - Flags:
C               ----------------
C               IFLAG = 51 : Some X(i)  is greater/lower than +/- 1.0D+8 (try to avoid huge values!)
C               IFLAG = 52 : Some XL(i) is greater/lower than +/- 1.0D+8 (try to avoid huge values!)
C               IFLAG = 53 : Some XU(i) is greater/lower than +/- 1.0D+8 (try to avoid huge values!)
C
C               IFLAG = 61 : Some X(i)  should be discrete (e.g. 1.0), but is continuous (e.g. 1.234)
C               IFLAG = 62 : Some XL(i) should be discrete (e.g. 1.0), but is continuous (e.g. 1.234)
C               IFLAG = 63 : Some XU(i) should be discrete (e.g. 1.0), but is continuous (e.g. 1.234)
C
C               IFLAG = 71 : Some XL(i) = XU(i) (fixed variable)
C
C               IFLAG = 81 : F(X) has value NaN for starting point X
C               IFLAG = 82 : Some G(X) has value NaN for starting point X
C
C               IFLAG = 91 : FSTOP is greater/lower than +/- 1.0D+8
C               IFLAG = 92 : ORACLE is greater/lower than +/- 1.0D+8
C
C               ERROR - Flags:
C               --------------
C               IFLAG = 101 :   P   <= 0   or   P  > 1.0D+6
C               IFLAG = 102 :   N   <= 0   or   N  > 1.0D+6
C               IFLAG = 103 :   NI  <  0
C               IFLAG = 104 :   NI  >  N
C               IFLAG = 105 :   M   <  0   or   M  > 1.0D+6
C               IFLAG = 106 :   ME  <  0
C               IFLAG = 107 :   ME  >  M
C
C               IFLAG = 201 :   some X(i)  has type NaN
C               IFLAG = 202 :   some XL(i) has type NaN
C               IFLAG = 203 :   some XU(i) has type NaN
C               IFLAG = 204 :   some X(i) < XL(i)
C               IFLAG = 205 :   some X(i) > XU(i)
C               IFLAG = 206 :   some XL(i) > XU(i)
C
C               IFLAG = 301 :   PARAM(1) < 0   or   PARAM(1) > 1.0D+6
C               IFLAG = 302 :   PARAM(2) < 0   or   PARAM(2) > 1.0D+12
C               IFLAG = 303 :   PARAM(3) greater/lower than +/- 1.0D+12
C               IFLAG = 304 :   PARAM(4) < 0   or   PARAM(4) > 1.0D+6
C               IFLAG = 305 :   PARAM(5) greater/lower than +/- 1.0D+12
C               IFLAG = 306 :   |PARAM(6)| < 1   or   PARAM(6) > 1.0D+12
C               IFLAG = 307 :   PARAM(7) < 0   or   PARAM(7) > 1.0D+8
C               IFLAG = 308 :   PARAM(8) < 0   or   PARAM(8) > 100
C               IFLAG = 309 :   PARAM(7) < PARAM(8)
C               IFLAG = 310 :   PARAM(7) > 0 but PARAM(8) = 0
C               IFLAG = 311 :   PARAM(8) > 0 but PARAM(7) = 0
C               IFLAG = 312 :   PARAM(9) < 0   or   PARAM(9) > 1000
C               IFLAG = 313 :   Some PARAM(i) has type NaN
C
C               IFLAG = 401 :   ISTOP < 0 or ISTOP > 1
C
C               IFLAG = 501 :   Double precision work space size LRW is too small.
C                               ---> RW must be at least of size LRW = 200*N+2*M+1000
C
C               IFLAG = 601 :   Integer work space size LIW is too small.
C                               ---> IW must be at least of size LIW = 2*N+P+1000
C
C               IFLAG = 701 :   Input check failed! MIDACO must be called initially with IFLAG = 0
C
C               IFLAG = 801 :   P > PMAX (user must increase PMAX in the MIDACO source code)
C               IFLAG = 802 :   P*M+1 > PXM (user must increase PXM in the MIDACO source code)
C
C               IFLAG = 900 :   Invalid or corrupted LICENSE-KEY
C
C               IFLAG = 999 :   N > 4. The free test version is limited up to 4 variables.
C
C    ISTOP :   Communication flag to stop MIDACO. If MIDACO is called with
C              ISTOP = 1, MIDACO returns the best found solution in 'X' with
C              corresponding 'F' and 'G'. As long as MIDACO should continue
C              its search, ISTOP must be equal to 0.
C
C    PARAM() :  Array containing 9 parameters that can be selected by the user to adjust MIDACO.
C               (See the user manual for a more detailed description of individual parameters)
C
C     PARAM(1) :  [ACCURACY] Accuracy for constraint violation (measured as the L1-Norm over G(X)).
C                 If PARAM(1) is set to 0, MIDACO uses a default accuracy of 0.001. If the user
C                 desires a more precise accuracy, set PARAM(1) = 0.000001 for example. If the
C                 user desires a less precise accuracy, set PARAM(1) = 0.05 for example.
C
C     PARAM(2) :  [RANDOM-SEED] Random seed used for MIDACO's internal pseudo random number
C                 generator. Value must be a positive discrete value, e.g. PARAM(2) = 1,2,3,...
C                 The default random seed for MIDACO is zero.
C
C     PARAM(3) :  [FSTOP] Stopping criteria for MIDACO. MIDACO will stop, if a feasible
c                 solution 'X' with F(X) <= FSTOP is found. FSTOP is disabled, when FSTOP = 0.
C                 In case the user wishes to use zero as FSTOP, a tiny value (e.g. 0.000001)
C                 should be used instead.
C
C     PARAM(4) :  [AUTOSTOP] Automatic stopping criteria within MIDACO. AUTOSTOP defines the
C                 number of successive internal MIDACO restarts, after which no progress was
C                 made on the objective function. AUTOSTOP must be a (discrete) value >= 0.
C                 Some (very rough) examples for AUTOSTOP values are:
C                 AUTOSTOP = 1 (lowest chance of global optimality, but fastest runtime)
C                 AUTOSTOP = 10 (medium chance of global optimality)
C                 AUTOSTOP = 50 (high chance of global optimality)
C                 AUTOSTOP = 500 (very high chance of global optimality, but very long runtime)
C                 AUTOSTOP = 0 disables this stopping criteria.
C
C     PARAM(5) :  [ORACLE] This parameter affects only constrained problems. If PARAM(4)=0.0,
C                 MIDACO will use its internal oracle strategy. If PARAM(4) is not equal to 0.0,
C                 MIDACO will use PARAM(5) as initial oracle for its in-build oracle penalty
C                 function. In case the user wishes to use zero as ORACLE, a tiny value
C                 (e.g. 0.000001) should be used as ORACLE. Note: For most problems it is
C                 better to overestimate an oracle, instead of underestimating it.
C
C     PARAM(6) :  [FOCUS] This parameter focuses the MIDACO search process around the
C                 current best solution (which might be the starting point).
C
C     PARAM(7) :  [ANTS] Number of ants (stochastically sampled iterates) used by MIDACO
C                 in every generation. ANTS must be greater or equal to KERNEL.
C
C     PARAM(8) :  [KERNEL] Solution archive size of MIDACO. KERNEL must be lower or equal to ANTS.
C
C     PARAM(9) :  [CHARACTER] This parameter allows the user to use highly customized MIDACO
C                 internal parameter settings for specific problems. This parameter is currently
C                 only available as a service. In case you are interested in this service, please
C                 contact info@midaco-solver.com.
C
C
C    RW(LRW) :  Real workarray (Type: double precision) of length 'LRW'
C    LRW :      Length of 'RW'. 'LRW' must be greater or equal to  200*N+2*M+1000
C
C    IW(LIW) :  Integer workarray (Type: long integer) of length 'LIW'
C    LIW :      Length of 'IW'. 'LIW' must be greater or equal to  2*N+P+1000
C
C    KEY :  License-Key for MIDACO. Note that any licensed copy of MIDACO comes with an
C           individual 'KEY' determining the license owner and its affiliation.
C
C
C    References:
C    -----------
C
C    [1] Schlueter, M., Egea, J. A., Banga, J. R.:
C        "Extended ant colony optimization for non-convex mixed integer nonlinear programming",
C        Computers & Operations Research (Elsevier), Vol. 36 , Issue 7, Page 2217-2229, 2009.
C
C    [2] Schlueter M., Gerdts M.: "The oracle penalty method",
C        Journal of Global Optimization (Springer), Vol. 47(2),pp 293-325, 2010.
C
C
C    Author (C) :   Dr. Martin Schlueter
C                   Information Initiative Center,
C                   Division of Large Scale Computing Systems,
C                   Hokkaido University, JAPAN.
C
C    Email :        info@midaco-solver.com
C    URL :          http://www.midaco-solver.com
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      SUBROUTINE MIDACO( P, N, NI, M, ME, X, F, G, XL, XU, IFLAG,
     &                   ISTOP, PARAM, RW, LRW, IW, LIW, KEY)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      INTEGER P, N, NI, M, ME, IFLAG, LRW, IW, LIW, ISTOP
      DOUBLE PRECISION X, F, G, XL, XU, PARAM, RW
      CHARACTER*60 KEY
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      DIMENSION X(P*N),F(P),G(P*M),XL(N),XU(N),RW(LRW),IW(LIW),PARAM(9)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C     Increase the size of PMAX and PXM if necessary (IFLAG: 801 or 802)
      INTEGER PMAX, PXM
      PARAMETER (PMAX =  100)
      PARAMETER (PXM  = 1000)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C     Scrambled Source Code Starts Here
C
      INTEGER I
      DOUBLE PRECISION A(PMAX),B(PMAX),GM(PXM)
      if(iflag.eq.-999) istop = 1
      if(iflag.eq.0)then
c         Check P <= PMAX
          if(P.gt.PMAX)then
              iflag = 801
              return
          endif
c         Check P*M+1 <= PXM
          if(P*M+1.gt.PXM)then
              iflag = 802
              return
          endif
      endif
      if(m.gt.0)then
          do I = 1,P*M
              GM(I) = g(I)
          enddo
      endif
      GM(P*M+1) = 0.0D0
      call jfk(P,n,ni,m,me,x,f,GM,xl,xu,iflag,istop,
     &                 param,rw,lrw,iw,liw,PMAX,A,B,KEY)
      end
                  subroutine i015(f,m,i8,g,i5,i2,i19,i49,
     &                           i4,i32,i6,i99,i42,i93,i36)
            implicit none
      integer f,m,i8,i19,i49,i32,i6,i99,i42,i93,i,j,i10,i43,
     &  i21,i39,i62,i44,i38,i20,i95,i36(*),
     &                        i71
                           double precision g,i5,i2,i4,i02,i35
            dimension g(m),i5(m),i2(m),i4(i32),i6(i99)
         data i10,i43,i21,i39,i62,i38,i20
     &    /0,0,0,0,0,0,0/
                            data i35 /0.0d0/
                                   if(i42.eq.-30)then
                i10 = 0
                                                     i43 = 31 + f + 1
                                                 do i = 1,m
                        j             = int(i*i02(i4(1))) + 1
          i6(i43+i-1) = i6(i43+j-1)
                 i6(i43+j-1) = i
                                                if(i.ge.-25-i42)then
                      do j = 1,i
             i4(j) = dble(i42)-i02(i4(1))
                                enddo
                                                            endif
                                   enddo
                                  i6(31) = 1
                                                 i38 = i43 + m
                                      do i = 1,m
                                        i6(i38+i-1) = 0
                                                         enddo
                                          if(i4(1).ge.0.9d0)then
         i44 = 92+m+f+m
                             i6(i44) = 0
                        do i=1,30
                             i6(i44) = i6(i44) + i6(f+40+m+i+m)
                                         enddo
                       if(1372-i6(i44).ne.0)then
                                     do i=1,i99
                                i4(i) = dble(i6(i))
                      enddo
          i42 = int(i4(1))*1000
                                                            goto 22
                       endif
        endif
                    endif
                                   i95 = 0
                                            if(i6(31).eq.0)then
                    i20 = i6(i43+i10-1)
                           i6(30) = i20
                                i39 = i39 + 1
                    i21 = - i21
                                       i35 = i35 / dble(2**i36(14))
                          if(i35.lt.1.0d0/dble(10*i36(15)))then
                             i35  = 1.0d0/dble(10*i36(15))
        endif
         if(i20.gt.m-i8.and.i39.gt.i62)then
                                            i6(i38+i20-1) = 1
                                                   if(i10.ge.m) goto 2
                                      i95 = 1
           endif
                        i71 = i36(16)
        if(i20.le.m-i8.and.i39.gt.i71)then
                                i6(i38+i20-1) = 1
                                                     if(i10.ge.m) goto 2
                        i95 = 1
                                           endif
                            if(dabs(i5(i20)-i2(i20)).le.1.0d-12)then
                                  i6(i38+i20-1) = 1
                                              if(i10.ge.m) goto 2
           i95 = 1
                     endif
                  endif
                     if(i6(31).eq.1.or.i95.eq.1)then
                                                          i10 = i10 + 1
                          if(i10.gt.m) goto 2
                                  i20 = i6(i43+i10-1)
        i6(30) = i20
                                                     i39 = 1
                              if(i20.gt.m-i8)then
         if(i4(i19+i20-1).eq.i5(i20).or.
     &      i4(i19+i20-1).eq.i2(i20))then
                                                  i62 = 1
          else
                                                i62 = 2
                                           endif
              endif
              if(i02(i4(1)).ge.0.5d0)then
                                                             i21 = 1
                        else
               i21 = -1
                               endif
                                     i35 = dsqrt(i4(i49+i20-1))
                                endif
                             do i = 1,m
                                       g(i) = i4(i19+i-1)
                                      enddo
                     if(i20.le.m-i8)then
               g(i20) = g(i20) + i21 * i35
            else
                                   g(i20) = g(i20) + i21
        if(g(i20).lt.i5(i20))then
                               g(i20) = i5(i20) + 1
                        endif
                               if(g(i20).gt.i2(i20))then
                            g(i20) = i2(i20) - 1
               endif
                      endif
                           if(g(i20).lt.i5(i20)) g(i20) = i5(i20)
             if(g(i20).gt.i2(i20)) g(i20) = i2(i20)
          if(i20.gt.m-i8) g(i20) = dnint(g(i20))
                   if(i10.eq.1.and.i39.eq.1)then
                                                     i42 = -30
                               else
                          i42 = -31
                                            endif
                             return
    2      i42 = -40
            do i = 1,m
                            if(i6(i38+i-1).eq.0) goto 22
                                                       enddo
                  i93 = 1
                         i42 = -99
   22                        return
                                                                   end
                                  subroutine i022(l,x,n)
                                              implicit none
             integer n,i
                    double precision l,x
               dimension x(n)
                                     if(l.ne.l)then
          l = 1.0d16
                                        endif
                   do i=1,n
                                if(x(i).ne.x(i))then
                        x(i) = - 1.0d16
                                                            endif
                                                        enddo
                                                       end
        subroutine i07(m,i8,g,i5,i2,i4,i32,i6,i99,i18,i19,i36)
                                  implicit none
                   integer m,i8,i32,i19,i,i6,i99,i18,i36(*)
                   double precision g,i5,i2,i4,i02,i34,i35,i04
                        dimension g(m),i5(m),i2(m),i4(i32),i6(i99)
                                                         do i = 1,m-i8
                         i34 = i02(i4(1))
                                         if(i34.le.0.25d0)then
             g(i) = i4(i19+i-1)
          else
       i35 = (i2(i)-i5(i)) / dble(i6(i18)**2)
                g(i) = i4(i19+i-1) + i35 *
     &           i04(i02(i4(1)),i02(i4(1)))
                                           endif
                  enddo
                                                  do i = m-i8+1,m
                                if(i02(i4(1)).le.0.75d0)then
                            g(i) = i4(i19+i-1)
                                                  else
                                       if(i02(i4(1)).le.0.5d0)then
                                   g(i) = i4(i19+i-1) + 1.0d0
       else
                   g(i) = i4(i19+i-1) - 1.0d0
           endif
                               endif
                                   if(g(i).lt.i5(i)) g(i) = i5(i)
                 if(g(i).gt.i2(i)) g(i) = i2(i)
                      enddo
                         i34 = i02(i4(1))
                                                         do i = 1,m-i8
                                  if(g(i).lt.i5(i))then
                            if(i34.ge.0.1d0*dble(i36(2)))then
      g(i) = i5(i) + (i5(i)-g(i)) / dble(3**i36(3))
         if(g(i).gt.i2(i)) g(i) = i2(i)
                                                        else
                                                      g(i) = i5(i)
                                    endif
                                                             goto 2
                         endif
                          if(g(i).gt.i2(i))then
                           if(i34.ge.0.1d0*dble(i36(2)))then
              g(i) = i2(i) - (g(i)-i2(i)) / dble(3**i36(3))
               if(g(i).lt.i5(i)) g(i) = i5(i)
                           else
      g(i) = i2(i)
                                          endif
                                                               endif
    2                continue
                                     enddo
                                                  end
           subroutine i05(l,i48,i17,i16,i42)
             implicit none
      double precision l,i48(*),i17,i16
                       integer i42
                        if(i48(3).eq.0.0d0) return
                                           if(l.le.i48(3))then
                                   if(i17.le.i16)then
                                         i42 = 5
                                                    endif
                                                 endif
                                                         end
                                               function i04(a,b)
                  implicit none
                   double precision i04, a, b, g(30), i(30)
                                                             data g /
     &        0.260390399999d0, 0.371464399999d0, 0.459043699999d0,
     &0.534978299999d0, 0.603856999999d0, 0.668047299999d0,
     &    0.728976299999d0, 0.787597599999d0, 0.844600499999d0,
     &            0.900516699999d0, 0.955780799999d0, 1.010767799999d0,
     &   1.065818099999d0, 1.121257099999d0, 1.177410099999d0,
     & 1.234617499999d0, 1.293250299999d0, 1.353728799999d0,
     &   1.416546699999d0, 1.482303899999d0, 1.551755799999d0,
     &      1.625888099999d0, 1.706040699999d0, 1.794122699999d0,
     &          1.893018599999d0, 2.007437799999d0, 2.145966099999d0,
     &       2.327251799999d0, 2.608140199999d0, 2.908140199999d0/
                                    data i /
     &      0.207911799999d0,  0.406736699999d0,  0.587785399999d0,
     &   0.743144899999d0,  0.866025499999d0,  0.951056599999d0,
     & 0.994521999999d0,  0.994521999999d0,  0.951056599999d0,
     &          0.866025499999d0,  0.743144899999d0,  0.587785399999d0,
     &    0.406736699999d0,  0.207911799999d0, -0.016538999999d0,
     & -0.207911799999d0, -0.406736699999d0, -0.587785399999d0,
     &   -0.743144899999d0, -0.866025499999d0, -0.951056599999d0,
     &        -0.994521999999d0, -0.994521999999d0, -0.951056599999d0,
     &  -0.866025499999d0, -0.743144899999d0, -0.587785399999d0,
     &-0.406736699999d0, -0.207911799999d0, -0.107911799999d0/
      i04 = g(max(nint(a*30.0d0),1)) *  i(max(nint(b*30.0d0),1))
             end
                  subroutine i09(m,i4,i32,i27,i66,i45,i16,
     /                 i19,i14,i40,i11,i6,i99,i12,o,i28,i24,i22,
     /                     i68,i48,i36)
                                               implicit none
      integer m,i6,i32,i99,i12,o,i28,i24,i22,i19,i14,i40,
     &                                i11,i27,i66,i45,i,j,i36(*)
                        double precision i4,i16,i68,i48(*),i02
                           dimension i4(i32),i6(i99)
                              integer i96,i74,i73,i72
                                                     data i96 /0/
               if(i6(i12).le.1)then
                                                        i6(10) = 0
          i96 = 0
                                                             else
                                         i96 = i96 + 1
                            i6(10) = i36(4)**i96
       endif
          if(i48(6).lt.0.0d0.and.i6(10).ne.0)then
                                 i6(10) = nint(dabs(i48(6)))
                                                              endif
                                i74 = i36(5)
                                                i73 = 5*i36(1)
         i72 = 2
         i6(i28) = i72 * nint( i02(i4(1)) * dble(m) )
        if(i48(7).ge.2.0d0) i6(i28) = nint(i48(7))
              if(i6(i28).lt.i74) i6(i28) = i74
                     if(i6(i28).gt.i73) i6(i28) = i73
        i6(o) = nint( i02(i4(1)) * dble(i6(i28)) )
                       if(i48(8).ge.2.0d0) i6(o) = nint(i48(8))
                           if(i6(o).lt.2) i6(o) = 2
         i6(i24) = i6(i28)
     &       + i36(6) * nint(i02(i4(1)) * dble(i6(i28)))
                               i6(i22) = nint( 1.5d0 * dble(i6(o)) )
              if(i27-i12.gt.i12**i12) i6(o) = i12-i27
                  i4(i45) = i68
      if(i4(i66).le.i16.and.i4(i27).lt.i68) i4(i45) = i4(i27)
      do j=1,i6(o)
                                                do i=1,m
            i4(i19+(j-1)*m+i-1) = 1.06660d90
                                            enddo
                        i4(i40+j-1)       = 1.07770d90
       i4(i14+j-1)         = 1.08880d90
                                     i4(i11+j-1)         = 1.09990d90
                                   enddo
               end
                 subroutine i01(f,m,i8,n,i0,g,l,x,i5,i2,
     &     i42,i41,i48,i4,i32,i6,i99,
     &                    i30,i52,i50,i100,
     &             i990)
                                          implicit none
                  integer f,m,i8,n,i0,i42,i32,i6,i99,i41
          double precision g,l,x,i5,i2,i48(*),i4
            dimension g( f*m ),l( f ),x( f*n+1 ),i5(m),i2(m)
                                            dimension i4(i32),i6(i99)
                      character*60 i990
                integer i30,i52,i50,i100
                                          integer i
                        if(m.gt.4)then
        i42 = 999
       goto 701
                              endif
                               i30 = 100
                                       if(f.le.0.or.f.gt.10**6)then
                                             i42 = 101
                                                 goto 701
         endif
               if(m.le.0.or.m.gt.10**6)then
        i42 = 102
                        goto 701
                                         endif
                                                    if(i8.lt.0)then
                                 i42 = 103
                                              goto 701
           endif
                                                         if(i8.gt.m)then
                                                i42 = 104
                                       goto 701
          endif
           if(n.lt.0.or.n.gt.10**6)then
                                               i42 = 105
                                                            goto 701
                                      endif
                                                     if(i0.lt.0)then
           i42 = 106
                     goto 701
                       endif
                                      if(i0.gt.n)then
               i42 = 107
                                                    goto 701
                                         endif
                     do i=1,m
                      if(g(i).ne.g(i))then
                i42 = 201
                                                   goto 701
                        endif
                                                 if(i5(i).ne.i5(i))then
         i42 = 202
         goto 701
           endif
              if(i2(i).ne.i2(i))then
                       i42 = 203
                                      goto 701
        endif
                              if(g(i).lt.i5(i)-1.0d-4)then
                  i42 = 204
                                         goto 701
                                 endif
                                           if(g(i).gt.i2(i)+1.0d-4)then
       i42 = 205
               goto 701
               endif
                  if(i5(i).gt.i2(i)+1.0d-4)then
                                         i42 = 206
             goto 701
                                                        endif
                                         enddo
           if(i48(1).lt.0.0d0.or.i48(1).gt.1.0d6)then
                         i42 = 301
               goto 701
                                                  endif
          if(i48(2).lt.0.0d0.or.i48(2).gt.1.0d12)then
                                                  i42 = 302
                                                goto 701
                                 endif
                             if(dabs(i48(3)).gt.1.0d12)then
                                                   i42 = 303
                                            goto 701
                                                           endif
       if(i48(4).lt.0.0d0.or.i48(4).gt.1.0d6)then
                       i42 = 304
                     goto 701
                                             endif
             if(dabs(i48(5)).gt.1.0d12)then
                                                  i42 = 305
                                           goto 701
                             endif
                if(i48(6).ne.0.0d0.and.dabs(i48(6)).lt.1.0d0
     &                  .or.dabs(i48(6)).gt.1.0d12)then
                i42 = 306
                                                               goto 701
                                                  endif
                  if(i48(7).lt.0.0d0.or.i48(7).gt.1.0d8)then
        i42 = 307
           goto 701
                                                  endif
                        if(i48(8).lt.0.0d0.or.i48(8).gt.i30)then
                           i42 = 308
                                             goto 701
                                                            endif
                 if(i48(7).gt.0.0d0.and.i48(7).lt.i48(8))then
                                                          i42 = 309
                  goto 701
                                      endif
                   if(i48(7).gt.0.0d0.and.i48(8).eq.0.0d0)then
            i42 = 310
       goto 701
                                                       endif
                        if(i48(7).eq.0.0d0.and.i48(8).gt.0.0d0)then
                                 i42 = 311
                                                             goto 701
                             endif
                 if(i48(9).lt.0.0d0.or.i48(9).gt.1.0d3)then
                                                      i42 = 312
                                   goto 701
                                            endif
                                                do i=1,9
                            if(i48(i).ne.i48(i))then
                                         i42 = 313
                                            goto 701
                                         endif
                                          enddo
                        if(i41.lt.0.or.i41.gt.1)then
                             i42 = 401
                     goto 701
                   endif
                                         i52 = 2*m+n+(m+5)*i30+8
                       i50 = 31+f+m+m
            if(i32.lt.i52+5+2*m+n)then
                           i42 = 501
                                                          goto 701
                                 endif
                        if(i99.lt.i50 + 61)then
                                                 i42 = 601
                                            goto 701
                                           endif
                                          do i=1,i52+5+m+n
                                             i4(i) = 0.0d0
                                                        enddo
         do i=1,i50
                          i6(i) = 0
                    enddo
                      call i023(i6(i50+1),i990)
                i6(i50+61) = 0
                                                   do i=1,60
                      i6(i50+61) =i6(i50+61) + i6(i50+i)
                          enddo
                       if(i6(i50+61).ne.2736)then
                            i42 = 900
                 goto 701
         endif
                                      i100 = 0
                                                         do i=1,m
              if(g(i).gt.1.0d+8.or.g(i).lt.-1.0d+8)then
                                     i42 = 51
                                                 goto 702
                          endif
              if(i5(i).gt.1.0d+8.or.i5(i).lt.-1.0d+8)then
                              i42 = 52
            goto 702
                                                           endif
               if(i2(i).gt.1.0d+8.or.i2(i).lt.-1.0d+8)then
           i42 = 53
                          goto 702
                                                         endif
                  if(i5(i).eq.i2(i))then
                i42 = 71
                      goto 702
                                                  endif
                 enddo
                                                    do i = m-i8+1,m
                           if(dabs(g(i)-dnint(g(i))).gt.1.0d-4)then
                                        i42 = 61
       goto 702
               endif
                        if(dabs(i5(i)-dnint(i5(i))).gt.1.0d-4)then
                                                          i42 = 62
                                       goto 702
                                                                  endif
           if(dabs(i2(i)-dnint(i2(i))).gt.1.0d-4)then
                                  i42 = 63
                                              goto 702
                              endif
                         enddo
                         if(l(1).ne.l(1))then
             i42 = 81
                                                       goto 702
                             endif
         do i = 1,n
                   if(x(i).ne.x(i))then
                                             i42 = 82
        goto 702
                                                   endif
                                   enddo
                                       if(dabs(i48(3)).gt.1.0d+8)then
               i42 = 91
                         goto 702
                                 endif
                                      if(dabs(i48(5)).gt.1.0d+8)then
                                             i42 = 92
                                                     goto 702
                endif
                             return
  701                     continue
                                                        i41 = 1
                                          return
  702                                continue
                                           i100 = 1
                                          return
                                           end
                          subroutine i023(i67s,i15)
                                   implicit none
                              character*60 i15
                        integer i67s(*),i
                             do i = 1,60
               call alphabet(i15(i:i),i67s(i))
            enddo
                                                          end
          subroutine jfk(f,m,i8,n,i0,g,l,x,i5,i2,
     &                              i42,i41,i48,i4,i32,i6,i99,
     &                                       i26,p,i17,i990)
                                            implicit none
           integer f,m,i8,n,i0,i42,i32,i6,i99,i41
                     double precision g,l,x,i5,i2,i48(*),i4
          dimension g(f*m),l(f),x(f*n+1),i5(m),i2(m),i4(i32),i6(i99)
         character*60 i990
                                                        integer i26
                       double precision p(i26),i17(i26)
                    logical i54
               integer i30,i97,i55,i53,j,i90,
     &    i50,i77,i64,i59,i56,i58,i75,
     &                 i101,i92,i94,i,c,i46,i68,
     &          i52,i63,i37,i100,i980,i36(20)
             double precision i02,i35,i60,i70,i79,i16
                                 double precision i04
                  data i30,i97,i63,i37,i55,i53,i68,
     &        i52,i50,i77,i64,i59,i56,i58,i75,
     &                    i101,i92,i94,i100,
     &      i980,i36
     &           /0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
     &            0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0/
       data i60,i70,i16  /0.0d0,0.0d0,0.0d0/
                                                   if(i42.ge.0)then
                               i63 = 0
             i37   = 0
                           if(i48(1).le.0.0d0)then
           i16 = 1.0d-3
                                  else
         i16 = i48(1)
                                           endif
          if(i42.gt.10.and.i42.lt.100)then
                                   i42 = -3
                  i100 = 0
                                         goto 79
                            endif
                                            call i024(i36,m,i8)
                                                    i97 = 0
                             call i01(f,m,i8,n,i0,g,l,x,i5,i2,
     &             i42,i41,i48,i4,i32,i6,i99,
     &                                          i30,i52,i50,
     &       i100,
     &      i990)
          if(i42.ge.100) goto 86
                    if(i100.eq.1)then
                                 i980 = i42
                                                               i42 = 0
                                            endif
                                                     i97 = 1
                i46 = int(i48(2))
                          i92 = int(i48(4))
                                   i60    = 1.0d16
                    i70  = 1.0d16
         i59      = 1
                                               i56      = i59 + m
                   i58      = i56 + 1
       i75    = i58 + n
                 i68     = i75 + 1
                                        do i=1,m
                                         i4(i52+i59+i-1) = g(i)
                                enddo
                                        do i=1,n
                                         i4(i52+i58+i-1) = x(i)
                                           enddo
                                         i4(i52+i56) = l(1)
                                 call i014(i4(i52+i75),x,n,i0,i16)
            i77 = 0
                                      i101 = 0
                                             i53 = i46
            if(i56-i59.gt.1+i97*3)then
                      i59      = 2
                 i56      = i59 + m
                                        i58      = i56 + 1
                    i75    = i58 + n
         i68     = i75 + 1
              do i=1,m
                           i4(i52+i59+i-1) = g(i)
                                                  enddo
             do i=1,n
          i4(i52+i58+i-1) = x(i)
                                                               enddo
             i48(6) = -1.0d16
                                                  i60    = 1.0d16
                                                      i70  = 1.0d16
                                                  endif
                       i55  = 0
                               if(i4(i52+i75).gt.i16)then
                                     if(i48(5).eq.0.0d0)then
                   i4(i52+i68) = 1.0d9 + i4(i52+i56)
                                                    else
                                    i4(i52+i68) = i48(5)
                                                     endif
       else
               i4(i52+i68) = i4(i52+i56)
         endif
                           else
                                   if(i97.ne.1)then
         i42 = 701
                                                   i41 = 1
                                     return
           endif
             do c = 1,f
                             call i022(l(c),x((c-1)*n+1),n)
                                                               enddo
                                     endif
   79                                                    continue
        if(i42.eq.-300)then
               i55 = 0
          i53 = i53 + 1
                        endif
                                if(i41.eq.0)then
                                                     i54 = .false.
                                   else
                                       i54 = .true.
       endif
                      call i08(f,m,i8,n,i0,g,l,x,i5,i2,i16,
     &                  i63,i37,i55,i54,i53,
     &   i4,i32,i6,i99,i30,i4(i52+i68),
     &       i48,p,i17,i26,i36)
                                              i42 = i55
      if(i42.eq.5) goto 1
                                   if(i42.eq.801)return
                            if(i54)then
                    if(i4(i52+i75).gt.i16. and .i4(2+m+1+n).lt.
     &                                 i4(i52+i75))    goto 1
                 if(i4(i52+i75).le.i16. and .i4(2+m+1+n).le.i16
     &      .and.i4(2+m).lt.i4(i52+i56)) goto 1
                                           goto 3
                                             endif
           if(i55.eq.-3) i77 = i77 + 1
          if(i32.gt.100.and.i56.gt.i59*5)then
                                          i4(i52+i56)     = i4(2+m)
                i4(i52+i75)   = i4(2+m+1+n)
             do i=1,m
                   i4(i52+i59+i-1) = i4(2+i-1)
                                                       enddo
                        do i=1,n
            i4(i52+i58+i-1) = i4(2+m+1+i-1)
                         enddo
                                  if(i4(i52+i75).le.i16)then
                i4(i52+i68) =  i4(i52+i56)
                 endif
                     l(1) = i4(i52+i56)
                                                       do i = 1,m
                      g(i) = i4(i52+i59+i-1)
                          enddo
                                       do i = 1,n
                   x(i) = i4(i52+i58+i-1)
                                                                enddo
                                                         goto 79
                                                          endif
                         i64 = i36(7)
                                       if(i77.ge.i64)then
        i101 = i101 + 1
                       if(i4(i52+i75).gt.i16. and .i4(2+m+1+n).lt.
     &                                 i4(i52+i75))   goto 11
            if(i4(i52+i75).le.i16. and .i4(2+m+1+n).le.i16.and.
     &       i4(2+m).lt.i4(i52+i56))goto 11
                               goto 12
   11          i4(i52+i56)     = i4(2+m)
                                   i4(i52+i75)   = i4(2+m+1+n)
                                      do i=1,m
                     i4(i52+i59+i-1) = i4(2+i-1)
                                         enddo
                                      do i=1,n
              i4(i52+i58+i-1) = i4(2+m+1+i-1)
              enddo
                                       if(i4(i52+i75).le.i16)then
                 i4(i52+i68) =  i4(i52+i56)
                endif
                                   goto 13
   12  continue
   13                                   do i = 2,i52
              i4(i) = 0.0d0
                                       enddo
                        do i = 1,i50
                                                  i6(i) = 0
                                              enddo
        if(i02(i4(1)).ge.0.1d0*dble(i36(9)).or.i48(6).lt.0.0d0)then
                                           do i=1,m
                               if(m.le.m-i8) i35 = (i2(i)-i5(i))
     &                     / (i02(i4(1))*dble(10**i36(17)))
                            i79 = ((i2(i)-i5(i))-(i2(i)-i5(i))
     &                                         / dsqrt(dble(i8)+0.1d0))
     &        / dble(i36(10))
                       if(m.gt.m-i8) i35 = (i2(i)-i5(i)) / 100
                            if(i.gt.m-i8.and.i35.lt.i79)then
                     i35 = i79
                                 endif
                           if(i48(6).ne.0.0d0)then
                  if(i.le.m-i8)then
                                i35 = (i2(i)-i5(i)) / dabs(i48(6))
                                              else
                     i35 = 1.0d0 / dsqrt(dabs(i48(6)))
                                                     endif
                                              endif
                          g(i) = i4(i52+i59+i-1) + i35 *
     &                         i04(i02(i4(1)),i02(i4(1)))
                                           if(g(i).lt.i5(i))then
                                        g(i)=i5(i)+(i5(i)-g(i)) / 2.0d0
                         endif
                                           if(g(i).gt.i2(i))then
          g(i)=i2(i)-(g(i)-i2(i)) / 2.0d0
                                                                   endif
                                                 if(g(i).lt.i5(i))then
                                                  g(i)=i5(i)
                                  endif
                                                if(g(i).gt.i2(i))then
                g(i)=i2(i)
                       endif
                        if(i.gt.m-i8)then
                                               g(i)=dnint(g(i))
                                   endif
                      if(i.gt.2*(i59*2))then
                                             do j=1,i32
                                 i4(j) = i02(i4(1))
                            enddo
                        do j=1,i99
                i6(j) = nint(i02(i4(1)))
                                                              enddo
                                   endif
                            enddo
                                                      else
                                          do i=1,m
                  g(i) = i5(i) + i02(i4(1)) * (i2(i)-i5(i))
                         if (i.gt.m-i8) g(i) = dnint(g(i))
                                             enddo
                              endif
                        i42 = -300
                          i77 = 0
          if(i92.gt.0)then
                                       if(i60.eq.1.0d16)then
                 i94 = 0
                                      i60   = i4(i52+i56)
                                    i70 = i4(i52+i75)
                                     else
                              if(i4(i52+i75).le.i70)then
                                         if(i70.le.i16)then
                                         if(i4(i52+i56).lt.
     &    i60-dabs(i60/1.0d6))then
                       i60   = i4(i52+i56)
                   i70 = i4(i52+i75)
                                                    i94 = 0
                                              else
                             i94 = i94 + 1
                                                     goto 76
                                                   endif
                                                          else
           i94 = 0
                     i60   = i4(i52+i56)
                  i70 = i4(i52+i75)
                                                     endif
                                 else
               i94 = i94 + 1
                               goto 76
           endif
                                                  endif
   76              continue
                                 if(i94.ge.i92)then
                             if(i4(i52+i75).le.i16)then
                          i42 = 3
                                    else
                                    i42 = 4
                            endif
            goto 3
                                                                  endif
                                  endif
                       endif
                                    if(i100.eq.1)then
                            i42 = i980
                  endif
    4                        return
    1   i4(i52+i56)     = i4(2+m)
                           i4(i52+i75)   = i4(2+m+1+n)
           do i=1,m
                        i4(i52+i59+i-1) = i4(2+i-1)
                    enddo
                                     do i=1,n
              i4(i52+i58+i-1) = i4(2+m+1+i-1)
                                                   enddo
                        if(i4(i52+i75).le.i16)then
             i4(i52+i68) =  i4(i52+i56)
                   endif
    3                                            l(1) = i4(i52+i56)
                              do i = 1,m
                   g(i) = i4(i52+i59+i-1)
                                                  enddo
                                 do i = 1,n
                x(i) = i4(i52+i58+i-1)
                                  enddo
          if(i42.ne.3.and.i42.ne.4.and.i42.ne.5)then
                                      if(i4(i52+i75).le.i16)then
                                        i42 = 1
                                  else
                                              i42 = 2
                                                    endif
                                            endif
                                                   i41 = 1
   86                                                     continue
       if(i42.eq.501.or.i42.eq.601) goto 4
                          i90 = i52+5+m+n
                           do i = 1,m
                                     if( g(i).gt.i2(i)+1.0d-6 )then
                                   i4(i90+i) = 91.0d0
                                           goto 87
                                                            endif
                                       if( g(i).lt.i5(i)-1.0d-6 )then
             i4(i90+i) = 92.0d0
                                              goto 87
                                                    endif
                   if( i5(i).gt.i2(i) )then
                          i4(i90+i) = 93.0d0
                                     goto 87
                        endif
      if( i5(i).eq.i2(i) )then
                                          i4(i90+i) = 90.0d0
               goto 87
                                                 endif
              if( dabs(g(i)-i5(i)) .le. (i2(i)-i5(i))/1000.0d0 )then
                  i4(i90+i) = 0.0d0
                             goto 87
                           endif
      if( dabs(g(i)-i2(i)) .le. (i2(i)-i5(i))/1000.0d0 )then
                                      i4(i90+i) = 22.0d0
            goto 87
                                          endif
                              do j = 1,21
                   if( g(i) .le. i5(i) + j * (i2(i)-i5(i))/21.0d0)then
                                           i4(i90+i) = dble(j)
                                                           goto 87
                            endif
                                                           enddo
   87         continue
            enddo
             goto 4
                             end
                    subroutine i03(m,i8,i4,i32,i49,i19,
     &   i6,i99,o,i18,i1,i36)
                                             implicit none
           integer m,i8,i6,i18,o,i32,i99,i49,i19,i1,i,j,i36(*)
                double precision i4,i33,i320,i76,i80,i79
          dimension i4(i32),i6(i99)
           i76  = dsqrt(dble(i6(i18)))
          i80 = (0.1d0*dble(i36(18))) / i76
       i79 = (1.0d0-1.0d0/dsqrt(dble(i8)+0.1d0)) / dble(i36(10))
                                                 do i=1,m
                      i33 = i4(i19+i-1)
                i320 = i4(i19+i-1)
                      do j=2,i6(o)
                                if(i4(i19+(j-1)*m+i-1).gt.i33)then
                     i33 = i4(i19+(j-1)*m+i-1)
                    endif
                    if(i4(i19+(j-1)*m+i-1).lt.i320)then
        i320 = i4(i19+(j-1)*m+i-1)
                                                  endif
                   enddo
                            i4(i49+i-1) = (i33-i320)/i76
                                 if(i1.gt.i49+i18)then
            i4(i19-(j-1)*m+i-1) = dble(i320)/dble(i33)
                      i33 = i4(i19+i-1)
                 i320 = i4(i19+i-1)
                                                   do j=2,i6(o)
                  if(i4(i19+(j-1)*m+i-1).gt.i33)then
              i33 = i4(i19+(j-1)*m+i-1)
                                endif
       if(i4(i19+(j-1)*m+i-1).lt.i320)then
                                i320 = i4(i19+(j-1)*m+i-1)
             endif
                                     enddo
                      endif
                                  if(i4(i49+i-1).lt.
     &   dabs(i4(i19+i-1))/(dble(10**i36(19))*dble(i6(i18))) )then
                                     i4(i49+i-1)  = dabs(i4(i19+i-1))
     &                             / (dble(10**i36(19))*dble(i6(i18)))
        endif
                                        if(i.gt.m-i8)then
                             if(i4(i49+i-1).lt.i80)then
            i4(i49+i-1) = i80
                    endif
        if(i4(i49+i-1).lt.i79)then
                        i4(i49+i-1) = i79
                                                 endif
                 endif
       enddo
                                                               end
           function i02(s)
                                              implicit none
                                 double precision i02,s,a,b,c,d
              data  a,b,c,d  /0.0d0,0.0d0,0.0d0,0.0d0/
                                if(s.eq.1.2d3)then
                                a =  0.241305836148d+00
              b =  0.502199827430d+00
        c =  0.150532162816d+00
                                          d =  0.278319003121d+00
          endif
                                      if(b.ge.0.5d0)then
                                   s = a + b + c
                                                                   else
                                                   s = a + b + c + d
                 endif
                                             if(s.ge.1.0d0)then
                                                 if(s.ge.2.0d0)then
                                      s = s - 2.0d0
       else
                  s = s - 1.0d0
                                                          endif
            endif
                                                            a = b
                                                        b = c
                                      c = s
                                                                 s = a
         i02 = s
                       end
        subroutine i011(m,i8,o,g,i5,i2,i4,i32,i19,i49,i9,i36)
                                       implicit none
          integer m,i8,t,o,i32,i19,i49,i9,i,j,i36(*)
                     double precision i4,g,i5,i2,i02,i34,i04
                            dimension i4(i32),g(m),i5(m),i2(m)
                                do i=1,m
                                        i34 = i02(i4(1))
                                                     do j=1,o
                               if(i34.le.i4(i9+j-1)) goto 1
                                                            enddo
    1      t = i-1
               g(i) = i4(i19+(j-2)*m+t) + i4(i49+t) *
     &            i04(i34,i02(i4(1)))
                                         if(g(i).lt.i5(i))then
                    if(i34.ge.0.1d0*dble(i36(2)))then
           g(i) = i5(i) + (i5(i)-g(i)) / dble(3**i36(3))
                        if(g(i).gt.i2(i)) g(i) = i2(i)
                   else
                                                          g(i) = i5(i)
                                               endif
                                         goto 2
                                                               endif
                                               if(g(i).gt.i2(i))then
           if(i34.ge.0.1d0*dble(i36(2)))then
                     g(i) = i2(i) - (g(i)-i2(i)) / dble(3**i36(3))
                                  if(g(i).lt.i5(i)) g(i) = i5(i)
                              else
                                        g(i) = i2(i)
                                                       endif
                                       endif
    2                                          continue
                                                      enddo
                                              do i=m-i8+1,m
                                                    g(i) = dnint(g(i))
                   enddo
                                                end
           subroutine i012(i6,i99,i18,i29,i28,i24,i22,i27,i31)
                              implicit none
                         integer i6,i99,i18,i29,i28,i24,i22
                                                      integer i27,i31
                                              dimension i6(i99)
           if(i6(i18).eq.1.and.i6(i22).eq.1)then
                                       i6(i29) = i6(i24)
                                                  else
                     i6(i29) = i6(i28)
                                endif
                   if(i27-i31.gt.i24-i29) i6(1) = 0
          if(i6(i18).le.i6(i22).and.i6(i22).gt.1)then
        i6(i29) = i6(i28) + (i6(i24)-i6(i28)) *
     &               int( dble((i6(i18)-1)) / dble((i6(i22)-1)) )
                                                                 endif
            if(i6(i18).gt.i6(i22).and.i6(i18).lt.2*i6(i22))then
                    i6(i29) = 2 * ( i6(i24) + (i6(i28)-i6(i24)) *
     &      int( dble(i6(i18)) / dble(2*i6(i22)) ) )
                       endif
                              end
                                 subroutine i013(m,i8,g,i5,i2,i4,i32)
                                     implicit none
                                   integer m,i8,i32,i
                                   double precision g,i5,i2,i4,i02
                             dimension g(m),i5(m),i2(m),i4(i32)
                                                       do i=1,m
        g(i) = i5(i) + i02(i4(1)) * (i2(i)-i5(i))
                  if (i.gt.m-i8) g(i) = dnint(g(i))
                         enddo
                                                      end
         subroutine i020(m,o,i4,i32,i19,i14,i40,i11,g,l,i17,p)
                  implicit none
                integer m,o,i32,i19,i14,i40,i11,pface,i,j
          double precision g,l,i17,p,i4
                                        dimension g(m),i4(i32)
                   if(p.ge.i4(i11+o-1))return
                                                            pface = 0
          do i=1,o
         if(p.le.i4(i11+o-i))then
          pface=o-i+1
                   else
                          goto 567
                                             endif
                  enddo
  567                                                do j=1,o-pface
                                                 do i=1,m
                   i4(i19+(o-j)*m+i-1) = i4(i19+(o-j-1)*m+i-1)
                                                                enddo
                                         i4(i14+o-j)     = i4(i14+o-j-1)
              i4(i40+o-j)   = i4(i40+o-j-1)
                                      i4(i11+o-j)     = i4(i11+o-j-1)
               enddo
       do i=1,m
         i4(i19+(pface-1)*m+i-1)  = g(i)
                     enddo
                   i4(i14+pface-1)          = l
                              i4(i40+pface-1)        = i17
                i4(i11+pface-1)          = p
                                                                    end
                  subroutine i08(f,m,i8,n,i0,g,l,x,i5,i2,i16,
     &                   i69,i25,i42,i41,i46,i4,i32,i6,i99,
     &               i30,i68,i48,p,i17,i26,i36)
                                                          implicit none
                       logical i41
      integer i31,i27,i23,i66,i45,i19,i14,i40,i11,i,j,c,i30,
     & i93,i9,w,i49,i1,i7,i170,o,i12,i18,i29,i13,i28,i24,i22,
     &         f,m,i8,n,i0,i69,i25,i42,i46,i32,i6,i99,i26
             double precision g,l,x,i5,i2,i16,i4,i68,i48(*),i02
        dimension g( f*m ),l( f ),x( f*n+1 ),i5(m),i2(m),i4(i32),i6(i99)
                             double precision p(i26),i17(i26)
                                  integer i65, i78, i36(*)
      data i93,i31,i27,i23,i66,i45,i19,i14,i40,i11,i22,
     &  i9,w,i49,i1,i7,i170,o,i12,i18,i29,i13,i28,i24,i78
     &   /0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0/
                                 if(i42.ge.0)then
       call i017(m,n,i31,i27,i23,i66,i45,i30,
     &              i19,i14,i40,i11,w,i49,i9,i1,i7,i170)
         call i016(i12,i29,o,i18,i13,i28,i24,i22)
                                i4(1) = 1.2d3
                 do i=1,i46
                                               i4(1) = i02(i4(1))
                                       enddo
                call i014(i17( 1 ),x( 1 ),n,i0,i16)
                           i4(i27)   = l( 1 )
                    i4(i66) = i17( 1 )
           do i=1,m
                                            i4(i31+i-1) = g( i )
                    enddo
                                                        do i=1,n
                                                i4(i23+i-1) = x( i )
                                                   enddo
            call i05(l( 1 ),i48,i17( 1) ,i16,i42)
                                  if(i42.eq.5) return
                      if(i6(f+2*m+32).ne.74-i6(m+m+f+33))then
                 goto 999
              endif
                                                           goto 101
                   endif
                                                             do c = 1,f
        call i014(i17( c ),x( (c-1)*n+1 ),n,i0,i16)
      if(n.gt.0) call i010(p( c ),l( c ),i17( c ),i4(i45),i16)
           if(n.eq.0) p( c ) = l( c )
             if(i42.gt.-30.or.i42.lt.-40)then
             call i020(m,i6(o),i4,i32,i19,i14,i40,i11,
     &                     g( (c-1)*m+1 ),l( c ),i17( c ),p( c ))
                              endif
             if(i42.le.-30.and.i42.ge.-40)then
           call i019(f,c,m,i4,i32,i6,i99,i19,i14,i40,i11,
     & g( (c-1)*m+1 ),l( c ),i17( c ),p( c ),i36)
                                                  endif
           if(i17( c ).lt.i4(i66)) goto 123
           if(i17( c ).eq.i4(i66).and.l( c ).lt.i4(i27)) goto 123
                                     goto 100
  123     i4(i27)   = l( c )
      i4(i66) = i17( c )
                                           do i=1,m
                                 i4(i31+i-1) = g( (c-1)*m+i )
           enddo
               do i=1,n
                                         i4(i23+i-1) = x( (c-1)*n+i )
                 enddo
      call i05(l( c ),i48,i17( c ),i16,i42)
             if(i42.eq.5) return
  100                      continue
                                       enddo
  101    if(i41)goto 999
                                     if(i42.le.-90)then
          if(i4(i170).gt.i16.and.i4(i40).lt.i4(i170))goto 81
                         if(i4(i170).le.i16.and.i4(i40).le.i16.
     &                         and.i4(i14).lt.i4(i7))goto 81
                              goto 82
   81                               i6(11) = 1
                                           goto 83
   82                  i6(11) = 0
   83            continue
             endif
                                 if(i42.eq.-10)then
            if(i4(i11).lt.i4(i1)-dabs(i4(i1))/dble(2*i36(11))) goto 84
                        i6(13) = 0
                                      goto 85
   84                                            i6(13) = 1
                        i4(i1) = i4(i11)
   85                        continue
                                      endif
              i65 = min(10+m*i36(8),5*i36(20))
                                             if(i6(i18).ge.i65)then
                                                i42 = -95
                                                  endif
              i69 = 0
                           i25 = 0
            if(i6(i18).eq.1.and.i42.eq.-10)then
                                             i78 = 0
                        endif
                               if(i42.eq.-10) i78 = i78 + 1
      if(i31+i31.lt.i27-i31) goto 101
      if(i41) goto 3
                                          if(i42.eq. -1)      goto 13
                                                   if(i42.eq. -2)then
                                            i42 = -1
                           goto 13
                                                   endif
                                         if(i42.eq. -3)then
                            if(i6(i13).ge.i6(i29))then
                                                 i42 = -30
                                                       goto 14
                                                              endif
                     i42 = -1
                                                 goto 13
                                                             endif
                                   if(i42.eq.-30)then
                                                            i42 = -31
                    goto 14
                                         endif
                                           if(i42.le.-31.and.
     &                  i42.ge.-39)then
                               i42 = i42
                                            goto 14
                                   endif
                    if(i42.eq.-40)then
                             i42 = -2
                 goto 12
                                                             endif
                             if(i42.eq.-10)then
                       i42 = -30
                                goto 14
                     endif
      if(i42.le.-90)then
                                             i42 = -3
                                                              goto 11
                            endif
                                       if(i42.eq.  0)then
        i42 = -3
                                                       goto 11
                                          endif
   11         i6(i12) = i6(i12)+1
                       i6(i18) = 0
       call i09(m,i4,i32,i27,i66,i45,i16,
     &          i19,i14,i40,i11,i6,i99,i12,o,i28,i24,i22,
     &                                         i68,i48,i36)
            call i021(i6(o),i4,i32,w)
              if(i23-3*i31.gt.i6(i18)+1)then
                    i4(i9)=0.0d0
                                          do j=1,i6(o)
            i4(i9+j) = i4(i9+j-1) + i4(w+j-1)
                                                        enddo
                                     i4(i7)   = i4(i27)
                       i4(i170) = i4(i66)
                                if(i6(i12).eq.1)then
           if(n.gt.0) call i010(p( 1 ),l( 1 ),i17( 1 ),i4(i45),i16)
                if(n.eq.0) p( 1 ) = l( 1 )
        call i020(m,i6(o),i4,i32,i19,i14,i40,i11,
     &                                   g( 1 ),l( 1 ),i17( 1 ),p( 1 ))
                           endif
               i6(o) = i42
                                  do c = 1,f
       i6(i13) = i6(i13) + 1
                                      if(i6(i18).eq.1)then
               if(i6(10).le.1)then
                                if(i48(6).ne.0.0d0)then
                    call i018(m,i8,g( (c-1)*m+1 ),i5,i2,
     &  i4,i32,i31,dabs(i48(6)),i36)
                                                              else
                call i013(m,i8,g( (c-1)*m+1 ),i5,i2,i4,i32)
              endif
                                    else
        call i06(m,i8,i4,i32,i6,i99,i31,
     &                 g( (c-1)*m+1 ),i5,i2,dabs(i48(6)),i36)
                  endif
                                             endif
         if(i6(i18).gt.1)call i011(m,i8,i6(o),g( (c-1)*m+1 ),i5,i2,
     &   i4,i32,i19,i49,i9,i36)
        if(i18.lt.i1-i49) call i013(m,i8,g( (c-1)*m+1 ),i5,i2,i4,i32)
                                     enddo
                    endif
                            i4(i9)=0.0d0
             do j=1,i6(o)
                               i4(i9+j) = i4(i9+j-1) + i4(w+j-1)
                         enddo
                                                  i4(i7)   = i4(i27)
       i4(i170) = i4(i66)
                                                           i93 = 0
                                                 if(i6(i12).eq.1)then
           if(n.gt.0) call i010(p( 1 ),l( 1 ),i17( 1 ),i4(i45),i16)
                 if(n.eq.0) p( 1 ) = l( 1 )
                 call i020(m,i6(o),i4,i32,i19,i14,i40,i11,
     &                   g( 1 ),l( 1 ),i17( 1 ),p( 1 ))
             endif
   12                                            i6(i18) = i6(i18) + 1
                                       i6(i13) = 0
         call i03(m,i8,i4,i32,i49,i19,i6,i99,o,i18,i1,i36)
                                 if(i48(7).gt.0.0d0)then
                                    i6(i29) = i6(i28)
                               else
           call i012(i6,i99,i18,i29,i28,i24,i22,i27,i31)
                                 endif
              if(i6(i18).eq.1) i4(i1) = 1.0d16
         if(i6(i18).gt.1) i4(i1) = i4(i11)
   13                                             do c = 1,f
           i6(i13) = i6(i13) + 1
                                       if(i6(i18).eq.1)then
                                   if(i6(10).le.1)then
             if(i48(6).ne.0.0d0)then
            call i018(m,i8,g( (c-1)*m+1 ),i5,i2,
     & i4,i32,i31,dabs(i48(6)),i36)
                                else
                     call i013(m,i8,g( (c-1)*m+1 ),i5,i2,i4,i32)
            endif
                                                              else
                          call i06(m,i8,i4,i32,i6,i99,i31,
     &          g( (c-1)*m+1 ),i5,i2,dabs(i48(6)),i36)
                                                 endif
                                                      endif
          if(i6(i18).gt.1)call i011(m,i8,i6(o),g( (c-1)*m+1 ),i5,i2,
     &                                   i4,i32,i19,i49,i9,i36)
        if(i18.lt.i1-i49) call i013(m,i8,g( (c-1)*m+1 ),i5,i2,i4,i32)
                 enddo
                 if(i6(i13).ge.i6(i29).and.i42.ne.-3) i42 = -10
    3                                                       return
   14                                continue
                                 if(i6(13).eq.1.or.i6(i18).eq.1)then
       i42 = -2
               goto 12
          else
             if(i42.lt.-30.and.i6(31).eq.1)then
                                i42 = -2
                                                   goto 12
                                  endif
                        if(i42.eq.-39)then
                                       i93 = 1
                             i42    = -99
               goto 101
                   endif
                                                    do c = 1,f
               if(f.gt.1) i6(31) = 0
                    call i015(f,m,i8,g( (c-1)*m+1 ),i5,i2,i19,i49,i4,
     &                                  i32,i6,i99,i42,i93,i36)
                            if(i42.eq.-30.and.f.gt.1) i42 = -31
                     if((i14-i19)/i30-i18.gt.0)then
                                                          i42 = -23
                  goto 14
                        endif
                     if(i93.eq.1.and.c.gt.1)then
                                call i07(m,i8,g( (c-1)*m+1 ),i5,i2,
     &                          i4,i32,i6,i99,i18,i19,i36)
                                                i93 = 0
            i42 = -39
                       endif
                 enddo
                             if(i93.eq.1) goto 101
        goto 3
                                     endif
  999      l( 1 )   = i4(i27)
                      i17( 1 ) = i4(i66)
                                                               do i=1,m
                           g(i) = i4(i31+i-1)
                                     enddo
                                                             do j=1,n
                                    x(j) = i4(i23+j-1)
                                                        enddo
           if(i17( 1 ).le.i16)then
                       i42 = 0
                    else
                                          i42 = 1
               endif
                                                   end
        subroutine i06(m,i8,i4,i32,i6,i99,i31,g,i5,i2,i47,i36)
          implicit none
              integer m,i8,i32,i31,i6,i99,i,i36(*)
                double precision i4,g,i5,i2,i47,i35,i02,i04
               double precision i34
       dimension i4(i32),i6(i99),g(m),i5(m),i2(m)
          do i=1,m
               i35 = (i2(i)-i5(i)) / dble(i6(10))
             if(i.gt.m-i8.and.i35.lt.0.1d0*dble(i36(13)))then
                                     i35 = 0.1d0*dble(i36(13))
                                     endif
                          if(i47.gt.0.0d0)then
                     if(i35.gt.(i2(i)-i5(i))/i47)then
                        i35 = (i2(i)-i5(i))/i47
                             endif
                     if(i.gt.m-i8)then
          if(i35.lt.1.0d0/dsqrt(i47))then
                                      i35 = 1.0d0 / dsqrt(i47)
                     endif
                                                              endif
                                                                  endif
                                                 i34 = i02(i4(1))
                         g(i) = i4(i31+i-1) + i35 *
     &  i04(i34,i02(i4(1)))
                         if(g(i).lt.i5(i))then
                                    if(i34.ge.0.1d0*dble(i36(2)))then
               g(i) = i5(i) + (i5(i)-g(i)) / dble(3**i36(3))
                            if(g(i).gt.i2(i)) g(i) = i2(i)
          else
                                                    g(i) = i5(i)
                                     endif
               goto 2
                     endif
                                               if(g(i).gt.i2(i))then
               if(i34.ge.0.1d0*dble(i36(2)))then
         g(i) = i2(i) - (g(i)-i2(i)) / dble(3**i36(3))
                   if(g(i).lt.i5(i)) g(i) = i5(i)
                                                        else
                        g(i) = i2(i)
                    endif
                             endif
    2                                 if(i.gt.m-i8) g(i) = dnint(g(i))
                                        enddo
                                           end
           subroutine i019(f,c,m,i4,i32,i6,i99,i19,i14,i40,i11,
     &                                  g,l,i17,p,i36)
                                                implicit none
                    integer f,c,m,i32,i99,i6,i19,i14,i40,i11,i,i36(*)
             double precision g,l,i17,p,i4
                               dimension g(m),i4(i32),i6(i99)
                          if(i17.le.0.0d0.and.i4(i40).le.0.0d0)then
         if(l.ge. i4(i14) - dabs(i4(i14)) / dble(10**i36(12)) )then
                                               i6(31 + c) = 0
                goto 1
                endif
             else
         if(p.ge. i4(i11) - dabs(i4(i11)) / dble(10**i36(12)) )then
                                          i6(31 + c) = 0
                   goto 1
         endif
                                                                endif
                   do i = 1,m
                                                    i4(i19+i-1) = g(i)
              enddo
        i4(i40) = i17
             i4(i14)   = l
                      i4(i11)   = p
                              i6(31 + c) = 1
   1                                                if(c.eq.f)then
                                           i6(31) = 0
                                               do i = 1,f
                                              i6(31) = i6(31) + i6(31+i)
             enddo
                            if(i6(31).gt.1) i6(31) = 1
                                               endif
         return
      end
                 subroutine i016(i12,i29,o,i18,i13,i28,i24,i22)
                                                implicit none
             integer i12,i29,o,i18,i13,i28,i24,i22
                                   o      = 1
                                      i12    = 2
          i29   = 3
                                                        i18    = 4
                                    i13    = 5
                  i28   = 6
                          i24   = 7
            i22   = 8
                                                  end
        subroutine i017(m,n,i31,i27,i23,i66,i45,i30,
     &                  i19,i14,i40,i11,w,i49,i9,i1,i7,i170)
                                   implicit none
             integer m,n
         integer i31,i27,i23,i66
                                               integer i45
                                                            integer i30
                              integer i19,i14,i40,i11
          integer w,i49
                           integer i9
        integer i1,i7,i170
                                i31        = 2
                                        i27        = i31        + m
             i23        = i27        + 1
                   i66      = i23        + n
                       i45       = i66      + 1
       i19         = i45       + 1
                               i14         = i19         + m * i30
                 i40       = i14         + i30
                                     i11         = i40       + i30
         i9          = i11         + i30
       w           = i9          + i30 + 1
       i49       = w           + i30
           i1          = i49       + m
             i7          = i1          + 1
                   i170        = i7          + 1
                                                end
        subroutine i018(m,i8,g,i5,i2,i4,i32,i31,i47,i36)
                                                        implicit none
                                  integer m,i8,i32,i31,i,i36(*)
           double precision g,i5,i2,i4,i47,i35,i02,i04,i34
           dimension g(m),i5(m),i2(m),i4(i32)
                            do i=1,m
              i35 = (i2(i)-i5(i)) / i47
                                          if(i.gt.m-i8)then
                if(i35.lt.1.0d0/dsqrt(i47))then
                                             i35 = 1.0d0 / dsqrt(i47)
          endif
             endif
              i34 = i02(i4(1))
                   g(i) = i4(i31+i-1) + i35 *
     &                        i04(i34,i02(i4(1)))
            if(g(i).lt.i5(i))then
                                if(i34.ge.0.1d0*dble(i36(2)))then
                  g(i) = i5(i) + (i5(i)-g(i)) / dble(3**i36(3))
                if(g(i).gt.i2(i)) g(i) = i2(i)
                                 else
                       g(i) = i5(i)
        endif
                                                      goto 2
               endif
                                      if(g(i).gt.i2(i))then
                          if(i34.ge.0.1d0*dble(i36(2)))then
           g(i) = i2(i) - (g(i)-i2(i)) / dble(3**i36(3))
                 if(g(i).lt.i5(i)) g(i) = i5(i)
                 else
                                            g(i) = i2(i)
                  endif
                                                       endif
    2                         if(i.gt.m-i8) g(i) = dnint(g(i))
                                 enddo
                                                  end
                                 subroutine i010(p,l,i17,i45,i16)
                                            implicit none
                              double precision p,l,i17,i45,i16,i61
                  i61 = l - i45
                           if (l.le.i45.and.i17.le.i16) then
                                                        p = i61
                                       return
                                                          else
       if (l.le.i45) then
             p = i17
                                                                 return
             else
                                              if (i17.le.i61) then
      p = i61 + i17**2/(2.0d0*i61) - i17/2.0d0
                                                       else
                     p = i17 + i61**2/(2.0d0*i17) - i61/2.0d0
                       endif
                                      endif
                                        endif
           end
                                       subroutine i021(o,i4,i32,w)
                                 implicit none
              integer o,i32,j,i57,w
                    double precision i4
                                                       dimension i4(i32)
                                 i57 = 0
                                   do j=1,o
                                                      i57 = i57 + j
                                                                enddo
               do j=1,o
                               i4(w+j-1) = dble(o-j+1)/dble(i57)
                                                              enddo
                                           end
       subroutine i014(i17,x,n,i0,i16)
                    implicit none
                                                    integer n,i0,i
             double precision x,i17,i16
         dimension x(n)
            i17 = 0.0d0
                       if(n.eq.0)return
                                    do i=1,n
                                     if (x(i).lt.-i16) i17 = i17 - x(i)
                                        enddo
                                                  do i=1,i0
       if (x(i).gt. i16) i17 = i17 + x(i)
          enddo
                                                           end
                                       subroutine i024(i36,m,i8)
                                    implicit none
         integer i36(*),m,i8
                   if(m-i8.gt.0)then
                                     i36(   1) =                 18
              i36(   2) =                  6
          i36(   3) =                  4
               i36(   4) =                  4
                            i36(   5) =                  5
      i36(   6) =                 18
                           i36(   7) =                 13
              i36(   8) =                  5
                         i36(   9) =                  2
                  i36(  10) =                  7
                                   i36(  11) =                 15
                              i36(  12) =                  8
                          i36(  13) =                  9
                                  i36(  14) =                  7
                                  i36(  15) =                 12
                                     i36(  16) =                  5
             i36(  17) =                  5
                           i36(  18) =                  5
                              i36(  19) =                  8
             i36(  20) =                 37
                 else
       i36(   1) =                  4
                    i36(   2) =                  4
                           i36(   3) =                  3
       i36(   4) =                  3
                                 i36(   5) =                  2
                    i36(   6) =                 49
                                       i36(   7) =                 37
          i36(   8) =                   7
            i36(   9) =                  0
              i36(  10) =                  8
                        i36(  11) =                  3
                                      i36(  12) =                  6
                   i36(  13) =                  4
                         i36(  14) =                  3
                                    i36(  15) =                  3
                        i36(  16) =                  2
                            i36(  17) =                 10
                              i36(  18) =                 12
                         i36(  19) =                  4
                          i36(  20) =                 58
                                     endif
                           end
      subroutine alphabet(a,b)
      implicit none
      character*1 a
      integer b
      b = 0
      if(a(1:1).eq.'A') b = 52
      if(a(1:1).eq.'B') b = 28
      if(a(1:1).eq.'C') b = 49
      if(a(1:1).eq.'D') b = 30
      if(a(1:1).eq.'E') b = 31
      if(a(1:1).eq.'F') b = 32
      if(a(1:1).eq.'G') b = 33
      if(a(1:1).eq.'H') b = 34
      if(a(1:1).eq.'I') b = 35
      if(a(1:1).eq.'J') b = 36
      if(a(1:1).eq.'K') b = 37
      if(a(1:1).eq.'L') b = 38
      if(a(1:1).eq.'M') b = 39
      if(a(1:1).eq.'N') b = 40
      if(a(1:1).eq.'O') b = 41
      if(a(1:1).eq.'P') b = 42
      if(a(1:1).eq.'Q') b = 43
      if(a(1:1).eq.'R') b = 44
      if(a(1:1).eq.'S') b = 45
      if(a(1:1).eq.'T') b = 46
      if(a(1:1).eq.'U') b = 47
      if(a(1:1).eq.'V') b = 48
      if(a(1:1).eq.'W') b = 29
      if(a(1:1).eq.'X') b = 50
      if(a(1:1).eq.'Y') b = 51
      if(a(1:1).eq.'Z') b = 27
      if(a(1:1).eq.'0') b = 53
      if(a(1:1).eq.'1') b = 54
      if(a(1:1).eq.'2') b = 55
      if(a(1:1).eq.'3') b = 56
      if(a(1:1).eq.'4') b = 57
      if(a(1:1).eq.'5') b = 58
      if(a(1:1).eq.'6') b = 59
      if(a(1:1).eq.'7') b = 60
      if(a(1:1).eq.'8') b = 61
      if(a(1:1).eq.'9') b = 62
      if(a(1:1).eq.'a') b = 23
      if(a(1:1).eq.'b') b = 2
      if(a(1:1).eq.'c') b = 3
      if(a(1:1).eq.'d') b = 16
      if(a(1:1).eq.'e') b = 5
      if(a(1:1).eq.'f') b = 13
      if(a(1:1).eq.'g') b = 7
      if(a(1:1).eq.'h') b = 8
      if(a(1:1).eq.'i') b = 9
      if(a(1:1).eq.'j') b = 10
      if(a(1:1).eq.'k') b = 11
      if(a(1:1).eq.'l') b = 12
      if(a(1:1).eq.'m') b = 6
      if(a(1:1).eq.'n') b = 14
      if(a(1:1).eq.'o') b = 15
      if(a(1:1).eq.'p') b = 4
      if(a(1:1).eq.'q') b = 17
      if(a(1:1).eq.'r') b = 18
      if(a(1:1).eq.'s') b = 19
      if(a(1:1).eq.'t') b = 20
      if(a(1:1).eq.'u') b = 21
      if(a(1:1).eq.'v') b = 22
      if(a(1:1).eq.'w') b = 1
      if(a(1:1).eq.'x') b = 24
      if(a(1:1).eq.'y') b = 25
      if(a(1:1).eq.'z') b = 26
      if(a(1:1).eq.'_') b = 64
      if(a(1:1).eq.'(') b = 65
      if(a(1:1).eq.')') b = 66
      if(a(1:1).eq.'+') b = 67
      if(a(1:1).eq.'-') b = 68
      if(a(1:1).eq.'&') b = 69
      if(a(1:1).eq.'.') b = 70
      if(a(1:1).eq.',') b = 71
      if(a(1:1).eq.':') b = 72
      if(a(1:1).eq.';') b = 73
      if(a(1:1).eq.'*') b = 74
      if(a(1:1).eq.'=') b = 75
      if(a(1:1).eq.'/') b = 76
      if(a(1:1).eq.'!') b = 80
      if(a(1:1).eq.'[') b = 83
      if(a(1:1).eq.']') b = 84
      end

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C     This subroutine handles all printing commands for MIDACO.
C     Note that this subroutine is called independently from MIDACO and
C     MIDACO itself does not include any print commands (due to
C     compiler portability and robustness).
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      SUBROUTINE MIDACO_PRINT(C,PRINTEVAL,SAVE2FILE,IFLAG,ISTOP,F,G,
     &           X,XL,XU,N,NI,M,ME,RW,LRW,MAXEVAL,MAXTIME,PARAM,P,KEY)
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      IMPLICIT NONE
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      INTEGER C,PRINTEVAL,SAVE2FILE,IFLAG,ISTOP,EVAL,N,NI,M,ME,LRW
      INTEGER TIC,Q,I,IOUT1,IOUT2,MAXEVAL,MAXTIME,KMAX,UPDATE
      INTEGER KF,KG,KRES,KX,WF,WG,WRES,WX,KBEST,WBEST,P
      DOUBLE PRECISION TNOW,TSTART,TMAX,RW(LRW),BESTF,BESTR,F,G(M),X(N)
      DOUBLE PRECISION XL(N),XU(N),ACC,PARAM(9),DUMMY_F,DUMMY_VIO
C     Increase size, if problems with N > 1000 are solved
      DOUBLE PRECISION BESTX(1000),BESTG(1000)
      CHARACTER*60 KEY
      DATA KF,KG,KRES,KX,WF,WG,WRES,WX /0,0,0,0,0,0,0,0/
      DATA KBEST,WBEST,TIC,Q,EVAL,IOUT1,IOUT2 /0,0,0,0,0,0,0/
      DATA BESTF,BESTR,DUMMY_F,DUMMY_VIO /0.0D0,0.0D0,0.0D0,0.0D0/
      DATA TNOW,TSTART,TMAX,ACC /0.0D0,0.0D0,0.0D0,0.0D0/
      IF(C.EQ.2)THEN
          CALL GET_TIME(TNOW)
          TNOW = TNOW - TSTART
          EVAL = EVAL + P
          IF(IFLAG.GE.10)THEN
              CALL WARNINGS_AND_ERRORS( IFLAG, 0 )
              IF(SAVE2FILE.GE.1)THEN
                  CALL WARNINGS_AND_ERRORS( IFLAG, IOUT1 )
                  CALL WARNINGS_AND_ERRORS( IFLAG, IOUT2 )
              ENDIF
              IF(IFLAG.GE.100) RETURN
          ENDIF
          IF(PRINTEVAL.GE.1)THEN
              TIC = TIC + P
              IF(TIC.GE.PRINTEVAL.OR.EVAL.EQ.P.OR.IFLAG.GE.1)THEN
                  IF(EVAL.GT.P) TIC = 0
                  IF(RW(KRES).EQ.RW(WRES))THEN
                    KBEST=KF
                    WBEST=WF
                  ELSE
                    KBEST=KRES
                    WBEST=WRES
                  ENDIF
                  IF(RW(WBEST).LT.RW(KBEST).OR.
     &              (IFLAG.GE.1.OR.IFLAG.EQ.-300))THEN
                      BESTF = RW(WF)
                      BESTR = RW(WRES)
                      DO I=1,M
                        BESTG(I) = RW(WG+I)
                      ENDDO
                      DO I=1,N
                        BESTX(I) = RW(WX+I)
                      ENDDO
                  ELSE
                      BESTF = RW(KF)
                      BESTR = RW(KRES)
                      DO I=1,M
                        BESTG(I) = RW(KG+I)
                      ENDDO
                      DO I=1,N
                        BESTX(I) = RW(KX+I)
                      ENDDO
                  ENDIF
                  CALL PRINT_LINE(EVAL,TNOW,BESTF,BESTR,0)
                  IF(SAVE2FILE.GE.1)THEN
                    CALL PRINT_LINE(EVAL,TNOW,BESTF,BESTR,IOUT1)
                  ENDIF
                  IF(SAVE2FILE.GE.1)THEN
                    UPDATE = 0
                    IF( (BESTR.LT.DUMMY_VIO) .OR.
     &                  (BESTR.EQ.DUMMY_VIO.AND.BESTF.LT.DUMMY_F) )THEN
                       DUMMY_F   = BESTF
                       DUMMY_VIO = BESTR
                       UPDATE    = 1
                    ENDIF

                    IF(UPDATE.EQ.1)THEN
                     WRITE(IOUT2,31)
   31                FORMAT(/,/,'            CURRENT BEST SOLUTION')
                     CALL PRINT_SOLUTION( N, M, ME, BESTX, BESTG, BESTF,
     &                     BESTR, XL, XU, ACC, EVAL, TNOW, IFLAG, IOUT2)
                     CONTINUE
                    ENDIF
                    CALL FORCE_OUTPUT( IOUT1 )
                    CALL FORCE_OUTPUT( IOUT2 )
                  ENDIF
              ENDIF
          ENDIF
          IF(ISTOP.EQ.0)THEN
              IF(TNOW.GE.TMAX     ) IFLAG = -999
              IF(EVAL.GE.MAXEVAL-1)THEN
                IF(MAXEVAL.LE.99999999) IFLAG = -999
              ENDIF
          ENDIF
          RETURN
      ENDIF
      IF(C.EQ.1)THEN
          IFLAG = 0
          ISTOP = 0
          TMAX = DBLE(MAXTIME)
          CALL GET_TIME(TSTART)
          EVAL = 0
          if(param(1).le.0.0D0)then
              acc = 1.0D-3
          else
              acc = param(1)
          endif
          KMAX = 100
          Q    = 2*N+M+(N+5)*KMAX+8
          KX   = 1
          KF   = 2+N
          KG   = 2+N
          KRES = 2+N+1+M
          WX   = Q
          WF   = Q+1+N
          WG   = Q+1+N
          WRES = Q+1+N+1+M
          IF(SAVE2FILE.GE.1)THEN
              IOUT1 = 36
              IOUT2 = 37
              OPEN(IOUT1,FILE='MIDACO_SCREEN.TXT'  ,STATUS='UNKNOWN')
              OPEN(IOUT2,FILE='MIDACO_SOLUTION.TXT',STATUS='UNKNOWN')
          ENDIF
          BESTF = 1.0D+32
          BESTR = 1.0D+32
          DUMMY_F   = 1.0D+32
          DUMMY_VIO = 1.0D+32
          TIC = 0
          IF(PRINTEVAL.GE.1)THEN
              CALL PRINT_HEAD( N, NI, M, ME, PARAM, MAXEVAL, MAXTIME,
     &                         PRINTEVAL, SAVE2FILE, KEY, 0)
              IF(SAVE2FILE.GE.1)THEN
                CALL PRINT_HEAD( N, NI, M, ME, PARAM, MAXEVAL, MAXTIME,
     &                        PRINTEVAL, SAVE2FILE, KEY, IOUT1)
              ENDIF
          ENDIF
          IF(SAVE2FILE.GE.1) WRITE(IOUT2,3)
    3     FORMAT('MIDACO - SOLUTION',/,
     &           '-----------------',/,
     &'This file saves the current best solution X found by MIDACO.',/,
     &'This file is updated after every PRINTEVAL function evaluation,',
     &/,'if X has been improved.',/,/)
          CALL FORCE_OUTPUT( IOUT1 )
          CALL FORCE_OUTPUT( IOUT2 )
      ENDIF
      IF(C.EQ.3.AND.PRINTEVAL.GE.1)THEN
          CALL PRINT_FINAL( IFLAG,TNOW,TMAX,EVAL,MAXEVAL,
     &                      N,M,ME,X,G,F,XL,XU,RW,ACC,WRES,PARAM,
     &                      0 )
          IF(SAVE2FILE.GT.0)THEN
             CALL PRINT_FINAL( IFLAG,TNOW,TMAX,EVAL,MAXEVAL,
     &                         N,M,ME,X,G,F,XL,XU,RW,ACC,WRES,PARAM,
     &                         IOUT1 )
             CALL PRINT_FINAL( IFLAG,TNOW,TMAX,EVAL,MAXEVAL,
     &                         N,M,ME,X,G,F,XL,XU,RW,ACC,WRES,PARAM,
     &                         IOUT2 )
          ENDIF
          CALL FORCE_OUTPUT( IOUT1 )
          CALL FORCE_OUTPUT( IOUT2 )
      ENDIF
      END
      SUBROUTINE PRINT_HEAD( N, NI, M, ME, PARAM, MAXEVAL, MAXTIME,
     &                       PRINTEVAL, SAVE2FILE, KEY, IOUT)
      IMPLICIT NONE
      INTEGER N,NI,M,ME,IOUT,MAXTIME,MAXEVAL,SAVE2FILE,PRINTEVAL,I
      DOUBLE PRECISION PARAM(*)
      CHARACTER*60 KEY
      WRITE(IOUT,1) KEY,N,MAXEVAL,NI,
     &              MAXTIME,M,PRINTEVAL,ME,SAVE2FILE
      DO I = 1,9
          IF(PARAM(I).NE.0.0D0) GOTO 66
      ENDDO
      IF(PRINTEVAL.GE.1) WRITE(IOUT,101)
      GOTO 67
   66  IF(PARAM(1).NE.0.0D0.AND.DABS(PARAM(1)).LT.0.000001D0)THEN
          WRITE(IOUT,202) PARAM(1)
       ELSE
          IF(PARAM(1).EQ.0.0D0) WRITE(IOUT,203) 0.001D0
          IF(PARAM(1).NE.0.0D0) WRITE(IOUT,203) PARAM(1)
       ENDIF
       WRITE(IOUT,104) PARAM(2)
       IF(PARAM(3).NE.0.0D0)THEN
          WRITE(IOUT,205) PARAM(3)
       ELSE
          WRITE(IOUT,204) PARAM(3)
       ENDIF
       WRITE(IOUT,105) PARAM(4)
       IF(PARAM(5).NE.0.0D0)THEN
          WRITE(IOUT,207) PARAM(5)
       ELSE
          WRITE(IOUT,206) PARAM(5)
       ENDIF
       WRITE(IOUT,106) PARAM(6),PARAM(7),PARAM(8),PARAM(9)
   67 CONTINUE
       WRITE(IOUT,103)
    1 FORMAT(/,
     &       ' MIDACO 4.0    (www.midaco-solver.com)',/,
     &       ' -------------------------------------',/,/,
     &       ' LICENSE-KEY:  ',A60,/,/,
     &       ' -------------------------------------',/,
     &       ' | N', I7, '    | MAXEVAL',    I12,' |',/,
     &       ' | NI',I6, '    | MAXTIME',    I12,' |',/,
     &       ' | M', I7, '    | PRINTEVAL',  I10,' |',/,
     &       ' | ME',I6 ,'    | SAVE2FILE',  I10,' |',/,
     &       ' |-----------------------------------|')
  101 FORMAT(' | PARAMETER:   All by default (0)   |',/,
     &       ' -------------------------------------')
  202 FORMAT(' | PARAM(1)', D10.1,'  ACCURACY G(X) |')
  203 FORMAT(' | PARAM(1)', F10.6,'  ACCURACY G(X) |')
  104 FORMAT(' |-----------------------------------|',/,
     &       ' | PARAM(2)  ', F10.1,'  RANDOM-SEED |')
  204 FORMAT(' | PARAM(3)  ', F10.1,'  FSTOP       |')
  205 FORMAT(' | PARAM(3)  ', D10.3,'  FSTOP       |')
  105 FORMAT(' | PARAM(4)  ', F10.1,'  AUTOSTOP    |')
  206 FORMAT(' | PARAM(5)  ', F10.1,'  ORACLE      |')
  207 FORMAT(' | PARAM(5)  ', D10.3,'  ORACLE      |')
  106 FORMAT(' | PARAM(6)  ', F10.1,'  FOCUS       |',/,
     &       ' | PARAM(7)  ', F10.1,'  ANTS        |',/,
     &       ' | PARAM(8)  ', F10.1,'  KERNEL      |',/,
     &       ' | PARAM(9)  ', F10.1,'  CHARACTER   |',/,
     &       ' -------------------------------------')
  103   FORMAT(/,/,
     &  ' [     EVAL,    TIME]        OBJECTIVE FUNCTION VALUE    ',
     &  '     VIOLATION OF G(X)',/,
     &  ' -----------------------------------------------',
     &  '-------------------------------')
      END
      SUBROUTINE PRINT_LINE( EVAL, TNOW, F, VIO, IOUT)
      IMPLICIT NONE
      INTEGER EVAL, IOUT
      DOUBLE PRECISION TNOW, F, VIO
      IF(DABS(F).LE.1.0D+10)THEN
        IF(VIO.LE.1.0D+5)THEN
           WRITE(IOUT,1) EVAL, TNOW, F, VIO
        ELSE
           WRITE(IOUT,2) EVAL, TNOW, F, VIO
        ENDIF
      ELSE
        IF(VIO.LE.1.0D+5)THEN
           WRITE(IOUT,3) EVAL, TNOW, F, VIO
        ELSE
           WRITE(IOUT,4) EVAL, TNOW, F, VIO
        ENDIF
      ENDIF
    1 FORMAT(
     &' [',I9,',',F8.2,']        F(X):',F19.8,'         VIO:',F13.6)
    2 FORMAT(
     &' [',I9,',',F8.2,']        F(X):',F19.8,'         VIO:',D13.6)
    3 FORMAT(
     &' [',I9,',',F8.2,']        F(X):',D19.10,'         VIO:',F13.6)
    4 FORMAT(
     &' [',I9,',',F8.2,']        F(X):',D19.10,'         VIO:',D13.6)
      END
      SUBROUTINE PRINT_SOLUTION( N, M, ME, X, G, F, VIO,
     &                           XL, XU, ACC, EVAL, TNOW, IFLAG, IOUT)
      IMPLICIT NONE
      INTEGER N,M,ME,EVAL,IFLAG,IOUT,I,J
      DOUBLE PRECISION X(*),G(*),F,VIO,XL(*),XU(*),ACC,TNOW,PROFIL
      WRITE(IOUT,4)
      WRITE(IOUT,41) EVAL,TNOW,IFLAG
      IF(DABS(F).LE.1.0D+14)THEN
          WRITE(IOUT,42) F
      ELSE
          WRITE(IOUT,82) F
      ENDIF
          IF(M.GT.0) WRITE(IOUT,142)
          IF(IFLAG.LT.100)THEN
              IF(VIO.LE.1.0D+12)THEN
                WRITE(IOUT,183) VIO
              ELSE
                WRITE(IOUT,184) VIO
              ENDIF
          ENDIF
          IF(M.GT.0) WRITE(IOUT,144)
          DO I = 1,ME
              IF(DABS(G(I)).LE.ACC)THEN
                  IF(DABS(G(I)).LE.1.0D+14)THEN
                      WRITE(IOUT,43) I,G(I)
                  ELSE
                      WRITE(IOUT,83) I,G(I)
                  ENDIF
              ELSE
                  IF(DABS(G(I)).LE.1.0D+14)THEN
                      WRITE(IOUT,431) I,G(I)
                  ELSE
                      WRITE(IOUT,831) I,G(I)
                  ENDIF
              ENDIF
          ENDDO
          DO I = ME+1,M
              IF(G(I).GE.-ACC)THEN
                  IF(DABS(G(I)).LE.1.0D+14)THEN
                      WRITE(IOUT,44) I,G(I)
                  ELSE
                      WRITE(IOUT,84) I,G(I)
                  ENDIF
              ELSE
                  IF(DABS(G(I)).LE.1.0D+14)THEN
                      WRITE(IOUT,441) I,G(I)
                  ELSE
                      WRITE(IOUT,841) I,G(I)
                  ENDIF
              ENDIF
          ENDDO
          WRITE(IOUT,145)
          DO I = 1,N

          PROFIL = -1.0D0

          IF( X(I).GT.XU(I)+1.0D-6 )THEN
              PROFIL = 91.0D0
              GOTO 88
          ENDIF
          IF( X(I).LT.XL(I)-1.0D-6 )THEN
              PROFIL = 92.0D0
              GOTO 88
          ENDIF
          IF( XL(I).GT.XU(I) )THEN
              PROFIL = 93.0D0
              GOTO 88
          ENDIF
          IF( XL(I).EQ.XU(I) )THEN
              PROFIL = 90.0D0
              GOTO 88
          ENDIF
          IF( DABS(X(I)-XL(I)) .LE. (XU(I)-XL(I))/1000.0D0 )THEN
              PROFIL = 0.0D0
              GOTO 88
          ENDIF
          IF( DABS(XU(I)-X(I)) .LE. (XU(I)-XL(I))/1000.0D0 )THEN
              PROFIL = 22.0D0
              GOTO 88
          ENDIF
          DO J = 1,21
              IF( X(I) .LE. XL(I) + DBLE(J) * (XU(I)-XL(I))/21.0D0)THEN
                  PROFIL = DBLE(J)
                  GOTO 88
              ENDIF
          ENDDO
   88    CONTINUE
         IF(DABS(X(I)).LE.1.0D+14)THEN
             IF(PROFIL.EQ. 0.0D0) WRITE(IOUT,400) I,X(I)
             IF(PROFIL.EQ. 1.0D0) WRITE(IOUT,401) I,X(I)
             IF(PROFIL.EQ. 2.0D0) WRITE(IOUT,402) I,X(I)
             IF(PROFIL.EQ. 3.0D0) WRITE(IOUT,403) I,X(I)
             IF(PROFIL.EQ. 4.0D0) WRITE(IOUT,404) I,X(I)
             IF(PROFIL.EQ. 5.0D0) WRITE(IOUT,405) I,X(I)
             IF(PROFIL.EQ. 6.0D0) WRITE(IOUT,406) I,X(I)
             IF(PROFIL.EQ. 7.0D0) WRITE(IOUT,407) I,X(I)
             IF(PROFIL.EQ. 8.0D0) WRITE(IOUT,408) I,X(I)
             IF(PROFIL.EQ. 9.0D0) WRITE(IOUT,409) I,X(I)
             IF(PROFIL.EQ.10.0D0) WRITE(IOUT,410) I,X(I)
             IF(PROFIL.EQ.11.0D0) WRITE(IOUT,411) I,X(I)
             IF(PROFIL.EQ.12.0D0) WRITE(IOUT,412) I,X(I)
             IF(PROFIL.EQ.13.0D0) WRITE(IOUT,413) I,X(I)
             IF(PROFIL.EQ.14.0D0) WRITE(IOUT,414) I,X(I)
             IF(PROFIL.EQ.15.0D0) WRITE(IOUT,415) I,X(I)
             IF(PROFIL.EQ.16.0D0) WRITE(IOUT,416) I,X(I)
             IF(PROFIL.EQ.17.0D0) WRITE(IOUT,417) I,X(I)
             IF(PROFIL.EQ.18.0D0) WRITE(IOUT,418) I,X(I)
             IF(PROFIL.EQ.19.0D0) WRITE(IOUT,419) I,X(I)
             IF(PROFIL.EQ.20.0D0) WRITE(IOUT,420) I,X(I)
             IF(PROFIL.EQ.21.0D0) WRITE(IOUT,421) I,X(I)
             IF(PROFIL.EQ.22.0D0) WRITE(IOUT,422) I,X(I)
             IF(PROFIL.EQ.90.0D0) WRITE(IOUT,490) I,X(I)
             IF(PROFIL.EQ.91.0D0) WRITE(IOUT,491) I,X(I)
             IF(PROFIL.EQ.92.0D0) WRITE(IOUT,492) I,X(I)
             IF(PROFIL.EQ.93.0D0) WRITE(IOUT,493) I,X(I)
             IF(PROFIL.LT.0.0D0) WRITE(IOUT,*) 'PROFIL-ERROR'
         ELSE
             IF(PROFIL.EQ. 0.0D0) WRITE(IOUT,800) I,X(I)
             IF(PROFIL.EQ. 1.0D0) WRITE(IOUT,801) I,X(I)
             IF(PROFIL.EQ. 2.0D0) WRITE(IOUT,802) I,X(I)
             IF(PROFIL.EQ. 3.0D0) WRITE(IOUT,803) I,X(I)
             IF(PROFIL.EQ. 4.0D0) WRITE(IOUT,804) I,X(I)
             IF(PROFIL.EQ. 5.0D0) WRITE(IOUT,805) I,X(I)
             IF(PROFIL.EQ. 6.0D0) WRITE(IOUT,806) I,X(I)
             IF(PROFIL.EQ. 7.0D0) WRITE(IOUT,807) I,X(I)
             IF(PROFIL.EQ. 8.0D0) WRITE(IOUT,808) I,X(I)
             IF(PROFIL.EQ. 9.0D0) WRITE(IOUT,809) I,X(I)
             IF(PROFIL.EQ.10.0D0) WRITE(IOUT,810) I,X(I)
             IF(PROFIL.EQ.11.0D0) WRITE(IOUT,811) I,X(I)
             IF(PROFIL.EQ.12.0D0) WRITE(IOUT,812) I,X(I)
             IF(PROFIL.EQ.13.0D0) WRITE(IOUT,813) I,X(I)
             IF(PROFIL.EQ.14.0D0) WRITE(IOUT,814) I,X(I)
             IF(PROFIL.EQ.15.0D0) WRITE(IOUT,815) I,X(I)
             IF(PROFIL.EQ.16.0D0) WRITE(IOUT,816) I,X(I)
             IF(PROFIL.EQ.17.0D0) WRITE(IOUT,817) I,X(I)
             IF(PROFIL.EQ.18.0D0) WRITE(IOUT,818) I,X(I)
             IF(PROFIL.EQ.19.0D0) WRITE(IOUT,819) I,X(I)
             IF(PROFIL.EQ.20.0D0) WRITE(IOUT,820) I,X(I)
             IF(PROFIL.EQ.21.0D0) WRITE(IOUT,821) I,X(I)
             IF(PROFIL.EQ.22.0D0) WRITE(IOUT,822) I,X(I)
             IF(PROFIL.EQ.90.0D0) WRITE(IOUT,890) I,X(I)
             IF(PROFIL.EQ.91.0D0) WRITE(IOUT,891) I,X(I)
             IF(PROFIL.EQ.92.0D0) WRITE(IOUT,892) I,X(I)
             IF(PROFIL.EQ.93.0D0) WRITE(IOUT,893) I,X(I)
             IF(PROFIL.LT.0.0D0) WRITE(IOUT,*) 'PROFIL-ERROR'
          ENDIF
      ENDDO
      WRITE(IOUT,47)
    4 FORMAT(' --------------------------------------------')
   41 FORMAT(' EVAL:',I9,',  TIME:',F9.2,',  IFLAG:',I4,/,
     &       ' --------------------------------------------')
   42 FORMAT(' F(X) =',F38.15)
   82 FORMAT(' F(X) =',D38.6)
  142 FORMAT(' --------------------------------------------')
  183 FORMAT(' VIOLATION OF G(X)',F27.12)
  184 FORMAT(' VIOLATION OF G(X)',D27.6)
  144 FORMAT(' --------------------------------------------')
   43 FORMAT(' G(',I4,') = ',F15.8,'  (EQUALITY CONSTR)')
   44 FORMAT(' G(',I4,') = ',F15.8,'  (IN-EQUAL CONSTR)')
  431 FORMAT(' G(',I4,') = ',F15.8,'  (EQUALITY CONSTR)',
     &       '  <---  INFEASIBLE  ( G NOT = 0 )')
  441 FORMAT(' G(',I4,') = ',F15.8,'  (IN-EQUAL CONSTR)',
     &       '  <---  INFEASIBLE  ( G < 0 )')
   83 FORMAT(' G(',I4,') = ',D15.2,'  (EQUALITY CONSTR)')
   84 FORMAT(' G(',I4,') = ',D15.2,'  (IN-EQUAL CONSTR)')
  831 FORMAT(' G(',I4,') = ',D15.2,'  (EQUALITY CONSTR)',
     &       '  <---  INFEASIBLE  ( G NOT = 0 )')
  841 FORMAT(' G(',I4,') = ',D15.2,'  (IN-EQUAL CONSTR)',
     &       '  <---  INFEASIBLE  ( G < 0 )')
  145 FORMAT(' --------------------------------------------',
     &                              '            BOUNDS-PROFIL    ')
  400 FORMAT(' X(',I4,') = ',F34.15,'    !   XL___________________')
  401 FORMAT(' X(',I4,') = ',F34.15,'    !   x____________________')
  402 FORMAT(' X(',I4,') = ',F34.15,'    !   _x___________________')
  403 FORMAT(' X(',I4,') = ',F34.15,'    !   __x__________________')
  404 FORMAT(' X(',I4,') = ',F34.15,'    !   ___x_________________')
  405 FORMAT(' X(',I4,') = ',F34.15,'    !   ____x________________')
  406 FORMAT(' X(',I4,') = ',F34.15,'    !   _____x_______________')
  407 FORMAT(' X(',I4,') = ',F34.15,'    !   ______x______________')
  408 FORMAT(' X(',I4,') = ',F34.15,'    !   _______x_____________')
  409 FORMAT(' X(',I4,') = ',F34.15,'    !   ________x____________')
  410 FORMAT(' X(',I4,') = ',F34.15,'    !   _________x___________')
  411 FORMAT(' X(',I4,') = ',F34.15,'    !   __________x__________')
  412 FORMAT(' X(',I4,') = ',F34.15,'    !   ___________x_________')
  413 FORMAT(' X(',I4,') = ',F34.15,'    !   ____________x________')
  414 FORMAT(' X(',I4,') = ',F34.15,'    !   _____________x_______')
  415 FORMAT(' X(',I4,') = ',F34.15,'    !   ______________x______')
  416 FORMAT(' X(',I4,') = ',F34.15,'    !   _______________x_____')
  417 FORMAT(' X(',I4,') = ',F34.15,'    !   ________________x____')
  418 FORMAT(' X(',I4,') = ',F34.15,'    !   _________________x___')
  419 FORMAT(' X(',I4,') = ',F34.15,'    !   __________________x__')
  420 FORMAT(' X(',I4,') = ',F34.15,'    !   ___________________x_')
  421 FORMAT(' X(',I4,') = ',F34.15,'    !   ____________________x')
  422 FORMAT(' X(',I4,') = ',F34.15,'    !   ___________________XU')
  490 FORMAT(' X(',I4,') = ',F34.15,'    !   WARNING: XL = XU     ')
  491 FORMAT(' X(',I4,') = ',F34.15,'  <---  *** ERROR *** (X > XU) ')
  492 FORMAT(' X(',I4,') = ',F34.15,'  <---  *** ERROR *** (X < XL) ')
  493 FORMAT(' X(',I4,') = ',F34.15,'  <---  *** ERROR *** (XL > XU)')
  800 FORMAT(' X(',I4,') = ',D34.1,'    !   XL___________________')
  801 FORMAT(' X(',I4,') = ',D34.1,'    !   x____________________')
  802 FORMAT(' X(',I4,') = ',D34.1,'    !   _x___________________')
  803 FORMAT(' X(',I4,') = ',D34.1,'    !   __x__________________')
  804 FORMAT(' X(',I4,') = ',D34.1,'    !   ___x_________________')
  805 FORMAT(' X(',I4,') = ',D34.1,'    !   ____x________________')
  806 FORMAT(' X(',I4,') = ',D34.1,'    !   _____x_______________')
  807 FORMAT(' X(',I4,') = ',D34.1,'    !   ______x______________')
  808 FORMAT(' X(',I4,') = ',D34.1,'    !   _______x_____________')
  809 FORMAT(' X(',I4,') = ',D34.1,'    !   ________x____________')
  810 FORMAT(' X(',I4,') = ',D34.1,'    !   _________x___________')
  811 FORMAT(' X(',I4,') = ',D34.1,'    !   __________x__________')
  812 FORMAT(' X(',I4,') = ',D34.1,'    !   ___________x_________')
  813 FORMAT(' X(',I4,') = ',D34.1,'    !   ____________x________')
  814 FORMAT(' X(',I4,') = ',D34.1,'    !   _____________x_______')
  815 FORMAT(' X(',I4,') = ',D34.1,'    !   ______________x______')
  816 FORMAT(' X(',I4,') = ',D34.1,'    !   _______________x_____')
  817 FORMAT(' X(',I4,') = ',D34.1,'    !   ________________x____')
  818 FORMAT(' X(',I4,') = ',D34.1,'    !   _________________x___')
  819 FORMAT(' X(',I4,') = ',D34.1,'    !   __________________x__')
  820 FORMAT(' X(',I4,') = ',D34.1,'    !   ___________________x_')
  821 FORMAT(' X(',I4,') = ',D34.1,'    !   ____________________x')
  822 FORMAT(' X(',I4,') = ',D34.1,'    !   ___________________XU')
  890 FORMAT(' X(',I4,') = ',D34.1,'    !   WARNING: XL = XU     ')
  891 FORMAT(' X(',I4,') = ',D34.1,'  <---  *** ERROR *** (X > XU) ')
  892 FORMAT(' X(',I4,') = ',D34.1,'  <---  *** ERROR *** (X < XL) ')
  893 FORMAT(' X(',I4,') = ',D34.1,'  <---  *** ERROR *** (XL > XU)')
   47 FORMAT(/,' ')
      END
      SUBROUTINE PRINT_FINAL(IFLAG,TNOW,TMAX,EVAL,MAXEVAL,
     &                       N,M,ME,X,G,F,XL,XU,RW,ACC,WRES,PARAM,IOUT)
      IMPLICIT NONE
      INTEGER IFLAG,EVAL,MAXEVAL,N,M,ME,WRES,IOUT
      DOUBLE PRECISION TNOW,TMAX,X(*),G(*),F,XL(*),XU(*),RW(*),ACC
      DOUBLE PRECISION PARAM(*)
      IF(IFLAG.EQ.1.OR.IFLAG.EQ.2)THEN
        IF(TNOW.GE.TMAX)    WRITE(IOUT,411)
        IF(EVAL.GE.MAXEVAL) WRITE(IOUT,412)
      ENDIF
      IF(IFLAG.EQ.3.OR.IFLAG.EQ.4) WRITE(IOUT,413) NINT(PARAM(4))
      IF(IFLAG.EQ.5.OR.IFLAG.EQ.6) WRITE(IOUT,414)
 411  FORMAT(/,' OPTIMIZATION FINISHED  --->  MAXTIME REACHED')
 412  FORMAT(/,' OPTIMIZATION FINISHED  --->  MAXEVAL REACHED')
 413  FORMAT(/,' OPTIMIZATION FINISHED  --->  AUTOSTOP (=',I3,')')
 414  FORMAT(/,' OPTIMIZATION FINISHED  --->  FSTOP REACHED')
      WRITE(IOUT,42)
   42 FORMAT(/,/,'         BEST SOLUTION FOUND BY MIDACO       ')
      CALL PRINT_SOLUTION( N, M, ME, X, G, F, RW(WRES),
     &                     XL, XU, ACC, EVAL, TNOW, IFLAG, IOUT)
      END
      SUBROUTINE WARNINGS_AND_ERRORS( IFLAG , IOUT )
      IMPLICIT NONE
      INTEGER IFLAG,IOUT
      IF(IFLAG.LT.100)THEN
      WRITE(IOUT,1) IFLAG
    1 FORMAT(/,' *** WARNING ***   ( IFLAG =',I6,' )',/)
      ELSE
      WRITE(IOUT,2) IFLAG
    2 FORMAT(/,/,/,' *** MIDACO INPUT ERROR ***   ( IFLAG =',I6,' )')
      ENDIF
      END
      SUBROUTINE GET_TIME( SECONDS )
      DOUBLE PRECISION SECONDS
      CALL CPU_TIME( SECONDS )
      END

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C     This subroutine forces to flush the MIDACO output to text files.
c     Some Fortran compiler (e.g. g77 and f77) do not accept the 'flush'
c     command. In case you have any problems with 'flush', you can
c     savely remove or comment the below 'flush' command.
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      SUBROUTINE FORCE_OUTPUT( IOUT )
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      INTEGER IOUT
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      FLUSH( IOUT )
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
      END
c     END OF FILE

max_stars_repo_name: liujiamingustc/phd

max_issues_repo_name: liujiamingustc/phd

max_forks_repo_name: liujiamingustc/phd