1#! /usr/bin/bc 2# 3# SPDX-License-Identifier: BSD-2-Clause 4# 5# Copyright (c) 2018-2023 Gavin D. Howard and contributors. 6# 7# Redistribution and use in source and binary forms, with or without 8# modification, are permitted provided that the following conditions are met: 9# 10# * Redistributions of source code must retain the above copyright notice, this 11# list of conditions and the following disclaimer. 12# 13# * Redistributions in binary form must reproduce the above copyright notice, 14# this list of conditions and the following disclaimer in the documentation 15# and/or other materials provided with the distribution. 16# 17# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 18# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 19# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 20# ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE 21# LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 22# CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 23# SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 24# INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 25# CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 26# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 27# POSSIBILITY OF SUCH DAMAGE. 28# 29 30scale = 20 31 32# Adjust this number to try ranges below different powers of 10. 33shift = 4 34 35# Adjust this to try extra digits. For example, a value of one means that one 36# digit is checked (such as 0.09 through 0.01), a value of two means that two 37# digits are checked (0.090 through 0.010), etc. 38max = shift + 2 39 40n = (9 >> shift) 41inc = (1 >> max) 42stop = (1 >> shift) 43 44# Uncomment this to test the high part of the ranges. 45#n += (1 - (1 >> max + 5)) >> shift 46 47for (i = n; i >= stop; i -= inc) 48{ 49 # This is the lower limit. 50 t1 = sqrt(1/(3*i)) 51 52 # Start with the inverse. 53 t2 = (1/i) 54 55 # And take half its length of course. 56 l = length(t2$)/2 57 58 temp = i 59 odd = 0 60 61 # We go by powers of 10 below, but there is a degenerate case: an exact 62 # power of 10, for which length() will return one digit more. So we check 63 # for that and fix it. 64 while (temp < 1) 65 { 66 temp <<= 1 67 odd = !odd 68 } 69 70 if (temp == 1) 71 { 72 odd = !odd 73 } 74 75 print "i: ", i, "\n" 76 print "t2: ", t2, "\n" 77 #print "l: ", l, "\n" 78 print "odd: ", odd, "\n" 79 80 if (odd) 81 { 82 # Limit between 6 and 7.5. 83 limit1 = 6.7 >> (l$ * 2 + 1) 84 85 # Limit between 1.5 and 1.83-ish. 86 limit2 = 1.7 >> (l$ * 2 + 1) 87 print "limit1: ", limit1, "\n" 88 print "limit2: ", limit2, "\n" 89 90 if (i >= limit1) 91 { 92 t2 = (t2 >> l$) 93 } 94 else if (i >= limit2) 95 { 96 t2 = (t2 >> l$) / 2 97 } 98 else 99 { 100 t2 = (t2 >> l$) / 4 101 } 102 } 103 else 104 { 105 # Limit between 2.4 and 3. 106 limit = 2.7 >> (l$ * 2) 107 print "limit: ", limit, "\n" 108 109 if (i >= limit) 110 { 111 t2 = (t2 >> l$) * 2 112 } 113 else 114 { 115 t2 = (t2 >> l$) 116 } 117 } 118 #t2 = 1 119 t3 = sqrt(5/(3*i)) 120 good = (t1 < t2 && t2 < t3) 121 122 print t1, " < ", t2, " < ", t3, ": ", good, "\n\n" 123 if (!good) sqrt(-1) 124} 125 126halt 127