Merge pull request #28 from UofUEpiBio/cpp-lecture

Adding the C++ lecture

05-cpp/hello-world.cpp

@@ -0,0 +1,6 @@
+#include<iostream>
+
+int main() {
+  std::cout << "Hello world" << std::endl;
+  return 0;
+}

05-cpp/lab.qmd

@@ -0,0 +1,115 @@
+---
+format: html
+---
+
+
+<div id="quarto-content"
+class="page-columns page-rows-contents page-layout-article">
+
+<div id="quarto-document-content" class="content" role="main">
+
+<div id="title-block-header" class="quarto-title-block default">
+
+<div class="quarto-title">
+
+# Lab 05 - Rcpp
+
+</div>
+
+<div class="quarto-title-meta">
+
+</div>
+
+</div>
+
+<div id="learning-goals" class="section level1">
+
+# Learning goals
+
+- Use the different data types in Rcpp.
+- Learn some fundamentals about C++ optimization.
+- Practice your GitHub skills.
+
+</div>
+
+<div id="lab-description" class="section level1">
+
+# Lab description
+
+For this lab, we will create a function for propensity score matching.
+The goal is simple: write out a C++ function with Rcpp and measure how
+faster it is compared to the following R implementation:
+
+<div id="cb1" class="sourceCode">
+
+``` sourceCode
+ps_matchR <- function(x) {
+  
+  match_expected <- as.matrix(dist(x))
+  diag(match_expected) <- .Machine$integer.max
+  indices <- apply(match_expected, 1, which.min)
+  
+  list(
+    match_id = as.integer(unname(indices)),
+    match_x  = x[indices]
+  )
+  
+}
+```
+
+</div>
+
+<div id="question-1-create-a-simple-function" class="section level2">
+
+## Question 1: Create a simple function
+
+Use the following pseudo-code template to get started:
+
+<div id="cb2" class="sourceCode">
+
+``` sourceCode
+#include <Rcpp.h>
+
+using namespace Rcpp;
+
+// [[Rcpp::export]]
+[output must be list] ps_match1(const NumericVector & x) {
+
+    ...prepare the output (save space)...
+    ...it should be an integer vector indicating the id of the match...
+    ...and a numeric vector with the value of `x` for the match...
+
+    for (...loop over i...) {
+
+        for (...loop over j and check if it is the optimum...) {
+            if (...the closests so far...) {
+                ...update the optimum...
+            }
+        }
+        
+    }
+
+    return [a list like the R function]
+
+}
+```
+
+</div>
+
+</div>
+
+<div id="question-2-things-can-be-done-faster" class="section level2">
+
+## Question 2: Things can be done faster
+
+In the previous question, we have a double loop running twice over the
+full set of observations. We need you to write the C++ so that the
+computational complexity goes below `n^2`. (hint: Distance is symmetric)
+
+</div>
+
+</div>
+
+</div>
+
+</div>

05-cpp/means.cpp

@@ -0,0 +1,27 @@
+#include<iostream> // To print
+#include<vector>   // To use vectors
+
+int main() {
+
+  // Defining the data
+  std::vector< double > dat = {1.0, 2.5, 4.4};
+
+  // Making room for the output
+  double ans = 0.0;
+
+  // For-loops have three components:
+  // - Starts in i = 0
+  // - Until i reaches dat.size() (stops)
+  // - Increments i + 1
+  for (int i = 0; i < dat.size(); ++i)
+    ans = ans + dat[i];
+
+  ans = ans/dat.size();
+
+  // Print out the value to the screen
+  std::cout << "The mean of dat is " << ans << std::endl;
+
+  // Returning
+  return 0;
+
+}

05-cpp/means2.cpp

@@ -0,0 +1,20 @@
+#include<iostream> // To print
+#include<vector>   // To use vectors
+#include<numeric>  // To use the accumulate function
+
+int main() {
+
+  // Defining the data
+  std::vector< double > dat = {1.0, 2.5, 4.4};
+
+  // Making room for the output
+  double ans = std::accumulate(dat.begin(), dat.end(), 0.0);
+  ans /= dat.size();
+
+  // Print out the value to the screen
+  printf("The mean of dat is %.2f\n", ans);
+
+  // Returning
+  return 0;
+
+}

05-cpp/person.cpp

@@ -0,0 +1,12 @@
+#include<iostream>
+#include "person.hpp"
+
+int main() {
+  Person p1; // Default constructor
+  Person p2("John", 30, 1.80); // Other constructor
+
+  std::cout << p1.get_name() << std::endl;
+  std::cout << p2.get_name() << std::endl;
+
+  return 0;
+}

05-cpp/person.hpp

@@ -0,0 +1,36 @@
+#ifndef PERSON_HPP
+#define PERSON_HPP
+
+#include<string>
+#include<iostream>
+
+class Person {
+private: // <1>
+    std::string name;
+    int age;
+    double height;
+
+public: // <2>
+  // Constructor
+  Person(std::string n, int a, double h) { // <3>
+    name = n;
+    age = a;
+    height = h;
+  };
+
+  // Default constructor
+  Person() : name("Unknown"), age(0), height(0.0) {}; // <3>
+
+  // Destructor
+  ~Person() { // <4>
+    std::cout <<
+      this->name + " destroyed" <<
+      std::endl;
+  };
+
+  // Getters and setters
+  std::string get_name() { return name; };    // <5>
+  void set_name(std::string n) { name = n; }; // <5>
+};
+
+#endif

05-cpp/pi.cpp

@@ -0,0 +1,37 @@
+#include <vector>
+#include <random> // <1>
+int main() {
+
+  // Setting the seed
+  std::mt19937 rng_engine; // <2>
+  rng_engine.seed(123); // <2>
+
+  std::uniform_real_distribution<double> dist(-1.0, 1.0); // <3>
+
+  // Number of simulations
+  size_t n_sims = 5e6;
+
+  // Defining the data
+  double pi_approx = 0.0;
+  for (size_t i = 0u; i < n_sims; ++i)
+  {
+
+    // Generating a point in the unit square
+    double x = dist(rng_engine);
+    double y = dist(rng_engine);
+
+    double dist = std::sqrt(
+        std::pow(x, 2.0) + std::pow(y, 2.0) // <2>
+        );
+
+    // Checking if the point is inside the unit circle 
+    if (dist <= 1.0)
+      pi_approx += 1.0;
+
+  }
+
+  printf("pi approx to %.4f\n", 4.0*pi_approx/n_sims);
+
+  return 0;
+
+}

05-cpp/pointers.cpp

@@ -0,0 +1,24 @@
+#include <cstdio>
+
+void set_x_copy(int x, int y) {x = y;};
+void set_x(int * x, int y) {*x = y;};
+void set_x_ref(int & x, int y) {x = y;};
+// This would generate an error
+// void set_x_ref(const int & x, int y) {x = y;};
+
+int main() {
+
+    int x = 0;
+
+    set_x_copy(x, 3);
+    std::printf("x = %d\n", x);
+
+    set_x(&x, 2);
+    std::printf("x = %d\n", x);
+
+    set_x_ref(x, 1);
+    std::printf("x = %d\n", x);
+
+    return 0;
+
+}