Computational Physics - Darren Walker - E-Book

Computational Physics E-Book

Darren Walker

0,0
29,99 €

-100%
Sammeln Sie Punkte in unserem Gutscheinprogramm und kaufen Sie E-Books und Hörbücher mit bis zu 100% Rabatt.

Mehr erfahren.
Beschreibung

This updated edition introduces computational physics for performing experiments on the computer. It provides a grounding in scientific programming with examples in the context of physics problems. Using C++ as the primary language, it covers interpolation, integration, and solving differential equations, from simple concepts to advanced topics. The book includes a chapter on high-performance computing and parallel programming.
Understanding computational physics is crucial for data manipulation and real-world simulations. The book starts with introductory topics and progresses to advanced material, including a C++ library and sample programs. Readers will learn problem-solving methods and constructing models.
This journey equips readers with skills for scientific tasks, demonstrating efficient problem-solving and model construction. The practical approach ensures users can apply these concepts in real-world scenarios, making it an essential resource for those interested in computational physics.

Das E-Book können Sie in Legimi-Apps oder einer beliebigen App lesen, die das folgende Format unterstützen:

EPUB
MOBI

Seitenzahl: 464

Veröffentlichungsjahr: 2024

Bewertungen
0,0
0
0
0
0
0
Mehr Informationen
Mehr Informationen
Legimi prüft nicht, ob Rezensionen von Nutzern stammen, die den betreffenden Titel tatsächlich gekauft oder gelesen/gehört haben. Wir entfernen aber gefälschte Rezensionen.



COMPUTATIONALPHYSICS

LICENSE, DISCLAIMER OF LIABILITY, AND LIMITED WARRANTY

By purchasing or using this book (the “Work”), you agree that this license grants permission to use the contents contained herein, but does not give you the right of ownership to any of the textual content in the book or ownership to any of the information or products contained in it. This license does not permit uploading of the Work onto the Internet or on a network (of any kind) without the written consent of the Publisher. Duplication or dissemination of any text, code, simulations, images, etc. contained herein is limited to and subject to licensing terms for the respective products, and permission must be obtained from the Publisher or the owner of the content, etc., in order to reproduce or network any portion of the textual material (in any media) that is contained in the Work.

MERCURY LEARNINGAND INFORMATION (“MLI” or “the Publisher”) and anyone involved in the creation, writing, or production of the companion disc, accompanying algorithms, code, or computer programs (“the software”), and any accompanying Web site or software of the Work, cannot and do not warrant the performance or results that might be obtained by using the contents of the Work. The author, developers, and the Publisher have used their best efforts to insure the accuracy and functionality of the textual material and/or programs contained in this package; we, however, make no warranty of any kind, express or implied, regarding the performance of these contents or programs. The Work is sold “as is” without warranty (except for defective materials used in manufacturing the book or due to faulty workmanship).

The author, developers, and the publisher of any accompanying content, and anyone involved in the composition, production, and manufacturing of this work will not be liable for damages of any kind arising out of the use of (or the inability to use) the algorithms, source code, computer programs, or textual material contained in this publication. This includes, but is not limited to, loss of revenue or profit, or other incidental, physical, or consequential damages arising out of the use of this Work.

The sole remedy in the event of a claim of any kind is expressly limited to replacement of the book, and only at the discretion of the Publisher. The use of “implied warranty” and certain “exclusions” vary from state to state, and might not apply to the purchaser of this product.

COMPUTATIONALPHYSICS

SECOND EDITION

DARREN J. WALKER

Copyright ©2022 by MERCURY LEARNING AND INFORMATION LLC. All rights reserved.

Original title and copyright: Computational Physics: An Undergraduate’s Guide, 2/E. Copyright ©2021 by D.J. Walker. All rights reserved. Published by Pantaneto Press.

This publication, portions of it, or any accompanying software may not be reproduced in any way, stored in a retrieval system of any type, or transmitted by any means, media, electronic display or mechanical display, including, but not limited to, photocopy, recording, Internet postings, or scanning, without prior permission in writing from the publisher.

Publisher: David PallaiMERCURY LEARNING AND INFORMATION22841 Quicksilver DriveDulles, VA [email protected](800) 232-0223

D. J. Walker. Computational Physics, Second Edition.ISBN: 978-1-68392-832-4

The publisher recognizes and respects all marks used by companies, manufacturers, and developers as a means to distinguish their products. All brand names and product names mentioned in this book are trademarks or service marks of their respective companies. Any omission or misuse (of any kind) of service marks or trademarks, etc. is not an attempt to infringe on the property of others.

Library of Congress Control Number: 2021953014

222324321 Printed on acid-free paper in the United States of America

Our titles are available for adoption, license, or bulk purchase by institutions, corporations, etc. For additional information, please contact the Customer Service Dept. at 800-232-0223(toll free).

All of our titles are available in digital format at academiccourseware.com and other digital vendors. The sole obligation of MERCURY LEARNING AND INFORMATION to the purchaser is to replace the book or disc, based on defective materials or faulty workmanship, but not based on the operation or functionality of the product.

To Charlotte

CONTENTS

Chapter 1:Introduction

1.1Getting Started with Coding

1.2Getting To Know The Linux Command Line

1.3Bonjour Tout Le Monde

1.4The Rest of the Book

Chapter 2:Getting Comfortable

2.1Computers: What You Should Know

2.1.1Hardware

2.1.2Software

2.1.3Number Representation and Precision

2.2Some Important Mathematics

2.2.1Taylor Series

2.2.2Matrices: A Brief Overview

Exercises

