*Please download the file homework8.R and write your
results there. Send the your answers to my
mailbox.*

# 1. Replicate plots from last class

In the last class we simulated the *water formation* system
under different rates and initial conditions. The system is represented
by this graph:

The code to make one simulation is this:

```
<- function(N, r1_rate, r2_rate, H_ini, O_ini, W_ini) {
water_formation <- d_W <- rep(NA, N) # Water, quantity and change on each time
W <- d_H <- rep(NA, N) # Hydrogen
H <- d_O <- rep(NA, N) # Oxygen
O 1] <- W_ini
W[1] <- H_ini
H[1] <- O_ini
O[1] <- d_H[1] <- d_O[1] <- 0 # the initial change is zero
d_W[
for(i in 2:N) {
<- r1_rate*H[i-1]*H[i-1]*O[i-1] - r2_rate*W[i-1]
d_W[i] <- -r1_rate*H[i-1]*H[i-1]*O[i-1] + r2_rate*W[i-1]
d_O[i] <- -2*r1_rate*H[i-1]*H[i-1]*O[i-1] + 2*r2_rate*W[i-1]
d_H[i]
<- W[i-1] + d_W[i]
W[i] <- O[i-1] + d_O[i]
O[i] <- H[i-1] + d_H[i]
H[i]
}return(data.frame(Step=1:N, W, H, O))
}
```

Please write the code to draw the following plot. You will need to simulate several times for different rates and initial conditions

# 2. Simulate a chaotic system

In the class we also saw a *chaotic system*, that can have
high sensitivity to initial conditions. The system is very simple

After simplification, this system can be simulated by the following code:

```
<- function(N, A, x_ini) {
quad_map <- rep(NA, N)
x 1] <- x_ini
x[for(i in 2:N) {
<- A * x[i-1] * (1 - x[i-1])
x[i]
}return(x)
}
```

As we said, the behavior depends on the value of `A`

.
Please write the code that simulates this system for different values of
`A`

from 2.9 to 4 by 0.001, with `x_ini=0.4`

and
`N=1500`

. Take the final 500 values of each simulation and
draw them using `pch="."`

. You should see this:

There are at least two ways to do this plot.

- Make a very long vector with 500 values from each simulation, and plot all in one command
- Draw first an empty plot using
`type="n"`

, and then use the function`points()`

to add the 500 values of each simulation separately

# 3. Many dice

In class we simulated throwing two dice at the same time. What if we throw 3 dice? Or 10? Or 100?

Please write the code to simulate `N`

times the experiment
**“throw M dice and count the total of
points”**. Something like

```
<- function(N, M) {
many_dice # repeat N times
# throw M dice
# sum them
# put the result in a vector `ans`
# return `ans` vector
}
```

If you do it right, you should be able to draw a plot like this: