[cig-commits] [commit] kommu/aspect-citcom-benchmarks: added benchmark files, Python scripts and workflow documentation to benchmarks folder (7d64cb0)

cig_noreply at geodynamics.org cig_noreply at geodynamics.org
Wed Nov 5 19:03:37 PST 2014


Repository : https://github.com/geodynamics/citcoms

On branch  : kommu/aspect-citcom-benchmarks
Link       : https://github.com/geodynamics/citcoms/compare/464e1b32299b15819f93efd98d969cddb84dfe51...f97ae655a50bdbd6dac1923a3471ee4dae178fbd

>---------------------------------------------------------------

commit 7d64cb0fdcaf86bd73986331b367a830e15dce71
Author: Rajesh Kommu <rajesh.kommu at gmail.com>
Date:   Fri Sep 5 16:19:57 2014 -0700

    added benchmark files, Python scripts and workflow documentation to benchmarks folder


>---------------------------------------------------------------

7d64cb0fdcaf86bd73986331b367a830e15dce71
 benchmarks/A1-benchmark                   |  11 +-
 benchmarks/A2-benchmark                   |  11 +-
 benchmarks/A3-benchmark                   |  25 ++-
 benchmarks/A4-benchmark                   |   9 +-
 benchmarks/{A3-benchmark => A5-benchmark} |  25 ++-
 benchmarks/{A3-benchmark => A6-benchmark} |  25 ++-
 benchmarks/{A3-benchmark => A7-benchmark} |  25 ++-
 benchmarks/A8-benchmark                   |  31 ++-
 benchmarks/{A3-benchmark => A9-benchmark} |  26 +--
 benchmarks/benchmark-workflow.pdf         | Bin 0 -> 133703 bytes
 benchmarks/benchmark-workflow.tex         | 336 ++++++++++++++++++++++++++++++
 benchmarks/citcoms-benchmark.py           |  70 +++++++
 benchmarks/table6.py                      |  99 +++++++++
 13 files changed, 595 insertions(+), 98 deletions(-)

diff --git a/benchmarks/A1-benchmark b/benchmarks/A1-benchmark
index 5302637..f14fc1a 100644
--- a/benchmarks/A1-benchmark
+++ b/benchmarks/A1-benchmark
@@ -1,15 +1,15 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=40001
+minstep=10001
 solver=full
 
 # CitcomS.controller
 checkpointFrequency=5000
-storage_spacing=100
+storage_spacing=10
 
 # CitcomS.solver
-datadir=A1-benchmark
-datafile=A1-benchmark
+datadir=/scratch/02204/rakuko/A1_09
+datafile=A1_09
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -48,7 +48,7 @@ perturbl=3
 perturblayer=16
 perturbm=2
 perturbmag=0.01
-tic_method=0
+tic_method=3
 
 # CitcomS.solver.output
 output_ll_max=40
@@ -64,4 +64,3 @@ rheol=1
 visc0=1,1,1,1
 viscE=0.0,0.0,0.0,0.0
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A2-benchmark b/benchmarks/A2-benchmark
index f33d0b0..8f699b5 100644
--- a/benchmarks/A2-benchmark
+++ b/benchmarks/A2-benchmark
@@ -1,15 +1,15 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=15001
 solver=full
 
 # CitcomS.controller
 checkpointFrequency=5000
-storage_spacing=100
+storage_spacing=50
 
 # CitcomS.solver
-datadir=A2-benchmark
-datafile=A2-benchmark
+datadir=/scratch/02204/rakuko/A2_02
+datafile=A2_02
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -48,7 +48,7 @@ perturbl=3
 perturblayer=16
 perturbm=2
 perturbmag=0.01
-tic_method=0
+tic_method=3
 
 # CitcomS.solver.output
 output_ll_max=40
@@ -64,4 +64,3 @@ rheol=1
 visc0=1,1,1,1
 viscE=2.3026,2.3026,2.3026,2.3026
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A3-benchmark b/benchmarks/A3-benchmark
index f1c0f16..dc5f171 100644
--- a/benchmarks/A3-benchmark
+++ b/benchmarks/A3-benchmark
@@ -1,15 +1,15 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=11001
 solver=full
 
 # CitcomS.controller
 checkpointFrequency=5000
-storage_spacing=100
+storage_spacing=50
 
 # CitcomS.solver
