## ../../testts/R/arma_Q0dotdotstats.R
## Do not edit this file manually.
## It has been automatically generated from *.org sources.
library(sarima)

## "Gardner1980"
arma_Q0Gardner <- function(phi, theta, tol = .Machine$double.eps){
    ## tol is not used here
    .Call(stats:::C_getQ0, phi, theta)
}

## "Rossignol2011"
arma_Q0bis <- function(phi, theta, tol = .Machine$double.eps){
    .Call(stats:::C_getQ0bis, phi, theta, tol)
}

## arma_Q0naive <- function(phi, theta, tol = .Machine$double.eps){
##     p <- length(phi)
##     q <- length(theta)
##     r <- max(p, q+1)
##
##     T <- cbind(numeric(r), rbind(diag(r - 1), numeric(r - 1)))
##     if(p > 0)
##         T[1:p, 1] <- phi
##
##     Sigma <- matrix(0, r, r)
##     Sigma[1:(q+1), 1:(q + 1)] <- c(1, theta) %o% c(1, theta)
##
##     A <- T %x% T
##     v <- solve(diag(nrow(A)) - A, as.vector(Sigma))
##     matrix(v, r)
## }
##
## arma_Q0gnbR <- function(phi, theta, tol = .Machine$double.eps){
##     p <- length(phi)
##     q <- length(theta)
##     r <- max(p, q+1)
##
##     if(r == 1) # q = 0;  p = 0 or 1 here
##         return( matrix(if(p == 1) 1/(1-phi^2) else 1.0, 1, 1) )
##
##     rphi <- if(p == r) phi  else c(phi, numeric(r - p))
##
##     Sigma <- matrix(0, r, r)
##     Sigma[1:(q+1), 1:(q + 1)] <- c(1, theta) %o% c(1, theta)
##
##     A <- diag(r)
##     b <- numeric(r)
##     for(k in seq_len(r)){
##         b[k] <- Sigma[1, k]
##         A[k,2] <- A[k,2] - rphi[k]
##             # ak1 <- rphi[r - k + 1] * rphi[k + (r - k + 1) - 1]
##         ak1 <- rphi[r - k + 1] * rphi[r]
##         if(k < r){
##             for(i in 1:(r - k)){
##                 b[k] <- b[k] + Sigma[1 + i, k + i]
##                 A[k, k + i] <- A[k, k + i] - rphi[i]
##                 if(2 + i <= r)
##                     A[k, 2 + i] <- A[k, 2 + i] - rphi[k + i]
##
##                 ak1 <- ak1 + rphi[i] * rphi[k + i - 1]
##             }
##         }
##         A[k, 1] <- A[k, 1] - ak1
##     }
##
##     g1 <- solve(A, b)
##
##     res <- matrix(0.0, r, r)
##     res[1, ] <- g1
##     res[r, r] <- rphi[r]^2 * g1[1] + Sigma[r, r]
##     if(r > 2){
##         for(i in (r-1):2){
##             res[i, r] <- rphi[i] * rphi[r] * g1[1] + rphi[r] * g1[i + 1] + Sigma[i, r]
##             for(j in (r-1):i){
##                 res[i, j] <- rphi[i] * rphi[j] * g1[1] +
##                     rphi[i] * g1[j + 1] + rphi[j] * g1[i + 1] +
##                     Sigma[i,j] + res[i + 1, j + 1]
##             }
##         }
##     }
##     ## fill the elements under the diagonal, symmetrically;
##     ## r > 1 here, so we don't check
##     res[2,1] <- res[1,2]
##     if(r > 2){
##         for(i in (1:(r - 1)))
##             for(j in (i + 1) : r)
##                 res[j,i] <- res[i, j]
##     }
##
##     res
## }