Add examples of parameter optimization

Author: gjbex

source-code/REEADME.md

@@ -0,0 +1,27 @@
+#| Source code
+
+This is source code that is either used in the presentation, or was developed
+to create it.  There is some material not covered in the presentation as well.
+
+## Requirements
+
+* Python version: at least 3.7
+* Packages (names listed that can be used with `pip` or `conda` to install):
+  * pandas
+  * xlrd
+  * seaborn
+  * matplotlib
+  * scipy
+  * jupyter
+  * scikit-learn
+  * keras
+
+## What is it?
+* [`keras`](keras): illustration of using Keras for machine learning.
+* [`parameter-optimization`](parameter-optimization): example of parameter
+  optimization kusing hyperopt, although the examples do not optimize
+  hyperparameters in machine learning, that would be very similar.
+* [`scikit-learn`](scikit-learn): examples of scikit-learn for machine learning,
+  examples are provided forsupervised (regression and classification) and
+  unsupervised (clustering) learnign, as well as dimensionality reduction for
+  visualization of high-demensional data.

source-code/parameter-optimization/.gitignore

@@ -0,0 +1,2 @@
+*.exe
+*.pbs.*

source-code/parameter-optimization/JuliaJobs/Makefile

@@ -0,0 +1,12 @@
+CXX = g++
+CXXFLAGS = -fopenmp -std=c++14 -O2 -g -Wall -Wextra
+LDLIBS = -lm
+
+all: julia_omp.exe
+
+%.exe: %.o
+	$(CXX) $(CXXFLAGS) -o $@ $^ $(LDLIBS)
+
+clean:
+	$(RM) $(wildcard *.exe) $(wildcard *.o)
+	$(RM) core $(wildcard core.*)

source-code/parameter-optimization/JuliaJobs/README.md

@@ -0,0 +1,17 @@
+# JuliaJobs
+
+Illustration of optimizing application parameters by running jobs.
+
+## What is it?
+
+  1. `julia_omp.cpp`: C++ application that compute the Julia set
+    fractal, using OpenMP for parallelization.  The parallel
+    for-loop schedule is `runtime`, and hence can be controled
+    using the `OMP_SCHEDULE` environment variable. Timing
+    information is writting to standard error.
+  1. `Makefile`: make file to build the Julia set application.
+  1. `julia.pbs`: PBS script to run the Julia set applciation
+    as a job.
+  1. `julia_optimize.py`: Python script to find optimal OpenMP
+    schedule, chunk size and number of threads for the Julia set
+    application.

source-code/parameter-optimization/JuliaJobs/julia.pbs

@@ -0,0 +1,52 @@
+#!/usr/bin/env bash
+#PBS -A lp_sys
+#PBS -l pmem=2gb
+#PBS -l walltime=00:05:00
+#PBS -j oe
+
+cd $PBS_O_WORKDIR
+
+source switch_to_2018a
+module load foss
+
+# verify that OpenMP environment variables have been passed
+# using -v options
+if [[ -z "$schedule" ]]
+then
+    (>&2 echo "### error: OpenMP schedule not set")
+    exit 1
+fi
+if [[ -z "$chunk" ]]
+then
+    (>&2 echo "### error: OpenMP chunk not set")
+    exit 1
+fi
+
+OMP_SCHEDULE="${schedule},${chunk}"
+(>&2 echo "### info: OMP_SCHEDULE=\'${OMP_SCHEDULE}\'")
+
+if [[ -z "$OMP_NUM_THREADS" ]]
+then
+    (>&2 echo "### error: OpenMP number of threads not set")
+    exit 2
+else
+    (>&2 echo "### info: OMP_NUM_THREADS=\'${OMP_NUM_THREADS}\'")
+fi
+
+# ensure the executable is op to date
+make all
+if [[ $? -ne 0 ]]
+then
+    (>&2 echo "### error: make failed")
+    exit 3
+fi
+
+# determine time file name
+time_file="julia.pbs.time${PBS_JOBID/\.*/}"
+
+# run application on a sufficiently large matrix, disregarding
+# output, timing information is writting to a file
+for i in $(seq 3)
+do
+    ./julia_omp.exe 4096 > /dev/null 2>> "${time_file}"
+done

source-code/parameter-optimization/JuliaJobs/julia_omp.cpp

