!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  String search

! Section 6.1, Exercise 10
  program Count_Me_In
    implicit none
    integer, parameter :: S_Len = 500, P_Len = 20
    character (len = S_Len) :: String
    character (len = P_Len) :: Pattern
    integer :: LS, LP, Me_In, Count, Loc
! start program Count_Me_In
    open (1, "str_500.txt")
    read (1, *) String
    LS = Len_Trim (String)
    print *, String(: LS)
    print *, " Number of characters: ", LS
    do
      print *, " How long is the search pattern? "
      read *, LP
      if (LP == 0) exit
      LP = Min (LP, S_Len)
      print *, " Please enter the search pattern: "
      read *, Pattern
      Me_In = 1 - LP
      Count = 0
      print *, ' Search for: "', Pattern(: LP), '"'
      do
        Loc = Index (String(Me_In + LP: LS), Pattern(: LP))
        if (Loc == 0) exit
        Count = Count + 1
        Me_In = Me_In + LP + Loc - 1
        print *, Count, Me_In, '"', String(Me_In: Me_In + LP - 1), '"'
      end do
      print *, " Number of occurrences: ", Count
    end do
    stop " OK. "
  end program Count_Me_In

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  "Old fashioned" (non-fast) Fourier transform -- with Complex exponential
!  See Hildebrand, "Intro to Numerical Analysis", McGraw Hill 1956

  program Periodic
    implicit none
    integer, parameter :: N = 6, L = 5
    real, parameter :: Pi = 3.1415926536
    real, dimension(-N: N) :: Data_In
    complex, dimension(0: N) :: Coeffs
    real, dimension(0: N) :: Pt_Sum
    real :: Theta
    integer :: K, R, NL
! start program Periodic
    open (unit = 1, file = "prob4.txt", position = "REWIND")
    read (1, *) (Theta, Data_In(R), R = -N, N)

    do R = -N, N
      print *, R, Real(R) * Pi / Real(N), Data_In(R)
    end do

    print *, " COEFFICIENTS: "

    Coeffs(0) = T_Sum(0) / Real(2 * N)
    print "(tr1, 2 f10.3)", Coeffs(0)

    do K = 1, N - 1
      Coeffs(K) = T_Sum(K) / Real(N)
      print "(tr1, 2 f10.3)", Coeffs(K)
    end do

    Coeffs(N) = T_Sum(N) / Real(2 * N)
    print "(tr1, 2 f10.3)", Coeffs(N)

    NL = L * N           ! L is number of output lines per input interval

    print *, " TIME PLOT: "

    do R = -NL,  NL
      Theta = Real(R) * Pi / Real(NL)
      Pt_Sum(0) = Real (Coeffs(0))
      do K = 1, N       ! Compute all partial sums, but don't print them all.
          ! Funny way to get A*C + B*D where (A, B) and (C, D) are complex:
        Pt_Sum(K) = Pt_Sum(K - 1) + Real (Coeffs(K) * Exp (Cmplx (0.0, - (K * Theta))) )
      end do
      if (Mod (R, L) == 0) print  *
      print "(tr1, 10 f10.4)",  Theta, (Pt_Sum(K), K = N / 2, N)
      if (Mod (R, L) == 0) print *
    end do

    stop " OK! "

contains

    function T_Sum (K) result (Sum)
      integer, intent (in) :: K
      complex :: Sum
      integer :: R
 !  start function T_Sum
      ! Data_In is periodic so -N and N are equal by def. So average them.
      Sum =  0.5 * ( Data_In(-N) * Exp (Cmplx (0.0, Real(- K) * Pi ))                &
                   + Data_In( N) * Exp (Cmplx (0.0, Real(  K) * Pi )) )
      do R = -N + 1, N - 1
        Sum = Sum + Data_In(R) * Exp (Cmplx (0.0, Real(K * R) * Pi / Real (N)))
      end do
      return
    end function T_Sum

end program Periodic

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Determine which Integer and Real KINDS are supported on this processor.
!    (Check for ranges and precisions up to 100, exponents up to 1000.)

program Check_Kinds
   implicit none
   integer, parameter :: MIR = 100, MRP = 100, MRE = 1000
   integer :: Range, Precision, KInt(0: MIR) = 0, KReal(2, 0: MRE) = 0

!start program Check_Kinds

   print *, " Integer Ranges and Kinds "
   KInt(MIR) = 0
   do Range = 1, MIR
      KInt(Range) = Selected_Int_Kind (R = Range)
      if (KInt(Range) /= KInt(Range - 1)) print *, " R: ",Range, " K: ", KInt(Range)
   end do

   print *, " Real Precisions, Ranges, and Kinds "
   do Precision = 1, MRP
      Kreal(1, :) = KReal(2, :)
      do Range = 1, MRE
         KReal (2, Range) = Selected_Real_Kind (P= Precision, R = Range)
         if ((KReal(2, Range) /= KReal(1, Range)) &
            .and. (KReal(2, Range) /= KReal(2, Range - 1))) &
               print *, " P: ", Precision, " R: ", Range, " K: ", KReal(2, Range)
      end do
   end do

   close (unit = 1)
   stop " OK! "

end program Check_Kinds

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 8.10. Binary search tree.
  program Binary_Tree
    implicit none
    integer, parameter :: S_LEN = 20
    type :: Box_Type
      character (len = S_LEN) :: Key
      type (Box_Type), pointer :: Left, Right
    end type Box_Type
    type (Box_Type), pointer :: Root
    character (len = S_LEN) :: New_Key = " "
    integer :: EoF
  ! start program Binary_Tree
    nullify (Root)
    open (unit = 1, file = "simply.txt")
  Input_Loop: &
    do
      read (1, *, iostat = EoF) New_Key
      if (EoF < 0) exit Input_Loop
      call Look_Up (Root)
    end do Input_Loop
  ! This point is reached at end of input file.
    call Print_Tree ( )
    stop " OK. "

  contains

    recursive subroutine Look_Up (Current_Ptr)
      type (Box_Type), pointer :: Current_Ptr, New_Item
  ! start subroutine Look_Up
      if (Associated (Current_Ptr)) then  ! Pointer is not NIL.
        if (New_Key == Current_Ptr % Key) then
          print *, New_Key, " Found. "
        else if (New_Key > Current_Ptr % Key) then
          call Look_Up (Current_Ptr % Right)  ! Recursion
        else  ! New_Key < Current_Ptr % Key
          call Look_Up (Current_Ptr % Left)  ! Recursion
        end if
      else  ! Insertion occurs only when Current_Ptr is NIL.
        print *, New_Key, " Not found. "
        allocate (New_Item)
        New_Item % Key = New_Key
        nullify (New_Item % Left)
        nullify (New_Item % Right)
        Current_Ptr => New_Item  ! Assign to argument pointer.
      end if
      return
    end subroutine Look_Up

    subroutine Print_Tree ( )
  ! start subroutine Print_Tree
      print *, " ====BINARY SEARCH TREE==== "  ! Print heading.
      call Print_Subtree (Root)
      return
    end subroutine Print_Tree

    recursive subroutine Print_Subtree (Current_Ptr)
      type (Box_Type), pointer :: Current_Ptr
  ! start subroutine Print_Subtree
      if (Associated (Current_Ptr)) then
        call Print_Subtree (Current_Ptr % Left)
        print *, Current_Ptr % Key
        call Print_Subtree (Current_Ptr % Right)
      end if
      return
    end subroutine Print_Subtree

  end program Binary_Tree

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 7.16. Generic vector maximum and minimum.
  module Max_Loc_Mod

  ! Specification part of module
    implicit none
    integer :: R(1)
    interface Max_Loc
      module procedure I_Max_Loc, R_Max_Loc
    end interface

  ! Subprogram part of module
  contains

    function I_Max_Loc (Vector, Mask) result (Ans)
      integer, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function I_Max_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MaxLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MaxLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function I_Max_Loc

    function R_Max_Loc (Vector, Mask) result (Ans)
      real, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function R_Max_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MaxLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MaxLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function R_Max_Loc

  end module Max_Loc_Mod

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
! Example 7.12. Selection sort.
! Uses Generic Max-Min Module
  module Sel_Sort_Mod
    use Max_Loc_Mod  ! Specification part of module
    implicit none

  contains  ! Subprogram part of module
    recursive subroutine Selection_Sort (Unsorted)  ! Module subprogram
      real, intent (in out) :: Unsorted (:)  ! Assumed-shape array
      integer :: N, Loc
  ! start subroutine Selection_Sort
      N = Size (Unsorted)
      if (N > 1) then  ! Skip if only one element remains.
        Loc = Max_Loc (Unsorted)  ! Location of largest element
        call Swap_It (N, Loc)
        call Selection_Sort (Unsorted(1: N - 1))  ! Recursion
      end if
      return

    contains  ! Internal procedure in Selection_Sort
      subroutine Swap_It (I, J)  ! Internal subprogram
        ! Obtains Unsorted by inheritance from host subprogram.
        integer, intent (in) :: I, J
        real :: Aux
    ! start subroutine Swap_It
        Aux = Unsorted(I)
        Unsorted(I) = Unsorted(J)
        Unsorted(J) = Aux
        return
      end subroutine Swap_It

    end subroutine Selection_Sort

  end module Sel_Sort_Mod

  program Variety
    use Sel_Sort_Mod
    implicit none
    integer, parameter :: MANY = 12
    real :: Vector (MANY)
