December 2003 - cypherpunks-legacy

Re: Mastercard, Visa, Access, Barclaycard, Amex, JCB ...
by jamiel＠sybase.com 17 Dec '03

17 Dec '03

At 11:20 AM 07/07/94 -0400, Paul J. Ste. Marie wrote: >The credit card swipers send the entire transaction in for both >approval and reporting. I'd imagine that the account number, merchant >number, and transaction amount are all there, since that's all needed. >Whether or not they encrypt the data I don't know. I believe modern card readers for at least MC/Visa use some form of encryption, but for backwards compatibility the central offices also work unencrypted. In the hotel I used to work in, the card reader certainly didn't encrypt. -j -- "Blah Blah Blah" ___________________________________________________________________ Jamie Lawrence <jamiel(a)sybase.com>

2 1

Re: Any news on the crypto export bill?
by jmdaluz＠kquest.com 17 Dec '03

17 Dec '03

Apologies if this appeared ealier on this list. ----- Begin Forwarded Message ----- From: "Shabbir J. Safdar" <shabbir(a)panix.com> Subject: URGENT: House Rules committee postpones vote until week of 7/11/94 Date: Wed, 6 Jul 1994 08:43:58 -0400 (EDT) Content-Length: 20437 Precedence: list To: eff-activists(a)eff.org (eff-activists mailing list) [updated July 6, 1994 with Rules Committee information, shabbir] [My apologies for the errors in the fax numbers. I still am unable to get the GPO to sell me an electronic copy of the Congressional Directory, so I end up playing "catchup". If you can help, let me know. I'll pay for the directory, I just can't navigate the gpo's terrible interface.. -Shabbir] [POTENTIAL ACTION SCHEDULED FOR THE WEEK OF JULY 11TH, 1994] [PLEASE CHECK THE "WHAT YOU CAN DO RIGHT NOW" SECTION!] ********************************************************************* DISTRIBUTE WIDELY ********************************************************************* Table of contents: Introduction & Alert Status of the bill What you can do right now List of legislators supporting HR 3937 (formerly HR 3627) List of legislators wavering on HR 3937 (formerly HR 3627) List of legislators opposing HR 3937 (formerly HR 3627) What is the Cantwell bill? ------------------------------------------------------------------------------- INTRODUCTION & ALERT Voters Telecomm Watch keeps scorecards on legislators' positions on legislation that affects telecommunications and civil liberties. If you have updates to a legislator's positions, from either: -public testimony, -reply letters from the legislator, -stated positions from their office, please send them to vtw(a)panix.com so they can be added to this list. General questions: vtw(a)panix.com Mailing List Requests: vtw-list-request(a)panix.com Press Contact: stc(a)panix.com Gopher URL: gopher://gopher.panix.com:70/1/1/vtw WWW URL: Be patient; we're working on it. :-) ------------------------------------------------------------------------------- STATUS OF THE BILL (updated 7/6/94) The Cantwell bill, that allows for fewer restrictions exports of cryptography, has an interesting history. It was rolled into the General Export Administration Act HR 3937. The House Foreign Affairs Committee passed the full strength version out of committee after open, public hearings. The House Intelligence Committee took the bill and gutted it after a day of closed, secret hearings. The gutted version is making its way to the House floor. There is a crucial stop-off point, however. The House Rules Committee planned to hold a hearing on 6/30/94 to determine if the bill can be amended on the House floor (an "open" bill) or not (a "closed" bill). *** The vote was put off, because in one member's words, "all the legislators *** went home early". This gives us more time to lobby! If they mark the bill as "open", then the Cantwell bill could be restored to its previous version, removing the language put in by the House Select Intelligence Committee which gutted it without a public hearing. YOUR LOBBYING HAS ALREADY WORKED! A constituent reported back to me that their legislator (Rep. Quillen) had taken a position on marking the bill as "open". In a land where legislators are loathe to offend, this is terrific! Call and urge your legislator to take a position in favor of "open"! This may be the last thing you can do for the cryptographic export legislation. Take the time to make a call! Chronology of the bill Jul 11, 94 House Rules Comm. reconvenes; possibly votes again this week Jun 30, 94 [*** vote postponed, perhaps till the week of 7/11/94] House Rules Comm. decides whether to allow amendments on the bill when it reaches the House floor Jun 14, 94 Gutted by the House Select Committee on Intelligence May 20, 94 Referred to the House Select Committee on Intelligence May 18, 94 Passed out of the House Foreign Affairs Committee on May 18 attached to HR 3937, the General Export Administration Act Dec 6, 93 Referred to the Subcommittee on Economic Policy, Trade and Nov 22, 93 Referred to the House Committee on Foreign Affairs. ------------------------------------------------------------------------------- WHAT YOU CAN DO RIGHT NOW Estimated time to do this good deed: Two minutes Show your support for HR 3937 (formerly HR 3627) by contacting a member of the House Rules Committee and ask them to mark the bill as "open" (allowing amendments) when it reaches the House floor. The phone numbers of the members of the House Rules Committee are listed below. Please pick one from your state and call them. If your state isn't listed please call the Chairman, Rep. Joe Moakley. Feel free to use the following sample communique: The Honorable ____________ address Washington DC, 20515 Dear Congressman or Congresswoman, Please mark the General Export Administration Act (HR 3937) as "open" (allowing amendments on the House floor). Recently the House Intelligence Committee removed several provisions of the General Export Administration Act, HR 3937, dealing with the export of cryptographic technology. The House Intelligence Committee did this in a closed, secret hearing which provided for no public input. The House Foreign Affairs Committee previously held an open hearing on the same issue and received a flood of people testifying in favor of the bill, which the committee then reported out in full. I urge you to allow the democratic process to take its course on the House floor and mark the bill as "open". Sincerely, _________________________________ Phone/Fax/Addresses of members of the House Rules Committee 103rd United States Congress, 1993-1994 All addresses are Washington, D.C. 20515 Dist ST Name, Address, and Party Phone Fax ==== == ======================== ============== ============== 9 MA Moakley, John Joseph (D) 1-202-225-8273 1-202-225-3984 235 Cannon UNSPECIFIED POSITION 3 SC Derrick, Butler (D) 1-202-225-5301 1-202-225-5383 221 Cannon UNSPECIFIED POSITION 24 CA Beilenson, Anthony (D) 1-202-225-5911 no reliable fax 2465 RHOB UNSPECIFIED POSITION 24 TX Frost, Martin (D) 1-202-225-3605 1-202-225-4951 2459 RHOB UNSPECIFIED POSITION 10 MI Bonior, David E. (D) 1-202-225-2106 1-202-226-1169 2207 RHOB UNSPECIFIED POSITION 3 OH Hall, Tony P. (D) 1-202-225-6465 1-202-225-9272 2264 RHOB UNSPECIFIED POSITION 5 MO Wheat, Alan (D) 1-202-225-4535 1-202-225-5990 2334 RHOB UNSPECIFIED POSITION 6 TN Gordon, Bart (R) 1-202-225-4231 1-202-225-6887 103 Cannon UNSPECIFIED POSITION 28 NY Slaughter, Louise M. (D) 1-202-225-3615 1-202-225-7822 2421 RHOB UNSPECIFIED POSITION 22 NY Solomon, Gerald B. (R) 1-202-225-5614 1-202-225-1168 2265 RHOB UNSPECIFIED POSITION 1 TN Quillen, James H. (R) 1-202-225-6356 1-202-225-7812 102 Cannon WILL VOTE FOR "OPEN" - CALL AND THANK HIM 28 CA Dreier, David (R) 1-202-225-2305 no reliable fax 411 Cannon UNSPECIFIED POSITION 14 FL Goss, Porter J. (R) 1-202-225-2536 1-202-225-6820 330 Cannon UNSPECIFIED POSITION ------------------------------------------------------------------------- LIST OF LEGISLATORS SUPPORTING CRYPTOGRAPHY EXPORT LEGISLATION The following legislators have formally registered support for cryptography export legislation. Call them with your cheers. All addresses are Washington, D.C. 20515 Dist ST Name, Address, and Party Phone Fax ==== == ======================== ============== ============== 1 WA Cantwell, Maria (D) 1-202-225-6311 1-202-225-2286 1520 LHOB HR 3627's sponsor; thank her for her work! 16 IL Manzullo, Donald (R) 1-202-225-5676 1-202-225-5284 506 Cannon Cosponsored HR 3627 on 11/22/93 3 UT Orton, William H. (D) 1-202-225-7751 1-202-226-1223 1122 LHOB Cosponsored HR 3627 on 03/22/94 3 OR Wyden, Ronald (D) 1-202-225-4811 1-202-225-8941 1111 LHOB Cosponsored HR 3627 on 03/22/94 16 CA Edwards, Donald (D) 1-202-225-3072 1-202-225-9460 2307 RHOB Cosponsored HR 3627 on 03/22/94 19 OH Fingerhut, Eric D. (D) 1-202-225-5731 1-202-225-9114 431 Cannon Cosponsored HR 3627 on 03/22/94 4 MA Frank, Barney (D) 1-202-225-5931 1-202-225-0182 2404 RHOB Cosponsored HR 3627 on 03/22/94 2 UT Shepherd, Karen (D) 1-202-225-3011 1-202-226-0354 414 Cannon Cosponsored HR 3627 on 03/22/94 3 WA Unsoeld, Jolene (D) 1-202-225-3536 1-202-225-9095 1527 LHOB Cosponsored HR 3627 on 03/22/94 19 FL Johnston II, Harry (D) 1-202-225-3001 1-202-225-8791 204 Cannon Cosponsored HR 3627 on 03/22/94 9 WA Kreidler, Mike (D) 1-202-225-8901 1-202-226-2361 1535 LHOB Cosponsored HR 3627 on 03/22/94 4 WA Inslee, Jay (D) 1-202-225-5816 1-202-226-1137 1431 LHOB Cosponsored HR 3627 on 03/22/94 7 WA McDermott, James A. (D) 1-202-225-3106 1-202-225-9212 1707 LHOB Cosponsored HR 3627 on 03/22/94 8 IN McCloskey, Frank (D) 1-202-225-4636 1-202-225-4688 306 Cannon Cosponsored HR 3627 on 03/22/94 14 CA Eshoo, Anna G. (D) 1-202-225-8104 1-202-225-8890 1505 LHOB Cosponsored HR 3627 on 03/22/94 10 NC Ballenger, Thomas C. (R) 1-202-225-2576 1-202-225-0316 2238 RHOB Cosponsored HR 3627 on 05/04/94 2 WA Swift, Al (D) 1-202-225-2605 1-202-225-2608 1502 LHOB Cosponsored HR 3627 on 05/04/94 ------------------------------------------------------------------------------- LIST OF LEGISLATORS WAVERING ON CRYPTOGRAPHY EXPORT LEGISLATION [Feel free to use the sample communique at the end of the FAQ when calling or writing a legislator.] 26 NY Hinchey, Maurice D. (D) 1-202-225-6335 1-202-226-0774 1313 LHOB Recently told a constituent that he is taking the Cantwell bill under consideration, but has "national security concerns" about allowing encryption to be exported outside the United States. 1 IA Leach, James (R) 1-202-225-6576 1-202-226-1278 2186 RHOB Has yet to answer a constituent letter with a stated position. 13 NY Molinari, Susan (D) 1-202-225-3371 1-202-226-1272 123 Cannon Has yet to answer a constituent letter with a stated position. (has taken inordinately long) 8 NY Nadler, Jerrold (D) 1-202-225-5635 1-202-225-6923 424 Cannon Met with lobbying constituent in April '94; no position taken yet 25 CA McKeon, Howard P. (R) 1-202-225-1956 1-202-226-0683 307 Cannon Responded to a constituent with a "non-position", May '94 Had a favorable meeting with a constituent and a VTW volunteer in May '94. ------------------------------------------------------------------------------- LIST OF LEGISLATORS OPPOSING CRYPTOGRAPHY EXPORT LEGISLATION [Feel free to use the sample communique at the end of the FAQ when calling or writing a legislator.] Dist ST Name, Address, and Party Phone Fax ==== == ======================== ============== ============== 5 AL Cramer Jr, Robert E. (D) 1-202-225-4801 1-202-225-4392 1318 LHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 8 CA Pelosi, Nancy (D) 1-202-225-4965 1-202-225-8259 240 Cannon FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 32 CA Dixon, Julian C. (D) 1-202-225-7084 1-202-225-4091 2400 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 40 CA Lewis, Jerry (R) 1-202-225-5861 1-202-225-6498 2312 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 46 CA Dornan, Robert K. (R) 1-202-225-2965 no reliable fax 2402 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 2 CO Skaggs, David E. (D) 1-202-225-2161 1-202-225-9127 1124 LHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 10 FL Young, C. W. (R) 1-202-225-5961 1-202-225-9764 2407 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 4 KS Glickman, Daniel (D) 1-202-225-6216 1-202-225-5398 2371 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 1 NE Bereuter, Douglas (R) 1-202-225-4806 1-202-226-1148 2348 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 9 NJ Torricelli, Robert (D) 1-202-224-5061 1-202-225-0843 2159 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 3 NM Richardson, William (D) 1-202-225-6190 no reliable fax 2349 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 1 NV Bilbray, James H. (D) 1-202-225-5965 1-202-225-8808 2431 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 17 PA Gekas, George W. (R) 1-202-225-4315 1-202-225-8440 2410 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 2 RI Reed, John F. (D) 1-202-225-2735 1-202-225-9580 1510 LHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 14 TX Laughlin, Gregory H. (D) 1-202-225-2831 1-202-225-1108 236 Cannon FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 16 TX Coleman, Ronald D. (D) 1-202-225-4831 None 440 Cannon FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 19 TX Combest, Larry (R) 1-202-225-4005 1-202-225-9615 1511 LHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 1 UT Hansen, James V. (R) 1-202-225-0453 1-202-225-5857 2466 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. 6 WA Dicks, Norman D. (D) 1-202-225-5916 1-202-226-1176 2467 RHOB FAILED Cryptography exports: Voted to kill Rep. Cantwell's export provisions in the House Intelligence Committee on 6/15/94. ------------------------------------------------------------------------------- What is the Cantwell bill? The Cantwell bill would permit companies to export products with encryption technology in them. US companies are currently not permitted to export products (hardware or software) with this technology in them. What is encryption technology? Encryption technology, or cryptography, is the art of scrambling a conversation so that only the people communicating can decode it. Other people (such as eavesdroppers) cannot learn about the conversation. Where is cryptography being used? Cryptography is used to encrypt electronic mail to protect its confidentiality in transit. It's used by bank automatic teller machines to protect sensitive data (such as your account number, your Personal Identification Number, and your bank balance). It can be implemented into software (such as electronic mail programs and word processors) as well as hardware (such as telephones and "walkie-talkies") to ensure your privacy. Why is there a restriction on exporting products with technology in them? For many years the United States was a leading researcher in cryptography. High quality cryptographic technology was available only within the United States. The US government thought that if they did not let this technology be exported, foreign individuals would not be able to obtain it and use it against us (by keeping US intelligence agencies from eavesdropping on their communications) Since then, cryptography research has been published in international journals. Companies have been created throughout the world who export cryptographic technology from countries that do not have these restrictions. You can now buy the same, high-quality cryptographic technology from many international firms. Although the marketplace has changed, the regulations have not. Why should the regulations be changed? US companies compete in a global marketplace. Because of the export regulations, they often compete alongside products with superior cryptographic capabilities built into them. The result is that US companies build their products with an inferior encryption technology. The result of this is that you, as an American consumer, have great difficulty obtaining products with strong encryption in them. Because US products cannot compete against products with better privacy features, and because the laws are outdated, the regulations should be changed. The Cantwell bill fixes these regulations to more accurately resemble the current situation of the world marketplace. How can I help encourage more privacy-enhanced products and pass the Cantwell bill? Call or write your representative and ask them to support or cosponsor Rep. Cantwell's export provisions (formerly HR 3627) in the General Export Administration Act, HR 3937. You can base your letter on the sample communication below. SAMPLE LETTER OR PHONE CALL The Honorable ____________ address Washington DC, 20515 Dear Congressman or Congresswoman, As a citizen concerned for my privacy, as well as a supporter of American business, I urge you to cosponsor the Rep. Cantwell's cryptographic export provisions (formerly HR 3627) in the General Export Administration Act, HR 3937. The bill would allow US companies to produce and export products with cryptographic privacy-enhancing technology in them. These products are already available from firms throughout the world. US companies lose nearly $100 million per year in exports to them. By encouraging this industry, ordinary citizens like you and me would be able to purchase products with better privacy features. Please support or co-sponsor HR 3937. Sincerely, ___________________________________ ------------------------------------------------------------------------------- ----- End Forwarded Message ----- Jose M. daLuz KnowledgeQuest Online Research jmdaluz(a)kquest.com (508) 996-6101 (vox) (508) 996-6215 (fax)

1 0

LAW STUDENTS?
by Sandy Sandfort 17 Dec '03

17 Dec '03

C'punks, If anyone on the Cypherpunks Channel is *currently* a law student, please contact me via private e-mail. I have an favor to ask of you. Thanks, S a n d y

1 0

No Subject
by Peter Capek (TL-863-6721) 17 Dec '03

17 Dec '03

Pardon me if this has already been mentioned here -- I haven't noticed it -- but an article of considerable interest to this group has just been published.. The May 1994 issue of the IBM Journal of Research and Development contains "The Data Encryption Standard (DES) and its strength against attacks" by Don Coppersmith of IBM Research, who participated in the design and testing of DES, particularly in the design of the famous S-boxes. Included in the paper is a discussion of differential cryptanalysis and what was done to thwart that kind of attack. Peter Capek

1 0

(fwd) Random Numbers - Request for feedback
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.stat.math,sci.math,sci.math.num-analysis Path: bga.com!news.sprintlink.net!news.onramp.net!convex!cs.utexas.edu!swrinde!ihnp4.ucsd.edu!agate!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Random Numbers - Request for feedback Message-ID: <deleydCsIB28.KEI(a)netcom.com> Organization: NETCOM On-line Communication Services (408 261-4700 guest) Date: Wed, 6 Jul 1994 06:51:43 GMT Lines: 43 Xref: bga.com sci.stat.math:1315 sci.math:15353 sci.math.num-analysis:3354 RE: Computer Generated Random Numbers A few closing comments and requests for further information: 1. All my tests on random number generators were performed on VAX/VMS computers. VAX uses a 32-bit architecture, so the random number generators I tested were ones which used a word size of 32 bits or less. I would be interested in anybody's test results of a random number generator utilizing a larger word size, such as xrand() using SIZE=63. 2. Anyone know of some good references on primitive polynomials mod 2 and their applications? They're used in additive congruential random number generators like the xrand() one tested here. They're also used by file transfer programs such as xmodem to insure error free transmission, and they're used in cryptography too. Anyone know of a good book on Abstract Algebra? (The ones I have just briefly touch the topic and then move on.) 3. Resolution: Usually the random number generator is set up to return a floating point value between 0 and 1. A typical floating point variable R can only represent a finite number of different values between 0 and 1. If you magnify the result too much the discreetness of the floating point datum will become obvious. For example, in VAX architecture the F-floating datum has a precision of approximately one part in 2**23. Multiplying R by a very large number N to create a random variable between 0 and N will fail if N is too large because some of the values between 0 and N have no corresponding R value which maps to them (i.e. the mapping is no longer a surjection or onto map). For an F_floating datum, N above 2**23 is obviously too large. But even below 2**23 there's still a problem of some bins having 2 R values which map to them while other bins have only 1. We need to get N small enough so that the number of R values which maps to any bin is about the same, close enough so that differences aren't noticed when we test the random number generator. -David Deley deleyd(a)netcom.com

