Many additions to parallel course, a few to Fortran

content/courses/fortran-introduction/codes/formats.f90

@@ -0,0 +1,73 @@
+program printit
+   implicit none
+
+   integer, parameter :: rk=kind(1.0)
+   integer, parameter :: dk=kind(1.d0)
+
+   complex  :: z
+   real(rk) :: rpi=4.0*atan(1.0)
+   real(rk) :: dpi=4.d0*atan(1.d0)
+   real(rk) :: res
+   integer  :: i,j
+   integer  :: i1=42, i2=99
+   logical  :: is_zero, is_finite
+   character(len=12) :: mess1, mess2
+
+   real(rk), dimension(4)   :: b
+   real(rk), dimension(4,4) :: arr
+
+   res=cos(rpi)
+
+   do i=1,size(arr,1)
+      b(i)=real(i)
+      do j=1,size(arr,2)
+          arr(i,j)=real(i+j)
+      enddo
+   enddo
+
+   z=cmplx(res,rpi)
+
+   mess1="First result"
+   mess2="An integer"
+
+   is_zero=.false.
+   is_finite=.true.
+
+   write(*,'(i5,2x,i6)') i1,i2
+   write(*,'(i5,a,i6)') i1,"     ",i2
+   write(*,'(a,f0.6)') "The result is  ",res
+   write(*,'(a,i4,es15.7)') "The answer",i1,dpi
+   write(*,'(4(i2,3x,f8.3))') (j,b(j),j=1,4)
+   write(*,'(2f8.2)') z !complex
+   write(*,'(2L)')is_zero,is_finite
+   write(*,'(2p,f8.2,0p,f8.2)') rpi, dpi
+   write(*,'(a,f8.2,/,a,i6)') mess1,res,mess2,i1
+   write(*,'(a)') ' '
+   write(*,'("first value ",f8.2,", second value ",i6)') res
+   write(*,'(a)') ' '
+   write(*,'("first value ",f8.2,:," second value ",i6)') res
+
+   write(*,*) '----------------------'
+
+   do i=1,size(arr,1)
+      write(*,'(*(g0))') arr(i,:)
+   enddo
+
+   write(*,*) '----------------------'
+
+   do i=1,size(arr,1)
+      write(*,'(*(f12.4))') arr(i,:)
+   enddo
+
+! Old style
+   WRITE ( *, 2 ) i1
+   WRITE ( *, 2 ) i1,i2
+2  FORMAT ( 1X 'RESULT = ', I2, :, 3X, 'AT INDEX = ', I2 )
+
+   write(*,*) 
+
+   write(*,'(f8.2)',advance='no') res
+   write(*,'(i6)',advance='yes') i2
+   write(*,*) 'All done'
+
+end program

content/courses/fortran-introduction/formatted_io.md