! start program Variety
    open (1, file = "file_in.txt", position = "REWIND")
    read (1, *) Vector
    call Selection_Sort (Vector)
    print *, Vector
    stop
  end program Variety

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 7.16. Generic vector maximum and minimum.
! With main program to call MinVal with character arg
  module Max_Min_Loc_Mod

  ! Specification part of module
    implicit none
    integer :: R(1)
    interface Max_Loc
      module procedure I_Max_Loc, R_Max_Loc, S_Max_Loc
    end interface
    interface MaxVal  ! Same as intrinsic function name
      module procedure S_Max_Val
    end interface
    interface Min_Loc
      module procedure I_Min_Loc, R_Min_Loc, S_Min_Loc
    end interface
    interface MinVal  ! Same as intrinsic function name
      module procedure S_Min_Val
    end interface

  ! Subprogram part of module
  contains

    function I_Max_Loc (Vector, Mask) result (Ans)
      integer, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function I_Max_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MaxLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MaxLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function I_Max_Loc

    function I_Min_Loc (Vector, Mask) result (Ans)
      integer, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function I_Min_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MinLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MinLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function I_Min_Loc

    function R_Max_Loc (Vector, Mask) result (Ans)
      real, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function R_Max_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MaxLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MaxLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function R_Max_Loc

    function R_Min_Loc (Vector, Mask) result (Ans)
      real, intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      integer :: Ans
  ! start function R_Min_Loc
      if (Size (Vector) > 0) then
        if (Present (Mask)) then
          if (All (Shape (Vector) == Shape (Mask)) .and. Any (Mask)) then
            R = MinLoc (Vector, Mask = Mask)
          else
            R = 0
          end if
        else
          R = MinLoc (Vector)
        end if
      else
        R = 0
      end if
      Ans = R(1)
      return
    end function R_Min_Loc

    function S_Max_Loc (Vector, Mask) result (Ans)
      character (len = *), intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      character (len = Len(Vector)) :: Max_String
      logical, dimension(Size(Vector)) :: M
      integer :: Ans
      integer :: Loop
  ! start function S_Max_Loc
      if (Present (Mask)) then
        M = Mask
      else
        M = .TRUE.
      end if
      Max_String = " "
      Ans = 0
      do Loop = 1, Size (Vector)
        if (M(Loop) .and. (Vector(Loop) > Max_String)) then
          Max_String = Vector(Loop)
          Ans = Loop
        end if
      end do
      return
    end function S_Max_Loc

    function S_Max_Val (Vector, Mask) result (Ans)
      character (len = *), intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7
      character (len = Len(Vector)) :: Ans
      integer :: Loc
  ! start function S_Max_Val
      if (Present (Mask)) then
        Loc = S_Max_Loc (Vector, Mask)
      else
        Loc = S_Max_Loc (Vector)
      end if
      if (Loc > 0) then
        Ans = Vector(Loc)
      else
        Ans = " "
      end if
      return
    end function S_Max_Val

    function S_Min_Loc (Vector, Mask) result (Ans)
      character (len = *), intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7.
      character (len = Len(Vector)) :: Min_String
      logical, dimension(Size(Vector)) :: M
      integer :: Ans
      integer :: Loop
  ! start function S_Min_Loc
      if (Present (Mask)) then
        M = Mask
      else
        M = .TRUE.
      end if
      Min_String = " "
      Ans = 0
      do Loop = 1, Size (Vector)
        if (M(Loop)) then
        Min_String = Vector(Loop)
        Ans = Loop
        end if
        end do
      do Loop = 1, Size (Vector)
        if (M(Loop) .and. (Vector(Loop) < Min_String)) then
          Min_String = Vector(Loop)
          Ans = Loop
        end if
      end do
      return
    end function S_Min_Loc

    function S_Min_Val (Vector, Mask) result (Ans)
      character (len = *), intent (in), dimension(:) :: Vector
      logical, intent (in), optional, dimension(:) :: Mask  ! See Section 8.7
      character (len = Len(Vector)) :: Ans
      integer :: Loc
  ! start function S_Min_Val
      if (Present (Mask)) then
        Loc = S_Min_Loc (Vector, Mask)
      else
        Loc = S_Min_Loc (Vector)
      end if
      if (Loc > 0) then
        Ans = Vector(Loc)
      else
        Ans = " "
      end if
      return
    end function S_Min_Val

  end module Max_Min_Loc_Mod

  program Variety
    use Max_Min_Loc_Mod
    implicit none
    integer, parameter :: MANY = 12
    character (len = 5), dimension(MANY) :: Vector
    character (len = 5) :: Small
    integer :: Loop
! start program Variety
    open (1, file = "char_in.txt", position = "REWIND")
    read (1, *) Vector
    do Loop = 1, MANY
      print *, Vector(Loop)
      end do
    Small = MinVal(Vector)
    print *, ' Smallest string is: "', Small, '"'
    stop " Vector Minimum found. "
  end program Variety

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 7.2. Complex numbers in Fortran style and in polar form.
  module Ursa

    implicit none
    type :: Polar
      real :: R
      real :: Theta
    end type Polar

  contains

    function P_Complex (Z) result (S)
      complex, intent (in) :: Z
      type (Polar) :: S
  ! start function P_Complex
      S % R = Abs (Z)
      S % Theta = ATan2 (AImag (Z), Real (Z))
      return
    end function P_Complex

    function F_Complex (P) result (Bear)
      type (Polar), intent (in) :: P
      complex :: Bear
      real :: X, Y
  ! start function F_Complex
      X = P % R * Cos (P % Theta)
      Y = P % R * Sin (P % Theta)
      Bear = Cmplx (X, Y)
      return
    end function F_Complex

    function Mult_Polar (P1, P2) result (P)
      type (Polar), intent (in) :: P1, P2
      type (Polar) :: P
  ! start function Mult_Polar
      P % R = P1 % R * P2 % R
      P % Theta = P1 % Theta + P2 % Theta
      return
    end function Mult_Polar

  end module Ursa

  program Driver
    use Ursa  ! Import type and operations from module.
    implicit none
    complex :: Z1, Z2, Z3
    type (Polar) :: P1, P2, P3  ! Derived type imported from module