1 0

(fwd) Random Numbers - Results of testing BSD random()
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.stat.math,sci.math,sci.math.num-analysis Path: bga.com!news.sprintlink.net!news.onramp.net!convex!cs.utexas.edu!swrinde!ihnp4.ucsd.edu!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Random Numbers - Results of testing BSD random() Message-ID: <deleydCsIAyv.K6n(a)netcom.com> Organization: NETCOM On-line Communication Services (408 261-4700 guest) Date: Wed, 6 Jul 1994 06:49:42 GMT Lines: 119 Xref: bga.com sci.stat.math:1314 sci.math:15352 sci.math.num-analysis:3353 BSD random() Here are the partial results. Further tests were not performed due to lack of time. So far the generator appears to be comparable to a shuffled linear congruential generator. DEFINITION: Generating polynomial: x^31 + x^3 + 1 (primitive polynomial) Initialize circular queue of 31 elements using ANSI C linear congruential generator. Recursion formula: a[i] = a[i] + a[i-3] RATING: 1-D FAILS above 800,000 bpd (bins per dimension) 2-D FAILS above 3000 bpd 3-D FAILS above 210 bpd 4-D PASSES at 50 bpd (highest tested so far) 5-D not tested 6-D not tested 7-D not tested 8-D not tested This is an additive congruential type random number generator. An array table[31] is initially filled with random numbers using the ANSI C linear congruential random number generator. Random numbers are then generated using the recursion formula: table[k] = (table[k-31] + table[k-3]) mod 32 (Note that x**31 + x**3 + 1 is a primitive polynomial mod 2, which is being used here as a generator.) Since we are using the array table[] as a circular queue with 31 elements then table[k-31] is just table[k] before it gets replaced with the new value. The recursion formula becomes: table[k] = table[k] + table[k-3] The generator works well in practice. Knuth claims the sequence will have period 2**31 - 1. Knuth also claims there is very little theory to prove that this generator does or does not have desirable random properties. I would be interested if anyone knows of any recent developments in this area. -David Deley deleyd(a)netcom.com (So sorry, I lost the name of the original person who posted this code below which was used in the tests. -D.D.) /*** Code to implement random() & srandom() of BSD Unix. It was taken (though coded somewhat differently) from the Gnu BSD implementation. ***/ #include <stdio.h> #include <stdlib.h> #define LONG31 #ifdef LONG31 /* x^31 + x^3 + 1 */ #define SIZE 31 #define SIZE1 30 #define P1 3 #define P2 0 #else /* LONG63: x^63 + x + 1 */ #define SIZE 63 #define SIZE1 62 #define P1 1 #define P2 0 #endif #define LONG_MAX 0x7fffffff int p1=P1, p2=P2; long table[SIZE]; /*** return a "random" number in range [0, LONG_MAX] */ long xrand () { int r; table[p1] = table[p1] + table[p2]; /* add two table elements */ r = (table[p1] >> 1) & LONG_MAX; /* throw least significant bit away */ if (p1 == SIZE1) { /* increment the table indexes */ p1 = 0; p2 = p2 + 1; } else if (p2 == SIZE1) { p1 = p1 + 1; p2 = 0; } else { p1 = p1 + 1; p2 = p2 + 1; } return (r); } /*** use a linear congruential type generator to seed the state table & cycle the entire table 10 times */ void sxrand (seed) long seed; { int i; table[0] = seed; for (i=1; i<SIZE; ++i) table[i] = (table[i-1] * 1103515145) + 12345; /* lousy */ for (i=0; i<10*SIZE; ++i) (void) xrand(); }

1 0

