! { dg-do compile }
! PR 62142 - this used to segfault
! Original test case by Ondřej Čertík .
program test_segfault
  implicit none
  real, allocatable :: X(:)
  allocate (x(1))
  x = 1.
  X = floor(X)
end program


text:
subroutine rom_GetRefinementCriteria(a,markel)
   use typre
   use Mod_Element
   use Mod_ZZErrorEstimator
   use Mod_PodRom

   implicit none
   class(PodRomProblem)         :: a
   class(FiniteElement), pointer :: e => NULL()
   integer(ip)                   :: markel(*)
   integer(ip)                   :: ndime,nelem,idofr
   real(rp), allocatable         :: error(:),errorROM(:,:)

   call a%Problem%Timer%Refine%Tic

   call a%Mesh%GetNelem(nelem)
   call a%Memor%alloc(nelem,error,'error','rom_GetRefinementCriteria')
   call a%Memor%alloc(nelem,a%ndofr,errorROM,'errorROM','rom_GetRefinementCriteria')
   call a%Mesh%GetNdime(ndime)

   if (a%RefinerErrorEstimator == 'ZZ   ') then
      do idofr = 1,a%ndofr
         call ZZErrorEstimator(a%Mesh,a%Memor,a%ndofn,a%Basis(:,:,idofr),errorROM(:,idofr))
      end do
      error = maxval(errorROM,2)
      call ApplyErrorCriteria(a%MPIcomm,nelem,a%RefinerErrorCriteria,a%RefinerErrorLimits,error,markel)
   elseif (a%RefinerErrorEstimator == 'GRADI') then
      do idofr = 1,a%ndofr
         call GradientErrorEstimator(a%Mesh,a%Memor,a%ndofn,a%Basis(:,:,1),errorROM(:,idofr))
      end do
      error = maxval(errorROM,2)
      call ApplyErrorCriteria(a%MPIcomm,nelem,a%RefinerErrorCriteria,a%RefinerErrorLimits,error,markel)
   elseif (a%RefinerErrorEstimator == 'PHYSI') then
      call a%Problem%GetRefinementCriteria(markel)
   endif

   call a%Memor%dealloc(nelem,a%ndofr,errorROM,'errorROM','rom_GetRefinementCriteria')
   call a%Memor%dealloc(nelem,error,'error','rom_GetRefinementCriteria')

   call a%Problem%Timer%Refine%Toc

end subroutine

max_stars_repo_name: ciaid-colombia/InsFEM

max_issues_repo_name: ciaid-colombia/InsFEM

max_forks_repo_name: ciaid-colombia/InsFEM


text:
Module Gauss_quad
!
  Real, Dimension(5,5) :: XMU_Gauss_Quad, WMU_Gauss_Quad
!
  Data XMU_Gauss_Quad(1,:)/0.5,0.,0.,0.,0./
  Data WMU_Gauss_Quad(1,:)/1. ,0.,0.,0.,0./
!
  Data XMU_Gauss_Quad(2,:)/0.1127017,0.7886751,0.,0.,0./
  Data WMU_Gauss_Quad(2,:)/0.5,0.5,0.,0.,0./
!
  Data XMU_Gauss_Quad(3,:)/0.1127017,0.5000000,0.8872983,0.,0./
  Data WMU_Gauss_Quad(3,:)/0.2777778,0.4444444,0.2777778,0.,0./
!
  Data XMU_Gauss_Quad(4,:)/6.9431841E-02,0.3300095,0.6699905,0.9305682,0./
  Data WMU_Gauss_Quad(4,:)/0.1739274,0.3260726,0.3260726,0.1739274,0./
!
  Data XMU_Gauss_Quad(5,:)/4.6910077E-02,0.2307653,0.5000000,.7692347,0.9530900/
  Data WMU_Gauss_Quad(5,:)/0.1184634,0.2393143,0.2844445,0.2393143,0.1184634/
!
End Module Gauss_quad

max_stars_repo_name: aasensio/DeepLearning

max_issues_repo_name: aasensio/DeepLearning

max_forks_repo_name: aasensio/DeepLearning


!************************************************************************************************************************
! Module containing common data for all methods and databases
!************************************************************************************************************************
MODULE commonData

IMPLICIT NONE

SAVE

!=========================================LBL COMMON DATA=================================================================
INTEGER,PARAMETER,PUBLIC :: iSOOT=0,CO2=1,H2O=2,CO=3,CH4=4,C2H4=5
TYPE GasInfo
  INTEGER :: Gas
  REAL(8) :: P,T,x
END TYPE GasInfo
REAL(8),PARAMETER,PUBLIC :: pwr= 0.1d0

LOGICAL,PUBLIC            :: Pressure_based= .FALSE.            ! linear or pressure-based Planck mean absco
REAL(8),PARAMETER,PUBLIC  :: wvnm_b=200.0d0,wvnm_e=15000.0d0    ! wavenumber range for nbk-dists(must adjust
INTEGER,PARAMETER,PUBLIC  :: n_pwrk= 5000 ! number of k-g points from lbl
REAL(8),PARAMETER,PUBLIC  :: wvnmst=1.0d-2 ! wvnmstep for user specification

!=========================================================================================================================

!======================================k-distribution method identifiers==================================================
INTEGER,PARAMETER,PUBLIC :: FSK=1,NB_FSK=2,MSFSK=3,NB_MSFSK=4,MGFSK=5,NB_MGFSK=6,MSMGFSK=7,NB_MSMGFSK=8
!========================================================================================================================

!=========================================kp file data====================================================================
INTEGER,PUBLIC :: ndata=126 ! number of data in Planck mean absco database
!=========================================================================================================================

!==================================GENERAL DB PARAMETERS==================================================================
! database identifiers
INTEGER,PARAMETER,PUBLIC :: FSMGDB=1,NBSGDB=2,NBMGDB=3

TYPE GasMixInfo
  REAL(8)  :: P,T  ! Gas properties
  REAL(8)  :: xCO2,xH2O,xCO,xCH4,xC2H4  ! species mole fraction (more can be added)
  REAL(8)  :: fv  ! soot volume fraction
END TYPE GasMixInfo

!!$ database k-g distribution species id
INTEGER,PARAMETER,PUBLIC  :: kg_CO2=1,kg_H2O=2,kg_CO=3,kg_CH4=4,kg_C2H4=5
INTEGER,PARAMETER,PUBLIC  :: spDB=5   ! number of species having k-g database

! number of points in k-g databases
INTEGER,PUBLIC,PARAMETER :: Ndb= 257, NoP=128, nfxk= 64, NNU=248

! number of groups in k-g distribution files
INTEGER,PUBLIC,PARAMETER :: Ng=4

INTEGER,PARAMETER,PUBLIC  :: imLinear= 0, im3DSpline= 7, im1Dspline= 2, im2Dspline=6, im4DLinear=0
! 0--B'000', 7--B'111', 2--B'010', 6--B'0110'
!bitwise description of Interp_Method (except im2Dspline): (B3 B2 B1)  1 -- spline, 0 -- linear. B3: P, B2: T, B1: x.
!bitwise description of im2Dspline, im4DLinear: (B4 B3 B2 B1)  1 -- spline, 0 -- linear. B4: P, B3: T, B2: Tp, B1: x.

REAL(8),PUBLIC :: kg_nb_Low(NNU),kg_nb_High(NNU)

!=========================================================================================================================

!================================NB SG DATA RETREIVAL PARAMETERS==========================================================
! data structure of index file record
! ****(MAY NEED MODIFICATION)###########
TYPE NBSGIndexRec
  INTEGER(KIND=4):: kgRecNum    ! k-g data record in which a certain k-g distribution is contained
  INTEGER(KIND=2):: offset      ! offset inside the k-g data record
  INTEGER(KIND=2):: NoP         ! the number of points included in a certain k-g distribution
  REAL(KIND=4)   :: kmax        ! the maximum value of the k-g distribution
END TYPE NBSGIndexRec

! data structure of memory index record
TYPE NBSGIndexMemRec
  INTEGER(KIND=4):: offset
  INTEGER(KIND=2):: NoP
  REAL(KIND=4)   :: kmax
END TYPE NBSGIndexMemRec

! data structure of a narrow band information
TYPE NBSGkgInfo
  INTEGER(KIND=2):: NB          ! narrow band number. 1 <= NB <= 248
  INTEGER(KIND=2):: Gas         ! species
  REAL(KIND=8)   :: P           ! total pressure
  REAL(KIND=8)   :: T           ! temperature
  REAL(KIND=8)   :: x           ! mole fraction
END TYPE NBSGkgInfo

!=========================================================================================================================

!============================NB MG DATA RETREIVAL PARAMETERS==============================================================
! data structure of index file record
TYPE NBMGIndexRec
  INTEGER(KIND=4):: kgRecNum    ! k-g data record in which a certain k-g distribution is contained
  INTEGER(KIND=2):: offset      ! offset inside the k-g data record
  INTEGER(KIND=2):: NoP         ! the number of points included in a certain k-g distribution
  REAL(KIND=4)   :: kmax        ! the maximum value of the k-g distribution
  REAL(KIND=4)   :: gmin        ! minimum g for each group within each NB
END TYPE NBMGIndexRec

! data structure of memory index record
TYPE NBMGIndexMemRec
  INTEGER(KIND=4):: offset
  INTEGER(KIND=2):: NoP
  REAL(KIND=4)   :: kmax
  REAL(KIND=4)   :: gmin
END TYPE NBMGIndexMemRec

! data structure of a narrow band information
TYPE NBMGkgInfo
  INTEGER(KIND=2):: NB          ! narrow band number. 1 <= NB <= 248
  INTEGER(KIND=2):: GP          ! groups (=4) for each NB
  INTEGER(KIND=2):: Gas         ! species
  REAL(KIND=8)   :: P           ! total pressure
  REAL(KIND=8)   :: T           ! temperature
  REAL(KIND=8)   :: x           ! mole fraction
END TYPE NBMGkgInfo
!===========================================================================================================================

!============================FS MG DATA RETREIVAL PARAMETERS==============================================================
! data structure of index file record
TYPE FSMGIndexRec
  INTEGER(KIND=4):: kgRecNum    ! k-g data record in which a certain k-g distribution is contained
  INTEGER(KIND=2):: offset      ! offset inside the k-g data record
  INTEGER(KIND=2):: NoP         ! the number of points included in a certain k-g distribution
  REAL(KIND=4)   :: kmax        ! the maximum value of the k-g distribution
  REAL(KIND=4)   :: gmin        ! minimum g for each group
END TYPE FSMGIndexRec

! data structure of memory index record
TYPE FSMGIndexMemRec
  INTEGER(KIND=4):: offset
  INTEGER(KIND=2):: NoP
  REAL(KIND=4)   :: kmax
  REAL(KIND=4)   :: gmin
END TYPE FSMGIndexMemRec

! data structure of a full-spectrum k-g information
TYPE FSMGkgInfo
  INTEGER(KIND=2):: GP          ! groups (=4) for each NB
  INTEGER(KIND=2):: Gas         ! species
  REAL(KIND=8)   :: P           ! total pressure
  REAL(KIND=8)   :: T           ! temperature
  REAL(KIND=8)   :: Tp          ! Planck function temperature
  REAL(KIND=8)   :: x           ! mole fraction
END TYPE FSMGkgInfo
!===========================================================================================================================

END MODULE commonData


text:
subroutine fsolve(N,M,fun,X,X0,Jac,options,info)
!
!   Copyright 2016, Viktor Reshiak, All rights reserved.
!
!   Purpose
!   =======
!   Find solution of a system of nonlinear equations.
!
!
!   Method
!   ======
!   Intel MKL Trust-Region algorithm.
!
!
!   IN
!   ==
!   1) N       - dimension of the input vector
!   2) M       - dimension of the objective function
!   3) fun     - objective function
!   5) X0      - (optional) initial guess as a separate vector
!   6) Jac     - (optional) N-by-N Jacobian of the objective function
!   7) options - (optional) struct with options of the solver
!
!
!   INOUT
!   =====
!   4) X       - initial guess / final solution
!
!
!   OUT
!   ===
!   8) info    - (optional) structure with information about solution process
!
	implicit none

    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    ! dummy arguments

    procedure(fsolve_objective_fun)                        :: fun
    procedure(fsolve_objective_Jac), optional              :: Jac

    integer(kind=4), intent(in)                            :: N, M
       real(kind=8), intent(inout), dimension(N)           :: X
       real(kind=8), intent(in),    dimension(N), optional :: X0

    type (fsolve_options), intent(in),  optional           :: options
    type (fsolve_info),    intent(out), optional           :: info

    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

	! parameters of the MKL solver

    ! dtrnlsp_init
    type(fsolve_options) :: op
    type(HANDLE_TR)      :: handle

    ! dtrnlsp_solve
       real(kind=8), dimension(N)   :: fvec
       real(kind=8), dimension(N,N) :: fjac
    integer(kind=4)                 :: RCI_Request

    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

	logical :: solved

    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    if ( present(options) ) then
        op = options
    else
        op = default_fsolve_options
    endif

    RCI_Request = 0
    solved      = .false.

    if ( present(X0) ) then
        X = X0
    endif

    call fun(M, N, X, fvec)
    if ( present(Jac) ) then
        call Jac(M, N, X, fjac)
    else
        if ( djacobi(fun,N,M,fjac,X,op%eps(1)) /= TR_SUCCESS ) then
            print*, '| error in djacobix'
            call MKL_FREE_BUFFERS
            stop 1;
        end if
    endif

    ! initialize solver
	if ( dtrnlsp_init(handle, N, M, X, op%eps, op%MaxIter, op%MaxIterTrial, op%rs) /= TR_SUCCESS ) then
        print*, '| error in dtrnlsp_init'
        call MKL_FREE_BUFFERS
        stop 1;
    endif
    !tmp = dtrnlsp_check(handle, n, m, fjac, fvec, op%eps, inf1)

    ! run solver
	do while (.not. solved)
        if (dtrnlsp_solve (handle, fvec, fjac, RCI_Request) /= TR_SUCCESS) then
            print*, '| error in dtrnlsp_solve'
            call MKL_FREE_BUFFERS
            stop 1;
        end if

        if (RCI_Request < 0) then
            solved = .true.
        end if

        ! recalculate function value
        if (RCI_Request == 1) then
            call fun(M, N, X, fvec)
        end if

        ! recalculate Jacobian
        if (RCI_Request == 2) then
            if ( present(Jac) ) then
        	   call Jac(M, N, X, fjac)
            else
                if ( djacobi(fun,N,M,fjac,X,op%eps(1)) /= TR_SUCCESS ) then
                    print*, '| error in djacobix'
                    call MKL_FREE_BUFFERS
                    stop 1;
                end if
            endif
        end if
	enddo

    ! retrieve information
    if ( present(info) ) then
        if (dtrnlsp_get(handle, info%iter, info%st_cr, info%init_residual, info%final_residual) /= TR_SUCCESS) then
            print*, '| error in dtrnlsp_get'
            call MKL_FREE_BUFFERS
            stop 1;
        end if
        !if ( info%iter > 10 ) then
        !    print*, info%st_cr, info%iter, info%init_residual, info%final_residual
        !endif
    endif

    if (dtrnlsp_delete(handle) /= TR_SUCCESS) then
        print*, '| error in dtrnlsp_delete'
        call MKL_FREE_BUFFERS
        stop 1;
    end if