Chapter 3:Interpolation and Data Fitting

3.1Interpolation

3.1.1Linear Interpolation

3.1.2Polynomial Interpolation

3.1.2Cubic Spline

3.2Data Fitting

3.2.1Regression: Illustrative Example

3.2.2Linear Least Squares: Matrix Form

3.2.3Realistic Example: Millikan’s Experiment

Exercises

Chapter 4:Searching for Roots

4.1Finding Roots

4.1.1Bisection

4.1.2Newton–Raphson

4.1.3Secant

4.2Hybrid Methods

4.2.1Bisection–Newton–Raphson

4.2.2Brute Force Search

4.3What’s The Point of Root Searching?

4.3.1The Infinite Square Well

4.3.2The Finite Square Well

4.3.3Programming the Root Finder

Exercises

Chapter 5:Numerical Quadrature

5.1Simple Quadrature

5.1.1The Mid-Ordinate Rule

5.1.2The Trapezoidal Rule

5.1.3Simpson’s Rule

5.2. Advanced Quadrature

5.2.1Euler–Maclaurin Integration

5.2.2Adaptive Quadrature

5.2.3Multidimensional Integration

Exercises

Chapter 6:Ordinary Differential Equations

6.1Classification of Differential Equations

6.1.1Types of Differential Equations

6.1.2Types of Solution and Initial Conditions

6.2Solving First-Order ODEs

6.2.1Simple Euler Method

6.2.2Modified and Improved Euler Methods

6.2.3The Runge–Kutta Method

6.2.4Adaptive Runge–Kutta

6.3Solving Second-Ordered ODEs

6.3.1Coupled 1st Order ODEs

6.3.2Oscillatory Motion

6.3.3More Than One Dimension

Exercises

Chapter 7:Fourier Analysis

7.1The Fourier Series

7.2Fourier Transforms

7.3The Discrete Fourier Transform

7.4The Fast Fourier Transform

7.4.1Brief History and Development

7.4.2Implementation and Sampling

Exercises

Chapter 8:Monte Carlo Methods

8.1Monte Carlo Integration

8.1.1Dart Throwing

8.1.2General Integration Using Monte Carlo

8.1.3Importance Sampling

8.2Monte Carlo Simulations

8.2.1Random Walk

8.2.2Radioactive Decay

Exercises

Chapter 9:Partial Differential Equations

9.1Classes, Boundary Values, and Initial Conditions

9.2Finite Difference Methods

9.2.1Difference Formulas

9.2.2Application of Difference Formulas

9.3Richardson Extrapolation

9.4Numerical Methods to Solve PDEs

9.4.1The Heat Equation with Dirichlet Boundaries

9.4.2The Heat Equation with Neumann Boundaries

9.4.3The Steady-State Heat Equation

9.4.4The Wave Equation

9.5Pointers To The Finite Element Method

Exercises

Chapter 10:Advanced Numerical Quadrature

10.1General Quadrature

10.2Orthogonal Polynomials

10.3Gauss–Legendre Quadrature

10.4Programming Gauss–Legendre

10.5Gauss–Laguerre Quadrature

Exercises

Chapter 11:Advanced ODE Solver and Applications

11.1Runge–Kutta–Fehlberg

11.2Phase Space

11.3Van Der Pol Oscillator

11.3.1Van der Pol in Phase Space

11.3.2Van der Pol FFT

11.4The “Simple” Pendulum

11.4.1Finite Amplitude

11.4.2Utter Chaos?

11.5Halley’s Comet

11.6To Infinity and Beyond

11.7To The Infinitesimal and Below

Exercises

Chapter 12:High-Performance Computing

12.1Indexing and Blocking

12.1.1Heap and Stack

12.1.2Computer Memory

12.1.3Loopy Indexing

12.1.4Blocking

12.1.5Loop Unrolling

12.2Parallel Programming

12.2.1Many (Hello) Worlds

12.2.2Vector Summation

12.2.3Overheads: Amdahl versus Gustafson

Exercises

Bibliography

Appendix:A Crash Course in C++ Programming

Index

CHAPTER 1

INTRODUCTION

Computational physics sits at the juncture of arguably three of the cornerstone subjects of modern times, physics, mathematics, and computer science. Many see it as sitting between theoretical physics, where there is a focus on mathematics and rigorous proof, and experimental physics, which is based on taking observations and quantitative measurements. The computational physicist performs numerical experimentation within the confines of the computer environment, applying mathematics to both simulate and examine complex models of physical systems. Just as the theoretician needs to master analytical mathematics, the experimentalist requires a working knowledge of laboratory apparatus, so does the computational physicist need to know about numerical analysis and computer programming. Any of these skills require (significant) practice to master but it is up to the physicist to know how to use them to interpret and, ultimately, understand the physical universe.

1.1 GETTING STARTED WITH CODING

You need two things to produce a computer program:

1.A text editor in which to write all your code in whatever language you choose.

2.A compiler to convert the code you have written into machine language (binary executable).

There are two methods by which you can write computer programs. The first method is via command-line control whereby you explicitly type in commands to compile a source code file written in a text editor. The second method uses what is called an Integrated Development Environment (IDE) that is essentially a compiler and text editor wrapped up into one neat application, for example, Microsoft’s Visual Studio. I would suggest trying out different text editors and IDEs to discover what suits you best. If your university uses Unix-based operating systems and you find it easier to code on those machines but do not want to splash out either on a Unix based machine (though the Raspberry Pi is reasonably priced) at home or make your Windows PC dual-booting (it can run either a Unix OS or Windows OS on one machine) an alternative is Cygwin. Cygwin creates a Unix type feel on a Windows PC and it’s free to download and install. Cygwin also comes with many different optional libraries and programs that are extremely useful to scientific programming, including the linear algebra package (LAPACK) library and Octave, a free alternative to MATLAB. If you can get your hands on a student version of MATLAB, I recommend you use it as it is a powerful programming tool and can be used to find quick programming solutions to problems, or as a first step towards a solution. A further alternative is to use a Virtual Machine.

