# Load the required packages
library(deSolve)
library(ggplot2)

# Define the function for the differential equation system
water_tank <- function(time, state, parameters) {
  water_level <- state[1]
  inflow_rate <- parameters[1]
  outflow_rate <- parameters[2]
  
  if (water_level < 90)
    inflow_rate <- 5
  else if (water_level >= 90)
    inflow_rate <- 0
  
  net_flow_rate <- inflow_rate - outflow_rate
  d_water_level <- net_flow_rate
  
  return(list(d_water_level))
}

# Set the initial conditions and parameters
initial_state <- c(water_level = 100)
parameters <- c(inflow_rate = 0, outflow_rate = 1)

# Define the time span of the simulation
times <- seq(0, 200, by = 0.1)

# Simulate the system using the ode() function from deSolve
output <- ode(y = initial_state, times = times, func = water_tank, parms = parameters)

tail(output)
# Convert the output data to a data frame
df <- as.data.frame(output)

# Plot the water level over time using ggplot
p <- ggplot(data = df, aes(x = time, y = water_level)) +
  geom_line() +
  xlab("Time (s)") +
  ylab("Water Level (liters)")

# Save the plot as a PDF file
ggsave(filename = "water_tank_plot.pdf", plot = p, device = "pdf")