Lecture notes on Computational Physics

computational physics matlab lecture notes pdf free download
FrankRoberts Profile Pic
FrankRoberts,France,Researcher
Published Date:11-07-2017
Your Website URL(Optional)
Comment
COMPUTATIONAL PHYSICS Morten Hjorth-Jensen University of Oslo, Fall 2009Contents I Introduction to Numerical Methods in Physics 1 1 Introduction 3 1.1 Choice of programming language . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 1.2 Designing programs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 2 Introduction to C++ and Fortran 9 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.2 Getting started . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.2.1 Scientific hello world . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.3 Representation of integer numbers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.3.1 Fortran codes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.3.2 Python codes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 2.4 Real numbers and numerical precision . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 2.4.1 Representation of real numbers . . . . . . . . . . . . . . . . . . . . . . . . . . 21 2.4.2 Machine numbers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.5 Programming examples on loss of precision and round-off errors . . . . . . . . . . . . . 25 −x 2.5.1 Algorithms fore . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 2.5.2 Fortran codes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.5.3 Further examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 2.6 Additional features of C++ and Fortran . . . . . . . . . . . . . . . . . . . . . . . . . . 34 2.6.1 Operators in C++ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 2.6.2 Pointers and arrays in C++. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 2.6.3 Macros in C++ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 2.6.4 Structures in C++ and TYPE in Fortran . . . . . . . . . . . . . . . . . . . . . . 39 2.7 Exercises and projects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 3 Numerical differentiation 45 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 3.2 Numerical differentiation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 x 3.2.1 The second derivative ofe . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 3.2.2 Error analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 3.3 Exercises and projects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 4 Linear algebra 63 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.2 Mathematical intermezzo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 ixContents 4.3 Programming details . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 4.3.1 Declaration of fixed-sized vectors and matrices . . . . . . . . . . . . . . . . . . 68 4.3.2 Runtime declarations of vectors and matrices in C++ . . . . . . . . . . . . . . . 69 4.3.3 Matrix operations and C++ and Fortran features of matrix handling . . . . . . . 74 4.4 Linear Systems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.4.1 Gaussian elimination . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 4.4.2 LU decomposition of a matrix . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.4.3 Solution of linear systems of equations . . . . . . . . . . . . . . . . . . . . . . 87 4.4.4 Inverse of a matrix and the determinant . . . . . . . . . . . . . . . . . . . . . . 89 4.4.5 Tridiagonal systems of linear equations . . . . . . . . . . . . . . . . . . . . . . 94 4.5 Exercises and projects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 5 Non-linear equations and roots of polynomials 105 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105 5.2 Iteration methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106 5.3 Bisection method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 5.4 Newton-Raphson’s method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110 5.5 The secant method and other methods . . . . . . . . . . . . . . . . . . . . . . . . . . . 113 5.6 Exercises and projects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 6 Numerical interpolation, extrapolation and fitting of data 117 6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 6.2 Interpolation and extrapolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 6.2.1 Polynomial interpolation and extrapolation . . . . . . . . . . . . . . . . . . . . 117 6.3 Richardson’s deferred extrapolation method . . . . . . . . . . . . . . . . . . . . . . . . 120 6.4 Cubic spline interpolation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121 7 Numerical integration 125 7.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125 7.2 Newton-Cotes quadrature: equal step methods . . . . . . . . . . . . . . . . . . . . . . . 125 7.3 Adaptive integration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 7.4 Gaussian quadrature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133 7.4.1 Orthogonal polynomials, Legendre . . . . . . . . . . . . . . . . . . . . . . . . 136 7.4.2 Integration points and weights with orthogonal polynomials . . . . . . . . . . . 138 7.4.3 Application to the caseN = 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . 140 7.4.4 General integration intervals for Gauss-Legendre . . . . . . . . . . . . . . . . . 142 7.4.5 Other orthogonal polynomials . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 7.4.6 Applications to selected integrals . . . . . . . . . . . . . . . . . . . . . . . . . 144 7.5 Treatment of singular Integrals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 146 7.6 Scattering equation and principal value integrals . . . . . . . . . . . . . . . . . . . . . . 149 7.7 Parallel computing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 151 7.7.1 Brief survey of supercomputing concepts and terminologies . . . . . . . . . . . 151 7.7.2 Parallelism . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 152 7.7.3 MPI with simple examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 154 7.7.4 Numerical integration with MPI . . . . . . . . . . . . . . . . . . . . . . . . . . 159 xChapter 1 Introduction ... Die Untersuchungsmethode, deren ich mich bedient habe und die auf ökonomis- che Probleme noch nicht angewandt wurde, macht die Lektüre der ersten Kapitel ziemlich schwierig, und es ist zu befürchten, daβ das französische Publikum, stets ungeduldig nach dem Ergebnis und begierig, den Zusammenhang zwischen den allgemeinen Grundsätzen und den Fragen zu erkennen, die es unmittelbar bewegen, sich abschrecken läβt, weil es nicht sofort weiter vordringen kann. Das ist ein Nachteil, gegen den ich nichts weiter unternehmen kann, als die nach Wahrheit strebenden Leser von vornherein darauf hinzuweisen und gefaβt zu machen. Es gibt keine Landstraβe für die Wissenschaft, und nur diejenigen haben, Aussicht, ihre lichten Höhen zu erreichen, die die Mühe nicht scheuen, ihre steilen Pfade zu erklimmen. Karl Marx, preface to the french edition of ’Das Kapital’, Vol. I In the physical sciences we often encounter problems of evaluating various properties of a given function f(x). Typical operations are differentiation, integration and finding the roots of f(x). In most cases we do not have an analytical expression for the function f(x) and we cannot derive explicit formulae for derivatives etc. Even if an analytical expression is available, the evaluation of certain operations on f(x) are so difficult that we need to resort to a numerical evaluation. More frequently,f(x) is the result of complicated numerical operations and is thus known only at a set of discrete points and needs to be approximated by some numerical methods in order to obtain derivatives, etc etc. The aim of these lecture notes is to give you an introduction to selected numerical methods which are encountered in the physical sciences. Several examples, with varying degrees of complexity, will be used in order to illustrate the application of these methods. The text gives a survey over some of the most used methods in computational physics and each chapter ends with one or more applications to realistic systems, from the structure of a neutron star to the description of quantum mechanical systems through Monte-Carlo methods. Among the algorithms we discuss, are some of the top algorithms in computational science. In recent surveys by Dongarra and Sullivan 1 and Cipra 2, the list over the ten top algorithms of the 20th century include 1. The Monte Carlo method or Metropolis algorithm, devised by John von Neumann, Stanislaw Ulam, and Nicholas Metropolis, discussed in chapters 8-11. 2. The simplex method of linear programming, developed by George Dantzig. 3. Krylov Subspace Iteration method for large eigenvalue problems in particular, developed by Mag- nus Hestenes, Eduard Stiefel, and Cornelius Lanczos, discussed in chapter 12. 3Introduction 4. The Householder matrix decomposition, developed by Alston Householder and discussed in chap- ter 12. 5. The Fortran compiler, developed by a team lead by John Backus, codes used throughout this text. 6. The QR algorithm for eigenvalue calculation, developed by Joe Francis, discussed in chapter 12 7. The Quicksort algorithm, developed by Anthony Hoare. 8. Fast Fourier Transform, developed by James Cooley and John Tukey, discussed in chapter 19 9. The Integer Relation Detection Algorithm, developed by Helaman Ferguson and Rodney 10. The fast Multipole algorithm, developed by Leslie Greengard and Vladimir Rokhlin; (to calculate 2 gravitational forces in an N-body problem normally requiresN calculations. The fast multipole method uses order N calculations, by approximating the effects of groups of distant particles using multipole expansions) The topics we cover start with an introduction to C++, Python and Fortran programming combining it with a discussion on numerical precision, a point we feel is often neglected in computational science. This chapter serves also as input to our discussion on numerical derivation in chapter 3. In that chapter we introduce several programming concepts such as dynamical memory allocation and call by reference and value. Several program examples are presented in this chapter. For those who choose to program in C++ we give also an introduction to the auxiliary library Blitz++, which contains several useful classes for numerical operations on vectors and matrices. The link to Blitz++, matrices and selected algorithms for linear algebra problems are dealt with in chapter 4. Chapters 5 and 6 deal with the solution of non- linear equations and the finding of roots of polynomials and numerical interpolation, extrapolation and data fitting. Therafter we switch to numerical integration for integrals with few dimensions, typically less than three, in chapter 7. The numerical integration chapter serves also to justify the introduction of Monte- Carlo methods discussed in chapters 8 and 9. There, a variety of applications are presented, from in- tegration of multidimensional integrals to problems in statistical physics such as random walks and the derivation of the diffusion equation from Brownian motion. Chapter 10 continues this discussion by ex- tending to studies of phase transitions in statistical physics. Chapter 11 deals with Monte-Carlo studies of quantal systems, with an emphasis on variational Monte Carlo methods and diffusion Monte Carlo meth- ods. In chapter 12 we deal with eigensystems and applications to e.g., the Schrödinger equation rewritten as a matrix diagonalization problem. Problems from scattering theory are also discussed, together with the most used solution methods for systems of linear equations. Finally, we discuss various methods for solving differential equations and partial differential equations in chapters 13-15 with examples ranging from harmonic oscillations, equations for heat conduction and the time dependent Schrödinger equation. The emphasis is on various finite difference methods. We assume that you have taken an introductory course in programming and have some familiarity with high-level or low-level and modern languages such as Java, Python, C++, Fortran 77/90/95, etc. 1 Fortran and C++ are examples of compiled low-level languages, in contrast to interpreted ones like Maple or Matlab. In such compiled languages the computer translates an entire subprogram into basic machine instructions all at one time. In an interpreted language the translation is done one statement at a time. This clearly increases the computational time expenditure. More detailed aspects of the above two programming languages will be discussed in the lab classes and various chapters of this text. 1 With Fortran we will consistently mean Fortran 2003. There are no programming examples in Fortran 77 in this text. 41.1 – Choice of programming language There are several texts on computational physics on the market, see for example Refs. 3–10, ranging from introductory ones to more advanced ones. Most of these texts treat however in a rather cavalier way the mathematics behind the various numerical methods. We’ve also succumbed to this approach, mainly due to the following reasons: several of the methods discussed are rather involved, and would thus require at least a one-semester course for an introduction. In so doing, little time would be left for problems and computation. This course is a compromise between three disciplines, numerical methods, problems from the physical sciences and computation. To achieve such a synthesis, we will have to relax our presentation in order to avoid lengthy and gory mathematical expositions. You should also keep in mind that computational physics and science in more general terms consist of the combination of several fields and crafts with the aim of finding solution strategies for complicated problems. However, where we do indulge in presenting more formalism, we have borrowed heavily from several texts on mathematical analysis. 1.1 Choice of programming language As programming language we have ended up with preferring C++, but all examples discussed in the text have their corresponding Fortran and Python programs on the webpage of this text. Fortran (FORmula TRANslation) was introduced in 1957 and remains in many scientific computing environments the language of choice. The latest standard, Fortran 11–14, includes extensions that are familiar to users of C++. Some of the most important features of Fortran include recursive subroutines, dynamic storage allocation and pointers, user defined data structures, modules, and the ability to manip- ulate entire arrays. However, there are several good reasons for choosing C++ as programming language for scientific and engineering problems. Here are some: – C++ is now the dominating language in Unix and Windows environments. It is widely available and is the language of choice for system programmers. It is very widespread for developments of non-numerical software – The C++ syntax has inspired lots of popular languages, such as Perl, Python and Java. – It is an extremely portable language, all Linux and Unix operated machines have a C++ compiler. – In the last years there has been an enormous effort towards developing numerical libraries for C++. Numerous tools (numerical libraries such as MPI 15–17) are written in C++ and interfacing them requires knowledge of C++. Most C++ and Fortran compilers compare fairly well when it comes to speed and numerical efficiency. Although Fortran 77 and C are regarded as slightly faster than C++ or Fortran, compiler improvements during the last few years have diminshed such differences. The Java numerics project has lost some of its steam recently, and Java is therefore normally slower than C++ or Fortran, see however the Java Numerics homepage for a discussion on numerical aspects of Java 18. – Complex variables, one of Fortran’s strongholds, can also be defined in the new ANSI C++ stan- dard. – C++ is a language which catches most of the errors as early as possible, typically at compilation time. Fortran has some of these features if one omits implicit variable declarations. 5Introduction – C++ is also an object-oriented language, to be contrasted with C and Fortran. This means that it supports three fundamental ideas, namely objects, class hierarchies and polymorphism. For- tran has, through theMODULE declaration the capability of defining classes, but lacks inheritance, although polymorphism is possible. Fortran is then considered as an object-based programming language, to be contrasted with C++ which has the capability of relating classes to each other in a hierarchical way. An important aspect of C++ is its richness with more than 60 keywords allowing for a good balance between object orientation and numerical efficiency. Furthermore, careful programming can results in an efficiency close to Fortran 77. The language is well-suited for large projects and has presently good standard libraries suitable for computational science projects, although many of these still lag behind the large body of libraries for numerics available to Fortran programmers. However, it is not difficult to interface libraries written in Fortran with C++ codes, if care is exercised. Other weak sides are the fact that it can be easy to write inefficient code and that there are many ways of writing the same things, adding to the confusion for beginners and professionals as well. The language is also under continuous development, which often causes portability problems. C++ is also a difficult language to learn. Grasping the basics is rather straightforward, but takes time to master. A specific problem which often causes unwanted or odd errors is dynamic memory management. The efficiency of C++ codes are close to those provided by Fortran. This means often that a code written in Fortran 77 can be faster, however for large numerical projects C++ and Fortran are to be preferred. If speed is an issue, one could port critical parts of the code to Fortran 77. Future plans Since our undergraduate curriculum has changed considerably from the beginning of Fall-2007, with the introduction of Python as programming language, the content of this course will change accordingly from the fall semester 2009. C++ and Fortran will then coexist with Python and students can choose between these three programming languages. The emphasis in the text will be on C++ programming, but how to interface C++ or Fortran programs with Python codes will also be discussed. 1.2 Designing programs Before we proceed with a discussion of numerical methods, we would like to remind you of some aspects of program writing. In writing a program for a specific algorithm (a set of rules for doing mathematics or a precise description of how to solve a problem), it is obvious that different programmers will apply different 2 styles, ranging from barely readable (even for the programmer) to well documented codes which can be used and extended upon by others in e.g., a project. The lack of readability of a program leads in many cases to credibility problems, difficulty in letting others extend the codes or remembering oneself what a certain statement means, problems in spotting errors, not always easy to implement on other machines, and so forth. Although you should feel free to follow your own rules, we would like to focus certain suggestions which may improve a program. What follows here is a list of our recommendations (or biases/prejudices). 2 As an example, a bad habit is to use variables with no specific meaning, like x1, x2 etc, or names for subprograms which go like routine1, routine2 etc. 61.2 – Designing programs First about designing a program. – Before writing a single line, have the algorithm clarified and understood. It is crucial to have a logical structure of e.g., the flow and organization of data before one starts writing. – Always try to choose the simplest algorithm. Computational speed can be improved upon later. – Try to write a as clear program as possible. Such programs are easier to debug, and although it may take more time, in the long run it may save you time. If you collaborate with other people, it reduces spending time on debugging and trying to understand what the codes do. A clear program will also allow you to remember better what the program really does – Implement a working code with emphasis on design for extensions, maintenance etc. Focus on the design of your code in the beginning and don’t think too much about efficiency before you have a thoroughly debugged and verified your program. A rule of thumb is the so-called 80−20 rule, 80 % of the CPU time is spent in 20 % of the code and you will experience that typically onlya small part of your code is responsible for most of the CPU expenditure. Therefore, spend most of your time in devising a good algorithm. – The planning of the program should be from top down to bottom, trying to keep the flow as linear as possible. Avoid jumping back and forth in the program. First you need to arrange the major tasks to be achieved. Then try to break the major tasks into subtasks. These can be represented by functions or subprograms. They should accomplish limited tasks and as far as possible be independent of each other. That will allow you to use them in other programs as well. – Try always to find some cases where an analytical solution exists or where simple test cases can be applied. If possible, devise different algorithms for solving the same problem. If you get the same answers, you may have coded things correctly or made the same error twice. – When you have a working code, you should start thinking of the efficiency. Analyze the efficiency with a tool (profiler) to predict the CPU-intensive parts. Attack then the CPU-intensive parts after the program reproduces benchmark results. However, although we stress that you should post-pone a discussion of the efficiency of your code to the stage when you are sure that it runs correctly, there are some simple guidelines to follow when you design the algorithm. – Avoid lists, sets etc., when arrays can be used without too much waste of memory. Avoid also calls to functions in the innermost loop since that produces an overhead in the call. – Heavy computation with small objects might be inefficient, e.g., vector of class complex objects – Avoid small virtual functions (unless they end up in more than (say) 5 multiplications) – Save object-oriented constructs for the top level of your code. – Use taylored library functions for various operations, if possible. – Reduce pointer-to-pointer-to....-pointer links inside loops. – Avoid implicit type conversion, use rather the explicit keyword when declaring constructors in C++. 7Introduction – Never return (copy) of an object from a function, since this normally implies a hidden allocation. Finally, here are some of our favoured approaches for writing a code. – Use always the standard ANSI version of the programming language. Avoid local dialects if you wish to port your code to other machines. – Add always comments to describe what a program or subprogram does. Comment lines help you remember what you did e.g., one month ago. – Declare all variables. Avoid totally theIMPLICIT statement in Fortran. The program will be more readable and help you find errors when compiling. – Do not useGOTO structures in Fortran. Although all varieties of spaghetti are great culinaric temp- tations, spaghetti-like Fortran with manyGOTO statements is to be avoided. Extensive amounts of time may be wasted on decoding other authors’ programs. – When you name variables, use easily understandable names. Avoidv1 when you can use speed_of_light . Associatives names make it easier to understand what a specific subprogram does. – Use compiler options to test program details and if possible also different compilers. They make errors too. – Writing codes in C++ and Fortran may often lead to segmentation faults. This means in most cases that we are trying to access elements of an array which are not available. When developing a code it is then useful to compile with debugging options. The use of debuggers like gdb is something we highly recommend during the development of a program. 8Chapter 2 Introduction to C++ and Fortran Computers in the future may weigh no more than 1.5 tons. Popular Mechanics, 1949 There is a world market for maybe five computers. Thomas Watson, IBM chairman, 1943 2.1 Introduction This chapters aims at catching two birds with a stone; to introduce to you essential features of the pro- gramming languages C++ and Fortran with a brief reminder on Python specific topics, and to stress problems like overflow, underflow, round off errors and eventually loss of precision due to the finite amount of numbers a computer can represent. The programs we discuss are taylored to these aims. 2.2 Getting started 1 In programming languages we encounter data entities such as constants, variables, results of evaluations of functions etc. Common to these objects is that they can be represented through the type concept. There are intrinsic types and derived types. Intrinsic types are provided by the programming language whereas derived types are provided by the programmer. If one specifies the type to be for example 2 INTEGER(KIND=2) for Fortran orshortint/int in C++, the programmer selects a particular date type with 2 bytes (16 bits) for every item of the classINTEGER(KIND=2) orint. Intrinsic types come in two classes, numerical (like integer, real or complex) and non-numeric (as logical and character). The general form for declaring variables isdatatypenameofvariable and Table 2.2 lists the standard variable declarations of C++ and Fortran (note well that there be may compiler and machine differences from the table below). An important aspect when declaring variables is their region of validity. 1 For more detailed texts on C++ programming in engineering and science are the books by Flowers 19 and Bar- ton and Nackman 20. The classic text on C++ programming is the book of Bjarne Stoustrup 21. See also the lec- ture notes on C++ athttp://heim.ifi.uio.no/hpl/INF-VERK4830. For Fortran we recommend the online lectures at http://folk.uio.no/gunnarw/INF-VERK4820. These web pages contain extensive references to other C++ and Fortran resources. Both web pages contain enough material, lecture notes and exercises, in order to serve as material for own studies. The Fortran 95 standard is well documented in Refs. 11–13 while the new details of Fortran 2003 can be found in Ref. 14. The reader should note that this is not a text on C++ or Fortran . It is therefore important than one tries to find additional literature on these programming languages. Good Python texts on scientific computing are 22, 23. 2 Our favoured display mode for Fortran statements will be capital letters for language statements and low key letters for user-defined statements. Note that Fortran does not distinguish between capital and low key letters while C++ does. 9Introduction to C++ and Fortran Inside a function we define a a variable through the expressionintvar orINTEGER::var . The question is whether this variable is available in other functions as well, moreover where is var initialized and finally, if we call the function where it is declared, is the value conserved from one call to the other? type in C++ and Fortran bits range int/INTEGER (2) 16 −32768 to 32767 unsigned int 16 0 to 65535 signed int 16 −32768 to 32767 short int 16 −32768 to 32767 unsigned short int 16 0 to 65535 signed short int 16 −32768 to 32767 int/long int/INTEGER(4) 32 −2147483648 to 2147483647 signed long int 32 −2147483648 to 2147483647 −44 +38 float/REAL(4) 32 10 to10 −322 +308 double/REAL(8) 64 10 to10e Table 2.1: Examples of variable declarations for C++ and Fortran . We reserve capital letters for Fortran declaration statements throughout this text, although Fortran is not sensitive to upper or lowercase letters. Note that there are machines which allow for more than 64 bits for doubles. The ranges listed here may therefore vary. Both C++ and Fortran operate with several types of variables and the answers to these questions depend on how we have defined for example an integer via the statementintvar. Python on the other hand does not use variable or function types (they are not explicitely written), allowing thereby for a better potential for reuse of the code. The following list may help in clarifying the above points: type of variable validity local variables defined within a function, only available within the scope of the function. formal parameter If it is defined within a function it is only available within that specific function. global variables Defined outside a given function, available for all functions from the point where it is defined. In Table 2.2 we show a list of some of the most used language statements in Fortran and C++. In addition, both C++ and Fortran allow for complex variables. In Fortran we would declare a com- plex variable asCOMPLEX(KIND=16)::x,y which refers to a double with word length of 16 bytes. In C++ we would need to include a complex library through the statements i n c l u d e complex complex double x , y ; We will discuss the above declaration complexdouble x,y; in more detail in appendix A. 102.2 – Getting started Fortran C++ Program structure PROGRAM something main () FUNCTION something(input) double (int) something(input) SUBROUTINE something(inout) Data type declarations REAL (4) x, y float x, y; REAL(8) :: x, y double x, y; INTEGER :: x, y int x,y; CHARACTER :: name char name; REAL(8), DIMENSION(dim1,dim2) :: x double xdim1dim2; INTEGER, DIMENSION(dim1,dim2) :: x int xdim1dim2; LOGICAL :: x TYPE name struct name declarations declarations; END TYPE name POINTER :: a double (int) a; ALLOCATE new; DEALLOCATE delete; Logical statements and control structure IF ( a == b) THEN if ( a == b) b=0 b=0; ENDIF DO WHILE (logical statement) while (logical statement) do something do something ENDDO IF ( a= b ) THEN if ( a= b) b=0 b=0; ELSE else a=0 a=0; ENDIF SELECT CASE (variable) switch(variable) CASE (variable=value1) do something case 1: CASE (... ) variable=value1; ... do something; break; END SELECT case 2: do something; break;... DO i=0, end, 1 for( i=0; i= end; i++) do something do something ; ENDDO Table 2.2: Elements of programming syntax. 11Introduction to C++ and Fortran 2.2.1 Scientific hello world Our first programming encounter is the ’classical’ one, found in almost every textbook on computer languages, the ’hello world’ code, here in a scientific disguise. We present first the C version. http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/ pp/program1. pp /∗ comments i n C b e g i n l i k e t h i s and end w i t h ∗ / i n c l u d e s t d l i b . h /∗ a t o f f u n c t i o n ∗ / i n c l u d e math . h /∗ s i n e f u n c t i o n ∗ / i n c l u d e s t d i o . h /∗ p r i n t f f u n c t i o n ∗ / i n t main ( i n t argc , char∗ arg v ) double r , s ; /∗ d e c l a r e v a r i a b l e s ∗ / r = a t o f ( arg v 1 ) ; /∗ c o n v e r t t h e t e x t argv 1 t o d o u b le ∗ / s = s i n ( r ) ; p r i n t f ("Hello,Worldsin(%g)=%g\n" , r , s ) ; return 0 ; /∗ s u c c e s s e x e c u t i o n o f t h e program ∗ / The compiler must see a declaration of a function before you can call it (the compiler checks the argument and return types). The declaration of library functions appears in so-called header files that must be included in the program, for example include stdlib .h. We call three functions atof , sin , printf and these are declared in three different header files. The main program is a function called main with a return value set to an integer, returning 0 if success. The operating system stores the return value, and other programs/utilities can check whether the execution was successful or not. The command-line arguments are transferred to the main function through the statement int main ( int argc , char∗ argv ) . The integer argc stands for the number of command-line arguments, set to one in our case, while argv is a vector of strings containing the command-line argu- ments with argv0 containing the name of the program and argv1, argv2, ... are the command-line args, i.e., the number of lines of input to the program. This means that we would run the programs as mhjensen ompphys:./myprogram.exe0.3 argv0 while the text string0.2 enters argv1. Here we define a floating point variable, see also below, through the keywords float for single pre- cision real numbers anddouble for double precision. The function atof transforms a text (argv 1) to a float. The sine function is declared in math.h, a library which is not automatically included and needs to be linked when computing an executable file. With the command printf we obtain a formatted printout. The printf syntax is used for formatting output in many C-inspired languages (Perl, Python, awk, partly C++). In C++ this program can be written as / / A comment l i n e b e g i n s l i k e t h i s i n C++ programs using namespace s t d ; i n c l u d e i o s t r e a m i n t main ( i n t argc , char∗ arg v ) / / c o n v e r t t h e t e x t argv 1 t o d o u b le u s i n g a t o f : double r = a t o f ( arg v 1 ) ; 122.2 – Getting started double s = s i n ( r ) ; c o u t "Hello,Worldsin(" r ")=" s e n d l ; / / s u c c e s s return 0 ; We have replaced the call to printf with the standard C++ function cout. The header file iostream is then needed. In addition, we don’t need to declare variables liker ands at the beginning of the program. I personally prefer however to declare all variables at the beginning of a function, as this gives me a feeling of greater readability. Note that we have used the declaration using namespace std;. Namespace is a way to collect all functions defined in C++ libraries. If we omit this declaration on top of the program we would have to add the declaration std in front of cout or cin. Our program would then read / / H e l l o world code w i t h o u t u s i n g namespace s t d i n c l u d e i o s t r e a m i n t main ( i n t argc , char∗ arg v ) / / c o n v e r t t h e t e x t argv 1 t o d o u b le u s i n g a t o f : double r = a t o f ( arg v 1 ) ; double s = s i n ( r ) ; s t d : : c o u t "Hello,Worldsin(" r ")=" s e n d l ; / / s u c c e s s return 0 ; Another feature which is worth noting is that we have skipped exception handlings here. In chapter 3 we discuss examples that test our input from the command line. But it is easy to add such a feature, as shown in our modified hello world program / / H e l l o world code w i t h e x c e p t i o n h a n d l i n g using namespace s t d ; i n c l u d e i o s t r e a m i n t main ( i n t argc , char∗ arg v ) / / Read i n o u t p u t f i l e , a b o r t i f t h e r e a re t o o few command−l i n e a rg u men ts i f ( a r g c = 1 ) c o u t "BadUsage:" arg v 0 "readalsoanumberonthesameline,e.g.,prog.exe0.2" e n d l ; e x i t ( 1 ) ; / / h ere t h e program s t o p s . / / c o n v e r t t h e t e x t argv 1 t o d o u b le u s i n g a t o f : double r = a t o f ( arg v 1 ) ; double s = s i n ( r ) ; c o u t "Hello,Worldsin(" r ")=" s e n d l ; / / s u c c e s s return 0 ; Here we test that we have more than one argument. If not, the program stops and writes to screen an error message. To run these programs, you need first to compile and link them in order to obtain an executable file under operating systems like e.g., UNIX or Linux. Before we proceed we give therefore examples on how to obtain an executable file under Linux/Unix. 13Introduction to C++ and Fortran In order to obtain an executable file for a C++ program, the following instructions under Linux/Unix can be used ++- -Wallmyprogram. ++-omyprogrammyprogram.o where the compiler is called through the command ++. The compiler option -Wall means that a warning is issued in case of non-standard language. The executable file is in this casemyprogram. The option- is for compilation only, where the program is translated into machine code, while the-o option links the produced object filemyprogram.o and produces the executablemyprogram . The corresponding Fortran code is http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/f90/program1.f90 PROGRAM shw IMPLICIT NONE REAL (KIND =8) : : r I n p u t number REAL (KIND=8) : : s R e s u l t Get a number from u s e r WRITE(∗ ,∗ ) ’ I n p u t a number : ’ READ(∗ ,∗ ) r C a l c u l a t e t h e s i n e o f t h e number s = SIN ( r ) W r i t e r e s u l t t o s c r e e n WRITE(∗ ,∗ ) ’ H e l l o World SINE o f ’ , r , ’ = ’ , s END PROGRAM shw The first statement must be a program statement; the last statement must have a corresponding end pro- gram statement. Integer numerical variables and floating point numerical variables are distinguished. The names of all variables must be between 1 and 31 alphanumeric characters of which the first must be a letter and the last must not be an underscore. Comments begin with a and can be included anywhere in the program. Statements are written on lines which may contain up to 132 characters. The asterisks (,) following WRITE represent the default format for output, i.e., the output is e.g., written on the screen. Similarly, the READ(,) statement means that the program is expecting a line input. Note also the IMPLICIT NONE statement which we strongly recommend the use of. In many Fortran 77 programs one can find statements like IMPLICIT REAL8(a-h,o-z), meaning that all variables beginning with any of the above letters are by default floating numbers. However, such a usage makes it hard to spot eventual errors due to misspelling of variable names. With IMPLICIT NONE you have to declare all variables and therefore detect possible errors already while compiling. I recommend strongly that you declare all variables when using Fortran. We call the Fortran compiler (using free format) through f90- -freemyprogram.f90 f90-omyprogram.xmyprogram.o Under Linux/Unix it is often convenient to create a so-called makefile, which is a script which includes possible compiling commands, in order to avoid retyping the above lines every once and then we have made modifcations to our program. A typical makefile for the abovecc compiling options is listed below Generalmakefilefor - hoosePROG=nameofgivenprogram 142.2 – Getting started Herewedefine ompileroption,librariesandthetarget CC= ++-Wall PROG=myprogram Herewemaketheexe utablefile PROG:PROG.o CCPROG.o-oPROG whereasherewe reatetheobje tfile PROG.o:PROG. pp CC- PROG. pp If you name your file for ’makefile’, simply type the command make and Linux/Unix executes all of the statements in the above makefile. Note that C++ files have the extension .cpp For Fortran, a similar makefile is GeneralmakefileforF90- hoosePROG=nameofgivenprogram Herewedefine ompileroptions,librariesandthetarget F90=f90 PROG=myprogram Herewemaketheexe utablefile PROG:PROG.o F90PROG.o-oPROG whereasherewe reatetheobje tfile PROG.o:PROG.f90 F90- PROG.f Finally, for the sake of completeness, we list the corresponding Python code http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/python/program1.py / u s r / b i n / env p yth o n import sys , math Read i n a s t r i n g a c o n v e r t i t t o a f l o a t r = f l o a t ( s y s . arg v 1 ) s = math . s i n ( r ) p r i n t"Hello,Worldsin(%g)=%12.6e" % ( r , s ) where we have used a formatted printout with scientific notation. 15Introduction to C++ and Fortran 2.3 Representation of integer numbers In Fortran a keyword for declaration of an integer is INTEGER (KIND=n) , n = 2 reserves 2 bytes (16 bits) of memory to store the integer variable wheras n = 4 reserves 4 bytes (32 bits). In Fortran, although it may be compiler dependent, just declaring a variable as INTEGER , reserves 4 bytes in memory as default. In C++ keywords areshort int , int , long int , long long int . The byte-length is compiler dependent within some limits. The GNU C++-compilers (called by gcc or g++) assign 4 bytes (32 bits) to variables declared by int and long int . Typical byte-lengths are 2, 4, 4 and 8 bytes, for the types given above. To see how many bytes are reserved for a specific variable, C++ has a library function called sizeof (type) which returns the number of bytes for type . An example of a program declaration is Fortran: INTEGER (KIND=2) :: age_of_participant C++: short int age_of_participant; Note that the (KIND=2) can be written as (2). Normally however, we will for Fortran programs just use the 4 bytes default assignment INTEGER . In the above examples one bit is used to store the sign of the variable age_of_participant and the other 15 15 bits are used to store the number, which then may range from zero to2 −1 = 32767. This should definitely suffice for human lifespans. On the other hand, if we were to classify known fossiles by age we may need Fortran: INTEGER (4) :: age_of_fossile C++: int age_of_fossile; Again one bit is used to store the sign of the variable age_of_fossile and the other 31 bits are used to 31 store the number which then may range from zero to 2 −1 = 2.147.483.647. In order to give you a feeling how integer numbers are represented in the computer, think first of the decimal representation of the number417 2 1 0 417 = 4×10 +1×10 +7×10 , which in binary representation becomes n n−1 n−2 0 417 =1×a 2 +a 2 +a 2 +···+a 2 , n n−1 n−2 0 where the coefficientsa withk = 0,...,n are zero or one. They can be calculated through successive k division by 2 and using the remainder in each division to determine the numbersa toa . A given integer n 0 in binary notation is then written as n n−1 n−2 0 a 2 +a 2 +a 2 +···+a 2 . n n−1 n−2 0 In binary notation we have thus (417) = (110100001) , 10 2 since we have 8 7 6 5 4 3 2 2 1 0 (110100001) = 1×2 +1×2 +0×2 +1×2 +0×2 +0×2 +0×2 +0×2 +0×2 +1×2 . 2 To see this, we have performed the following divisions by 2 162.3 – Representation of integer numbers 0 417/2=208 remainder 1 coefficient of2 is 1 1 208/2=104 remainder 0 coefficient of2 is 0 2 104/2=52 remainder 0 coefficient of2 is 0 3 52/2=26 remainder 0 coefficient of2 is 0 4 26/2=13 remainder 0 coefficient of2 is 0 5 13/2= 6 remainder 1 coefficient of2 is 1 6 6/2= 3 remainder 0 coefficient of2 is 0 7 3/2= 1 remainder 1 coefficient of2 is 1 8 1/2= 0 remainder 1 coefficient of2 is 1 We see that nine bits are sufficient to represent 417. Normally we end up using 32 bits as default for integers, meaning that our number reads (417) = (00000000000000000000000110100001) , 10 2 A simple program which performs these operations is listed below. Here we employ the modulus operation (with division by 2), which in C++ is given by the a%2 operator. In Fortran we would call the function MOD(a,2) in order to obtain the remainder of a division by2. http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/ pp/program2. pp using namespace s t d ; i n c l u d e i o s t r e a m i n t main ( i n t argc , char∗ arg v ) i n t i ; i n t t e r m s 3 2 ; / / s t o r a g e o f a0 , a1 , e t c , up t o 32 b i t s i n t number = a t o i ( arg v 1 ) ; / / i n i t i a l i s e t h e term a0 , a1 e t c f o r ( i =0 ; i 32 ; i ++) t e r m s i = 0 ; f o r ( i =0 ; i 32 ; i ++) t e r m s i = number%2; number /= 2 ; / / w r i t e o u t r e s u l t s c o u t ‘ ‘ Number o f b y t e s used =’’ s i z e o f ( number ) e n d l ; f o r ( i =0 ; i 32 ; i ++) c o u t ‘ ‘ Term n r : ‘ ‘ i ‘ ‘ Value = ‘ ‘ t e r m s i ; c o u t e n d l ; return 0 ; The C++ function sizeof yields the number of bytes reserved for a specific variable. Note also the for construct. We have reserved a fixed array which contains the values ofa being 0 or1, the remainder of i a division by two. We have enforced the integer to be represented by 32 bits, or four bytes, which is the default integer representation. Note that for417 we need 9 bits in order to represent it in a binary notation, while a number like the number 3 is given in an 32 bits word as (3) = (00000000000000000000000000000011) . 10 2 17Introduction to C++ and Fortran For this number 2 significant bits would be enough. With these prerequesites in mind, it is rather obvious that if a given integer variable is beyond the range assigned by the declaration statement we may encounter problems. If we multiply two large integersn ×n and the product is too large for the bit size allocated for that 1 2 specific integer assignement, we run into an overflow problem. The most significant bits are lost and the least significant kept. Using 4 bytes for integer variables the result becomes 20 20 2 ×2 = 0. However, there are compilers or compiler options that preprocess the program in such a way that an error message like ’integer overflow’ is produced when running the program. Here is a small program which may cause overflow problems when running (try to test your own compiler in order to be sure how such problems need to be handled). http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/ pp/program3. pp / / Program t o c a l c u l a t e 2∗∗n using namespace s t d ; i n c l u d e i o s t r e a m i n t main ( ) i n t i n t 1 , i n t 2 , i n t 3 ; / / p r i n t t o s c r e e n c o u t "ReadintheexponentialNfor2N=\n" ; / / read from s c r e e n c i n i n t 2 ; i n t 1 = ( i n t ) pow ( 2 . , ( double ) i n t 2 ) ; c o u t "2N2N=" i n t 1∗ i n t 1 "\n" ; i n t 3 = i n t 1 − 1 ; c o u t "2N(2N-1)=" i n t 1 ∗ i n t 3 "\n" ; c o u t "2N-1=" i n t 3 "\n" ; return 0 ; / / End : program main ( ) If we run this code with an exponentN = 32, we obtain the following output 2N2N=0 2N(2N-1)=-2147483648 2N-1=2147483647 64 We notice that2 exceeds the limit for integer numbers with 32 bits. The program returns0. This can be N N dangerous, since the results from the operation2 (2 −1) is obviously wrong. One possibility to avoid such cases is to add compilation option which flag if an overflow or underflow is reached. 2.3.1 Fortran codes The corresponding Fortran code is http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/f90/program2.f90 PROGRAM b i n a r y _ i n t e g e r 182.3 – Representation of integer numbers IMPLICIT NONE INTEGER i , number , t e r m s ( 0 : 3 1 ) s t o r a g e o f a0 , a1 , e t c , up t o 32 b i t s , n o t e a r r a y l e n g t h r u n n i n g from 0 : 3 1 . F o r t r a n a l l o w s n e g a t i v e i n d e x e s as w e l l . WRITE(∗ ,∗ ) ’ Give a number to t r a n s f o r m to b i n a r y n o t a t i o n ’ READ(∗ ,∗ ) number I n i t i a l i s e t h e t e r m s a0 , a1 e t c t e r m s = 0 F o r t r a n t a k e s o n l y i n t e g e r lo o p v a r i a b l e s DO i =0 , 31 t e r m s ( i ) = MOD( number , 2 ) Modulus f u n c t i o n i n F o r t r a n number = number / 2 ENDDO w r i t e o u t r e s u l t s WRITE(∗ ,∗ ) ’ Bin ary r e p r e s e n t a t i o n ’ DO i =0 , 31 WRITE(∗ ,∗ ) ’ Term n r and v alu e ’ , i , t e r m s ( i ) ENDDO END PROGRAM b i n a r y _ i n t e g e r and http://www.fys.uio.no/ ompphys/ p/programs/FYS3150/ hapter02/f90/program3.f90 PROGRAM i n t e g e r _ e x p IMPLICIT NONE INTEGER : : i n t 1 , i n t 2 , i n t 3 T h i s i s t h e b e g i n o f a comment l i n e i n F o r t r a n 90 Now we read from s c r e e n t h e v a r i a b l e i n t 2 WRITE(∗ ,∗ ) ’Read in t h e number to be e x p o n e n t i a t e d ’ READ(∗ ,∗ ) i n t 2 i n t 1 =2∗∗ i n t 2 WRITE(∗ ,∗ ) ’2N∗2N’ , i n t 1∗ i n t 1 i n t 3 = i n t 1−1 WRITE(∗ ,∗ ) ’2N∗ (2 N−1) ’ , i n t 1∗ i n t 3 WRITE(∗ ,∗ ) ’2N−1 ’ , i n t 3 END PROGRAM i n t e g e r _ e x p In F o r t r a n t h e modulus d i v i s i o n i s p erfo rm ed by t h e i n t r i n s i c f u n c t i o n \ l s t i n l i n e MOD( number , 2 ) in ca se o f a d i v i s i o n by 2 . The e x p o n e n t a t i o n o f a number i s g i v e n by f o r example \ l s t i n l i n e 2∗∗N i n s t e a d o f t h e c a l l to t h e \ l s t i n l i n e pow f u n c t i o n in C++. 2.3.2 Python codes In preparation for fall 2009 19

Advise: Why You Wasting Money in Costly SEO Tools, Use World's Best Free SEO Tool Ubersuggest.