For a list of freely available text editors just use your favorite search engine. Emacs is a popular programming text editor and is the default editor on most Unix-based machines; Cygwin also contains the GNU version of Emacs. On Windows you could use Notepad, however, it does not have any of the functionality of text editors specifically designed for coding. For example, programming languages have certain keywords reserved that have special meaning, for example, if, for, and while to name but a few. Once written these keywords are automatically distinguished from the rest of the text in some way, different color, different font, bolded, and so on. In Notepad all you will get is the same black text on a white background, which is not useful for reading and debugging the code you have written. Notepad++ is a good (and free) programming text editor for Windows that supports multiple languages.

If you prefer to use IDEs, there are a number available that are free to use. Some of these only support one language, for example, Dev C++, whereas others support multiple languages, for example, NetBeans, Code::Blocks or Eclipse. Microsoft do a “Community” version of their Visual Studio IDE which is free to use, as is the arguably more powerful Visual Studio Code.

Most of the code that accompanies this book has been written using the C++ programming language (C++ 11 onwards), using the Eclipse IDE for C/C++ Developers. I will not review the merits of the different programming languages here as the differences only really come into their own once you start to consider high-performance computing, Web applications, game programming, or other more specific applications. The basics of programming are sufficiently covered using just one language. That said, please be aware of different programming languages and how they can be used to produce different applications. For a challenge, you could convert the programs in this book into another language.

The next section gives a crash course in using the Linux command line.

1.2 GETTING TO KNOW THE LINUX COMMAND LINE

On modern operating systems a terminal emulator is a program that allows the use of the terminal in a graphical interface. In a Linux system, the shell is a command-line interface that interprets the user’s commands and passes them on to the underlying operating system. There are several shells in current use, such as the Bourne-Again shell (bash) or The C shell (tcsh), and each has its own set of features and operations, but all provide a means to interact with the machine.

When you open a new terminal emulator window the command prompt will be at the home directory (synonymous with “Folder” on Windows) of the current user. The information displayed in the command prompt is customizable by the user but typically consist of the user’s username, the host machine name, the current directory, and is ended by the prompt symbol. For an example of what this looks like please see Figure 1.1 that shows a macOS terminal.

FIGURE 1.1: Example of a Linux terminal emulator.

Commands can be issued after the command prompt by typing the name of an executable file, which can be a binary program or a script. There are many standard commands that are installed as default with the operating system that allows for system configuration, file system navigation, creation of new directories and files, installing third party programs, among other operations.

A useful command to start off with is pwd. It displays the full path to the current, working directory and can be useful if we ever get lost in the directory structure. The ls command will list, on the terminal, all the files and subdirectories of the current directory. Commands can also take arguments and options (or flags) that can affect their behavior. For instance, ls -l will nicely format the files and subdirectories with additional information such as attributes, permissions, sizes, and modification dates. The cd command is typically passed an argument of the directory to which we would like to navigate. For example, cd foo/bar will navigate to the subdirectory bar of the directory foo, assuming foo is a subdirectory of the current directory. The command cd alone will navigate us back to the user’s home directory. The Linux file system has two symbols reserved to represent the current directory and the parent directory namely “.” and “..,” respectively. Issuing cd . does not do much but cd .. will take us up one level into the parent directory. Using our example if we were in the bar subdirectory and we issued cd .. our current directory would now be foo. All these commands and more can be found on the infamous Linux Man pages (linux.die.net/man/ or man7.org/linux/man-pages/index.html) that outline all the possible options and arguments for these commands take. The man pages are daunting at first but once you learn how to read them offer a particularly useful resource when discovering or reusing various Linux commands.

The following table summarizes some of the more common commands you will likely use:

Command

Examples of use

mkdir

mkdir foo

creates a directory called foo in the current directory

mkdir foo/bar

creates a directory called bar in the directory foo

cd

cd

changes the current directory to the home directory

cd bar

changes the current directory to subdirectory bar

cd ..

changes the current directory to one level up

rmdir

rmdir foo/bar

removes directory bar from directory foo, if bar is empty

ls

ls

lists directories and files in current directory

ls foo

lists the directories and files in subdirectory foo

pwd

pwd

displays current location within the tree

touch

touch foo.log

if foo.log does not exist creates file “foo.log” in the current directory, else modifies the file’s timestamps

rm

rm foo.log

deletes the file “foo.log” in the current directory

1.3 BONJOUR TOUT LE MONDE

This section explains the basics of C++ syntax and structure. If you are already familiar with the C++ language, then please skip this section. A more in-depth introduction to the C++ language can be found in the Appendix.

The Hello World program is typically the first one that anyone learning a programming language gets to write. It gives us the basic syntax of a particular language and how to output something to the terminal or console. To begin let us first create a suitable directory structure to hold our code. Start by creating a “root” directory for our project called CompPhys in your home directory. Change into our newly created directory and create a subdirectory called Hello-World. Change into that directory and create the file helloworld.cpp. We now want to edit that file to fill it with the ground-breaking code for our hello world program.

Open the helloworld.cpp file in your text editor (or IDE) of choice and type the following:

// Helloworld program – displays message on stdout

#include <iostream>

int main () {

 std::cout << “Hello World” << std::endl;

 return 0;

}

Let us examine this line-by-line. At the very top, we have a comment line. These are either started using a double forward slash for a single-line comment or for a multiline (block) comment anything between “/*” and “*/” is treated as a comment. Comments are important. They should be used to explain the intention of code where this is not obviously apparent from the code itself. I urge you to use comments liberally; can you remember what you were doing yesterday, last week, last month, last year?

The next line down is how we include header files in source files. The hash symbol “#” is used to send instructions or directives to the pre-processor that is run before the compiler. In this case, we are using the directive include to insert the contents of the file iostream.h at this point in the source file, and that’s all it does. Note that iostream.h is a standard library header and as such we can drop the dot “h” extension from the filename when including it in source code. The angled brackets tell the preprocessor to look for the file in the standard external locations and implementation specified include directories only. A pair of double quotes around an included file tells the pre-processor to check the local directory of the source code first before checking the other locations. Generally, angled brackets are used for standard and system headers, and double quotes are used for programmer-defined headers.

White space is typically ignored by the compiler but is useful to humans by writing code that is easier to read. Having at least one blank line between the header include statement(s) and the start of the main function definition nicely separates the different parts of the source code.

The main function is the entry point to our program. It is this function that is called by the shell to process the code therein. As such it must be defined as a function that returns an integer value and the keyword int is that which represents the integer type in C and C++. The value of the returned integer indicates the “exit status” of the program; it is a convention that zero indicates success and that any other value indicates failure. The empty parenthesis tells the compiler that the main is a function declaration. Without them main would just be an integer variable declaration. We will discuss the different types of declarations later on.

The first statement in our Hello World program uses the C++ standard output stream object std::cout in conjunction with the output stream operator “<<” to send to the standard output (typically the terminal or console) the string sequence Hello World. Here I have introduced many new concepts namely streams, objects, and operators which we will get to in due course. All you need to know for the present is that this is how we can output data from a C++ program. Note that output from a program is stored in a data buffer and will only be printed on the terminal once that buffer is full, or the output is “flushed.” A flush means that a program will produce a line of output immediately. The std::endl is a stream “manipulator” that inserts a newline character into the output sequence and flushes the sequence.

As an aside the :: (double colon) is called the scope resolution operator. It is used to resolve the names of code elements (classes, functions, variables) that are contained in different namespaces. For example, cout and endl are scoped to the namespace “std,” short for standard. If we had omitted the namespace and scope resolution operator in the Hello World program, we would have received a compilation error on building the binary; remove the namespace and scope resolution operator from cout and endl to see the specific compilation error. The point of the namespace feature in C++ is to avoid naming clashes between various code elements, which is especially useful when developing large projects that might involve several third-party libraries. A more detailed discussion of the use of namespaces can be found in the Appendix.

The return statement defines the exit point of the function. In this case, we are returning the integer value zero to the calling environment (the shell) to indicate the successful execution of the program. At this point a program will flush all of its output buffers meaning any remaining data will be printed to the terminal; in our case, the buffer is already empty due to the use of the std::endl stream manipulator.

The code now has to be compiled into binary so that the computer can read and execute the instructions. In order to do this, we need a compiler. Most Linux distributions will come with Gnu’s C++ compiler installed as standard. To check that it is installed, and discover what version you have, type on the command line:

g++ -v

If the terminal output reads something along the lines of “g++ command not found” you will need to install the compiler; refer to your specific OS manual. Once you have confirmed the installation of the compiler, ensure you are in the same directory as your “helloworld” source file and type the following at the prompt:

g++ helloworld.cpp -o helloworld

All being well the compiler will have translated your source code contained within the helloworld.cpp file and produced a binary or executable file called “helloworld.” The “-o” flag tells the compiler to name the executable as the text you specify after the flag. If you do not give an output name, the file gets saved as “a.out” by default. To run this program type at the command prompt

./helloworld

and it will print to the command terminal the text “Hello World.” Note that if you are doing this in Windows system, binary executables are given a “.exe” extension.

We can now edit our source code to make the program more sophisticated, though only a little. Open the helloworld.cpp source file and modify the code as follows:

// Helloworld program – displays message on stdout

#include <iostream>

int main (int argc, char * argv[]) {

 if (argc < 2) {

std::cerr << “usage: ” << argv[0] << “ <name>”;

  return -1;

}

 std::cout << “Hello ” << argv[1] << std::endl;

 return 0;

}

Our program can now access the so-called command-line arguments, the parameters that are given on the command line after the program name when it is run. The argc variable contains the total number of command-line arguments, and the individual values of those arguments are available as strings accessed via the array argv. Note that the program name is counted as a command-line argument and is accessed as the first entry in the array, argv[0]; array indexing starts at zero for both C and C++. Compile this modified code and have a play with the command line arguments to gain an understanding of how it works. Note the use of a new stream object std::cerr to print an error message for when we do not pass sufficient arguments to the program. When a program is started by the shell it normally gains three open file descriptors: descriptor 0 is standard input, descriptor 1 is standard output, and descriptor 2 is the standard error. These descriptors are usually connected to the shell terminal that started the program but can be redirected to separate locations. As you may have guessed std::cout uses the standard output descriptor and std::cerr uses the standard error descriptor. The standard (input) stream object std::cin uses the standard input descriptor and can be used as an interactive means to obtain input from the user. Feel free to modify the Hello World program to use std::cin in some way; I find cppreference.com a useful resource.

