search.f

92 lines · 1 definitions · on GitHub (7be4133)

Python version: search.py, 42 lines, on the right (switch to Python above the code).

Variables on this page (1)

json

Purpose

A binary search in an increasing table, part of the shape-preserving quadratic splines of McAllister and Roulier (coded by S. L. Dodd and M. Roulier, N.C. State University, as the header says).

How it fits

Its only caller is MEVAL, which evaluates the spline and is called by OSPLIN. No COMMON blocks.

The method

Given XTAB(1..NUM) and a value S with XTAB(1) ≤ S ≤ XTAB(NUM), it returns LCN, the index of the last entry below S, or the index of the entry equal to S (and then FND = 1). The interval is halved until it has length one. FND relies on exact equality of reals (line 55 ↗, line 60 ↗, line 77 ↗), so a value equal to a node up to rounding is reported as not found; MEVAL then evaluates the spline of the interval next to the node, which passes through the node's value, so the result differs only by about the same rounding (settled from the code, Nov 9, 2026: spline.f, lines 66–73 ↗).

Definitions

subroutineSEARCH

lines 4–92 · GitHub

Uses modules DIMENSIONS_MOD

Called by MEVAL

Compiler warnings (3)

gfortran -Wall -Wextra, compiling only; nothing was run.

line 55: Inequality comparison for REAL(8) [-Wcompare-reals]
line 60: Inequality comparison for REAL(8) [-Wcompare-reals]
line 77: Equality comparison for REAL(8) [-Wcompare-reals]

Change log

The changes since rp_1

Changed since rp_1. Declarations only: the array size JSON from DIMENSIONS_MOD, IMPLICIT NONE.

1C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 2C SEARCH 3C$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ 4 SUBROUTINE SEARCH(XTAB,NUM,S,LCN,FND) 5C 10/15/26 MHP via Claude 6C code changed to take JSON from DIMENSIONS_MOD instead of declaring 7C PARAMETER(JSON=5000) here. 8 USE DIMENSIONS_MOD 9C 10/20/26 MHP via Claude 10C code changed to replace implicit typing with IMPLICIT NONE and 11C explicit declarations. 12 IMPLICIT NONE 13 INTEGER LAST,LCN,MIDDLE,NUM 14C 15C SHAPE PRESERVING QUADRATIC SPLINES 16C BY D.F.MCALLISTER & J.A.ROULIER 17C CODED BY S.L.DODD & M.ROULIER 18C N.C. STATE UNIVERSITY 19C 20 REAL*8 XTAB(JSON),S 21 INTEGER FND,FIRST 22C 23C SEARCH CONDUCTS A BINARY SEARCH FOR S. SEARCH IS CALLED ONLY IF S IS 24C BETWEEN XTAB(1) AND XTAB(NUM). 25C 26C ON INPUT-- 27C 28C XTAB CONTAINS THE ABSCISSAS OF THE DATA POINTS OF INTERPOLATION. 29C 30C NUM IS THE DIMENSION OF XTAB 31C 32C S IS THE VALUE WHOSE RELATIVE POSITION IN XTAB IS LOCATED BY SEARCH. 33C 34C 35C ON OUTPUT-- 36C 37C FND IS SET EQUAL TO 1 IF S IS FOUND IN XTAB AND IS SET EQUAL TO 0 38C OTHERWISE. 39C 40C LCN IS THE INDEX OF THE LARGEST VALUE IN XTAB FOR WHICH XTAB(I) 41C .LT. S. 42C 43C AND 44C 45C SEARCH DOES NOT ALTER XTAB,NUM,S. 46C 47C---------------------------------------------------------------------- 48C 49 save 50 51 FIRST=1 52 LAST=NUM 53 FND=0 54C 55 IF (XTAB(1) .NE. S) GO TO 10 56 LCN=1 57 FND=1 58 RETURN 59C 60 10 10 IF (XTAB(NUM) .NE. S) GO TO 20 61 LCN=NUM 62 FND=1 63 RETURN 64C 65C IF (LAST-FIRST) .EQ. 1, S IS NOT IN XTAB. SET POSITION EQUAL TO 66C FIRST. 67 20 20 IF ((LAST-FIRST) .EQ. 1) GO TO 30 68C 69 MIDDLE=(FIRST+LAST)/2 70C 71C CHECK IF S .EQ. XTAB(MIDDLE). IF NOT, CONTINUE THE SEARCH IN THE 72C APPROPRIATE HALF OF THE VECTOR XTAB. 73C KC 2025-05-30 fixed "Arithmetic IF statement" 74C IF (XTAB(MIDDLE) - S) 40,50,60 75 IF (XTAB(MIDDLE) .LT. S) THEN 76 GOTO 40 77 ELSE IF (XTAB(MIDDLE) .EQ. S) THEN 78 GOTO 50 79 ELSE 80 GOTO 60 81 END IF 82C 83 30 30 LCN=FIRST 84 RETURN 85 40 40 FIRST=MIDDLE 86 GO TO 20 87 50 50 LCN=MIDDLE 88 FND =1 89 RETURN 90 60 60 LAST=MIDDLE 91 GO TO 20 92 END

1"""Python version of search.f (YREC 7be4133): a binary search in an 2increasing table, for the shape-preserving quadratic splines of McAllister 3and Roulier.""" 4 5import numpy as np 6 7 8def search(xtab: np.ndarray, num: int, s: float, lcn: int, fnd: int) -> tuple[int, int]: 9 """Locate s in xtab[1..num], with xtab[1] <= s <= xtab[num]: lcn is the 10 index of the last entry below s, or of the entry equal to s (then 11 fnd = 1; otherwise fnd = 0). 12 13 Translates: search.f#SEARCH (7be4133) 14 15 Arguments: xtab[1..num] (increasing); num; s; lcn and fnd, set only. 16 Returns (lcn, fnd). 17 18 """ 19 20 # Fortran lines 51–53. 21 first = 1 22 last = num 23 fnd = 0 24 25 # The ends of the table (Fortran lines 55–63); equality of reals is 26 # exact (see the page). 27 if xtab[1] == s: 28 return 1, 1 29 if xtab[num] == s: # label 10 30 return num, 1 31 32 # Bisection (labels 20–60, Fortran lines 67–91). 33 while last - first != 1: # label 20 34 middle = (first + last) // 2 35 if xtab[middle] < s: 36 first = middle # label 40 37 elif xtab[middle] == s: 38 return middle, 1 # label 50 39 else: 40 last = middle # label 60 41 lcn = first # label 30 42 return lcn, fnd