Degenerate Conic

Binary Search

Binary search is a search algorithm that finds the position of a target value within a sorted array. The following binary search routine is taken from NASTRAN, a venerable finite element analysis program from NASA:

    SUBROUTINE BISLOC (*,ID,ARR,LEN,KN,JLOC)
C-----
C BINARY SEARCH - LOCATE KEY WORD 'ID' IN ARRAY 'ARR', 1ST ENTRY
C IF FOUND, 'JLOC' IS THE MATCHED POSITION IN 'ARR'
C IF NOT FOUND, NON-STANDARD RETURN
C I.E.
C ID = KEY WORD TO MATCH IN ARR. MATCH AGAINST 1ST COL OF ARR
C ARR = ARRAY TO SEARCH. ARR(ROW,COL)
C LEN = LENGTH OF EACH ENTRY IN ARRAY. LEN=ROW
C KN = NUMBER OF ENTRIES IN THE ARR. KN =COL
C JLOC= POINTER RETURNED - FIRST WORD OF ENTRY. MATCHED ROW
C-----
C
    INTEGER ARR(1)
    DATA ISWTCH / 16 /
C
    JJ = LEN - 1
    IF (KN .LT. ISWTCH) GO TO 120
    KLO = 1
    KHI = KN
10  K = (KLO+KHI+1)/2
20  J = K*LEN - JJ
    IF (ID-ARR(J)) 30,90,40
30  KHI = K
    GO TO 50
40  KLO = K
50  IF (KHI-KLO -1) 100,60,10
60  IF (K .EQ. KLO) GO TO 70
    K = KLO
    GO TO 80
70  K = KHI
80  KLO = KHI
    GO TO 20
90  JLOC = J
    RETURN
100 JLOC = KHI*LEN - JJ
    J = KN *LEN - JJ
    IF (ID .GT.ARR(J)) JLOC = JLOC + LEN
110 RETURN 1
C
C SEQUENTIAL SEARCH MORE EFFICIENT
C
120 KHI = KN*LEN - JJ
    DO 130 J = 1,KHI,LEN
    IF (ARR(J)-ID) 130,90,140
130 CONTINUE
    JLOC = KHI + LEN
    GO TO 110
140 JLOC = J
    GO TO 110
    END

As you can see, it is a scene of FORTRAN 66 horrors! However, it can be refactored into a nice modern Fortran routine like so:

pure function bisloc(id,arr) result(jloc)

!! binary search of a sorted array.

implicit none

integer,intent(in) :: id
    !! key word to match in `arr`
integer,dimension(:),intent(in) :: arr
    !! array to search (it is
    !! assumed to be sorted)
integer :: jloc
    !! the first matched index in 'arr'
    !! (if not found, 0 is returned)

integer :: j,k,khi,klo,n
integer,parameter :: iswtch = 16

n = size(arr)
jloc = 0

if ( n<iswtch ) then

    ! sequential search more efficient

    do j = 1 , n
        if ( arr(j)==id ) then
            jloc = j
            return
        else if (arr(j)>id) then
            return ! error
        end if
    end do
    return ! error

else

    klo = 1
    khi = n
    k = (klo+khi+1)/2
    do
        j = k
        if ( id<arr(j) ) then
            khi = k
        else if ( id==arr(j) ) then
            jloc = j
            return
        else
            klo = k
        end if
        if ( khi-klo<1 ) then
            return ! error
        else if ( khi-klo==1 ) then
            if ( k==klo ) then
                k = khi
            else
                k = klo
            end if
            klo = khi
        else
            k = (klo+khi+1)/2
        end if
    end do

end if

end function bisloc

I did change a couple of things:

The original routine was using an obsolete alternate return (note the * argument) for an error condition. The new one just returns 0 if there is an error.
The original routine was passing in the array as a vector and treating it like a matrix (it's either searching a row or column of this matrix, I admit I had a hard time telling). The new version assumes the input is a rank-1 array. You can pass in whatever matrix slices you want (say a column: a(:,1) or a row: a(1,:)).
I also made it a pure function, because that's how I roll.

nasaLogo

Note that NASTRAN was released under the NASA Open Source Agreement version 1.3, a somewhat unusual open source license. It isn't considered a "free" license by the FSF since it includes a statement requiring changes to be your "original creation".

Unfortunately, the trend at NASA (and elsewhere) is to declare old Fortran code like this to be totally unmaintainable. So, instead of refactoring them into modern Fortran standards, they are thrown away and rewritten from scratch using a more "modern" language like C++ (whose roots go back to 1972!).

Algorithms • Modern Fortran Programming • Orbital Mechanics

Mar 11, 2017

Binary Search

See also