What you have just done is a simple development cycle. Software development generally goes through three main phases:

1.Edit,

2.Compile, and

3.Execute.

Editing means writing or modifying the program contained in a text file using a text editor or IDE.

Compiling means translating the program to executable code, at this phase the code is checked for syntax errors and if found these are flagged for additional editing. We may also get linker errors at this stage.

Executing the program means running the code on your machine. At this phase, we may get logical errors or errors with the semantics (the meaning) of our code. These are known as runtime errors, and it is the process of debugging that removes these errors or bugs.

As you write computer programs in C++ and other languages, you will repeatedly edit, compile, and run programs. Sometimes the compiler will give you error messages. Often the messages can be quite cryptic. Just as trainee doctors learn best by giving known diseases to patients it is a good idea to make some deliberate errors, to get a feel for what the errors that you will encounter in the future might mean. This knowledge of deliberately inflicted “diseases” and observing their “symptoms” should help you diagnose future errors more effectively.

Errors in the semantics of your code can be quite insidious, especially if they cause undefined behavior; it is these runtime errors that we refer to as bugs. For instance, some piece of code complies with no syntax errors, linker errors, or any warnings. The code also runs with no logical errors and exits successfully. But the results you get are garbage. This will be highly likely due to a misunderstanding of the programmer as to how a particular feature of the language works, or a mistake in the use of a third-party library. Whatever the problem is, to fix the bug, you will first need to find its location within your code. To do this we step through the code line-by-line and observe how the data contained in the program variables change with each instruction. It is this stepping and observing a process that is called debugging. Most Linux distributions will come with the GNU debugger tool installed, invoked using gdb on the command line, and I highly recommend you familiarize yourself with its use. Note that this is an extremely brief exposition of debugging and to go into its details is beyond the scope of this book.

Just a quick word on nomenclature. A “program” can refer either to the source code contained in a text file or the binary executable file itself, they are semantically the same thing. The source code is readable by humans, the binary is readable by machines. An instance of a running program is called a “process” such that you can have several processes running the same program.

1.4 THE REST OF THE BOOK

Each chapter describes an individual topic within the general subject area of computational physics. Where there is a cross over between topics this has been explicitly referred to in the text. Throughout the remaining chapters, there are frequent references to C++ files that contain example programs for your study and use. These can be found on GitHub, an online repository for all sorts of different coding projects and applications, at the following URL: github.com/DJWalker42/laserRacoon.

These source code files come with GNU Makefiles such that compiling the code can be done by just typing “make” in the appropriate directory. These were developed on Mac OSX so contain variables specific to that operating system. You will have to modify some of the variables if you have a different OS. For more information on the GNU Makefile framework go to: gnu.org/software/make/manual/make.html.

The laserRacoon library makes use of OpenCV for a “visualization” module. If you do not have OpenCV installed on your system you can either use your install manager to get a copy (the library has been tested with OpenCV3) or visit their official site, opencv.org, for more options. If you cannot get OpenCV or would prefer not to use the visualization module then you will have to remove the related header and source files from the library (Visualise.h and Visualise.cpp) and remove any use of that module from the programs provided (anything using namespace phys::visual). Note that OpenCV is not really plotting software; OpenCV is an open-source library that performs image processing, video analysis, object and feature detection, camera calibration, 3D reconstruction, among other functions. At the time of writing the laserRacoon library, I needed a built-in way of visualizing the data being produced by the C++ programs. I had some experience of using OpenCV so challenged myself to make it able to plot data. The Viewer class does just that but know that it is not a fully optimized class and may contain bugs that have yet to be discovered (but that’s part of the fun of coding, right?).

Please note that the laserRacoon code has been designed by me and, as such, it should NOT be taken as gospel. I have taken the utmost care to make the classes, functions, and algorithms perform correctly but they have not been rigorously tested. It certainly could be redesigned to be optimized for performance or made more user friendly, I am not precious about it. Use it, abuse it, change it, that is how you learn. The code in the repository will only get updated if major bugs are found.

The code for the Fortran version of this book can also be found on the GitHub repository and may lend additional insight into the topics we discuss in this book. Indeed, I have not (directly) converted the Fortran code written for Chapter 9 on partial differential equations into C++ code; this is left as a challenge for the reader (and partly because we all have time constraints).

The chapters are arranged to provide some logical flow to the exploration of computational physics, starting out with the basic topics such as data fitting and root finding, and building to more advanced techniques, such as performing Fourier transforms and solving partial differential equations. At the end of each chapter are some exercises for the reader to do. These are designed to test you and to get you thinking like a physicist so do not be put off if you find them overly difficult at first. Use the resources available to you to find solutions, which includes fellow students, tutors, and professors, as well as that repository of all knowledge, the Internet—do not forget the library also.

In the Bibliography, you will find a guide to more general reading around each of the topics discussed, including pointers to other introductory texts in computational physics, the C++ programming language, and the Linux operating system. The Appendix contains a more thorough crash course in the C++ language.

CHAPTER 2

GETTING COMFORTABLE

2.1. COMPUTERS: WHAT YOU SHOULD KNOW

Computers are machines that help solve complex or tedious numerical problems. To make the hardware perform such tasks it must be programmed; in other words, told what to do. Remember a computer program cannot think by itself and is only as clever as the programmer who wrote the code. Understanding the underlying structure of a computer can help the programmer write smart code that takes advantage of that structure. For a comparison think about driving a car. You do not need to know how the car works at a component level to drive one. However, should you wish to improve the performance of the car, for racing, or rallying, or off-roading say, then you will have to know about the engine, the suspension, gearing, different types of tires and fuels, streamlining the bodywork, and so forth. This is no different for computers. Anyone can use a computer, but you really need to understand the details in order to get the most out of it.

