!************************************************************************
!*
!*
!*     * Test problems for large nonsmooth unconstrained minimization 
!*       including also partially separable problems
!*
!*     by Napsu Karmitsa (2014, last modified 2015)
!*
!*
!*     * Problems *
!*      
!*     1.  Generalization of MAXQ (convex).
!*     2.  Generalization of MXHILB (convex).
!*     3.  Chained LQ (convex), also a partially separable version is available.
!*     4.  Chained CB3 I (convex), also a partially separable version is available.
!*     5.  Chained CB3 II (convex).
!*     6.  Number of active faces (nonconvex).
!*     7.  Nonsmooth generalization of Brown function 2 (nonconvex), also a partially 
!*         separable version is available.
!*     8.  Chained Mifflin 2 (nonconvex), also a partially separable version is available.
!*     9.  Chained crescent I (nonconvex). 
!*     10. Chained crescent II (nonconvex), also a partially separable version is available.

!*
!*     * References*
!*
!*     Haarala M., Miettinen K. and Mäkelä M.M.: New Limited Memory
!*     Bundle Method for Large-Scale Nonsmooth Optimization, Optimization
!*     Methods and Software, Vol. 19, No. 6, 2004, 673-692.
!*
!*     Bagirov A., Karmitsa N. and Mäkelä M.M.: Introduction to nonsmooth
!*     optimization: theory, practice and software, Springer, 2014.
!*
!************************************************************************
!*
!*
!*     File testps.f95 includes three modules: MODULE problem_data, 
!*     MODULE large_problems and MODULE ps_problems. 
!*
!*     MODULE problem_data is just for storing some global data needed in
!*     calculations (e.g. the number of the problem). 
!*
!*     MODULE large_problems includes the general large-scale nonsmooth 
!*     problems. It consists of the following subroutines:
!*
!*     S   startx          Initiation of variables.
!*     S   func            Computation of the value of the objective 
!*                           function.
!*     S   subgra          Computation of the subgradient of the 
!*                           objective function.
!*     S   fungra          Computation of the value and the subgradient 
!*                           of the objective function.
!*
!*     MODULE ps_problems includes the partially separable problems. It 
!*     consists of the following subroutines:
!*
!*     S   startx_ps       Initiation of variables.
!*     S   startpartial    Initiation of partially separable structure.
!*     S   func_ps         Computation of values of the objective 
!*                           function and partial functions.
!*     S   subgra_ps       Computation of subgradients of the 
!*                           objective function and partial functions.
!*
!*     For more details of these subroutines see below.
!*
!*     In addition to these modules, you need a MODULE r_precision to run 
!*     these problems (you may uncomment MODULE r_precision from below if 
!*     you do not have it elsewhere).
!* 

!MODULE r_precision ! 
!  IMPLICIT NONE
!  INTEGER, PARAMETER, PUBLIC :: prec = SELECTED_REAL_KIND(2*PRECISION(1.0))  ! old DOUBLE PRECISION
!!  INTEGER, PARAMETER, PUBLIC :: prec = SELECTED_REAL_KIND(12) 
!!  INTEGER, PARAMETER, PUBLIC :: prec = SELECTED_REAL_KIND(6) 
!END MODULE r_precision

MODULE problem_data
  USE r_precision, ONLY : prec 
  IMPLICIT NONE
  INTEGER, SAVE :: next, &                   ! Number of the problem
                   hit                       ! Number of the active partial objective
  REAL(KIND=prec), SAVE :: f_sign            ! sign of the active partial objective
  REAL(KIND=prec), DIMENSION(3), SAVE :: f_k ! partial values of the objective
END MODULE problem_data

MODULE large_problems 
  USE r_precision, ONLY : prec 
  USE problem_data, ONLY : next

  IMPLICIT NONE
  PRIVATE
  PUBLIC :: startx,func,subgra,fungra
  
CONTAINS

!************************************************************************
!*      
!*     * SUBROUTINE startx *
!*
!*      
!*     * Purpose *
!*
!*     Initiation of x.
!*
!*      
!*     * Calling sequence *
!*
!*     CALL startx(n,x)
!*
!*
!*     * Parameters *
!*      
!*     II  n          Number of variables.
!*     RO  x(n)       Vector of variables.
!*
!*
!*     * Modules *
!*     
!*     problem_data, ONLY : next
!*     r_precision, ONLY : prec 
!*
      
  SUBROUTINE startx(n,x)

!    USE problem_data, ONLY : next ! no need, if inside MODULE large_problems
!    USE r_precision, ONLY : prec ! no need, if inside MODULE large_problems
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: n                         ! Number of variables

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(OUT) :: x  ! Vector of variables

!   Local Arguments
    INTEGER :: i

