program TEST_CST
!     Test of sine, cosine & tangent
! 1. sin(2.x)  = 2.sin(x).cos(x)
! 2. cos(2.x)  = (cos(x) + sin(x)).(cos(x) - sin(x))
! 3. sec(x)**2 = 1 + tan(x)**2

use QUADRUPLE_PRECISION
implicit none

type (QUAD)                        :: X, TWOX, SINE, COSINE, TANGENT, LHS, &
                                      RHS, DIFF
real (kind=DP), parameter          :: HALF = 0.5_DP, SMALL = epsilon(HALF)
real (kind=DP)                     :: R
integer, allocatable, dimension(:) :: SEED
integer                            :: K, I

!     Set the random number seed.

call random_seed(size=K)
allocate (SEED(K))
call random_seed(get=SEED)
write(unit=*, fmt=*)"Old random number seeds: ", SEED

write(unit=*, fmt="(a, i4, a)") " Enter ", K, " integers as random number seeds: "
read(unit=*, fmt=*) SEED
call random_seed(put=SEED)

do I = 1, 100
  call random_number(X%HI)
  call random_number(R)
  X%HI = (X%HI - HALF) / R
  X%LO = X%HI * SMALL
                             ! sin(2.x)  = 2.sin(x).cos(x)
  TWOX = 2.0_DP * X
  LHS = sin(TWOX)
  SINE = sin(X)
  COSINE = cos(X)
  RHS = SINE * COSINE
  RHS = 2.0_DP * RHS
  DIFF = LHS - RHS
  write(unit=*, fmt="(a, g13.5, a, g12.4)")  &
                    " sin(2x)   lhs =", LHS%HI, "  Diff. =", DIFF%HI

                             ! cos(2.x)  = (cos(x) + sin(x)).(cos(x) - sin(x))
  LHS = cos(TWOX)
  RHS = (COSINE + SINE) * (COSINE - SINE)
  DIFF = LHS - RHS
  write(unit=*, fmt="(a, g13.5, a, g12.4)")  &
                    " cos(2x)   lhs =", LHS%HI, "  Diff. =", DIFF%HI

                             ! sec(x)**2 = 1 + tan(x)**2
  DIFF%HI = 1.0_DP
  DIFF%LO = 0.0_DP
  LHS = DIFF / COSINE
  LHS = LHS * LHS
  TANGENT = tan(X)
  RHS = DIFF + TANGENT * TANGENT
  DIFF = LHS - RHS
  write(unit=*, fmt="(a, g13.5, a, g12.4)")  &
                    " sec(x)^2  lhs =", LHS%HI, "  Diff. =", DIFF%HI
end do

stop
end program TEST_CST