2.1.1 Hardware

Due to the rapid advancements in computer technology, quantifying statements made in this section may well be out-of-date. However, the general qualifying remarks should still hold true (unless some paradigm-shifting technology has been invented since).

The physical elements that make up your computer is called hardware and consists of several components. The motherboard is the large, printed circuit board that contains all the ports, plugs, and electronics required to make the required components talk to each other. The central processing unit (CPU) handles most tasks in the computer. The speed at which the CPU handles these tasks is dependent on its clock frequency measured in Hertz. A 2 GHz single-core processor can handle at most 2 billion operations per second; operations may include additions, logic comparisons, and memory calls among others. Before 2004, clock frequencies were roughly doubling every 18 months. This followed the prediction made by Moore in the 1960s that the transistor density on silicon chips would double every 18 months. However, as the power consumed by the CPU goes up as the clock frequency squared, and with global concerns over energy usage, the frequency of the CPU is now capped at or around 4 GHz. The performance of computers has continued to increase according to Moore’s prediction using multiple core machines. At the date of writing the current commercially available state-of-the-art is 16 cores, with most “standard” computers having 4 cores, though that is rapidly changing to 8 cores. Multiple cores allow for parallel operation, whereby tasks can be handled simultaneously rather than having to be performed serially. For an introduction to parallel programming see Chapter 12 in this book.

For the CPU to be useful it must have a place to store information. Generally, there are several places for this information to be stored namely cache levels I and II (some CPUs have an additional third level of cache), random access memory (RAM), and storage either on a hard disk drive or in more modern systems on a solid state drive (SSD); from here onwards we will just refer to the storage device as such, or simply storage. This memory system has a hierarchical structure whereby the caches are the fastest but smallest memory levels and the storage is the largest but slowest memory level. Level I cache typically has a size of several tens of kilobytes (if not hundreds of kilobytes in modern CPUs) and can be accessed at the full processor speed. It is split into two separate areas, one for data and one for instructions, both required by the CPU to function. Level II cache has a typical size of several megabytes and can also be accessed at the full processor speed. Level II cache acts as a fast storage area for program code or variables required by that code. If the level II cache is filled by program code, then the overflow is put into RAM. RAM typically has a size of several gigabytes but is accessed at slower speeds than the caches. If the RAM is filled, then the CPU can store information on the storage device in a place called virtual memory. Storage devices today are immense coming in at relatively conservative 100 GB all the way up to 1 TB and beyond. However, the communication between CPU and virtual memory is limited by the speed at which data can be read from and written to the storage device. This speed of access can be a bottleneck for programs requiring large portions of memory; this is less true for modern PCIe/M.2 SSDs that can achieve around 1.5 GB/s read/write speed. Typically, memory considerations only come into play if it is dealt with high-definition images, video, or 3D graphics. However, some numerical methods can produce matrices of extremely large size that must be dealt with efficiently for a computer to produce timely (and accurate) results. The rise of the graphics card sometimes referred to as a graphical processing unit (GPU) has allowed for the development of some very sophisticated software without the need to use up CPU resources.

Other parts of a computer consist of input and output devices. Input devices are the things with which you communicate with the computer, for example, the keyboard and mouse. Output devices are how the computer communicates with the user, for example, the monitor and printer. Other devices can be considered as “slaves” being both controlled by the computer and relaying data back to the computer on command, for instance, a thermostat used to keep the room temperature constant.

2.1.2 Software

How do you make all that hardware do something? Computers are controlled using programs, referred to as software. The main program that is run on your computer is the operating system or OS. Mostly, Microsoft Windows OS of some version is used in the past; the latest version at the time of writing is Windows 10. Another widely available OS is UNIX which comes in various flavors. An undergraduate will almost certainly encounter a UNIX OS called Linux in their computer lab.

Programs are written in what is known as high-level languages, for example, Fortran and C++, and are compiled into machine language or binary via another program called a compiler. Note that programs such as MATLAB and Python are interpreted languages that are designed to run effortlessly on multiplatform machines (different OSs). Java is interesting in that it compiles sources into what is called “Byte Code,” files identified with the jar extension. Byte code can be interpreted by any machine that has the Java Runtime Environment installed.

Before you attempt any programming, please have a look at the following guidelines that may make your life easier:

1.Use your universities’ or work’s resources, which includes those sat next to you should you be in a computer lab or in your office. Failing an actual person who can communicate at least on some level, use the Internet. If you’ve got a complicated problem to solve it is very likely someone else has solved it already, and elegantly too (though never believe they managed it in one go without scratching their head at least once, drinking a lot of caffeine-based beverages, and swearing on several occasions). Try not to treat their solution as a black box that takes your inputs and gives the desired outputs without at least trying to understand what the code is doing. There is a practical limit to everyone’s knowledge and if it really makes no sense to accept that it works and that, out there, somewhere, is someone much cleverer than you and you’re ok with that.

2.Design your program first. Sit down, go to the old school with a pencil and paper, and write down the problem you are going to solve. What do you want to get as the output and what are going to be your inputs? Draw a flow chart if it helps. Write it out as pseudo-code; English phrases that mimic actual code and describe the program’s intended function line-by-line. This will, in the long run, save you time. Probably not straight away but practice makes perfect, allegedly. Now once this is done open your favorite text editor/IDE and start tapping away, but be aware (or beware) …