-datadir=A3-benchmark
-datafile=A3-benchmark
+datadir=/scratch/02204/rakuko/A3_07
+datafile=A3_07
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -29,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -50,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -64,4 +64,3 @@ rheol=1
 visc0=1,1,1,1
 viscE=3.0,3.0,3.0,3.0
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A4-benchmark b/benchmarks/A4-benchmark
index f7e23df..8002f02 100644
--- a/benchmarks/A4-benchmark
+++ b/benchmarks/A4-benchmark
@@ -1,6 +1,6 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=75001
+minstep=30001
 solver=full
 
 # CitcomS.controller
@@ -8,8 +8,8 @@ checkpointFrequency=5000
 storage_spacing=100
 
 # CitcomS.solver
-datadir=A4-benchmark
-datafile=A4-benchmark
+datadir=/scratch/02204/rakuko/A4_02
+datafile=A4_02
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -48,7 +48,7 @@ perturbl=3
 perturblayer=16
 perturbm=2
 perturbmag=0.01
-tic_method=0
+tic_method=3
 
 # CitcomS.solver.output
 output_ll_max=40
@@ -64,4 +64,3 @@ rheol=1
 visc0=1,1,1,1
 viscE=4.6052,4.6052,4.6052,4.6052
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A3-benchmark b/benchmarks/A5-benchmark
similarity index 88%
copy from benchmarks/A3-benchmark
copy to benchmarks/A5-benchmark
index f1c0f16..52ea7d7 100644
--- a/benchmarks/A3-benchmark
+++ b/benchmarks/A5-benchmark
@@ -1,6 +1,6 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=31001
 solver=full
 
 # CitcomS.controller
@@ -8,8 +8,8 @@ checkpointFrequency=5000
 storage_spacing=100
 
 # CitcomS.solver
-datadir=A3-benchmark
-datafile=A3-benchmark
+datadir=/scratch/02204/rakuko/A5_01
+datafile=A5_01
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -29,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -50,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -62,6 +62,5 @@ VISC_UPDATE=on
 num_mat=4
 rheol=1
 visc0=1,1,1,1
-viscE=3.0,3.0,3.0,3.0
+viscE=6.9078,6.9078,6.9078,6.9078
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A3-benchmark b/benchmarks/A6-benchmark
similarity index 88%
copy from benchmarks/A3-benchmark
copy to benchmarks/A6-benchmark
index f1c0f16..4d4ad81 100644
--- a/benchmarks/A3-benchmark
+++ b/benchmarks/A6-benchmark
@@ -1,6 +1,6 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=110001
 solver=full
 
 # CitcomS.controller
@@ -8,8 +8,8 @@ checkpointFrequency=5000
 storage_spacing=100
 
 # CitcomS.solver
-datadir=A3-benchmark
-datafile=A3-benchmark
+datadir=/scratch/02204/rakuko/A6_01
+datafile=A6_01
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -29,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -50,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -62,6 +62,5 @@ VISC_UPDATE=on
 num_mat=4
 rheol=1
 visc0=1,1,1,1
-viscE=3.0,3.0,3.0,3.0
+viscE=9.2103,9.2103,9.2103,9.2103
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A3-benchmark b/benchmarks/A7-benchmark
similarity index 88%
copy from benchmarks/A3-benchmark
copy to benchmarks/A7-benchmark
index f1c0f16..8b6ae17 100644
--- a/benchmarks/A3-benchmark
+++ b/benchmarks/A7-benchmark
@@ -1,6 +1,6 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=70001
 solver=full
 
 # CitcomS.controller
@@ -8,8 +8,8 @@ checkpointFrequency=5000
 storage_spacing=100
 
 # CitcomS.solver
-datadir=A3-benchmark
-datafile=A3-benchmark
+datadir=/scratch/02204/rakuko/A7_01
+datafile=A7_01
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -29,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -50,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -62,6 +62,5 @@ VISC_UPDATE=on
 num_mat=4
 rheol=1
 visc0=1,1,1,1
-viscE=3.0,3.0,3.0,3.0
+viscE=11.5129,11.5129,11.5129
 viscT=0.5,0.5,0.5,0.5
-
diff --git a/benchmarks/A8-benchmark b/benchmarks/A8-benchmark
index c6016db..d24dbb1 100644
--- a/benchmarks/A8-benchmark
+++ b/benchmarks/A8-benchmark
@@ -1,17 +1,16 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=76000
+minstep=76001
 solver=full
 
 # CitcomS.controller
 checkpointFrequency=5000
-storage_spacing=500
+storage_spacing=50
 
 # CitcomS.solver