@@ -0,0 +1,94 @@
+#include <algorithm>
+#include <chrono>
+#include <cmath>
+#include <complex>
+#include <iostream>
+#include <valarray>
+
+using namespace std;
+using my_time_t = chrono::nanoseconds;
+
+using cmplx = complex<float>;
+
+valarray<float> coordinates(float min_coord, float max_coord,
+                             size_t steps);
+valarray<cmplx> z_values(const valarray<float>& x_coords,
+                         const valarray<float>& y_coords);
+valarray<int> iterate_zs(valarray<cmplx>& zs, const complex<float>& c,
+                         size_t max_iters);
+void print_results(const valarray<int>& ns);
+
+int main(int argc, char *argv[]) {
+    const cmplx c(-0.62772, - 0.42193);
+    const float x1 {-1.8};
+    const float x2 {1.8};
+    const float y1 {-1.8};
+    const float y2 {1.8};
+    const size_t max_iters {255};
+    size_t steps {100};
+    if (argc > 1)
+        steps = stoi(argv[1]);
+    valarray<float> x_coords = coordinates(x1, x2, steps);
+    valarray<float> y_coords = coordinates(y1, y2, steps);
+    valarray<cmplx> zs = z_values(x_coords, y_coords);
+    auto start_time = chrono::steady_clock::now();
+    valarray<int> ns = iterate_zs(zs, c, max_iters);
+    auto end_time = chrono::steady_clock::now();
+    auto duration = chrono::duration_cast<my_time_t>(end_time - start_time);
+    cerr << duration.count()*1.0e-9 << endl;
+    print_results(ns);
+    return 0;
+}
+
+valarray<float> coordinates(float min_coord, float max_coord,
+                             size_t steps) {
+    valarray<float> coords(steps);
+    const float step {(max_coord - min_coord)/steps};
+    float value {min_coord};
+    for (size_t i = 0; i < steps; i++) {
+        coords[i] = value;
+        value += step;
+    }
+    return coords;
+}
+
+valarray<cmplx> z_values(const valarray<float>& x_coords,
+                         const valarray<float>& y_coords) {
+    valarray<cmplx> zs(x_coords.size()*y_coords.size());
+    size_t i {0};
+    for (auto y: y_coords)
+        for (auto x: x_coords) {
+            complex<float> z(x, y);
+            zs[i++] = z;
+        }
+    return zs;
+}
+
+int iterate_z(cmplx z, const cmplx& c, size_t max_iters) {
+    size_t n {0};
+    const complex<float> z_in {z};
+    while (abs(z) < 2.0 && n++ < max_iters)
+        z = z*z + c;
+    return n;
+}
+
+valarray<int> iterate_zs(valarray<cmplx>& zs, const complex<float>& c,
+                         size_t max_iters) {
+    valarray<int> ns(zs.size());
+#pragma omp parallel for schedule(runtime)
+    for (size_t i = 0; i < zs.size(); i++)
+        ns[i] = iterate_z(zs[i], c, max_iters);
+    return ns;
+}
+
+void print_results(const valarray<int>& ns) {
+    size_t steps = ((size_t) sqrt(ns.size()) + 0.1);
+    size_t count {0};
+    for (auto n: ns) {
+        cout << n;
+        if (++count % steps == 0)
+            cout << endl;
+        else
+            cout << " ";
+    }
+}

source-code/parameter-optimization/JuliaJobs/julia_optimize.py