@@ -19,15 +19,19 @@ Formatted output in Fortran is similar to other languages (the general layout de
 
 The edit descriptor modifies how to output the variables.  They are combined into forms like
 ```
-RaF.w
+RdF.w
 ```
-where `R` is a repeat count, `a` is the descriptor, `F` is the total field width _including_ space for `+-`, and if requested `+-e` and exponent, and `w` is the number of digits to the right of the decimal point. If you are willing to let the compiler calculate the number of characters to use, then `Ra0.w` alone works.  For floating-point numbers, `Ra.0` will print the integer part.
+where `R` is a repeat count, `d` is the descriptor, `F` is the total field width _including_ space for `+-`, and if requested `+-e` and exponent, and `w` is the number of digits to the right of the decimal point. If you are willing to let the compiler calculate the number of characters to use, use `Rd0.w`. You may also omit the number of decimal places with `Rd0` and the compiler will use its default for the type.
+
+For floating-point numbers, `Rd.0` will print the integer part. 
 
 Strings take only `RaF` and do not usually require the `F` since the length will be known to the compiler.
 
 Integers can be written as `iF` and any of the `F` spaces not needed will be blank filled, with the digits right justified.  When written as `iF.m` they will be printed with at least `m` digits and the rest of the field zero-filled on the left if all of `F` is not needed.
 
-If the field width is specified and the requested literal does not fit, the compiler will output a string of asterisks, e.g. `********`.
+If the field width is specified and the requested literal does not fit, the compiler will output a string of asterisks, e.g. `********`.  
+
+If you allow the compiler to compute the total field width, note that it will not include spaces before or after the item.
 
 ### Common Edit Descriptors
 
@@ -42,8 +46,10 @@ D  !double precision (prints D rather than E for exponent)
 A  !character (does not require a field width in most cases)
 X  !space
 /  !write an EOL and go to the next line (record) within the format
+:  !terminate the output if there are no more variables to write
 ```
 The real descriptors `F`, `E`, `G`, and `D` all work for both single and double precision. `G` allows the compiler to choose whether to use decimal or exponential format.
+
 The default exponential format writes in machine normalization, with the leading digit between 0 and 1. `ES` causes it to write with the leading digit between 1 and 9, which is what most humans can read most easily.  `ES` ignores `p` on output.
 
 ### Modifiers
@@ -68,28 +74,49 @@ For most purposes it is best to put the format string into the write statement.
 
 **Examples**
 ```fortran
-write(*,'(i5,2x,i6)') i1,i2
-write(*,'(i5,a,i6)')) i1,"  ",i2
-write(*,'(a,f0.6)')) "The result is  ",res
-write(*,'(a,i4,es15.7)') "row",n,var
-write(*,'(3(i2,3x,f8.3)') (r(j),var(j),j=1,3)
-write(*,'(2f8.2)') z !complex
-write(*,'(2L)')is_zero,is_finite
-write(*,'(2p,f8.2,0p,f8.2)') var1, var2
-write(*,'(a,f8.2,/,a,i6)') mess1,x,mess2,i
-```
-A format string may be a variable
+   write(*,'(i5,2x,i6)') i1,i2
+   write(*,'(i5,a,i6)') i1,"     ",i2
+   write(*,'(a,f0.6)') "The result is  ",res
+   write(*,'(a,i4,es15.7)') "The answer",i1,dpi
+   write(*,'(2p,f8.2,0p,f8.2)') rpi, dpi
+   write(*,'(a,f8.2,/,a,i6)') mess1,res,mess2,i1
+   write(*,'(a)') ' '
+   write(*,'("first value ",f8.2,", second value ",i6)') res
+   write(*,'(a)') ' '
+   write(*,'("first value ",f8.2,:," second value ",i6)') res
+```
+
+A format string may be a variable.
 ```fortran
 character(len=32) :: formatstr
    code
    formt='(f8.3,es15.7)'
    write(*,formatstr) A, B
 ```
 
+### Repetition
+
+Format strings can be repeated for multiple variables. If more than one descriptor is present, the format to be repeated should be enclosed in parentheses.
+```fortran
+   write(*,'(2L)')is_zero,is_finite
+   write(*,'(2f8.2)') z !complex
+   write(*,'4(f0.6)') a(1,:)
+   write(*,'(4(i2,3x,f8.3))') (j,b(j),j=1,4)
+```
+
+Especially when an array is allocatable, it may be awkward to specify the repeat count if it is unknown at compile time. List-directed I/O allows the compiler to choose to add end-of-line markers to line up columns, so the output can differ between different compilers and may not be what is desired for later processing. In principle a variable format string can be constructed with internal writes, but this can be complicated.  The Fortran 2008 standard introduced the `*` repetition count.  The compiler will repeat until it runs out of items. 
+
+```fortran
+   do i=1,size(arr,1)
+      write(*,'(*(f12.4))') arr(i,:)
+   enddo
+```
+
 ## Format Statements
 
 Format statements are abundant in older code, before the strings could be inserted into writes.
 FORMAT is non-executable but can appear anywhere in the source.  It is the only non-executable statement that can do so.
+
 It can still be useful for a particularly complex format (to keep the write statement short and readable) or for formats that are repeated in many write statements.
 The second parameter to the write is then an integer statement label.  The label marks the format statement.
 
@@ -134,7 +161,7 @@ read(quantity_char,*) quantity
 
 ## Fortran Non-Advancing IO
 