-datadir=A8-benchmark
-datafile=A8-benchmark
-#rayleigh=7.681755e4
-rayleigh=7e3
+datadir=/scratch/02204/rakuko/A8_04
+datafile=A8_04
+rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
 fi_max=1
@@ -30,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -51,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -63,6 +62,6 @@ VISC_UPDATE=on
 num_mat=4
 rheol=1
 visc0=1,1,1,1
-viscE=13.82,13.82,13.82,13.82 
+viscE=13.816,13.816,13.816,13.816
 viscT=0.5,0.5,0.5,0.5
-
+CDEPV=off
diff --git a/benchmarks/A3-benchmark b/benchmarks/A9-benchmark
similarity index 86%
copy from benchmarks/A3-benchmark
copy to benchmarks/A9-benchmark
index f1c0f16..cc1299c 100644
--- a/benchmarks/A3-benchmark
+++ b/benchmarks/A9-benchmark
@@ -1,6 +1,6 @@
 # CitcomS
 cpu_limits_in_seconds=360000000
-minstep=50001
+minstep=65001
 solver=full
 
 # CitcomS.controller
@@ -8,8 +8,8 @@ checkpointFrequency=5000
 storage_spacing=100
 
 # CitcomS.solver
-datadir=A3-benchmark
-datafile=A3-benchmark
+datadir=/scratch/02204/rakuko/A9_01
+datafile=A9_01
 rayleigh=7.681755e4
 
 # CitcomS.solver.mesher
@@ -29,14 +29,6 @@ mgunitx=2
 mgunity=2
 mgunitz=2
 
-# CitcomS.solver.ic
-num_perturbations=1
-perturbl=3
-perturblayer=16
-perturbm=2
-perturbmag=0.01
-tic_method=0
-
 # CitcomS.solver.vsolver
 Solver=multigrid
 precond=on
@@ -50,6 +42,14 @@ vhighstep=3
 max_mg_cycles=200
 piterations=1500
 
+# CitcomS.solver.ic
+num_perturbations=1
+perturbl=3
+perturblayer=16
+perturbm=2
+perturbmag=0.01
+tic_method=3
+
 # CitcomS.solver.output
 output_ll_max=40
 output_optional=geoid,surf,botm,horiz_avg,volume_avg
@@ -62,6 +62,6 @@ VISC_UPDATE=on
 num_mat=4
 rheol=1
 visc0=1,1,1,1
-viscE=3.0,3.0,3.0,3.0
+viscE=16.1181,16.1181,16.1181,16.1181
 viscT=0.5,0.5,0.5,0.5