3.Make certain to comment on your code. As previously described, comments not only let others know what is intended with the code but also tell you what to be done. Comments should be clear, concise, descriptive, and written as understandable by the user. Do not worry if you have more comments than code.

4.Make names descriptive. This includes programs, classes, functions, and variables. If the program spanning is large as many hundreds of lines and/or several source files then it would be grateful to have given names that mean something. Additionally, many software companies will have their own naming convention for the various data types, structures, classes, and so forth that can be defined and declared in a program.

5.Do not be afraid to try something out. The worse thing that can happen is that your program crashes at runtime. Control-C starts over. Nowadays, it is very unfortunate to crash the entire computer but just turn it off and on again.

Some suggest that before anything else you should check that the problem you want to solve is suited to the use of a computer to avoid wasting your time and computer resources. While this is a helpful tip for experienced programmers (and who have several higher degrees in mathematics and physics), and it is arguable that only after getting into this habit, it will be comfortable to write a computer code. Sometimes the simple problems allow to explore writing novel and occasionally elegant or clever code that one may have missed trying to tackle a more complex problem.

2.1.3 Number Representation and Precision

As we are scientists, we will be dealing with real numbers obtained from measurements. During A-level physics course (or equivalent) teachers will likely have banged on about significant figures, rounding off, and the difference between precision and accuracy, when taking measurements from experiments. They would have found it bemusing that you quoted every figure on your calculator when figuring out, say, the strength of gravity at the Earth’s surface using a free-fall technique. Using a simple stopwatch to determine the time in free-fall and distanced traveled measured with a ruler the best you could hope to achieve is around three significant figures, limited by human reactions on the stopwatch. This precision could be increased using more precise equipment; for example, a computer control timing circuit measuring the distance using a laser. The point is that the precision of the results depends on the equipment used.

Computers can only express integer values exactly and are limited to a maximum integer value that can be expressed. Numbers in computers are stored as bits in binary format. A bit can have a logical value of 0 or 1, and strings of bits can be used to express integer numbers. A byte generally means a string of 8 bits, and 4 bytes, that is, composed of 32 bits is referred to as one word. Here, the binary format is referred as a big-endian, that is with the most significant bit written first at the left, as one would write decimal numbers. In contrast, little-endian puts the most significant bit at the end on the right, which is a natural format when performing binary addition, and the bits are in arithmetic order.

Take into consideration a byte or 8 bits. Each bit represents a power of two, starting at seven and ending with zero:

For example, the decimal number 6 would be represented by 0000 0110 in binary format; the equation governing this is

where k represents the bit location and s represents the bit value. Note that binary is easily read in 4-bit strings that the astute reader may notice leads naturally to the hexadecimal format − honest. The maximum integer can be expressed with 8 bits, that is . Note that numbers can be represented, one of which is zero hence the −1. This data type is known as an 8-bit, unsigned integer. Note that color images tend to be saved in this format with 8-bit unsigned integer values defining the three color channels (RGB) leading to the statement that color images have 16 million colors (256 × 256 × 256). A single channel image is generally referred to as a grayscale image, zero representing black, 255 representing white. Thus, there are 256 shades of grey but this makes less of a snappy title for a book.

Negative integers may also be expressed using this binary format. To do this the first bit (most significant bit) is taken as the sign bit. This now leaves us with only 7 bits to represent a number which gives a maximum number of. However, negative numbers can now be formed by taking the two’s compliment of a positive number. To do this you first form the one’s compliment by swapping the ones and zeros in the number, then add one to the result. For instance,

.

Zero is still represented by all zeros in the bit locations (take the two’s compliment of zero and you should still get zero). Given this conversion what is the largest negative number we can represent in this format? Take the largest positive number we can represent and take its two’s compliment:

However, note that the one’s compliment of +127 is available for use and, by definition, it is one less than the two’s compliment. Hence the largest negative number we can represent is −128. Note we have not lost any depth of numbers, we can still represent numbers; 128 negative numbers plus 127 positive numbers plus 1 for the zero.

Larger numbers can be represented using larger bit stings. A 32-bit word length can represent a maximum unsigned integer of or the signed integers in the range [].

Computational physics would be somewhat limited if computers could only use integer numbers. We need a way of representing floating-point decimal numbers. To do this, we take our 32-bit word length and split it into three blocks. Figure 2.1 illustrates this representation. The first block is one bit long and represents the sign of the number, 0 and 1 representing positive and negative values, respectively. The second block, typically 8 bits long represents the exponent, and the third block, containing the remaining 23 bits is the mantissa.

FIGURE 2.1: 32-bit representation of a floating-point number.

The most significant bit in the mantissa is on the left and represents The next bit represents and so forth. To calculate the floating-point decimal from the 32-bit representation we use the following equation

where s is the value of the sign bit, and the mantissa and exponent are the decimal values obtained from their respective binary format blocks. The bias is an implicit value that is included for the following reason. The 8-bit exponent does not contain an explicit sign bit and so can only represent positive numbers up to the maximum of 255. To circumvent this drastic limitation on floating-point number representation an implicit bias of 127 is included in the floating-point calculation. The range of exponents hence becomes []. The largest positive or negative number that can be represented is then approximately , and the smallest, not considering zero, is approximately . However, do not confuse this number as the computer’s precision. The computer’s precision is governed by the bit length of the mantissa; the exponent just defines the range of representable numbers.

