# Title: Parallel simulation on four processors
#
# Description:
#
# The simulation is run in parallel on four processors. Load balancing
# is used to dynamically redistribute the elements across processors
# in order to maintain roughly equal mesh sizes on each processor as
# the resolution varies due to adaptive mesh refinement.
#
# \begin{figure}[htbp]
# \caption{\label{pid}MPEG movie of the processor number
# assigned to each cell (colour field) together with isolines of
# vorticity.}
# \begin{center}
# \htmladdnormallinkfoot{\includegraphics[width=\hsize]{pid.eps}}{pid.mpg}
# \end{center}
# \end{figure}
#
# \begin{figure}[htbp]
# \caption{\label{balance}Number of elements per processor as a
# function of time.}
# \begin{center}
# \includegraphics[width=0.8\hsize]{balance.eps}
# \end{center}
# \end{figure}
#
# Author: St\'ephane Popinet
# Command: sh parallel.sh
# Version: 090713
# Required files: parallel.sh pid.gfv
# Running time: 32 minutes
# Generated files: pid.mpg pid.eps balance.eps
#
8 7 GfsSimulation GfsBox GfsGEdge {} {

  # Stop the simulation at t = 15
  Time { end = 15 }

  # Insert the solid boundary defined as x*x + y*y - 0.0625*0.0625 = 0
  # (i.e. a cylinder of radius 0.0625 centered on the origin)
  Solid (x*x + y*y - 0.0625*0.0625)

  # Use an initial refinement of 6 levels (i.e. 2^6=64x64 for each
  # box) only around the solid boundary
  RefineSolid 6

  # Add a passive tracer called T
  VariableTracer {} T

  # Set the initial x-component of the velocity to 1
  Init {} { U = 1 }

  # Adapt the mesh using the vorticity criterion at every timestep
  # down to a maximum level of 6 and with a maximum tolerance of 1e-2
  AdaptVorticity { istep = 1 } { maxlevel = 6 cmax = 1e-2 }

  # Adapt the mesh using the gradient criterion on variable T at
  # every timestep, down to a maximum level of 6 and with a maximum tolerance of 1e-2
  AdaptGradient { istep = 1 } { maxlevel = 6 cmax = 1e-2 } T

  # Set a viscosity source term on the velocity vector
  # The Reynolds number is Re = D*U/Nu = 0.125*1/0.00078125 = 160
  # where D is the cylinder diameter (as defined in cylinder.gts)
  SourceViscosity 0.00078125

  # Balance the number of elements across parallel subdomains at every
  # timestep if the imbalance is larger than 0.1 (i.e. 10% difference
  # between the largest and smallest subdomains).
  EventBalance { istep = 1 } 0.1

  # Writes the time and timestep every 10 timesteps on standard error
  OutputTime { istep = 10 } stderr

  # Writes the time and simulation balance every timestep in 'balance'
  OutputTime { istep = 1 } balance
  OutputBalance { istep = 1 } balance

  # Writes info about the convergence of the Poisson solver on standard error
  OutputProjectionStats { istep = 10 } stderr

  # Save MPEG movie using GfsView module
  GModule gfsview
  OutputView { step = 0.05 } { 
      ppm2mpeg -s 800x100 > pid.mpg 
  } { width = 1600 height = 200 } pid.gfv

  # Outputs profiling information at the end of the simulation to standard error
  OutputTiming { start = end } stderr

  # Generate graphics
  OutputSimulation { start = end } end.gfs
  EventScript { start = end } {
      echo "Save pid.eps { format = EPS width = 800 height = 100 line_width = 0.2 }" | \
	  gfsview-batch2D pid.gfv end.gfs
      awk '{
        if ($1 == "step:")
          t = $4;
        else if ($1 == "domain")
          print t, 100.*($9/$3 - 1.), $3, $5, $9;
      }' < balance > balance1
      cat <<EOF | gnuplot
      set term postscript eps lw 3 solid 20 colour
      set output 'balance.eps'
      set xlabel 'Time'
      set ylabel 'Number of elements per processor'
      set key bottom right
      plot 'balance1' u 1:3 w l t 'minimum', '' u 1:4 w l t 'average', '' u 1:5 w l t 'maximum'
EOF
  }
}
GfsBox {
  # Left boundary on the leftmost box is:
  #   Dirichlet U=1 for the x-component of the velocity
  #   Dirichlet T = 1 if y < 0, 0 otherwise
  left = Boundary {
    BcDirichlet U 1
    BcDirichlet T { return y < 0. ? 1. : 0.; }
  }
}
GfsBox {}
GfsBox {}
GfsBox {}
GfsBox {}
GfsBox {}
GfsBox {}
# Right boundary on the rightmost box is outflow
GfsBox { right = BoundaryOutflow }
# All the boxes are linked by left to right links
1 2 right
2 3 right
3 4 right
4 5 right
5 6 right
6 7 right
7 8 right