CODEWITHSUNDEEP
rmatlabwolfram-mathematicaequation
14thTimeLord
14thTimeLord

Reputation: 373

How to plot an equation with no analytical solution in R?

I have been trying to plot an equation in R. The variables are fixed at (for this instance) a=0.06, b=-0.01, i=0, t=0.005, r=0.0025, s=0.015 while I intend to vary variable e. Some of the functions in the main equations are,

functions

The equation looks like this,

enter image description here

I'd ideally like a plot between Q and x as e increases.

t = 0.0005
A = 0.06
B = -0.01
r = 0.00025
s = 0.015
i=0


M= function(e) e*A - i*B
N= function(e) t/(1+t*e-t*i)
P= function(e) (s/t)*(1+N(e)*(r+t))

I don't know how to proceed. I'd like to maybe create a list of (e, Q, x) where x is the root of the equation for a given e and then maybe use interpolation and then plot (Q, x) in R.

Is there a more direct way of plotting x vs Q? If not, what else could I try? Also, if it's easier to do with Mathematica or MATLAB, I'd be keen to know.

Upvotes: 1

Views: 160

Answers (1)

Onyambu
Onyambu

Reputation: 79338

t <- 0.0005
A <- 0.06
B <- -0.01
r <- 0.00025
s <- 0.015
i<-0

Q <- function(e){
  N <-  t/(1+t*e-t*i)
  M <- e*A - i*B
  P <- (s/t)*(1+N*(r+t))
  M/P
}


RHS <- function(x, e){
  N <-  t/(1+t*e-t*i)
  M <- e*A - i*B
  P <- (s/t)*(1+N*(r+t))
  Q <- M/P
  V <- exp(-(r+2*N)/x/r/N)
  exp(-1/(x*N)) +  V* (0.5-Q) + Q*V^2
}

LHS <- function(e){
  1 - Q(e)
}


x <- function(e){
  suppressWarnings(sapply(e, \(u)
         optimise(\(x) (RHS(x,u) - LHS(u))^2, 
                  c(-100000, 100000))$min))
}


# assume 
e <- 0.1
LHS(e)
RHS(x(e), e)


#plot:

e <- seq(0, 1000)
plot(x(e), Q(e), ty='l')

Upvotes: 2

Related Questions