end subroutine fsolve

max_stars_repo_name: revitmt/ChemicalKinetics

max_issues_repo_name: revitmt/ChemicalKinetics

max_forks_repo_name: revitmt/ChemicalKinetics


! ##################################################################################################################################
! Begin MIT license text.
! _______________________________________________________________________________________________________

! Copyright 2019 Dr William R Case, Jr (dbcase29@gmail.com)

! Permission is hereby granted, free of charge, to any person obtaining a copy of this software and
! associated documentation files (the "Software"), to deal in the Software without restriction, including
! without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
! copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to
! the following conditions:

! The above copyright notice and this permission notice shall be included in all copies or substantial
! portions of the Software and documentation.

! THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
! OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
! FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
! AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
! LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
! OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
! THE SOFTWARE.
! _______________________________________________________________________________________________________

! End MIT license text.

      SUBROUTINE WRITE_GRD_PCH_OUTPUTS ( JSUB, NUM, WHAT )

! Writes "punch" output for grid point related quantities (accels, displacements, eigenvectors, applied loads and SPC, MPC forces)

      USE PENTIUM_II_KIND, ONLY       :  BYTE, LONG, DOUBLE
      USE IOUNT1, ONLY                :  WRT_LOG, ERR, F04, F06, PCH
      USE SCONTR, ONLY                :  BLNK_SUB_NAM, FATAL_ERR, INT_SC_NUM, PCH_LINE_NUM, SOL_NAME
      USE TIMDAT, ONLY                :  TSEC
      USE SUBR_BEGEND_LEVELS, ONLY    :  WRITE_GRD_PCH_OUTPUTS_BEGEND
      USE NONLINEAR_PARAMS, ONLY      :  LOAD_ISTEP
      USE LINK9_STUFF, ONLY           :  GID_OUT_ARRAY, OGEL
      USE MODEL_STUF, ONLY            :  GRID, LABEL, SCNUM, SUBLOD, STITLE, TITLE
      USE EIGEN_MATRICES_1 , ONLY     :  EIGEN_VAL

      USE WRITE_GRD_PCH_OUTPUTS_USE_IFs

      IMPLICIT NONE
      CHARACTER(LEN=LEN(BLNK_SUB_NAM)):: SUBR_NAME = 'WRITE_GRD_PCH_OUTPUTS'
      CHARACTER(LEN=*) , INTENT(IN)   :: WHAT              ! Indicator whether to process displ or force output requests
      CHARACTER(LEN=1)                :: G_OR_S            ! 'G' if a grid point or 'S' if a scalar point
      CHARACTER(LEN=19)               :: OUTNAM            ! An output name for a header for the PCH file

      INTEGER(LONG), INTENT(IN)       :: JSUB              ! Solution vector number
      INTEGER(LONG), INTENT(IN)       :: NUM               ! The number of rows of OGEL to write out
      INTEGER(LONG)                   :: I,J               ! DO loop indices
      INTEGER(LONG), PARAMETER        :: SUBR_BEGEND = WRITE_GRD_PCH_OUTPUTS_BEGEND

! **********************************************************************************************************************************
      IF (WRT_LOG >= SUBR_BEGEND) THEN
         CALL OURTIM
         WRITE(F04,9001) SUBR_NAME,TSEC
 9001    FORMAT(1X,A,' BEGN ',F10.3)
      ENDIF

! **********************************************************************************************************************************
!  Make sure that WHAT is a valid value
      IF ((WHAT == 'ACCE') .OR. (WHAT == 'DISP') .OR. (WHAT == 'OLOAD') .OR. (WHAT == 'SPCF') .OR. (WHAT(1:4) == 'MPCF')) THEN
         CONTINUE
      ELSE
         WRITE(ERR,9100) WHAT
         WRITE(F06,9100) WHAT
         FATAL_ERR = FATAL_ERR + 1
         CALL OUTA_HERE ( 'Y' )
      ENDIF

! Write output headers.

      PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9003)  TITLE(INT_SC_NUM)(1:61), PCH_LINE_NUM
      PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9004) STITLE(INT_SC_NUM)(1:61), PCH_LINE_NUM
      PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9005)  LABEL(INT_SC_NUM)(1:61), PCH_LINE_NUM

      IF      (WHAT == 'ACCE' ) THEN
         OUTNAM = 'ACCELERATION       '
      ELSE IF (WHAT == 'DISP' ) THEN
         IF      ((SOL_NAME(1:7) == 'STATICS') .OR. (SOL_NAME(1:8) == 'NLSTATIC')) THEN
            OUTNAM = 'DISPLACEMENT OUTPUT'
         ELSE IF ((SOL_NAME(1:8) == 'BUCKLING') .AND. (LOAD_ISTEP == 1)) THEN
            OUTNAM = 'DISPLACEMENT OUTPUT'
         ELSE IF ((SOL_NAME(1:8) == 'BUCKLING') .AND. (LOAD_ISTEP == 2)) THEN
            OUTNAM = 'EIGENVECTOR OUTPUT '
         ELSE IF  (SOL_NAME(1:5) == 'MODES') THEN
            OUTNAM = 'EIGENVECTOR OUTPUT '
         ENDIF
      ELSE IF (WHAT == 'OLOAD') THEN
         OUTNAM = 'OLOAD OUTPUT       '
      ELSE IF (WHAT == 'SPCF' ) THEN
         OUTNAM = 'SPCF OUTPUT        '
      ELSE IF (WHAT == 'MPCF' ) THEN
         OUTNAM = 'MPCF OUTPUT        '
      ENDIF
      PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9006) OUTNAM, PCH_LINE_NUM

      PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9006) 'REAL OUTPUT        ', PCH_LINE_NUM

      IF    ((SOL_NAME(1:7) == 'STATICS') .OR. (SOL_NAME(1:8) == 'NLSTATIC')) THEN

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9011) SCNUM(JSUB), PCH_LINE_NUM

      ELSE IF ((SOL_NAME(1:8) == 'BUCKLING') .AND. (LOAD_ISTEP == 1)) THEN

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9011) SCNUM(LOAD_ISTEP), PCH_LINE_NUM

      ELSE IF ((SOL_NAME(1:8) == 'BUCKLING') .AND. (LOAD_ISTEP == 2)) THEN

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9011) SCNUM(LOAD_ISTEP), PCH_LINE_NUM
         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9013) EIGEN_VAL(JSUB), JSUB, PCH_LINE_NUM

      ELSE IF (SOL_NAME(1:5) == 'MODES') THEN

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9011) SCNUM(LOAD_ISTEP), PCH_LINE_NUM
         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9013) EIGEN_VAL(JSUB), JSUB, PCH_LINE_NUM

      ENDIF

! Write accels, displ's, applied forces or SPC forces (also calc TOTALS for forces if that is being output)
! TOTALS(J) is summation of G.P. values of applied forces, SPC forces, or MFC forces, for each of the J=1,6 components.

      DO I=1,NUM

         IF      (GRID(I,6) == 1) THEN
            G_OR_S = 'S'
         ELSE IF (GRID(I,6) == 6) THEN
            G_OR_S = 'G'
         ELSE
! ERROR
         ENDIF

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9902) GID_OUT_ARRAY(I,1), G_OR_S, (OGEL(I,J),J=1,3), PCH_LINE_NUM

         PCH_LINE_NUM = PCH_LINE_NUM + 1  ;  WRITE(PCH,9903) (OGEL(I,J),J=4,6), PCH_LINE_NUM

      ENDDO

! **********************************************************************************************************************************
      IF (WRT_LOG >= SUBR_BEGEND) THEN
         CALL OURTIM
         WRITE(F04,9002) SUBR_NAME,TSEC
 9002    FORMAT(1X,A,' END  ',F10.3)
      ENDIF

      RETURN

! **********************************************************************************************************************************
 9003 FORMAT('$TITLE   = ',A,I8)

 9004 FORMAT('$SUBTITLE= ',A,I8)

 9005 FORMAT('$LABEL   = ',A,I8)

 9006 FORMAT('$',A,52X,I8)

 9011 FORMAT('$SUBCASE ID = ',I11,47X,I8)

 9012 FORMAT('$EIGENVECTOR ',I8,51X,I8)

 9013 FORMAT('$EIGENVALUE = ',1ES14.7,2X,'MODE =',I6,30X,I8)

 9100 FORMAT(' *ERROR  9100: PROGRAMMING ERROR IN SUBROUTINE ',A                                                                   &
                    ,/,14X,' ILLEGAL INPUT FOR VARIABLE "WHAT" = ',A)

 9902 FORMAT(I10,7X,A1,3(1ES18.6),I8)

 9903 FORMAT('-CONT-',12X,3(1ES18.6),I8)

! **********************************************************************************************************************************

      END SUBROUTINE WRITE_GRD_PCH_OUTPUTS


text:
C     ( Last modified on 14 Jan 2001 at 19:03:33 )
C-----------------------------------------------------------------------------
C
      PROGRAM SELECT
C
C-----------------------------------------------------------------------------
C
C  The purpose of this program is to interrogate the file containing problem
C  classifications (by default, $MASTSIF/CLASSF.DB) for obtaining a list
C  of problems matching interactively defined characteristics.
C
C  The dialog with the user is on the standard input/output.
C
C  Programming: I. Bongartz, A.R. Conn and Ph. Toint for CGT Productions.
C               Salford fortran by Kristjan Jonasson
C
C--------- THE FOLLOWING SPECIFICATIONS MAY BE MODIFIED BY THE USER ----------
C
C  Standard input default definition is device 5.  Standard output default
C  definition is device 6. Change them to whatever values are appropriate
C  on your system.
C
      INTEGER      STDIN,     STDOUT
      PARAMETER  ( STDIN = 5, STDOUT = 6 )
C
C  Name of the classification database
C
      CHARACTER*32 NAME
      PARAMETER  ( NAME = 'CLASSF.DB' )
C
C  Device number for reading the classification database
C
      INTEGER      CLSDVC, FLSDVC
      PARAMETER  ( CLSDVC = 55, FLSDVC = 56 )
C
C    default directory for classification file.
C
C  Device number and file name for file containing default directory
C
      INTEGER      DATDVC
      PARAMETER  ( DATDVC = 57 )
      CHARACTER*32 DATNAM
      PARAMETER  ( DATNAM = 'SLCT.DAT' )
C
C---------------- END OF THE USER MODIFIABLE SPECIFICATION ------------------
C
C  Classification constants
C
      INTEGER      OBJ, CON, REG, DER, INTRST, INTVAR, VARN, VARM
      PARAMETER  ( OBJ    = 1, CON  = 2, REG  = 3, DER = 4, INTRST = 5,
     *             INTVAR = 6, VARN = 7, VARM = 8 )
C
C  Maximum number of simultaneous targets in search
C
      INTEGER      MAXTRG
      PARAMETER  ( MAXTRG = 7 )
C
C  Variable definitions
C
      CHARACTER*80 LINE
C    names for classification file
      CHARACTER*72 FILEN
C    default directory for classification file
      CHARACTER*256 DFTDIR
