Update value in function using uniroot-CodePudding

I wanna update a specific value in a function based on the result from uniroot in R.

For example, if I'm trying to solve for x, where s=60000 and t=19.95 I have:

(s-x-0.05*(0.04*x 1810.726 - mu(40,t)*(s-x)))=59789

For next iteration, I wanna update the value of s so that s=59789, as well as the value of t so that t=19.90. Repeatedly, this should be updated all the way down to t=0, with 0.05 steps downwards for t. Thus, 399 iterations, since 19.95/0.05 = 399.

(mu is just a predefined function, e.g. mu(40, 19.95) = 0.003204.)

Here is some sample code:

s <- 60000
t <- 19.95

f <- function(x) (s - x - 0.05*(0.04*x   1810.726 - mu(40, t)*(s - x)))

uniroot(f, lower=0.1, upper=100000000)$root

Could anyone give me some advice on how to implement?

The mu() is given as follows:

mu <- function(x, t) {
  A <- .00022
  B <- 2.7*10^(-6)
  c <- 1.124
  mutemp <- A   B*c^(x   t)
  out <- ifelse(t <= 2, 0.9^(2 - t)*mutemp, mutemp)
  out
}

The expected result should look s follows:

t x
0 0.0000
1 1853.8638 11 26882.9244
2 3817.7860 12 30070.8515
3 5894.9409 13 33384.7327
4 8088.4838 14 36823.9198
5 10400.9021 15 40387.3491
6 12834.5166 16 44073.5260
7 15391.4745 17 47880.5110
8 18073.7445 18 51805.9074
9 20883.1160 19 55846.8507
10 23821.2011 20 60000.0000

The value of 59789 corrsponds to V_{19.95} = 59789 and V_{20} = 60000 is a given starting value.

CodePudding user response：

This for loop may be helpful.

1. Run all of your codes

s <- 60000
t <- 20

mu <- function(x, t) {
  A <- .00022
  B <- 2.7*10^(-6)
  c <- 1.124
  mutemp <- A   B*c^(x   t)
  out <- ifelse(t <= 2, 0.9^(2 - t)*mutemp, mutemp)
  out}

f <- function(x) (s - x - 0.05*(0.04*x   1810.726 - mu(40, t)*(s - x)))

2. Run the for loop below for iteration

2.1 Predefine the length of the outcome. In your case is 400 (t/0.05 = 400).

output <- vector(mode = "numeric", length = t/0.05)

2.2 Run through the for loop from 1 to 400. Save each uniroot result to step 2.1, and then reassign both s and t accordingly.

for (i in 1:400) {
  output[i] <- uniroot(f, lower=0.1, upper=100000000)$root
  s <- output[i]
  t <- 20 - i * 0.05
}

3. Inspect the result

output

Hope this is helpful.

CodePudding user response：

You could use vapply on a defined t sequence.

s <- 6e4
tseq <- seq.int(19.95, 0, -.05)

x <- vapply(tseq, \(t) {
  s <<- uniroot(\(x) (s - x - 0.05*(0.04*x   1810.726 - mu(40, t)*(s - x))), lower=0.1, upper=100000000)$root
}, numeric(1L))

Note, that <<- changes s in the global environment, and at the end gets the last value.

s
# [1] 2072.275

res <- cbind(t=tseq, x)

head(res)
#          t        x
# [1,] 19.95 59789.92
# [2,] 19.90 59580.25
# [3,] 19.85 59371.01
# [4,] 19.80 59162.18
# [5,] 19.75 58953.77
# [6,] 19.70 58745.77