annotate src/fftw-3.3.3/genfft/util.ml @ 169:223a55898ab9 tip default

Add null config files
author Chris Cannam <cannam@all-day-breakfast.com>
date Mon, 02 Mar 2020 14:03:47 +0000
parents 89f5e221ed7b
children
rev   line source
cannam@95 1 (*
cannam@95 2 * Copyright (c) 1997-1999 Massachusetts Institute of Technology
cannam@95 3 * Copyright (c) 2003, 2007-11 Matteo Frigo
cannam@95 4 * Copyright (c) 2003, 2007-11 Massachusetts Institute of Technology
cannam@95 5 *
cannam@95 6 * This program is free software; you can redistribute it and/or modify
cannam@95 7 * it under the terms of the GNU General Public License as published by
cannam@95 8 * the Free Software Foundation; either version 2 of the License, or
cannam@95 9 * (at your option) any later version.
cannam@95 10 *
cannam@95 11 * This program is distributed in the hope that it will be useful,
cannam@95 12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
cannam@95 13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
cannam@95 14 * GNU General Public License for more details.
cannam@95 15 *
cannam@95 16 * You should have received a copy of the GNU General Public License
cannam@95 17 * along with this program; if not, write to the Free Software
cannam@95 18 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
cannam@95 19 *
cannam@95 20 *)
cannam@95 21
cannam@95 22 (* various utility functions *)
cannam@95 23 open List
cannam@95 24 open Unix
cannam@95 25
cannam@95 26 (*****************************************
cannam@95 27 * Integer operations
cannam@95 28 *****************************************)
cannam@95 29 (* fint the inverse of n modulo m *)
cannam@95 30 let invmod n m =
cannam@95 31 let rec loop i =
cannam@95 32 if ((i * n) mod m == 1) then i
cannam@95 33 else loop (i + 1)
cannam@95 34 in
cannam@95 35 loop 1
cannam@95 36
cannam@95 37 (* Yooklid's algorithm *)
cannam@95 38 let rec gcd n m =
cannam@95 39 if (n > m)
cannam@95 40 then gcd m n
cannam@95 41 else
cannam@95 42 let r = m mod n
cannam@95 43 in
cannam@95 44 if (r == 0) then n
cannam@95 45 else gcd r n
cannam@95 46
cannam@95 47 (* reduce the fraction m/n to lowest terms, modulo factors of n/n *)
cannam@95 48 let lowest_terms n m =
cannam@95 49 if (m mod n == 0) then
cannam@95 50 (1,0)
cannam@95 51 else
cannam@95 52 let nn = (abs n) in let mm = m * (n / nn)
cannam@95 53 in let mpos =
cannam@95 54 if (mm > 0) then (mm mod nn)
cannam@95 55 else (mm + (1 + (abs mm) / nn) * nn) mod nn
cannam@95 56 and d = gcd nn (abs mm)
cannam@95 57 in (nn / d, mpos / d)
cannam@95 58
cannam@95 59 (* find a generator for the multiplicative group mod p
cannam@95 60 (where p must be prime for a generator to exist!!) *)
cannam@95 61
cannam@95 62 exception No_Generator
cannam@95 63
cannam@95 64 let find_generator p =
cannam@95 65 let rec period x prod =
cannam@95 66 if (prod == 1) then 1
cannam@95 67 else 1 + (period x (prod * x mod p))
cannam@95 68 in let rec findgen x =
cannam@95 69 if (x == 0) then raise No_Generator
cannam@95 70 else if ((period x x) == (p - 1)) then x
cannam@95 71 else findgen ((x + 1) mod p)
cannam@95 72 in findgen 1
cannam@95 73
cannam@95 74 (* raise x to a power n modulo p (requires n > 0) (in principle,
cannam@95 75 negative powers would be fine, provided that x and p are relatively
cannam@95 76 prime...we don't need this functionality, though) *)
cannam@95 77
cannam@95 78 exception Negative_Power
cannam@95 79
cannam@95 80 let rec pow_mod x n p =
cannam@95 81 if (n == 0) then 1
cannam@95 82 else if (n < 0) then raise Negative_Power
cannam@95 83 else if (n mod 2 == 0) then pow_mod (x * x mod p) (n / 2) p
cannam@95 84 else x * (pow_mod x (n - 1) p) mod p
cannam@95 85
cannam@95 86 (******************************************
cannam@95 87 * auxiliary functions
cannam@95 88 ******************************************)
cannam@95 89 let rec forall id combiner a b f =
cannam@95 90 if (a >= b) then id
cannam@95 91 else combiner (f a) (forall id combiner (a + 1) b f)
cannam@95 92
cannam@95 93 let sum_list l = fold_right (+) l 0
cannam@95 94 let max_list l = fold_right (max) l (-999999)
cannam@95 95 let min_list l = fold_right (min) l 999999
cannam@95 96 let count pred = fold_left
cannam@95 97 (fun a elem -> if (pred elem) then 1 + a else a) 0
cannam@95 98 let remove elem = List.filter (fun e -> (e != elem))
cannam@95 99 let cons a b = a :: b
cannam@95 100 let null = function
cannam@95 101 [] -> true
cannam@95 102 | _ -> false
cannam@95 103 let for_list l f = List.iter f l
cannam@95 104 let rmap l f = List.map f l
cannam@95 105
cannam@95 106 (* functional composition *)
cannam@95 107 let (@@) f g x = f (g x)
cannam@95 108
cannam@95 109 let forall_flat a b = forall [] (@) a b
cannam@95 110
cannam@95 111 let identity x = x
cannam@95 112
cannam@95 113 let rec minimize f = function
cannam@95 114 [] -> None
cannam@95 115 | elem :: rest ->
cannam@95 116 match minimize f rest with
cannam@95 117 None -> Some elem
cannam@95 118 | Some x -> if (f x) >= (f elem) then Some elem else Some x
cannam@95 119
cannam@95 120
cannam@95 121 let rec find_elem condition = function
cannam@95 122 [] -> None
cannam@95 123 | elem :: rest ->
cannam@95 124 if condition elem then
cannam@95 125 Some elem
cannam@95 126 else
cannam@95 127 find_elem condition rest
cannam@95 128
cannam@95 129
cannam@95 130 (* find x, x >= a, such that (p x) is true *)
cannam@95 131 let rec suchthat a pred =
cannam@95 132 if (pred a) then a else suchthat (a + 1) pred
cannam@95 133
cannam@95 134 (* print an information message *)
cannam@95 135 let info string =
cannam@95 136 if !Magic.verbose then begin
cannam@95 137 let now = Unix.times ()
cannam@95 138 and pid = Unix.getpid () in
cannam@95 139 prerr_string ((string_of_int pid) ^ ": " ^
cannam@95 140 "at t = " ^ (string_of_float now.tms_utime) ^ " : ");
cannam@95 141 prerr_string (string ^ "\n");
cannam@95 142 flush Pervasives.stderr;
cannam@95 143 end
cannam@95 144
cannam@95 145 (* iota n produces the list [0; 1; ...; n - 1] *)
cannam@95 146 let iota n = forall [] cons 0 n identity
cannam@95 147
cannam@95 148 (* interval a b produces the list [a; 1; ...; b - 1] *)
cannam@95 149 let interval a b = List.map ((+) a) (iota (b - a))
cannam@95 150
cannam@95 151 (*
cannam@95 152 * freeze a function, i.e., compute it only once on demand, and
cannam@95 153 * cache it into an array.
cannam@95 154 *)
cannam@95 155 let array n f =
cannam@95 156 let a = Array.init n (fun i -> lazy (f i))
cannam@95 157 in fun i -> Lazy.force a.(i)
cannam@95 158
cannam@95 159
cannam@95 160 let rec take n l =
cannam@95 161 match (n, l) with
cannam@95 162 (0, _) -> []
cannam@95 163 | (n, (a :: b)) -> a :: (take (n - 1) b)
cannam@95 164 | _ -> failwith "take"
cannam@95 165
cannam@95 166 let rec drop n l =
cannam@95 167 match (n, l) with
cannam@95 168 (0, _) -> l
cannam@95 169 | (n, (_ :: b)) -> drop (n - 1) b
cannam@95 170 | _ -> failwith "drop"
cannam@95 171
cannam@95 172
cannam@95 173 let either a b =
cannam@95 174 match a with
cannam@95 175 Some x -> x
cannam@95 176 | _ -> b