Debugging, Profiling and Optimization

# Debugging, Profiling and Optimization
## RaukR 2021 • Advanced R for Bioinformatics
### <b>Ashfaq Ali (Slides Credit Marcin Kierczak)</b>
### NBIS, SciLifeLab

---

---

## Contents

* [Overview](#overview)
* [Types of bugs](#types-of-bugs)
* [Handling Errors](#error-handling-1)
* [Debugging -- Errors and Warnings](#errors-and-warnings)
* [Debugging -- What are my Options?](#debugging-options)
* [Profiling](#profiling_proc_time)
* [Advanced profiling - **`Rprof`**](#Rprof)
* [Profiling **`profr`** package](#profr_package)
* [Profiling with **`profvis`** package](#profviz_package)
* [Optimizing your code](#optimizing_code1)
* [Copy-on-modify and memory allocation](#copy-on-modify)
* [Allocating memory](#memory-benchmark)
* [GPU](#gpu-1)
* [Parallelization using package **`parallel`**](#prallelization)

---
name: overview
## Topics of This Presentation
<br><br>
Code:
<br>

* **Debugging** -- my code does not run.

* **Profiling** -- now it does run but... out of memory!

* **Optimization** -- making things better.

---

### There are different types of bugs we can introduce:
* Syntax -- `prin(var1), mean(sum(seq((x + 2) * (y - 9 * b)))`

* Arithmetic -- `x/0` (not in R, though!) `Inf/Inf`

* Type -- `mean('a')`

* Logic -- everything works and produces seemingly valid output that is WRONG!

### To avoid bugs:
* Encapsulate your code in smaller units (functions), you can test.

* Use classes and type checking.

* Test at the boundaries, e.g. loops at min and max value.

* Feed your functions with test data that should result with a known output.

* Use *antibugging*: `stopifnot(y <= 75)`

---
name: arithmetic-bugs
## Arithmetic bugs

```r
(vec <- seq(0.1, 0.9, by=0.1))
vec == 0.7 
vec == 0.5
(0.5 + 0.1) - 0.6
(0.7 + 0.1) - 0.8 # round((0.7 + 0.1) , digits = 2) - 0.8
```

```
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
## [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [1] FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE
## [1] 0
## [1] -1.110223e-16
```
Beware of floating point arithmetics!
.small[

```r
head(unlist(.Machine))
head(unlist(.Platform))
```

```
##     double.eps double.neg.eps    double.xmin    double.xmax    double.base  double.digits 
##   2.220446e-16   1.110223e-16  2.225074e-308  1.797693e+308   2.000000e+00   5.300000e+01 
##      OS.type     file.sep   dynlib.ext          GUI       endian      pkgType 
##       "unix"          "/"        ".so"    "RStudio"     "little" "mac.binary"
```
]

---
name: error-handling-1

## Handling Errors

```r
input <- c(1, 10, -7, -2/5, 0, 'char', 100, pi)
for (val in input) {
  (paste0('Log of ', val, 'is ', log10(val)))
}
```

```
## Error in log10(val): non-numeric argument to mathematical function
```

One option is to use the `try` block:

```r
for (val in input) {
  val <- as.numeric(val)
  try(print(paste0('Log of ', val, ' is ', log10(val))))
}
```

```
## [1] "Log of 1 is 0"
## [1] "Log of 10 is 1"
## [1] "Log of -7 is NaN"
## [1] "Log of -0.4 is NaN"
## [1] "Log of 0 is -Inf"
## [1] "Log of NA is NA"
## [1] "Log of 100 is 2"
## [1] "Log of 3.14159265358979 is 0.497149872694133"
```

---
name: error-handling-2

## Handling Errors with `tryCatch`

```r
for (val in input) {
  val <- as.numeric(val)
  result <- tryCatch(log10(val), 
              warning = function(w) { print('Negative argument supplied. Negating.'); log10(-val) }, 
              error = function(e) { print('Not a number!'); NaN })
  print(paste0('Log of ', val, ' is ', result))
}
```

```
## [1] "Log of 1 is 0"
## [1] "Log of 10 is 1"
## [1] "Negative argument supplied. Negating."
## [1] "Log of -7 is 0.845098040014257"
## [1] "Negative argument supplied. Negating."
## [1] "Log of -0.4 is -0.397940008672038"
## [1] "Log of 0 is -Inf"
## [1] "Log of NA is NA"
## [1] "Log of 100 is 2"
## [1] "Log of 3.14159265358979 is 0.497149872694133"
```

---
name: errors-and-warnings

## Debugging -- Errors and Warnings
* An error in your code will result in a call to the `stop()` function that:

+ breaks the execution of the program (loop, if-statemetnt, etc.),
 
  + performs the action defined by the global parameter `error`.
  
* A warning just prints out the warning message (or reports it in another way).

* Global parameter `error` defines what R should do when an error occurs.

```r
options(error = )
```

* You can use `simpleError()` and `simpleWarning()` to generate errors and warnings in your code:

```r
f <- function(x) {
  if (x < 0) {
    x <- abs(x)
    w <- simpleWarning("Value less than 0. Taking abs(x)")
    w
  }
}
```

---
name: debugging-options

## Debugging -- What are my Options?

* Old-school debugging: a lot of `print` statements
  + print values of your variables at some checkpoints,
  + sometimes fine but often laborious,
  + need to remove/comment out manually after debugging.
  
* Dumping frames
  + on error, R state will be saved to a file,
  + file can be read into debugger,
  + values of all variables can be checked,
  + can debug on another machine, e.g. send dump to your colleague!

* Traceback
  + a list of the recent function calls with values of their params,

* Step-by-step debugging
  + execute code line by line within the debugger
  
---
name: debugging_dump_frames
## Option 1: Dumping Frames

```r
options(error = quote(dump.frames("testdump", TRUE)))

f <- function(x) {
    sin(x)
}
f('test')
```

```
## Error in sin(x): non-numeric argument to mathematical function
```

```r
options(error = NULL)
load("testdump.rda")
# debugger(testdump)
```

.smaller[<tt>Message:  Error in sin(x) : non-numeric argument to mathematical function <br>
Available environments had calls: <br>
1: f("test") <br>
 <br>
Enter an environment number, or 0 to exit   <br>
Selection: 1 <br>
Browsing in the environment with call: <br>
   f("test") <br>
Called from: debugger.look(ind) <br>
Browse[1]> x <br>
[1] "test" <br>
Browse[1]>  <br>
[1] "test" <br>
Browse[1]> 
</tt>
]
Last empty line brings you back to the environments menu.

---
name: debugging_traceback
## Option 2: Traceback

```r
f <- function(x) {
  log10(x)  
}

g <- function(x) {
  f(x)
}
g('test')
```

```
## Error in log10(x): non-numeric argument to mathematical function
```
<tt>
> traceback()<br>
2: f(x) at #2<br>
1: g("test")<br>
</tt>

`traceback()` shows what were the function calls and what parameters were passed to them when the error occured.

---
name: step_by_step-1
## Option 3: Debug step-by-step

Let us define a new function `h(x, y)`:

```r
h <- function(x, y) { 
  f(x) 
  f(y) 
  }
```
Now, we can use `debug()` to debug the function in a step-by-step manner:

```r
debug(h)
h('text', 7)
undebug(h)
```

---
name: step_by_step-2
## Option 3: Debug step-by-step cted.

`n` -- execute next line, `c` -- execute whole function, `q` -- quit debugger mode.

<tt>
> debug(h)<br>
> h('text', 7)<br>
debugging in: h("text", 7)<br>
debug at #1: {<br>
    f(x)<br>
    f(y)<br>
}<br>
</tt>
--
<tt>Browse[2]> x<br>
[1] "text"<br></tt>
--
<tt>Browse[2]> y<br>
[1] 7<br></tt>
--
<tt>Browse[2]> n<br>
debug at #2: f(x)<br></tt>
--
<tt>Browse[2]> x<br>
[1] "text"<br></tt>
--
<tt>Browse[2]> n<br>
Error in log10(x) : non-numeric argument to mathematical function<br>
</tt>

---

## Profiling -- `proc.time()`
Profiling is the process of identifying memory and time bottlenecks in your code.

```r
proc.time()
```

```
##       user     system    elapsed 
##   4826.843   2731.461 565646.913
```
* `user time` -- CPU time charged for the execution of user instructions of the calling process, 
* `system time` -- CPU time charged for execution by the system on behalf of the calling process,
* `elapsed time` -- total CPU time elapsed for the currently running R process.

```r
pt1 <- proc.time()
tmp <- runif(n =  10e5)
pt2 <- proc.time()
pt2 - pt1
```

```
##    user  system elapsed 
##   0.037   0.004   0.035
```
---
name: profiling_system_time

## Profiling -- `system.time()`

```r
system.time(runif(n = 10e6))
system.time(rnorm(n = 10e6))
```

```
##    user  system elapsed 
##   0.277   0.020   0.298 
##    user  system elapsed 
##   0.577   0.015   0.594
```

An alternative approach is to use `tic` and `toc` statements from the `tictoc` package.

```r
library(tictoc)
tic()
tmp1 <-runif(n = 10e6)
toc()
```

```
## 0.347 sec elapsed
```

---
name: profiling_system_time

## Profiling in Action

Let's see profiling in action! We will define four functions that fill a large vector in two different ways:
--

```r
fun_fill_loop1 <- function(n = 10e6, f) {
  result <- NULL
  for (i in 1:n) {
    result <- c(result, eval(call(f, 1)))
  }
  return(result)
}
```
--

```r
fun_fill_loop2 <- function(n = 10e6, f) {
  result <- vector(length = n)
  for (i in 1:n) {
    result[i] <- eval(call(f, 1))
  }
  return(result)
}
```
 
---
name: profiling_in_action

## Profiling in Action cted.

It is maybe better to use...
--
vectorization!
--

```r
fun_fill_vec1 <- function(n = 10e6, f) {
  result <- NULL
  result <- eval(call(f, n))
  return(result)
}
```
--

```r
fun_fill_vec2 <- function(n = 10e6, f) {
  result <- vector(length = n)
  result <- eval(call(f, n))
  return(result)
}
```

---
name: compare_loop_vec_sys_time

## Profiling our functions

```r
system.time(fun_fill_loop1(n = 10e4, "runif")) # Loop 1
system.time(fun_fill_loop2(n = 10e4, "runif")) # Loop 2
system.time(fun_fill_vec1(n = 10e4, "runif"))  # Vectorized 1
system.time(fun_fill_vec2(n = 10e4, "runif"))  # Vectorized 2
```

```
##    user  system elapsed 
##  17.745  11.612  29.469 
##    user  system elapsed 
##   0.368   0.053   0.424 
##    user  system elapsed 
##   0.004   0.000   0.004 
##    user  system elapsed 
##   0.003   0.000   0.003
```

The `system.time()` function is not the most accurate though. During the lab, we will experiment with package `microbenchmark`.

---
name: Rprof

## More advanced profiling
We can also do a bit more advanced profiling, including the memory profiling, using, e.g. `Rprof()` function.

And let us summarise:
.small[

```r
summary <- summaryRprof('profiler_test.out', memory='both')
datatable(summary$by.self, options=list(pageLength = 10, searching = F, info = F))
#knitr::kable(summary$by.self)
```

<div id="htmlwidget-a5ac6d90465682c41b01" style="width:100%;height:auto;" class="datatables html-widget"></div>
<script type="application/json" data-for="htmlwidget-a5ac6d90465682c41b01">{"x":{"filter":"none","data":[["\"print.default\"","\"runif\"","\"eval\"","\"fun_fill_loop2\"","\"paste0\"","\"parent.frame\"","\"is.list\""],[102.82,0.89,0.54,0.13,0.08,0.05,0.02],[98.36,0.85,0.52,0.12,0.08,0.05,0.02],[102.82,0.89,104.45,1.63,0.08,0.05,0.02],[98.36,0.85,99.92,1.56,0.08,0.05,0.02],[[36382.2],[981.9],[38299.1],[1916.9],[9.4],[42.9],[26.7]]],"container":"<table class=\"display\">\n  <thead>\n    <tr>\n      <th> <\/th>\n      <th>self.time<\/th>\n      <th>self.pct<\/th>\n      <th>total.time<\/th>\n      <th>total.pct<\/th>\n      <th>mem.total<\/th>\n    <\/tr>\n  <\/thead>\n<\/table>","options":{"pageLength":10,"searching":false,"info":false,"columnDefs":[{"className":"dt-right","targets":[1,2,3,4,5]},{"orderable":false,"targets":0}],"order":[],"autoWidth":false,"orderClasses":false}},"evals":[],"jsHooks":[]}</script>
]

---
name: profr_package

## Profiling -- `profr` package
There are also packages available that enable even more advanced profiling:

```r
library(profr)
Rprof("profiler_test2.out", interval = 0.01)
tmp <- table(sort(rnorm(1e5)))
Rprof(NULL)
profile_df <- parse_rprof('profiler_test2.out')
```
This returns a table that can be visualised:

<div id="htmlwidget-7f27c5b3112aef697b8f" style="width:100%;height:auto;" class="datatables html-widget"></div>
<script type="application/json" data-for="htmlwidget-7f27c5b3112aef697b8f">{"x":{"filter":"none","data":[["1","2","3","4","5","6","7","8","9","10","11","12"],[1,2,2,2,3,3,3,3,4,4,4,4],[1,1,2,3,1,2,3,4,1,2,3,4],[1,1,1,1,1,1,1,1,1,1,1,1],["table","sort","factor","as.integer","rnorm","sort.default","unique","as.character","sort.int","unique.default","as.character","unique.default"],[0,0,0.04,0.6,0,0.02,0.04,0.4,0.02,0.04,0.06,0.38],[0.62,0.04,0.6,0.62,0.02,0.04,0.4,0.58,0.04,0.06,0.38,0.4],[1,1,1,1,1,1,1,1,1,1,1,1],[false,false,false,true,true,false,false,true,true,true,true,true],[0.62,0.04,0.56,0.02,0.02,0.02,0.36,0.18,0.02,0.02,0.32,0.02],["base","base","base","base","stats","base","base","base","base","base","base","base"]],"container":"<table class=\"display\">\n  <thead>\n    <tr>\n      <th> <\/th>\n      <th>level<\/th>\n      <th>g_id<\/th>\n      <th>t_id<\/th>\n      <th>f<\/th>\n      <th>start<\/th>\n      <th>end<\/th>\n      <th>n<\/th>\n      <th>leaf<\/th>\n      <th>time<\/th>\n      <th>source<\/th>\n    <\/tr>\n  <\/thead>\n<\/table>","options":{"pageLength":5,"searching":false,"info":false,"columnDefs":[{"className":"dt-right","targets":[1,2,3,5,6,7,9]},{"orderable":false,"targets":0}],"order":[],"autoWidth":false,"orderClasses":false,"lengthMenu":[5,10,25,50,100]}},"evals":[],"jsHooks":[]}</script>

---
name: profr_package_cted

## Profiling -- `profr` package cted.

We can also plot the results using -- `proftools` package-

```r
library("proftools") # depends on "graph" and "Rgraphviz" packages
profile_df2 <- readProfileData("profiler_test2.out")
plotProfileCallGraph(profile_df2, style = google.style, score = "total")
```

---
name: profviz_package

## Profiling with `profvis`

Yet another nice way to profile your code is by using Hadley Wickham's `profvis` package:

```r
library(profvis)
profvis({fun_fill_loop2(1e4, 'runif')
  fun_fill_vec2(1e4, 'runif')
  })
```

---
name: profvis_run

## Profiling with `profvis` cted.
<div id="htmlwidget-b37886f13bf7ed8933a8" style="width:100%;height:600px;" class="profvis html-widget"></div>
<script type="application/json" data-for="htmlwidget-b37886f13bf7ed8933a8">{"x":{"message":{"prof":{"time":[1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,5,5,5,5,5,6,6,7,7,7,7,7,8,8,8,8,8,9,10,10,10,10,10,11,11,11,11,11],"depth":[10,9,8,7,6,5,4,3,2,1,14,13,12,11,10,9,8,7,6,5,4,3,2,1,14,13,12,11,10,9,8,7,6,5,4,3,2,1,20,19,18,17,16,15,14,13,12,11,10,9,8,7,6,5,4,3,2,1,5,4,3,2,1,2,1,5,4,3,2,1,5,4,3,2,1,1,5,4,3,2,1,5,4,3,2,1],"label":["unlist","Filter","findLocalsList","funEnv","make.functionContext","cmpfun","compiler:::tryCmpfun","evalq","evalq","tryCatch","findCenvVar","getInlineInfo","isBaseVar","FUN","lapply","Filter","findLocalsList","funEnv","make.functionContext","cmpfun","compiler:::tryCmpfun","evalq","evalq","tryCatch","findCenvVar","getInlineInfo","isBaseVar","FUN","lapply","Filter","findLocalsList","funEnv","make.functionContext","cmpfun","compiler:::tryCmpfun","evalq","evalq","tryCatch","paste0","cb$makelabel","tryInline","cmpCall","cmp","cmpSymbolAssign","h","tryInline","cmpCall","cmp","h","tryInline","cmpCall","cmp","genCode","cmpfun","compiler:::tryCmpfun","evalq","evalq","tryCatch","runif","fun_fill_loop2","evalq","evalq","tryCatch","is.list","tryCatch","runif","fun_fill_loop2","evalq","evalq","tryCatch","runif","fun_fill_loop2","evalq","evalq","tryCatch","tryCatch","runif","fun_fill_loop2","evalq","evalq","tryCatch","eval","fun_fill_loop2","evalq","evalq","tryCatch"],"filenum":[null,null,null,null,null,null,1,null,null,null,null,null,null,null,null,null,null,null,null,null,1,null,null,null,null,null,null,null,null,null,null,null,null,null,1,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,1,null,null,null,null,1,null,null,null,null,null,null,1,null,null,null,null,1,null,null,null,null,null,1,null,null,null,null,1,null,null,null],"linenum":[null,null,null,null,null,null,2,null,null,null,null,null,null,null,null,null,null,null,null,null,2,null,null,null,null,null,null,null,null,null,null,null,null,null,2,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,2,null,null,null,null,2,null,null,null,null,null,null,2,null,null,null,null,2,null,null,null,null,null,2,null,null,null,null,2,null,null,null],"memalloc":[107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621612548828,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621726989746,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.621940612793,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.648765563965,107.828941345215,107.828941345215,107.828941345215,107.828941345215,107.828941345215,110.63037109375,110.63037109375,113.841697692871,113.841697692871,113.841697692871,113.841697692871,113.841697692871,117.655723571777,117.655723571777,117.655723571777,117.655723571777,117.655723571777,121.995300292969,126.272239685059,126.272239685059,126.272239685059,126.272239685059,126.272239685059,131.096405029297,131.096405029297,131.096405029297,131.096405029297,131.096405029297],"meminc":[0,0,0,0,0,0,0,0,0,0,0.00011444091796875,0,0,0,0,0,0,0,0,0,0,0,0,0,0.000213623046875,0,0,0,0,0,0,0,0,0,0,0,0,0,0.026824951171875,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0.18017578125,0,0,0,0,2.80142974853516,0,3.21132659912109,0,0,0,0,3.81402587890625,0,0,0,0,0,4.27693939208984,0,0,0,0,0,0,0,0,0],"filename":[null,null,null,null,null,null,"<expr>",null,null,null,null,null,null,null,null,null,null,null,null,null,"<expr>",null,null,null,null,null,null,null,null,null,null,null,null,null,"<expr>",null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,null,"<expr>",null,null,null,null,"<expr>",null,null,null,null,null,null,"<expr>",null,null,null,null,"<expr>",null,null,null,null,null,"<expr>",null,null,null,null,"<expr>",null,null,null]},"interval":10,"files":[{"filename":"<expr>","content":"library(profvis)\nprofvis({fun_fill_loop2(1e4, 'runif')\n  c()\n  })","normpath":"<expr>"}],"prof_output":"/Users/cob-aaf/Documents/GitHub/RaukR-2021/docs/debugging_profiling_Ashfaq/.Rproj.user/77B5F86E/profiles-cache/file34b854aad517.Rprof","highlight":{"output":["^output\\$"],"gc":["^<GC>$"],"stacktrace":["^\\.\\.stacktraceo(n|ff)\\.\\.$"]},"split":"h"}},"evals":[],"jsHooks":[]}</script>

---
name: optimizing_code1

## Optimizing your code

*We should forget about small efficiencies, say about 97% of the   time: premature optimization is the root of all evil. Yet we should not pass up our opportunities in that critical 3%. A good programmer will not be delulled into complacency by such reasoning, he will be wise to look carefully at the critical code; but only after that code has been identified.*
<br><br>
*-- Donald Knuth*

<div class="pull-left">
<img src="./assets/xkcd_automation.png" style="height:300px;">
<br>
.vsmall[source: http://www.xkcd/com/1319]
</div>

<div class="pull-right">
<img src="./assets/xkcd_is_it_worth_the_time_2x.png" style="height:300px;">
<br>
.vsmall[source: http://www.xkcd/com/1205]
</div>

---
name: optimization_types

## Ways to optimize the code

* write it in a more efficient way, e.g. use vectorization or `*apply` family instead of loops etc.,
* allocating memory to avoid copy-on-modify,
* use package `BLAS` for linear algebra,
* use `bigmemory` package,
* GPU computations,
* multicore support, e.g. `multicore`, `snow`
* use `data.table` or `tibble` instead of `data.frame`

---
name: copy-on-modify

## Copy-on-modify

```r
library(pryr)
order <- 1024
matrix_A <- matrix(rnorm(order^2), nrow = order)
matrix_B <- matrix_A
```
--
Check where the objects are in the memory:
--

```r
address(matrix_A)
address(matrix_B)
```

```
## [1] "0x7ffb52b9b000"
## [1] "0x7ffb52b9b000"
```
--
What happens if we modify a value in one of the matrices?
--

```r
matrix_B[1,1] <- 1
address(matrix_A)
address(matrix_B)
```

```
## [1] "0x7ffb52b9b000"
## [1] "0x7ffb5339c000"
```

---
name: avoid-copy-on-modify-1

## Avoid copying by allocating memory

### No memory allocation

```r
f1 <- function(to = 3, silent=F) {
  tmp <- c()
  for (i in 1:to) {
    a1 <- address(tmp)
    tmp <- c(tmp, i)
    a2 <- address(tmp)
    if(!silent) { print(paste0(a1, " --> ", a2)) } 
  }
}
f1()
```

```
## [1] "0x7ffbb100aee0 --> 0x7ffb9e10f158"
## [1] "0x7ffb9e10f158 --> 0x7ffb9e10f4d8"
## [1] "0x7ffb9e10f4d8 --> 0x7ffb9f17f208"
```

---
name: avoid-copy-on-modify-1
## Avoid copying by allocating memory cted.

### With allocation

```r
f2 <- function(to = 3, silent = FALSE) {
  tmp <- vector(length = to, mode='numeric')
  for (i in 1:to) {
    a1 <- address(tmp)
    tmp[i] <- i
    a2 <- address(tmp)
    if(!silent) { print(paste0(a1, " --> ", a2)) }
  }
}
f2()
```

```
## [1] "0x7ffb9b21ba18 --> 0x7ffb9b21ba18"
## [1] "0x7ffb9b21ba18 --> 0x7ffb9b21ba18"
## [1] "0x7ffb9b21ba18 --> 0x7ffb9b21ba18"
```

---
name: memory-benchmark

## Allocating memory -- benchmark.

```r
library(microbenchmark)
benchmrk <- microbenchmark(f1(to = 1e3, silent = T), 
                           f2(to = 1e3, silent = T), 
                           times = 100L)
autoplot(benchmrk)
```

---
name: gpu-1
## GPU

```r
library(gpuR)
library(microbenchmark)
A = matrix(rnorm(1000^2), nrow=1000) # stored: RAM, computed: CPU
B = matrix(rnorm(1000^2), nrow=1000) 
gpuA = gpuMatrix(A, type = "double") # stored: RAM, computed: GPU
gpuB = gpuMatrix(B, type = "double")
vclA = vclMatrix(A, type = "double") # stored: GPU, computed: GPU
vclB = vclMatrix(B, type = "double")
bch <- microbenchmark(
  cpuC = A %*% B,
  gpuC = gpuA %*% gpuB,
  vclC = vclA %*% vclB, 
  times = 10L) 
```

```
## Number of platforms: 1
## - platform: Apple: OpenCL 1.2 (May  8 2021 03:14:28)
##   - context device index: 0
##     - Intel(R) Core(TM) i7-8750H CPU @ 2.20GHz
##   - context device index: 1
##     - Intel(R) UHD Graphics 630
##   - context device index: 2
##     - AMD Radeon Pro 555X Compute Engine
## checked all devices
## completed initialization
```
.small[
More on [Charles Determan's Blog](https://www.r-bloggers.com/r-gpu-programming-for-all-with-gpur/).
]

---
name:gpu-2

## GPU cted.

```r
autoplot(bch)
```

---
name: prallelization

## Parallelization using package `parallel`
Easiest to paralllelize is `lapply`:

```r
result <- lapply(1:2, function(x) { c(x, x^2, x^3) } )
result
```

```
## [[1]]
## [1] 1 1 1
## 
## [[2]]
## [1] 2 4 8
```

```r
library(parallel)
num_cores <- detectCores() - 1
cl <- makeCluster(num_cores) # Init cluster
parLapply(cl, 1:2, function(x) { c(x, x^2, x^3)} )
stopCluster(cl)
```

```
## [[1]]
## [1] 1 1 1
## 
## [[2]]
## [1] 2 4 8
```

---
name: end-slide
class: end-slide, middle
count: false

# Thank you. Questions?

<p>R version 4.1.0 (2021-05-18)<br><p>Platform: x86_64-apple-darwin17.0 (64-bit)</p><p>OS: macOS Big Sur 11.4</p><br>

Built on : <i class='fa fa-calendar' aria-hidden='true'></i> 15-Jun-2021 at <i class='fa fa-clock-o' aria-hidden='true'></i> 08:11:49

<b>2021</b> • [SciLifeLab](https://www.scilifelab.se/) • [NBIS](https://nbis.se/) • [RaukR](https://nbisweden.github.io/workshop-RaukR-2106/)