(fwd) Random Numbers - CHIKSN.FOR
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.stat.math,sci.math,sci.math.num-analysis Path: bga.com!news.sprintlink.net!news.onramp.net!convex!cs.utexas.edu!swrinde!ihnp4.ucsd.edu!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Random Numbers - CHIKSN.FOR Message-ID: <deleydCsIBAp.KuF(a)netcom.com> Organization: NETCOM On-line Communication Services (408 261-4700 guest) Date: Wed, 6 Jul 1994 06:56:49 GMT Lines: 526 Xref: bga.com sci.stat.math:1316 sci.math:15354 sci.math.num-analysis:3355 {Approx. 520 lines} C CHIKSN.FOR C C This is the program which impliments the chi-square test used to test C random number generators. Presented here if you would wish to play C with it yourself, maybe do some testing of your own. See the paper C "Computer Generated Random Numbers" sections 4 and 6 for an explanation C of what this program does. C C This is not polished code you put on a shelf and admire, this is code C you dig your hands into and work with to make it do what you want it to. C C The main routine is the one to tinker with. This code is meant to be C modified to suit your needs. The goal is to fill up the bins with C balls. Here's a brief outline: C C 1. ASK USER INPUT: C a. Number of dimensions? NDIM C b. Number of bins per dimension? NBINSPD C c. Total number of balls to throw at bins? NBALLS C d. Number of tests to run? NCHITESTS C e. Random number generator to use (if more than one defined) C f. SEED value to initialize generator with C C 2. CREATE ARRAY BINS and ZERO ARRAY C C 3. THROW THE BALLS AT THE BINS and CALCULATE PROBABILITY: C C LOOP (do CHITEST=1 to NCHITESTS) C zero BIN array C LOOP (do J=1 to NBALLS) C get random numbers r1,r2,...,rn) C increment BIN(r1,r2,...,rn) C ENDLOOP C CALL CHSONE to calculate chi-square probability C ENDLOOP C CALL KSONE to calculate Kolmogorov-Smirnov probability C C The main routine here is where you put a call to your random number C generator, or for speed you can attempt a direct implimentation of your C random number generator to save the overhead of a call. (It can make a C difference when you call the random number generator 100 million times C in one test). C C The main routine defines a large one-dimensional array called BINS, the C maximum size of which would depend on your account quotas and machine C specific limitations. The array BINS keeps track of how many balls have C fallen into each bin. The size of array BINS determines the maximum C number of bins a user may select for a test. (10,000,000 bins is a C typical number you may want to use if possible.) C C So steps are: C C 1. Check definition of one-dimensional array NBINS C that it's not too large for your account quotas C or system limitations. C C 2. Place your random number generator to be tested where it C bluntly says "PLACE YOUR RANDOM NUMBER GENERATOR HERE". C The output is an integer IRANDOM between 0 and NBINS-1 C (NBINS is the number of bins chosen by the user). C C Currently the program is set up to use subroutine RAN1, a portable C random number generator from the book "NUMERICAL RECIPES: The Art C of Scientific Computing". I've had trouble on our UNIX system C not making array R(97) static even though the code says to. C Compiling with the -static qualifier works. C C The random number generator being tested is used to "randomly" select a C bin for the ball to fall in, and the counter for that bin is C incremented. Note for a multi-dimensional test we calculate the C appropriate index into the linear array BINS by hand. After all the C balls are thrown we call the subroutines to do the heavy math. C C All the subroutines should be fairly standard FORTRAN-77 modified C versions of routines from the book "NUMERICAL RECIPES: The Art of C Scientific Computing" by William H. Press, Brian P. Flannery, Saul A. C Teukolsky, and William T. Vetterling, and you should look there for C further reference as to what the routines are doing. (Note: the book C comes in several programming language forms including C, PASCAL, BASIC, C as well as FORTRAN, so you can take your pick and rewrite this code in C any language you please.) C C Note: C CHIKSN.FOR is currently set up to be run by a process with a very C large page file quota (pgflquo). If you get a 'exceed quota' error C attempting to run this then all you need to do is change the line C which reads: C C INTEGER*2 BINS(20 000 000) !The bins. C C to something smaller like: C C INTEGER*2 BINS(1 000 000) !The bins. C C C To compile: C $ FORTRAN CHIKSN C $ LINK CHIKSN C or C % f77 chiksn.f !(Some UNIX F77 compilers require -save option) C C Sample run: C Test MTH$RANDOM in 3-D with 10 bins per dimension and 10 balls per bin: C C $ RUN CHIKSN C Input number of dimensions NDIM: 3 C Input number of bins per dimension NBINSPD: 10 C Total number of bins = NBINSPD**NDIM = 1000 C Minimum number of balls = 5*NBINS = 5000 C Input total number of balls NBALLS: 10000 C Input number of Chi-Square tests NCHITESTS (min=2) : 2 C Choose random number generator to test C (1) MTH$RANDOM, C (2) RANDU, C (3) ANSI C, C (4) Microsoft C C (5) Turbo Pascal C (6) DES C : 1 C Input starting SEED value: 1 C C BALLS= 10000 CHISQ= 993.0002441 PROB= 0.4524292 C BALLS= 10000 CHISQ= 974.0001831 PROB= 0.2915459 C KS D= 0.5475708 PROB= 0.5863269 C C----------------------------------------------------------------------- PROGRAM CHIKSN C Perform a CHI-SQUARE test on a sequence of sets of N random numbers C NDIM = number of dimensions C NBINSPD = number of bins per dimension C NBINS = total number of bins. NBINS = NBINSPD**NDIM C NBALLS = total number of balls. Should be at least 5*(NBINS**NDIM) C NCHITESTS = Number of chi-square tests to do. Must be 2 or more. C EBINS = Expected value for each bin. EBINS = NBALLS/NBINS C SEED = Initial seed value for random number generator C C Note 1: The maximum size of array NBINS may be determined by the users C page file quota (pgflquo in AUTHORIZE). Also, it is recommended C the user have a very large working set quota (wsquo,wsextent) C to reduce page faulting. This can greatly improve speed. C We use INTEGER*2 array here to save space. C C Note 2: The maximum number of Chi-Square tests that can be saved is C arbitrary (array SAVEPROB). The user may choose any value. C C Note 3: The user may choose any starting seed value. Some restrictions C may apply depending upon the particular random number generator C being used. For example, RANDU should always be started with C an odd value of SEED. MTH$RANDOM may be started with any value C of SEED. C C Note 4: The MTH$RANDOM generator is used by the VAX FORTRAN intrinsic C function RAN and the VAX BASIC function RND. It is defined as: C C SEED = 69069*SEED + 1 mod 2**32 C X = SEED/2**32 C C Note 5: The RANDU generator is obsolite due to very strong correlation C in 3d space. ( Prove to yourself using 65539 = 2**16 + 3 that C SEED[i+2] = 6*SEED[i+1] - 9*SEED[i] ). It is defined as: C C SEED = 65539*SEED mod 2**31 C X = SEED/2**31 C C The RANDU generator should be started with an odd value of SEED. C C Note 6: The C standard library function rand() is defined as: C C SEED = 1103515245*SEED + 12345 mod 2**32 C IX = SEED mod 2**31 C C This standard random number generator is defined in the book: C The C Programming Language C Brian W. Kernighan and Dennis M. Ritchie C Prentice Hall, 1978 C C The same generator is defined in the ANSI C version by the same C authors above, and the same generator is used in VAX C. C C Note 7: The Microsoft C version 4.0 library function rand() impliments C the following: C C SEED = 214013*SEED + 2531011 mod 2**32 C IX = bits 16-31 of SEED C C Note 8: The Turbo Pascal version 6.0 function impliments the following: C C SEED = 134775813*SEED + 1 mod 2**32 C IX = bits 16-32 of SEED C IMPLICIT NONE INTEGER NDIM !Number of dimensions INTEGER NBINS !Number of bins INTEGER NBINSPD !Number of bins per dimension INTEGER NBALLS !Number of random numbers per !chi-square test. INTEGER NCHITESTS !Number of chi-square tests to do C INTEGER*2 BINS(20 000 000) !The bins. (see note 1) INTEGER*2 BINS(200 000) !Less bins. (see note 1) REAL EBINS !Expected number of balls per bin REAL SAVEPROB(100) !Array to save results of !chi-square tests (see note 2) INTEGER*4 SEED(2) !Only SEED(1) result is ever used. INTEGER*2 W(4) !Seeds for RANDU EQUIVALENCE(SEED,W) !for RANDU COMMON / SEEDSTORE / SEED INTEGER I,J,K,MRANDO,NBYTES,CLEAR,CHITEST,INDEX,IRANDOM,NCLEAR REAL*4 FRANDOM, RRANDOM CHARACTER*8 TIMEBUF EQUIVALENCE (IRANDOM,FRANDOM) REAL FOR$IRAN !The RANDU random number generator REAL RAND !UNIX rand() INTEGER*4 xrand !BSD random() REAL RAN1 !test generator supplied REAL D INTEGER JISHFT,IRANDOM2,COUNT REAL RANDES !DES FUNCTION (not supplied) INTEGER KEY(2) REAL CHSQ,PROB !Chi-square value, !chi-square probability REAL*4 FNBINSPD !float(NBINSPD) REAL*4 TWO31F REAL*4 TWO16F REAL*4 TWO15F TWO31F = 2.0**31.0 TWO16F = 2.0**16.0 TWO15F = 2.0**15.0 C*DES KEY(1) = 12345 !Choose any number you want C*DES KEY(2) = 678901 !to initialize DES with C*DES CALL DES_INIT(KEY) !DES code not included. 104 FORMAT(' Input number of dimensions NDIM: ',$) 100 FORMAT(' Input number of bins per dimension NBINSPD: ',$) 105 FORMAT(' Total number of bins = NBINSPD**NDIM = ',I) 106 FORMAT(' Minimum number of balls = 5*NBINS = ',I) 101 FORMAT(' Input total number of balls NBALLS: ',$) 103 FORMAT(' Input number of Chi-Square tests NCHITESTS (min=2) : ',$) 102 FORMAT(' Choose random number generator to test'/, 1 ' /*(1)*/ xrand(),'/ 1 ' /*(2)*/ UNIX rand(),'/ 1 ' /*(3)*/ MTH$RANDOM,'/ 2 ' /*(4)*/ RANDU,'/ 3 ' /*(5)*/ ANSI C,'/ 4 ' /*(6)*/ Microsoft C'/ 5 ' /*(7)*/ Turbo Pascal'/ 7 ' /*(8)*/ DES'/ 8 ' (9) another random number generator (choose this one)'/ 6 ' : ',$) 107 FORMAT(' Input starting SEED value: ',$) 200 FORMAT(BN,I) C ***GET USER INPUT*** 10 WRITE(6,104) !Input number of dimensions READ(5,200) NDIM WRITE(6,100) !Input number of bins per dimension READ(5,200) NBINSPD FNBINSPD = FLOAT(NBINSPD) NBINS = NBINSPD**NDIM !Calculate total number of bins WRITE(6,105) NBINS !Total number of bins is... WRITE(6,106) 5*NBINS !Minimum number of balls is... WRITE(6,101) !Input total number of balls READ(5,200) NBALLS WRITE(6,103) !Input number of chi-square tests to do READ(5,200) NCHITESTS WRITE(6,102) !Choose random number generator to test READ(5,200) MRANDO WRITE(6,107) !Starting SEED value READ(5,200) SEED(1) SEED(2) = 1 !Used only if random number generator !uses bigger than 32 bits C INITIALIZE GENERATOR IF NEEDED C*XRAND CALL sxrand(SEED(1)) !Initialize xrand() CALL RAN1(-SEED(1)) !Initialize RAN1 generator C Calculate expected average number of balls for each bin EBINS = FLOAT(NBALLS)/FLOAT(NBINS) C CALL TIME(TIMEBUF) C WRITE(6,201) TIMEBUF C201 FORMAT(1X,A8) DO CHITEST=1,NCHITESTS C *** ZERO BIN ARRAY *** DO I=1,NBINS BINS(I) = 0 ENDDO C*VMS !Quickly set BINS(k) = 0, k=1,...NBINS C*VMS !Does the equivalent of above C*VMS !but a lot faster. C*VMS K = 1 C*VMS NBYTES = NBINS*2 !total number of bytes to zero C*VMS DO WHILE (NBYTES .GT. 0) C*VMS IF (NBYTES .LE. 65534) THEN !maximum number of bytes we can clear C*VMS NCLEAR = NBYTES !in one call to LIB$MOVC5 is 65535 C*VMS ELSE C*VMS NCLEAR = 65534 !max that LIB$MOVC3 can do in one call C*VMS ENDIF !(make nclear an even number so we can divide evenly by 2) C*VMS CALL LIB$MOVC5(0,0,0,NCLEAR,BINS(K)) !Clear a block of memory C*VMS NBYTES = NBYTES - NCLEAR !Number of bytes still left to clear C*VMS K = K + NCLEAR/2 !Number of bytes cleared so far + 1 C*VMS ENDDO C Main Loop DO J=1,NBALLS INDEX = 1 DO I=0,NDIM-1 C ***PLACE YOUR RANDOM NUMBER GENERATOR HERE*** C Set IRANDOM using whatever random number generator you choose C IRANDOM = integer between 0 and NBINS-1 c IF (MRANDO .EQ. 1) THEN c IRANDOM = INT( ( float( xrand() ) /TWO31F ) *FNBINSPD) c ELSEIF (MRANDO .EQ. 2) THEN c IRANDOM = INT( RAND(SEED(1)) *FNBINSPD) !UNIX rand() c ELSEIF (MRANDO .EQ. 3) THEN c IRANDOM = INT( RAN(SEED(1)) *FNBINSPD) !VMS mth$random c ELSEIF (MRANDO .EQ. 4) THEN c IRANDOM = INT( FOR$IRAN(W(2),W(1)) *FNBINSPD) !Infamous randu c ELSEIF (MRANDO .EQ. 5) THEN c CALL LIB$EMUL(1103515245,SEED,12345,SEED) !ANSI C c IRANDOM = SEED(1) .AND. '7FFFFFFF'X c IRANDOM = INT( FLOAT(IRANDOM)/(TWO31F) *FNBINSPD) c ELSEIF (MRANDO .EQ. 6) THEN c CALL LIB$EMUL(214013,SEED,2531011,SEED) !Microsoft C 4.0 c IRANDOM = W(2) .AND. '7FFF'X c IRANDOM = INT( FLOAT(IRANDOM)/(TWO15F) *FNBINSPD) c ELSEIF (MRANDO .EQ. 7) THEN c CALL LIB$EMUL(134775813,SEED,1,SEED) !Turbo Pascal 6.0 c IRANDOM = SEED(1) .AND. 'FFFF0000'X c IRANDOM = JISHFT(IRANDOM,-16) c IRANDOM = INT( FLOAT(IRANDOM)/(TWO16F) * FNBINSPD) c ELSEIF (MRANDO .EQ. 8) THEN c IRANDOM = INT( RANDES() * FNBINSPD ) !DES (not supplied) c ELSEIF (MRANDO .EQ. 9) THEN IRANDOM = INT( RAN1(SEED(1)) * FNBINSPD ) c ENDIF C Calculate index by hand. INDEX = INDEX + IRANDOM*(NBINSPD**I) ENDDO BINS(INDEX) = BINS(INDEX) + 1 !ball fell in this bin C IF ( MOD(J, 1 000 000) .EQ. 0 ) THEN C CALL TIME(TIMEBUF) C WRITE(6,302) J, TIMEBUF 302 FORMAT(1X,'AT BALL:',I,3X,'TIME=',A8) C WRITE(6,303) SEED(2), SEED(1) 303 FORMAT(1X,'HEX: SEED(2)= ',Z,' SEED(1)= ',Z) C WRITE(6,304) SEED(2), SEED(1) 304 FORMAT(1X,'DEC: SEED(2)= ',I,' SEED(1)= ',I) C ENDIF ENDDO 400 CALL CHSONE(BINS,EBINS,NBINS,CHSQ,PROB) SAVEPROB(CHITEST) = PROB WRITE(6,1) NBALLS,CHSQ,PROB 1 FORMAT(' BALLS=',I,' CHISQ=',F,' PROB=',F) ENDDO C Now see if all the chi-square values are chi-square distributed: IF (NCHITESTS .GT. 1) THEN CALL KSONE(SAVEPROB,NCHITESTS,D,PROB) WRITE(6,2) D,PROB 2 FORMAT(1X,'KS D=',F,' PROB=',F) ENDIF END C============================================================================ C From book NUMERICAL RECIPES: The Art of Scientific Computing C Here for demonstration purposes C Replace this with whatever random number generator you want to test C Initialize with negative number FUNCTION RAN1(IDUM) REAL R(97) SAVE R !(Some UNIX F77 compilers require -save option on compile) PARAMETER (M1=259200,IA1=7141,IC1=54773,RM1=3.8580247E-6) PARAMETER (M2=134456,IA2=8121,IC2=28411,RM2=7.4373773E-6) PARAMETER (M3=243000,IA3=4561,IC3=51349) DATA IFF /0/ IF (IDUM.LT.0.OR.IFF.EQ.0) THEN IFF=1 IX1=MOD(IC1-IDUM,M1) IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IX1,M2) IX1=MOD(IA1*IX1+IC1,M1) IX3=MOD(IX1,M3) DO 11 J=1,97 IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IA2*IX2+IC2,M2) R(J)=(FLOAT(IX1)+FLOAT(IX2)*RM2)*RM1 11 CONTINUE IDUM=1 ENDIF IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IA2*IX2+IC2,M2) IX3=MOD(IA3*IX3+IC3,M3) J=1+(97*IX3)/M3 IF(J.GT.97.OR.J.LT.1)PAUSE write(1,100) R write(1,102) R(J) 100 format(f) 102 format(1x,'RAN1 = ', F) RAN1=R(J) R(J)=(FLOAT(IX1)+FLOAT(IX2)*RM2)*RM1 RETURN END C---------------------------------------------------------------------------- C CALCULATE THE CHI-SQUARE PROBABILITY. SINCE NBINS IS LARGE, IT IS JUST C THE CUMULATIVE GAUSSIAN DISTRIBUTION AFTER WE NORMALIZE THE VARIABLES. C OR ERROR FUNCTION. FUNCTION CHIPROB(NBINS,CHISQ) C Formula is the inverse of one given in Knuth for going the other way. INTEGER NBINS,DF REAL*4 CHISQ,Z DF = NBINS-1 Z = ( SQRT(24.0*CHISQ - 6.0*DF + 16.0) - 3*SQRT(2.0*DF) ) / 4.0 CHIPROB = ERF(Z) RETURN END FUNCTION ERF(X) C Return approximation to the complimentary error function erfc(X). C Return is not normalized err function. See book for details. C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C Modified to return normalized error function erf(X) C (It's a polynomial approximation) REAL ERFCC,Z,T Z=ABS(X/1.414213) !Normalize T=1./(1.+0.5*Z) ERFCC=T*EXP(-Z*Z-1.26551223+T*(1.00002368+T*(.37409196+ * T*(.09678418+T*(-.18628806+T*(.27886807+T*(-1.13520398+ * T*(1.48851587+T*(-.82215223+T*.17087277))))))))) IF (X.LT.0.) ERFCC=2.-ERFCC ERF = 1.0 - ERFCC/2.0 !Normalize and compliment RETURN END C---------------------------------------------------------------------------- C THE FOLLOWING SUBROUTINES CALCULATE THE CHI-SQUARE VALUE: SUBROUTINE CHSONE(BINS,EBINS,NBINS,CHSQ,PROB) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing INTEGER NBINS INTEGER*2 BINS(NBINS) REAL EBINS,CHSQ,PROB CHSQ=0. IF(EBINS.LE.0.) PAUSE 'CHSONE: EBINS must be > 0' DO 11 J=1,NBINS CHSQ=CHSQ+(BINS(J)-EBINS)**2/EBINS 11 CONTINUE PROB=CHIPROB(NBINS,CHSQ) RETURN END C============================================================================ C THE FOLLOWING SUBROUTINES CALCULATE THE KOLMOGOROV-SMIRNOV PROBABILITY SUBROUTINE KSONE(DATA,N,D,PROB) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C DF - degrees of freedom. Passsed to FUNC INTEGER N REAL DATA(N) REAL D,PROB CALL PIKSRT(N,DATA) EN=N D=0. FO=0. DO 11 J=1,N FN=J/EN FF=DATA(J) DT=AMAX1(ABS(FO-FF),ABS(FN-FF)) IF(DT.GT.D)D=DT FO=FN 11 CONTINUE PROB=PROBKS(SQRT(EN)*D) RETURN END C---------------------------------------------------------------------------- FUNCTION PROBKS(ALAM) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C Note the routine in the Numerical Recipes book erronously returns C 1 instead of 0 for large values of ALAM. PARAMETER (EPS1=0.001, EPS2=1.E-8) A2=-2.*ALAM**2 FAC=2. PROBKS=0. TERMBF=0. DO 11 J=1,100 TERM=FAC*EXP(A2*J**2) PROBKS=PROBKS+TERM C Error in Numerical Recipes book. Terminate if TERM underflows. C** IF(ABS(TERM).LT.EPS1*TERMBF.OR.ABS(TERM).LT.EPS2*PROBKS)RETURN IF(ABS(TERM).LE.EPS1*TERMBF.OR.ABS(TERM).LE.EPS2*PROBKS)RETURN FAC=-FAC TERMBF=ABS(TERM) 11 CONTINUE PROBKS=1.0 RETURN END C---------------------------------------------------------------------------- SUBROUTINE PIKSRT(N,ARR) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C See book for details. INTEGER N REAL ARR(N) DO 12 J=2,N A=ARR(J) DO 11 I=J-1,1,-1 IF(ARR(I).LE.A)GO TO 10 ARR(I+1)=ARR(I) 11 CONTINUE I=0 10 ARR(I+1)=A 12 CONTINUE RETURN END