-* If we’d like to write to and read from standard input on the same line we can use non-advancing IO:
+Unlike most languages, Fortran `print` and `write` by default add an end-of-line marker at the end of the output. If we’d like to suppress this so that we can write multiple groups of output on the same line, or we would like to write some output to the console and read something from the console, we can use _non-advancing IO_.
 ```fortran
 write(*,'(a)',advance='no') "Enter input value:"
 read(*,*) value
@@ -146,7 +173,15 @@ Non-advancing IO _must_ be formatted
 
 **Exercises**
 
-1. Write a program that computes pi using a trig identity such as `pi=4*atan(1)`.
+1. Examine this example code:
+
+{{< spoiler text="A variety of formatting examples" >}}
+{{< code file="courses/fortran-introduction/code/formats.f90" lang="fortran" >}}
+{{< /spoiler >}}
+
+Make sure you understand all the formats. Correct unattractive outputs.  Experiment with changing the formats.
+
+2. Write a program that computes pi using a trig identity such as `pi=4*atan(1)`.
    * Use kind to switch between real and double precision
      * integer, parameter ::rk=kind(1.0)  or (1.0d0)
    * Using single precision, print pi in
@@ -160,7 +195,7 @@ Repeat for double precision.
 {{< /spoiler >}}
 
 
-2. In an “infinite” while loop:
+3. In an “infinite” while loop:
  Request an integer from the user with non-advancing input/output, e.g.
 ```
 "Please enter an integer:" <then read integer>

content/courses/fortran-introduction/project2.md

@@ -25,19 +25,7 @@ Print the result in Python format,
 
 A consultant was given a program that has code to generate a [format string](/courses/fortran-introduction/formatted_io) dynamically.  In particular, it can print items that might be arrays, with the repeat value generated automatically.  
 For example, given n1=5, n2=1, n3=3 and a pattern '(#e15.8,#i5,#f8.2)' 
-the result would be '(5e15.8,1i5,3f8.2)'. However, it didn’t work if any of the n1, n2, n3 variables were two digits.  
-The algorithm was convoluted and hard to understand.  It did work for two 
-digits if the programmer
-used two hash marks, e.g. ##e15.8, but that required hand-editing the several 
-files with output routines to find every place she wanted to write out more 
-than 9 array elements.  The author of the original code didn’t use any 
-character or string 
-functions other than substrings. This would surely be implemented more 
-generally with better use of strings.  Your task is to come up with a way to 
-do this.  
-If you have time, come up with a way to handle a 
-0 (i.e. skip printing).
-Test your program carefully.
+the result would be '(5e15.8,1i5,3f8.2)'. However, it didn’t work if any of the n1, n2, n3 variables were two digits. The algorithm was convoluted and hard to understand.  It did work for two digits if the programmer used two hash marks, e.g. ##e15.8, but that required hand-editing the several files with output routines to find every place she wanted to write out more than 9 array elements.  The author of the original code didn’t use any character or string functions other than substrings. This would surely be implemented more generally with better use of strings.  Your task is to come up with a way to do this.  If you have time, come up with a way to handle a 0 (i.e. skip printing).  Test your program carefully.
 
 _Hints_.  You do not need to use a variable-length string, but if not, be sure to declare a generously-sized set of strings (32 or even 64 characters, for example).  If using a fixed-length string, remember that you will need to remove blank space.  Test your program for n1=10, n2=2, n3=3.  Try another pattern like 
 '(#es12.4,#i2,#f15.7,#i4)'.  Suggestion: use an allocatable array for the coefficients (both numerical and character).  Use array size to make sure they match.    

content/courses/parallel-computing-introduction/codes/contour.py

@@ -0,0 +1,21 @@
+import sys
+import argparse
+import glob
+import numpy as np
+import pylab as plt
+
+parser = argparse.ArgumentParser()
+parser.add_argument("-f", "--fortran", help="Fortran ordering", action="store_true")
+parser.add_argument("filename", help="output file name")
+args = parser.parse_args()
+base = args.filename
+
+data=np.loadtxt(args.filename,unpack=False)
+
+if args.fortran:
+    data.T
+
+fig=plt.figure()
+plt.contourf(data)
+plt.colorbar()
+plt.show()

content/courses/parallel-computing-introduction/codes/contour_ranks.py

