Skip to content

Commit

Permalink
working in rc76
Browse files Browse the repository at this point in the history
  • Loading branch information
abarbour committed Apr 7, 2021
1 parent 72d0787 commit 86288f8
Show file tree
Hide file tree
Showing 2 changed files with 49 additions and 45 deletions.
73 changes: 37 additions & 36 deletions R/ricecleary76.R
Original file line number Diff line number Diff line change
Expand Up @@ -6,34 +6,35 @@
#' @param Times numeric; times to calculate at \code{depth} [s]
#' @param slip numeric; relative slip along the edge dislocation [m]
#' @param radial_dist numeric; distance from xxx to xxx [m]
#' @param angFromEast
#' @param diffusiv
#' @param nuu
#' @param nu
#' @param B
#' @param mu
#' @param ...
#' @param flt_angle numeric; angle counter clockwise from fault strike
#' @param diffusiv numeric; in plane hydraulic diffusivity [m^2/s]
#' @param nuu numeric; undrained Poisson's ratio [-]
#' @param nu numeric; drained Poisson's ratio [-]
#' @param B numeric; Skempton's coefficient [-]
#' @param mu numeric; elastic shear modulus [Pa]
#' @param ... additional arguments
#'
#' @return
#' @return tibble
#' @export
#'
#' @examples
rc76_edge_disloc_poro <- function(Times, slip, radial_dist, angFromEast=30,
#' NULL
rc76_edge_disloc_poro <- function(Times, slip, radial_dist, flt_angle=30,
diffusiv=1,
nuu=0.33, nu=0.25, B = 0.6,
nuu=0.33, nu=0.25, B=0.6,
mu=30e9,
...){
# Equation X of Rice and Cleary 1976, or equations 7.149 - 7.152 of Wang (2000)
data.table::inrange(radial_dist, 0, Inf)
stopifnot(length(slip) == 1)
stopifnot(length(radial_dist) == 1)
stopifnot(length(angFromEast) == 1)
stopifnot(length(flt_angle) == 1)

zer <- Times <= 0
calc_times <- Times
calc_times[zer] <- NA

theta <- angFromEast * pi / 180
theta <- flt_angle * pi / 180
r <- radial_dist
geom_fac <- slip / r
nu_fac <- (nuu - nu) / (1 - nu) / (1 - nuu)
Expand All @@ -60,30 +61,30 @@ rc76_edge_disloc_poro <- function(Times, slip, radial_dist, angFromEast=30,

examp <- function(){

library(dplyr)
ti <- c(-100, -1, 0, 10**seq(-4, log10(6598), length.out=301))
Mw <- 4.33
M0 <- Mw2M0(Mw) #3.935501e+15 N*m for Mw 4.33
bx <- 10^(-3.22 + 0.69 * Mw) #0.5
rd <- 250
mu <- 2e8
Diffu <- 60
rc <- rc76_edge_disloc_poro(ti, bx, rd, diffusiv=Diffu, mu=mu)
nrc <- rc76_edge_disloc_poro(ti, bx, -rd, diffusiv=Diffu, mu=mu)
all.equal(rc$pp, -nrc$pp)

plot(pp / 1e3 + 5 ~ I(Times / 3600), rc, type='l', lwd=2, ylab='kPa',
ylim=c(-8.8,36.9), yaxs='i', xaxs='i')
axis(4)
abline(v=1*60/3600, lty=3)
rect(0,0,100,100)
#lines(sigma_tt / 1e3 ~ Times, rc, lty=2)
#lines(sigma_rt / 1e3 ~ Times, rc, lty=3)
#lines(sigma_rr / 1e3 ~ Times, rc)
#plot((sigma_rt - dplyr::first(na.omit(sigma_rt)))/1e3 ~ Times, rc, type='l')
SrcArea <- M0 / (mu * bx)
SrcRad <- sqrt(SrcArea / pi)
print(c(Area=SrcArea, Rad=SrcRad))
library(dplyr)
ti <- c(-100, -1, 0, 10**seq(-4, log10(6598), length.out=301))
Mw <- 4.33
M0 <- Mw2M0(Mw) #3.935501e+15 N*m for Mw 4.33
bx <- 10^(-3.22 + 0.69 * Mw) #0.5
rd <- 250
mu <- 2e8
Diffu <- 60
rc <- rc76_edge_disloc_poro(ti, bx, rd, diffusiv=Diffu, mu=mu)
nrc <- rc76_edge_disloc_poro(ti, bx, -rd, diffusiv=Diffu, mu=mu)
all.equal(rc$pp, -nrc$pp)
plot(pp / 1e3 + 5 ~ I(Times / 3600), rc, type='l', lwd=2, ylab='kPa',
ylim=c(-8.8,36.9), yaxs='i', xaxs='i')
axis(4)
abline(v=1*60/3600, lty=3)
rect(0,0,100,100)
#lines(sigma_tt / 1e3 ~ Times, rc, lty=2)
#lines(sigma_rt / 1e3 ~ Times, rc, lty=3)
#lines(sigma_rr / 1e3 ~ Times, rc)
#plot((sigma_rt - dplyr::first(na.omit(sigma_rt)))/1e3 ~ Times, rc, type='l')
SrcArea <- M0 / (mu * bx)
SrcRad <- sqrt(SrcArea / pi)
print(c(Area=SrcArea, Rad=SrcRad))

}

Expand Down
21 changes: 12 additions & 9 deletions man/rc76_edge_disloc_poro.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

0 comments on commit 86288f8

Please sign in to comment.