2 1

(fwd) Random Numbers - CORELA.FOR
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.stat.math,sci.math,sci.math.num-analysis Path: bga.com!news.sprintlink.net!news.onramp.net!convex!cs.utexas.edu!swrinde!ihnp4.ucsd.edu!agate!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Random Numbers - CORELA.FOR Message-ID: <deleydCsIBCu.KzB(a)netcom.com> Organization: NETCOM On-line Communication Services (408 261-4700 guest) Date: Wed, 6 Jul 1994 06:58:06 GMT Lines: 211 Xref: bga.com sci.stat.math:1317 sci.math:15355 sci.math.num-analysis:3356 {Approx. 200 lines} PROGRAM CORELA C Perform a KS test comparing the first 100 elements of a random C number generator, starting with SEED values of 1..10 C C This is not polished code you put on a shelf and admire, this is code C you dig your hands into and work with to make it do what you want it to. C C Place your random number generator to be tested where it C bluntly says "PLACE YOUR RANDOM NUMBER GENERATOR HERE". C The output is a floating point between 0 (inclusive) and 1 (exclusive). C C Currently the program is set up to use subroutine RAN1, a portable C random number generator from the book "NUMERICAL RECIPES: The Art C of Scientific Computing". I've had trouble on our UNIX system C not making array R(97) static even though the code says to. C Compiling with the -static qualifier works. C IMPLICIT NONE INTEGER SINC,I,J REAL SEQ(100,10) REAL AR(10) REAL D,PROB INTEGER MRANDO, SEEDINIT, IRANDOM INTEGER*4 SEED(2) !Only SEED(1) result is ever used. INTEGER*2 W(4) !Seeds for RANDU EQUIVALENCE(SEED,W) !for RANDU COMMON / SEEDSTORE / SEED REAL*4 FRANDOM REAL FOR$IRAN !The RANDU random number generator REAL RAN1 !test generator supplied INTEGER JISHFT,IRANDOM2,COUNT REAL*4 FNBINSPD !float(NBINSPD) REAL*4 TWO31F REAL*4 TWO16F REAL*4 TWO15F TWO31F = 2.0**31.0 TWO16F = 2.0**16.0 TWO15F = 2.0**15.0 102 FORMAT(' Choose random number generator to test'/, 1 ' /*(1)*/ MTH$RANDOM,'/ 2 ' /*(2)*/ RANDU,'/ 3 ' /*(3)*/ ANSI C,'/ 4 ' /*(4)*/ Microsoft C'/ 5 ' /*(5)*/ Turbo Pascal'/ 8 ' (9) another random number generator (choose this one)'/ 6 ' : ',$) 107 FORMAT(' Input starting SEED value: ',$) 108 FORMAT(' Input increment between SEED values: ',$) 200 FORMAT(BN,I) 10 CONTINUE WRITE(6,102) !Choose random number generator to test READ(5,200) MRANDO WRITE(6,107) !Starting SEED value READ(5,200) SEED(1) SEEDINIT = SEED(1) SEED(2) = 1 WRITE(6,108) !INCREMENT VALUE READ(5,200) SINC C Main Loop DO J=1,10 !10 sequences DO I=1,100 !sequence length of first 100 numbers C ***PLACE YOUR RANDOM NUMBER GENERATOR HERE*** C Set FRANDOM using whatever random number generator you choose C to a floating point value in the range [0,1) FRANDOM = RAN1(SEED(1)) C IF (MRANDO .EQ. 1) THEN C FRANDOM = RAN(SEED(1)) !mth$random C ELSEIF (MRANDO .EQ. 2) THEN C FRANDOM = FOR$IRAN(W(2),W(1)) !randu C ELSEIF (MRANDO .EQ. 3) THEN C CALL LIB$EMUL(1103515245,SEED,12345,SEED) !VAX C C IRANDOM = SEED(1) .AND. '7FFFFFFF'X C FRANDOM = FLOAT(IRANDOM)/(TWO31F) C ELSEIF (MRANDO .EQ. 4) THEN C CALL LIB$EMUL(214013,SEED,2531011,SEED) !Microsoft C 4.0 C IRANDOM = W(2) .AND. '7FFF'X C FRANDOM = FLOAT(IRANDOM)/(TWO15F) C ELSEIF (MRANDO .EQ. 5) THEN C CALL LIB$EMUL(134775813,SEED,1,SEED) !Turbo Pascal 6.0 C IRANDOM = SEED(1) .AND. 'FFFF0000'X C IRANDOM = JISHFT(IRANDOM,-16) C FRANDOM = FLOAT(IRANDOM)/(TWO16F) C ENDIF SEQ(I,J) = FRANDOM ENDDO SEEDINIT = SEEDINIT + SINC !calculate new initial seed SEED(1) = SEEDINIT !set new initial seed ENDDO C Do a KS test on each edlement comparing the 10 sequences DO I=1,100 DO J=1,10 AR(J) = SEQ(I,J) !Transfer to short array ENDDO CALL KSONE(AR,10,D,PROB) WRITE(6,2) I,PROB 2 FORMAT(1X,'I=',I4,' KS PROB=',F) ENDDO END C============================================================================ C From book NUMERICAL RECIPES: The Art of Scientific Computing C Here for demonstration purposes C Replace this with whatever random number generator you want to test C Initialize with negative number FUNCTION RAN1(IDUM) DIMENSION R(97) SAVE R !(Some UNIX F77 compilers require -save option on compile) PARAMETER (M1=259200,IA1=7141,IC1=54773,RM1=3.8580247E-6) PARAMETER (M2=134456,IA2=8121,IC2=28411,RM2=7.4373773E-6) PARAMETER (M3=243000,IA3=4561,IC3=51349) DATA IFF /0/ IF (IDUM.LT.0.OR.IFF.EQ.0) THEN IFF=1 IX1=MOD(IC1-IDUM,M1) IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IX1,M2) IX1=MOD(IA1*IX1+IC1,M1) IX3=MOD(IX1,M3) DO 11 J=1,97 IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IA2*IX2+IC2,M2) R(J)=(FLOAT(IX1)+FLOAT(IX2)*RM2)*RM1 11 CONTINUE IDUM=1 ENDIF IX1=MOD(IA1*IX1+IC1,M1) IX2=MOD(IA2*IX2+IC2,M2) IX3=MOD(IA3*IX3+IC3,M3) J=1+(97*IX3)/M3 IF(J.GT.97.OR.J.LT.1)PAUSE RAN1=R(J) R(J)=(FLOAT(IX1)+FLOAT(IX2)*RM2)*RM1 RETURN END C============================================================================== C THE FOLLOWING SUBROUTINES CALCULATE THE KOLMOGOROV-SMIRNOV PROBABILITY SUBROUTINE KSONE(DATA,N,D,PROB) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C DF - degrees of freedom. Passsed to FUNC INTEGER N REAL DATA(N) REAL D,PROB CALL PIKSRT(N,DATA) EN=N D=0. FO=0. DO 11 J=1,N FN=J/EN FF=DATA(J) DT=AMAX1(ABS(FO-FF),ABS(FN-FF)) IF(DT.GT.D)D=DT FO=FN 11 CONTINUE PROB=PROBKS(SQRT(EN)*D) RETURN END C------------------------------------------------------------------------------ FUNCTION PROBKS(ALAM) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C Note the routine in the Numerical Recipes book erronously returns C 1 instead of 0 for large values of ALAM. PARAMETER (EPS1=0.001, EPS2=1.E-8) A2=-2.*ALAM**2 FAC=2. PROBKS=0. TERMBF=0. DO 11 J=1,100 TERM=FAC*EXP(A2*J**2) PROBKS=PROBKS+TERM C Error in Numerical Recipes book. Terminate if TERM underflows. C** IF(ABS(TERM).LT.EPS1*TERMBF.OR.ABS(TERM).LT.EPS2*PROBKS)RETURN IF(ABS(TERM).LE.EPS1*TERMBF.OR.ABS(TERM).LE.EPS2*PROBKS)RETURN FAC=-FAC TERMBF=ABS(TERM) 11 CONTINUE PROBKS=1.0 RETURN END C------------------------------------------------------------------------------ SUBROUTINE PIKSRT(N,ARR) C Adapted from book NUMERICAL RECIPES: The Art of Scientific Computing C See book for details. INTEGER N REAL ARR(N) DO 12 J=2,N A=ARR(J) DO 11 I=J-1,1,-1 IF(ARR(I).LE.A)GO TO 10 ARR(I+1)=ARR(I) 11 CONTINUE I=0 10 ARR(I+1)=A 12 CONTINUE RETURN END

1 0

(fwd) Random Numbers - SPECTRAL.FOR
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.stat.math,sci.math,sci.math.num-analysis Path: bga.com!news.sprintlink.net!news.onramp.net!convex!cs.utexas.edu!swrinde!ihnp4.ucsd.edu!agate!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Random Numbers - SPECTRAL.FOR Message-ID: <deleydCsIBEG.L2B(a)netcom.com> Organization: NETCOM On-line Communication Services (408 261-4700 guest) Date: Wed, 6 Jul 1994 06:59:04 GMT Lines: 426 Xref: bga.com sci.stat.math:1318 sci.math:15356 sci.math.num-analysis:3357 {Approx. 420 lines} PROGRAM SPECTRAL ! Performs the spectral test for a linear congruential random number generator. ! ! This program adapted from: ! ! ALGORYTHM AS 193 APPLIED STATISTICS, (1983) VOL. 32, NO.3 PG. 328-335 ! T. R. Hopkins ! Modified to run on VAX/VMS systems using REAL*16 (64 bit) variables. ! The original is FORTRAN-66 compliant. ! ! Consider linear congruential generators of the form: ! SEED = (A*SEED + C) mod M ! ! Given A, and M, the spectral test calculates NUSQ (NU**2), LOGNU (base 2), ! and MU. As a guide, Knuth suggests a multiplier A may be considered ! adequate if the values of MU returned by the spectral test are > 0.1 . ! For an exceptionally good multiplier, these values will all be greater ! than unity. ! ! The spectral test may be applied if: ! 1. The sequence has maximal period, or ! 2. M is prime and C = 0 and the period length is M-1, or ! 3. M = 2**e and A mod 8 = 5 or A mod 8 = 3. ! In this third case the spectral test is applied using ! A = A and M = 2**(e-2). For example, in analyzing RANDU, ! use A = 65539 and M = 536870912 (2**29). ! ! Further information on the spectral test is in: ! ! Knuth, Donald E. "The Art of Computer Programming Vol. 2: Seminumerical ! algorithms, 2nd edition. Reading, Mass.: Addison-Wesley. 1981 ! ! The value of parameter BIGT determines how many dimensions are calculated. ! Higher dimensions may be obtained by changing this parameter and recompiling. ! Note that 12 is about the highest feasible. Above 12 the program may take ! days to complete. C Example: C MTH$RANDOM is defined as C C SEED = (69069*SEED + 1) MOD 2**32 C C Here A = 69069 C and M = 2**32 = 4294967296 C C $ RUN SPECTRAL C INPUT A: 69069 C INPUT M: 4294967296 C C A= 69069.0 C M= 4294967296.0 C BIGT= 6 C NUSQ= C NUSQ ( 2)= 4243209856.000000 C NUSQ ( 3)= 2072544.000000 C NUSQ ( 4)= 52804.000000 C NUSQ ( 5)= 6990.000000 C NUSQ ( 6)= 242.000000 C LOGNU= C LOGNU( 2)= 15.991254 C LOGNU( 3)= 10.491486 C LOGNU( 4)= 7.844180 C LOGNU( 5)= 6.385538 C LOGNU( 6)= 3.959432 C MU= C MU= ( 2)= 3.103734 C MU= ( 3)= 2.909942 C MU= ( 4)= 3.203639 C MU= ( 5)= 5.006469 C MU= ( 6)= 0.017052 C C Now examine the MU values. All values are above 1 except the very last C value MU(6) is 0.01, indicating MTH$RANDOM may not perform as well in a C 6-D test. C C C Run spectral again this time trying the values for the bad RANDU generator: C C MTH$RANDOM is defined as C C SEED = (65539*SEED) MOD 2**31 C C Here A = 65539 C and M = 2**31 but we use M = 2**29 for reasons discussed above C C C $ RUN SPECTRAL C INPUT A: 65539 C INPUT M: 536870912 !(2**29) C C A= 65539.0 C M= 536870912.0 C BIGT= 6 C NUSQ= C NUSQ ( 2)= 536936458.000000 C NUSQ ( 3)= 118.000000 C NUSQ ( 4)= 116.000000 C NUSQ ( 5)= 116.000000 C NUSQ ( 6)= 116.000000 C LOGNU= C LOGNU( 2)= 14.500088 C LOGNU( 3)= 3.441322 C LOGNU( 4)= 3.428990 C LOGNU( 5)= 3.428990 C LOGNU( 6)= 3.428990 C MU= C MU= ( 2)= 3.141976 C MU= ( 3)= 0.000010 C MU= ( 4)= 0.000124 C MU= ( 5)= 0.001421 C MU= ( 6)= 0.015025 C C Notice here the MU values for dimensions 2 through 6 are all extremely C small. This generator does horribly on these dimensions. The spectral C test noticed it right away. PARAMETER BIGT = 6 !Number of dimensions to go up to. Max is 12. PARAMETER IU = BIGT !(Beyond 12 program may take days to run.) PARAMETER IV = BIGT INTEGER*4 IFAULT REAL*16 A, M, MU(BIGT), NUSQ(BIGT), LOGNU(BIGT), U(IU,BIGT), 2 V(IV,BIGT), Z(BIGT) 100 FORMAT(' INPUT A: ',$) 101 FORMAT(' INPUT M: ',$) 200 FORMAT(BN,G33.0) 201 WRITE(6,100) READ(5,200) A !MTH$RANDOM example: A = 69069.0 WRITE(6,101) READ(5,200) M !MTH$RANDOM example: M = 4294967296.0 (2**32) CALL SPECT(A,M,BIGT,MU,NUSQ,LOGNU,U,IU,V,IV,Z,IFAULT) IF (IFAULT .GT. 0) THEN IF (IFAULT .EQ. 1) THEN PRINT*, ' BIGT < 2' ELSEIF (IFAULT .EQ. 2) THEN PRINT*, ' A .GE. M .OR. A .LE. 0 .OR. M .LE. 0' ELSEIF (IFAULT .EQ. 3) THEN PRINT*, ' M > Mmax' ELSEIF (IFAULT .EQ. 4) THEN PRINT*, ' A and M not relatively prime' ELSEIF (IFAULT .EQ. 5) THEN PRINT*, ' Intermediate result > Mmax * Mmax' ELSE PRINT*, ' IFAULT .GT. 5' ENDIF STOP ENDIF WRITE(6,1) A WRITE(6,2) M WRITE(6,3) BIGT WRITE(6,41) DO I=2,BIGT WRITE(6,4) I,NUSQ(I) ENDDO WRITE(6,51) DO I=2,BIGT WRITE(6,5) I,LOGNU(I) ENDDO WRITE(6,61) DO I=2,BIGT WRITE(6,6) I,MU(I) ENDDO 1 FORMAT(' A=',F33.1) 2 FORMAT(' M=',F33.1) 3 FORMAT(' BIGT=',I) 41 FORMAT(' NUSQ=') 4 FORMAT(' NUSQ (',I,')=',F33.6) 51 FORMAT(' LOGNU=') 5 FORMAT(' LOGNU(',I,')=',F33.6) 61 FORMAT(' MU=') 6 FORMAT(' MU= (',I,')=',F33.6) C GOTO 201 END SUBROUTINE SPECT(A, M, BIGT, MU, NUSQ, LOGNU, U, IU, V, IV, Z, * IFAULT) C C ALGORYTHM AS 193 APPLIED STATISTICS, (1983) VOL. 32, NO.3 PG. 328-335 C T. R. Hopkins C C A REVISED ALGORITHM FOR THE SPECTRAL TEST C Modified to use REAL*16 variables for VAX/VMS C IMPLICIT NONE INTEGER*4 I, I2, J, K INTEGER*4 BIGT, IU, IV, T, T1, IFAULT REAL*16 A, M, MU(BIGT), NUSQ(BIGT), LOGNU(BIGT), * U(IU, BIGT), V(IV, BIGT), Z(BIGT), * H, HPRIME, MMAX, MMAX2, MSQ, P, PI, PPRIME, Q, * QTEMP, R, S, SIGN, UC, VC, VIJ, VJJ, W, ZERO, ONE, TWO, FOUR, * DINT, DNINT, VPROD DATA ZERO /0.0Q0/, ONE /1.0Q0/, TWO /2.0Q0/, FOUR /4.0Q0/ C C SUITABLE VALUES FOR C 1) IBM REAL*8 C DATA MMAX/33554432.0D0/ C 2) IBM REAL*16 C 3) CDC 7600 DOUBLE PRECISION C DATA MMAX/35184372088832.0D0/ C DATA MMAX /9007199254740992.0D0/ C C A VAX/VMS REAL*16 has precision approximately one part in 2**112 C Knuth claims values rarely if ever exceed M**2 C So Hopkins takes maxval = 8*m**2 and solves 2**112 = 8*m**2 for M C giving Mmax = 2**(112/2)/8 DATA MMAX /9.0Q15/ C C TEST THE VALIDITY OF THE INPUT PARAMETERS C MMAX2 = MMAX * MMAX IFAULT = 0 IF (BIGT .LT. 2) IFAULT = 1 IF (A .GE. M .OR. A .LE. ZERO .OR. M .LE. ZERO) IFAULT = 2 IF (M .GT. MMAX) IFAULT = 3 IF (IFAULT .GT. 0) RETURN C C CHECK A AND M ARE RELATIVELY PRIME C NEED VALID A AND M C USE EUCLIDS ALGORITHM C H = A HPRIME = M 10 R = QMOD(HPRIME, H) IF (R .EQ. ZERO) GOTO 20 HPRIME = H H = R GOTO 10 20 IF (H .NE. ONE) IFAULT = 4 ! A and M not relatively prime IF (IFAULT .NE. 0) RETURN MSQ = M * M C C ALL STEPS REFER TO THOSE IN KNUTHS ALGORITHM C STEP 1 - INITIALIZATION C H = A HPRIME = M P = ONE PPRIME = ZERO R = A S = ONE + A * A C C STEP 2 - EUCLIDEAN STEP C 30 Q = QINT(HPRIME / H) UC = HPRIME - Q * H VC = PPRIME - Q * P W = UC * UC + VC * VC IF (W .GE. S) GOTO 40 S = W HPRIME = H H = UC PPRIME = P P = VC GOTO 30 C C STEP 3 - COMPUTE NU(2) C 40 UC = UC - H VC = VC - P W = UC * UC + VC * VC IF (W .GE. S) GOTO 50 S = W HPRIME = UC PPRIME = VC 50 NUSQ(2) = S C C INITIALIZE U AND V MATRICES C NOTE WE STORE BY COLUMNS WHEREAS KNUTH STORES BY ROWS C T = 2 U(1, 1) = -H U(1, 2) = -HPRIME U(2, 1) = P U(2, 2) = PPRIME SIGN = ONE IF (PPRIME .GT. ZERO) SIGN = -ONE V(1, 1) = SIGN * PPRIME V(1, 2) = -SIGN * P V(2, 1) = SIGN * HPRIME V(2, 2) = -SIGN * H C C STEP 4 - ADVANCE T C 60 IF (T .EQ. BIGT) GOTO 200 T1 = T T = T + 1 R = QMOD(A * R, M) U(1, T) = -R U(T, T) = ONE U(T, 1) = ZERO V(1, T) = ZERO V(T, T) = M DO 70 I = 2, T1 U(I, T) = ZERO U(T, I) = ZERO V(I, T) = ZERO 70 CONTINUE DO 90 I = 1, T1 QTEMP = V(1, I) * R Q = QNINT(QTEMP / M) V(T, I) = QTEMP - Q * M DO 80 I2 = 1, T 80 U(I2, T) = U(I2, T) + Q * U(I2, I) 90 CONTINUE S = QMIN1(S, VPROD(U(1, T), U(1, T), T)) K = T J = 1 C C STEP 5 - TRANSFORM C 100 DO 120 I = 1, T IF (I .EQ. J) GOTO 120 VIJ = VPROD(V(1, I), V(1, J), T) VJJ = VPROD(V(1, J), V(1, J), T) IF (TWO * QABS(VIJ) .LE. VJJ) GOTO 120 Q = QNINT(VIJ / VJJ) DO 110 I2 = 1, T V(I2, I) = V(I2, I) - Q * V(I2, J) U(I2, J) = U(I2, J) + Q * U(I2, I) 110 CONTINUE K = J 120 CONTINUE C C STEP 6 - EXAMINE NEW BOUND C IF (K .EQ. J) S = QMIN1(S, VPROD(U(1, J), U(1, J), T)) C C STEP 7 - ADVANCE J C J = J + 1 IF (J .EQ. T + 1) J = 1 IF (J .NE. K) GOTO 100 C C STEP 8 - PREPARE FOR SEARCH C C MU AND LOGNU ARE USED TO STORE KNUTHS X AND Y RESPECTIVELY C DO 130 I = 1, T MU(I) = ZERO LOGNU(I) = ZERO QTEMP = VPROD(V(1, I), V(1, I), T) IF (QTEMP .GT. MMAX2) GOTO 240 !Intermediate result > Mmax * Mmax QTEMP = QTEMP / MSQ Z(I) = QINT(QSQRT(QINT(QTEMP * S))) 130 CONTINUE K = T C C STEP 9 - ADVANCE XK C 140 IF (MU(K) .EQ. Z(K)) GOTO 190 MU(K) = MU(K) + ONE DO 150 I = 1, T 150 LOGNU(I) = LOGNU(I) + U(I, K) C C STEP 10 - ADVANCE K C 160 K = K + 1 IF (K .GT. T) GOTO 180 MU(K) = -Z(K) DO 170 I = 1, T 170 LOGNU(I) = LOGNU(I) - TWO * Z(K) * U(I, K) GOTO 160 180 S = QMIN1(S, VPROD(LOGNU, LOGNU, T)) C C STEP 11 - DECREASE K C 190 K = K - 1 IF (K .GE. 1) GOTO 140 NUSQ(T) = S GOTO 60 C C CALCULATE NU AND LOG(NU) C 200 DO 210 I = 2, BIGT MU(I) = QSQRT(NUSQ(I)) LOGNU(I) = QLOG(MU(I)) / QLOG(TWO) 210 CONTINUE C C CALCULATE TRANSFORMED MU VALUES C PI = 3.14159 26535 89793 23846 26433 83279 50288 41971 69399 37511 Q = ONE DO 220 T = 2, BIGT, 2 Q = Q * PI * TWO / QEXT(T) MU(T) = Q * MU(T) ** T / M 220 CONTINUE IF (BIGT .EQ. 2) RETURN Q = TWO DO 230 T = 3, BIGT, 2 Q = Q * PI * TWO / QEXT(T) MU(T) = Q * MU(T) ** T / M 230 CONTINUE RETURN 240 IFAULT = 5 !Intermediate result > Mmax * Mmax RETURN END REAL*16 FUNCTION VPROD(U, V, T) C C ALGORYTHM AS 193 APPLIED STATISTICS, (1983) VOL. 32, NO.3 PG. 328-335 C C AUXILIARY FUNCTION TO CALCULATE THE INNER PRODUCT OF C THE TWO VECTORS U AND V OF LENGTH T. C Modified to REAL*16 C INTEGER T REAL*16 U(T), V(T), SUM, ZERO DATA ZERO /0.0Q0/ C SUM = ZERO DO 10 I = 1, T 10 SUM = SUM + U(I) * V(I) VPROD = SUM RETURN END