C    names for output listing file
C                  Addition by Kristjan Jonasson
      CHARACTER*72 FILES
      CHARACTER*28 PBCLS
      CHARACTER*8  TARGET( MAXTRG ), LIST(5)
      CHARACTER*1  CHOICE, CHAR, UPPER
      INTEGER      I, IM1, J, NVAR( MAXTRG ), NCON( MAXTRG ), L, K, NUM,
     *             NMATCH, CONVERT, NBT,  NBI, LN, UN, LM, SIZE, UM
      LOGICAL      REJECT, ANYFNV,  ANYFNC, MATCH
      INTRINSIC    MIN
C
C  Banner
C
      WRITE ( STDOUT, 1000 )
C
C    in order that this program can give the full path name for the
C    default classification file.
C
C  Open the file containing the default directory name.
C
      OPEN ( UNIT = DATDVC, FILE = DATNAM, STATUS = 'OLD' )
C
C  Read in the default directory
C
      READ ( DATDVC, 8000 ) DFTDIR
      CLOSE ( DATDVC )
      DO 910 I = 1, 256
         IF ( DFTDIR( I : I ) .EQ. ' ' ) THEN
            SIZE = I - 1
            GO TO 920
         END IF
  910 CONTINUE
      SIZE = 256
  920 CONTINUE
      FILEN = DFTDIR( 1 : SIZE ) // '/' // NAME
      SIZE  = LEN( FILEN )
C
C  MAIN LOOP
C
  1   CONTINUE
C
C  Target initialization
C
      DO 4 I = 1, MAXTRG
        NVAR(I)   = -2
        NCON(I)   = -2
        TARGET(I) = 'XXXXXXXX'
  4   CONTINUE
      TARGET(1) = '********'
      ANYFNV = .FALSE.
      ANYFNC = .FALSE.
C
C  Bounds on the number of variables and constraints are initialized
C  to be inactive.
C
      LN     = 0
      LM     = 0
      UN     = 99999999
      UM     = 99999999
C
C  Verify the name of the classification database file
C
 77   CONTINUE
      WRITE ( STDOUT, 4950 ) FILEN(1:SIZE)
      WRITE ( STDOUT, 4957 )
      READ  ( STDIN, '( A1 )' ) CHOICE
      IF ( UPPER( CHOICE ) .EQ. 'Y' ) THEN
        WRITE ( STDOUT, 4955 )
        READ  ( STDIN , FMT = '( A )', ERR = 78 ) FILEN
        SIZE = LEN( FILEN )
        GO TO 77
 78     WRITE ( STDOUT, 1201 )
        GO TO 77
      ENDIF
C
C  Writes the current specification
C
 5    CONTINUE
      WRITE ( STDOUT, 5002 )
      NUM = NBT ( OBJ, TARGET, MAXTRG )
      WRITE ( STDOUT, 5003 ) ( TARGET(L)(OBJ:OBJ), L = 1, NUM )
      NUM = NBT ( CON, TARGET, MAXTRG )
      WRITE ( STDOUT, 5004 ) ( TARGET(L)(CON:CON), L = 1, NUM )
      WRITE ( STDOUT, 5005 ) TARGET(1)(REG:REG)
      NUM = NBT ( DER, TARGET, MAXTRG )
      WRITE ( STDOUT, 5006 ) ( TARGET(L)(DER:DER), L = 1, NUM )
      NUM = NBT ( INTRST, TARGET, MAXTRG )
      WRITE ( STDOUT, 5007 ) ( TARGET(L)(INTRST:INTRST), L = 1, NUM )
      WRITE ( STDOUT, 5008 ) TARGET(1)(INTVAR:INTVAR)
      IF ( TARGET(1)(VARN:VARN) .EQ. 'V' ) THEN
        WRITE ( STDOUT, 5009 )
      ELSE IF ( TARGET(1)(VARN:VARN) .EQ. '*' ) THEN
        WRITE ( STDOUT, 5010 )
      ELSE IF ( TARGET(1)(VARN:VARN) .EQ. 'I' ) THEN
        WRITE ( STDOUT, 4008 ) LN, UN
      ELSE
        IF ( ANYFNV ) THEN
          WRITE ( STDOUT, 5011 )
        ELSE
          NUM = NBI ( NVAR, MAXTRG )
          WRITE ( STDOUT, 5012 ) ( NVAR(L), l = 1, NUM )
        ENDIF
      ENDIF
      IF ( TARGET(1)(VARM:VARM) .EQ. 'V' ) THEN
        WRITE ( STDOUT, 5013 )
      ELSE IF ( TARGET(1)(VARM:VARM) .EQ. '*' ) THEN
        WRITE ( STDOUT, 5014 )
      ELSE IF ( TARGET(1)(VARM:VARM) .EQ. 'I' ) THEN
        WRITE ( STDOUT, 6008 ) LM, UM
      ELSE
        IF ( ANYFNV ) THEN
          WRITE ( STDOUT, 5015 )
        ELSE
          NUM = NBI ( NCON, MAXTRG )
          WRITE ( STDOUT, 5016 ) ( NCON(L), l = 1, NUM )
        ENDIF
      ENDIF
C
C  Open the classification file.
C
      OPEN( UNIT = CLSDVC, FILE = FILEN, STATUS = 'OLD' )
C
C  Read in the problem characteristic one wishes to specify
C
      WRITE ( STDOUT, 5000 )
 6    CONTINUE
      WRITE ( STDOUT, 5001 )
      READ  ( STDIN, '( A1 )', ERR = 2 ) CHAR
      CHAR = UPPER( CHAR )
C
C  Get objective function target types
C
      IF ( CHAR .EQ. 'O' ) THEN
 41     CONTINUE
        NUM = 0
        WRITE ( STDOUT, 1002 )
        DO 30 I = 1, MIN( MAXTRG, 6 )
 10       CONTINUE
          IF ( NUM .EQ. 0) WRITE ( STDOUT, 1103 )
          IF ( NUM .GT. 0) WRITE ( STDOUT, 1003 )
          READ  ( STDIN , FMT = '( A1 )', ERR = 20 ) CHOICE
          CHOICE = UPPER( CHOICE )
          IF ( REJECT( CHOICE, OBJ, NUM ) ) GO TO 20
          IF ( I .GT. 1 ) THEN
            IM1 = I - 1
            DO 35 J = 1, IM1
              IF ( TARGET(J)(OBJ:OBJ) .EQ. CHOICE ) GO TO 25
 35         CONTINUE
          ENDIF
          TARGET(I)(OBJ:OBJ) = CHOICE
          IF ( CHOICE .NE. ' ' ) NUM = NUM + 1
          IF ( CHOICE .EQ. ' ' .OR. CHOICE .EQ. '*' ) GO TO 40
          GO TO 30
 20       CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 10
 25       CONTINUE
          WRITE ( STDOUT, 1101 )
          GO TO 10
 30     CONTINUE
C
C   Verify the type of objective function
C
 40     CONTINUE
        IF ( NUM .EQ. 0 ) THEN
          WRITE ( STDOUT, 2102 )
          GO TO 41
        ELSE
          WRITE ( STDOUT, 2002 ) ( TARGET(K)(OBJ:OBJ), K = 1, NUM )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 41
C
C  Get constraints target types
C
      ELSE IF ( CHAR .EQ. 'C' ) THEN
 111    CONTINUE
        NUM = 0
        WRITE ( STDOUT, 1004 )
        DO 130 I = 1, MIN( MAXTRG, 6 )
 110      CONTINUE
          IF ( NUM .EQ. 0) WRITE ( STDOUT, 1105 )
          IF ( NUM .GT. 0) WRITE ( STDOUT, 1005 )
          READ  ( STDIN , FMT = '( A1 )', ERR = 120 ) CHOICE
          CHOICE = UPPER( CHOICE )
          IF ( REJECT( CHOICE, CON, NUM ) ) GO TO 120
          IF ( I .GT. 1 ) THEN
            IM1 = I - 1
            DO 135 J = 1, IM1
              IF ( TARGET(J)(CON:CON) .EQ. CHOICE ) GO TO 125
 135        CONTINUE
          ENDIF
          TARGET(I)(CON:CON) = CHOICE
          IF ( CHOICE .NE. ' ' ) NUM = NUM + 1
          IF ( CHOICE .EQ. ' ' .OR. CHOICE .EQ. '*' ) GO TO 140
          GO TO 130
 120      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 110
 125      CONTINUE
          WRITE ( STDOUT, 1101 )
          GO TO 110
 130    CONTINUE
C
C   Verify the type of constraints
C
 140    CONTINUE
        IF ( NUM .EQ. 0 ) THEN
          WRITE ( STDOUT, 2104 )
          GO TO 111
        ELSE
          WRITE ( STDOUT, 2004 ) ( TARGET(K)(CON:CON), K = 1, NUM )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 111
C
C  Get regularity target types
C
      ELSE IF ( CHAR .EQ. 'R' ) THEN
 211    CONTINUE
        WRITE ( STDOUT, 1006 )
 210    CONTINUE
        WRITE ( STDOUT, 1007 )
        READ  ( STDIN , FMT = '( A1 )', ERR = 220 ) CHOICE
        CHOICE = UPPER( CHOICE )
        IF ( REJECT( CHOICE, REG, 0 ) ) GO TO 220
        TARGET(1)(REG:REG) = CHOICE
        GO TO 240
 220    CONTINUE
        WRITE ( STDOUT, 1001 )
        GO TO 210
 240    CONTINUE
C
C   Verify the problem's regularity
C
        WRITE ( STDOUT, 2005 ) TARGET(1)(REG:REG)
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 211
C
C  Get derivative degree
C
      ELSE IF ( CHAR .EQ. 'D' ) THEN
 311    CONTINUE
        NUM = 0
        WRITE ( STDOUT, 1008 )
        DO 330 I = 1, MIN( MAXTRG, 3 )
 310      CONTINUE
          IF ( NUM .EQ. 0 ) WRITE ( STDOUT, 1109 )
          IF ( NUM .GT. 0 ) WRITE ( STDOUT, 1009 )
          READ  ( STDIN , '( A1 )', ERR = 320 ) CHOICE
          CHOICE = UPPER( CHOICE )
          IF ( REJECT( CHOICE, DER, NUM ) ) GO TO 320
          IF ( I .GT. 1 ) THEN
            IM1 = I - 1
            DO 335 J = 1, IM1
              IF ( TARGET(J)(DER:DER) .EQ. CHOICE ) GO TO 325
 335         CONTINUE
          ENDIF
          TARGET(I)(DER:DER) = CHOICE
          IF ( CHOICE .NE. ' ' ) NUM = NUM + 1
          IF ( CHOICE .EQ. ' ' .OR. CHOICE .EQ. '*' ) GO TO 340
          GO TO 330
 320      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 310
 325      CONTINUE
          WRITE ( STDOUT, 1101 )
          GO TO 310
 330    CONTINUE
C
C   Verify the degree of available derivatives
C
 340    CONTINUE
        IF ( NUM .EQ. 0 ) THEN
          WRITE ( STDOUT, 2106 )
          GO TO 311
        ELSE
          WRITE ( STDOUT, 2006 ) ( TARGET(K)(DER:DER), K = 1, NUM )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 311
C
C  Get problem interest
C
      ELSE IF ( CHAR .EQ. 'I' ) THEN
 411    CONTINUE
        NUM = 0
        WRITE ( STDOUT, 1010 )
        DO 430 I = 1, MIN( MAXTRG, 3 )
 410      CONTINUE
          IF ( NUM .EQ. 0 ) WRITE ( STDOUT, 1111 )
          IF ( NUM .GT. 0 ) WRITE ( STDOUT, 1011 )
          READ  ( STDIN , FMT = '( A1 )', ERR = 420 ) CHOICE
          CHOICE = UPPER( CHOICE )
          IF ( REJECT( CHOICE, INTRST, NUM ) ) GO TO 420
          IF ( I .GT. 1 ) THEN
            IM1 = I - 1
            DO 435 J = 1, IM1
              IF ( TARGET(J)(INTRST:INTRST) .EQ. CHOICE ) GO TO 425
 435        CONTINUE
          ENDIF
          TARGET(I)(INTRST:INTRST) = CHOICE
          IF ( CHOICE .NE. ' ' ) NUM = NUM + 1
          IF ( CHOICE .EQ. ' ' .OR. CHOICE .EQ. '*' ) GO TO 440
          GO TO 430
 420      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 410
 425      CONTINUE
          WRITE ( STDOUT, 1101 )
          GO TO 410
 430    CONTINUE
C
C   Verify the problem's interest
C
 440    CONTINUE
        IF ( NUM .EQ. 0 ) THEN
          WRITE ( STDOUT, 2107 )
          GO TO 411
        ELSE
          WRITE ( STDOUT, 2007 ) ( TARGET(K)(INTRST:INTRST), K = 1, NUM)
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 411
C
C  Get the internal variables indicator
C
      ELSE IF ( CHAR .EQ. 'S' ) THEN
 511    CONTINUE
        WRITE ( STDOUT, 1012 )
 510    CONTINUE
        WRITE ( STDOUT, 1013 )
        READ  ( STDIN , '( A1 )', ERR = 520 ) CHOICE
        CHOICE = UPPER( CHOICE )
        IF ( REJECT( CHOICE, INTVAR, 0 ) ) GO TO 520
        TARGET(1)(INTVAR:INTVAR) = CHOICE
        GO TO 540
 520    CONTINUE
        WRITE ( STDOUT, 1001 )
        GO TO 510
