diff --git a/NEWS.md b/NEWS.md
index 056cc06..f83e5a1 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -39,6 +39,9 @@
    Because there was no recorded direct usage for any of these, I am opting to just rip the band-aid
    off and un-export them in this release as opposed to waiting a full cycle more to do so.
 
+1. A replacement in an integer64 vector using `[<-` or `[[<-` with a character leads to an R consistent 
+   coercion of the integer64 object to a character object.
+   
 ## NEW FEATURES
 
 1. `anyNA` gets an `integer64` method. Thanks @hcirellu.
@@ -46,6 +49,7 @@
 ## BUG FIXES
 
 1. `min.integer64`, `max.integer64` and `range.integer64` now support `na.rm=TRUE` correctly when combining across mutliple inputs like `min(x, NA_integer64_, na.rm=TRUE)` (#142).
+1. `[.integer64` now runs faster and correctly regarding `NA` and arrays. (#176)
 
 ## NOTES
 
diff --git a/R/integer64.R b/R/integer64.R
index 46c70a1..a065dfa 100644
--- a/R/integer64.R
+++ b/R/integer64.R
@@ -104,14 +104,15 @@ NULL
 #' Methods to extract and replace parts of an integer64 vector.
 #'
 #' @param x an atomic vector
-#' @param i indices specifying elements to extract
+#' @param i,j indices specifying elements to extract
+#' @param drop relevant for matrices and arrays. If TRUE the result is coerced to the lowest possible dimension.
 #' @param value an atomic vector with values to be assigned
 #' @param ... further arguments to the [NextMethod()]
 #'
 #' @note
 #'   You should not subscript non-existing elements and not use `NA`s as subscripts.
 #'   The current implementation returns `9218868437227407266` instead of `NA`.
-#' @returns A vector or scalar of class 'integer64'
+#' @returns A vector, matrix, array or scalar of class 'integer64'
 #' @keywords classes manip
 #' @seealso [`[`][base::Extract] [integer64()]
 #' @examples
@@ -846,90 +847,122 @@ str.integer64 <- function(object,
 
 #' @rdname extract.replace.integer64
 #' @export
-`[.integer64` <- function(x, i, ...) {
-    cl <- oldClass(x)
-    ret <- NextMethod()
-    # Begin NA-handling from Leonardo Silvestri
-    if (!missing(i)) {
-        if (inherits(i, "character")) {
-          na_idx <- union(which(!(i %in% names(x))), which(is.na(i)))
-          if (length(na_idx))
-                ret[na_idx] <- NA_integer64_
-        } else {
-      ni <- length(i)
-      nx <- length(x)
-      if (inherits(i, "logical")) {
-            if (ni>nx) {
-              na_idx <- is.na(i) | (i & seq_along(i)>nx)
-              na_idx <- na_idx[is.na(i) | i]
-            } else {
-          i <- i[is.na(i) | i]
-          na_idx <- rep_len(is.na(i), length(ret))
-            }
-          } else if (ni && min(i, na.rm=TRUE)>=0L) {
-            i <- i[is.na(i) | i>0L]
-            na_idx <- is.na(i) | i>length(x)
-          } else {
-            na_idx <- FALSE
-          }
-          if (any(na_idx))
-                ret[na_idx] <- NA_integer64_
-        }
+`[.integer64` = function(x, i, j, ..., drop=TRUE) {
+  args = lapply(as.list(sys.call())[-(1:2)], {function(el) {
+    if(is.symbol(el) && el == substitute()) return(el)
+    el = eval(el, parent.frame(3L))   
+    if (is.integer64(el))
+      el = as.integer(el)
+    el
+  }})
+  args$drop = FALSE
+  if (length(args) == 1L) return(x)
+  oldClass(x) = NULL
+  ret = do.call("[", c(list(x=x), args))
+  NA_integer64_real = NA_integer64_
+  oldClass(NA_integer64_real) = NULL
+
+  # NA handling
+  if (length(dim(ret)) <= 1L) {
+    # vector mode
+    if (!is.symbol(args[[1L]]) || args[[1L]] != substitute()) {
+      arg1Value = args[[1L]]
+      if (is.logical(arg1Value)) {
+        ret[is.na(arg1Value[arg1Value])] = NA_integer64_real
+      } else if (is.character(arg1Value)) {
+        ret[is.na(arg1Value) | arg1Value == "" | !arg1Value %in% names(x)] = NA_integer64_real
+      } else if (anyNA(arg1Value) || suppressWarnings(max(arg1Value, na.rm=TRUE)) > length(x)) {
+        arg1Value = arg1Value[arg1Value != 0]
+        ret[which(is.na(arg1Value) | arg1Value > length(x))] = NA_integer64_real
+      }
     }
-    # End NA-handling from Leonardo Silvestri
-    oldClass(ret) <- cl
-    remcache(ret)
-    ret
-}
-
-
-`[.integer64` <- function(x, i, ...) {
-  cl <- oldClass(x)
-  ret <- NextMethod()
-  # Begin NA-handling from Leonardo Silvestri
-  if (!missing(i)) {
-    if (inherits(i, "character")) {
-      na_idx <- union(which(!(i %in% names(x))), which(is.na(i)))
-      if (length(na_idx))
-        ret[na_idx] <- NA_integer64_
-    } else {
-      na_idx <- is.na(rep(TRUE, length(x))[i])
-      if (any(na_idx))
-        ret[na_idx] <- NA_integer64_
+  } else {
+    # array/matrix mode
+    dimSelect = args[seq_along(dim(x))]
+    for (ii in seq_along(dimSelect)) {
+      if (is.symbol(dimSelect[[ii]]) && dimSelect[[ii]] == substitute()) next
+      dsValue = dimSelect[[ii]]
+      if (is.logical(dsValue) && anyNA(dsValue)) {
+        naIndex = which(is.na(seq_len(dim(x)[ii])[dsValue]))
+      } else {
+        naIndex = which(is.na(dsValue[dsValue != 0L]))
+      }
+      if (length(naIndex)) {
+        setArgs = rep(list(substitute()), length(dimSelect))
+        setArgs[[ii]] = naIndex
+        ret = do.call("[<-", c(list(x=ret), setArgs, list(value=NA_integer64_real)))
+      }
     }
   }
-  # End NA-handling from Leonardo Silvestri
-  oldClass(ret) <- cl
-  remcache(ret)
+
+  # dimension handling
+  if (!isFALSE(drop) && !is.null(dim(ret))) {
+    newDim = dim(ret)[dim(ret) != 1L]
+    dim(ret) = {if (length(newDim)) newDim else NULL}
+    if(length(dim(ret)) <= 1L)
+      dim(ret) = NULL
+  }
+
+  oldClass(ret) = "integer64"
   ret
 }
 
 #' @rdname extract.replace.integer64
 #' @export
-`[<-.integer64` <- function(x, ..., value) {
-  cl <- oldClass(x)
-  value <- as.integer64(value)
-  ret <- NextMethod()
-  oldClass(ret) <- cl
+`[<-.integer64` = function(x, ..., value) {
+  sc = as.list(sys.call())
+  args = lapply(sc[-c(1:2, length(sc))], {function(el) {
+    if(is.symbol(el) && el == substitute()) return(el)
+    el = eval(el, parent.frame(3L))   
+    if (is.integer64(el))
+      el = as.integer(el)
+    el
+  }})
+  if (is.character(value) || is.complex(value) || (is.double(value) && class(value)[1L] != "numeric")) {
+    args$value = value
+    x = structure(as(x, class(value)[1L]), dim = dim(x), dimnames = dimnames(x))
+    ret = do.call("[<-", c(list(x=x), args))
+  } else {
+    args$value = as.integer64(value)
+    oldClass(x) = NULL
+    ret = do.call("[<-", c(list(x=x), args))
+    oldClass(ret) = "integer64"
+  }
   ret
 }
 
 #' @rdname extract.replace.integer64
 #' @export
-`[[.integer64` <- function(x, ...) {
-  cl <- oldClass(x)
-  ret <- NextMethod()
-  oldClass(ret) <- cl
+`[[.integer64` = function(x, ...) {
+  args = lapply(list(...), {function(el) {
+    if (is.integer64(el))
+      el = as.integer(el)
+    el
+  }})
+  oldClass(x) = NULL
+  ret = do.call("[[", c(list(x=x), args))
+  oldClass(ret) = "integer64"
   ret
 }
 
 #' @rdname extract.replace.integer64
 #' @export
-`[[<-.integer64` <- function(x, ..., value) {
-  cl <- oldClass(x)
-  value <- as.integer64(value)
-  ret <- NextMethod()
-  oldClass(ret) <- cl
+`[[<-.integer64` = function(x, ..., value) {
+  args = lapply(list(...), {function(el) {
+    if (is.integer64(el))
+      el = as.integer(el)
+    el
+  }})
+  if (is.character(value) || is.complex(value) || (is.double(value) && class(value)[1L] != "numeric")) {
+    args$value = value
+    x = structure(as(x, class(value)[1L]), dim = dim(x), dimnames = dimnames(x))
+    ret = do.call("[[<-", c(list(x=x), args))
+  } else {
+    value = as.integer64(value)
+    oldClass(x) = NULL
+    ret = do.call("[[<-", c(list(x=x), args, list(value=value)))
+    oldClass(ret) = "integer64"
+  }
   ret
 }
 
diff --git a/man/extract.replace.integer64.Rd b/man/extract.replace.integer64.Rd
index 78b0166..c88b650 100644
--- a/man/extract.replace.integer64.Rd
+++ b/man/extract.replace.integer64.Rd
@@ -8,7 +8,7 @@
 \alias{[[<-.integer64}
 \title{Extract or Replace Parts of an integer64 vector}
 \usage{
-\method{[}{integer64}(x, i, ...)
+\method{[}{integer64}(x, i, j, ..., drop = TRUE)
 
 \method{[}{integer64}(x, ...) <- value
 
@@ -19,14 +19,16 @@
 \arguments{
 \item{x}{an atomic vector}
 
-\item{i}{indices specifying elements to extract}
+\item{i, j}{indices specifying elements to extract}
 
 \item{...}{further arguments to the \code{\link[=NextMethod]{NextMethod()}}}
 
+\item{drop}{relevant for matrices and arrays. If TRUE the result is coerced to the lowest possible dimension.}
+
 \item{value}{an atomic vector with values to be assigned}
 }
 \value{
-A vector or scalar of class 'integer64'
+A vector, matrix, array or scalar of class 'integer64'
 }
 \description{
 Methods to extract and replace parts of an integer64 vector.
@@ -38,10 +40,10 @@ The current implementation returns \code{9218868437227407266} instead of \code{N
 \examples{
   as.integer64(1:12)[1:3]
   x <- as.integer64(1:12)
-  dim(x) <- c(3,4)
+  dim(x) <- c(3, 4)
   x
   x[]
-  x[,2:3]
+  x[, 2:3]
 }
 \seealso{
 \code{\link[base:Extract]{[}} \code{\link[=integer64]{integer64()}}
diff --git a/tests/testthat/helper.R b/tests/testthat/helper.R
index 04a0d39..e32fe8d 100644
--- a/tests/testthat/helper.R
+++ b/tests/testthat/helper.R
@@ -1,13 +1,15 @@
 # TODO(#45): use matrix() directly
-matrix64 = function(x, nrow=1L, ncol=1L, byrow=FALSE) {
+matrix64 = function(x, nrow=1L, ncol=1L, byrow=FALSE, dimnames=NULL) {
   x = as.integer64(x)
   if (byrow) {
     dim(x) = c(ncol, nrow)
-    t(x)
+    x = t(x)
+    dimnames(x) = dimnames
   } else {
     dim(x) = c(nrow, ncol)
-    x
+    dimnames(x) = dimnames
   }
+  x
 }
 
 array64 = function(x, dim) {
diff --git a/tests/testthat/test-bit64-package.R b/tests/testthat/test-bit64-package.R
index 2db72f3..423905a 100644
--- a/tests/testthat/test-bit64-package.R
+++ b/tests/testthat/test-bit64-package.R
@@ -55,28 +55,28 @@ test_that("Minus and plus edge cases and 'rev'", {
 
 test_that("'range.integer64', multiplication, integer division, sqrt, power, and log", {
   i64 = integer64(63L)
-  i64[1L] = 1.0
+  i64[1L] = 1L
   for (i in 2:63)
-    i64[i] = 2.0 * i64[i-1L]
+    i64[i] = 2L * i64[i-1L]
   expect_true(identical.integer64(i64 * rev(i64), rep(i64[63L], 63L)))
   for (i in 63:2)
-    i64[i-1L] = i64[i] %/% 2.0
+    i64[i-1L] = i64[i] %/% 2L
   expect_true(identical.integer64(i64 * rev(i64), rep(i64[63L], 63L)))
   for (i in 63:2)
-    i64[i-1L] = i64[i] / 2.0
+    i64[i-1L] = as.integer64(i64[i] / 2L)
   expect_true(identical.integer64(i64 * rev(i64), rep(i64[63L], 63L)))
   expect_true(identical.integer64(
     c(
-      -i64[63L] - (i64[63L] - 1.0),
-      i64[63L] + (i64[63L] - 1.0)
+      -i64[63L] - (i64[63L] - 1L),
+      i64[63L] + (i64[63L] - 1L)
     ),
     lim.integer64()
   ))
 
-  expect_true(identical.integer64(i64[-1L] %/%2.0 * as.integer64(2L), i64[-1L]))
-  expect_true(identical.integer64(i64[-1L] %/%2L * as.integer64(2L), i64[-1L]))
-  expect_true(identical.integer64(i64[-1L] / 2.0 * as.integer64(2L), i64[-1L]))
+  expect_true(identical.integer64(i64[-1L] %/% 2.0 * as.integer64(2L), i64[-1L]))
+  expect_true(identical.integer64(i64[-1L] %/% 2L * as.integer64(2L), i64[-1L]))
   expect_true(identical.integer64(i64[-1L] / 2.0 * as.integer64(2L), i64[-1L]))
+  expect_true(identical.integer64(i64[-1L] / 2L * as.integer64(2L), i64[-1L]))
 
   expect_true(identical.integer64(i64[-63L] * 2.0 %/% 2.0, i64[-63L]))
   expect_true(identical.integer64(i64[-63L] * 2L %/% 2L, i64[-63L]))
diff --git a/tests/testthat/test-integer64.R b/tests/testthat/test-integer64.R
index c66d8d3..a0b5801 100644
--- a/tests/testthat/test-integer64.R
+++ b/tests/testthat/test-integer64.R
@@ -40,22 +40,62 @@ test_that("S3 class basics work", {
 })
 
 test_that("indexing works", {
-  x = as.integer64(1:10)
 
-  x[1.0] = 2.0
+  x = as.integer64(1:10)
+  x[1.0] = 2L
   x[2L] = 3L
   expect_identical(x, as.integer64(c(2:3, 3:10)))
 
-  x[[1.0]] = 3.0
+  x = as.integer64(1:10)
+  x[1.0] = 2L
+  x[2L] = 3.0
+  expect_identical(x, as.integer64(c(2:3, 3:10)))
+
+  # has to be commented until as.complex.integer64 exists
+  # x = as.integer64(1:10)
+  # x[1.0] = 2L
+  # x[2L] = 3+0i
+  # expect_identical(x, as.complex(c(2:3, 3:10)))
+  
+  x = as.integer64(1:10)
+  x[1.0] = 2L
+  x[2L] = "3"
+  expect_identical(x, as.character(c(2:3, 3:10)))
+  
+  x = as.integer64(1:10)
+  x[[1.0]] = 3L
   x[[2L]] = 4L
   expect_identical(x, as.integer64(c(3:4, 3:10)))
 
+  x = as.integer64(1:10)
+  x[[1.0]] = 3L
+  x[[2L]] = 4.0
+  expect_identical(x, as.integer64(c(3:4, 3:10)))
+
+  # has to be commented until as.complex.integer64 exists
+  # x = as.integer64(1:10)
+  # x[[1.0]] = 3L
+  # x[[2L]] = 4+0i
+  # expect_identical(x, as.complex(c(3:4, 3:10)))
+  
+  x = as.integer64(1:10)
+  x[[1.0]] = 3L
+  x[[2L]] = "4"
+  expect_identical(x, as.character(c(3:4, 3:10)))
+  
+  x = as.integer64(1:10)
   expect_identical(x[3L], as.integer64(3L))
   expect_identical(x[[4L]], as.integer64(4L))
 
   names(x) = letters[1:10]
   expect_identical(x[c("b", "c")], x[2:3])
   expect_identical(x[["d"]], x[[4L]])
+  
+  expect_no_warning(expect_identical(integer64()[integer()], integer64()))
+  expect_no_warning(expect_identical(structure(as.integer64(1L), dim=c(1))[1], as.integer64(1L)))
+
+  expect_no_warning(expect_identical(as.integer64(1L)[NA_integer_], NA_integer64_))
+  expect_no_warning(expect_identical(as.integer64(1L)[NA_integer64_], NA_integer64_))
 })
 
 test_that("arithmetic & basic math works", {
@@ -493,3 +533,236 @@ test_that("match works with zero length input", {
   expect_identical(match(x64, integer(), nomatch=10L), match(x32, integer(), nomatch=10L))
   expect_identical(match(integer(), x64), match(integer(), x32))
 })
+
+
+test_that("extraction and replacement works consistent to integer (except for double)", {
+
+  skip_if_not_r_version("4.0.0")
+  # extraction with `[`
+  x = as.integer(1:10)
+  names(x) = letters[seq_along(x)]
+  y = as.integer64(x)
+  names(y) = letters[seq_along(y)]
+  sel = c(TRUE, FALSE, NA, TRUE)
+  expect_identical(y[sel], structure(as.integer64(x[sel]), names=names(x)[sel]))
+  sel = c(1, NA, 3, 11)
+  expect_identical(y[sel], structure(as.integer64(x[sel]), names=names(x)[sel]))
+  expect_identical(y[as.integer64(sel)], structure(as.integer64(x[sel]), names=names(x)[sel]))
+  sel = c(-1, -3, 0, -11)
+  expect_identical(y[sel], structure(as.integer64(x[sel]), names=names(x)[sel]))
+  sel = c(-1, -3, 0, -11, NA)
+  expect_error(x[sel], "only 0's may be mixed with negative subscripts", fixed=TRUE)
+  expect_error(y[sel], "only 0's may be mixed with negative subscripts", fixed=TRUE)
+
+  expect_identical(as.integer64(c("9218868437227407266", "1"))[c(1,NA,3,4)], structure(as.integer64(c("9218868437227407266", NA_character_, NA_character_, NA_character_))))
+
+  sel = c("d", "", "b", NA_character_)
+  expect_identical(y[sel], structure(as.integer64(x[sel]), names=names(x)[match(sel, names(x))]))
+
+  m32 = matrix(1:10, nrow=2L)
+  m64 = matrix64(as.integer64(m32), nrow=dim(m32)[1L], ncol=dim(m32)[2L])
+  expect_identical(m32[integer(), 1:2, drop=TRUE], structure(integer(), dim = c(0L, 2L)))
+  expect_identical(m64[integer(), 1:2, drop=TRUE], structure(integer64(), dim = c(0L, 2L)))
+
+  expect_identical(m32[1:2, integer(), drop=TRUE], structure(integer(), dim = c(2L, 0L)))
+  expect_identical(m64[1:2, integer(), drop=TRUE], structure(integer64(), dim = c(2L, 0L)))
+
+  expect_identical(m32[integer(), 1:2, drop=FALSE], structure(integer(), dim = c(0L, 2L)))
+  expect_identical(m64[integer(), 1:2, drop=FALSE], structure(integer64(), dim = c(0L, 2L)))
+
+  expect_identical(m32[1:2, integer(), drop=FALSE], structure(integer(), dim = c(2L, 0L)))
+  expect_identical(m64[1:2, integer(), drop=FALSE], structure(integer64(), dim = c(2L, 0L)))
+
+  expect_identical(m32[1:2, 1:3, drop=TRUE], structure(as.integer(1:6), dim = c(2L, 3L)))
+  expect_identical(m64[1:2, 1:3, drop=TRUE], structure(as.integer64(1:6), dim = c(2L, 3L)))
+
+  expect_identical(m32[1:2], structure(as.integer(1:2)))
+  expect_identical(m64[1:2], structure(as.integer64(1:2)))
+
+  expect_identical(m32[1:2, drop=TRUE], structure(as.integer(1:2)))
+  expect_identical(m64[1:2, drop=TRUE], structure(as.integer64(1:2)))
+
+  expect_identical(m32[j = 1:3, drop=TRUE], structure(as.integer(1:3)))
+  expect_identical(m64[j = 1:3, drop=TRUE], structure(as.integer64(1:3)))
+
+  expect_identical(m32[1:2, , drop=TRUE], structure(as.integer(1:10), dim = c(2L, 5L)))
+  expect_identical(m64[1:2, , drop=TRUE], structure(as.integer64(1:10), dim = c(2L, 5L)))
+
+  expect_identical(m32[, 1:3, drop=TRUE], structure(as.integer(1:6), dim = c(2L, 3L)))
+  expect_identical(m64[, 1:3, drop=TRUE], structure(as.integer64(1:6), dim = c(2L, 3L)))
+
+  expect_identical(m32[1, , drop=TRUE], structure(as.integer(c(1L, 3L, 5L, 7L, 9L))))
+  expect_identical(m64[1, , drop=TRUE], structure(as.integer64(c(1L, 3L, 5L, 7L, 9L))))
+
+  expect_identical(m32[1, , drop=FALSE], structure(as.integer(c(1L, 3L, 5L, 7L, 9L)), dim = c(1L, 5L)))
+  expect_identical(m64[1, , drop=FALSE], structure(as.integer64(c(1L, 3L, 5L, 7L, 9L)), dim = c(1L, 5L)))
+
+  expect_identical(m32[, 1, drop=TRUE], structure(as.integer(1:2)))
+  expect_identical(m64[, 1, drop=TRUE], structure(as.integer64(1:2)))
+
+  expect_identical(m32[, 1, drop=FALSE], structure(as.integer(1:2), dim = 2:1))
+  expect_identical(m64[, 1, drop=FALSE], structure(as.integer64(1:2), dim = 2:1))
+
+  expect_identical(m32[c(9, NA, 11, 12), drop=FALSE], structure(as.integer(c(9L, NA, NA, NA))))
+  expect_identical(m64[c(9, NA, 11, 12), drop=FALSE], structure(as.integer64(c(9L, NA, NA, NA))))
+
+  expect_identical(m32[integer(), c(1:2, 0, NA), drop=TRUE], structure(integer(), dim = c(0L, 3L)))
+  expect_identical(m64[integer(), c(1:2, 0, NA), drop=TRUE], structure(integer64(), dim = c(0L, 3L)))
+  expect_identical(m64[integer64(), c(1:2, 0, NA), drop=TRUE], structure(integer64(), dim = c(0L, 3L)))
+
+  expect_identical(m32[, c(1:2, 0, NA), drop=TRUE], structure(as.integer(c(1:4, NA, NA)), dim = c(2L, 3L)))
+  expect_identical(m64[, c(1:2, 0, NA), drop=TRUE], structure(as.integer64(c(1:4, NA, NA)), dim = c(2L, 3L)))
+
+  expect_identical(m32[c(1, NA, 2), 1:3, drop=TRUE], structure(as.integer(c(1L, NA, 2L, 3L, NA, 4L, 5L, NA, 6L)), dim = c(3L, 3L)))
+  expect_identical(m64[c(1, NA, 2), 1:3, drop=TRUE], structure(as.integer64(c(1L, NA, 2L, 3L, NA, 4L, 5L, NA, 6L)), dim = c(3L, 3L)))
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+
+  expect_error(m32[c("B", "D", "A"), c("d", "a")], "subscript out of bounds", fixed=TRUE)
+  expect_error(m64[c("B", "D", "A"), c("d", "a")], "subscript out of bounds", fixed=TRUE)
+
+  expect_identical(m32[c("B", "D", "A")], rep(NA_integer_, 3L))
+  expect_identical(m64[c("B", "D", "A")], rep(NA_integer64_, 3L))
+
+  a32 = array(as.integer(1:27), c(3,3,3))
+  a64 = array64(as.integer64(1:27), c(3,3,3))
+
+  expect_identical(a32[2, , 3, drop=FALSE], structure(as.integer(c(20L, 23L, 26L)), dim = c(1L, 3L, 1L)))
+  expect_identical(a64[2, , 3, drop=FALSE], structure(as.integer64(c(20L, 23L, 26L)), dim = c(1L, 3L, 1L)))
+
+  expect_identical(a32[2, , 3, drop=TRUE], structure(as.integer(c(20L, 23L, 26L))))
+  expect_identical(a64[2, , 3, drop=TRUE], structure(as.integer64(c(20L, 23L, 26L))))
+
+  expect_identical(a32[1, c(1, 3, 2), 2:3, drop=TRUE], structure(as.integer(c(10L, 16L, 13L, 19L, 25L, 22L)), dim = 3:2))
+  expect_identical(a64[1, c(1, 3, 2), 2:3, drop=TRUE], structure(as.integer64(c(10L, 16L, 13L, 19L, 25L, 22L)), dim = 3:2))
+
+  expect_identical(a32[, c(1, 2, 0, 3, NA, 1), c(TRUE, FALSE, NA), drop=FALSE], structure(as.integer(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, NA, NA, NA, 1L, 2L, 3L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), dim = c(3L, 5L, 2L)))
+  expect_identical(a64[, c(1, 2, 0, 3, NA, 1), c(TRUE, FALSE, NA), drop=FALSE], structure(as.integer64(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, NA, NA, NA, 1L, 2L, 3L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), dim = c(3L, 5L, 2L)))
+
+  expect_identical(a32[, c(1, 2, 0, 3, NA, 1), c(TRUE, FALSE, NA), drop=TRUE], structure(as.integer(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, NA, NA, NA, 1L, 2L, 3L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), dim = c(3L, 5L, 2L)))
+  expect_identical(a64[, c(1, 2, 0, 3, NA, 1), c(TRUE, FALSE, NA), drop=TRUE], structure(as.integer64(c(1L, 2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, NA, NA, NA, 1L, 2L, 3L, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), dim = c(3L, 5L, 2L)))
+
+  expect_identical(a32[c(1, 0, 7, NA, 27, 28), drop=FALSE], structure(as.integer(c(1L, 7L, NA, 27L, NA))))
+  expect_identical(a64[c(1, 0, 7, NA, 27, 28), drop=FALSE], structure(as.integer64(c(1L, 7L, NA, 27L, NA))))
+
+  expect_identical(a32[c(TRUE, FALSE, NA, TRUE), drop=FALSE], structure(as.integer(c(1L, NA, 4L, 5L, NA, 8L, 9L, NA, 12L, 13L, NA, 16L, 17L, NA, 20L, 21L, NA, 24L, 25L, NA))))
+  expect_identical(a64[c(TRUE, FALSE, NA, TRUE), drop=FALSE], structure(as.integer64(c(1L, NA, 4L, 5L, NA, 8L, 9L, NA, 12L, 13L, NA, 16L, 17L, NA, 20L, 21L, NA, 24L, 25L, NA))))
+
+  expect_identical(a32[-1, , -c(0, 2:3), drop=FALSE], structure(as.integer(c(2L, 3L, 5L, 6L, 8L, 9L)), dim = c(2L, 3L, 1L)))
+  expect_identical(a64[-1, , -c(0, 2:3), drop=FALSE], structure(as.integer64(c(2L, 3L, 5L, 6L, 8L, 9L)), dim = c(2L, 3L, 1L)))
+
+  expect_identical(a32[-1, 2, -c(0, 2:3), drop=FALSE], structure(as.integer(5:6), dim = c(2L, 1L, 1L)))
+  expect_identical(a64[-1, 2, -c(0, 2:3), drop=FALSE], structure(as.integer64(5:6), dim = c(2L, 1L, 1L)))
+
+  expect_identical(a32[-1, 2, -c(0, 2:3), drop=TRUE], structure(as.integer(5:6)))
+  expect_identical(a64[-1, 2, -c(0, 2:3), drop=TRUE], structure(as.integer64(5:6)))
+
+  # replacement with `[<-`
+  x = as.integer(1:10)
+  names(x) = letters[seq_along(x)]
+  y = as.integer64(x)
+  names(y) = letters[seq_along(y)]
+
+  sel = c("d", "", "b", NA_character_)
+  x[sel] = 100L
+  y[sel] = 100L
+  expect_identical(y, structure(as.integer64(x), names = names(x)))
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+
+  m32[1, c(1, 3, NA)] = 100L
+  m64[1, c(1, 3, NA)] = as.integer64(100L)
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  m32[1, c(1, 4, NA)] = 101L
+  m64[1, c(1, 4, NA)] = 101L
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  m32[1, c(1, 5, NA)] = 102
+  m64[1, c(1, 5, NA)] = 102
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  # has to be commented until as.complex.integer64 exists
+  # m32[1, c(1, 5, NA)] = 102+0i
+  # m64[1, c(1, 5, NA)] = 102+0i
+  # expect_identical(m64, m32)
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m32[1, c(1, 3, NA)] = "103"
+  m64[1, c(1, 3, NA)] = "103"
+  expect_identical(m64, m32)
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m32[1, c(1, 3, NA)] = 101L
+  m64[1, as.integer64(c(1, 3, NA))] = 101L
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  m32[, -(1:3)] = 102L
+  m64[, -(1:3)] = 102L
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  # extraction with `[[`  
+  x = as.integer(1:10)
+  names(x) = letters[seq_along(x)]
+  y = as.integer64(x)
+  names(y) = letters[seq_along(y)]
+  expect_identical(y[[3]], as.integer64(x[[3]]))
+  expect_identical(y[["d"]], as.integer64(x[["d"]]))
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  expect_identical(m64[[1, 2]], as.integer64(m32[[1, 2]]))
+  expect_identical(m64[[as.integer64(1L), as.integer64(2L)]], as.integer64(m32[[1, 2]]))
+  expect_identical(m64[["A", "d"]], as.integer64(m32[["A", "d"]]))
+
+  expect_identical(m64[[1]], as.integer64(m32[[1]]))
+  expect_identical(m64[[as.integer64(1L)]], as.integer64(m32[[1]]))
+
+  expect_error(m32[[NA]], "subscript out of bounds", fixed=TRUE)
+  expect_error(m64[[NA]], "subscript out of bounds", fixed=TRUE)
+  expect_error(m64[[as.integer64(NA)]], "subscript out of bounds", fixed=TRUE)
+
+  expect_error(m32[[0L]], "attempt to select less than one element in integerOneIndex", fixed=TRUE)
+  expect_error(m64[[0L]], "attempt to select less than one element in integerOneIndex", fixed=TRUE)
+  expect_error(m64[[as.integer64(0L)]], "attempt to select less than one element in integerOneIndex", fixed=TRUE)
+
+  expect_error(m32[[integer()]], "attempt to select less than one element in get1index", fixed=TRUE)
+  expect_error(m64[[integer()]], "attempt to select less than one element in get1index", fixed=TRUE)
+  expect_error(m64[[as.integer64()]], "attempt to select less than one element in get1index", fixed=TRUE)
+
+  # replacement with `[[<-`
+  x[["e"]] = 100L
+  y[["e"]] = 100L
+  expect_identical(y, structure(as.integer64(x), names = names(x)))
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+
+  m32[[1, 3]] = 110L
+  m64[[1, 3]] = 110L
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  m32[["A", "e"]] = 112L
+  m64[["A", "e"]] = 112L
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  m32[[1, 3]] = 111
+  m64[[1, 3]] = 111
+  expect_identical(m64, structure(as.integer64(m32), dim = dim(m32), dimnames = dimnames(m32)))
+
+  # has to be commented until as.complex.integer64 exists
+  # m32[[1, 3]] = 111+0i
+  # m64[[1, 3]] = 111+0i
+  # expect_identical(m64, m32)
+
+  m32 = matrix(1:10, 2L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m64 = matrix64(as.integer64(1:10), nrow=2L, ncol=5L, dimnames = list(LETTERS[1:2], letters[1:5]))
+  m32[[1, 4]] = "112"
+  m64[[1, 4]] = "112"
+  expect_identical(m64, m32)
+
+})