1 0

(fwd) Re: BSD random() - any good (source included)?
by Jim choate 17 Dec '03

17 Dec '03

Newsgroups: sci.math,sci.stat.math Path: bga.com!news.sprintlink.net!hookup!europa.eng.gtefsd.com!howland.reston.ans.net!spool.mu.edu!agate!library.ucla.edu!csulb.edu!csus.edu!netcom.com!deleyd From: deleyd(a)netcom.com Subject: Re: BSD random() - any good (source included)? Message-ID: <deleydCs1L9A.2L6(a)netcom.com> Followup-To: sci.math,sci.stat.math Organization: NETCOM On-line Communication Services (408 261-4700 guest) X-Newsreader: TIN [version 1.2 PL1] References: <2ua9ln$4lv(a)news.tamu.edu> Date: Mon, 27 Jun 1994 06:12:45 GMT Lines: 13 Xref: bga.com sci.math:14916 sci.stat.math:1219 I did a research paper on Computer Generated Random Number Sequences in 1991. Included are the results of testing numerous popular generators. The code used for testing the generators is also available if one is so inclined to do some testing of a particular generator. (The only thing is a thorough test to determine the limits of the generator can take many hours of CPU time). Perhaps later this week I'll post the paper and see what the response is. I'm always a bit apprehensive to post. Never sure what the response will be. Maybe someone will think it's interesting. David Deley deleyd(a)netcom.com

1 0