C
C   Verify the indicator for explicit internal variables
C
 540    CONTINUE
        IF ( CHOICE .EQ. 'Y' ) THEN
          WRITE ( STDOUT, 2008 )
        ELSE IF ( CHOICE .EQ. 'N' ) THEN
          WRITE ( STDOUT, 2018 )
        ELSE
          WRITE ( STDOUT, 2017 )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 511
C
C  Get the number of variables
C
      ELSE IF ( CHAR .EQ. 'N' ) THEN
 611    CONTINUE
        NUM = 0
        WRITE ( STDOUT, 1014 )
 610    CONTINUE
        WRITE ( STDOUT, 1015 )
        READ  ( STDIN , FMT = '( A1 )', ERR = 620 ) CHOICE
        CHOICE = UPPER( CHOICE )
        IF ( REJECT( CHOICE, VARN, 0 ) ) GO TO 620
        TARGET(1)(VARN:VARN) = CHOICE
        GO TO 640
 620    CONTINUE
        WRITE ( STDOUT, 1001 )
        GO TO 610
C
C  Verify the number of variables type
C
 640    CONTINUE
        IF ( CHOICE .EQ. 'F' ) THEN
          WRITE ( STDOUT, 2009 )
        ELSE IF ( CHOICE .EQ. 'V' ) THEN
          WRITE ( STDOUT, 2019 )
        ELSE IF ( CHOICE .EQ. 'I' ) THEN
          WRITE ( STDOUT, 4009 )
        ELSE
          WRITE ( STDOUT, 2029 )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 611
C
C  Get the fixed number of variables
C
        IF ( TARGET(1)(VARN:VARN) .EQ. 'F' ) THEN
 641      CONTINUE
          NUM = 0
          DO 670 I = 1, MAXTRG
            WRITE( STDOUT, 1018 )
 650        CONTINUE
            IF ( NUM .EQ. 0) WRITE( STDOUT, 1119 )
            IF ( NUM .GT. 0) WRITE( STDOUT, 1019 )
            READ ( STDIN , '( A80 )', ERR = 660 ) LINE
            CHOICE = UPPER( LINE(1:1) )
            IF ( NUM .EQ. 0 ) ANYFNV =
     1           CHOICE .EQ. '*' .OR. CHOICE .EQ. ' '
            IF ( NUM .NE. 0 ) ANYFNV = CHOICE .EQ. '*'
            IF ( CHOICE .EQ. ' ' .OR. ANYFNV ) GO TO 695
            NVAR( I ) = CONVERT( LINE(1:5) )
            IF ( I .GT. 1 ) THEN
              IM1 = I - 1
              DO 655 J = 1, IM1
                IF ( NVAR(J) .EQ. NVAR(I) ) GO TO 665
 655          CONTINUE
            ENDIF
            IF ( NVAR(I) .LT. 0 .OR. NVAR(I) .GT. 99999999 ) GO TO 660
            NUM = NUM + 1
            GO TO 670
 660        CONTINUE
            WRITE ( STDOUT, 1001 )
            GO TO 650
 665        CONTINUE
            WRITE ( STDOUT, 1101 )
            GO TO 650
 670      CONTINUE
C
C   Verify the number of variables
C
 695      CONTINUE
          IF ( ANYFNV ) THEN
            WRITE ( STDOUT, 2020 )
          ELSE
            IF ( NUM .EQ. 0 ) THEN
              WRITE ( STDOUT, 2110 )
              GO TO 641
            ELSE
              WRITE ( STDOUT, 2010 ) ( NVAR(K), K = 1, NUM )
            ENDIF
          ENDIF
          WRITE ( STDOUT, 2003 )
          READ  ( STDIN, '( A1 )' ) CHOICE
          IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 641
C
C  Get an interval for the number of variables
C  a) lower bound
C
        ELSE IF ( TARGET(1)(VARN:VARN) .EQ. 'I' ) THEN
 800      CONTINUE
          LN = 0
          UN = 99999999
          WRITE ( STDOUT, 4003 )
 801      CONTINUE
          WRITE ( STDOUT, 4004 )
          READ ( STDIN, '( A80 ) ', ERR = 802 ) LINE
          CHOICE = UPPER( LINE(1:1) )
          IF ( CHOICE .NE. ' ' ) THEN
            LN = CONVERT( LINE(1:5) )
            IF ( LN .LT. 0 .OR. LN. GT. 99999999 ) GO TO 802
          ENDIF
          GO TO 805
 802      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 801
C
C  b) upper bound
C
 805      CONTINUE
          WRITE ( STDOUT, 4005 )
 803      CONTINUE
          WRITE ( STDOUT, 4006 )
          READ ( STDIN, '( A80 ) ', ERR = 804 ) LINE
          CHOICE = UPPER( LINE(1:1) )
          IF ( CHOICE .NE. ' ' ) THEN
            UN = CONVERT( LINE(1:5) )
            IF ( UN .LT. 0 .OR. UN. GT. 99999999 .OR. LN. GT. UN )
     1         GO TO 804
          ENDIF
          GO TO 806
 804      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 803
C
C  Verify the bounds on the number of variables
C
 806      CONTINUE
          WRITE ( STDOUT, 4007 ) LN, UN
          WRITE ( STDOUT, 2003 )
          READ  ( STDIN, '( A1 )' ) CHOICE
          IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 800
        ENDIF
C
C  Get the number of constraints
C
      ELSE IF ( CHAR .EQ. 'M' ) THEN
        NUM = 0
 711    CONTINUE
        WRITE ( STDOUT, 1016 )
 710    CONTINUE
        WRITE ( STDOUT, 1017 )
        READ  ( STDIN , FMT = '( A1 )', ERR = 720 ) CHOICE
        CHOICE = UPPER( CHOICE )
        IF ( REJECT( CHOICE, VARM, 0 ) ) GO TO 720
        TARGET(1)(VARM:VARM) = CHOICE
        GO TO 740
 720    CONTINUE
        WRITE ( STDOUT, 1001 )
        GO TO 710
C
C  Verify the constraints number type
C
 740    CONTINUE
        IF ( CHOICE .EQ. 'F' ) THEN
          WRITE ( STDOUT, 2011 )
        ELSE IF ( CHOICE .EQ. 'V' ) THEN
          WRITE ( STDOUT, 2021 )
        ELSE IF ( CHOICE .EQ. 'I' ) THEN
          WRITE ( STDOUT, 6009 )
        ELSE
          WRITE ( STDOUT, 2031 )
        ENDIF
        WRITE ( STDOUT, 2003 )
        READ  ( STDIN, '( A1 )' ) CHOICE
        IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 711
C
C  Get the fixed number of constraints
C
        IF ( TARGET(1)(VARM:VARM) .EQ. 'F' ) THEN
 741      CONTINUE
          NUM = 0
          DO 770 I = 1, MAXTRG
            WRITE( STDOUT, 1020 )
 750        CONTINUE
            IF ( NUM .EQ. 0) WRITE( STDOUT, 1121 )
            IF ( NUM .GT. 0) WRITE( STDOUT, 1021 )
            READ ( STDIN , FMT = '( A80 )', ERR = 760 ) LINE
            CHOICE = UPPER( LINE(1:1) )
            IF ( NUM .EQ. 0 ) ANYFNC =
     1           CHOICE .EQ. '*' .OR. CHOICE .EQ. ' '
            IF ( NUM .NE. 0 ) ANYFNC = CHOICE .EQ. '*'
            IF ( CHOICE .EQ. ' ' .OR. ANYFNC ) GO TO 795
            NCON(I) = CONVERT( LINE(1:5) )
            IF ( I .GT. 1 ) THEN
              IM1 = I - 1
              DO 755 J = 1, IM1
                IF ( NCON(J) .EQ. NCON(I) ) GO TO 765
 755          CONTINUE
            ENDIF
            IF ( NCON(I) .LT. 0 .OR. NCON(I) .GT. 99999999 ) GO TO 760
            NUM = NUM + 1
            GO TO 770
 760        CONTINUE
            WRITE ( STDOUT, 1001 )
            GO TO 750
 765        CONTINUE
            WRITE ( STDOUT, 1101 )
            GO TO 750
 770      CONTINUE
C
C  Verify the number of constraints
C
 795      CONTINUE
          IF ( ANYFNC ) THEN
            WRITE ( STDOUT, 2022 )
          ELSE
            IF ( NUM .EQ. 0 ) THEN
              WRITE ( STDOUT, 2120 )
              GO TO 741
            ELSE
              WRITE ( STDOUT, 2012 ) ( NCON(K), K = 1, NUM )
            ENDIF
          ENDIF
          WRITE ( STDOUT, 2003 )
          READ  ( STDIN, '( A1 )' ) CHOICE
          IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 741
C
C  Get an interval for the number of constraints
C  a) lower bound
C
        ELSE IF ( TARGET(1)(VARM:VARM) .EQ. 'I' ) THEN
 900      CONTINUE
          LM = 0
          UM = 99999999
          WRITE ( STDOUT, 6003 )
 901      CONTINUE
          WRITE ( STDOUT, 6004 )
          READ ( STDIN, '( A80 ) ', ERR = 902 ) LINE
          CHOICE = UPPER( LINE(1:1) )
          IF ( CHOICE .NE. ' ' ) THEN
            LM = CONVERT( LINE(1:5) )
            IF ( LM .LT. 0 .OR. LM. GT. 99999999 ) GO TO 902
          ENDIF
          GO TO 905
 902      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 901
C
C  b) upper bound
C
 905      CONTINUE
          WRITE ( STDOUT, 6005 )
 903      CONTINUE
          WRITE ( STDOUT, 6006 )
          READ ( STDIN, '( A80 ) ', ERR = 904 ) LINE
          CHOICE = UPPER( LINE(1:1) )
          IF ( CHOICE .NE. ' ' ) THEN
            UM = CONVERT( LINE(1:5) )
            IF ( UM .LT. 0 .OR. UM. GT. 99999999 .OR. LM. GT. UM )
     1          GO TO 904
          ENDIF
          GO TO 906
 904      CONTINUE
          WRITE ( STDOUT, 1001 )
          GO TO 903
C
C  Verify the bounds on the number of constraints
C
 906      CONTINUE
          WRITE ( STDOUT, 6007 ) LM, UM
          WRITE ( STDOUT, 2003 )
          READ  ( STDIN, '( A1 )' ) CHOICE
          IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 900
        ENDIF
C
C  All characteristics have been recorded. Quit
C
      ELSE IF ( CHAR .EQ. ' ' ) THEN
        GO TO 7
C
C  Error in the choice of characteristic
C
      ELSE
        GO TO 2
      ENDIF
C
C  Loop for another characteristic
C
      CLOSE ( CLSDVC )
      GO TO 5
C
C  Handle the error in characteristic choice
C
 2    CONTINUE
      WRITE ( STDOUT, 1001 )
      GO TO 6
C
C  The target(s) are now defined.
C
 7    CONTINUE
      WRITE (STDOUT, 4001 )
C
C  Loop on the problem classification and print names of matching ones
C  in lots of 5.  Also accumulate the number of matching problems found.
C
      NMATCH = 0
      L      = 0
 3000 CONTINUE
      READ ( CLSDVC, '( A28 )', END = 3010 ) PBCLS
      IF ( MATCH( PBCLS(10:28), TARGET, MAXTRG, ANYFNV, ANYFNC,
     1             NVAR, NCON, LN, UN, LM, UM ) ) THEN
        NMATCH = NMATCH + 1
        L = L + 1
        LIST(L) = PBCLS(1:8)
        IF ( L .EQ. 5 ) THEN
          WRITE ( STDOUT, 4002 ) ( LIST(I), I = 1, 5 )
          L = 0
        ENDIF
      ENDIF
      GO TO 3000
C
C  End of the database processing for the main loop.
C
C  Print selected problem names still
C  in waiting list for output.
C
 3010 CONTINUE
      IF ( L .GE. 1 ) WRITE ( STDOUT, 4002 ) ( LIST(I), I = 1, L )
      WRITE ( STDOUT, 4000 ) NMATCH
      CLOSE ( CLSDVC )
      WRITE ( STDOUT, 7001 )
      READ  ( STDIN, '( A1 )' ) CHOICE
      IF ( UPPER( CHOICE ) .EQ. 'Y' ) THEN
        OPEN( UNIT = CLSDVC, FILE = FILEN, STATUS = 'OLD' )
        WRITE ( STDOUT, 4955 )
        READ  ( STDIN , FMT = '( A )' ) FILES
        SIZE = LEN( FILES )
        OPEN( UNIT = FLSDVC, FILE = FILES, STATUS = 'UNKNOWN' )
 3020   CONTINUE
        READ ( CLSDVC, '( A28 )', END = 3050 ) PBCLS
        IF ( MATCH( PBCLS(10:28), TARGET, MAXTRG, ANYFNV, ANYFNC,
     1               NVAR, NCON, LN, UN, LM, UM ) ) THEN
          NBLK = 8
          DO 3030 NLBK = 8, 2, -1
            IF ( PBCLS(NBLK:NBLK) .NE. ' ' ) GO TO 3040
 3030     CONTINUE
 3040     CONTINUE
          WRITE ( FLSDVC, '(A)' ) PBCLS(1:NBLK)
        ENDIF
        GO TO 3020
 3050   CLOSE ( FLSDVC )
        CLOSE ( CLSDVC )
      ENDIF
      WRITE ( STDOUT, 7000 )
      READ  ( STDIN, '( A1 )' ) CHOICE
      IF ( UPPER( CHOICE ) .EQ. 'Y' ) GO TO 1