@@ -0,0 +1,37 @@
+import sys
+import argparse
+import glob
+import numpy as np
+import pylab as plt
+
+parser = argparse.ArgumentParser()
+parser.add_argument("-f", "--fortran", help="Fortran ordering", action="store_true")
+parser.add_argument("filename", help="base name of the output files")
+args = parser.parse_args()
+base = args.filename
+endpos=len(base)
+
+files=glob.glob(base+"*")
+numbers=[]
+suffixes=[]
+for file in files:
+    suffix=file[endpos:]
+    numbers.append(int(suffix))
+    suffixes.append(suffix)
+
+pvec = np.argsort(numbers)
+
+arrays=[]
+for n in pvec:
+    fname=files[n]
+    arrays.append(np.loadtxt(fname,unpack=False))
+
+if args.fortran:
+    data=np.concatenate(arrays,axis=1)
+else:
+    data=np.concatenate(arrays,axis=0)
+
+fig=plt.figure()
+plt.contourf(data)
+plt.colorbar()
+plt.show()

content/courses/parallel-computing-introduction/codes/exchange_deadlock.cxx

@@ -0,0 +1,39 @@
+#include <iostream>
+#include <mpi.h>
+
+using namespace std;
+
+int main(int argc, char **argv) {
+
+    int rank, nprocs, message, neighbor;
+    MPI_Status status;
+
+    MPI_Init(&argc, &argv);
+    MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
+    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
+
+    if (nprocs < 2) {
+        cout<<"This program works only for at least two processes\n";
+	MPI_Finalize();
+        return 1;
+    }
+    else if (nprocs%2 != 0) {
+        cout<<"This program works only for an even number of processes\n";
+	MPI_Finalize();
+        return 2;
+    }
+
+    if (rank%2==0) {
+        neighbor = rank+1;
+    }
+    else {
+	neighbor = rank-1;
+    }
+
+    MPI_Recv(&message, 1, MPI_INT, neighbor, 0, MPI_COMM_WORLD, &status);
+    MPI_Send(&rank, 1, MPI_INT, neighbor, 0, MPI_COMM_WORLD);
+    cout<<rank<<" "<<message<<endl;
+    MPI_Finalize();
+
+    return 0;
+}

content/courses/parallel-computing-introduction/codes/exchange_deadlock.f90

@@ -0,0 +1,37 @@
+program exchange
+   use mpi
+   implicit none
+
+   integer :: rank, nprocs, neighbor, ierr
+   integer :: status(MPI_STATUS_SIZE)
+   integer :: message
+
+   call MPI_Init(ierr)
+   call MPI_Comm_rank(MPI_COMM_WORLD, rank, ierr)
+   call MPI_Comm_size(MPI_COMM_WORLD, nprocs, ierr)
+
+   if (nprocs < 2) then
+        write(6,*) "This program works only for at least two processes"
+        call MPI_Finalize(ierr)
+        stop
+   else if ( mod(nprocs,2) /= 0 ) then
+        write(6,*) "This program works only for an even number of processes"
+        call MPI_Finalize(ierr)
+        stop
+   end if
+
+   if ( mod(rank,2) == 0) then
+       neighbor = rank+1
+   else
+       neighbor = rank-1
+   end if
+
+   call MPI_Recv(message, 1, MPI_INTEGER, neighbor, 0, MPI_COMM_WORLD, status, ierr)
+   call MPI_Send(rank,1, MPI_INTEGER, neighbor, 0, MPI_COMM_WORLD, ierr)
+
+   write(*,*) rank, message
+
+   call MPI_Finalize(ierr)
+
+end program
+

content/courses/parallel-computing-introduction/codes/exchange_deadlock.py

@@ -0,0 +1,27 @@
+import sys
+import numpy as np
+from mpi4py import MPI
+
+comm=MPI.COMM_WORLD
+nprocs=comm.Get_size()
+rank=comm.Get_rank()
+
+if nprocs<2:
+    print("This program works only for at least two processes.")
+    sys.exit()
+elif nprocs%2!=0:
+    print("This program works only for an even number of processes.")
+    sys.exit()
+
+message=np.zeros(1,dtype='int')
+rank_val=rank*np.ones(1,dtype='int')
+
+if rank%2==0:
+    neighbor=rank+1
+else:
+    neighbor=rank-1
+
+comm.Recv([message,MPI.INT],source=neighbor)
+comm.Send([rank_val,MPI.INT],dest=neighbor)
+
+print(rank,message)