!   Intrinsic Functions
    INTRINSIC MOD,REAL


    SELECT CASE (next)  ! Number of the problem

    CASE(1)  ! Generalization of MAXQ (convex)
       DO i=1,n/2
          x(i) = REAL(i, KIND=prec)
       END DO
       DO i=n/2+1,n
          x(i) = -REAL(i, KIND=prec)
       END DO

    CASE(2)  ! Generalization of MXHILB (convex)
       DO i=1,n
          x(i) = 1.0_prec
       END DO

    CASE(3)  ! Chained LQ (convex)
       DO i=1,n
          x(i) = -0.50_prec
       END DO

    CASE(4,5)  ! Chained CB3 I and II (convex)
       DO i=1,n
          x(i) = 2.0_prec
       END DO

    CASE(6)  ! Number of active faces (nonconvex)
       DO i=1,n
          x(i) = 1.0_prec
       END DO
      
    CASE(7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       DO i=1,n
          IF (MOD(i,2) == 1) THEN
             x(i) = -1.0_prec
          ELSE
             x(i) = 1.0_prec
          END IF
       END DO
      
    CASE(8)  ! Chained Mifflin 2 (nonconvex)
       DO i=1,n
          x(i) = -1.0_prec
       END DO
      
    CASE(9,10)  ! Chained crescent I and II (nonconvex)
       DO i=1,n
          IF (MOD(i,2) == 1) THEN
             x(i) = -1.50_prec
          ELSE
             x(i) = 2.0_prec
          END IF
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE startx
      
!************************************************************************
!*      
!*     * SUBROUTINE func *
!*
!*      
!*     * Purpose *
!*
!*     Computation of the value of the objective function.
!*
!*      
!*     * Calling sequence *
!*
!*     CALL func(n,x,f)
!*
!*
!*     * Parameters *
!*      
!*     II  n          Number of variables.
!*     RI  x(n)       Vector of variables.
!*     RO  f          Value of the objective function.
!*
!*
!*     * Modules *
!*     
!*     problem_data, ONLY : next,hit,f_sign,f_k
!*     r_precision, ONLY : prec 
!*
      
  SUBROUTINE func(n,x,f)

!    USE problem_data, ONLY : next  ! no need, if inside MODULE large_problems
!    USE precision, ONLY : prec  ! no need, if inside MODULE large_problems
    USE problem_data, ONLY : hit,f_sign,f_k
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: n                        ! Number of variables
    REAL(KIND=prec), INTENT(OUT) :: f               ! Value of the objective function

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(IN) :: x  ! Vector of variables

!   Local Arguments
    INTEGER :: i,j
    REAL(KIND=prec), PARAMETER :: zero=0.0, one=1.0, two=2.0
    REAL(KIND=prec) :: y,z,v,w

!   Intrinsic Functions
    INTRINSIC ABS,MAX,SIGN,LOG,EXP,REAL


    SELECT CASE (next)  ! Number of the problem
       
    CASE (1)  ! Generalization of MAXQ (convex)
       f=x(1)*x(1)
       hit=1
       DO i=2,n
          y=x(i)*x(i)
          IF (y > f) THEN
             f=y
             hit=i
          END IF
       END DO

    CASE (2)  ! Generalization of MXHILB (convex)
       f = zero
       hit = 1
       DO j = 1,n
          f = f + x(j)/REAL(j, KIND=prec)
       END DO
       f_sign=SIGN(one,f)
       f = ABS(f)
       DO i = 2,n
          y = zero
          DO j = 1,n
             y = y + x(j)/REAL(i+j-1, KIND=prec)
          END DO
          z = ABS(y)
          IF (z > f) THEN
             f=z
             hit=i
             f_sign=SIGN(one,y)
          END IF
       END DO

    CASE (3)  ! Chained LQ (convex)
       f=zero
       DO i=1,n-1
          y = -x(i)-x(i+1)
          z = -x(i)-x(i+1)+(x(i)*x(i)+x(i+1)*x(i+1)-one)
          IF (y >= z) THEN
             f=f+y
          ELSE
             f=f+z
          END IF
       END DO

    CASE (4)  ! Chained CB3 I (convex)
       f=zero
       DO i=1,n-1
          z = x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          v = (two-x(i))*(two-x(i))+ &
               (two-x(i+1))*(two-x(i+1))
          w = two*EXP(-x(i)+x(i+1))
          y = MAX(z,v,w)
          f = f+y
       END DO
      
    CASE (5)  ! Chained CB3 II (convex)
       f=zero
       f_k=zero 

       DO i=1,n-1
          f_k(1)=f_k(1)+x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          f_k(2)=f_k(2)+(two-x(i))*(two-x(i))+ &
               (two-x(i+1))*(two-x(i+1))
          f_k(3)=f_k(3)+two*EXP(-x(i)+x(i+1))
       END DO
       f=MAX(f_k(1),f_k(2),f_k(3))

    CASE (6)  ! Number of active faces (nonconvex)
       f_sign = one
       f_k(1) = -x(1)
       f = LOG(ABS(x(1))+one)
       hit = 1
       DO i=2,n
          f_k(1) = f_k(1) - x(i)
          y = LOG(ABS(x(i))+one) 
          IF(f < y) THEN
             hit = i
             f = y
          END IF
       END DO
       y = LOG(ABS(f_k(1))+one) 
       IF(f < y) THEN
          hit=n+1
          f = y
          IF (f_k(1) >= zero) f_sign = -one
       ELSE
          IF (x(hit) < zero) f_sign = -one
       END IF
      
    CASE (7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       f=zero
       DO i=1,n-1
          y = ABS(x(i))
          z = ABS(x(i+1))
          v = x(i)*x(i)+one
          w = x(i+1)*x(i+1)+one
          f = f+z**v+y**w
       END DO

    CASE (8) ! Chained mifflin 2 (nonconvex)
       f=zero
       DO i=1,n-1
          y = x(i)*x(i) + x(i+1)*x(i+1) - one
          f = f -x(i) + two*y + 1.75_prec*ABS(y)
       END DO

    CASE (9)  ! Chained crescent I (nonconvex)
       f_k(1)=zero
       f_k(2)=zero
       DO i=1,n-1
          f_k(1) = f_k(1) + x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          f_k(2) = f_k(2) - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
       END DO
       f = MAX(f_k(1),f_k(2))
       
    CASE (10)  ! Chained crescent II (nonconvex)
       f=zero
       DO i=1,n-1
          y =  x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          z =  - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
          IF (y >= z) THEN
             f=f+y
          ELSE
             f=f+z
          END IF
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE func

!************************************************************************
!*      
!*     * SUBROUTINE subgra *
!*
!*      
!*     * Purpose *
!*
!*     Computation of the subgradient of the objective function.
!*
!*      
!*     * Calling sequence *
!*
!*     CALL subgra(n,x,g)
!*
!*
!*     * Parameters *
!*      
!*     II  n          Number of variables.
!*     RI  x(n)       Vector of variables.
!*     RO  g(n)       Subgradient of the objective function.
!*
!*
!*     * Modules *
!*     
!*     problem_data, ONLY : next,hit,f_sign,f_k
!*     r_precision, ONLY : prec 
!*
      
  SUBROUTINE subgra(n,x,g)

!    USE problem_data, ONLY : next  ! no need, if inside MODULE large_problems
!    USE precision, ONLY : prec  ! no need, if inside MODULE large_problems
    USE problem_data, ONLY : hit,f_sign,f_k 
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: n                          ! Number of variables

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(IN) :: x    ! Vector of variables
    REAL(KIND=prec), DIMENSION(n), INTENT(OUT)  :: g  ! Subgradient of the objective

!   Local Arguments
    INTEGER :: i,j
    REAL(KIND=prec), PARAMETER :: zero=0.0, one=1.0, two=2.0
    REAL(KIND=prec) :: y,z,v,w,u

!   Intrinsic Functions
    INTRINSIC ABS,MAX,SIGN,LOG,EXP,REAL


    SELECT CASE (next)  ! Number of the problem
       
    CASE (1)  ! Generalization of MAXQ (convex)
       g = zero
       g(hit) = two*x(hit)

    CASE (2)  ! Generalization of MXHILB (convex)
       DO i = 1,n
          g(i) = f_sign/REAL(hit+i-1, KIND=prec)
       END DO
      
    CASE (3)  ! Chained LQ (convex)
       g(1) = zero
       DO i=1,n-1
          g(i+1) = zero
          y = -x(i)-x(i+1)
          z = -x(i)-x(i+1)+(x(i)*x(i)+x(i+1)*x(i+1)-one)
          IF (y >= z) THEN
             g(i) = g(i)-one
             g(i+1) = -one
          ELSE
             g(i) = g(i)-one+two*x(i)
             g(i+1) = -one+two*x(i+1)
          ENDIF
       END DO

    CASE (4)  ! Chained CB3 I (convex)
       g(1) = zero
       DO i=1,n-1
          g(i+1) = zero
          z = x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          v = (two-x(i))*(two-x(i))+ &
               (two-x(i+1))*(two-x(i+1))
          w = two*EXP(-x(i)+x(i+1))
          y = MAX(z,v,w) 
          IF (y == z) THEN
             g(i) = g(i)+4.0_prec*x(i)*x(i)*x(i)
             g(i+1) = two*x(i+1)
          ELSE IF (y == v) THEN
             g(i) = g(i)+two*x(i)-4.0_prec
             g(i+1) = two*x(i+1)-4.0_prec
          ELSE
             g(i) = g(i) - w
             g(i+1) = w
          END IF
       END DO
      
    CASE (5)  ! Chained CB3 II (convex)
       g(1)=zero 
       y = MAX(f_k(1),f_k(2),f_k(3))
       IF (y == f_k(1)) THEN
          DO i=1,n-1
             g(i) = g(i)+4.0_prec*x(i)*x(i)*x(i)
             g(i+1) = two*x(i+1)
          END DO
       ELSE IF (y == f_k(2)) THEN
          DO i=1,n-1
             g(i) = g(i)+two*x(i)-4.0_prec
             g(i+1) = two*x(i+1)-4.0_prec
          END DO
       ELSE
          DO i=1,n-1
             g(i) = g(i) - two*EXP(-x(i)+x(i+1))
             g(i+1) = two*EXP(-x(i)+x(i+1))
          END DO
       END IF

    CASE (6)  ! Number of active faces (nonconvex)
       IF (hit == n+1) THEN
          DO i=1,n
             g(i) = f_sign/(ABS(f_k(1))+one)
          END DO
       ELSE
          g = zero
          g(hit) = f_sign/(ABS(x(hit))+one)
       END IF
      
    CASE (7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       g(1)=zero
       DO i=1,n-1
          y = ABS(x(i))
          z = ABS(x(i+1))
          v = x(i)*x(i)+one
          w = x(i+1)*x(i+1)+one
          u = zero
          IF (z > u) u = LOG(z)
          IF (x(i) < zero) THEN
             g(i) = g(i)-w*y**(w-one)+two*x(i)*u*z**v
          ELSE
             g(i) = g(i)+w*y**(w-one)+two*x(i)*u*z**v
          END IF

          u = zero
          IF (x(i+1) == zero) THEN
             g(i+1) = zero
          ELSE IF (x(i+1) < zero) THEN
             IF (y > u) u = LOG(y)
             g(i+1)=-v*z**(v-one)+two*x(i+1)*u*y**w
          ELSE
             IF (y > u) u = LOG(y)
             g(i+1)=v*z**(v-one)+two*x(i+1)*u*y**w
          END IF
       END DO
      
    CASE (8)  ! Chained mifflin 2 (nonconvex)
       g(1)=zero
       DO i=1,n-1
          y = x(i)*x(i) + x(i+1)*x(i+1) - one
          f_sign = SIGN(3.5_prec,y) + 4.0_prec
          g(i) = g(i) + f_sign*x(i) - one
          g(i+1) = f_sign*x(i+1)
       END DO

    CASE (9)  ! Chained crescent I (nonconvex)
       g(1)=zero
       IF (f_k(1) >= f_k(2)) THEN
          DO i=1,n-1
             g(i)=g(i)+two*x(i)
             g(i+1)=two*(x(i+1)-one) + one
          END DO
       ELSE
          DO i=1,n-1
             g(i)=g(i)-two*x(i)
             g(i+1)=-two*(x(i+1)-one) + one
          END DO
       END IF
         
    CASE (10)  ! Chained crescent II (nonconvex)
       g(1)=zero
       DO i=1,n-1
          y =  x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          z =  - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
          IF (y >= z) THEN
             g(i)=g(i)+two*x(i)
             g(i+1)=two*(x(i+1)-one) + one
          ELSE
             g(i)=g(i)-two*x(i)
             g(i+1)=-two*(x(i+1)-one) + one
          END IF
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE subgra

!************************************************************************
!*      
!*     * SUBROUTINE fungra *
!*
!*      
!*     * Purpose *
!*
!*     Computation of the value and the subgradient of the objective
!*     function.
!*
!*      
!*     * Calling sequence *
!*
!*     CALL fungra(n,x,f,g)
!*
!*
!*     * Parameters *
!*      
!*     II  n          Number of variables.
!*     RI  x(n)       Vector of variables.
!*     RO  f          Value of the objective function.
!*     RO  g(n)       Subgradient of the objective function.
!*
!*
!*     * Modules *
!*     
!*     problem_data, ONLY : next,hit
!*     r_precision, ONLY : prec 
!*
      
  SUBROUTINE fungra(n,x,f,g)

!    USE problem_data, ONLY : next  ! no need, if inside MODULE large_problems
!    USE precision, ONLY : prec  ! no need, if inside MODULE large_problems
!    USE problem_data, ONLY : hit ! no need, we can have local hit "pip"
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: n                          ! Number of variables
    REAL(KIND=prec), INTENT(OUT) :: f                 ! Value of the objective

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(IN) :: x    ! Vector of variables
    REAL(KIND=prec), DIMENSION(n), INTENT(OUT)  :: g  ! Subgradient

!   Local Arguments
    INTEGER :: i,j,pip
    REAL(KIND=prec), PARAMETER :: zero=0.0, one=1.0, two=2.0
    REAL(KIND=prec) :: y,z,v,w,u

!   Intrinsic Functions
    INTRINSIC ABS,MAX,SIGN,LOG,EXP,REAL

    SELECT CASE (next)
       
    CASE (1)  ! Generalization of MAXQ (convex)
       f = x(1)*x(1)
       g(1) = zero
       pip = 1
       DO i=2,n
          y = x(i)*x(i)
          IF (y > f) THEN
             f = y
             pip = i
          END IF
          g(i) = zero
       END DO
       g(pip) = two*x(pip)

    CASE (2)  ! Generalization of MXHILB (convex)
       f = zero
       pip = 1
       DO j=1,n
          f = f + x(j)/REAL(j, KIND=prec)
       END DO
       g(1) = SIGN(one,f)
       f = ABS(f)
       DO i=2,n
          y = zero
          DO j=1,n
             y = y + x(j)/REAL(i+j-1, KIND=prec)
          END DO
          g(i) = SIGN(one,y)
          y = ABS(y)
          IF (y > f) THEN
             f = y
             pip = i
          END IF
       END DO
       z = g(pip)
       DO j=1,n
          g(j) = z/REAL(pip+j-1, KIND=prec)
       END DO
      
    CASE (3)  ! Chained LQ (convex)
       f = zero
       g(1) = zero
       DO i=1,n-1
          g(i+1) = zero
          y = -x(i)-x(i+1)
          z = -x(i)-x(i+1)+(x(i)*x(i)+x(i+1)*x(i+1)-one)
          IF (y >= z) THEN
             f = f+y
             g(i) = g(i)-one
             g(i+1) = -one
          ELSE
             f = f+z
             g(i) = g(i)-one+two*x(i)
             g(i+1) = -one+two*x(i+1)
          ENDIF
       END DO

    CASE (4)  ! Chained CB3 I (convex)
       f = zero
       g(1) = zero
       DO i=1,n-1
          g(i+1) = zero
          z = x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          v = (two-x(i))*(two-x(i))+(two-x(i+1))*(two-x(i+1))
          w = two*EXP(-x(i)+x(i+1))
          y = MAX(z,v,w)
          IF (y == z) THEN
             g(i) = g(i)+4.0_prec*x(i)*x(i)*x(i)
             g(i+1) = two*x(i+1)
          ELSE IF (y == v) THEN
             g(i) = g(i)+two*x(i)-4.0_prec
             g(i+1) = two*x(i+1)-4.0_prec
          ELSE
             g(i) = g(i) - w
             g(i+1) = w
          END IF
          f = f+y
       END DO
      
    CASE (5)  ! Chained CB3 II (convex)
       f = zero
       g(1) = zero
       y = zero
       z = zero
       v = zero

       DO i=1,n-1
          g(i+1) = zero
          y = y+x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          z = z+(two-x(i))*(two-x(i))+(two-x(i+1))*(two-x(i+1))
          v = v+two*EXP(-x(i)+x(i+1))
       END DO
       f=MAX(y,z,v)
       IF (f == y) THEN
          DO i=1,n-1
             g(i) = g(i)+4.0_prec*x(i)*x(i)*x(i)
             g(i+1) = two*x(i+1)
          END DO
       ELSE IF (f == z) THEN
          DO i=1,n-1
             g(i) = g(i)+two*x(i)-4.0_prec
             g(i+1) = two*x(i+1)-4.0_prec
          END DO
       ELSE
          DO i=1,n-1
             g(i) = g(i) - two*EXP(-x(i)+x(i+1))
             g(i+1) = two*EXP(-x(i)+x(i+1))
          END DO
       END IF

    CASE (6)  ! Number of active faces (nonconvex)
       w = one
       y = -x(1)
       f = LOG(ABS(x(1))+one)
       pip = 1
       DO i=2,n
          y = y - x(i)
          z = LOG(ABS(x(i))+one) 
          IF(f < z) THEN
             pip = i
             f = z
          END IF
       END DO
       z = LOG(ABS(y)+one) 
       IF(f < z) THEN
          f = z
          IF (y >= zero) w = -one
          DO i=1,n
             g(i) = w/(ABS(y)+one)
          END DO
       ELSE
          IF (x(pip) < zero) w = -one
          g = zero
          g(pip) = w/(ABS(x(pip))+one)
       END IF

    CASE (7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       f = zero
       g(1) = zero
       DO i=1,n-1
          y =ABS(x(i))
          z =ABS(x(i+1))
          v =x(i)*x(i)+one
          w =x(i+1)*x(i+1)+one
          f = f+z**v+y**w

          u = zero
          IF (z > u) u = LOG(z)
          IF (x(i) < zero) THEN
             g(i) = g(i)-w*y**(w-one)+two*x(i)*u*z**v
          ELSE
             g(i) = g(i)+w*y**(w-one)+two*x(i)*u*z**v
          END IF

          u = zero
          IF (x(i+1) == zero) THEN
             g(i+1)=zero
          ELSE IF (x(i+1) < zero) THEN
             IF (y > u) u = LOG(y)
             g(i+1)=-v*z**(v-one)+two*x(i+1)*u*y**w
          ELSE
             IF (y > u) u = LOG(y)
             g(i+1)=v*z**(v-one)+two*x(i+1)*u*y**w
          END IF
       END DO
      
    CASE (8)  ! Chained mifflin 2 (nonconvex)
       f = zero
       g(1) = zero
       DO i=1,n-1
          y = x(i)*x(i) + x(i+1)*x(i+1) - one
          f = f -x(i) + two*y + 1.75_prec*ABS(y)
          y = SIGN(3.5_prec,y) + 4.0_prec
          g(i) = g(i) + y*x(i) - one
          g(i+1) = y*x(i+1)
       END DO

    CASE (9)  ! Chained crescent I (nonconvex)
       y = zero
       z = zero
       DO i=1,n-1
          y = y + x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          z = z - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
       END DO
       f = MAX(y,z)
       
       g(1) = zero
       IF (y >= z) THEN
          DO i=1,n-1
             g(i) = g(i)+two*x(i)
             g(i+1) = two*(x(i+1)-one) + one
          END DO
       ELSE
          DO i=1,n-1
             g(i) = g(i)-two*x(i)
             g(i+1) = -two*(x(i+1)-one) + one
          END DO
       END IF
         
    CASE (10)  ! Chained crescent II (nonconvex)
       f = zero
       g(1) = zero
      
       DO i=1,n-1
          y = x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          z = - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
          IF (y >= z) THEN
             f = f+y
             g(i) = g(i)+two*x(i)
             g(i+1) = two*(x(i+1)-one) + one
          ELSE
             f = f+z
             g(i) = g(i)-two*x(i)
             g(i+1) = -two*(x(i+1)-one) + one
          END IF
       END DO

    CASE DEFAULT
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE fungra

END MODULE large_problems

MODULE ps_problems 
  USE r_precision, ONLY : prec 
  USE problem_data, ONLY : next

  IMPLICIT NONE
  PRIVATE
  PUBLIC :: startx_ps,startpartial,func_ps,subgra_ps
  
CONTAINS

!************************************************************************
!*      
!*     * SUBROUTINE startx_ps *
!*
!*      
!*     * Purpose *
!*
!*     Initiation of x.
!*
      
  SUBROUTINE startx_ps(n,x)

!    USE problem_data, ONLY : next ! no need, if inside MODULE ps_problems
!    USE r_precision, ONLY : prec ! no need, if inside MODULE ps_problems
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: n                         ! Number of variables

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(OUT) :: x  ! Vector of variables

!   Local Arguments
    INTEGER :: i

!   Intrinsic Functions
    INTRINSIC MOD,REAL


    SELECT CASE (next)  ! Number of the problem

    CASE(1)  ! Generalization of MAXQ (convex)
       DO i=1,n/2
          x(i) = REAL(i, KIND=prec)
       END DO
       DO i=n/2+1,n
          x(i) = -REAL(i, KIND=prec)
       END DO

    CASE(2)  ! Generalization of MXHILB (convex)
       DO i=1,n
          x(i) = 1.0_prec
       END DO

    CASE(3)  ! Chained LQ (convex)
       DO i=1,n
          x(i) = -0.50_prec
       END DO

    CASE(4,5)  ! Chained CB3 I and II (convex)
       DO i=1,n
          x(i) = 2.0_prec
       END DO

    CASE(6)  ! Number of active faces (nonconvex)
       DO i=1,n
          x(i) = 1.0_prec
       END DO
      
    CASE(7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       DO i=1,n
          IF (MOD(i,2) == 1) THEN
             x(i) = -1.0_prec
          ELSE
             x(i) = 1.0_prec
          END IF
       END DO
      
    CASE(8)  ! Chained Mifflin 2 (nonconvex)
       DO i=1,n
          x(i) = -1.0_prec
       END DO
      
    CASE(9,10)  ! Chained crescent I and II (nonconvex)
       DO i=1,n
          IF (MOD(i,2) == 1) THEN
             x(i) = -1.50_prec
          ELSE
             x(i) = 2.0_prec
          END IF
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE startx_ps

!************************************************************************
!*      
!*     * SUBROUTINE startpartial *
!*
!*      
!*     * Purpose *
!*
!*     Initialization of the number of variables n_i(n_fi) in each partial function f_i and
!*     the packed variables x_i for each partial function f_i in a form of nonzero indices 
!*     given in the array i_xi(n_imax,n_fi). E.g. if the function f_2 depends on n_2=3 
!*     variables x(1), x(3) and x(4), then i_xi(1,2)=1, i_xi(2,2)=3, and i_xi(3,2)=4.
!*
      
  SUBROUTINE startpartial(n,n_fi,n_imax,n_i,i_xi)

!    USE problem_data, ONLY : next ! no need, if inside MODULE ps_problems
!    USE r_precision, ONLY : prec ! no need, if inside MODULE ps_problems
    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: &
      n, &                  ! Number of variables
      n_fi, &               ! Number of partial functions f_i.
      n_imax                ! The maximum number of variables in partial functions f_i.

!   Array Arguments
    INTEGER, DIMENSION(n_fi), INTENT(OUT) :: n_i  ! Number of variables in each partial 
                                                  ! function f_i.
    INTEGER, DIMENSION(n_imax,n_fi), INTENT(OUT) :: i_xi  ! Array of packed variable indices.

!   Local Arguments
    INTEGER :: i

!   Intrinsic Functions
    INTRINSIC MAXVAL

    SELECT CASE (next)  ! Number of the problem

    CASE(1,2,5,6,9)     ! Generalization of MAXQ (convex)
                        ! Generalization of MXHILB (convex)
                        ! Chained CB3 II (convex)
                        ! Number of active faces (nonconvex)
                        ! Chained crescent I (nonconvex)
       next = -1

    CASE(3,4,7,8,10)    ! Chained LQ (convex) 
                        ! Chained CB3 I (convex)
                        ! Nonsmooth generalization of Brown function 2 (nonconvex)
                        ! Chained Mifflin 2 (nonconvex)
                        ! Chained crescent II (nonconvex)
       n_i = 2 
       DO i=1,n_fi
          i_xi(1,i) = i
          i_xi(2,i) = i+1
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

    IF (n_imax < MAXVAL(n_i)) THEN
       next=-2
    ELSE IF (n_imax > MAXVAL(n_i)) THEN
       next=-3
    END IF

  END SUBROUTINE startpartial
      
!************************************************************************
!*      
!*     * SUBROUTINE func_ps *
!*
!*      
!*     * Purpose *
!*
!*     Computation of the value of the objective function and the values of component functions.
!*
      
  SUBROUTINE func_ps(n,n_fi,x,f,f_i)

!    USE problem_data, ONLY : next  ! no need, if inside MODULE ps_problems
!    USE precision, ONLY : prec     ! no need, if inside MODULE ps_problems
    IMPLICIT NONE

!   Scalar Arguments
    REAL(KIND=prec), INTENT(OUT) :: f  ! Value of the objective function
    INTEGER, INTENT(IN) :: &
      n, &                  ! Number of variables
      n_fi                  ! Number of partial functions f_i.
!      n_imax                ! The maximum number of variables in partial functions f_i.

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n_fi), INTENT(OUT) :: f_i  ! Values of partial functions.
    REAL(KIND=prec), DIMENSION(n), INTENT(IN) :: x        ! Vector of variables
!    INTEGER, DIMENSION(n_imax,n_fi), INTENT(IN) :: ipx    ! Array of packed variable indices.
!    INTEGER, DIMENSION(n_fi), INTENT(in) :: n_i           ! Number of variables in each partial 
!                                                          ! function f_i.


!   Local Arguments
    REAL(KIND=prec), PARAMETER :: zero=0.0, one=1.0, two=2.0
    REAL(KIND=prec) :: y,z,v,w
    INTEGER :: i

!   Intrinsic Functions
    INTRINSIC ABS,MAX,EXP


    SELECT CASE (next)  ! Number of the problem
       
    CASE (1,2,5,6,9)  ! STOP
       next=-1

    CASE (3)  ! Chained LQ (convex)
       f=zero
       DO i=1,n_fi
          f_i(i) = -x(i)-x(i+1)
          z = -x(i)-x(i+1)+(x(i)*x(i)+x(i+1)*x(i+1)-one)
          IF (z >= f_i(i)) f_i(i)=z
          f = f + f_i(i)
       END DO

    CASE (4)  ! Chained CB3 I (convex)
       f=zero
       DO i=1,n_fi
          z = x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          v = (two-x(i))*(two-x(i))+ &
               (two-x(i+1))*(two-x(i+1))
          w = two*EXP(-x(i)+x(i+1))
          f_i(i) = MAX(z,v,w)
          f = f + f_i(i)
       END DO
      
    CASE (7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       f=zero
       DO i=1,n_fi
          y = ABS(x(i))
          z = ABS(x(i+1))
          v = x(i)*x(i)+one
          w = x(i+1)*x(i+1)+one
          f_i(i) = z**v+y**w
          f = f + f_i(i) 
       END DO

    CASE (8) ! Chained mifflin 2 (nonconvex)
       f=zero
       DO i=1,n_fi
          y = x(i)*x(i) + x(i+1)*x(i+1) - one
          f_i(i) = - x(i) + two*y + 1.75_prec*ABS(y)
          f = f + f_i(i)
       END DO

    CASE (10)  ! Chained crescent II (nonconvex)
       f=zero
       DO i=1,n_fi
          f_i(i) =  x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
                    + x(i+1) - one
          z =  - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
          IF (z >= f_i(i)) f_i(i)=z
          f = f + f_i(i)
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE func_ps

!************************************************************************
!*      
!*     * SUBROUTINE subgra_ps *
!*
!*      
!*     * Purpose *
!*
!*     Computation of the subgradient of the objective function and
!*     packed subgradients of partial functions.
!*
      
  SUBROUTINE subgra_ps(n,n_fi,n_imax,x,g,g_i)

!    USE problem_data, ONLY : next  ! no need, if inside MODULE ps_problems
!    USE precision, ONLY : prec  ! no need, if inside MODULE ps_problems

    IMPLICIT NONE

!   Scalar Arguments
    INTEGER, INTENT(IN) :: &
      n, &                  ! Number of variables
      n_fi, &               ! Number of partial functions f_i.
      n_imax                ! The maximum number of variables in partial functions f_i.

!   Array Arguments
    REAL(KIND=prec), DIMENSION(n), INTENT(OUT) :: g        ! Subgradient of the objective
    REAL(KIND=prec), DIMENSION(n_imax,n_fi), INTENT(OUT) :: g_i  ! Subgradient of partial functions
    REAL(KIND=prec), DIMENSION(n), INTENT(IN) :: x         ! Vector of variables
!    REAL(KIND=prec), DIMENSION(n_fi), INTENT(IN) :: f_i   ! Values of partial functions.
!    INTEGER, DIMENSION(n_imax,n_fi), INTENT(IN) :: ipx    ! Array of packed variable indices.
!    INTEGER, DIMENSION(n_fi), INTENT(in) :: n_i           ! Number of variables in each partial 
!                                                          ! function f_i.

!   Local Arguments
    INTEGER :: i
    REAL(KIND=prec), PARAMETER :: zero=0.0, one=1.0, two=2.0
    REAL(KIND=prec) :: y,z,v,w,u,f_sign

!   Intrinsic Functions
    INTRINSIC ABS,MAX,SIGN,LOG,EXP


    SELECT CASE (next)  ! Number of the problem
       
    CASE (1,2,5,6,9)  ! STOP
       next=-1

      
    CASE (3)  ! Chained LQ (convex)
       g(1) = zero
       DO i=1,n_fi
          y = -x(i)-x(i+1)
          z = -x(i)-x(i+1)+(x(i)*x(i)+x(i+1)*x(i+1)-one)
          IF (y >= z) THEN
             g_i(1,i) = -one
             g_i(2,i) = -one
          ELSE
             g_i(1,i) = -one+two*x(i)
             g_i(2,i) = -one+two*x(i+1)
          ENDIF
          g(i)=g(i)+g_i(1,i)
          g(i+1)=g_i(2,i)
       END DO


    CASE (4)  ! Chained CB3 I (convex)
       g(1) = zero
       DO i=1,n-1
          z = x(i)*x(i)*x(i)*x(i)+x(i+1)*x(i+1)
          v = (two-x(i))*(two-x(i))+ &
               (two-x(i+1))*(two-x(i+1))
          w = two*EXP(-x(i)+x(i+1))
          y = MAX(z,v,w) 
          IF (y == z) THEN
             g_i(1,i) = 4.0_prec*x(i)*x(i)*x(i)
             g_i(2,i) = two*x(i+1)
          ELSE IF (y == v) THEN
             g_i(1,i) = two*x(i)-4.0_prec
             g_i(2,i) = two*x(i+1)-4.0_prec
          ELSE
             g_i(1,i) = - w
             g_i(2,i) = w
          END IF
          g(i)=g(i)+g_i(1,i)
          g(i+1)=g_i(2,i)
       END DO
      
      
    CASE (7)  ! Nonsmooth generalization of Brown function 2 (nonconvex)
       g(1)=zero
       DO i=1,n_fi
          y = ABS(x(i))
          z = ABS(x(i+1))
          v = x(i)*x(i)+one
          w = x(i+1)*x(i+1)+one
          u = zero
          IF (z > u) u = LOG(z)
          IF (x(i) < zero) THEN
             g_i(1,i) = -w*y**(w-one)+two*x(i)*u*z**v
          ELSE
             g_i(1,i) = +w*y**(w-one)+two*x(i)*u*z**v
          END IF

          u = zero
          IF (x(i+1) == zero) THEN
             g_i(2,i) = zero
          ELSE IF (x(i+1) < zero) THEN
             IF (y > u) u = LOG(y)
             g_i(2,i)=-v*z**(v-one)+two*x(i+1)*u*y**w
          ELSE
             IF (y > u) u = LOG(y)
             g_i(2,i)=v*z**(v-one)+two*x(i+1)*u*y**w
          END IF
          g(i)=g(i)+g_i(1,i)
          g(i+1)=g_i(2,i)
       END DO
      
    CASE (8)  ! Chained mifflin 2 (nonconvex)
       g(1)=zero
       DO i=1,n_fi
          y = x(i)*x(i) + x(i+1)*x(i+1) - one
          f_sign = SIGN(3.5_prec,y) + 4.0_prec
          g_i(1,i) = f_sign*x(i) - one
          g_i(2,i) = f_sign*x(i+1)
          g(i) = g(i) + g_i(1,i)
          g(i+1) = g_i(2,i)
       END DO


    CASE (10)  ! Chained crescent II (nonconvex)
       g(1)=zero
       DO i=1,n_fi
          y =  x(i)*x(i) + (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) - one
          z =  - x(i)*x(i) - (x(i+1)-one)*(x(i+1)-one) &
               + x(i+1) + one
          IF (y >= z) THEN
             g_i(1,i)=two*x(i)
             g_i(2,i)=two*(x(i+1)-one) + one
          ELSE
             g_i(1,i)=-two*x(i)
             g_i(2,i)=-two*(x(i+1)-one) + one
          END IF
          g(i) = g(i) + g_i(1,i)
          g(i+1) = g_i(2,i)
       END DO

    CASE DEFAULT  ! Error
       PRINT*,'Error: Not such a problem.'
       next=-1

    END SELECT

  END SUBROUTINE subgra_ps

END MODULE ps_problems