slumping model

[This article was first published on Dan Kelley Blog/R, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Introduction

I got interested in layered sedimentation from viewing a video and decided it would be interesting to code this into R. More on that in due course, but my first step was to code a syatem with one sediment “type”.

Procedure

The following code drops sediment particles at x=1, and lets them roll downhill until they reach the bottom or a ledge. It draws numbers at the sedimented particles’ final positions. Since the numbers start at 1, the values are like inverse ages.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
m <- 51  # number of particles
n <- 10  # grid width
debug <- FALSE  # put TRUE for debugging
info <- function(...) if (debug) cat(...)
pch <- 20
cex <- 4/log2(n)
type <- "t"
set.seed(1)
rollDownhill <- function(X, Z) {
    info("rollDownhill(", X, ",", Z, ")\n", sep = "")
    if (Z == 1) 
        return(list(x = X, z = Z))
    ## Particles roll down-slope until they hit the bottom...  ... or a ledge
    ## comprising two particles.
    XX <- X
    ZZ <- Z
    while (0 == S[XX + 1, ZZ - 1]) {
        # move down and to right
        info("  XX:", XX, " ZZ:", ZZ, "\n")
        XX <- XX + 1
        ZZ <- which(0 == S[XX, ])[1]
        if (ZZ == 1) 
            break
        if (XX == n) 
            break
    }
    return(list(x = XX, z = ZZ))
}

S <- matrix(0, nrow = n, ncol = n)  # 'S' means 'space'
par(mar = c(3, 3, 1, 1), mgp = c(2, 0.7, 0))
plot(1:n, 1:n, type = "n", xlab = "", ylab = "")
xDrop <- 1  # location of drop; everything goes here or to right
for (i in 1:m) {
    # 'p' means partcle
    while (0 == length(zDrop <- which(0 == S[xDrop, ])[1])) {
        info("in while line 72\n")
        xDrop <- xDrop + 1
        if (xDrop == n) {
            message("RHS")
            break
        }
    }
    info("particle:", i, " ")
    p <- rollDownhill(xDrop, zDrop)
    S[p$x, p$z] <- 1
    if (type == "p") {
        points(p$x, p$z, col = "gray", pch = pch, cex = cex)
    } else {
        text(p$x, p$z, i, col = "gray")
    }
}

center

Discussion and conclusions

Reading the numbers on the graph as inverse age, one can see an interesting age structure.

Viewed along diagonals, ages increase by 1 time unit with every lateral step away from the source.

Viewed along Z levels, though, the time step is more interesting. You can see this at a glance, by first-differencing the values along z=1, and then at z=2, etc.

I suppose that if something came along and sliced the sediment mound along z levels, we’d see this more interesting pattern of time variation in the lateral.

I wonder if these patterns (or the code) are of interest to geologists?

Resources

To leave a comment for the author, please follow the link and comment on their blog: Dan Kelley Blog/R.

R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)