! start program Driver
    open (1, file = "h02.dat", position = "rewind")
    read (1, *) Z1, Z2
    print *, " Complex numbers in Fortran style: "
    print *, Z1, Z2
    P1 = P_Complex (Z1)
    P2 = P_Complex (Z2)
    print *, " Same complex numbers in polar form: "
    print *, P1, P2
    P3 = Mult_Polar (P1, P2)
    print *, " Complex product in polar form: "
    print *, P3
    Z3 = F_Complex (P3)
    print *, " Complex product in Fortran style: "
    print *, Z3
    stop
  end program Driver

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 6.17. Convert uppercase letters to lowercase letters.
  program Decapitation
    implicit none
    integer, parameter :: STRING_LENGTH = 80
    character (len = STRING_LENGTH) :: Text_1, Text_2
! start program Decapitation
    open (1, file = "g17.dat", position = "rewind")
    read (1, *) Text_1
    print *, Text_1
    Text_2 = Decap (Text_1)
    print *, Text_2
    stop
  contains
    function Decap (String) result (New_String)
      character (len = *), intent(in) :: String  ! Assumed length (see Section 8.3)
      character (len = Len(String)) :: New_String ! Automatic length result
      integer, parameter :: ALPHABET = 26
      character (len = ALPHABET) :: UCAlph = "ABCDEFGHIJKLMNOPQRSTUVWXYZ", &
          LCAlph = "abcdefghijklmnopqrstuvwxyz"
      integer :: Next, S
  ! start function Decap
      New_String = String
      do S = 1, Len (String)
        Next = Index (UCAlph, String(S: S))
        if (Next /= 0) New_String(S: S) = LCAlph (Next: Next)
      end do
      return
    end function Decap
  end program Decapitation

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 5.45. Ising transition model.
  program Window
    implicit none
    logical, allocatable, dimension(:, :, :) :: Lodge
    integer, allocatable, dimension(:, :, :) :: Ising, Counter
    real, allocatable, dimension(:, :, :) :: P, Ram
    real :: Changes, RL3
    integer :: L, Time, More, T, I
! start program Window
    open (unit = 1, file = "ising.dat", position = "REWIND")
    read (1, *) L
    allocate (Lodge(L, L, L), Ising(L, L, L), Counter(L, L, L), P(L, L, L), Ram(L, L, L))
    P = 1.0
    read (1, *) Lodge
    Ising = Merge (1, 0, Mask = Lodge)  ! Merge elementwise according to Mask.
    print '(" ", I6, 64 A1)', L**3, Merge ("+", "-", Ising(: 4, : 4, : 4)  == 1)
    RL3 = L**3
    T = 1
    do
      Time = T - 1
      print *, " How many steps? [", Time, "] .. (enter 0 to stop.)"
      read *, More
      if (More <= 0) exit
      Changes = 0.0
      do T = Time + 1, Time + More
        call Transition ( )
      end do
      print '(" ", F6.4, 64 A1)', Changes / More, Merge ("+", "-", Ising(: 4, : 4, : 4)  == 1)
    end do
    stop

  contains

    subroutine Transition ( )
      real, parameter :: PROB(4: 6) = (/ 0.2, 0.1, 0.0 /)
  ! start subroutine Transition
      Counter = CShift (Ising, Shift = -1, Dim = 1) + CShift (Ising, Shift = 1, Dim = 1) &
        + CShift (Ising, Shift = -1, Dim = 2) + CShift (Ising, Shift = 1, Dim = 2) &
        + CShift (Ising, Shift = -1, Dim = 3) + CShift (Ising, Shift = 1, Dim = 3)   ! Circular shift
      where (Ising == 0) Counter = 6 - Counter  ! Assignment with "where"
      do I = 4, 6
        where (Counter == I) P = PROB(I)  ! Assignment with "where"
      end do
      call Random_Number (Harvest = Ram)
      where (Ram < P) Ising = 1 - Ising  ! Assignment with "where"
      Changes = Changes + Count (Ram < P) / RL3
      return
    end subroutine Transition

  end program Window

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 5.39. Compute exponential of a matrix.
  program Matress
    implicit none
    integer, parameter :: P = 4, K = Selected_Real_Kind (P)
    integer :: Dim
    real (kind = K), allocatable, dimension(:, :) :: Array, X
! start program Matress
    open (1, file = "f39.dat", position = "rewind")
    read (1, *) Dim
    allocate (Array(Dim, Dim), X(Dim, Dim))
    read (1, *) Array
    X = Mat_Exp (Array)
    print *, X
    stop
  contains
    function Mat_Exp (M) result (E)
      real (kind = K), dimension(:, :), intent (in) :: M
      real (kind = K), dimension(Size(M, Dim = 1), Size(M, Dim = 2)) :: E
      integer :: I, N, M_Shape(2)
      real (kind = K) :: Sign = 1.0
      real (kind = K), allocatable, dimension(:, :) :: Term
    ! start function Mat_Exp
      E = 0.0_K      ! Initialize power series.
      M_Shape = Shape(M)
      if (Sum (M ** 2) > 0.5_K) return  ! M will probably not converge fast enough.
      if (M_Shape(1) /= M_Shape(2)) return ! M must be a square matrix
      allocate (Term(M_Shape(1), M_Shape(2)))
      Term = 0.0_K  ! Initialize Term to identity matrix.
      do I = 1, Dim
        Term(I, I) = 1.0_K
      end do
      do N = 0, 100  ! Exit when Term is small or after Term is 100.
        E = E + Term  ! Exp (M) power series
        if (Sum(Term ** 2) < 2.0_K * Spacing (1.0_K)) return  ! Converged OK.
        Term = Sign * MatMul (Term, M) / Real (N + 1, Kind = K)  ! Calculate next Term.
        Sign = - Sign
      end do
      return
    end function Mat_Exp
  end program Matress

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 6.29. Expanding list as an array of structures.
  program Array_List
    implicit none
    integer, parameter ::  NAME_LENGTH = 24, TBL_SIZE = 50
    type :: Compound
      integer :: Count
      character (len = NAME_LENGTH) :: Name
    end type Compound
    type :: Table  ! Structure to hold list and current list size
      type (Compound) :: Data (TBL_SIZE)  ! Array of compounds
      integer :: L_Size  ! Current list size
    end type Table
    type (Compound) :: Element = Compound (1, " ")  ! Structure constructor
    type (Table) :: C_Table
    integer :: EoF = 0
! start program Array_List
    C_Table % L_Size = 0
    open (unit = 1, file = "compound.txt")
    do
      read (1, *, iostat = EoF) Element % Name
        if (EoF < 0) exit
      if (.not. Find (C_Table, TBL_SIZE, Element)) exit
    end do
  ! At this point EoF < 0 or attempted to insert in full list.
    call Print_List ( )
    stop

  contains

    function Find (C_Tbl, T_Size, Element) result (Fit)
  ! Result is .FALSE. on attempt to insert when Tbl is full.
      integer, intent (in):: T_Size  ! Maximum list size
      type (Table), intent (in out) :: C_Tbl  ! List
      type (Compound), intent (in) :: Element  ! Element to find or insert in list
      logical :: Fit
      integer, dimension(2) :: Loc  ! To hold result from Search
!  start function Find
      Loc = Search &
        (C_Tbl % Data(: C_Tbl % L_Size) % Name, Element % Name)
      if (Loc(1) == Loc(2)) then
        print *, " Found at  ", Loc(2), " of ", C_Tbl % L_Size, ": ", &
          Element % Name, C_Tbl % Data(Loc(2)) % Count
        C_Tbl % Data(Loc(2)) % Count = C_Tbl % Data(Loc(2)) % Count + 1
        Fit = .TRUE.  ! Current list size has not changed.
      else  ! No matching element; insert.
        if (C_Tbl % L_Size < T_Size) then
          Fit = .TRUE.  ! There is room to insert element.
          C_Tbl % Data(Loc(2) + 1 : C_Tbl % L_Size + 1) = C_Tbl % Data(Loc(2) : C_Tbl % L_Size)
  ! Move elements up to make space at insertion point.
          C_Tbl % Data(Loc(2)) = Element  ! Insert.
  ! Adjust current size.
          C_Tbl % L_Size = C_Tbl % L_Size + 1
          print *, " Insert at ", Loc(2), " of ", C_Tbl % L_Size, ": ", Element % Name
        else
          Fit = .FALSE.  ! There is no room to insert element.
          print *, " List full. ", Loc(2), " of ", C_Tbl % L_Size, ": ", Element % Name
        end if
      end if
      return
    end function Find

    subroutine Print_List ( )
      integer :: Loop
  ! start subroutine Print_List
      print *, " Final size: ", C_Table % L_Size
      do Loop = 1, C_Table % L_Size
        print *, C_Table % Data(Loop) % Count, " ", C_Table % Data(Loop) % Name
      end do
      return
    end subroutine Print_List

  ! Example 5.36. Binary search in ordered array.
    function Search (List, Key) result (Location)
      integer :: Location(2), Middle, Last
      character (len = *), intent(in) :: List(:), Key
  ! start function Search
      Last = Size(List)
      Location(1) = 0
      Location(2) = Last + 1
      do while (Location(2) - Location(1) > 1)
        Middle = (Location(1) + Location(2)) / 2
        if (List(Middle) == Key) then
          Location = Middle  ! Assign scalar to both vector elements.
          exit
        else if (List(Middle) > Key) then
          Location(2) = Middle
        else
          Location(1) = Middle
        end if
      end do
      return
    end function Search

  end program Array_List

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 8.15. Optional arguments and the intrinsic function Present.
  program Discretion
    implicit none
! start program Discretion
  ! All of the following function references are equivalent.
    print *, Series_Sum (0, 700, 100.0)  ! D assumes default value.
    print *, Series_Sum (0, 700, S = 100.0)  ! D assumes default value.
    print *, Series_Sum (N = 700, S = 100.0)  ! M and D assume default values.
    print *, Series_Sum (D = 0.1, S = 100.0, N = 700)  ! M assumes default value.
    stop " OK. "
  contains
    real function Series_Sum (M, N, S, D) result (Series_Sum_Result)
      integer, intent (in), optional :: M
      integer, intent (in) :: N
      real, intent (in) :: S
      real, intent (in), optional :: D
      real :: Series_Sum_Result
      integer :: Temp_M, Loop
      real :: Temp_D
  ! start function Series_Sum
      if (Present (M)) then
        Temp_M = M
      else
        Temp_M = 0
      end if
      if (Present (D)) then
        Temp_D = D
        else
          Temp_D = 0.1
        end if
      Series_Sum_Result = 0.0
      do Loop = Temp_M, N
        Series_Sum_Result = Series_Sum_Result + S + Temp_D * real (Loop)
      end do
      return
    end function Series_Sum
  end program Discretion

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Earlier versions of ELF failed because INTENT was not declared for pointers.

  program Pointer_Arg_Intent
    implicit none
    integer, parameter :: Str_Len = 15
    character (len = Str_Len), pointer :: Root

! start program Pointer_Arg_Intent
    nullify (Root)
    call Create (Root)
    print *, Get_Target (Root)
    stop

  contains

    subroutine Create (Temp_Ptr)
      character (len = Str_Len), pointer :: Temp_Ptr   ! Dummy Arg is a pointer; can't declare intent
      character (len = Str_Len), pointer :: New_Item   ! Local
  ! start subroutine Create
      allocate (New_Item)
      New_Item = " A NEW STRING "
      Temp_Ptr => New_Item  ! Assign to argument pointer.
      return
    end subroutine Create

    function Get_Target (Temp_Ptr) result (String)
      character (len = Str_Len), pointer :: Temp_Ptr   ! Dummy Arg is a pointer; can't declare intent
      character (len = Str_Len) :: String              ! Function result

  ! start subroutine Get_Target
      String = Temp_Ptr                                ! Copy target
      return
    end function Get_Target

  end program Pointer_Arg_Intent

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! Example 4.30. Statistics with Linear Regression, by updating method.
  program Statistics_4
    implicit none
    real :: Count = 0.0, X, Y, X_Sum = 0.0, Y_Sum = 0.0, &
      XX_Res = 0.0, YY_Res = 0.0, XY_Res = 0.0, &
      X_Dev, Y_Dev, Covar, Correlation, A, B
    integer :: EoF
! start program Statistics_4
    open (2, file = "numeroso.txt")
    do
      read (2, *, iostat = EoF) X, Y
      if (EoF < 0) exit
      Count = Count + 1.0
      X_Sum = X_Sum + X
      Y_Sum = Y_Sum + Y
      if (Count > 1.0) then
        XX_Res = XX_Res + (Count * X - X_Sum) ** 2 / (Count * (Count - 1.0))
        YY_Res = YY_Res + (Count * Y - Y_Sum) ** 2 / (Count * (Count - 1.0))
        XY_Res = XY_Res + (Count * X - X_Sum) * (Count * Y - Y_Sum) &
          / (Count * (Count - 1.0))
      end if
    end do
    if (Count > 1.0) then
      X_Dev = SqRt (XX_Res / Count)
      Y_Dev = SqRt (YY_Res / Count)
      Covar = XY_Res / Count
      Correlation = Covar / (X_Dev * Y_Dev)
      print *, X_Dev, Y_Dev, Covar, Correlation
      B = XY_Res/XX_Res
      A = (Y_Sum - B * X_Sum) / Count
      print *, " Regression: ", A, B
    else
      print *, " Insufficient data in this file. "
    end if
    stop
  end program Statistics_4

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Gauss elimination with vector subscripts
!  Example 5.17a
  module Gauss_Elimination

  ! specification part
    implicit none
    integer, parameter :: LOW = Selected_Real_Kind (6),   HIGH = Selected_Real_Kind (12)
    real (kind = LOW) :: R_Norm
    real (kind = LOW), allocatable, private :: Save_A(:, :), LU(:, :), R(:), W(:)
    integer, private :: N
    integer, private, allocatable :: P(:)
    logical, private :: Saved = .FALSE., Exist = .FALSE.
    private :: Reduce   ! Private procedure in module

! procedure part
  contains

    subroutine Factor (A, B, X, E, X_Norm, E_Norm, Flag)   ! Module procedure (public)
      ! This subroutine solves the system of simultaneous linear algebraic equations
      !   A x = b, where A is an n by n matrix and x and b are vectors of length n.
      !   A and B are assumed-shape arrays supplied by the calling program.
      !   X, E, X_Norm, E_Norm, and Flag are computed and returned by Factor.
      !   E is the error vector computed as A e = r, where r = B - A x(calculated).
      !   Factor does not change A or B.
      ! On return from Factor, X contains the solution vector and E contains the error
      !   vector.  X_Norm and E_Norm contain the (Euclidean vector) norms of X and E,
      !   respectively.  However, if Flag is false, the contents of X, E, X_Norm, and
      !   E_Norm are unpredictable, since the algorithm may have short-circuited
      !   somewhere in the middle.
      ! The problem size N, a copy of A, and the LU decomposition matrix (with
      !   permutation vector P) are saved in the module as private arrays, and are not
      !   changed except when Factor is called.
      ! The subroutine Solve may be called to solve additional systems with the same
      !   matrix A, between successful (i.e., Flag = true) calls to Factor.
      !   Call Solve with a new right-hand side B whose length equals the original N.
      !   Solve will compute the solution vector X, the corresponding error vector E,
      !   and the norms.
      ! DUMMY ARGUMENTS
      real (kind = LOW), intent (in) :: A(:, :), B(:)
      real (kind = LOW), intent (out) :: X(:), E(:), X_Norm, E_Norm
      logical, intent(out) :: Flag
      ! LOCAL DATA
      integer :: I, J, M, I_Temp
      real (kind = LOW) :: Temp
  ! start subroutine Factor

      ! Set Flag to false for quick error exit
      Flag = .FALSE.

      ! Determine problem size; allocate and initialize private arrays
      N = Size (A, dim = 1)
      if ((Size (A, dim = 2) /= N) .or. (Size (B) /= N) .or. &
        (Size (X) /= N) .or. (Size (E) /= N))  stop " Size incompatibility."   ! Error
      allocate (Save_A(N, N), LU(N, N), R(N), W(N), P(N))
      print *, " Arrays Save_A, LU, R, W, P have been allocated: ", N
      Exist = .TRUE.
      Save_A = A
      LU = A
      P = (/ (I, I = 1, N) /)
      ! Store row norms in W
      W = (/ (MaxVal(Abs (LU(I, :))), I = 1, N) /)

      if (Any (W == 0)) then                                                   ! Error
        print *, " Row norms: ", W
        stop " One of the row norms is zero. "
      end if

      ! Perform factorization A = L * U with scaled partial pivoting.
      do M = 1, N
        ! Reduce column M and choose Pivot
        Temp = 0.0
        I_Temp = 0

        do I = M, N
          ! LU(P(...), ...) has a Vector Subscript
          LU(P(I), M) = Reduce (LU(P(I), M), LU(P(I), 1: M - 1), LU(P(1: M - 1), M))
          if (Abs (LU(P(I), M) / W(P(I))) > Temp) then
            Temp = Abs (LU(P(I), M) / W(P(I)))
            I_Temp = I
          end if
        end do

        if (Temp <= 0.0) stop " All pivot candidates are zero. "               ! Error
        call Permute (M, I_Temp)

        ! Reduce row M
        do J = M + 1, N
          ! LU(P(...), ...) has a Vector Subscript
          LU(P(M), J) = Reduce (LU(P(M), J), LU(P(M), 1 : M - 1), LU(P(1: M - 1), J)) / LU(P(M), M)
        end do

      end do

      print *, " LU Factors"
      do I = 1, N
        print *, LU(I, 1: I - 1), "|", LU(I, I: N)
      end do

      ! Apply LU with P to right-hand side B; compute residual (private) and error.
      Saved = .TRUE.                                  ! Module has stored the Factors.
      call Solve (B, X, E, X_Norm, E_Norm)
      if (E_Norm < X_Norm / 2.0) Flag =.TRUE.
      Saved = Flag                ! Module has stored the Factors ready to Solve again
      return

    contains

      subroutine Permute (I, J)       ! Internal procedure in Factor
        integer, intent(in) :: I, J
        integer :: Temp
    ! start subroutine Permute
        Temp = P(I)
        P(I) = P(J)
        P(J) = Temp
        return
      end subroutine Permute

    end subroutine Factor

    subroutine Unsolve ( )
  ! start subroutine Unsolve
      Saved = .FALSE.
      if (Exist) then
        print *, " Time to deallocate Save_A, LU, R, W, P: ", N
        deallocate (Save_A, LU, R, W, P)
      end if
      Exist = .FALSE.
      return
    end subroutine Unsolve

    function Reduce (A, Row, Col) result (Sum)      ! Module procedure (private)
      ! DUMMY ARGUMENTS
      real (kind = LOW), intent (in) :: A, Row(:), Col(:)
      real (kind = LOW) :: Sum
  ! start function Reduce
      ! Change A, Row, and Col to HIGH precision;
      ! Compute A - Row times Col;
      ! Change result to LOW precision.
      Sum = Real (Real (A, Kind = HIGH) - Dot_Product (Real (Row, Kind = HIGH), Real (Col, Kind = HIGH)), Kind = LOW)
      return
    end function Reduce

    subroutine Solve (B, X, E, X_Norm, E_Norm)          ! Module procedure (private)
      ! DUMMY ARGUMENTS
      real (kind = LOW), intent (in) :: B(:)
      real (kind = LOW), intent (out) :: X(:), X_Norm, E(:), E_Norm
  ! start subroutine Solve
      if ((Size (B) /= N) .or. (Size (X) /= N) .or. (Size (E) /= N) .or. .not. Saved) then
        print *, " B, X, or E is the wrong size. "
        print *, Size(B), Size(X), Size(E), " ", Saved
        X = 0.0_LOW
        E = 0.0_LOW
        X_Norm = 0.0_LOW
        E_Norm = Huge (1.0_LOW)
        return
      end if
      call For_Bak (B, X, X_Norm)     ! Destroys W
      call Residual (B, R, R_Norm)
        ! The two following print statements are optional and can be deleted:
        print *, " Residual Norm and Residual Vector: ", R_Norm
        print *, R
      call For_Bak (R, E, E_Norm)     ! Destroys W
      return

    contains

      subroutine For_Bak (B, X, X_Norm)            ! Internal procedure in Solve
        ! This subroutine performs forward and backward substitution,
        !   using LU and P stored at the latest call to Factor.
        ! B is a right-hand side vector.
        ! The result will be returned in X, and its Euclidean norm in X_Norm.
        ! DUMMY ARGUMENTS
        real (kind = LOW), intent (in) :: B(:)
        real (kind = LOW), intent (out) :: X(:), X_Norm
        ! LOCAL DATA
        integer :: K
        real (kind = LOW) :: SXX
    ! start subroutine For_Bak
        SXX = 0.0_LOW
        do K = 1, N
          W(K) = Reduce (B(P(K)), LU(P(K), 1: K - 1), W(1: K - 1)) / LU(P(K), K)
        end do
        do K = N, 1, -1
          X(K) = Reduce (W(K), LU(P(K), K + 1: N), X(K + 1: N))
          SXX = SXX + X(K) ** 2
        end do
        X_Norm = SqRt (SXX)
        return
      end subroutine For_Bak

      subroutine Residual (B, R, R_Norm)               ! Internal procedure in Solve
        ! Uses Save_A from latest call to Factor
        ! DUMMY ARGUMENTS
        real (kind = LOW), intent (in) :: B(:)
        real (kind = LOW), intent (out) :: R(:), R_Norm
        ! LOCAL DATA
        integer :: I
        real (kind = LOW) :: SRR
   !  start subroutine Residual
        SRR = 0.0_LOW
        do I = 1, N
          R(I) = Reduce (B(I), Save_A(I, :), X)
          SRR = SRR + R(I) ** 2
        end do
        R_Norm = SqRt (SRR)
        return
      end subroutine Residual

    end subroutine Solve

  end module Gauss_Elimination

  program Driver
    use Gauss_Elimination
    implicit none
    integer :: EoF, N
    real (kind = LOW), allocatable :: A(:, :), B(:), X(:), E(:)
    real (kind = LOW) :: X_Norm, E_Norm
    integer :: I
    logical :: Flag
! start program Driver
    open (11, file = "GauData.txt", position = "REWIND")
    do
      read (11, *, iostat = EoF) N
      if (EoF <  0) stop " End of data file reached.  "
      print *
      print *, " Allocating A, B, X, E: ", N
      print *
      allocate (A(N, N), B(N), X(N), E(N))

      do I = 1, N       ! Read input matrix A by rows and echo print it.
        read (11, *) A(I, :)
        print *, A(I, :)
      end do

      do I = 1, N       ! Generate columns of Identity matrix, one at a time, and solve.
        print *, " Column ", I
        B = 0.0
        B(I) = 1.0

        if (I == 1) then
          call Factor (A, B, X, E, X_Norm, E_Norm, Flag)

          if (Flag) then
            print *, " X_Norm, E_Norm, X, E: ", X_Norm, E_Norm
            print *, X
            print *, E
            print *, " OK: ", Flag
          else
            print *, " Factor cannot solve this one."
            exit
          end if

        else
          call Solve (B, X, E, X_Norm, E_Norm)
          print *, " X_Norm, E_Norm, X, E: ", X_Norm, E_Norm
          print *, X
          print *, E
        end if

      end do

      call Unsolve ()
      print *, " Time to deallocate A, B, X, E: ", N
      deallocate (A, B, X, E)
    end do

    stop
  end program Driver

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Example 8.12a. Bisection procedure for finding an approximate
!  zero of any real function of one real argument.
!  No interface blocks except for dummy procedures.
  module Zero_M
    implicit none

  contains

    function Zero (Fun, A_In, B_In) result (C)

      interface
        function Fun (X) result (Y)
          real, intent (in) :: X
          real :: Y
        end function Fun
      end interface

      real, intent (in) :: A_In, B_In
      integer, parameter :: Max_Div = 20
      real :: A, B, C, FA, FB, FC
      integer :: Loop

  ! start function Zero
      A = A_In
      FA = Fun (A)
      B = B_In
      FB = Fun (B)
      if (FA * FB > 0.0) then  ! Check for error in initial values of A and B.
        print *, " Improper A and B values. "
      else        ! Subdivide the (A, B) interval as many as Max_Div times.
        do Loop = 1, Max_Div
          C = 0.5 * (A + B)
          FC = Fun (C)
          if (Abs (A - B) <= 100.0 * Spacing (A)) exit
          if (FC * FA > 0.0) then
            A = C
          else
            B = C
          end if
        end do
      end if
      return
    end function Zero

  end module Zero_M

! Example 8.12b.
  module Zero_Tests
    implicit none

  contains

  ! First test polynomial for Zero.
    function Poly_A (Z) result (A)  ! Module subprogram
      real, intent (in) :: Z
      real :: A
  ! start function Poly_A
      A = 17.3 + Z * (23.4 + Z * (137.0 + Z * (64.2 - Z)))
      return
    end function Poly_A

  ! Second test polynomial for Zero.
    function Poly_B (Z) result (B)  ! Module subprogram
      real, intent (in) :: Z
      real :: B
  ! start function Poly_B
      B = -102.4 + Z * (51.2 + Z * (25.6 + Z))
    return
    end function Poly_B

  ! Intrinsic function imbedded in module subprogram
    function My_Sin (Z) result (C)
      real, intent(in) :: Z
      real :: C
  ! start function My_Sin
      C = Sin (Z)  ! Selects Real version of generic intrinsic
      return
    end function My_Sin

  end module Zero_Tests

! Example 8.12c. Find zeros of Poly_A and Poly_B.
  program George
    use Zero_M
    use Zero_Tests  ! Procedures that appear as actual arguments
    implicit none

! start program George
    print * , Zero (Poly_A, 0.0, 100.0), Zero (Poly_B, 0.0, 100.0), Zero (My_Sin, 0.0, 100.0)
    stop

  end program George

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Array imiplementation of STACK

! Loren Meissner, CS 112, 4 Nov 93
  module Stack_Array_Module
    !**SPECIFICATION PART OF MODULE
    implicit none
    type :: Element_Type
      integer :: Key
      integer :: Sequence
    end type Element_Type
    type :: Stack_Type
      private
      integer :: Current_Size
      type (Element_Type), pointer, dimension (:) :: Space
    end type Stack_Type
    type :: Van
      logical :: L
      type (Element_Type) :: E
    end type Van

    integer, private :: Max_Size

  !**PROCEDURE PART OF MODULE
  contains

    subroutine Create_Stack(S_Ptr, N)
      type (Stack_Type), pointer :: S_Ptr
      integer, intent(in) :: N
  ! start subroutine Create_Stack
      allocate (S_Ptr)
      allocate (S_Ptr % Space(N))
      Max_Size = N
      S_Ptr % Current_Size = 0
      return
    end subroutine Create_Stack

    subroutine Destroy_Stack(S_Ptr)
      type (Stack_Type), pointer :: S_Ptr
  ! start subroutine Destroy_Stack
      S_Ptr % Current_Size = 0
      deallocate (S_Ptr % Space)
      deallocate (S_Ptr)
      return
    end subroutine Destroy_Stack

    function Push (S_Ptr, Element) result (LE)
      type (Stack_Type), pointer :: S_Ptr
      type (Element_Type), intent(in) :: Element
      type (Van) :: LE
      integer :: N
  ! start function Push
      N = S_Ptr % Current_Size + 1
      LE % E = Element
      LE % L = (N <= Max_Size)
      if (.not. LE % L) return
      S_Ptr % Space(N) = Element
      S_Ptr % Current_Size = N
      return
    end function Push

    function Pop (S_Ptr) result (LE)
      type (Stack_Type), pointer :: S_Ptr
      type (Van) :: LE
      integer :: N
  ! start function Pop
      N = S_Ptr % Current_Size
      LE % L = (N > 0)
      if (.not. LE % L) return
      LE % E = S_Ptr % Space(N)
      S_Ptr % Current_Size = N - 1
      return
    end function Pop

    function Peek (S_Ptr) result (LE)
      type (Stack_Type), pointer :: S_Ptr
      type (Van) :: LE
      integer :: N
  ! start function Peek
      N = S_Ptr % Current_Size
      LE % L = (N > 0)
      if (.not. LE % L) return
      LE % E = S_Ptr % Space(N)
      return
    end function Peek

  end module Stack_Array_Module

!**MAIN PROGRAM
  program Lab_Prog
    use Stack_Array_Module
    implicit none
    integer, parameter :: MANY = 6
    integer :: Op, Push_Sequence = 1
    type (Stack_Type), pointer :: Stack
    type (Van) :: Truck = Van (.TRUE., Element_Type (99, 0))

! start program Lab_Prog
    call Create_Stack (Stack, MANY)

    do
      print *, " Please enter Operation: "
      print *, " 0 to quit; 1 to push; 2 to pop, 3 to peek."
      read *, Op

      select case (Op)

      case (0)
        print *, " Operation 0 means quit."
        exit

      case (1)
        print *, " Please enter value to push."
        read *, Truck % E % Key
        Push_Sequence = Push_Sequence + 1
        Truck = Push (Stack, Truck % E)
        if (.not. Truck % L) then
          print *, " Error: Push failed."
          exit
        end if

      case (2)
        Truck = Pop (Stack)
        if (.not. Truck % L) then
          print *, " Error: Pop failed."
          exit
        else
          print *, " Popped: ", Truck % E % Key, Truck % E % Sequence
        end if

      case (3)
        Truck = Peek (Stack)
        if (.not. Truck % L) then
          print *, " Error: Peek failed."
          exit
        else
          print *, " Peeked: ", Truck % E % Key, Truck % E % Sequence
        end if

      case default
        print *, Op, " is not a valid operation."

      end select

    end do

    do
      Truck = Pop (Stack)
      if (.not. Truck % L) exit
      print *, " Left on stack: ", Truck% E % Key, Truck % E % Sequence
    end do

    call Destroy_Stack (Stack)

    stop

  end program Lab_Prog

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Example 8.11a. Function to perform quadrature by Simpsons rule.
!  No interface block except for dummy procedures.

   module Quadrature_M
      implicit none

   contains

      function Quadrature( Integrand, A, B, N ) result (Simpson)
         interface
            function Integrand( X ) result ( Y )
               real, intent (in) :: X
               real :: Y
            end function Integrand
         end interface
         real, intent(in) :: A, B
         integer, intent (in) :: N
         real :: Simpson
   ! start function Quadrature
         call Simp( Integrand, A, B, N, Simpson )
         return
      end function Quadrature

      subroutine Simp( Integrand, A, B, N_In, Simpson )
         interface
            function Integrand ( X ) result ( Y )
               real, intent(in) :: X
               real :: Y
            end function Integrand
         end interface
         real, intent (in) :: A, B
         integer, intent (in) :: N_In
         real, intent(out) :: Simpson
         real :: Evens, Odds, X, Delta_X
         integer :: I, N
   !  start subroutine Simp
         N = N_In + Modulo( N_In, 2 )   ! Add 1 to N if it is not even.
         Delta_X = (B - A) / Real( N )
         Evens = 0.0
         Odds = 0.0
         do I = 1, N - 1
            X = A + Real( I ) * Delta_X
            if (Modulo( I, 2 ) == 0) then
               Evens = Evens + Integrand( X )
            else
               Odds = Odds + Integrand( X )
            end if
         end do
         Simpson = Delta_X * &
            (Integrand( A ) + Integrand( B ) + 4.0 * Odds + 2.0 * Evens) / 3.0
         return
      end subroutine Simp

   end module Quadrature_M

  module Quad_Tests
  implicit none
  contains

  !  Fourth-degree polynomial test for Quadrature.
    real function Poly_4( Z ) result ( A )  ! External subprogram
      real, intent (in) :: Z
  !  start function Poly_4
      A = 17.3 + Z * (23.4 + Z * (137.0 + Z * (64.2 - Z)))
      return
    end function Poly_4

  !  Gaussian distribution test for Quadrature.
    real function Distribution( Z ) result ( B )  ! External subprogram
      real, intent (in) :: Z
  !  start function Distribution
     B = Exp( -Z ** 2 )
     return
    end function Distribution

  end module Quad_Tests

  program Ronald
    use Quadrature_M
    use Quad_Tests  ! Procedures that appear as actual arguments
    implicit none
! start program Ronald
    print *, Quadrature( Poly_4, -2.0, 2.0, 100 ), Quadrature( Distribution, -2.0, 2.0, 100 )
    stop
  end program Ronald

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
! I have found some relaxation of the restrictions in the Elf90 Spec:
!  Character dummy args can have fixed length; the spec says that assumed length is required.
!  The "nX" edit descriptor is supported; the spec says that it has to be replaced by "TRn".
! The following program runs OK under Elf90 v 1.00d -- Loren Meissner
! [Some of these may be rejected by in later versions.]

program Relax
implicit none
character (len = 14) :: String = " Hello world. "
! start program Relax
call Fix_It (String)
stop
contains
subroutine Fix_It (S)
character (len = 14), intent (in) :: S
! start subroutine Fix_It
write (*, "(' ...', 5x, A)") S
return
end subroutine Fix_It
end program Relax

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Linked List operations with Pointer to Pointer

   module Linked_List_Ops
      implicit none

      integer, parameter :: S_LEN = 20, LIMIT = Huge( 1 )

      type :: Info_Type
         character (len = S_LEN) :: Key
         integer :: Count(2)
      end type Info_Type

      type :: Box_Type
         private
         type (Node_Type), pointer :: NextP
      end type Box_Type

      type :: Node_Type
         private
         type (Info_Type) :: Info
         type (Box_Type) :: Box
      end type Node_Type

      type (Box_Type), pointer, private :: BTP  ! Local "travelling" box pointer

   contains

      subroutine Create_List( BR )
         type (Box_Type), target, intent(in out) :: BR

   !  start subroutine Create_List

         nullify (BR % NextP)
         return

      end subroutine Create_List

      subroutine Look_Up( BR, Item )
         type (Box_Type), target, intent(in) :: BR
         type (Info_Type), intent(in) :: Item

   !  start subroutine Look_Up

         BTP => BR  ! Make dummy arg BR the target of BTP

         do while (Associated ( BTP % NextP ))

            if (Item % Key == BTP % NextP % Info % Key) then
               print *, Item % Key, " Found. "
               call Modify_Target( BTP % NextP, Item % Count(1) )
               return
            else if (Item % Key > BTP % NextP % Info % Key) then   ! Keep looking.
               BTP => BTP % NextP % Box   ! New target of BTP is (next node) % Box
            else   ! Item % Key < BTP % NextP % Info % Key.
               print *, Item % Key, " Insert. "   ! Insert a new item.
               call Insert_Target( BTP % NextP, Item )
               return
            end if

         end do

         ! Pointer BTP % NextP is NIL.
         print *, Item % Key, " Insert. "   ! Insert a new item.
         call Insert_Target( BTP % NextP, Item )
         return

      contains

         subroutine Modify_Target( NQP, I )
            type (Node_Type), pointer :: NQP
            integer, intent(in) :: I

      !  start subroutine Modify_Target

            if (Associated( NQP )) then
               NQP % Info % Count(1) = I
               NQP % Info % Count(2) = NQP % Info % Count(2) + 1
            end if

            return

         end subroutine Modify_Target

         subroutine Insert_Target( NQP, Information )
            type (Node_Type), pointer :: NQP, N_TempP
            type (Info_Type), intent(in) :: Information
            integer :: AS

      !  start subroutine Insert_Target

            allocate (N_TempP, stat = AS)

            if (AS /= 0) then
               print *, Information
               stop " Allocation failure. "
            end if

            N_TempP % Info = Information
            N_TempP % Box % NextP => NQP
            NQP => N_TempP
            nullify (N_TempP)   ! Optional
            return

         end subroutine Insert_Target

      end subroutine Look_Up

      subroutine Print_List( BR )
         type (Box_Type), target, intent(in) :: BR
         type (Node_Type), pointer :: NTP

   !  start subroutine Print_List

         NTP => BR % NextP

         do while (Associated ( NTP ))
            call Print_Target( NTP )
            NTP => NTP % Box % NextP   ! Advance to next node
         end do

         return

      contains

         subroutine Print_Target( NQP )
            type (Node_Type), pointer :: NQP

      !  start subroutine Print_Target

            if (Associated( NQP )) print *, " Print: ", NQP % Info
            return

         end subroutine Print_Target

      end subroutine Print_List

      subroutine Delete_List( BR )
         type (Box_Type), target, intent(in) :: BR

   !  start subroutine Delete_List

         do while (Associated ( BR % NextP ))
            print *, " Delete: ", Delete_Target( BR % NextP )
         end do

         return

      contains

         function Delete_Target( NQP ) result (Information)
            type (Node_Type), pointer :: NQP, Temp
            type (Info_Type) :: Information

      !  start function Delete_Target

            Temp => NQP

            if (Associated( NQP )) then
              Information = NQP % Info
              NQP => NQP % Box % NextP  ! New target of NQP is (next node) % Box
              deallocate( Temp )
            end if

            return

         end function Delete_Target

      end subroutine Delete_List

   end module Linked_List_Ops

   program Literate
      use Linked_List_Ops
      implicit none
      type (Box_Type) :: Root_Box
      integer :: Loop, Eof
      type (Info_Type) :: New_Item

!  start program Literate

      call Create_List( Root_Box )
      open (1, file = "simply.txt", position = "REWIND")

      do Loop = 1, LIMIT
         read (1, *, iostat = EoF) New_Item % Key   ! Read to Key component of new structure.
         if (EoF < 0) exit
         New_Item % Count(1) = Loop
         New_Item % Count(2) = 1
         call Look_Up( Root_Box, New_Item )
      end do

      print *, " Contents of Linked List: "
      call Print_List( Root_Box )
      call Delete_List( Root_Box )
      stop " OK. "

   end program Literate

!==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ==== CUT HERE ====
!  copyright 1996 Loren P. Meissner -- May be distributed if this line is included.
!  Linked List operations with Recursive search

   module Linked_List_Ops
      implicit none

      integer, parameter :: S_LEN = 20, LIMIT = Huge( 1 )

      type :: Info_Type
         character (len = S_LEN) :: Key
         integer :: Count(2)
      end type Info_Type

      type :: Node_Type
         private
         type (Info_Type) :: Info
         type (Node_Type), pointer :: NextP
      end type Node_Type

   contains

      subroutine Create_List( RootP )
         type (Node_Type), pointer :: RootP

   !  start subroutine Create_List

         nullify (RootP)
         return

      end subroutine Create_List

      recursive subroutine Look_Up( P, Item )
         type (Node_Type), pointer :: P
         type (Info_Type), intent(in) :: Item

   !  start subroutine Look_Up

         if (.not. Associated ( P )) then
            call Insert_Target( P, Item )   ! Insert at end of linked list.
         else

            if (Item % Key == P % Info % Key) then
               call Modify_Target( P, Item % Count(1) )
            else if (Item % Key < P % Info % Key) then
               call Insert_Target( P, Item )   ! Insert ahead of larger node.
            else
               call Look_Up( P % NextP, Item )   ! Keep looking.
            end if
         end if

         return

       contains


         subroutine Modify_Target( NQP, I )
            type (Node_Type), pointer :: NQP
            integer, intent(in) :: I

      !  start subroutine Modify_Target

            if (Associated( NQP )) then
               NQP % Info % Count(1) = I
               NQP % Info % Count(2) = NQP % Info % Count(2) + 1
            end if

            return

         end subroutine Modify_Target

         subroutine Insert_Target( NQP, Information )
            type (Node_Type), pointer :: NQP, N_TempP
            type (Info_Type), intent(in) :: Information
            integer :: AS

      !  start subroutine Insert_Target

            allocate (N_TempP, stat = AS)

            if (AS /= 0) then
               print *, Information
               stop " Allocation failure. "
            end if

            N_TempP % Info = Information
            N_TempP % NextP => NQP
            NQP => N_TempP
            nullify (N_TempP)   ! Optional
            return

         end subroutine Insert_Target

      end subroutine Look_Up

      subroutine Print_List( RootP )
         type (Node_Type), pointer :: RootP, NTP

   !  start subroutine Print_List

         NTP => RootP

         do while (Associated ( NTP ))
            call Print_Target( NTP )
            NTP => NTP % NextP   ! Advance to next node
         end do

         return

      contains

         subroutine Print_Target( NQP )
            type (Node_Type), pointer :: NQP

      !  start subroutine Print_Target

            if (Associated( NQP )) print *, " Print: ", NQP % Info
            return

         end subroutine Print_Target

      end subroutine Print_List

      subroutine Delete_List( RootP )
         type (Node_Type), pointer :: RootP

   !  start subroutine Delete_List

         call Recursive_Delete_List( RootP )
         return

      contains

         recursive subroutine Recursive_Delete_List( P )
            type (Node_Type), pointer :: P

      !  start subroutine Recursive_Delete_List

            if (Associated ( P )) then
              call Recursive_Delete_List( P % NextP )
              print *, " Deleted: ", Delete_Target( P )
            end if
            return

         end subroutine Recursive_Delete_List

         function Delete_Target( NQP ) result (Information)
            type (Node_Type), pointer :: NQP, Temp
            type (Info_Type) :: Information

      !  start function Delete_Target

            Temp => NQP

            if (Associated( NQP )) then
              Information = NQP % Info
              NQP => NQP % NextP  ! New target of NQP is (next node) % Box
              deallocate( Temp )
            end if

            return

         end function Delete_Target

      end subroutine Delete_List

   end module Linked_List_Ops

   program Literate
      use Linked_List_Ops
      implicit none
      type (Node_Type), pointer :: RootP
      integer :: Loop, Eof
      type (Info_Type) :: New_Item

!  start program Literate

      call Create_List( RootP )
      open (1, file = "simply.txt", position = "REWIND")

      do Loop = 1, LIMIT
         read (1, *, iostat = EoF) New_Item % Key   ! Read to Key component of new structure.
         if (EoF < 0) exit
         New_Item % Count(1) = Loop
         New_Item % Count(2) = 1
         call Look_Up( RootP, New_Item )
      end do

      print *, " Contents of Linked List: "
      call Print_List( RootP )
      call Delete_List( RootP )
      stop " OK. "

   end program Literate


== f39.dat =====================================================================
3
0.11 0.12 0.12
0.11 0.10 0.14
0.11 0.14 0.11
0.10 0.10 0.11

== g17.dat =====================================================================
"Every Word In This Sentence Begins With A Capital Letter."

== h02.dat =====================================================================
(0.8660254, 0.5) (0.5, 0.8660254)

== ising.dat ===================================================================
(0.8660254, 0.5) (0.5, 0.8660254)

== char_in.txt =================================================================
X
Aaaaa
B
C
Z
" "
FUNNY_STUFF
8
9
10
a
Aaaa
Aaa
Aa

== compound.txt ================================================================
"zinc oxide"
"sulfuric acid"
"carbon monoxide"
"calcium carbonate"
"hydrochloric acid"
"calcium chloride"
"sodium chloride"
"carbon dioxide"
"iron (III) chloride"
"sodium hydroxide"
"calcium phosphate"
"sulfuric acid"
"calcium sulfate"
"zinc oxide"
"sodium chloride"
"phosphoric acid"
"silver nitrate"
"sodium chloride"
"aluminum hydroxide"
"sodium chloride"
"sulfuric acid"
"water"
"copper (II) sulfate"
"sodium chloride"
"silver nitrate"
"sodium chloride"
"copper (II) nitrate"
"sodium hydrogen sulfate"
"hydrochloric acid"
"sodium peroxide"
"hydrogen peroxide"
"barium chloride"
"sodium sulfate"
"sodium chloride"
"barium sulfate"
"sodium hydroxide"
"lead (II) acetate"
"sodium chloride"
"sodium chloride"
"iron (II) sulfide"
"lead (II) nitrate"
"iron (III) chloride"
"sodium bromide"

== file_in.txt =================================================================
11
9
10
6
8
7
5
1
3
2
0
2

== gaudata.txt =================================================================
4
1 2 2 0
1 0 4 0
1 4 1 0
0 0 1 1
1
1
0

== numbers.txt =================================================================
11
9
10
6
8
7
5
1
3
2
0
2

== numeroso.txt ================================================================
1001 3009
1002 3008
1003 3007
1004 3006
1005 3005
1006 3004
1007 3003
1008 3002
1009 3001

== prob4.txt ===================================================================
   -3.14159        22.5000
   -2.61799        17.3000
   -2.09440        15.3000
   -1.57080        14.7000
   -1.04720        15.3000
  -0.523599        17.3000
   0.000000        22.5000
   0.523599        35.8000
    1.04720        71.6000
    1.57080        117.900
    2.09440        71.6000
    2.61799        35.8000
    3.14159        22.5000

== simply.txt ==================================================================
SIMPLY
STATED
A
COMPILER
IS
A
PROGRAM
THAT
READS
A
PROGRAM
WRITTEN
IN
ONE
LANGUAGE
THE
SOURCE
LANGUAGE
AND
TRANSLATES
IT
INTO
AN
EQUIVALENT
PROGRAM
IN
ANOTHER
LANGUAGE
THE
TARGET
LANGUAGE
AS
AN
IMPORTANT
PART
OF
THIS
TRANSLATION
PROCESS
THE
COMPILER
REPORTS
TO
ITS
USER
THE
PRESENCE
OF
ERRORS
IN
THE
SOURCE
PROGRAM
AT
FIRST
GLANCE
THE
VARIETY
OF
COMPILERS
MAY
APPEAR
OVERWHELMING
THERE
ARE
THOUSANDS
OF
SOURCE
LANGUAGES
RANGING
FROM
TRADITIONAL
PROGRAMMING
LANGUAGES
SUCH
AS
FORTRAN
AND
PASCAL
TO
SPECIALIZED
LANGUAGES
THAT
HAVE
ARISEN
IN
VIRTUALLY
EVERY
AREA
OF
COMPUTER
APPLICATION

== str_500.txt == COMBINE THE LINES OF THIS FILE INTO A SINGLE RECORD!!! =======
"A character string consists of zero or more characters. Even though it is made 
up of individual characters, a character string is considered to be a scalar. As 
 with any data type, it is possible to declare an array of character strings, al
l the same length. A substring is a contiguous portion of a character string tha
t has a starting point and an ending point within the character string. It is po
ssible to reference a substring of a character scalar variable or constant."