-
+CDEPV=off
diff --git a/benchmarks/benchmark-workflow.pdf b/benchmarks/benchmark-workflow.pdf
new file mode 100644
index 0000000..a095d3b
Binary files /dev/null and b/benchmarks/benchmark-workflow.pdf differ
diff --git a/benchmarks/benchmark-workflow.tex b/benchmarks/benchmark-workflow.tex
new file mode 100644
index 0000000..700e540
--- /dev/null
+++ b/benchmarks/benchmark-workflow.tex
@@ -0,0 +1,336 @@
+\documentclass[10pt]{article}
+%\usepackage{fourier}
+\usepackage[margin=1in]{geometry}
+\usepackage{minted}
+\usemintedstyle{friendly}
+
+\title{Workflow for generating CitcomS benchmark results}
+\author{}
+
+\begin{document}
+\maketitle
+\section*{Code changes to CitcomS}
+\subsection*{Calculating the volume average quantities}
+The following function has been added to \texttt{\textbf{Process\_buoyancy.c}} to compute
+$\langle T\rangle$ and $\langle V_{rms}\rangle$
+\begin{minted}{c}
+/*
+ * compute the volume average of temperature and RMS velocity 
+ */
+void compute_volume_avg(struct All_variables *E, double *T_avg, double *Vrms_avg)
+{
+  int m, n, i;
+  double *S1[NCS], *S2[NCS];
+
+  for(m=1; m<=E->sphere.caps_per_proc; m++)
+  {
+    S1[1] = (double *)malloc((E->lmesh.nno+1)*sizeof(double));
+    S2[1] = (double *)malloc((E->lmesh.nno+1)*sizeof(double));
+  }
+
+  for(m=1; m<=E->sphere.caps_per_proc; m++)
+  {
+    for(i=1; i<=E->lmesh.nno+1; i++)
+    {
+      S1[m][i] = E->T[m][i];
+      S2[m][i] = 
+        E->sphere.cap[m].V[1][i]*E->sphere.cap[m].V[1][i] +
+        E->sphere.cap[m].V[2][i]*E->sphere.cap[m].V[2][i] +
+        E->sphere.cap[m].V[3][i]*E->sphere.cap[m].V[3][i];
+    }
+  }
+  
+  *T_avg = return_bulk_value_d(E, S1, 1); /* 1 => need volume average */
+  *Vrms_avg = sqrt(return_bulk_value_d(E, S2, 1)); /* 1 => need volume average */
+
+  for(m=1;m<=E->sphere.caps_per_proc;m++) {
+    free((void *)S1[m]);
+    free((void *)S2[m]);
+  }
+}
+\end{minted}
+
+The following function has been added to \textbf{\texttt{Output.c}}
+
+\begin{minted}{c}
+void output_volume_avg(struct All_variables *E, int cycles)
+{
+  /* volume average output of temperature and rms velocity */
+  void compute_volume_avg();
+
+  int j;
+  char output_file[255];
+  FILE *fp1;
+  double T_avg=0.0, V_rms_avg=0.0;
+
+  /* compute horizontal average here.... */
+  compute_volume_avg(E, &T_avg, &V_rms_avg);
+
+  if (E->parallel.me == 0)  
+  {
+    sprintf(output_file,"%s.volume_avg", E->control.data_file);
+    fp1=fopen(output_file,"a+");
+    fprintf(fp1,"%d %.4e %.4e %.4e\n", 
+	    cycles, E->monitor.elapsed_time, T_avg, V_rms_avg);
+    fclose(fp1);
+  }
+}
+\end{minted}
+
+\subsection*{Temperature Initial Conditions for the B Benchmarks}
+The function \textbf{\texttt{construct\_tic\_from\_input}} has been modified to set up the temperature initial 
+condition for the B benchmarks. This modified code is called by setting \textbf{\texttt{tic\_method=6}} in the 
+parameter file
+
+\begin{minted}{c}
+static void construct_tic_from_input(struct All_variables *E)
+{
+    double mantle_temperature;
+
+    switch (E->convection.tic_method){
+    .
+    .
+    .
+    case 6:
+      /* a conductive temperature profile + perturbations at all layers */
+      conductive_temperature_profile(E);
+      add_perturbations_at_all_layers_B(E);
+      break;
+    .
+    .
+    .
+    }
+}
+\end{minted}
+
+The new function \texttt{\textbf{add\_perturbations\_at\_all\_layers\_B}} is defined as 
+\begin{minted}{c}
+static void add_perturbations_at_all_layers_B(struct All_variables *E)
+{
+  /*
+   * This method generates the initial temperature profile according to 
+   * equation (47) of the Zhong et. al. 2008 paper. (4,0)+(4,4) I.C. is
+   * hard coded
+   */
+
+    int m, i, j, k, node;
+    int p;
+    int nox, noy, noz, gnoz;
+    double r1, t1, f1, tlen, flen, rlen, con;
+
+    nox = E->lmesh.nox;
+    noy = E->lmesh.noy;
+    noz = E->lmesh.noz;
+    gnoz = E->mesh.noz;
+
+    rlen = M_PI / (E->sphere.ro - E->sphere.ri);
+
+    for (p=0; p<E->convection.number_of_perturbations; p++) {
+        con = E->convection.perturb_mag[p];
+
+        if (E->parallel.me_loc[1] == 0 && E->parallel.me_loc[2] == 0
+            && E->sphere.capid[1] == 1 )
+            fprintf(stderr,"Initial temperature perturbation:  mag=%g\n", con);
+
+        if(E->sphere.caps == 1) {
+          myerror(E, "add_pertubrbations_at_all_layers_B is not implemented for the Regional model");
+        }
+        else {
+            /* global mode, add spherical harmonics perturbation */
+
+            for(m=1; m<=E->sphere.caps_per_proc; m++)
+                for(i=1; i<=noy; i++)
+                    for(j=1; j<=nox;j ++)
+                        for(k=1; k<=noz; k++) {
+                            node = k + (j-1)*noz + (i-1)*nox*noz;
+                            t1 = E->sx[m][1][node];
+                            f1 = E->sx[m][2][node];
+                            r1 = E->sx[m][3][node];
+
+                            E->T[m][node] += con * (modified_plgndr_a(4,0,t1) + 
+                                (5.0/7.0)*cos(4*f1)*modified_plgndr_a(4,4,t1))
+                                * sin((r1-E->sphere.ri) * rlen);
+                        }
+        } /* end if */
+    } /* end for p */
+
+    return;
+}
+\end{minted}
+
+\section*{Files needed for generating the plots and tables}
+The following four files are needed for generating the plots and tables. The ``A1\_10'' prefix is the value of the 
+\textbf{\texttt{datafile}} parameter from the configuration file. In this particular example, we have
+\textbf{\texttt{nprocz=2}}, hence there are only 2 \textbf{\texttt{*.horiz\_avg.*}} files. Also, 10000 is the 
+\textbf{last time step} for which we have the recorded data.
+
+\begin{description}
+\item[A1\_10.o3993348] The output logfile generated when CitcomS is run on Stampede (or some other cluster). This file
+has information for the surface and bottom heat flux.
+\item[A1\_10.horiz\_avg.0.10000] The horizontal average file for the 0th z process, at the 10000th time step. 
+\item[A1\_10.horiz\_avg.1.10000] The horizontal average file for the 1st z process, at the 10000th time step.
+\item[A1\_10.volume\_avg] The volume average file that has the numbers for $\langle T\rangle$ and 
+$\langle V_{rms}\rangle$
+\end{description}
+
+\subsection*{Generating the data files}
+\subsubsection*{Surface heat flux}
+\begin{minted}[mathescape]{bash}
+cat A1_10.o3993348 | sed -n 's/surface heat flux= //p' | awk 'NR%2==0' >> surface-heat-flux-A1_10
+\end{minted}
+\subsubsection*{Extract the horizontally averaged $r$ values}
+\begin{minted}[mathescape]{bash}
+awk '{print $1}' A1_10.horiz_avg.0.10000 >> R_vals_0
+\end{minted}
+\begin{minted}[mathescape]{bash}
+awk '{print $1}' A1_10.horiz_avg.1.10000 >> R_vals_1
+\end{minted}
+\begin{minted}[mathescape]{bash}
+cat R_vals_0 R_vals_1 >> A1_10_r
+\end{minted}
+\textbf{\texttt{A1\_10\_r}} has duplicate entries; the last line from \textbf{\texttt{R\_vals\_0}} and the first line
+from \textbf{\texttt{R\_vals\_1}}. This needs to be fixed manually.
+\subsubsection*{Extract the horizontally averaged $T$ values}
+\begin{minted}[mathescape]{bash}
+awk '{print $2}' A1_10.horiz_avg.0.10000 >> T_vals_0
+\end{minted}
+\begin{minted}[mathescape]{bash}
+awk '{print $2}' A1_10.horiz_avg.1.10000 >> T_vals_1
+\end{minted}
+\begin{minted}[mathescape]{bash}
+cat T_vals_0 T_vals_1 >> A1_10_T
+\end{minted}
+\textbf{\texttt{A1\_10\_T}} has duplicate entries; the last line from \textbf{\texttt{T\_vals\_0}} and the first line
+from \textbf{\texttt{T\_vals\_1}}. This needs to be fixed manually.
+
+\subsection*{Generating the Plots}
+The Python script \textbf{\texttt{citcoms-benchmark.py}} can be used for generating figures 5(a), 5(b), 5(c) and
+figure 7(a) of the Zhong et. al. 2008 paper. You need to generate the files for the A3 and A8 benchmarks, using the
+procedure for A1 described above, and correctly specify the path to the data files.
+
+\begin{minted}[mathescape]{python}
+'''
+citcoms-benchmark.py
+'''
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib.ticker as ticker
+
+d=0.45
+kappa=1.0
+time_rescale = (d*d)/kappa # t_citcom / time_rescale = t_nondim (or t_paper) 
+vel_rescale = kappa / d    # v_citcom / vel_rescale = v_nondim (i.e. v_paper)
+r_t = 1.0
+r_b = 0.55
+top_prefac = r_t*(r_t-r_b)/r_b;
+
+dataA1 = np.loadtxt("/home/rkk/A1_10/A1_10.volume_avg")
+dataA3 = np.loadtxt("/home/rkk/A3_08/A3_08.volume_avg")
+dataA8 = np.loadtxt("/home/rkk/A8_05/A8_05.volume_avg")
+surfA1 = top_prefac*np.loadtxt("/home/rkk/A1_10/surface-heat-flux-A1_10")
+surfA3 = top_prefac*np.loadtxt("/home/rkk/A3_08/surface-heat-flux-A3_08")
+surfA8 = top_prefac*np.loadtxt("/home/rkk/A8_05/surface-heat-flux-A8_05")
+timeA1 = (1.0/time_rescale)*dataA1[:,1]
+timeA3 = (1.0/time_rescale)*dataA3[:,1]
+timeA8 = (1.0/time_rescale)*dataA8[:,1]
+T_avgA1 = dataA1[:,2]
+T_avgA3 = dataA3[:,2]
+T_avgA8 = dataA8[:,2]
+Vrms_avgA1 = (1.0/vel_rescale)*dataA1[:,3]
+Vrms_avgA3 = (1.0/vel_rescale)*dataA3[:,3]
+Vrms_avgA8 = (1.0/vel_rescale)*dataA8[:,3]
+
+plt.figure(1)
+# plt.subplot(3,1,1)
+plt.plot(timeA1,T_avgA1, 'g-', label='A1')
+plt.plot(timeA3,T_avgA3, 'r-', label='A3')
+plt.plot(timeA8,T_avgA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('<T>')
+plt.title('Time dependence of Volume-averaged temperature')
+plt.ylim((0,0.8))
+plt.xlim((0,1.0))
+plt.yticks([0.0,0.2,0.4,0.6,0.8])
+
+plt.figure(2)
+# plt.subplot(3,1,2)
+plt.plot(timeA1,Vrms_avgA1, 'g-', label='A1')
+plt.plot(timeA3,Vrms_avgA3, 'r-', label='A3')
+plt.plot(timeA8,Vrms_avgA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('<V_rms>')
+plt.title('Time dependence of Root-mean squared velocity')
+plt.ylim((0,120.0))
+plt.xlim((0,1.0))
+plt.yticks(np.arange(0.0,121.0,20.0))
+
+plt.figure(3)
+# plt.subplot(3,1,3)
+plt.plot(timeA1,surfA1, 'g-', label='A1')
+plt.plot(timeA3,surfA3, 'r-', label='A3')
+plt.plot(timeA8,surfA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('Nu_t')
+plt.title('Time dependence of surface Nusselt number')
+plt.ylim((0,6.0))
+plt.xlim((0,1.0))
+
+plt.show()
+\end{minted}
+
+\subsection*{Calculating mean and standard deviations}
+The Python script \textbf{\texttt{table6.py}} can be used for computing the mean and standard deviations of the
+values reported in table 6 and 7 of the paper. Again, the path to the data files needs to be set appropriately
+\begin{minted}[mathescape]{python}
+#-----------------------------------------------------------------------------------------------------------------------
+# table6.py : Generate the entries for Table 6 of the Zhong et. al. paper
+#-----------------------------------------------------------------------------------------------------------------------
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib.ticker as ticker
+
+d=0.45
+kappa=1.0
+time_rescale = (d*d)/kappa # t_citcom / time_rescale = t_nondim (or t_paper) 
+vel_rescale = kappa / d    # v_citcom / vel_rescale = v_nondim (i.e. v_paper)
+r_t = 1.0
+r_b = 0.55
+top_prefac = r_t*(r_t-r_b)/r_b;
+
+def get_array_slice(arr, low, high):
+  '''
+  returns a boolean array with True for arr in[low,high]
+  '''
+
+  lower_limit = np.greater_equal(arr, low)
+  upper_limit = np.less_equal(arr, high)
+  combined = np.logical_and(lower_limit, upper_limit)
+  return combined
+
+def compute_and_print_mean_and_std(benchmark, volume_avg_file, heat_flux_file, tlow, thigh):
+  data = np.loadtxt(volume_avg_file)
+  surf = top_prefac * np.loadtxt(heat_flux_file)
+  time = (1.0/time_rescale) * data[:,1]
+  T_avg = data[:,2]
+  Vrms_avg = (1.0/vel_rescale)*data[:,3]
+  combined = get_array_slice(time, tlow, thigh)
+  sub_Tavg = T_avg[combined]
+  sub_Vrms = Vrms_avg[combined]
+  sub_Nu_t = surf[combined]
+
+  print("----------------- %s Benchmark Results -----------------------" % benchmark)
+  print("<T> = %8.4f std = %e" % (np.mean(sub_Tavg), np.std(sub_Tavg)))
+  print("<V_rms> = %8.4f std = %e" % (np.mean(sub_Vrms), np.std(sub_Vrms)))
+  print("<Nu_t> = %8.4f std = %e" % (np.mean(sub_Nu_t), np.std(sub_Nu_t)))
+  print("")
+
+compute_and_print_mean_and_std("A1", 
+  "/home/rkk/A1_10/A1_10.volume_avg", 
+  "/home/rkk/A1_10/surface-heat-flux-A1_10", 
+  0.7, 1.0)
+\end{minted}
+
+\end{document}
diff --git a/benchmarks/citcoms-benchmark.py b/benchmarks/citcoms-benchmark.py
new file mode 100644
index 0000000..4e4d10a
--- /dev/null
+++ b/benchmarks/citcoms-benchmark.py
@@ -0,0 +1,70 @@
+'''
+citcoms-benchmark.py
+'''
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib.ticker as ticker
+
+d=0.45
+kappa=1.0
+time_rescale = (d*d)/kappa # t_citcom / time_rescale = t_nondim (or t_paper) 
+vel_rescale = kappa / d    # v_citcom / vel_rescale = v_nondim (i.e. v_paper)
+r_t = 1.0
+r_b = 0.55
+top_prefac = r_t*(r_t-r_b)/r_b;
+
+dataA1 = np.loadtxt("/home/rkk/A1_10/A1_10.volume_avg")
+dataA3 = np.loadtxt("/home/rkk/A3_08/A3_08.volume_avg")
+dataA8 = np.loadtxt("/home/rkk/A8_05/A8_05.volume_avg")
+surfA1 = top_prefac*np.loadtxt("/home/rkk/A1_10/surface-heat-flux-A1_10")
+surfA3 = top_prefac*np.loadtxt("/home/rkk/A3_08/surface-heat-flux-A3_08")
+surfA8 = top_prefac*np.loadtxt("/home/rkk/A8_05/surface-heat-flux-A8_05")
+timeA1 = (1.0/time_rescale)*dataA1[:,1]
+timeA3 = (1.0/time_rescale)*dataA3[:,1]
+timeA8 = (1.0/time_rescale)*dataA8[:,1]
+T_avgA1 = dataA1[:,2]
+T_avgA3 = dataA3[:,2]
+T_avgA8 = dataA8[:,2]
+Vrms_avgA1 = (1.0/vel_rescale)*dataA1[:,3]
+Vrms_avgA3 = (1.0/vel_rescale)*dataA3[:,3]
+Vrms_avgA8 = (1.0/vel_rescale)*dataA8[:,3]
+
+plt.figure(1)
+# plt.subplot(3,1,1)
+plt.plot(timeA1,T_avgA1, 'g-', label='A1')
+plt.plot(timeA3,T_avgA3, 'r-', label='A3')
+plt.plot(timeA8,T_avgA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('<T>')
+plt.title('Time dependence of Volume-averaged temperature')
+plt.ylim((0,0.8))
+plt.xlim((0,1.0))
+plt.yticks([0.0,0.2,0.4,0.6,0.8])
+
+plt.figure(2)
+# plt.subplot(3,1,2)
+plt.plot(timeA1,Vrms_avgA1, 'g-', label='A1')
+plt.plot(timeA3,Vrms_avgA3, 'r-', label='A3')
+plt.plot(timeA8,Vrms_avgA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('<V_rms>')
+plt.title('Time dependence of Root-mean squared velocity')
+plt.ylim((0,120.0))
+plt.xlim((0,1.0))
+plt.yticks(np.arange(0.0,121.0,20.0))
+
+plt.figure(3)
+# plt.subplot(3,1,3)
+plt.plot(timeA1,surfA1, 'g-', label='A1')
+plt.plot(timeA3,surfA3, 'r-', label='A3')
+plt.plot(timeA8,surfA8, 'b-', label='A8')
+plt.legend(title='Benchmark', loc='best')
+plt.xlabel('time')
+plt.ylabel('Nu_t')
+plt.title('Time dependence of surface Nusselt number')
+plt.ylim((0,6.0))
+plt.xlim((0,1.0))
+
+plt.show()
diff --git a/benchmarks/table6.py b/benchmarks/table6.py
new file mode 100644
index 0000000..ee36f69
--- /dev/null
+++ b/benchmarks/table6.py
@@ -0,0 +1,99 @@
+#-----------------------------------------------------------------------------------------------------------------------
+# table6.py : Generate the entries for Table 6 of the Zhong et. al. paper
+#-----------------------------------------------------------------------------------------------------------------------
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib.ticker as ticker
+
+d=0.45
+kappa=1.0
+time_rescale = (d*d)/kappa # t_citcom / time_rescale = t_nondim (or t_paper) 
+vel_rescale = kappa / d    # v_citcom / vel_rescale = v_nondim (i.e. v_paper)
+r_t = 1.0
+r_b = 0.55
+top_prefac = r_t*(r_t-r_b)/r_b;
+
+def get_array_slice(arr, low, high):
+  '''
+  returns a boolean array with True for arr in[low,high]
+  '''
+
+  lower_limit = np.greater_equal(arr, low)
+  upper_limit = np.less_equal(arr, high)
+  combined = np.logical_and(lower_limit, upper_limit)
+  return combined
+
+def compute_and_print_mean_and_std(benchmark, volume_avg_file, heat_flux_file, tlow, thigh):
+  data = np.loadtxt(volume_avg_file)
+  surf = top_prefac * np.loadtxt(heat_flux_file)
+  time = (1.0/time_rescale) * data[:,1]
+  T_avg = data[:,2]
+  Vrms_avg = (1.0/vel_rescale)*data[:,3]
+  combined = get_array_slice(time, tlow, thigh)
+  sub_Tavg = T_avg[combined]
+  sub_Vrms = Vrms_avg[combined]
+  sub_Nu_t = surf[combined]
+
+  print("----------------- %s Benchmark Results -----------------------" % benchmark)
+  print("<T> = %8.4f std = %e" % (np.mean(sub_Tavg), np.std(sub_Tavg)))
+  print("<V_rms> = %8.4f std = %e" % (np.mean(sub_Vrms), np.std(sub_Vrms)))
+  print("<Nu_t> = %8.4f std = %e" % (np.mean(sub_Nu_t), np.std(sub_Nu_t)))
+  print("")
+
+compute_and_print_mean_and_std("A1", 
+  "/home/rkk/A1_10/A1_10.volume_avg", 
+  "/home/rkk/A1_10/surface-heat-flux-A1_10", 
+  0.7, 1.0)
+
+'''
+compute_and_print_mean_and_std("A2", 
+  "/home/rkk/A2_02/A2_02.volume_avg", 
+  "/home/rkk/A2_02/surface-heat-flux-A2_02", 
+  1.0, 1.3)
+'''
+compute_and_print_mean_and_std("A3", 
+  "/home/rkk/A3_08/A3_08.volume_avg", 
+  "/home/rkk/A3_08/surface-heat-flux-A3_08", 
+  0.6, 0.9)
+'''
+compute_and_print_mean_and_std("A4", 
+  "/home/rkk/A4_02/A4_02.volume_avg", 
+  "/home/rkk/A4_02/surface-heat-flux-A4_02", 
+  1.5, 2.0)
+
+compute_and_print_mean_and_std("A5", 
+  "/home/rkk/A5_01/A5_01.volume_avg", 
+  "/home/rkk/A5_01/surface-heat-flux-A5_01", 
+  1.0, 1.5)
+
+compute_and_print_mean_and_std("A7", 
+  "/home/rkk/A7_01/A7_01.volume_avg", 
+  "/home/rkk/A7_01/surface-heat-flux-A7_01", 
+  1.2, 1.7)
+'''
+compute_and_print_mean_and_std("A8", 
+  "/home/rkk/A8_05/A8_05.volume_avg", 
+  "/home/rkk/A8_05/surface-heat-flux-A8_05", 
+  0.8, 1.0)
+'''
+compute_and_print_mean_and_std("B1", 
+  "/home/rkk/B1_01/B1_01.volume_avg", 
+  "/home/rkk/B1_01/surface-heat-flux-B1_01", 
+  1.2, 1.7)
+
+compute_and_print_mean_and_std("B2", 
+  "/home/rkk/B2_01/B2_01.volume_avg", 
+  "/home/rkk/B2_01/surface-heat-flux-B2_01", 
+  0.8, 1.1)
+
+compute_and_print_mean_and_std("B3", 
+  "/home/rkk/B3_01/B3_01.volume_avg", 
+  "/home/rkk/B3_01/surface-heat-flux-B3_01", 
+  0.8, 1.2)
+
+compute_and_print_mean_and_std("B4", 
+  "/home/rkk/B4_01/B4_01.volume_avg", 
+  "/home/rkk/B4_01/surface-heat-flux-B4_01", 
+  1.0, 1.3)
+'''
+



More information about the CIG-COMMITS mailing list