~ubuntu-branches/ubuntu/quantal/ruby1.9.1/quantal

« back to all changes in this revision

Viewing changes to lib/prime.rb

  • Committer: Bazaar Package Importer
  • Author(s): Lucas Nussbaum
  • Date: 2011-09-24 19:16:17 UTC
  • mfrom: (1.1.8 upstream) (13.1.7 experimental)
  • Revision ID: james.westby@ubuntu.com-20110924191617-o1qz4rcmqjot8zuy
Tags: 1.9.3~rc1-1
* New upstream release: 1.9.3 RC1.
  + Includes load.c fixes. Closes: #639959.
* Upload to unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
46
46
# The set of all prime numbers.
47
47
#
48
48
# == Example
49
 
#  Prime.each(100) do |prime|
50
 
#    p prime  #=> 2, 3, 5, 7, 11, ...., 97
51
 
#  end
 
49
#
 
50
#   Prime.each(100) do |prime|
 
51
#     p prime  #=> 2, 3, 5, 7, 11, ...., 97
 
52
#   end
 
53
#
 
54
# Prime is Enumerable:
 
55
#
 
56
#   Prime.first 5 # => [2, 3, 5, 7, 11]
52
57
#
53
58
# == Retrieving the instance
 
59
#
54
60
# +Prime+.new is obsolete. Now +Prime+ has the default instance and you can
55
61
# access it as +Prime+.instance.
56
62
#
58
64
# as a class method of +Prime+.
59
65
#
60
66
# e.g.
61
 
#  Prime.instance.prime?(2)  #=> true
62
 
#  Prime.prime?(2)           #=> true
 
67
#   Prime.instance.prime?(2)  #=> true
 
68
#   Prime.prime?(2)           #=> true
63
69
#
64
70
# == Generators
 
71
#
65
72
# A "generator" provides an implementation of enumerating pseudo-prime
66
73
# numbers and it remembers the position of enumeration and upper bound.
67
74
# Futhermore, it is a external iterator of prime enumeration which is
80
87
#   is faster and uses much less memory than other generators. So,
81
88
#   it is suitable for factorizing an integer which is not large but
82
89
#   has many prime factors. e.g. for Prime#prime? .
 
90
 
83
91
class Prime
84
92
  include Enumerable
85
93
  @the_instance = Prime.new
105
113
  # Iterates the given block over all prime numbers.
106
114
  #
107
115
  # == Parameters
 
116
  #
108
117
  # +ubound+::
109
118
  #   Optional. An arbitrary positive number.
110
119
  #   The upper bound of enumeration. The method enumerates
113
122
  #   Optional. An implementation of pseudo-prime generator.
114
123
  #
115
124
  # == Return value
 
125
  #
116
126
  # An evaluated value of the given block at the last time.
117
127
  # Or an enumerator which is compatible to an +Enumerator+
118
128
  # if no block given.
119
129
  #
120
130
  # == Description
 
131
  #
121
132
  # Calls +block+ once for each prime number, passing the prime as
122
133
  # a parameter.
123
134
  #
126
137
  #   yields all prime numbers p <= +ubound+.
127
138
  #
128
139
  # == Note
 
140
  #
129
141
  # +Prime+.+new+ returns a object extended by +Prime+::+OldCompatibility+
130
142
  # in order to compatibility to Ruby 1.8, and +Prime+#each is overwritten
131
143
  # by +Prime+::+OldCompatibility+#+each+.
141
153
  # Returns true if +value+ is prime, false for a composite.
142
154
  #
143
155
  # == Parameters
 
156
  #
144
157
  # +value+:: an arbitrary integer to be checked.
145
158
  # +generator+:: optional. A pseudo-prime generator.
146
159
  def prime?(value, generator = Prime::Generator23.new)
161
174
  #        and a natural number -- an exponent.
162
175
  #
163
176
  # == Example
164
 
  # For [[p_1, e_1], [p_2, e_2], ...., [p_n, e_n]], it returns
165
 
  # p_1**e_1 * p_2**e_2 * .... * p_n**e_n.
166
 
  #
167
 
  #  Prime.int_from_prime_division([[2,2], [3,1]])  #=> 12
 
177
  # For <tt>[[p_1, e_1], [p_2, e_2], ...., [p_n, e_n]]</tt>, it returns:
 
178
  #
 
179
  #   p_1**e_1 * p_2**e_2 * .... * p_n**e_n.
 
180
  #
 
181
  #   Prime.int_from_prime_division([[2,2], [3,1]])  #=> 12
168
182
  def int_from_prime_division(pd)
169
183
    pd.inject(1){|value, (prime, index)|
170
184
      value *= prime**index
185
199
  # +ZeroDivisionError+:: when +value+ is zero.
186
200
  #
187
201
  # == Example
188
 
  # For an arbitrary integer
189
 
  # n = p_1**e_1 * p_2**e_2 * .... * p_n**e_n,
190
 
  # prime_division(n) returns
191
 
  # [[p_1, e_1], [p_2, e_2], ...., [p_n, e_n]].
192
 
  #
193
 
  #  Prime.prime_division(12) #=> [[2,2], [3,1]]
 
202
  # For an arbitrary integer:
 
203
  #
 
204
  #   n = p_1**e_1 * p_2**e_2 * .... * p_n**e_n,
 
205
  #
 
206
  # prime_division(n) returns:
 
207
  #
 
208
  #   [[p_1, e_1], [p_2, e_2], ...., [p_n, e_n]].
 
209
  #
 
210
  #   Prime.prime_division(12) #=> [[2,2], [3,1]]
194
211
  #
195
212
  def prime_division(value, generator= Prime::Generator23.new)
196
213
    raise ZeroDivisionError if value == 0
203
220
    for prime in generator
204
221
      count = 0
205
222
      while (value1, mod = value.divmod(prime)
206
 
             mod) == 0
207
 
        value = value1
208
 
        count += 1
 
223
             mod) == 0
 
224
        value = value1
 
225
        count += 1
209
226
      end
210
227
      if count != 0
211
 
        pv.push [prime, count]
 
228
        pv.push [prime, count]
212
229
      end
213
230
      break if value1 <= prime
214
231
    end
259
276
    def each(&block)
260
277
      return self.dup unless block
261
278
      if @ubound
262
 
        last_value = nil
263
 
        loop do
264
 
          prime = succ
265
 
          break last_value if prime > @ubound
266
 
          last_value = block.call(prime)
267
 
        end
 
279
        last_value = nil
 
280
        loop do
 
281
          prime = succ
 
282
          break last_value if prime > @ubound
 
283
          last_value = block.call(prime)
 
284
        end
268
285
      else
269
 
        loop do
270
 
          block.call(succ)
271
 
        end
 
286
        loop do
 
287
          block.call(succ)
 
288
        end
272
289
      end
273
290
    end
274
291
 
279
296
    def with_object(obj)
280
297
      return enum_for(:with_object) unless block_given?
281
298
      each do |prime|
282
 
        yield prime, obj
 
299
        yield prime, obj
283
300
      end
284
301
    end
285
302
  end
334
351
 
335
352
    def succ
336
353
      loop do
337
 
        if (@step)
338
 
          @prime += @step
339
 
          @step = 6 - @step
340
 
        else
341
 
          case @prime
342
 
          when 1; @prime = 2
343
 
          when 2; @prime = 3
344
 
          when 3; @prime = 5; @step = 2
345
 
          end
346
 
        end
347
 
        return @prime
 
354
        if (@step)
 
355
          @prime += @step
 
356
          @step = 6 - @step
 
357
        else
 
358
          case @prime
 
359
          when 1; @prime = 2
 
360
          when 2; @prime = 3
 
361
          when 3; @prime = 5; @step = 2
 
362
          end
 
363
        end
 
364
        return @prime
348
365
      end
349
366
    end
350
367
    alias next succ
353
370
    end
354
371
  end
355
372
 
356
 
 
357
 
 
358
 
 
359
373
  # Internal use. An implementation of prime table by trial division method.
360
374
  class TrialDivision
361
375
    include Singleton
385
399
    # +index+ is a 0-based index.
386
400
    def [](index)
387
401
      while index >= @primes.length
388
 
        # Only check for prime factors up to the square root of the potential primes,
389
 
        #   but without the performance hit of an actual square root calculation.
390
 
        if @next_to_check + 4 > @ulticheck_next_squared
391
 
          @ulticheck_index += 1
392
 
          @ulticheck_next_squared = @primes.at(@ulticheck_index + 1) ** 2
393
 
        end
394
 
        # Only check numbers congruent to one and five, modulo six. All others
 
402
        # Only check for prime factors up to the square root of the potential primes,
 
403
        #   but without the performance hit of an actual square root calculation.
 
404
        if @next_to_check + 4 > @ulticheck_next_squared
 
405
          @ulticheck_index += 1
 
406
          @ulticheck_next_squared = @primes.at(@ulticheck_index + 1) ** 2
 
407
        end
 
408
        # Only check numbers congruent to one and five, modulo six. All others
395
409
 
396
 
        #   are divisible by two or three.  This also allows us to skip checking against
397
 
        #   two and three.
398
 
        @primes.push @next_to_check if @primes[2..@ulticheck_index].find {|prime| @next_to_check % prime == 0 }.nil?
399
 
        @next_to_check += 4
400
 
        @primes.push @next_to_check if @primes[2..@ulticheck_index].find {|prime| @next_to_check % prime == 0 }.nil?
401
 
        @next_to_check += 2
 
410
        #   are divisible by two or three.  This also allows us to skip checking against
 
411
        #   two and three.
 
412
        @primes.push @next_to_check if @primes[2..@ulticheck_index].find {|prime| @next_to_check % prime == 0 }.nil?
 
413
        @next_to_check += 4
 
414
        @primes.push @next_to_check if @primes[2..@ulticheck_index].find {|prime| @next_to_check % prime == 0 }.nil?
 
415
        @next_to_check += 2
402
416
      end
403
417
      return @primes[index]
404
418
    end
428
442
      n = (n-1).div(2)*2+3 # the next odd number to given n
429
443
      table_index, integer_index, bit_index = indices(n)
430
444
      loop do
431
 
        extend_table until @tables.length > table_index
432
 
        for j in integer_index...ENTRIES_PER_TABLE
433
 
          if !@tables[table_index][j].zero?
434
 
            for k in bit_index...BITS_PER_ENTRY
435
 
              return NUMS_PER_TABLE*table_index + NUMS_PER_ENTRY*j + 2*k+1 if !@tables[table_index][j][k].zero?
436
 
            end
437
 
          end
438
 
          bit_index = 0
439
 
        end
440
 
        table_index += 1; integer_index = 0
 
445
        extend_table until @tables.length > table_index
 
446
        for j in integer_index...ENTRIES_PER_TABLE
 
447
          if !@tables[table_index][j].zero?
 
448
            for k in bit_index...BITS_PER_ENTRY
 
449
              return NUMS_PER_TABLE*table_index + NUMS_PER_ENTRY*j + 2*k+1 if !@tables[table_index][j][k].zero?
 
450
            end
 
451
          end
 
452
          bit_index = 0
 
453
        end
 
454
        table_index += 1; integer_index = 0
441
455
      end
442
456
    end
443
457
 
459
473
      ubound = lbound + NUMS_PER_TABLE
460
474
      new_table = [FILLED_ENTRY] * ENTRIES_PER_TABLE # which represents primarity in lbound...ubound
461
475
      (3..Integer(Math.sqrt(ubound))).step(2) do |p|
462
 
        i, j, k = indices(p)
463
 
        next if @tables[i][j][k].zero?
 
476
        i, j, k = indices(p)
 
477
        next if @tables[i][j][k].zero?
464
478
 
465
 
        start = (lbound.div(p)+1)*p  # least multiple of p which is >= lbound
466
 
        start += p if start.even?
467
 
        (start...ubound).step(2*p) do |n|
468
 
          _, j, k = indices(n)
469
 
          new_table[j] &= FILLED_ENTRY^(1<<k)
470
 
        end
 
479
        start = (lbound.div(p)+1)*p  # least multiple of p which is >= lbound
 
480
        start += p if start.even?
 
481
        (start...ubound).step(2*p) do |n|
 
482
          _, j, k = indices(n)
 
483
          new_table[j] &= FILLED_ENTRY^(1<<k)
 
484
        end
471
485
      end
472
486
      @tables << new_table.freeze
473
487
    end
483
497
 
484
498
    # Overwrites Prime#each.
485
499
    #
486
 
    # Iterates the given block over all prime numbers. Note that enumeration starts from
487
 
    # the current position of internal pointer, not rewound.
 
500
    # Iterates the given block over all prime numbers. Note that enumeration
 
501
    # starts from the current position of internal pointer, not rewound.
488
502
    def each(&block)
489
503
      return @generator.dup unless block_given?
490
504
      loop do
491
 
        yield succ
 
505
        yield succ
492
506
      end
493
507
    end
494
508
  end