R: Variance of ARIMA as a Starting Point for Simulating from ARIMA Model

Question

I want to simulate time series data from ARIMA(1, 1, 1) manually (without using arima.sim() function) using R. I can do that for ARIMA(1, 0, 0), ARIMA(2, 0, 0), ARIMA(0, 0, 1), ARIMA(0, 0, 2) and ARIMA(1, 0, 1) as follows:

# simulate ARIMA(1, 0, 0)
n <- 50
phi <- 0.5
set.seed(1)
wn <- rnorm(n, mean=0, sd=1)
v.ar1 <- sqrt((wn[1])^2/(1-phi^2))
for(i in 2:n){
  ar1[i] <- v.ar1[i - 1] * phi + wn[i]
}

# simulate ARIMA(2, 0, 0)
n <- 50
phi1 <- 0.5
phi2 <- -0.7
set.seed(1)
wn <- rnorm(n, mean=0, sd=1)
x1 <- sqrt(1/(1-phi1^2))
x2 <- sqrt((1-phi2^2)/((1-phi2)*((1-phi2)^2-phi1^2)))
ar1 <- c(x1,x2)
for(i in 3:n){
  ar1[i] <- ar1[i - 1] * phi1 + ar1[i - 2]* phi2 + wn[i]
}

# simulate ARIMA(0, 0, 1)
n=50
theta <- 0.5
set.seed(1)
wn <- rnorm(n, mean=0, sd=1)
y1 <- sqrt(1*(1+theta^2))
ma1 <- y1
for(i in 2:n){
  ma1[i] <- wn[i - 1] * theta + wn[i]
}

# simulate ARIMA(0, 0, 2)
n=50
theta1 <- 0.5
theta2 <- 0.2
set.seed(1)
wn <- rnorm(n, mean=0, sd=1)
y1 <- sqrt(1*(1+theta1^2))
y2 <- sqrt(1*(1+theta1^2+theta2^2))
ma2 <- c(y1,y2)
for(i in 3:n){
  ma2[i] <- ma2[i - 1] * theta1 + ma2[i - 2]* theta2 + wn[i]
}

# simulate ARIMA(1, 0, 1)
n <- 500
phi1 <- 0.9
theta1 <- 0.5
set.seed(5)
wn <- rnorm(n,mean=0,sd=1)
v.ARMA11 <- sqrt((1+2*phi1*theta1+theta1^2)*1/(1-phi1^2))
ARMA11 <- v.ARMA11
for(i in 2:n){
  ARMA11[i] <- phi1*z[i-1] +theta1*wn[i - 1] + wn[i]
}

The variance of the ARIMA model one is simulating from seems to be the starting point of the randomly generated number.

WHAT I WANT

What is the variance of the ARIMA(1, 1, 1) so I can use it as my starting point?

EDITED

The reason I am looking away from arima.sim() function I that so I can simulate a non-stationary series by making coefficient of AR to be grater than 1 $\phi > 1$