Re: busting sp datatypes

On Thu, 28 Aug 2008 03:15:32 +0000 (UTC), Steven G. Kargl posted:

In article <otfdd3cv5rxp.dlg@xxxxxxxxxxxxxxx>,
Ron Ford <ron@xxxxxxxxxxxxxxx> writes:
! inner loop rolls die till all values attained
do
call random_number(harvest)
b = ceiling(harvest*real(sides))

b = min(max(1,ceiling(harvest*sides)),sides)

I thought the problem was that I recycled b later in the program. I cut
out that as well as the inferior selection sort and added a print staement
to diagnose:

! quicksort for largish ints

module qsort_c_module

! Recursive Fortran 95 quicksort routine
! Author: Juli Rew

implicit none
public :: QsortC
private :: Partition

contains

recursive subroutine QsortC(A)
integer, parameter :: i13 = selected_int_kind(13)
integer(i13), intent(in out), dimension(:) :: A
integer (i13):: iq

if(size(A) > 1) then
call Partition(A, iq)
call QsortC(A(:iq-1))
call QsortC(A(iq:))
endif
end subroutine QsortC

subroutine Partition(A, marker)
integer, parameter :: i13 = selected_int_kind(13)
integer(i13), intent(in out), dimension(:) :: A
integer(i13), intent(out) :: marker
integer(i13) :: i, j
integer(i13) :: temp
integer(i13) :: x ! pivot point
x = A(1)
i= 0
j= size(A) + 1

do
j = j-1
do
if (A(j) <= x) exit
j = j-1
end do
i = i+1
do
if (A(i) >= x) exit
i = i+1
end do
if (i < j) then
! exchange A(i) and A(j)
temp = A(i)
A(i) = A(j)
A(j) = temp
elseif (i == j) then
marker = i+1
return
else
marker = i
return
endif
end do

end subroutine Partition

end module qsort_c_module

use qsort_c_module

implicit none
integer, parameter :: i13 = selected_int_kind(13)

integer(i13), parameter:: sides = 8
integer(i13), parameter:: trials = 300
integer(i13), parameter:: array_size = trials
integer(i13), parameter:: bins = 50
integer(i13), parameter:: percentile = 95
integer(i13), dimension(sides)::A
integer(i13), dimension(trials)::C
integer(i13), dimension(bins)::D
integer(i13), dimension(trials)::F
integer(i13):: b, ii, i, counter, &
tab, maxput, diff2, c3, c4
real:: harvest, tot, diff, t3, t4

! seed random num generator

CALL init_random_seed

! prime the pump
call random_number(harvest)
b = 3 + nint(10*harvest)
do i=1,b
call random_number(harvest)
print *, i, harvest
end do


! main control

! outer loop is the number of trials
C = 0
do ii=1,trials
tab = 0
A = 0
! inner loop rolls die till all values attained
do
call random_number(harvest)
b = ceiling(harvest*real(sides))
print *, "b is ", b, " A(b) is", A(b)
A(b) = A(b) + 1

! count until all outcomes non-zero
counter = 0
do i = 1, sides
if (real(A(i)) > .5) then
counter = counter + 1
end if
end do
tab = tab + 1
if (counter == sides) then
C(ii) = tab
exit
end if
end do !inner
end do !outer

! end main control
! print *, "C is", C

! sort into bins
D=0
maxput = bins
do ii =1, trials
if (C(ii) > maxput) then
C(ii) = maxput
end if
D(C(ii))=D(C(ii)) + 1
end do
tot= sum(D)
print *, "total trials=",tot

! determine 95th percentile
tot = 0.01 *trials*real(percentile)
print *, "95th % is", tot

counter = 0
do ii = 1, bins
counter = counter + D(ii)

if (real(counter).ge. tot) then
exit
end if
end do
print *, "95th percentile was in bin number", ii

! time the sort
F=C
call system_clock (c3)
call cpu_time (t3)
call QsortC(F)
call cpu_time (t4)
call system_clock (c4)
diff=t4-t3
diff2=c4-c3

! output
print *, "system clock for qsort was ", diff
print *, "cpu time for qsort was ", diff2
! end output
contains
SUBROUTINE init_random_seed()
INTEGER :: i, n, clock
INTEGER, DIMENSION(:), ALLOCATABLE :: seed

CALL RANDOM_SEED(size = n)
print *, "n=", n
ALLOCATE(seed(n))

CALL SYSTEM_CLOCK(COUNT=clock)
print *, "clock=", clock

seed = clock + 37 * (/ (i - 1, i = 1, n) /)
CALL RANDOM_SEED(PUT = seed)
print *, "seed= ", seed

DEALLOCATE(seed)
END SUBROUTINE
end program
! gfortran -o sort -fbounds-check freeformat52.f95

output is:
b is 7 A(b) is 0
b is 6 A(b) is 0
b is 2 A(b) is 0
b is 2 A(b) is 1
b is 8 A(b) is 0
b is 8 A(b) is 1
b is 1 A(b) is 0
b is 2 A(b) is 2
b is 3 A(b) is 0
b is 2 A(b) is 3
b is 1 A(b) is 1
b is 4 A(b) is 0
b is 4 A(b) is 1
b is 4 A(b) is 2
b is 5 A(b) is 0
b is 7 A(b) is 0
total trials= 300.000
95th % is 285.000
95th percentile was in bin number 39
system clock for qsort was 0.00000
cpu time for qsort was 3

Press RETURN to close window . . .


A(b) = A(b) + 1 !this is line 158

F:\gfortran\source>sort
n= 8
clock= 108897905
seed= 108897905 108897942 108897979 108898016 108898053
108898090
108898127 108898164
1 0.42305011
2 8.57861638E-02
3 0.40694624
At line 158 of file freeformat50.f95
Fortran runtime error: Array reference out of bounds for array 'a', lower
bound of dimension 1 exceeded (0 < 1)

Read the error message.

The error message leads me to think that numbers get so large that A(b) is
A(0). b doesn't get bigger than 8. A(b) doesn't exceed one hundred.

Then read the description of random_number().

I thought I worked these kinks out earlier.
--
Unquestionably, there is progress. The average American now pays out twice
as much in taxes as he formerly got in wages. 1
H. L. Mencken
.

Relevant Pages

  • Re: ordering integer array
    ... recursive subroutine selection_sort ... tab, maxput ... real:: harvest, tot ... sort into bins ...
    (comp.lang.fortran)
  • Re: ordering integer array
    ... recursive subroutine selection_sort ... tab, maxput ... real:: harvest, tot ... sort into bins ...
    (comp.lang.fortran)
  • Re: adapting a quicksort
    ... recursive subroutine selection_sort ... type:: temp ... tab, maxput ... real:: harvest, tot ...
    (comp.lang.fortran)
  • Re: ordering integer array
    ... 50 bins would suffice well for this sort. ... I use maxput to deal with the ... tab, maxput, penultimate, goal ... real:: harvest, tot ...
    (comp.lang.fortran)
  • Re: ordering integer array
    ...   implicit none ...   recursive subroutine selection_sort ... real:: harvest, tot ... sort into bins ...
    (comp.lang.fortran)