C
C  End of the database processing.
C
      STOP
C
C  Non excutable statements
C
 1000 FORMAT( /'      *************************************************'
     1        /'      *                                               *'
     1        /'      *         CONSTRAINED AND UNCONSTRAINED         *'
     1        /'      *              TESTING ENVIRONMENT              *'
     1        /'      *                                               *'
     1        /'      *                   ( CUTE )                    *'
     1        /'      *                                               *'
     1        /'      *         INTERACTIVE PROBLEM SELECTION         *'
     1        /'      *                                               *'
     1        /'      *                CGT PRODUCTIONS                *'
     1        /'      *                                               *'
     1        /'      *************************************************'
     1        / )
 1101 FORMAT( '   *** THIS SELECTION IS REPETITION.  Please choose ',
     1             'again.'
     1        / )
 1201 FORMAT( '   *** YOU USED MORE THAN 32 CHARACTERS.  Please choose',
     1             ' again.'
     1        / )
 1001 FORMAT( '   *** YOUR ANSWER IS NOT ALLOWED.  Please choose again.'
     1        / )
 4950 FORMAT( /'   Your current classification file is : ',
     1                A )
 4957 FORMAT( /'   Do you wish to change this [<CR> = N] ?',
     1             ' (N/Y)')
 4955 FORMAT(  '   Input the filename you want (up to 32 characters): ')
 5000 FORMAT( /'   CHOOSE A PROBLEM CHARACTERISTIC THAT YOU WANT',
     1         ' TO SPECIFY :'
     1        /'   ---------------------------------------------',
     1         '-------------' )
 5001 FORMAT( /'     O   : Objective type          C : Constraint',
     1               ' type'
     1        /'     R   : Regularity              I : Problem',
     1               ' interest'
     1        /'     N   : Number of variables     M : Number of',
     1               ' constraints'
     1        /'     D   : Degree of available analytic derivatives'
     1        /'     S   : Presence of explicit internal variables'
     1        /'    <CR> : No further characteristic, perform',
     1               ' selection'
     1        /' '
     1        /'   Your choice :' )
 5002 FORMAT( /'   Your current problem selection key is:'
     1        /'     ( * = anything goes )' )
 5003 FORMAT( /'     Objective function type         : ',
     1               10 ( A1, 1X ) )
 5004 FORMAT(  '     Constraints type                : ',
     1               10 ( A1, 1X ) )
 5005 FORMAT(  '     Regularity                      : ',
     1               10 ( A1, 1X ) )
 5006 FORMAT(  '     Degree of available derivatives : ',
     1               10 ( A1, 1X ) )
 5007 FORMAT(  '     Problem interest                : ',
     1               10 ( A1, 1X ) )
 5008 FORMAT(  '     Explicit internal variables     : ',
     1               10 ( A1, 1X ) )
 5009 FORMAT(  '     Number of variables             : v' )
 5010 FORMAT(  '     Number of variables             : *' )
 5011 FORMAT(  '     Number of variables             : any fixed number
     1               ')
 5012 FORMAT(  '     Number of variables             : ',
     1               10 ( I5, 1X ) )
 5013 FORMAT(  '     Number of constraints           : v' )
 5014 FORMAT(  '     Number of constraints           : *' )
 5015 FORMAT(  '     Number of constraints           : any fixed number
     1               ')
 5016 FORMAT(  '     Number of constraints           : ',
     1               10 ( I5, 1X ) )
 1002 FORMAT( /'   OBJECTIVE FUNCTION TYPE :'
     1        /'   -------------------------' )
 1103 FORMAT( /'     C   : Constant                L : Linear'
     1        /'     Q   : Quadratic               S : Sum of squares'
     1        /'     N   : No objective'
     1        /'     O   : Other (that is none of the above)'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1003 FORMAT( /'     C   : Constant                L : Linear'
     1        /'     Q   : Quadratic               S : Sum of squares'
     1        /'     N   : No objective'
     1        /'     O   : Other (that is none of the above)'
     1        /'    <CR> : No further type'
     1        /' '
     1        /'   Your choice :' )
 1004 FORMAT( /'   CONSTRAINTS TYPE :'
     1        /'   ------------------' )
 1105 FORMAT( /'     U   : No constraint           X : Fixed variables',
     1                                                 ' only'
     1        /'     B   : Bounds only             N : Linear network'
     1        /'     L   : Linear                  Q : Quadratic'
     1        /'     O   : Other (that is more general than any of the',
     1                     ' above alone)'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1005 FORMAT( /'     U   : No constraint           X : Fixed variables',
     1                                                 ' only'
     1        /'     B   : Bounds only             N : Linear network'
     1        /'     L   : Linear                  Q : Quadratic'
     1        /'     O   : Other (that is more general than any of the',
     1                     ' above alone)'
     1        /'    <CR> : No further type'
     1        /' '
     1        /'   Your choice :' )
 1006 FORMAT( /'   PROBLEM REGULARITY TYPE :'
     1        /'   -------------------------' )
 1007 FORMAT( /'     R   : Twice continuously differentiable'
     1        /'     I   : Other'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1008 FORMAT( /'   DEGREE OF AVAILABLE ANALYTICAL DERIVATIVES :'
     1        /'   --------------------------------------------' )
 1109 FORMAT( /'     0   : No analytical der.      1 : Analytical first'
     1        /'     2   : Analytical second'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1009 FORMAT( /'     0   : No analytical der.      1 : Analytical first'
     1        /'     2   : Analytical second'
     1        /'    <CR> : No further degree'
     1        /' '
     1        /'   Your choice :' )
 1010 FORMAT( /'   PROBLEM INTEREST TYPE :'
     1        /'   -----------------------' )
 1111 FORMAT( /'     A   : Academic                R : Real application'
     1        /'     M   : Modelling'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1011 FORMAT( /'     A   : Academic                R : Real application'
     1        /'     M   : Modelling'
     1        /'    <CR> : No further type'
     1        /' '
     1        /'   Your choice :' )
 1012 FORMAT( /'   PRESENCE OF EXPLICIT INTERNAL VARIABLES :'
     1        /'   -----------------------------------------' )
 1013 FORMAT( /'     Y   : Yes                     N : No'
     1        /'    <CR> : Any of the above (*)'
     1        /' '
     1        /'   Your choice :' )
 1014 FORMAT( /'   NUMBER OF VARIABLES :'
     1        /'   ---------------------' )
 1015 FORMAT( /'     F   : Fixed                   V : Variable'
     1        /'     I   : In an interval'
     1        /'    <CR> : Any number of variables (*)'
     1        /' '
     1        /'   Your choice :' )
 1016 FORMAT( /'   NUMBER OF CONSTRAINTS :'
     1        /'   -----------------------' )
 1017 FORMAT( /'     F   : Fixed                   V : Variable'
     1        /'     I   : In an interval'
     1        /'    <CR> : Any number of constraints (*)'
     1        /' '
     1        /'   Your choice :' )
 1018 FORMAT( /'   SELECT A NUMBER OF VARIABLES:'
     1        /'   -----------------------------' )
 1119 FORMAT( /'    (INT) : Select only problems with (INT) variables'
     1        /'            (minimum 0, maximum 99999999, multiple ',
     1               'choices are allowed)'
     1        /'    <CR>  : Any fixed number of variables (*)'
     1        /' '
     1        /'   Your choice :' )
 1019 FORMAT( /'    (INT) : Select only problems with (INT) variables'
     1        /'            (minimum 0, maximum 99999999, multiple ',
     1               'choices are allowed)'
     1        /'     *    : Any fixed number of variables'
     1        /'    <CR>  : No further selection',
     1        /' '
     1        /'   Your choice :' )
 1020 FORMAT( /'   SELECT A NUMBER OF CONSTRAINTS:'
     1        /'   -------------------------------' )
 1121 FORMAT( /'    (INT) : Select only problems with (INT) constraints'
     1        /'            (minimum 0, maximum 99999999, multiple ',
     1               'choices are allowed)'
     1        /'    <CR>  : Any fixed number of variables (*)'
     1        /' '
     1        /'   Your choice :' )
 1021 FORMAT( /'    (INT) : Select only problems with (INT) constraints'
     1        /'            (minimum 0, maximum 99999999, multiple ',
     1               'choices are allowed)'
     1        /'     *    : Any fixed number of variables '
     1        /'    <CR>  : No further selection',
     1        /' '
     1        /'   Your choice :' )
 2002 FORMAT( /'   You have specified objective of type(s): ',
     1             10( A1, 1X ) )
 2102 FORMAT( /'   *** Please choose an objective function type.' )
 2003 FORMAT(  '   Do you wish to reconsider your choice [<CR> = N] ?',
     1             ' (N/Y)')
 2004 FORMAT( /'   You have specified constraints of type(s): ',
     1             10( A1, 1X ) )
 2104 FORMAT( /'   *** Please choose a constraint type.' )
 2005 FORMAT( /'   You have specified regularity of type: ', A1 )
 2006 FORMAT( /'   You have specified derivatives of degree: ',
     1             10( A1, 1X ) )
 2106 FORMAT(/'   *** Please choose a degree of available derivatives.')
 2007 FORMAT( /'   You have specified problem interest of type: ',
     1             10( A1, 1X ) )
 2107 FORMAT( /'   *** Please choose a problem interest type.' )
 2008 FORMAT( /'   You have specified problems with explicit',
     1         ' internal variables.' )
 2018 FORMAT( /'   You have specified problems without explicit',
     1         ' internal variables.' )
 2017 FORMAT( /'   You have specified problems with or without',
     1         ' explicit internal variables.' )
 2009 FORMAT( /'   You have specified a fixed number of variables.' )
 2019 FORMAT( /'   You have specified a variable number of variables.' )
 2029 FORMAT( /'   You have specified any number of variables.' )
 2010 FORMAT( /'   You have specified a number of variables',
     1         ' in the set: ',
     1        /'       ', 10( I5, 1X ) )
 2020 FORMAT( /'   You have specified any fixed number of variables.' )
 2110 FORMAT( /'   *** Please choose a number of variables.' )
 2011 FORMAT( /'   You have specified a fixed number of constraints.' )
 2021 FORMAT(/'   You have specified a variable number of constraints.')
 2031 FORMAT( /'   You have specified any number of constraints.' )
 2012 FORMAT( /'   You have specified a number of constraints',
     1         ' in the set: ',
     1        /'       ', 10( I5, 1X ) )
 2022 FORMAT( /'   You have specified any fixed number of constraints.')
 2120 FORMAT( /'   *** Please choose a number of constraints.' )
 4000 FORMAT( /'   ', I5, ' Problem(s) match(es) the specification.' / )
 4001 FORMAT( /'   MATCHING PROBLEMS :'
     1        /'   -------------------' / )
 4002 FORMAT(  '     ', 5( A8, 3X ) )
 4003 FORMAT( /'   LOWER BOUND ON THE NUMBER OF VARIABLES :'
     1        /'   ----------------------------------------' )
 4004 FORMAT( /'    (INT) : Problems with at least (INT) variables'
     1        /'    <CR>  : No lower bound on the number of variables'
     1        /' '
     1        /'   Your choice : ' )
 4005 FORMAT( /'   UPPER BOUND ON THE NUMBER OF VARIABLES :'
     1        /'   ----------------------------------------' )
 4006 FORMAT( /'    (INT) : Problems with at most (INT) variables'
     1        /'    <CR>  : No upper bound on the number of variables'
     1        /' '
     1        /'   Your choice : ' )
 4007 FORMAT( /'   You have specified a number of variables in [ ',
     1        I5,', ',I5,' ]' )
 4008 FORMAT(  '     Number of variables             : in [ ',
     1        I5,', ',I5,' ]' )
 4009 FORMAT( /'   You have specified an interval for the number of',
     1         ' variables.' )
 6003 FORMAT( /'   LOWER BOUND ON THE NUMBER OF CONSTRAINTS :'
     1        /'   -----------------------------------------' )
 6004 FORMAT( /'    (INT) : Problems with at least (INT) constraints'
     1        /'    <CR>  : No lower bound on the number of constraints'
     1        /' '
     1        /'   Your choice : ' )
 6005 FORMAT( /'   UPPER BOUND ON THE NUMBER OF CONSTRAINTS :'
     1        /'   -----------------------------------------' )
 6006 FORMAT( /'    (INT) : Problems with at most (INT) constraints'
     1        /'    <CR>  : No upper bound on the number of constraints'
     1        /' '
     1        /'   Your choice : ' )
 6007 FORMAT( /'   You have specified a number of constraints in [ ',
     1        I5,', ',I5,' ]' )
 6008 FORMAT(  '     Number of constraints           : in [ ',
     1        I5,', ',I5,' ]' )
 6009 FORMAT( /'   You have specified an interval for the number of',
     1         ' constraints.' )
C                            Added by Kristjan Jonasson.
 7001 FORMAT( /'   Do you wish to save the problem names to a file',
     1             ' [<CR> = N] ? (N/Y)')
 7000 FORMAT( /'   Do you wish to make another selection [<CR> = N] ?',
     1             ' (N/Y)')
 8000 FORMAT( A256 )
      END
C
C-----------------------------------------------------------------------------
C
      LOGICAL FUNCTION REJECT( CHOICE, ITEM, NUM)
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to verify if the character CHOICE
C  is a valid specification for item ITEM, where ITEM is one of
C         1 = objective type,
C         2 = constraints type,
C         3 = problem regularity,
C         4 = degree of available derivatives,
C         5 = problem interest,
C         6 = presence of explicit internal variables,
C         7 = number of variables,
C         8 = number of constraints.
C
C  Programming: A. R. Conn and Ph. Toint for CGT Productions.
C
C-----------------------------------------------------------------------------
C
C
C  Classification constants
C
      INTEGER      OBJ, CON, REG, DER, INTRST, INTVAR, VARN, VARM
      PARAMETER  ( OBJ    = 1, CON  = 2, REG  = 3, DER = 4, INTRST = 5,
     1             INTVAR = 6, VARN = 7, VARM = 8 )
C
C  Arguments
C
      CHARACTER*1  CHOICE
      INTEGER      ITEM, NUM
C
C  Other variables
C
      LOGICAL      ADMIT
C
C  Control choices
C
      IF ( ( CHOICE .EQ. ' ' .OR. CHOICE .EQ. '*' )
     1                              .AND. NUM .EQ. 0 ) THEN
        CHOICE = '*'
        REJECT = .FALSE.
        RETURN
C
C  Objective function type
C
      ELSE IF ( ITEM .EQ. OBJ ) THEN
        ADMIT = CHOICE .EQ. 'N' .OR. CHOICE .EQ. 'C' .OR.
     1          CHOICE .EQ. 'L' .OR. CHOICE .EQ. 'Q' .OR.
     1          CHOICE .EQ. 'S' .OR. CHOICE .EQ. 'O' .OR.
     1          CHOICE .EQ. ' '
C
C  Constraint type
C
      ELSE IF ( ITEM .EQ. CON ) THEN
        ADMIT = CHOICE .EQ. 'U' .OR. CHOICE .EQ. 'B' .OR.
     1          CHOICE .EQ. 'N' .OR. CHOICE .EQ. 'L' .OR.
     1          CHOICE .EQ. 'Q' .OR. CHOICE .EQ. 'O' .OR.
     1          CHOICE .EQ. 'X' .OR. CHOICE .EQ. ' '
C
C  Problem regularity
C
      ELSE IF ( ITEM .EQ. REG ) THEN
        ADMIT = CHOICE .EQ. 'R' .OR. CHOICE .EQ. 'I'
C
C  Degree of analytical derivatives
C
      ELSE IF ( ITEM .EQ. DER ) THEN
        ADMIT = CHOICE .EQ. '0' .OR. CHOICE .EQ. '1' .OR.
     1          CHOICE .EQ. '2' .OR. CHOICE .EQ. ' '
C
C  Problem interest
C
      ELSE IF ( ITEM .EQ. INTRST ) THEN
        ADMIT = CHOICE .EQ. 'A' .OR. CHOICE .EQ. 'R' .OR.
     1          CHOICE .EQ. 'M' .OR. CHOICE .EQ. ' '
C
C  Presence of explicit internal variables
C
      ELSE IF ( ITEM .EQ. INTVAR ) THEN
         ADMIT = CHOICE .EQ. 'Y' .OR. CHOICE .EQ. 'N'
C
C  Number of variables and constraints
C
      ELSE IF ( ITEM .EQ. VARN .OR. ITEM .EQ. VARM ) THEN
        ADMIT = CHOICE .EQ. 'F' .OR. CHOICE .EQ. 'V' .OR.
     1          CHOICE .EQ. 'I'
      ENDIF
      REJECT = .NOT. ADMIT
      RETURN
      END
C
C-----------------------------------------------------------------------------
C
      LOGICAL FUNCTION MATCH ( C, T, MAXTRG, ANYFNV, ANYFNC, NV, NC,
     1                         LN, UN, LM, UM )
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to check if the classification C
C  matches one of the targets specified by T, ANYFNV, ANYFNC, NV, NC, LN
C  UN, LM and UM.
C
C  Programming: Ph. Toint, Dec 1991, for CGT Productions.
C
C-----------------------------------------------------------------------------
C
C  Arguments
C
      CHARACTER*8  T( MAXTRG )
      CHARACTER*19 C
      LOGICAL      ANYFNV, ANYFNC
      INTEGER      NV( MAXTRG ), NC( MAXTRG ), MAXTRG, LN, UN, LM, UM
C
C  Other variables
C
      INTEGER      I, J, CONVERT
      CHARACTER*1  CH
C
C  Match objective function type
C
      MATCH = .FALSE.
      DO 10 I = 1, MAXTRG
        CH = T(I)(1:1)
        IF ( CH .EQ. ' ' ) GO TO 20
        MATCH = MATCH .OR. CH .EQ. C(1:1) .OR. CH .EQ. '*'
 10   CONTINUE
 20   CONTINUE
      IF ( .NOT. MATCH ) RETURN
C
C  Match constraint type
C
      MATCH = .FALSE.
      DO 30 I = 1, MAXTRG
        CH = T(I)(2:2)
        IF ( CH .EQ. ' ' ) GO TO 40
        MATCH = MATCH .OR. CH .EQ. C(2:2) .OR. CH .EQ. '*'
 30   CONTINUE
 40   CONTINUE
      IF ( .NOT. MATCH ) RETURN
C
C  Match regularity type
C
      CH = T(1)(3:3)
      MATCH = CH .EQ. C(3:3) .OR. CH .EQ. '*'
      IF ( .NOT. MATCH ) RETURN
C
C  Match degree of available derivatives
C
      MATCH = .FALSE.
      DO 50 I = 1, MAXTRG
        CH = T(I)(4:4)
        IF ( CH .EQ. ' ' ) GO TO 60
        MATCH = MATCH .OR. CH .EQ. C(4:4) .OR. CH .EQ. '*'
 50   CONTINUE
 60   CONTINUE
      IF ( .NOT. MATCH ) RETURN
C
C  Match interest of the problem
C
      MATCH = .FALSE.
      DO 70 I = 1, MAXTRG
        CH = T(I)(5:5)
        IF ( CH .EQ. ' ' ) GO TO 80
        MATCH = MATCH .OR. CH .EQ. C(6:6) .OR. CH .EQ. '*'
 70   CONTINUE
 80   CONTINUE
      IF ( .NOT. MATCH ) RETURN
C
C  Match for explicit internal variables
C
      CH = T(1)(6:6)
      MATCH = CH .EQ. C(7:7) .OR. CH .EQ. '*'
      IF ( .NOT. MATCH ) RETURN
C
C  Match the number of variables
C
      MATCH = .FALSE.
      CH = T(1)(7:7)
      IF ( CH .EQ. '*' ) THEN
        MATCH = .TRUE.
      ELSE IF ( CH .EQ. 'V' ) THEN
        MATCH = C(13:13) .EQ. 'V'
      ELSE IF ( CH .EQ. 'I' ) THEN
C    interval.  A variable number of variables is no longer considered to
C    match a number in an interval.
        IF ( C(13:13) .NE. 'V' ) THEN
          I = CONVERT( C(9:13) )
          MATCH = I. GE. LN .AND. I .LE. UN
        ENDIF
      ELSE
        IF ( ANYFNV ) THEN
          MATCH = .TRUE.
        ELSE
C    A variable number of variables is no longer considered to match a fixed
C    number.
          IF ( C(13:13) .NE. 'V' ) THEN
            J = CONVERT( C(9:13) )
            DO 90 I = 1, MAXTRG
              IF ( NV(I) .LT. 0 ) GO TO 90
              MATCH = NV(I) .EQ. J
              IF ( MATCH ) GO TO 100
 90         CONTINUE
 100        CONTINUE
          ENDIF
        ENDIF
      ENDIF
      IF ( .NOT. MATCH ) RETURN
C
C  Match the number of constraints
C
      MATCH = .FALSE.
      CH = T(1)(8:8)
      IF ( CH .EQ. '*' ) THEN
        MATCH = .TRUE.
      ELSE IF ( CH .EQ. 'V' ) THEN
        MATCH = C(19:19) .EQ. 'V'
      ELSE IF ( CH .EQ. 'I' ) THEN
C    interval.  A variable number of constraints is no longer considered to
C    match a number in an interval.
        IF ( C(19:19) .NE. 'V' ) THEN
          I = CONVERT( C(15:19) )
          MATCH = I. GE. LM .AND. I .LE. UM
        ENDIF
      ELSE
        IF ( ANYFNC ) THEN
          MATCH = .TRUE.
        ELSE
C    A variable number of constraints is no longer considered to match a fixed
C    number.
          IF ( C(19:19) .NE. 'V' ) THEN
            J = CONVERT( C(15:19) )
            DO 110 I = 1, MAXTRG
              IF ( NC(I) .LT. 0 ) GO TO 110
              MATCH = NC(I) .EQ. J
              IF ( MATCH ) GO TO 120
 110        CONTINUE
 120        CONTINUE
          ENDIF
        ENDIF
      ENDIF
      RETURN
      END
C
C-----------------------------------------------------------------------------
C
      INTEGER FUNCTION CONVERT( LINE )
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to convert the nonnegative integer contained
C  in the first 5 characters of the string LINE into a proper integer.
C  -1 is returned in case of error.
C
C  Programming: Ph. Toint, Dec 1991, for CGT Productions.
C
C-----------------------------------------------------------------------------
C
C  Arguments
C
      CHARACTER*5    LINE
C
C  Other variables
C
      INTEGER L
C
C  Remove trailing blanks
C
      DO 10 L= 5, 1, -1
        IF ( LINE(L:L) .NE. ' ' ) GO TO 20
 10   CONTINUE
      GO TO 30
C
C  Read the integer
C
 20   CONTINUE
      IF ( L .EQ. 1 ) THEN
        READ ( LINE(1:1), '( I1 )', ERR = 30 ) CONVERT
      ELSE IF ( L .EQ. 2 ) THEN
        READ ( LINE(1:2), '( I2 )', ERR = 30 ) CONVERT
      ELSE IF ( L .EQ. 3 ) THEN
        READ ( LINE(1:3), '( I3 )', ERR = 30 ) CONVERT
      ELSE IF ( L .EQ. 4 ) THEN
        READ ( LINE(1:4), '( I4 )', ERR = 30 ) CONVERT
      ELSE
        READ ( LINE(1:5), '( I5 )', ERR = 30 ) CONVERT
      ENDIF
      RETURN
C
C  Error
C
 30   CONTINUE
      CONVERT = -1
      RETURN
      END
C
C-----------------------------------------------------------------------------
C
      INTEGER FUNCTION NBT ( I, T, MX )
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to return the number of non trivial
C  choices of item I specified in the target list T.
C
C  Programming: Ph. Toint, Dec 1991, for CGT Productions.
C
C-----------------------------------------------------------------------------
C
C  Arguments
C
      CHARACTER*8 T( MX )
      INTEGER     I, MX
C
C  Other variables
C
      INTEGER L
C
      NBT = 1
      DO 10 L = 2, MX
        IF ( T(L)(I:I) .EQ. 'X' ) RETURN
        NBT = NBT + 1
 10   CONTINUE
      RETURN
      END
C
C-----------------------------------------------------------------------------
C
      INTEGER FUNCTION NBI ( N, MX )
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to return the number of non trivial
C  dimensions specified in the integer dimension vector N.
C
C  Programming: Ph. Toint, Dec 1991, for CGT Productions.
C
C-----------------------------------------------------------------------------
C
C  Arguments
C
      INTEGER N( MX ), MX
C
C  Other variables
C
      INTEGER L
C
      NBI = 0
      DO 10 L = 1, MX
        IF ( N(L) .LT. 0 ) RETURN
        NBI = NBI + 1
 10   CONTINUE
      RETURN
      END
C
C-----------------------------------------------------------------------------
C
      CHARACTER*1 FUNCTION UPPER( CH )
C
C-----------------------------------------------------------------------------
C
C  The purpose of this function is to transform the character CH
C  to upper case, if it is not already.
C
C  Programming: Ph. Toint, Dec 1991, for CGT Productions.
C
C-----------------------------------------------------------------------------
C
      CHARACTER*1  CH
C
      INTRINSIC    ICHAR, CHAR
C
      INTEGER      ICH, LSTART
C
      LSTART = ICHAR( 'a' )
      ICH    = ICHAR( CH )
      IF ( ICH .GE. LSTART .AND. ICH .LE. ICHAR( 'z' ) ) THEN
        UPPER = CHAR( ICHAR( 'A' ) + ICH - LSTART )
      ELSE
        UPPER = CH
      ENDIF
      RETURN
      END

C
C-----------------------------------------------------------------------------
C

ext: f

max_forks_repo_name: dynamics-of-stellar-systems/dynamite


text:
#include "1.inc"

module same_file_mod
    use xmod_only_mod
    use other_file_mod
end module same_file_mod

program p1
    use inc_mod
    use same_file_mod
    use iso_c_binding
end program p1

max_stars_repo_name: clementval/claw-compiler

max_issues_repo_name: skosukhin/claw-compiler

max_forks_repo_name: skosukhin/claw-compiler


text:
!===============================================================================
program EC_grid_rotate_bin
!===============================================================================
! Rotate all the constants on an EC grid by some angles.

use EC_grid
use anisotropy_ajn

   implicit none

   ! Precision selector
   integer, parameter :: rs = 8

   ! Types of grid: assumed-int, explicit-real and explicit-int
   type(ECgrid) :: g
   type(ECgridr) :: gr
   type(ECgridi) :: gi
   integer :: grid_type,ix,iy,iz,ios
   real(rs) :: Cr(6,6)

   ! Angles
   real(rs) :: a,b,c

   ! In/outfiles
   character(len=250) :: infile, outfile, arg

   ! Verbosity
   logical :: quiet = .false.

   ! Check arguments
   if (command_argument_count() /= 5) then
      write(0,'(a)') 'Usage: EC_grid_rotate_bin [infile] [a] [b] [c] [outfile]',&
                     '  a, b and c are respectively rotations about the x1, x2 and x3 axes,',&
                     '  applied clockwise when looking down the axis'
      stop
   endif

   ! Get arguments
   call get_command_argument(1,infile)
   call get_command_argument(2,arg); read(arg,*,iostat=ios) a
   call bad_angle_read
   call get_command_argument(3,arg); read(arg,*,iostat=ios) b
   call bad_angle_read
   call get_command_argument(4,arg); read(arg,*,iostat=ios) c
   call bad_angle_read
   call get_command_argument(5,outfile)

   ! Determine the coordinate type of the grid
   call EC_grid_check_type_bin(infile,grid_type,quiet=quiet)

   ! Read in file, rotate constants and write out to new file
   ! Assumed-int
   if (grid_type == EC_grid_nver_assumed_i) then
      call EC_grid_load_bin(infile,g,quiet=quiet)
      do iz=1,g%nz
         do iy=1,g%ny
            do ix=1,g%nz
               call CIJ_rot3(g%ecs(ix,iy,iz,:,:),a,b,c,CR)
               g%ecs(ix,iy,iz,:,:) = CR
            enddo
         enddo
      enddo
      call EC_grid_write_bin(outfile,g)

   ! Int
   else if (grid_type == EC_grid_nver_i) then
      call EC_gridi_load_bin(infile,gi,quiet=quiet)
      do iz=1,gi%nz
         do iy=1,gi%ny
            do ix=1,gi%nz
               call CIJ_rot3(gi%ecs(ix,iy,iz,:,:),a,b,c,CR)
               gi%ecs(ix,iy,iz,:,:) = CR
            enddo
         enddo
      enddo
      call EC_gridi_write_bin(outfile,gi)

   ! Real
   else if (grid_type == EC_grid_nver_r) then
      call EC_gridr_load_bin(infile,gr,quiet=quiet)
      do iz=1,gr%nz
         do iy=1,gr%ny
            do ix=1,gr%nz
               call CIJ_rot3(gr%ecs(ix,iy,iz,:,:),a,b,c,CR)
               gr%ecs(ix,iy,iz,:,:) = CR
            enddo
         enddo
      enddo
      call EC_gridr_write_bin(outfile,gr)

   ! Bad file
   else
      write(0,'(a)') 'EC_grid_rotate_bin: input file is not a recognised EC_grid binary file.'
      stop
   endif

   contains

   !----------------------------------------------------------------------------
   subroutine bad_angle_read
   !----------------------------------------------------------------------------
   implicit none

      if (ios /= 0) then
         write(0,'(a)') 'EC_grid_rotate_bin: cannot read rotation angle from command line'
         stop
      endif
   end subroutine bad_angle_read
   !----------------------------------------------------------------------------

end program EC_grid_rotate_bin
!-------------------------------------------------------------------------------

max_stars_repo_name: anowacki/elasticity

max_issues_repo_name: anowacki/elasticity

max_forks_repo_name: anowacki/elasticity


text:
SUBROUTINE change_sec (j1,upstream)
! Code converted using TO_F90 by Alan Miller
! Date: 2005-08-23  Time: 12:57:15

use common_h
use time_h
use reservoir_h

IMPLICIT NONE

INTEGER, INTENT(IN OUT)                  :: upstream

INTEGER :: j1,g !,npt

!REAL :: max_dheight,min_dheight,dummy,dist_bank1,dist_bank2
!REAL :: max_adjust
REAL ::elev
!REAL :: predom_grain,perc_grain

!REAL :: dummy1,dummy2
real :: dummy3,dummy4 !,damelev, dummy,xmin_tw,
real :: wfactor(200),dummy5
! -----------------------------------------------------------------------

! calculation of reservoir bed elevation changes

! identification of the sediment balance process (erosion/deposition)

! 1) erosion
!dummy1=0

IF (dvol_sed(j1,res_index(upstream)) < 0.) THEN
  elev=erosion_level(j1,res_index(upstream))
!  elev=maxelev_sec(j1,upstream)
!  elev=watelev_sec(j1,upstream)
!  elev=min(minelev_sec(j1,upstream)+(depth_sec(j1,upstream)*1.5),maxelev_sec(j1,upstream))
!  elev=watelev_sec(j1,upstream)+((maxelev_sec(j1,upstream)-watelev_sec(j1,upstream))/20.)
!if (id_sec_extern(j1,upstream)==43) elev=max(elev,y_sec(11,j1,upstream))
!if (id_sec_extern(j1,upstream)==47) elev=max(elev,y_sec(18,j1,upstream))
!if (id_sec_extern(j1,upstream)==47)elev=max(elev,435.)
!  elev=min(watelev_sec(j1,upstream)+5.,maxelev_sec(j1,upstream))
!  elev=min(watelev_sec(j1,upstream)+1.5,watelev_sec(j1,upstream)+((maxelev_sec(j1,upstream)-watelev_sec(j1,upstream))/5.))
  dummy5=0.
  DO m=1,npoints(j1,res_index(upstream))
    if(y_sec(m,j1,res_index(upstream))<elev) then
	  dummy5=dummy5+(y_sec(m,j1,res_index(upstream))-y_actlay(m,j1,res_index(upstream)))
	endif
  enddo

  DO m=1,npoints(j1,res_index(upstream))
    if(y_sec(m,j1,res_index(upstream))<elev .and. dummy5/=0.) then
	  wfactor(m)=(y_sec(m,j1,res_index(upstream))-y_actlay(m,j1,res_index(upstream)))/dummy5
	else
	  wfactor(m)=0.
    endif
  enddo

  dummy3=0.
  DO m=2,npoints(j1,res_index(upstream))-1
    IF (geom(m,j1) == 1)dummy3=dummy3+(((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(1./2.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m+1)*(1./8.)))
    IF (geom(m,j1) == 2)dummy3=dummy3+((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m-1)*(1./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m+1)*(1./8.))
    IF (geom(m,j1) == 3)dummy3=dummy3+((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m-1)*(1./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(1./2.))
    IF (geom(m,j1) == 4)dummy3=dummy3+((x_sec(m+1,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(1./2.))
!if(j1<30)write(*,'(3I4,4F10.2)')j1,m,geom(m,j1),dummy3
  enddo

  dummy4=darea_sed(j1,res_index(upstream))/dummy3
  DO m=2,npoints(j1,res_index(upstream))-1
	y_laststep(m,j1,res_index(upstream))=y_sec(m,j1,res_index(upstream))
	IF (geom(m,j1) /= 0)y_sec(m,j1,res_index(upstream))=y_sec(m,j1,res_index(upstream))-dummy4*wfactor(m)
!if(j1<20)write(*,'(2I4,3I3,6F12.6)')t,d,j1,m,geom(m,j1),y_sec(m,j1,upstream),y_laststep(m,j1,upstream),dummy4*wfactor(m),dummy4,wfactor(m)
!write(*,'(2I4,3I3,6F12.6)')t,d,j1,m,geom(m,j1),y_sec(m,j1,upstream),y_laststep(m,j1,upstream),dummy4*wfactor(m),dummy4,wfactor(m)
  ENDDO

!if(d==2)stop

! 2) deposition
ELSE IF (dvol_sed(j1,res_index(upstream)) > 0.) THEN

  DO m=1,npoints(j1,res_index(upstream))
    if(y_sec(m,j1,res_index(upstream))<watelev_sec(j1,res_index(upstream))) then
!	  dummy5=(watelev_sec(j1,upstream)-y_sec(m,j1,upstream))/(watelev_sec(j1,upstream)-minelev_sec(j1,upstream))
!	  wfactor(m)=1.-((1.-(dummy5))**2.9)
	  wfactor(m)=(watelev_sec(j1,res_index(upstream))-y_sec(m,j1,res_index(upstream)))/(watelev_sec(j1,res_index(upstream))-minelev_sec(j1,res_index(upstream)))
	else
!	  dummy5=0.
	  wfactor(m)=0.
    endif
!if(j1<20)write(*,'(2I4,3I3,2F15.6)')t,d,j1,m,geom(m,j1),wfactor(m),dummy5
  enddo

  dummy3=0.
  DO m=2,npoints(j1,res_index(upstream))-1
    IF (geom(m,j1) == 1)dummy3=dummy3+(((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(1./2.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m+1)*(1./8.)))
    IF (geom(m,j1) == 2)dummy3=dummy3+((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m-1)*(1./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m+1)*(1./8.))
    IF (geom(m,j1) == 3)dummy3=dummy3+((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(3./8.))+ &
							((x_sec(m,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m-1)*(1./8.))+ &
							((x_sec(m+1,j1,res_index(upstream))-x_sec(m,j1,res_index(upstream)))*wfactor(m)*(1./2.))
    IF (geom(m,j1) == 4)dummy3=dummy3+((x_sec(m+1,j1,res_index(upstream))-x_sec(m-1,j1,res_index(upstream)))*wfactor(m)*(1./2.))
!if(j1<30)write(*,'(3I4,4F10.2)')j1,m,geom(m,j1),dummy3
  enddo

  dummy4=darea_sed(j1,res_index(upstream))/dummy3
  DO m=2,npoints(j1,res_index(upstream))-1
	y_laststep(m,j1,res_index(upstream))=y_sec(m,j1,res_index(upstream))
	IF (geom(m,j1) /= 0)y_sec(m,j1,res_index(upstream))=y_sec(m,j1,res_index(upstream))+dummy4*wfactor(m)
!if(j1<20)write(*,'(2I4,3I3,6F12.6)')t,d,j1,m,geom(m,j1),y_sec(m,j1,upstream),y_laststep(m,j1,upstream),dummy4*wfactor(m),dummy4,wfactor(m)
  ENDDO

!if(d==2)stop

END IF

!if(d==41)stop

DO g=1,n_sed_class
  if (frac_actlay(g,j1,res_index(upstream)) /= 0.) then
!	dummy1=diam(g)
  endif
enddo

!Ge User-specified positive multiplication factor for defining the thickness of active layer
!Ge constant=10/20/30
!dummy2=20.*dummy1*24.

!DO m=1,npoints(j1,upstream)
!  if (x_sec(m,j1,upstream) == x_minelev(j1,upstream)) then
!    y_actlay(m,j1,upstream)=max(y_sec(m,j1,upstream)-dummy,y_original(m,j1,upstream))
!	dummy=y_actlay(m,j1,upstream)
!  endif
!enddo
DO m=1,npoints(j1,res_index(upstream))
!  if (x_sec(m,j1,upstream) /= x_minelev(j1,upstream)) then
!    y_actlay(m,j1,upstream)=max(y_sec(m,j1,upstream)-dummy2,y_original(m,j1,upstream))
!  endif
!if(j1==11)write(*,'(2I4,3F15.8)')j1,m,y_original(m,j1,upstream),y_actlay(m,j1,upstream),y_sec(m,j1,upstream)
ENDDO

IF (dvol_sed(j1,res_index(upstream)) < 0.) THEN

if (j1/=nbrsec(res_index(upstream))) then
!write(*,'(I4,4F10.2,A10)')j1,minelev_sec(j1,upstream),dheight_sed(j1,upstream),minelev_sec(j1,upstream)-dheight_sed(j1,upstream),minelev_sec(j1+1,upstream)+(mean_slope*dist_sec(j1,upstream)),'eros'
else
!write(*,'(I4,4F10.2,A10)')j1,minelev_sec(j1,upstream),dheight_sed(j1,upstream),minelev_sec(j1,upstream)-dheight_sed(j1,upstream),minelev_sec(j1-1,upstream)-(mean_slope*dist_sec(j1-1,upstream)),'eros'
endif

else

if (j1/=nbrsec(res_index(upstream))) then
!write(*,'(I4,4F10.2,A10)')j1,minelev_sec(j1,upstream),dheight_sed(j1,upstream),minelev_sec(j1,upstream)+dheight_sed(j1,upstream),minelev_sec(j1+1,upstream)+(mean_slope*dist_sec(j1,upstream)),'depo'
else
!write(*,'(I4,4F10.2,A10)')j1,minelev_sec(j1,upstream),dheight_sed(j1,upstream),minelev_sec(j1,upstream)+dheight_sed(j1,upstream),minelev_sec(j1-1,upstream)-(mean_slope*dist_sec(j1-1,upstream)),'depo'
endif

endif

DO m=1,npoints(j1,res_index(upstream))
!write(*,*)j1,y_sec(m,j1,upstream),y_laststep(m,j1,upstream)
ENDDO

!if (d==7)stop

RETURN
END SUBROUTINE change_sec

max_stars_repo_name: plvoit/WASA-SED

max_issues_repo_name: TillF/WASA-SED

max_forks_repo_name: ERottler/WASA-SED


text:
C    Copyright(C) 2008-2017 National Technology & Engineering Solutions
C    of Sandia, LLC (NTESS).  Under the terms of Contract DE-NA0003525 with
C    NTESS, the U.S. Government retains certain rights in this software.
C
C    Redistribution and use in source and binary forms, with or without
C    modification, are permitted provided that the following conditions are
C    met:
C
C    * Redistributions of source code must retain the above copyright
C       notice, this list of conditions and the following disclaimer.
C
C    * Redistributions in binary form must reproduce the above
C      copyright notice, this list of conditions and the following
C      disclaimer in the documentation and/or other materials provided
C      with the distribution.
C
C    * Neither the name of NTESS nor the names of its
C      contributors may be used to endorse or promote products derived
C      from this software without specific prior written permission.
C
C    THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
C    "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
C    LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
C    A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
C    OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
C    SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
C    LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
C    DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
C    THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
C    (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
C    OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
C
      SUBROUTINE MYGET (MYLOC, MNGET, VOID, LVOID, NVOIDS,
     *   CHRCOL, MAXSIZ, LASTER, VROW)
C
      IMPLICIT INTEGER (A-Z)
      INCLUDE 'params.inc'
C
C     This subroutine returns the location (row number) of a void with
C     sufficient space for the memory request.  If necessary, memory is
C     requested from the system.  The memory is contiguous.
C     This routine is to be used only if CHRCOL = 2.
C
C
C***********************************************************************
C
C     MYLOC    Address of internal reference array
C     MNGET    Memory request in character storage units
C     VOID     Void table
C     LVOID    Dimension of void table
C     NVOIDS   Number of voids
               DIMENSION VOID(LVOID,CHRCOL,2), NVOIDS(2)
C     CHRCOL   Column for character tables (must be 2)
C     MAXSIZ   Dimension of static character array.
C     LASTER   Error return
C     VROW     Row number of void which satisfies the memory request
C
C***********************************************************************
C
C     IS THE MEMORY REQUEST SENSIBLE?
C
      IF (MNGET .LT. 0) THEN
         LASTER = BADLEN
         RETURN
      ELSE IF (MNGET .EQ. 0) THEN
         LASTER = SUCESS
         RETURN
      END IF
C
      CALL MXLOOK (MNGET, VOID(1,CHRCOL,1), CHRCOL*LVOID,
     *   NVOIDS(CHRCOL), VROW, LASTER)
      IF (LASTER .EQ. SUCESS) RETURN
C
C     CALL EXTENSION LIBRARY ROUTINE TO GET SPACE FROM SYSTEM.
C
      CALL MYMEMY (MNGET, LOC, MEMRET, MAXSIZ)
      LOC = LOC - MYLOC + 1
C
      IF (MEMRET .LT. 0) THEN
C
C        ILLEGAL MEMORY BLOCK SIZE.
C
         LASTER = ILBLK
         RETURN
C
      END IF
C
C     UPDATE VOID TABLE.
C
      CALL VTABLE (LOC, MEMRET, VOID(1,CHRCOL,1), LVOID,
     *   NVOIDS(CHRCOL), CHRCOL, LASTER)
      IF (LASTER .NE. SUCESS) RETURN
C
      CALL MXLOOK (MNGET, VOID(1,CHRCOL,1), CHRCOL*LVOID,
     *   NVOIDS(CHRCOL), VROW, LASTER)
C
      RETURN
      END

max_stars_repo_name: mathstuf/seacas

max_stars_repo_licenses: BSD-3-Clause NetCDF Zlib MIT

max_issues_repo_name: mathstuf/seacas

max_issues_repo_licenses: BSD-3-Clause NetCDF Zlib MIT

max_forks_repo_name: mathstuf/seacas

max_forks_repo_licenses: BSD-3-Clause NetCDF Zlib MIT


subroutine  VEGVEL(init_vel, PSLCALINIT)
  use simsphere_mod
  implicit none

  real :: SGMA
  integer(kind=1) :: init_vel
  real :: PES, RHOCP, SIGALF, CDL, PI, SDL, rhg, RMRATIODIF, RHA, QSTA
  real :: XLEFN, TFN

  integer :: PSLCALINIT

!      INCLUDE 'modvars.h'

  PES = ( XLAI / 2.0 ) + 1
  RHOCP = CP * PS1 * 100 / ( R * TAF )
  SIGALF = 1 - 0.5 / ( 0.5 + XLAI ) * EXP( -XLAI**2 / 8 )

!     CDL IS THE DRAG COEFFICIENT FOR THE LEAF
!     KEL'S FORMULATION

  CDL = 0.08
  PI = 3.1416

  SDL = XLAI/VEGHEIGHT
  UAF = USTAR * SQRT (PES/(CDL * SDL))

  chf = 0.011 * SQRT(UAF/WIDTH)*(XLAI/PES)

  IF (CHF < 0.001) CHF = 0.001

  RAF = 1 / ( CHF )

!      CHG = ( 1 - SIGALF ) * USTAR**2 / UAF
!        Kell's chg
        chg = Karman**2*uaf/(alog(vegheight*100/0.7))**2
!       Add surface layer resistance to internal air resistance and recompute conductance
  rhg = 1/chg + rtranw
  chg = 1/rhg

!     FOLLOWING IS KEL'S FORMULATION FROM GOUDRIAAN
!     TO GET CHA

!      ALEN = 2 * ((3 * WIDTH**2) / (4 * PI * XLAI/(VEGHEIGHT/2)))
!     /        ** (0.3333)
!     AKCAN = ALEN * 0.5 * UAF
!      RHA = VEGHEIGHT / (2 * AKCAN)
!      CHA = 1 / RHA
!       Above is outdated formulae from Kell. Use old method

  cha = ustar**2 / (uten - uaf )

  IF (CHA <= 0.001) CHA = 0.001

! Convert Q to ELTL and EA.

  QSTF = 10**( 6.1989 - ( 2353 / TF ) )
  SGMA = RHOCP / 0.666
  RMRATIODIF = QSTF - QAF
  IF (RMRATIODIF <= 0) RMRATIODIF = .0001
  VFL = (RMRATIODIF) * PS1 / 0.622

! Call Deardorff (D) or Carlson/Lynn formulation for RST.
! Note RST is actually RF because it contains a cuticular resistance (RCUT).

  IF (STMTYPE == 'D') THEN  ! Deardorff Formulation
    RS = RMIN * (800.0 / (1.0+SOL) + (1.2 * WILT/                       &
         (0.9*W2G+0.1*WGG))**2) ! Stomatal Resistance

    RST = RS * RCUT / (RS + RCUT) * PES / XLAI ! Total Leaf/Canopy Resistance

  ELSE IF (STMTYPE == 'L') THEN ! Lynn and Carlson
    CALL PSLCAL (SGMA, PES, PSLCALINIT)

    RL = RS * RCUT / ( RS + RCUT )   ! Leaf Resistance
    RST = RL * PES / XLAI    ! Total Leaf/Canopy Resistance

  ELSE IF (STMTYPE == 'B') THEN ! Ball-Berry
    CALL BALL (PES, RHA)

    RL = RS * RCUT / ( RS + RCUT )   ! Leaf Resistance
    RST = RL * PES / XLAI    ! Total Leaf/Canopy Resistance

  END IF

  QSTA = 10**( 6.1989 - ( 2353 / TA ) )
  XLEFN = DENS * LE * RMRATIODIF / ( RST + RAF )

!    Average TF and XLEF after initial time step.

  IF (init_vel == 1) THEN
    TF = TAF + ( RNETF - XLEFN ) / (CHF * DENS * CP)
    XLEF = XLEFN
    init_vel = 2
  ELSE
    XLEF = ( XLEF + XLEFN ) / 2
    TFN = TAF + ( RNETF - XLEFN ) / (CHF * DENS * CP)
    TF = (TF + TFN) / 2
  END IF

  RETURN
END


text:
C    @(#)kmpsrt.f	20.3 2/13/96
        function kmpsrt (m,n)
c
        include 'ipfinc/parametr.inc'

        include 'ipfinc/komps.inc'
c
        kmpsrt = ksort(1,m) - ksort(1,n)
        if (kmpsrt .eq. 0) kmpsrt = ksort(2,m) - ksort(2,n)
        if (kmpsrt .eq. 0) kmpsrt = ksort(3,m) - ksort(3,n)
        return
        end

max_stars_repo_name: mbheinen/bpa-ipf-tsp

max_issues_repo_name: cuihantao/bpa-ipf-tsp

max_forks_repo_name: mbheinen/bpa-ipf-tsp


text:
	SUBROUTINE GINDEX ( i, j, ij, ndir, ffi1, ffj1, lin1,
     +			    ffi2, ffj2, lin2, iret )
C************************************************************************
C* GINDEX                                                               *
C*                                                                      *
C* This subroutine returns intersection information given an (i,j)	*
C* location and a nominal direction.  If no intersection in direction	*
C* NDIR, locations will be set to zero with values set to RMISSD.	*
C* Note: On a compass, NDIR is as follows:				*
C*	1 - NE		2 - E		3 - SE		4 - S		*
C*	5 - SW		6 - W		7 - NW		8 - N		*
C*                                                                      *
C* GINDEX ( I, J, IJ, NDIR, FFI1, FFJ1, LIN1, FFI2, FFJ2, LIN2, IRET )	*
C*                                                                      *
C* Input parameters:                                                    *
C*      I		INTEGER		I location on grid		*
C*      J		INTEGER		J location on grid		*
C*      IJ		INTEGER		Index pointer into int arrays	*
C*      NDIR		INTEGER		Direction (1 thru 8)		*
C*                                                                      *
C* Output parameters:                                                   *
C*      FFI1		REAL		I location of first int		*
C*      FFJ1		REAL		J location of first int		*
C*      LIN1		INTEGER		Line number at first int	*
C*      FFI2		REAL		I location of second int	*
C*      FFJ2		REAL		J location of second int	*
C*      LIN2		INTEGER		Line number at second int	*
C*      IRET            INTEGER         Return code                     *
C**                                                                     *
C* Log:                                                                 *
C* D.W.Plummer/NCEP      9/98                                           *
C* D.W.Plummer/NCEP      7/03	Chgs for intersection efficiency	*
C* D.W.Plummer/NCEP      9/03	Return line number instead of line val	*
C************************************************************************
C
        INCLUDE         'GEMPRM.PRM'
        INCLUDE         'grphgd.cmn'
C
	iret = 0
C
	ffi1 = RMISSD
	ffj1 = RMISSD
	lin1 = IMISSD
	ffi2 = RMISSD
	ffj2 = RMISSD
	lin2 = IMISSD
C
	IF ( ndir .eq. 1 )  THEN
C
		IF ( intptrs(ij,3,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,3,STPTR)+intptrs(ij,3,NINTS)-1,
     &		        intptrs(ij,3,STPTR), -1
C
			IF ( i .lt. intsct(n,1) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 2 )  THEN
C
		IF ( intptrs(ij,2,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,2,STPTR)+intptrs(ij,2,NINTS)-1,
     &                  intptrs(ij,2,STPTR), -1
C
			IF ( i .lt. intsct(n,1) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 3 )  THEN
C
		IF ( intptrs(ij,4,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,4,STPTR)+intptrs(ij,4,NINTS)-1,
     &                  intptrs(ij,4,STPTR), -1
C
			IF ( i .lt. intsct(n,1) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 4 )  THEN
C
		IF ( intptrs(ij,1,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,1,STPTR),
     &                  intptrs(ij,1,STPTR)+intptrs(ij,1,NINTS)-1
C
			IF ( j .gt. intsct(n,2) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 5 )  THEN
C
		IF ( intptrs(ij,3,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,3,STPTR),
     &                  intptrs(ij,3,STPTR)+intptrs(ij,3,NINTS)-1
C
			IF ( i .gt. intsct(n,1) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 6 )  THEN
C
		IF ( intptrs(ij,2,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,2,STPTR),
     &                  intptrs(ij,2,STPTR)+intptrs(ij,2,NINTS)-1
C
			IF ( i .gt. intsct(n,1) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 7 )  THEN
C
		IF ( intptrs(ij,4,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,4,STPTR),
     &                  intptrs(ij,4,STPTR)+intptrs(ij,4,NINTS)-1
C
			IF ( i .gt. intsct(n,INT_X) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	ELSE IF ( ndir .eq. 8 )  THEN
C
		IF ( intptrs(ij,1,NINTS) .eq. 0 )  RETURN
C
		DO  n = intptrs(ij,1,STPTR)+intptrs(ij,1,NINTS)-1,
     &                  intptrs(ij,1,STPTR), -1
C
			IF ( j .lt. intsct(n,INT_Y) )  THEN
C
				ffi2 = ffi1
				ffj2 = ffj1
				lin2 = lin1
				ffi1 = intsct(n,INT_X)
				ffj1 = intsct(n,INT_Y)
				lin1 = intinfo(n,INT_LINE)
C
			END IF
C
		END DO
C
	END IF
C
	RETURN
	END

max_stars_repo_name: oxelson/gempak

max_issues_repo_name: oxelson/gempak

max_forks_repo_name: oxelson/gempak