The machine precision is best described in terms of how the computer performs floating-point arithmetic. Say you have the number 5 and wanted to add . Both numbers can be represented by the computer in floating-point notation, so far so good. To add them together the computer must match their exponents meaning that the bits in the mantissa of the smaller number get shifted to the right. By the time, the bits have been shifted to represent with the same exponent as 5 they have all gone past the least significant place and have been lost, in essence making equal to zero. The result of the addition would be 5.

The number has not just been plucked out of thin air. The least significant bit in the mantissa has a value of (2s.f.). This value represents a kind of number resolution; it is the smallest discernible difference between two numbers on a computer using a 32-bit word length. Note that it is a relative value; if you take the number then the next discernible number as far as the computer is concerned is 2000000.1. The machine epsilon or precision is the unit round-off error, essentially half the number resolution. For example, numbers in the range 2000000.000 to 2000000.049 would round down to 2000000.0, whereas numbers in the range 2000000.050 to 2000000.099 would round up to 2000000.1. Any result quoted from the computer should really include this rounding error, for example, . Because of the machine epsilon you should always consider whether the precision you are using is fit for purpose. If your calculations involve extreme differences between variable values, then unit round off may lead to large errors.

Clearly, the precision can be improved by adding more bits to the mantissa. This can be done by taking bits from the exponent but at the expense of the range of representable numbers. The other way of increasing the bit length of the mantissa is to double the word length from 32 to 64 bits. The standard format of a double data type is an 11-bit exponent and 52-bit mantissa plus the sign bit. What should the machine epsilon be using a double precision data type? To check you can write a few short lines of code to calculate the machine epsilon for both single and double precision variables. The pseudo-code for this task is written as follows:

Be wary that some compilers have been written to be smart and will try to “help” when producing the binary (executable) output. For instance, a C++ program written using the pseudo-code above gave a result that the single and double precision were both equal to , a precision of 64 bits. This clearly is incorrect. Changing the optimization flags one could managed to recover the expected result of for single precision and for double precision. The incorrect result is probably due to the compiler “helpfully’” converting the variables to extended precision that has a length of 80 bits, with a 64-bit mantissa. In any case, the initial result was clearly incorrect, and that brings us to an important point. Do not blindly accept what the computer outputs. If the answer looks wrong, then it most probably is wrong. When that guy in (A-level) physics class stated boldly that the strength of Earth’s gravity is two orders of magnitude larger than it is because that is what the calculator outputted, only to later realized that it would been using centimeters rather than meters.

2.2 SOME IMPORTANT MATHEMATICS

Physics describes the universe from tiniest sub-atomic particle to the shape of the universe itself. The language of physics is mathematics. However, do not confuse the two; physics is not the study of mathematics (and vice versa) but uses mathematics as a tool to describe and interpret the observation that we make of the universe. Nowhere is this truer than when dealing with computers that are, at the most basic level, efficient number crunching machines.

In this section, we will briefly review some fundamental mathematical concepts that are vital to any computational scientist performing numerical analysis.

2.2.1 Taylor Series

Brook Taylor was an English mathematician born in 1685 who devised an extremely useful way of approximating a function, the Taylor series expansion. This series expansion is arguably one of the most useful in mathematics and certainly within numerical analysis and will play a major role in much of the subject matter contained in this book.

The Taylor series is a mathematical technique for expressing a (potentially) complicated function in the form of a polynomial. The polynomial will have a similar value to the approximated function at least in some small neighborhood of a particular point. More precisely, a Taylor series is an infinite sum of power terms that represent a function at a single point. The summation terms are calculated from the values of the function’s derivatives at that point. Mathematically we write

(2.1)

where a is some point on x and we have used the notation that

.(2.2)

The last term in Equation (2.1) is the remainder or the error in the approximation where .

Usually, functions are approximated by using a finite number of terms of its Taylor series. Any finite number of initial terms of the Taylor series of a function is called a Taylor polynomial, the order of the polynomial governed by the highest power left in the approximation. For instance, a first ordered Taylor polynomial has the form

.(2.3)

Note the use of the approximately equals to sign as we have not included the remainder term here.

Please recollect the first ordered Taylor polynomial approximation as attempting to match the local neighborhood of the function at the point , using the function’s slope at that point. The second order polynomial, then, includes the curvature of the function at the point of interest. As more terms are added, higher ordered derivatives become utilized leading to a more accurate approximation of the function around the point of interest. Typically, the approximation is only usefully accurate over a closed interval about the point. A function that is equal to its Taylor series in an open interval is known as an analytic function. For instance, a straight-line function with some non-zero gradient would be given exactly by Equation (2.3) and is thus analytic.

The upper bound to the error in a Taylor polynomial can be estimated by analyzing the next term in the series from where we truncated the approximation. For example, take the Taylor series for the sine function taken about zero and truncated so that it is a seventh ordered polynomial approximation

. (2.4)

The upper bound to the error is then calculated by the next term in the series thus

. (2.5)

This upper bound is ignored further, higher ordered terms, which tend to improve the accuracy of the approximation, that is, to reduce the error.

The sine function and its seventh ordered Taylor polynomial are plotted in Figure 2.2. Here we can see the approximation is only reasonably (to the eye) accurate on the interval .

FIGURE 2.2: Seventh ordered Taylor polynomial approximation of the sine function.

2.2.2 Matrices: A Brief Overview

Matrices are incredibly important structures within mathematics, and thus within physics also. A very brief overview of their form and function were described in this section.

A matrix is an array of numbers. The dimensions of a matrix specify the number of rows and the number of columns the matrix has, in that order. Hence, when we say an n-by-m matrix we imply it has n rows and m columns. Vectors are essentially matrices of dimension n