! begin: keyword_control_lib.f
! copyright 2005,  J. E. Akin, all rights reserved.

  subroutine apply_key (key)
!---------------------------------------------------------
! free format input of application data
!---------------------------------------------------------
   use system_constants  ! for data definable here
   use keyword_buffer    ! for MAX_KEY, WORD_SIZE, ECHO_KEY
   implicit none

   character(len=MAX_KEY), intent (in) :: key     
   character(len=WORD_SIZE) :: on_off   
   character(len=WORD_SIZE) :: word    
   integer, parameter :: limit = 10 ! max allowed bad words
   integer, save      :: count = 0  ! number of bad words
   integer :: J, K  ! loops
   !b real(dp), allocatable :: body_force (:)

    word   = 'null'
    on_off = 'off'
    select case ( key ) 

! *** add any application specific control keyword case here ***
!  (they must be previously defined in module System_Constants)
!     Case ('
!       *** end of application dependent case actions ***

! Check allowed dummy keys of lines to skip
      Case (' ' )        ! ignore blank line
      Case ('#' )        ! ignore comment line
      Case ('!' )        ! ignore comment line
      Case ('?' )        ! ignore comment line
          
! problem title
      Case ('title' )         ; call get_string (title)    

! program logical controls
      Case ('average_mass')    ; AVERAGE_MASS   = .TRUE.
      Case ('axisymmetric')    ; AXISYMMETRIC   = .TRUE.
      Case ('bar_chart'    )   ; BAR_CHART      = .TRUE.
                                 NO_DIST        = 1
      Case ('buckling'     )   ; BUCKLING       = .TRUE.
                                 GEOMETRIC_K    = .TRUE.
      Case ('debug_adapt'  )   ; DEBUG_ADAPT    = .TRUE.
      Case ('debug_all'    )   ; DEBUG_ADAPT    = .TRUE.
                                 DEBUG_ASSEMBLY = .TRUE.
                                 DEBUG_B        = .TRUE.
                                 DEBUG_E        = .TRUE.
                                 DEBUG_EL_COL   = .TRUE.
                                 DEBUG_EL_SQ    = .TRUE.
                                 DEBUG_EL_POST  = .TRUE.
                                 DEBUG_EL_TYPE  = .TRUE.
                                 DEBUG_ERROR    = .TRUE.
                                 DEBUG_EXACT    = .TRUE.
                                 DEBUG_INCLUDE  = .TRUE.
                                 DEBUG_INERTIA  = .TRUE.
                                 DEBUG_MIX_SQ   = .TRUE.
                                 DEBUG_POST_EL  = .TRUE.
                                 DEBUG_PROPERTY = .TRUE.
                                 DEBUG_SCP      = .TRUE.
                                 DEBUG_SEG_COL  = .TRUE.
                                 I_BUG          = 1
      Case ('debug_assembly')  ; DEBUG_ASSEMBLY = .TRUE.
      Case ('debug_b'      )   ; DEBUG_B        = .TRUE.
      Case ('debug_e'      )   ; DEBUG_E        = .TRUE.
      Case ('debug_el_col' )   ; DEBUG_EL_COL   = .TRUE.
      Case ('debug_el_sq'  )   ; DEBUG_EL_SQ    = .TRUE.
      Case ('debug_el_post')   ; DEBUG_EL_POST  = .TRUE.
      Case ('debug_el_type')   ; DEBUG_EL_TYPE  = .TRUE.
      Case ('debug_error'  )   ; DEBUG_ERROR    = .TRUE.
      Case ('debug_exact'  )   ; DEBUG_EXACT    = .TRUE.
      Case ('debug_include')   ; DEBUG_INCLUDE  = .TRUE.
      Case ('debug_inertia')   ; DEBUG_INERTIA  = .TRUE.
      Case ('debug_mixed'  )   ; DEBUG_MIX_SQ   = .TRUE.
      Case ('debug_post_el')   ; DEBUG_POST_EL  = .TRUE.
      Case ('debug_property')  ; DEBUG_PROPERTY = .TRUE.
      Case ('debug_scp'    )   ; DEBUG_SCP      = .TRUE.
      Case ('debug_segment')   ; DEBUG_SEG_COL  = .TRUE.
      Case ('debug_units'  )   ; DEBUG_UNITS    = .TRUE.
      Case ('diagonal_mass')   ; DIAGONAL_MASS  = .TRUE.
      Case ('dof_vary')        ; DOF_VARY       = .TRUE.
      Case ('dynamic'     )    ; DYNAMIC        = .TRUE.
                                 EIGEN          = .FALSE.
                                 JACOBI         = .FALSE.
                                 STATIC         = .FALSE.
                                 TRANSIENT      = .FALSE.
      Case ('echo_key'    )    ; ECHO_KEY       = .TRUE.
      Case ('eigen'       )    ; EIGEN          = .TRUE.
                                 DYNAMIC        = .FALSE.
                                 STATIC         = .FALSE.
                                 TRANSIENT      = .FALSE.
      Case ('eigen_root'  )    ; EIGEN_ROOT     = .TRUE.
      Case ('el_damping'  )    ; EL_DAMPING     = .TRUE.
      Case ('fourier_law')     ; FOURIER_LAW    = .TRUE. 
      Case ('Fourier_law')     ; FOURIER_LAW    = .TRUE. 
      Case ('from_rest'  )     ; FROM_REST      = .TRUE. 
      Case ('f95'        )     ; F95            = .TRUE.
      Case ('geom_iter'  )     ; GEOM_ITER      = .TRUE.
      Case ('geometric_k')     ; GEOMETRIC_K    = .TRUE.
      Case ('get_deleted')     ; GET_DELETED_EL = .TRUE.
      Case ('grad_base_error') ; GRAD_BASE_ERROR= .TRUE.
      Case ('hughes_method'  ) ; HUGHES_METHOD  = .TRUE.
                                 TIME_METHOD    = 3
                                 NEWMARK_METHOD = .FALSE.
                                 WILSON_METHOD  = .FALSE.
      Case ('incremental')     ; INCREMENTAL    = .TRUE.
      Case ('initial_function'); INITIAL_FUNCTION = .TRUE.
      Case ('input_only' )     ; INPUT_ONLY     = .TRUE.
      Case ('is_null'    )     ; IS_NULL        = .TRUE.
      Case ('jacobi'     )     ; JACOBI         = .TRUE.
                                 DYNAMIC        = .FALSE.
                                 STATIC         = .FALSE.
                                 TRANSIENT      = .FALSE.
      Case ('list_el_to_el')   ; LIST_EL_TO_EL  = .TRUE.
      Case ('list_exact')      ; USE_EXACT      = .TRUE.
      Case ('list_exact_flux') ; USE_EXACT_FLUX = .TRUE.
      Case ('list_qp_flux')    ; LIST_QP_FLUX   = .TRUE.
      Case ('masses')          ; PT_MASSES      = .TRUE.  
      Case ('newmark_method')  ; NEWMARK_METHOD = .TRUE.  
                                 TIME_METHOD     = 2
                                 HUGHES_METHOD  = .FALSE.
                                 WILSON_METHOD  = .FALSE.
      Case ('nodal_tau')       ; NODAL_TAU      = .TRUE.  
      Case ('no_bc_echo')      ; ECHO_BC        = .FALSE. 
      Case ('no_el_echo')      ; ECHO_EL        = .FALSE. 
      Case ('no_flux_echo')    ; ECHO_FLUX      = .FALSE. 
      Case ('no_input_echo')   ; ECHO_BC        = .FALSE. 
                                 ECHO_EL        = .FALSE. 
                                 ECHO_FLUX      = .FALSE.
                                 ECHO_INPUT     = .FALSE. 
                                 ECHO_MIXED     = .FALSE.
                                 ECHO_PROP      = .FALSE. 
                                 ECHO_PTS       = .FALSE.
                                 ECHO_RHS       = .FALSE.
                                 ECHO_START     = .FALSE.
      Case ('no_mixed_echo')   ; ECHO_MIXED     = .FALSE. 
      Case ('no_printing')     ; NO_PRINTING    = .TRUE.  
      Case ('no_prop_echo')    ; ECHO_PROP      = .FALSE. 
      Case ('no_pt_echo')      ; ECHO_PTS       = .FALSE.
      Case ('no_rhs_echo')     ; ECHO_RHS       = .FALSE.
      Case ('no_start_echo')   ; ECHO_START     = .FALSE.
      Case ('omit_prints')     ; NO_PRINTING    = .TRUE.  
      Case ('plane_strain')    ; PLANE_STRAIN   = .TRUE.
      Case ('polar_coord')     ; POLAR_COORD    = .TRUE.
      Case ('property_first')  ; PROPERTY_FIRST = .TRUE.
      Case ('pt_el_sums')      ; PT_EL_SUMS     = .TRUE.
      Case ('pt_el_list')      ; PT_EL_LIST     = .TRUE.
      Case ('pt_masses')       ; PT_MASSES      = .TRUE.  
      Case ('pt_damping')      ; PT_DAMPING     = .TRUE.  
      Case ('ramp_iter' )      ; RAMP_ITER      = .TRUE.
      Case ('save_new_mesh')   ; SAVE_NEW_MESH  = .TRUE.
      Case ('save_1248')       ; SAVE_1248      = .TRUE.
      Case ('scp_center_only') ; SCP_CENTER_ONLY= .TRUE.
      Case ('scp_center_no')   ; SCP_CENTER_ONLY= .FALSE.
      Case ('scp_dof')         ; SCP_DOF        = .TRUE.
      Case ('scp_e_at_np')     ; SCP_E_AT_NP    = .TRUE.
      Case ('scp_e_at_qp')     ; SCP_E_AT_QP    = .TRUE.
      Case ('scp_e_smooth')    ; SCP_E_SMOOTH   = .TRUE.
      Case ('scp_neigh_el')    ; SCP_NEIGH_EL   = .TRUE.
                                 SCP_NEIGH_FACE = .FALSE.
                                 SCP_NEIGH_PT   = .FALSE.
                                 N_PATCH        = N_ELEMS
      Case ('scp_neigh_face')  ; SCP_NEIGH_FACE = .TRUE.
                                 SCP_NEIGH_EL   = .FALSE.
                                 SCP_NEIGH_PT   = .FALSE.
                                 N_PATCH        = N_ELEMS
      Case ('scp_neigh_pt')    ; SCP_NEIGH_PT   = .TRUE.
                                 SCP_NEIGH_FACE = .FALSE.
                                 SCP_NEIGH_EL   = .FALSE.
                                 N_PATCH        = MAX_NP 
      Case ('scp_only_once')   ; SCP_ONLY_ONCE  = .TRUE.
      Case ('scp_not_once')    ; SCP_ONLY_ONCE  = .FALSE.
      Case ('scp_2nd_deriv')   ; SCP_2ND_DERIV  = .TRUE.
                                 CALL LIST_2ND_DERIV_ORDER
      Case ('scp_3rd_deriv')   ; SCP_3RD_DERIV  = .TRUE.
      Case ('shock_capture')   ; TAU_RGN        = .TRUE.
      Case ('skip_error')      ; SKIP_ERROR     = .TRUE.
      Case ('space_time')      ; SPACE_TIME     = .TRUE.
      Case ('s_t_continuous')  ; S_T_CONTINUOUS = .TRUE.
                                 SPACE_TIME     = .TRUE.
      Case ('static')          ; STATIC         = .TRUE.
                                 DYNAMIC        = .FALSE.
                                 EIGEN          = .FALSE.
                                 TRANSIENT      = .FALSE.
      Case ('supg')            ; SUPG           = .TRUE.
                               ! TAU_S1         = .TRUE.  ! default
      Case ('supg_2nd_deriv')  ; SUPG           = .TRUE.
                                 GET_2ND_DERIV  = .TRUE.
                                 CALL LIST_2ND_DERIV_ORDER
      Case ('symmetric')       ; SYMMETRIC      = .TRUE.  ; MODE = 0
      Case ('tau_box' )        ; TAU_BOX        = .TRUE.
      Case ('tau_norm')        ; TAU_S1         = .TRUE.
      Case ('tau_geom')        ; TAU_GEOM       = .TRUE.
      Case ('tau_pec')         ; TAU_PEC        = .TRUE.
      Case ('tau_qp')          ; TAU_QP         = .TRUE.
      Case ('tau_rgn')         ; TAU_RGN        = .TRUE.
      Case ('tau_show')        ; TAU_SHOW       = .TRUE.
      Case ('tau_s1')          ; TAU_S1         = .TRUE.
      Case ('tau_ugn')         ; TAU_UGN        = .TRUE.
      Case ('tau_vol')         ; TAU_VOL        = .TRUE.
      !b Case ('time_slab')       ; TIME_SLAB      = .TRUE.
      Case ('transient')       ; TRANSIENT      = .TRUE.
                                 DYNAMIC        = .FALSE.
                                 EIGEN          = .FALSE.
                                 STATIC         = .FALSE.
      Case ('unsymmetric')     ; SYMMETRIC      = .FALSE. ; MODE = 1
      Case ('upwind')          ; UPWIND         = .TRUE.
      Case ('use_exact_bc')    ; USE_EXACT_BC   = .TRUE.
      Case ('use_exact_flux')  ; USE_EXACT_FLUX = .TRUE.
      Case ('use_exact_mixed') ; USE_EXACT_ROBIN = .TRUE.
      Case ('use_exact_robin') ; USE_EXACT_ROBIN = .TRUE.
      Case ('use_exact_source'); USE_EXACT_SOURCE = .TRUE.
      Case ('user_logic')      ; USER_LOGIC     = .TRUE.
      Case ('wilson_method')   ; WILSON_METHOD  = .TRUE.
                                 TIME_METHOD    = 1
                                 NEWMARK_METHOD = .FALSE.
                                 HUGHES_METHOD  = .FALSE.

! program control integers
      Case ('data_set')   ; call get_int (DATA_SET)
      Case ('DEG_HOMO')   ; call get_int (DEG_HOMO)
      Case ('eigen_post') ; call get_int (EIGEN_POST)
      Case ('eigen_scp')  ; call get_int (EIGEN_SCP) 
      Case ('eigen_show') ; call get_int (EIGEN_SHOW)
      Case ('EXACT_FL')   ; call get_int (EXACT_FL)
      Case ('exact_reals'); call get_int (EXACT_FL)
                            EXACT_REALS = EXACT_FL
      Case ('EXACT_FX')   ; call get_int (EXACT_FX)
      Case ('exact_int')  ; call get_int (EXACT_FX)
                            EXACT_INT = EXACT_FX
      Case ('example')    ; call get_int (EXAMPLE)
        PRINT *, 'NOTE: USING SOURCE EXAMPLE LIBRARY NUMBER ', EXAMPLE
      Case ('exact_case') ; call get_int (EXACT_CASE)
        PRINT *, 'NOTE: USING EXACT CASE SOLUTION NUMBER ', EXACT_CASE
      Case ('GET_DAT')    ; call get_int (GET_DAT)
      Case ('IN_RHS')     ; call get_int (IN_RHS)
      Case ('loads')      ; IN_RHS = 1
      Case ('I_BUG')      ; call get_int (I_BUG)
      Case ('debug' )     ; I_BUG = 1
      Case ('I_SAY')      ; call get_int (I_SAY)
      Case ('remarks')    ; call get_int (I_SAY)
      Case ('LINE_WARN')  ; call get_int (LINE_WARN)
      Case ('LEM_WRT')    ; call get_int (LEM_WRT)
      Case ('el_list')    ; LEM_WRT = 1
      Case ('LIST_LR')    ; call get_int (LIST_LR)
      !b Case ('el_react')   ; LIST_LR = 1 ; L_REACT = 21
      Case ('el_react')      ; L_REACT = 21
      Case ('list_el_react') ; LIST_LR        = 1      !b
      Case ('L_B_N')      ; call get_int (L_B_N)
      Case ('el_segment') ; call get_int (L_B_N) 
        IF ( L_B_N > M_B_N )  M_B_N = L_B_N
      Case ('continuity') ; call get_int (L_CONT)
      Case ('L_CONT')     ; call get_int (L_CONT)
      Case ('L_HOMO')     ; call get_int (L_HOMO)
      Case ('el_homo')    ; L_HOMO = 1
      Case ('L_REACT')    ; call get_int (L_REACT)
      Case ('L_SHAPE')    ; call get_int (L_SHAPE)
      Case ('shape'  )    ; call get_int (L_SHAPE)
      Case ('line_el')        ; L_SHAPE = 1
      Case ('lines')          ; L_SHAPE = 1
      Case ('tri_el')         ; L_SHAPE = 2
      Case ('triangles')      ; L_SHAPE = 2
      Case ('quad_el')        ; L_SHAPE = 3
      Case ('quads')          ; L_SHAPE = 3
      Case ('quadrilaterals') ; L_SHAPE = 3
      Case ('hex_el')         ; L_SHAPE = 4
      Case ('hexes')          ; L_SHAPE = 4
      Case ('hexahedrons')    ; L_SHAPE = 4
      Case ('tet_el')         ; L_SHAPE = 5
      Case ('tets')           ; L_SHAPE = 5
      Case ('tetrahedrons')   ; L_SHAPE = 5
      Case ('vary_e')         ; VARY_E = .TRUE.
      Case ('wedge_el')       ; L_SHAPE = 6
      Case ('wedges')         ; L_SHAPE = 6
      Case ('mat_real')   ; call get_int (MAT_FLO)
                            mat_real = MAT_FLO ; IP_TEST = 1
      !b Case ('MAX_DAT')    ; call get_int (MAX_DAT) ! >= 8
      Case ('MAX_NP')     ; call get_int (MAX_NP)
      Case ('nodes' )     ; call get_int (MAX_NP)
      Case ('MAX_TYP')    ; call get_int (MAX_TYP)
      Case ('bc_types')   ; call get_int (MAX_TYP)
      Case ('MISC_FL')    ; call get_int (MISC_FL)
      Case ('reals'  )    ; call get_int (MISC_FL)
                            reals = MISC_FL
      Case ('MISC_FX')    ; call get_int (MISC_FX)
      Case ('integers')   ; call get_int (MISC_FX)
                            integers = MISC_FX
      Case ('MODE')       ; call get_int (MODE)
      Case ('M_B_N')      ; call get_int (M_B_N)
      Case ('M_COEFF')    ; call get_int (M_COEFF)
      Case ('M_SHOW')     ; call get_int (M_SHOW)
      Case ('NEEDS' )     ; call get_int (NEEDS) 
      Case ('face_nodes') ; call get_int (NEEDS) 
      Case ('all_el_to_el'); NEEDS = 1
      !b Case ('NEIGH_L')    ; call get_int (NEIGH_L)
      !b Case ('NEIGH_N')    ; call get_int (NEIGH_N)
      Case ('NOD_PER_EL') ; call get_int (NOD_PER_EL)
      Case ('el_nodes'  ) ; call get_int (NOD_PER_EL)
      !b Case ('NOT_SYM')    ; call get_int (NOT_SYM)
      Case ('NO_DIST')    ; call get_int (NO_DIST)
      Case ('bar_short')  ; NO_DIST   = 1
                            BAR_CHART = .TRUE.
      Case ('bar_long' )  ; NO_DIST   = 0
                            BAR_CHART = .TRUE.
      Case ('NPT_WRT')    ; call get_int (NPT_WRT)
      Case ('pt_list')    ; NPT_WRT = 1
      Case ('NUL_COL')    ; call get_int (NUL_COL)
      Case ('el_no_col')  ; NUL_COL = 1
      Case ('el_analysis'); call get_int (N_ANAL)
      Case ('N_ANAL')     ; call get_int (N_ANAL)
      Case ('N_BS_FIX')   ; call get_int (N_BS_FIX)
      Case ('seg_int' )   ; call get_int (N_BS_FIX)
                            seg_int = N_BS_FIX ; IP_TEST = 1
      Case ('N_BS_FLO')   ; call get_int (N_BS_FLO)
      Case ('seg_real')   ; call get_int (N_BS_FLO)
                            seg_real = N_BS_FLO ; IP_TEST = 1
      Case ('N_BUG')      ; call get_int (N_BUG)
      Case ('checker')    ; call get_int (N_BUG)
      !b Case ('N_COEFF')    ; call get_int (N_COEFF)
      Case ('N_CORNER')   ; call get_int (N_CORNER)
      Case ('el_corners') ; call get_int (N_CORNER)
      Case ('N_CRD')      ; call get_int (N_CRD)
      Case ('reader')     ; call get_int (N_CRD)
      !b Case ('N_D_FRE')    ; call get_int (N_D_FRE)
      Case ('N_ELEMS')    ; call get_int (N_ELEMS)
      Case ('elems'  )    ; call get_int (N_ELEMS)
      Case ('N_EL_FRE')   ; call get_int (N_EL_FRE)
      Case ('N_GEOM')     ; call get_int (N_GEOM)
      Case ('el_geom')    ; call get_int (N_GEOM)
      Case ('N_G_DOF')    ; call get_int (N_G_DOF)
      Case ('dof'    )    ; call get_int (N_G_DOF)
      Case ('N_G_FLUX')   ; call get_int (N_G_FLUX)
      Case ('seg_pt_flux'); call get_int (N_G_FLUX)
      Case ('N_HOMO')     ; call get_int (N_HOMO)
      Case ('pt_homo')    ; N_HOMO = 1
      Case ('N_ITER')     ; call get_int (N_ITER)
      Case ('iters' )     ; call get_int (N_ITER)
      Case ('N_LP_FIX')   ; call get_int (N_LP_FIX)
      Case ('el_int'  )   ; call get_int (N_LP_FIX)
                            el_int = N_LP_FIX ; IP_TEST = 1
                            LP_TEST = 1
      Case ('N_LP_FLO')   ; call get_int (N_LP_FLO)
      Case ('el_real' )   ; call get_int (N_LP_FLO)
                            el_real = N_LP_FLO ; IP_TEST = 1
                            LP_TEST = 1
      !b Case ('N_L_DEL')    ; call get_int (N_L_DEL)
      Case ('N_L_TYPE')   ; call get_int (N_L_TYPE)
      Case ('el_types')   ; call get_int (N_L_TYPE)
      Case ('N_MAT')      ; call get_int (N_MAT)
      Case ('materials')  ; call get_int (N_MAT)
      Case ('N_MIXED')    ; call get_int (N_MIXED)
      Case ('mixed_segs') ; call get_int (N_MIXED)
      Case ('N_MX_FIX')   ; call get_int (N_MX_FIX)
      Case ('mixed_int')  ; call get_int (N_MX_FIX)
                            mixed_int = N_MX_FIX ; IP_TEST = 1
      Case ('N_MX_FLO')   ; call get_int (N_MX_FLO)
      Case ('mixed_real') ; call get_int (N_MX_FLO)
                            mixed_real = N_MX_FLO ; IP_TEST = 1
      Case ('N_NP_FIX')   ; call get_int (N_NP_FIX)
      Case ('pt_int'  )   ; call get_int (N_NP_FIX)
                            pt_int = N_NP_FIX ; IP_TEST = 1
      Case ('N_NP_FLO')   ; call get_int (N_NP_FLO)
      Case ('pt_real' )   ; call get_int (N_NP_FLO)
                            pt_real = N_NP_FLO ; IP_TEST = 1
      Case ('N_PARM')     ; call get_int (N_PARM)
      Case ('el_space')   ; call get_int (N_PARM)
      Case ('N_PRT')      ; call get_int (N_PRT)
      Case ('printer')    ; call get_int (N_PRT)
      Case ('N_QP')       ; call get_int (N_QP)
      Case ('gauss')      ; call get_int (N_QP)
      Case ('Gauss')      ; call get_int (N_QP)
      Case ('N_QP_C')     ; call get_int (N_QP_C)
      Case ('N_QP_R')     ; call get_int (N_QP_R)
      Case ('N_R_B')      ; call get_int (N_R_B)
      Case ('b_rows')     ; call get_int (N_R_B)
      Case ('N_F_SEG')    ; call get_int (N_F_SEG)
      Case ('segments')   ; call get_int (N_F_SEG)
      Case ('N_SPACE')    ; call get_int (N_SPACE)
      Case ('space'  )    ; call get_int (N_SPACE)
      !b Case ('N_STORE')    ; call get_int (N_STORE)
      Case ('N_FILE1')    ; call get_int (N_FILE1)
      Case ('post  ' )    ; N_FILE1 = 32
      Case ('post_el')    ; N_FILE1 = 32
                            POST_EL = .TRUE.
      Case ('post_1' )    ; N_FILE1 = 32
                            POST_EL = .TRUE.
      Case ('N_FILE2')    ; call get_int (N_FILE2)
      Case ('post_mixed') ; N_FILE2 = 33
                            POST_MIXED = .TRUE.
      Case ('post_2' )    ; N_FILE2 = 33
                          !bPOST_MIXED = .TRUE.
      Case ('N_FILE3')    ; call get_int (N_FILE3)
      Case ('post_3' )    ; N_FILE3 = 34
      Case ('N_FILE4')    ; call get_int (N_FILE4)
      Case ('post_4' )    ; N_FILE4 = 35
      Case ('N_FILE5')    ; call get_int (N_FILE5)
      Case ('post_5' )    ; N_FILE5 = 36
      Case ('N_2_DER')    ; call get_int (N_2_DER)
      Case ('N_3_DER')    ; call get_int (N_3_DER)
      Case ('restart_step') ; call get_int (RESTART_STEP)
      Case ('SCP_DEG')    ; call get_int (SCP_DEG)
      Case ('SCP_FIT')    ; call get_int (SCP_FIT)
      Case ('SCP_INC')    ; call get_int (SCP_INC)
      Case ('scp_deg_inc') ; call get_int (SCP_INC)
      Case ('SCP_PAD')    ; call get_int (SCP_PAD)
      Case ('jacobi_sweeps'); call get_int (SWEEPS)
! test multiple inputs per line (two integers)          !b
      Case ('test_2_i')   ; call get_int (T_SETS)       !b
                            print *,'T_SETS ', T_SETS   !b
                            call get_int (T_STEPS)      !b
                            print *,'T_STEPS ', T_STEPS !b
! test element types data
      Case ('type_nodes')  ! LT_DATA (1, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_NODES (J))
          END DO
          !b IF ( NOD_PER_EL == 0 )           &
          IF ( MAXVAL (TYPE_NODES) > NOD_PER_EL ) &
               NOD_PER_EL = MAXVAL (TYPE_NODES)
      Case ('type_gauss', 'type_Gauss')  ! LT_DATA (2, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_GAUSS (J))
          END DO
          IF ( N_QP == 0 )           &
            N_QP = MAXVAL (TYPE_GAUSS)
      Case ('type_geom')  ! LT_DATA (3, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_GEOM (J))
          END DO
          IF ( N_GEOM == 0 )          &
            N_GEOM = MAXVAL (TYPE_GEOM)
      Case ('type_parm' )  ! LT_DATA (4, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_PARM  (J))
          END DO
          IF (  N_PARM == 0 )         &
            N_PARM = MAXVAL (TYPE_PARM)
      Case ('type_shape')  ! LT_DATA (5, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_SHAPE (J))
          END DO
          IF ( L_SHAPE == 0 )           &
            L_SHAPE = MAXVAL (TYPE_SHAPE)
      Case ('type_faces')  ! LT_DATA (6, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_FACES (J))
          END DO
          IF ( MAX_FACES == 0 )           &
            MAX_FACES = MAXVAL (TYPE_FACES)
      Case ('type_sons' )  ! LT_DATA (7, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_SONS  (J))
          END DO
      Case ('type_user' )  ! LT_DATA (8, J)
        IF ( .NOT. TYPES_ALLOCATED ) CALL ALLOCATE_TYPE_CONTROLS
          DO J = 1, N_L_TYPE
            call get_int ( TYPE_USER  (J))
          END DO

! time marching data
      Case ('history_dof')     ; call get_int (HISTORY_DOF)
      Case ('history_node')    ; call get_int (HISTORY_NODE)
      Case ('inc_save')        ; call get_int (INC_SAVE)    
      Case ('time_groups')     ; call get_int (T_SETS) 
      Case ('T_SETS' )         ; call get_int (T_SETS) 
      Case ('T_STEPS')         ; call get_int (T_STEPS)
      Case ('time_steps')      ; call get_int (T_STEPS)

      Case ('initial_value')   ; call get_real (INITIAL_VALUE)
      Case ('start_value'  )   ; call get_real (INITIAL_VALUE)
      Case ('ramp_down')       ; call get_real (RAMP_DOWN_TIME)
      Case ('ramp_up')         ; call get_real (RAMP_UP_TIME)
      Case ('START_TIME'   )   ; call get_real (START_TIME)
      Case ('start_time'   )   ; call get_real (START_TIME)
      Case ('time_method')     ; call get_int (TIME_METHOD)
      Case ('crank-nicolson')  ; TIME_METHOD = 2
      Case ('euler_forward')   ; TIME_METHOD = 1
      Case ('galerkin_time')   ; TIME_METHOD = 3
      Case ('least_sq_time')   ; TIME_METHOD = 4
      Case ('trapezoidal_coeff')
             call get_real (TRAPEZOIDAL_COEFF)
             TIME_METHOD = 5

! unit for XYZ, E, B for SCP flux and/or post-processing flux
      Case ('U_FLUX')     ; call get_int (U_FLUX)
! unit for element integral of H
      Case ('U_INTGR')    ; call get_int (U_INTGR)

! "plot" units, for ascii tmp files
      Case ('U_PLT1')          ; call get_int (U_PLT1)
      Case ('save_bc_xyz')     ; U_PLT1 = 22
      Case ('omit_bc_xyz')     ; U_PLT1 = 0
      Case ('U_PLT2')          ; call get_int (U_PLT2)
      Case ('save_el_topo')    ; U_PLT2 = 23
      Case ('omit_el_topo')    ; U_PLT2 = 0
      Case ('U_PLT3')          ; call get_int (U_PLT3)
      Case ('save_pt_ans')     ; U_PLT3 = 24
      Case ('omit_pt_ans')     ; U_PLT3 = 0
      Case ('U_PLT4')          ; call get_int (U_PLT4)
      Case ('save_qp_flux')    ; U_PLT4 = 25
      Case ('omit_qp_flux')    ; U_PLT4 = 0 
      Case ('U_PLT5')          ; call get_int (U_PLT5)
      Case ('save_pt_aves')    ; U_PLT5 = 26
      Case ('omit_pt_aves')    ; U_PLT5 = 0
      Case ('U_PLT6')          ; call get_int (U_PLT6)
      Case ('save_el_err')     ; U_PLT6 = 27
      Case ('omit_el_err')     ; U_PLT6 = 0
      Case ('U_PLT7')          ; call get_int (U_PLT7)
      Case ('save_pt_err')     ; U_PLT7 = 28
      Case ('omit_pt_err')     ; U_PLT7 = 0
      Case ('U_PLT8')          ; call get_int (U_PLT8)
      Case ('save_exact')      ; U_PLT8 = 34
                                USE_EXACT = .TRUE.
      Case ('omit_exact')      ; U_PLT8 = 0
      Case ('U_PLT9')          ; call get_int (U_PLT9)
      Case ('save_exact_flux') ; U_PLT9 = 35
      Case ('omit_exact_flux') ; U_PLT9 = 0
      Case ('U2_PLT1')         ; call get_int (U2_PLT1)
      Case ('save_pt_flux_grad'); U2_PLT1 = 51
      Case ('omit_pt_flux_grad'); U2_PLT1 = 0
      Case ('U2_PLT2')         ; call get_int (U2_PLT2)
      Case ('save_ex_flux_grad'); U2_PLT2 = 52
      Case ('omit_ex_flux_grad'); U2_PLT2 = 0
      Case ('U2_PLT3')         ; call get_int (U2_PLT3)
   !  Case ('')     ; U2_PLT3 = 53
      Case ('U2_PLT4')         ; call get_int (U2_PLT4)
   !  Case ('')    ; U_PLT4 = 54
      !b XXX Case ('omit_saves')  ! tmp files for plots
        !b XXX U_PLT1 = 0 ; U_PLT2 = 0 ; U_PLT3 = 0 ; U_PLT4 = 0 ;
        !b XXX U_PLT5 = 0 ; U_PLT6 = 0 ; U_PLT7 = 0 ; U_PLT8 = 0 ;
        !b XXX U_PLT9 = 0
        !b XXX U2_PLT1 = 0 ; U2_PLT2 = 0 ; U2_PLT3 = 0 ; U2_PLT4 = 0 
      Case ('save_all')  ! tmp files for plots
        U_PLT1 = 22; U_PLT2 = 23; U_PLT3 = 24; U_PLT4 = 25;
        U_PLT5 = 26; U_PLT6 = 27; U_PLT7 = 28; U_PLT8 = 34;
        U_PLT9 = 35
        U2_PLT1 = 51 ; U2_PLT2 = 52 ; U2_PLT3 = 53 ; U2_PLT4 = 54

! super convergent patch units
      Case ('U_2_DER'     )   ; call get_int (U_2_DER) ! sec deriv
      Case ('U_3_DER'     )   ; call get_int (U_3_DER) ! 3rd deriv
      Case ('U_SCPR'      )   ; call get_int (U_SCPR)  ! scp unit
      Case ('no_error_est')   ; U_SCPR = 29  ! just average fluxes
                                NO_ERROR_EST = .TRUE.
      Case ('scp_aves'   )    ; U_SCPR = 29  ! just average fluxes
      Case ('turn_on_scp')    ; U_SCPR = 29  ! ave flux, error est
      Case ('no_scp'      )   ; U_SCPR = 0 ; U_FLUX = 0 !b
                                NO_SCP_AVE = .TRUE.
      Case ('no_scp_ave'  )   ; U_SCPR = 0 ; U_FLUX = 0 !b
                                NO_SCP_AVE = .TRUE.
      Case ('turn_off_scp')   ; U_SCPR = 0 ; U_FLUX = 0 !b
                                NO_SCP_AVE = .TRUE.

! program control reals   
      Case ('area_thick')      ; call get_real (AREA_THICK)
                                 GLOBAL_PROPERTY = .true.
      Case ('AVE_H_WT')        ; call get_real (AVE_H_WT)
      Case ('scp_h_wt')        ; call get_real (AVE_H_WT)
      Case ('body_force')
        HAS_BODY_FORCE = .true. ; GLOBAL_PROPERTY = .true.
        IF ( N_SPACE > 4 ) STOP 'BODY FORCE > 4, Re-dimension'
        DO J = 1, N_SPACE
          call get_real (BODY_FORCE (J))
        END DO
      Case ('constant_j')      ; CONSTANT_J = .TRUE.
      Case ('convect_coef')    ; call get_real (CONVECT_COEF) 
        CONVECTION = .TRUE.    ! constant convection on mixed bc's
        IF ( CONVECT_COEF <= 0.d0 ) THEN ! bad data
          PRINT *,'WARNING, KEYWORD convect_coef VALUE <= 0'   
          N_WARN = N_WARN + 1 ; END IF
      Case ('convect_temp')    ; call get_real (CONVECT_TEMP)
      Case ('convect_thick')   ; call get_real (CONVECT_THICK)
      Case ('convect_vary')    ; CONVECT_VARY = .TRUE.
      Case ('CUT_OFF')         ; call get_real (CUT_OFF)
      Case ('tolerance')       ; call get_real (CUT_OFF)
      Case ('density')         ; call get_real (DENSITY)  !b 
      Case ('flux_thick')      ; call get_real (FLUX_THICK)
      Case ('freq_cos')        ; call get_real (FREQ_COS)
                                 FREQ_CONTROL = .TRUE.
                                 SCALE_LOADS  = .TRUE.
      Case ('freq_sin')        ; call get_real (FREQ_SIN)
                                 FREQ_CONTROL = .TRUE.
                                 SCALE_LOADS  = .TRUE.
      Case ('hughes_alpha')    ; call get_real (HUGHES_ALPHA)
      Case ('mass_damping')    ; call get_real (MASS_DAMPING)
      Case ('MAX_DERIVATIVE')  ; call get_real (MAX_DERIVATIVE)
      Case ('newmark_beta')    ; call get_real (NEWMARK_BETA)
      Case ('newmark_gamma')   ; call get_real (NEWMARK_GAMMA)
      Case ('max_error_%')     ; call get_real (PERCENT_ERR_MAX)
      Case ('PERCENT_ERR_MAX') ; call get_real (PERCENT_ERR_MAX)
      Case ('mixed_thick')     ; call get_real (ROBIN_THICK)
        IF ( CONVECT_THICK == 1.d0 ) CONVECT_THICK = ROBIN_THICK
      Case ('robin_thick')     ; call get_real (ROBIN_THICK)
      Case ('robin_square')    ; call get_real (ROBIN_1_SEG)
      Case ('ROBIN_1_SEG')     ; call get_real (ROBIN_1_SEG)
      Case ('robin_column')    ; call get_real (ROBIN_2_SEG)
      Case ('ROBIN_2_SEG')     ; call get_real (ROBIN_2_SEG)
      Case ('robin_u_free')    ; call get_real (ROBIN_1_SEG)
                                 call get_real (ROBIN_2_SEG)
      Case ('normal_flux')     ; call get_real (Q_NORMAL_SEG)
                                 FLUX_NORMAL = .TRUE.
      Case ('Q_NORMAL_SEG')    ; call get_real (Q_NORMAL_SEG)
      Case ('scalar_source')   ; call get_real (SCALAR_SOURCE)
      Case ('scp_max_error')   ; call get_real (PERCENT_ERR_MAX)
      Case ('relaxation'   )   ; call get_real (RELAXATION)
      Case ('stif_damping')    ; call get_real (STIF_DAMPING)
! test multiple inputs per line (two reals)                  !b
      Case ('test_2_r')        ; call get_real (TS_A)        !b
                                 print *, 'TS_A ', TS_A      !b
                                 call get_real (TS_B)        !b
                                 print *, 'TS_B ', TS_B      !b
! test multiple inputs per line (integer then a real)        !b
      Case ('test_i_r')        ; call get_int  (T_SETS)      !b
                                 print *, 'T_SETS ', T_SETS  !b
                                 call get_real (TS_B)        !b
                                 print *, 'TS_B ', TS_B      !b
      Case ('THICKNESS_SEG')   ; call get_real (THICKNESS_SEG)
      Case ('seg_thick'    )   ; call get_real (THICKNESS_SEG)
        IF ( CONVECT_THICK == 1.d0 ) CONVECT_THICK = THICKNESS_SEG
        IF ( FLUX_THICK    == 1.d0 ) FLUX_THICK    = THICKNESS_SEG
        IF ( ROBIN_THICK   == 1.d0 ) ROBIN_THICK   = THICKNESS_SEG
      Case ('TIME_STEP'    )   ; call get_real (TIME_STEP)
      Case ('time_step'    )   ; call get_real (TIME_STEP)
      Case ('gen_trap_beta', 'TS_B')            
        call get_real (TS_B) ! time step ratio
        TS_A = 1.d0 - TS_B   ! time step ratio
      Case ('trapezoidal')     ; call get_real (TS_B) ! time step ratio
                                 TS_A = 1.d0 - TS_B   ! time step ratio
! all other keywords
      Case default  ! all other words
        Select Case (key(1:1)) ! on first character only
          Case (' ' ) ! ignore blank line
          Case ('#' ) ! ignore comment line
          Case ('!' ) ! ignore comment line
          Case ('?' ) ! ignore comment line
          Case default  ! all other words
            print *, 'WARNING, apply_key: unknown keyword ', key
            n_warn = n_warn + 1
            count  = count  + 1
            if ( count >= limit ) then
              print *, 'STOP: reached limit on unknown words'
              print *, 'Expecting to find "end" or "quit"'
              stop 'No end or quit keyword found in control'
            end if ! likely user error
        end select ! on first character only
    end select ! from key 
  end subroutine apply_key 

SUBROUTINE OLD_F77_CONTROL_TO_KEYWORDS
! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
!      CONVERT OLD F77 FIXED CONTROL DATA TO KEYWORD FORM
!  See Appendix IV, "FE for Analysis & Design", J. E. Akin, 1994
!               STORED IN FILE old_control.key
! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
Use System_Constants
  implicit none
  integer :: unit, last
  integer :: N_CURVE ! now inactive

  interface 
    function get_next_io_unit () result (next)
      integer :: next
    end function get_next_io_unit
  end interface

!    get a free unit and open the file to save (over write)
  unit = get_next_io_unit ()
  !b print *, 'unit = ', unit !b
  open (unit, file = "old_control.key", status = "unknown")

!    read the old old controls (again)
   REWIND (5) !b
   READ  (5, *  )   TITLE
   CALL INCREMENT_INPUT_LINE ! Update position count
   READ  (5, *  ) MAX_NP, N_ELEMS, N_G_DOF, NOD_PER_EL, N_SPACE,  &
     N_F_SEG, L_B_N, N_ITER, N_CURVE, IN_RHS, I_SAY, N_R_B, N_QP, &
     L_SHAPE, N_L_TYPE, MODE
   CALL INCREMENT_INPUT_LINE ! Update position count
   READ (5, *  ) N_NP_FIX, N_NP_FLO, N_LP_FIX, N_LP_FLO, MISC_FX, &
     MISC_FL, N_HOMO, L_HOMO, NPT_WRT, LEM_WRT, N_FILE1, N_FILE2, &
     N_G_FLUX, NUL_COL, N_BS_FIX, N_BS_FLO
   CALL INCREMENT_INPUT_LINE ! Update position count

!    save the old controls and title in keyword style
   last = len_trim (title)
   write (unit, *) 'title "', title (1:last), '" ! begin keywords'

!   first old control line
   write (unit, *) 'nodes      ', MAX_NP       , &
     ' ! Number of nodes in the mesh'
   write (unit, *) 'elems      ', N_ELEMS      , &
     ' ! Number of elements in the system'
   write (unit, *) 'dof        ', N_G_DOF      , &
     ' ! Number of unknowns per node'
   write (unit, *) 'el_nodes   ', NOD_PER_EL   , &
     ' ! Maximum number of nodes per element'
   write (unit, *) 'space      ', N_SPACE      , &
     ' ! Solution space dimension'
   write (unit, *) 'b_rows     ', N_R_B        , &
     ' ! Number of rows in the B (operator) matrix'
   write (unit, *) 'shape      ', L_SHAPE      , &
     ' ! Element shape, 1=line, 2=tri, 3=quad, 4=hex'
   if ( N_F_SEG  > 0 ) write (unit, *) 'segments   ', N_F_SEG     , &
     ' ! Number of element segments with flux input'
   if ( L_B_N    > 0 ) write (unit, *) 'el_segment ', L_B_N       , &
     ' ! Maximum nodes on element boundary segmen'
   if ( N_ITER   > 1 ) write (unit, *) 'iters      ', N_ITER      , &
    ' ! Number of iterations to run'
   if ( N_ITER   < 0 ) then
                       write (unit, *) 'iters      ', IABS(N_ITER), &
                       ' ! Number of iterations to run'
                       write (unit, *) 'debug      '              , &
                       ' ! Turn on debug prints'
   end if
   !b if ( LIST_LR  > 0 ) write (unit, *) 'el_react   ' , &
   if ( LIST_LR  > 0 ) write (unit, *) 'list_el_react ' , & !b
     ' ! Compute & list element reactions' ! LIST_LR
   if ( IN_RHS   > 0 ) write (unit, *) 'loads      ' , &
     ' ! An initial source vector is input' ! IN_RHS
   if ( I_SAY    > 0 ) write (unit, *) 'remarks    ', I_SAY     , &
     ' ! Number of user remarks'
   if ( N_QP     > 0 ) write (unit, *) 'gauss      ', N_QP      , &
     ' ! Maximum number of quadrature point'
   if ( N_L_TYPE > 0 ) write (unit, *) 'el_types   ', N_L_TYPE  , &
     ' ! Number of different types of elements'
   if ( MODE == 0 ) then
     write (unit, *) 'symmetric '     , &
       ' ! Symmetric skyline storage'
   else
     write (unit, *) 'unsymmetric'    , &
       ' ! Unsymmetric skyline storage'
   end if

!   second old control line
   if ( N_NP_FIX > 0 ) write (unit, *) 'pt_int     ', N_NP_FIX , &
     ' ! Number of integer properties per node'
   if ( N_NP_FLO > 0 ) write (unit, *) 'pt_real    ', N_NP_FLO , &
     ' ! Number of real properties per node'
   if ( N_LP_FIX > 0 ) write (unit, *) 'el_int     ', N_LP_FIX , &
     ' ! Number of integer properties per element'
   if ( N_LP_FLO > 0 ) write (unit, *) 'el_real    ', N_LP_FLO , &
     ' ! Number of real properties per element'
   if ( MISC_FX  > 0 ) write (unit, *) 'integers   ', MISC_FX  , &
     ' ! Number of miscellaneous integer properties'
   if ( MISC_FL  > 0 ) write (unit, *) 'reals      ', MISC_FL  , &
     ' ! Number of miscellaneous real properties'
   if ( N_HOMO   > 0 ) write (unit, *) 'pt_homo    ', & ! N_HOMO   
     ' ! Nodal properties are homogeneous'
   if ( L_HOMO   > 0 ) write (unit, *) 'el_homo    ', & ! L_HOMO  
     ' ! Element properties are homogeneous'
   if ( NPT_WRT == 0 ) write (unit, *) 'pt_list    ', & ! NPT_WRT
     ' ! List the answers at each node point'
   if ( LEM_WRT == 0 ) write (unit, *) 'el_list    ', & ! LEM_WRT  
     ' ! List results at each node of each element'
   if ( N_FILE1  > 0 ) write (unit, *) 'post_1     ', & ! N_FILE1 
     ' ! Require post-processing, create n_file1'
   if ( N_FILE2  > 0 ) write (unit, *) 'post_2     ', & ! N_FILE2
     ' ! Require post-processing, create n_file2'
   if ( N_G_FLUX > 0 ) write (unit, *) 'seg_pt_flux', N_G_FLUX , &
     ' ! Segment flux components input at flux nodes'
   if ( NUL_COL  > 0 ) write (unit, *) 'el_no_col  ', & ! NUL_COL  
     ' ! All element column matrice are null'
   if ( N_BS_FIX > 0 ) write (unit, *) 'seg_int    ', N_BS_FIX , &
     ' ! Number of integer properties per segment'
   if ( N_BS_FLO > 0 ) write (unit, *) 'seg_real   ', N_BS_FLO , &
     ' ! Number of real properties per segment'

!    terminate keyword control input
   write (unit, *) 'quit ! keyword input'
   close (unit)
END SUBROUTINE OLD_F77_CONTROL_TO_KEYWORDS

! end: keyword_control_lib.f