@@ -0,0 +1,81 @@
+#!/usr/bin/env python
+
+from argparse import ArgumentParser
+from pathlib import Path
+import subprocess
+import sys
+import time
+from hyperopt import fmin, hp, STATUS_OK, tpe, Trials
+
+
+def function(params):
+    schedule, chunk, ppn = params
+    # chunk is given as a floating point number
+    chunk = int(chunk)
+    # ppn ranges from 0 to 35 (inclusive)
+    ppn = 1 + int(ppn)
+    omp_env = (f'schedule={schedule},chunk={int(chunk)},' +
+               f'OMP_NUM_THREADS={ppn}')
+    cmd = ['qsub', '-l', f'nodes=1:ppn={ppn}:haswell',
+           '-v', omp_env, 'julia.pbs']
+    process = subprocess.run(cmd, stdout=subprocess.PIPE,
+                             encoding='utf8')
+    job_id, *_ = process.stdout.split('.')
+    print(f'### info: submitted job {job_id}', file=sys.stderr)
+    output_file = Path(f'julia.pbs.o{job_id}')
+    while not output_file.exists():
+        time.sleep(3)
+    print(f'### info: job {job_id} finished', file=sys.stderr)
+    runtimes = list()
+    with open(f'julia.pbs.time{job_id}', 'r') as time_file:
+        for line in time_file:
+            runtimes.append(float(time_file.readline()))
+    runtime = sum(runtimes)/len(runtimes)
+    return {
+        'loss': runtime, 'schedule': schedule, 'chunk': chunk,
+        'ppn': ppn, 'job_id': job_id, 'status': STATUS_OK,
+        'time': time.strftime('%Y-%m-%d %H:%M:%S'),
+    }
+
+
+def optimize(max_evals, max_ppn):
+    space = hp.choice('schedule', [
+        ('static', hp.qloguniform('chunk_s', 2, 11, 10),
+         hp.randint('ppn_s', max_ppn)),
+        ('dynamic', hp.qloguniform('chunk_d', 2, 11, 10),
+         hp.randint('ppn_d', max_ppn)),
+        ('guided', hp.qloguniform('chunk_g', 2, 11, 10),
+         hp.randint('ppn_g', max_ppn)),
+    ])
+    trials = Trials()
+    best = fmin(function, space=space, algo=tpe.suggest,
+                max_evals=max_evals, trials=trials)
+    return best, trials
+
+
+def main():
+    arg_parser = ArgumentParser(description='optimize external '
+                                            'process')
+    arg_parser.add_argument('--max-ppn', type=int, default=20,
+                            help='maximum number of cores to use')
+    arg_parser.add_argument('--max-evals', type=int,
+                            default=100, help='maximum evals')
+    arg_parser.add_argument('--trials', required=True,
+                            help='file to save trials')
+    options = arg_parser.parse_args()
+    _, trials = optimize(options.max_evals, options.max_ppn)
+    with open(options.trials, 'w') as trials_file:
+        print('schedule,chunk,ppn,job_id,runtime',
+              file=trials_file)
+        for trial in trials.results:
+            schedule = trial['schedule']
+            chunk = trial['chunk']
+            ppn = trial['ppn']
+            job_id = trial['job_id']
+            runtime = trial['loss']
+            print(f'{schedule},{chunk:d},{ppn:d},'
+                  f'({job_id}),{runtime}', file=trials_file)
+
+
+if __name__ == '__main__':
+    sys.exit(main())

source-code/parameter-optimization/JuliaProcesses/Makefile

@@ -0,0 +1,12 @@
+CXX = g++
+CXXFLAGS = -fopenmp -std=c++14 -O2 -g -Wall -Wextra
+LDLIBS = -lm
+
+all: julia_omp.exe
+
+%.exe: %.o
+	$(CXX) $(CXXFLAGS) -o $@ $^ $(LDLIBS)
+
+clean:
+	$(RM) $(wildcard *.exe) $(wildcard *.o)
+	$(RM) core $(wildcard core.*)

source-code/parameter-optimization/JuliaProcesses/README.md

@@ -0,0 +1,20 @@
+# JuliaProcesses
+
+Illustration of optimizing application parameters by running
+processes.
+
+## What is it?
+
+  1. `julia_omp.cpp`: C++ application that compute the Julia set
+    fractal, using OpenMP for parallelization.  The parallel
+    for-loop schedule is `runtime`, and hence can be controled
+    using the `OMP_SCHEDULE` environment variable. Timing
+    information is writting to standard error.
+  1. `Makefile`: make file to build the Julia set application.
+  1. `julia_optimize.py`: Python script to find optimal OpenMP
+    schedule, chunk size and number of threads for the Julia set
+    application.
+  1. `julia_gcc.pbs`: PBS script to optimize the Julia set
+    applciation parameters as a job (GCC compiler).
+  1. `julia_intel.pbs`: PBS script to optimize the Julia set
+    applciation parameters as a job (Intel compiler).

source-code/parameter-optimization/JuliaProcesses/julia_gcc.pbs

@@ -0,0 +1,27 @@
+#!/usr/bin/env bash
+#PBS -A lp_sys
+#PBS -l nodes=1:ppn=24:haswell
+#PBS -l pmem=2gb
+#PBS -l walltime=03:00:00
+#PBS -j oe
+#PBS -m ae
+#PBS -M geertjan.bex@uhasselt.be
+
+cd $PBS_O_WORKDIR
+
+source switch_to_2018a
+module load foss
+
+# ensure the executable is op to date
+make CXX=g++ CXXFLAGS='-fopenmp -std=c++14 -O2' clean all
+if [[ $? -ne 0 ]]
+then
+    (>&2 echo "### error: make failed")
+    exit 1
+fi
+
+module load conda
+source activate hyperopt
+
+./julia_optimize.py --max-ppn 20 --max-evals 100 \
+                    --trials trials_gcc.txt