c7-1 (779515)
Текст из файла
7.1 Uniform Deviates275As for references on this subject, the one to turn to first is Knuth [1]. Thentry [2]. Only a few of the standard books on numerical methods [3-4] treat topicsrelating to random numbers.Bratley, P., Fox, B.L., and Schrage, E.L. 1983, A Guide to Simulation (New York: SpringerVerlag). [2]Dahlquist, G., and Bjorck, A.
1974, Numerical Methods (Englewood Cliffs, NJ: Prentice-Hall),Chapter 11. [3]Forsythe, G.E., Malcolm, M.A., and Moler, C.B. 1977, Computer Methods for MathematicalComputations (Englewood Cliffs, NJ: Prentice-Hall), Chapter 10. [4]7.1 Uniform DeviatesUniform deviates are just random numbers that lie within a specified range(typically 0 to 1), with any one number in the range just as likely as any other. Theyare, in other words, what you probably think “random numbers” are. However,we want to distinguish uniform deviates from other sorts of random numbers, forexample numbers drawn from a normal (Gaussian) distribution of specified meanand standard deviation. These other sorts of deviates are almost always generated byperforming appropriate operations on one or more uniform deviates, as we will seein subsequent sections.
So, a reliable source of random uniform deviates, the subjectof this section, is an essential building block for any sort of stochastic modelingor Monte Carlo computer work.System-Supplied Random Number GeneratorsMost C implementations have, lurking within, a pair of library routines forinitializing, and then generating, “random numbers.” In ANSI C, the synopsis is:#include <stdlib.h>#define RAND_MAX ...void srand(unsigned seed);int rand(void);You initialize the random number generator by invoking srand(seed) withsome arbitrary seed.
Each initializing value will typically result in a differentrandom sequence, or a least a different starting point in some one enormously longsequence. The same initializing value of seed will always return the same randomsequence, however.You obtain successive random numbers in the sequence by successive calls torand(). That function returns an integer that is typically in the range 0 to thelargest representable positive value of type int (inclusive).
Usually, as in ANSI C,this largest value is available as RAND_MAX, but sometimes you have to figure it outfor yourself. If you want a random float value between 0.0 (inclusive) and 1.0(exclusive), you get it by an expression likeSample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).CITED REFERENCES AND FURTHER READING:Knuth, D.E.
1981, Seminumerical Algorithms, 2nd ed., vol. 2 of The Art of Computer Programming(Reading, MA: Addison-Wesley), Chapter 3, especially §3.5. [1]276Chapter 7.Random Numbersx = rand()/(RAND_MAX+1.0);Ij+1 = aIj + c(mod m)(7.1.1)Here m is called the modulus, and a and c are positive integers called the multiplierand the increment respectively. The recurrence (7.1.1) will eventually repeat itself,with a period that is obviously no greater than m.
If m, a, and c are properly chosen,then the period will be of maximal length, i.e., of length m. In that case, all possibleintegers between 0 and m − 1 occur at some point, so any initial “seed” choice of I0is as good as any other: the sequence just takes off from that point.Although this general framework is powerful enough to provide quite decentrandom numbers, its implementation in many, if not most, ANSI C libraries is quiteflawed; quite a number of implementations are in the category “totally botched.”Blame should be apportioned about equally between the ANSI C committee andthe implementors.
The typical problems are these: First, since the ANSI standardspecifies that rand() return a value of type int — which is only a two-byte quantityon many machines — RAND_MAX is often not very large. The ANSI C standardrequires only that it be at least 32767. This can be disastrous in many circumstances:for a Monte Carlo integration (§7.6 and §7.8), you might well want to evaluate 106different points, but actually be evaluating the same 32767 points 30 times each, notat all the same thing! You should categorically reject any library random numberroutine with a two-byte returned value.Second, the ANSI committee’s published rationale includes the followingmischievous passage: “The committee decided that an implementation should beallowed to provide a rand function which generates the best random sequencepossible in that implementation, and therefore mandated no standard algorithm. Itrecognized the value, however, of being able to generate the same pseudo-randomsequence in different implementations, and so it has published an example.
. . .[emphasis added]” The “example” isunsigned long next=1;int rand(void) /* NOT RECOMMENDED (see text) */{next = next*1103515245 + 12345;return (unsigned int)(next/65536) % 32768;}void srand(unsigned int seed){next=seed;}Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).Now our first, and perhaps most important, lesson in this chapter is: be very,very suspicious of a system-supplied rand() that resembles the one just described.If all scientific papers whose results are in doubt because of bad rand()s wereto disappear from library shelves, there would be a gap on each shelf about asbig as your fist.
System-supplied rand()s are almost always linear congruentialgenerators, which generate a sequence of integers I1 , I2 , I3 , . . ., each between 0 andm − 1 (e.g., RAND_MAX) by the recurrence relation7.1 Uniform Deviates277Correlation in k-space is not the only weakness of linear congruential generators.Such generators often have their low-order (least significant) bits much less randomthan their high-order bits. If you want to generate a random integer between 1 and10, you should always do it using high-order bits, as inj=1+(int) (10.0*rand()/(RAND_MAX+1.0));and never by anything resemblingj=1+(rand() % 10);(which uses lower-order bits). Similarly you should never try to take apart a“rand()” number into several supposedly random pieces.
Instead use separatecalls for every piece.Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited.
To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).This corresponds to equation (7.1.1) with a = 1103515245, c = 12345, and m = 232(since arithmetic done on unsigned long quantities is guaranteed to return thecorrect low-order bits). These are not particularly good choices for a and c, thoughthey are not gross embarrassments by themselves.
The real botches occur whenimplementors, taking the committee’s statement above as license, try to “improve”on the published example. For example, one popular 32-bit PC-compatible compilerprovides a long generator that uses the above congruence, but swaps the high-orderand low-order 16 bits of the returned value. Somebody probably thought that thisextra flourish added randomness; in fact it ruins the generator. While these kinds ofblunders can, of course, be fixed, there remains a fundamental flaw in simple linearcongruential generators, which we now discuss.The linear congruential method has the advantage of being very fast, requiringonly a few operations per call, hence its almost universal use. It has the disadvantagethat it is not free of sequential correlation on successive calls. If k random numbers ata time are used to plot points in k dimensional space (with each coordinate between0 and 1), then the points will not tend to “fill up” the k-dimensional space, but ratherwill lie on (k − 1)-dimensional “planes.” There will be at most about m1/k suchplanes.
If the constants m, a, and c are not very carefully chosen, there will be manyfewer than that. If m is as bad as 32768, then the number of planes on which triplesof points lie in three-dimensional space will be no greater than about the cube rootof 32768, or 32. Even if m is close to the machine’s largest representable integer,e.g., ∼ 232 , the number of planes on which triples of points lie in three-dimensionalspace is usually no greater than about the cube root of 232, about 1600. You mightwell be focusing attention on a physical process that occurs in a small fraction of thetotal volume, so that the discreteness of the planes can be very pronounced.Even worse, you might be using a generator whose choices of m, a, and c havebeen botched.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.















