diff --git a/CHANGELOG.md b/CHANGELOG.md
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -5,6 +5,23 @@
 The format is based on [Keep a Changelog](http://keepachangelog.com/en/1.0.0/)
 and this project adheres to [Semantic Versioning](http://semver.org/spec/v2.0.0.html).
 
+## [0.25.19]
+
+### Added
+
+* The compiler now does slightly less aggressive inlining. Use the
+  `#[inline]` attribute if you want to force inlining of some
+  function.
+
+* Arrays of opaque types now support indexing through the C API.
+  Arrays of records can also be constructed. (#2082)
+
+### Fixed
+
+* The `opencl` backend now always passes
+  `-cl-fp32-correctly-rounded-divide-sqrt` to the kernel compiler, in
+  order to match CUDA and HIP behaviour.
+
 ## [0.25.18]
 
 ### Added
diff --git a/docs/c-api.rst b/docs/c-api.rst
--- a/docs/c-api.rst
+++ b/docs/c-api.rst
@@ -64,7 +64,7 @@
 changes to the configuration must be made *before* calling
 :c:func:`futhark_context_new`.  A configuration object must not be
 freed before any context objects for which it is used.  The same
-configuration may *not* be used for multiple concurrent contexts.
+configuration must *not* be used for multiple concurrent contexts.
 Configuration objects are cheap to create and destroy.
 
 .. c:struct:: futhark_context_config
@@ -162,7 +162,7 @@
    :c:func:`futhark_context_get_error`, which will return non-``NULL``
    if initialisation failed.  If initialisation has failed, then you
    still need to call :c:func:`futhark_context_free` to release
-   resources used for the context object, but you may not use the
+   resources used for the context object, but you must not use the
    context object for anything else.
 
 .. c:function:: void futhark_context_free(struct futhark_context *ctx)
@@ -230,6 +230,11 @@
 contains the bitwise representation of the ``f16`` value in the IEEE
 754 binary16 format.
 
+.. _array-values:
+
+Arrays of Primitive Values
+~~~~~~~~~~~~~~~~~~~~~~~~~~
+
 For each distinct array type of primitives (ignoring sizes), an opaque
 C struct is defined.  Arrays of ``f16`` are presented as containing
 ``uint16_t`` elements.  For types that do not map cleanly to C,
@@ -272,7 +277,7 @@
 .. c:function:: int futhark_free_i32_1d(struct futhark_context *ctx, struct futhark_i32_1d *arr)
 
    Free the value.  In practice, this merely decrements the reference
-   count by one.  The value (or at least this reference) may not be
+   count by one.  The value (or at least this reference) must not be
    used again after this function returns.
 
 .. c:function:: int futhark_values_i32_1d(struct futhark_context *ctx, struct futhark_i32_1d *arr, int32_t *data)
@@ -282,6 +287,14 @@
    space to store the full array.  Multi-dimensional arrays are
    written in row-major form.
 
+.. c:function:: int futhark_index_i32_1d(struct futhark_context *ctx, int32_t *out, struct futhark_i32_1d *arr, int64_t i0);
+
+   Asynchronously copy a single element from the array and store it in
+   ``*out``. Returns a nonzero value if the index is out of bounds.
+   **Note:** if you need to read many elements, it is much faster to
+   retrieve the entire array with the ``values`` function,
+   particularly when using a GPU backend.
+
 .. c:function:: const int64_t *futhark_shape_i32_1d(struct futhark_context *ctx, struct futhark_i32_1d *arr)
 
    Return a pointer to the shape of the array, with one element per
@@ -326,7 +339,7 @@
 .. c:function:: int futhark_free_opaque_foo(struct futhark_context *ctx, struct futhark_opaque_foo *obj)
 
    Free the value.  In practice, this merely decrements the reference
-   count by one.  The value (or at least this reference) may not be
+   count by one.  The value (or at least this reference) must not be
    used again after this function returns.
 
 .. c:function:: int futhark_store_opaque_foo(struct futhark_context *ctx, const struct futhark_opaque_foo *obj, void **p, size_t *n)
@@ -474,7 +487,56 @@
    **Precondition:** ``t`` must be an instance of the ``foo`` variant,
    which can be determined with :c:func:`futhark_variant_opaque_t`.
 
+.. _arrays_of_opaques:
 
+Arrays of Non-Primitive Values
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+An array that contains a non-primitive type is considered an opaque
+value. However, it also supports a subset of the API documented in
+:ref:`array-values`.
+
+For an opaque array type ``[]t``, the following functions are always
+generated (assuming the generated C type is ``arr_t``):
+
+.. c:function:: int futhark_index_opaque_arr_t(struct futhark_context *ctx, struct futhark_opaque_t **out, struct futhark_opaque_arr_t *arr, int64_t i0);
+
+   Asynchronously copy a single element from the array and store it in
+   ``*out``. Returns a nonzero value if the index is out of bounds.
+
+.. c:function:: const int64_t *futhark_shape_opaque_arr_t(struct futhark_context *ctx, struct futhark_opaque_arr_t *arr);
+
+   Return a pointer to the shape of the array, with one element per
+   dimension. The lifetime of the shape is the same as ``arr``, and
+   must *not* be manually freed. Assuming ``arr`` is a valid object,
+   this function cannot fail.
+
+Additionally, if the element type is a record (or equivalently a
+tuple), for example if the array type is ``[](f32,f32)``, the
+following functions are also available:
+
+.. c:function:: int futhark_zip_opaque_arr1d_tup2_f32_f32(struct futhark_context *ctx, struct futhark_opaque_arr1d_tup2_f32_f32 **out, const struct futhark_f32_1d *f_0, const struct futhark_f32_1d *f_1);
+
+   Construct an array of records from arrays of the component values.
+   This is analogous to ``zip`` in the source language. The provided
+   arrays must have compatible shapes, and the function returns
+   nonzero if they do not.
+
+   **Note:** This is a cheap operation, as it does not copy array
+   elements.
+
+   **Note:** The resulting array aliases the original arrays.
+
+.. c:function:: int futhark_project_opaque_arr1d_tup2_f32_f32_0(struct futhark_context *ctx, struct futhark_f32_1d **out, const struct futhark_opaque_arr1d_tup2_f32_f32 *obj);
+
+   Retrieve an array of all the ``.0`` fields of the array elements. A
+   similar function is provided for each field.
+
+   **Note:** This is a cheap operation, as it does not copy array
+   elements.
+
+   **Note:** The resulting array aliases the original array.
+
 Entry points
 ------------
 
@@ -738,7 +800,7 @@
     functions are as documented above.  The following operations are
     listed:
 
-    * For arrays: ``free``, ``shape``, ``values``, ``new``.
+    * For arrays: ``free``, ``shape``, ``values``, ``new``, ``index``.
 
     * For opaques: ``free``, ``store``, ``restore``.
 
@@ -750,6 +812,23 @@
 
     * The name of the C *new* function for creating a record from
       field values.
+
+   * For opaques that are actually arrays of records:
+
+     * The element type and rank.
+
+     * The operations ``index``, ``shape``, ``zip``.
+
+     * The fields, which will be the fields of the element type, but
+       with the dimensions preprended. These are the types of the
+       arrays that should be passed to the ``zip`` function.
+
+   * For other opaques that are actually arrays:
+
+     * The element type and rank.
+
+     * The operations ``index`` and ``shape``.
+
 
 Manifests are defined by the following JSON Schema:
 
diff --git a/futhark.cabal b/futhark.cabal
--- a/futhark.cabal
+++ b/futhark.cabal
@@ -1,6 +1,6 @@
 cabal-version: 2.4
 name:           futhark
-version:        0.25.18
+version:        0.25.19
 synopsis:       An optimising compiler for a functional, array-oriented language.
 
 description:    Futhark is a small programming language designed to be compiled to
@@ -459,7 +459,7 @@
     , free >=5.1.10
     , futhark-data >= 1.1.0.0
     , futhark-server >= 1.2.2.1
-    , futhark-manifest >= 1.4.0.0
+    , futhark-manifest >= 1.5.0.0
     , githash >=0.1.6.1
     , half >= 0.3
     , haskeline
diff --git a/rts/c/backends/opencl.h b/rts/c/backends/opencl.h
--- a/rts/c/backends/opencl.h
+++ b/rts/c/backends/opencl.h
@@ -657,6 +657,12 @@
     w += snprintf(compile_opts+w, compile_opts_size-w, "-DEMULATE_F16 ");
   }
 
+  // By default, OpenCL allows imprecise (but faster) division and
+  // square root operations. For equivalence with other backends, ask
+  // for correctly rounded ones here.
+  w += snprintf(compile_opts+w, compile_opts_size-w,
+                "-cl-fp32-correctly-rounded-divide-sqrt");
+
   free(macro_names);
   free(macro_vals);
 
diff --git a/src/Futhark/CodeGen/Backends/GenericC/CLI.hs b/src/Futhark/CodeGen/Backends/GenericC/CLI.hs
--- a/src/Futhark/CodeGen/Backends/GenericC/CLI.hs
+++ b/src/Futhark/CodeGen/Backends/GenericC/CLI.hs
@@ -190,7 +190,7 @@
             [C.cstm|;|],
             [C.cexp|$id:dest|]
           )
-    Just (TypeOpaque desc _ _ _) ->
+    Just (TypeOpaque desc _ _) ->
       ( [C.citems|futhark_panic(1, "Cannot read input #%d of type %s\n", $int:i, $string:(T.unpack desc));|],
         [C.cstm|;|],
         [C.cstm|;|],
@@ -256,7 +256,7 @@
             [C.cexp|$id:result|],
             [C.cstm|assert($id:(arrayFree ops)(ctx, $id:result) == 0);|]
           )
-        Just (TypeOpaque t ops _ _) ->
+        Just (TypeOpaque t ops _) ->
           ( [C.citem|typename $id:t $id:result;|],
             [C.cexp|$id:result|],
             [C.cstm|assert($id:(opaqueFree ops)(ctx, $id:result) == 0);|]
@@ -269,7 +269,7 @@
     Nothing ->
       let info = tname <> "_info"
        in [C.cstm|write_scalar(stdout, binary_output, &$id:info, &$exp:e);|]
-    Just (TypeOpaque desc _ _ _) ->
+    Just (TypeOpaque desc _ _) ->
       [C.cstm|{
          fprintf(stderr, "Values of type \"%s\" have no external representation.\n", $string:(T.unpack desc));
          retval = 1;
diff --git a/src/Futhark/CodeGen/Backends/GenericC/EntryPoints.hs b/src/Futhark/CodeGen/Backends/GenericC/EntryPoints.hs
--- a/src/Futhark/CodeGen/Backends/GenericC/EntryPoints.hs
+++ b/src/Futhark/CodeGen/Backends/GenericC/EntryPoints.hs
@@ -23,11 +23,6 @@
 valueDescToType (ArrayValue _ _ pt signed shape) =
   ValueType signed (Rank (length shape)) pt
 
-allTrue :: [C.Exp] -> C.Exp
-allTrue [] = [C.cexp|true|]
-allTrue [x] = x
-allTrue (x : xs) = [C.cexp|$exp:x && $exp:(allTrue xs)|]
-
 prepareEntryInputs ::
   [ExternalValue] ->
   CompilerM op s ([(C.Param, Maybe C.Exp)], [C.BlockItem])
diff --git a/src/Futhark/CodeGen/Backends/GenericC/Monad.hs b/src/Futhark/CodeGen/Backends/GenericC/Monad.hs
--- a/src/Futhark/CodeGen/Backends/GenericC/Monad.hs
+++ b/src/Futhark/CodeGen/Backends/GenericC/Monad.hs
@@ -282,6 +282,7 @@
   mapM_
     earlyDecl
     [C.cunit|struct program {
+               int dummy;
                $sdecls:fields
              };
              static void setup_program(struct futhark_context* ctx) {
diff --git a/src/Futhark/CodeGen/Backends/GenericC/Server.hs b/src/Futhark/CodeGen/Backends/GenericC/Server.hs
--- a/src/Futhark/CodeGen/Backends/GenericC/Server.hs
+++ b/src/Futhark/CodeGen/Backends/GenericC/Server.hs
@@ -119,7 +119,7 @@
 cType manifest tname =
   case M.lookup tname $ manifestTypes manifest of
     Just (TypeArray ctype _ _ _) -> [C.cty|typename $id:(T.unpack ctype)|]
-    Just (TypeOpaque ctype _ _ _) -> [C.cty|typename $id:(T.unpack ctype)|]
+    Just (TypeOpaque ctype _ _) -> [C.cty|typename $id:(T.unpack ctype)|]
     Nothing -> uncurry primAPIType $ scalarToPrim tname
 
 -- First component is forward declaration so we don't have to worry
@@ -156,10 +156,10 @@
                 .aux = &$id:aux_name
               };|]
       )
-typeBoilerplate manifest (tname, TypeOpaque c_type_name ops record _sumops) =
+typeBoilerplate manifest (tname, TypeOpaque c_type_name ops extra_ops) =
   let type_name = typeStructName tname
       aux_name = type_name <> "_aux"
-      (record_edecls, record_init) = recordDefs type_name record
+      (record_edecls, record_init) = recordDefs type_name extra_ops
    in ( [C.cedecl|const struct type $id:type_name;|],
         [C.cinit|&$id:type_name|],
         record_edecls
@@ -179,8 +179,7 @@
               };|]
       )
   where
-    recordDefs _ Nothing = ([], [C.cinit|NULL|])
-    recordDefs type_name (Just (RecordOps fields new)) =
+    recordDefs type_name (Just (OpaqueRecord (RecordOps fields new))) =
       let new_wrap = new <> "_wrap"
           record_name = type_name <> "_record"
           fields_name = type_name <> "_fields"
@@ -212,6 +211,7 @@
              };|],
             [C.cinit|&$id:record_name|]
           )
+    recordDefs _ _ = ([], [C.cinit|NULL|])
 
 entryTypeBoilerplate :: Manifest -> ([C.Definition], [C.Initializer], [C.Definition])
 entryTypeBoilerplate manifest =
diff --git a/src/Futhark/CodeGen/Backends/GenericC/Types.hs b/src/Futhark/CodeGen/Backends/GenericC/Types.hs
--- a/src/Futhark/CodeGen/Backends/GenericC/Types.hs
+++ b/src/Futhark/CodeGen/Backends/GenericC/Types.hs
@@ -11,7 +11,7 @@
 import Control.Monad
 import Control.Monad.Reader (asks)
 import Control.Monad.State (gets, modify)
-import Data.Char (isDigit)
+import Data.List qualified as L
 import Data.Map.Strict qualified as M
 import Data.Maybe
 import Data.Text qualified as T
@@ -37,6 +37,25 @@
   modify $ \s -> s {compArrayTypes = add $ compArrayTypes s}
   pure [C.cty|struct $id:name|]
 
+prepareNewMem ::
+  (C.ToExp arr, C.ToExp dim) =>
+  arr ->
+  Space ->
+  [dim] ->
+  PrimType ->
+  CompilerM op s ()
+prepareNewMem arr space shape pt = do
+  let rank = length shape
+      arr_size = cproduct [[C.cexp|$exp:k|] | k <- shape]
+  resetMem [C.cexp|$exp:arr->mem|] space
+  allocMem
+    [C.cexp|$exp:arr->mem|]
+    [C.cexp|$exp:arr_size * $int:(primByteSize pt::Int)|]
+    space
+    [C.cstm|err = 1;|]
+  forM_ (zip [0 .. rank - 1] shape) $ \(i, dim_s) ->
+    stm [C.cstm|$exp:arr->shape[$int:i] = $exp:dim_s;|]
+
 arrayLibraryFunctions ::
   Publicness ->
   Space ->
@@ -56,28 +75,22 @@
   values_array <- publicName $ "values_" <> name
   values_raw_array <- publicName $ "values_raw_" <> name
   shape_array <- publicName $ "shape_" <> name
+  index_array <- publicName $ "index_" <> name
 
   let shape_names = ["dim" <> prettyText i | i <- [0 .. rank - 1]]
       shape_params = [[C.cparam|typename int64_t $id:k|] | k <- shape_names]
-      arr_size = cproduct [[C.cexp|$id:k|] | k <- shape_names]
+      shape = [[C.cexp|$id:k|] | k <- shape_names]
+      index_names = ["i" <> prettyText i | i <- [0 .. rank - 1]]
+      index_params = [[C.cparam|typename int64_t $id:k|] | k <- index_names]
+      arr_size = cproduct shape
       arr_size_array = cproduct [[C.cexp|arr->shape[$int:i]|] | i <- [0 .. rank - 1]]
+
   copy <- asks $ opsCopy . envOperations
 
   memty <- rawMemCType space
 
-  let prepare_new = do
-        resetMem [C.cexp|arr->mem|] space
-        allocMem
-          [C.cexp|arr->mem|]
-          [C.cexp|$exp:arr_size * $int:(primByteSize pt::Int)|]
-          space
-          [C.cstm|return NULL;|]
-        forM_ [0 .. rank - 1] $ \i ->
-          let dim_s = "dim" ++ show i
-           in stm [C.cstm|arr->shape[$int:i] = $id:dim_s;|]
-
   new_body <- collect $ do
-    prepare_new
+    prepareNewMem [C.cexp|arr|] space shape pt
     copy
       CopyNoBarrier
       [C.cexp|arr->mem.mem|]
@@ -109,6 +122,31 @@
         space
         [C.cexp|((size_t)$exp:arr_size_array) * $int:(primByteSize pt::Int)|]
 
+  let arr_strides = do
+        r <- [0 .. rank - 1]
+        pure $ cproduct [[C.cexp|arr->shape[$int:i]|] | i <- [r + 1 .. rank - 1]]
+      index_exp =
+        cproduct
+          [ [C.cexp|$int:(primByteSize pt::Int)|],
+            csum (zipWith (\x y -> [C.cexp|$id:x * $exp:y|]) index_names arr_strides)
+          ]
+      in_bounds =
+        allTrue
+          [ [C.cexp|$id:p >= 0 && $id:p < arr->shape[$int:i]|]
+            | (p, i) <- zip index_names [0 .. rank - 1]
+          ]
+  index_body <-
+    collect $
+      copy
+        CopyNoBarrier
+        [C.cexp|(unsigned char*)out|]
+        [C.cexp|0|]
+        DefaultSpace
+        [C.cexp|arr->mem.mem|]
+        index_exp
+        space
+        [C.cexp|$int:(primByteSize pt::Int)|]
+
   ctx_ty <- contextType
   ops <- asks envOperations
 
@@ -127,6 +165,9 @@
   proto
     [C.cedecl|int $id:values_array($ty:ctx_ty *ctx, $ty:array_type *arr, $ty:pt' *data);|]
   proto
+    [C.cedecl|int $id:index_array($ty:ctx_ty *ctx, $ty:pt' *out, $ty:array_type *arr,
+                                  $params:index_params);|]
+  proto
     [C.cedecl|$ty:memty $id:values_raw_array($ty:ctx_ty *ctx, $ty:array_type *arr);|]
   proto
     [C.cedecl|const typename int64_t* $id:shape_array($ty:ctx_ty *ctx, $ty:array_type *arr);|]
@@ -172,6 +213,18 @@
             return err;
           }
 
+          int $id:index_array($ty:ctx_ty *ctx, $ty:pt' *out, $ty:array_type *arr,
+                              $params:index_params) {
+            int err = 0;
+            if ($exp:in_bounds) {
+              $items:(criticalSection ops index_body)
+            } else {
+              err = 1;
+              set_error(ctx, strdup("Index out of bounds."));
+            }
+            return err;
+          }
+
           $ty:memty $id:values_raw_array($ty:ctx_ty *ctx, $ty:array_type *arr) {
             (void)ctx;
             return arr->mem.mem;
@@ -190,7 +243,8 @@
         Manifest.arrayValues = values_array,
         Manifest.arrayNew = new_array,
         Manifest.arrayNewRaw = new_raw_array,
-        Manifest.arrayValuesRaw = values_raw_array
+        Manifest.arrayValuesRaw = values_raw_array,
+        Manifest.arrayIndex = index_array
       }
 
 lookupOpaqueType :: Name -> OpaqueTypes -> OpaqueType
@@ -202,10 +256,15 @@
 opaquePayload :: OpaqueTypes -> OpaqueType -> [ValueType]
 opaquePayload _ (OpaqueType ts) = ts
 opaquePayload _ (OpaqueSum ts _) = ts
+opaquePayload _ (OpaqueArray _ _ ts) = ts
 opaquePayload types (OpaqueRecord fs) = concatMap f fs
   where
     f (_, TypeOpaque s) = opaquePayload types $ lookupOpaqueType s types
     f (_, TypeTransparent v) = [v]
+opaquePayload types (OpaqueRecordArray _ _ fs) = concatMap f fs
+  where
+    f (_, TypeOpaque s) = opaquePayload types $ lookupOpaqueType s types
+    f (_, TypeTransparent v) = [v]
 
 entryPointTypeToCType :: Publicness -> EntryPointType -> CompilerM op s C.Type
 entryPointTypeToCType _ (TypeOpaque desc) = opaqueToCType desc
@@ -305,48 +364,26 @@
                *v->$id:(tupleField i) = *$exp:e;
                (void)(*(v->$id:(tupleField i)->mem.references))++;}|]
 
-recordNewFunctions ::
+recordNewSetFields ::
   OpaqueTypes ->
-  Name ->
   [(Name, EntryPointType)] ->
   [ValueType] ->
-  CompilerM op s Manifest.CFuncName
-recordNewFunctions types desc fs vds = do
-  opaque_type <- opaqueToCType desc
-  ctx_ty <- contextType
-  ops <- asks envOperations
-  new <- publicName $ "new_" <> opaqueName desc
-
-  (params, new_stms) <-
-    fmap (unzip . snd)
-      . mapAccumLM onField 0
-      . zip fs
-      . recordFieldPayloads types (map snd fs)
-      $ vds
-
-  headerDecl
-    (OpaqueDecl desc)
-    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params);|]
-  libDecl
-    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params) {
-                $ty:opaque_type* v = malloc(sizeof($ty:opaque_type));
-                $items:(criticalSection ops new_stms)
-                *out = v;
-                return FUTHARK_SUCCESS;
-              }|]
-  pure new
+  CompilerM op s ([C.Id], [C.Param], [C.BlockItem])
+recordNewSetFields types fs =
+  fmap (L.unzip3 . snd)
+    . mapAccumLM onField 0
+    . zip fs
+    . recordFieldPayloads types (map snd fs)
   where
     onField offset ((f, et), f_vts) = do
-      let param_name =
-            if T.all isDigit (nameToText f)
-              then C.toIdent ("v" <> f) mempty
-              else C.toIdent ("f_" <> f) mempty
+      let param_name = C.toIdent ("f_" <> f) mempty
       case et of
         TypeTransparent (ValueType sign (Rank 0) pt) -> do
           let ct = primAPIType sign pt
           pure
             ( offset + 1,
-              ( [C.cparam|const $ty:ct $id:param_name|],
+              ( param_name,
+                [C.cparam|const $ty:ct $id:param_name|],
                 [C.citem|v->$id:(tupleField offset) = $id:param_name;|]
               )
             )
@@ -354,7 +391,8 @@
           ct <- valueTypeToCType Public vt
           pure
             ( offset + 1,
-              ( [C.cparam|const $ty:ct* $id:param_name|],
+              ( param_name,
+                [C.cparam|const $ty:ct* $id:param_name|],
                 [C.citem|{v->$id:(tupleField offset) = malloc(sizeof($ty:ct));
                           *v->$id:(tupleField offset) = *$id:param_name;
                           (void)(*(v->$id:(tupleField offset)->mem.references))++;}|]
@@ -367,11 +405,211 @@
                 pure [C.cexp|$id:param_name->$id:(tupleField i)|]
           pure
             ( offset + length f_vts,
-              ( [C.cparam|const $ty:ct* $id:param_name|],
+              ( param_name,
+                [C.cparam|const $ty:ct* $id:param_name|],
                 [C.citem|{$stms:(zipWith3 setFieldField [offset ..] param_fields f_vts)}|]
               )
             )
 
+recordNewFunctions ::
+  OpaqueTypes ->
+  Name ->
+  [(Name, EntryPointType)] ->
+  [ValueType] ->
+  CompilerM op s Manifest.CFuncName
+recordNewFunctions types desc fs vds = do
+  opaque_type <- opaqueToCType desc
+  ctx_ty <- contextType
+  ops <- asks envOperations
+  new <- publicName $ "new_" <> opaqueName desc
+
+  (_, params, new_stms) <- recordNewSetFields types fs vds
+
+  headerDecl
+    (OpaqueDecl desc)
+    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params);|]
+  libDecl
+    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params) {
+                $ty:opaque_type* v = malloc(sizeof($ty:opaque_type));
+                $items:(criticalSection ops new_stms)
+                *out = v;
+                return FUTHARK_SUCCESS;
+              }|]
+  pure new
+
+-- Because records and arrays-of-records are very similar in their
+-- actual representation, we can reuse most of the code. Only indexing
+-- requires something special.
+
+recordArrayProjectFunctions ::
+  OpaqueTypes ->
+  Name ->
+  [(Name, EntryPointType)] ->
+  [ValueType] ->
+  CompilerM op s [Manifest.RecordField]
+recordArrayProjectFunctions = recordProjectFunctions
+
+recordArrayZipFunctions ::
+  OpaqueTypes ->
+  Name ->
+  [(Name, EntryPointType)] ->
+  [ValueType] ->
+  Int ->
+  CompilerM op s Manifest.CFuncName
+recordArrayZipFunctions types desc fs vds rank = do
+  opaque_type <- opaqueToCType desc
+  ctx_ty <- contextType
+  ops <- asks envOperations
+  new <- publicName $ "zip_" <> opaqueName desc
+
+  (param_names, params, new_stms) <- recordNewSetFields types fs vds
+
+  headerDecl
+    (OpaqueDecl desc)
+    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params);|]
+  libDecl
+    [C.cedecl|int $id:new($ty:ctx_ty *ctx, $ty:opaque_type** out, $params:params) {
+                if (!$exp:(sameShape param_names)) {
+                  set_error(ctx, strdup("Cannot zip arrays with different shapes."));
+                  return 1;
+                }
+                $ty:opaque_type* v = malloc(sizeof($ty:opaque_type));
+                $items:(criticalSection ops new_stms)
+                *out = v;
+                return FUTHARK_SUCCESS;
+              }|]
+  pure new
+  where
+    valueShape TypeTransparent {} p =
+      [[C.cexp|$id:p->shape[$int:i]|] | i <- [0 .. rank - 1]]
+    -- We know that the opaque value must contain arrays.
+    valueShape TypeOpaque {} p =
+      [[C.cexp|$id:p->$id:(tupleField 0)->shape[$int:i]|] | i <- [0 .. rank - 1]]
+    sameShape param_names =
+      allTrue $ map allEqual $ L.transpose $ zipWith valueShape (map snd fs) param_names
+
+recordArrayIndexFunctions ::
+  Space ->
+  OpaqueTypes ->
+  Name ->
+  Int ->
+  Name ->
+  [ValueType] ->
+  CompilerM op s Manifest.CFuncName
+recordArrayIndexFunctions space _types desc rank elemtype vds = do
+  index_f <- publicName $ "index_" <> opaqueName desc
+  ctx_ty <- contextType
+  array_ct <- opaqueToCType desc
+  obj_ct <- opaqueToCType elemtype
+  copy <- asks $ opsCopy . envOperations
+
+  index_items <- collect $ zipWithM_ (setField copy) [0 ..] vds
+
+  headerDecl
+    (OpaqueDecl desc)
+    [C.cedecl|int $id:index_f($ty:ctx_ty *ctx, $ty:obj_ct **out, $ty:array_ct *arr,
+                              $params:index_params);|]
+  libDecl
+    [C.cedecl|int $id:index_f($ty:ctx_ty *ctx, $ty:obj_ct **out, $ty:array_ct *arr,
+                              $params:index_params) {
+                int err = 0;
+                if ($exp:in_bounds) {
+                  $ty:obj_ct* v = malloc(sizeof($ty:obj_ct));
+                  $items:index_items
+                  if (err == 0) {
+                    *out = v;
+                  }
+                } else {
+                  err = 1;
+                  set_error(ctx, strdup("Index out of bounds."));
+                }
+                return err;
+              }|]
+
+  pure index_f
+  where
+    index_names = ["i" <> prettyText i | i <- [0 .. rank - 1]]
+    index_params = [[C.cparam|typename int64_t $id:k|] | k <- index_names]
+    indexExp pt r shape =
+      cproduct
+        [ [C.cexp|$int:(primByteSize pt::Int)|],
+          csum (zipWith (\x y -> [C.cexp|$id:x * $exp:y|]) index_names strides)
+        ]
+      where
+        strides = do
+          d <- [0 .. r - 1]
+          pure $ cproduct [[C.cexp|$exp:shape[$int:i]|] | i <- [d + 1 .. r - 1]]
+
+    in_bounds =
+      allTrue
+        [ [C.cexp|$id:p >= 0 && $id:p < arr->$id:(tupleField 0)->shape[$int:i]|]
+          | (p, i) <- zip index_names [0 .. rank - 1]
+        ]
+
+    setField copy j (ValueType _ (Rank r) pt)
+      | r == rank =
+          -- Easy case: just copy the scalar from the array into the
+          -- variable.
+          copy
+            CopyNoBarrier
+            [C.cexp|&v->$id:(tupleField j)|]
+            [C.cexp|0|]
+            DefaultSpace
+            [C.cexp|arr->$id:(tupleField j)->mem.mem|]
+            (indexExp pt rank [C.cexp|arr->$id:(tupleField j)->shape|])
+            space
+            [C.cexp|$int:(primByteSize pt::Int)|]
+      | otherwise = do
+          -- Tricky case, where we first have to allocate memory.
+          let shape = do
+                i <- [rank .. r - 1]
+                pure [C.cexp|arr->$id:(tupleField j)->shape[$int:i]|]
+          stm [C.cstm|v->$id:(tupleField j) = malloc(sizeof(*v->$id:(tupleField j)));|]
+          prepareNewMem [C.cexp|v->$id:(tupleField j)|] space shape pt
+          -- Now we can copy into the freshly allocated memory.
+          copy
+            CopyNoBarrier
+            [C.cexp|v->$id:(tupleField j)->mem.mem|]
+            [C.cexp|0|]
+            DefaultSpace
+            [C.cexp|arr->$id:(tupleField j)->mem.mem|]
+            (indexExp pt r [C.cexp|arr->$id:(tupleField j)->shape|])
+            space
+            $ cproduct ([C.cexp|$int:(primByteSize pt::Int)|] : shape)
+
+recordArrayShapeFunctions :: Name -> CompilerM op s Manifest.CFuncName
+recordArrayShapeFunctions desc = do
+  shape_f <- publicName $ "shape_" <> opaqueName desc
+  ctx_ty <- contextType
+  array_ct <- opaqueToCType desc
+
+  -- We know that the opaque value consists of arrays of at least the
+  -- expected rank, and which have the same outer shape, so we just
+  -- return the shape of the first one.
+  headerDecl
+    (OpaqueDecl desc)
+    [C.cedecl|const typename int64_t* $id:shape_f($ty:ctx_ty *ctx, $ty:array_ct *arr);|]
+  libDecl
+    [C.cedecl|const typename int64_t* $id:shape_f($ty:ctx_ty *ctx, $ty:array_ct *arr) {
+                (void)ctx;
+                return arr->$id:(tupleField 0)->shape;
+              }|]
+
+  pure shape_f
+
+opaqueArrayIndexFunctions ::
+  Space ->
+  OpaqueTypes ->
+  Name ->
+  Int ->
+  Name ->
+  [ValueType] ->
+  CompilerM op s Manifest.CFuncName
+opaqueArrayIndexFunctions = recordArrayIndexFunctions
+
+opaqueArrayShapeFunctions :: Name -> CompilerM op s Manifest.CFuncName
+opaqueArrayShapeFunctions = recordArrayShapeFunctions
+
 sumVariants ::
   Name ->
   [(Name, [(EntryPointType, [Int])])] ->
@@ -513,33 +751,49 @@
               }|]
   pure variant
 
-processOpaqueRecord ::
+opaqueExtraOps ::
+  Space ->
   OpaqueTypes ->
   Name ->
   OpaqueType ->
   [ValueType] ->
-  CompilerM op s (Maybe Manifest.RecordOps, Maybe Manifest.SumOps)
-processOpaqueRecord _ _ (OpaqueType _) _ =
-  pure (Nothing, Nothing)
-processOpaqueRecord _types desc (OpaqueSum _ cs) vds =
-  (Nothing,) . Just
+  CompilerM op s (Maybe Manifest.OpaqueExtraOps)
+opaqueExtraOps _ _ _ (OpaqueType _) _ =
+  pure Nothing
+opaqueExtraOps _ _types desc (OpaqueSum _ cs) vds =
+  Just . Manifest.OpaqueSum
     <$> ( Manifest.SumOps
             <$> sumVariants desc cs vds
             <*> sumVariantFunction desc
         )
-processOpaqueRecord types desc (OpaqueRecord fs) vds =
-  (,Nothing) . Just
+opaqueExtraOps _ types desc (OpaqueRecord fs) vds =
+  Just . Manifest.OpaqueRecord
     <$> ( Manifest.RecordOps
             <$> recordProjectFunctions types desc fs vds
             <*> recordNewFunctions types desc fs vds
         )
+opaqueExtraOps space types desc (OpaqueRecordArray rank elemtype fs) vds =
+  Just . Manifest.OpaqueRecordArray
+    <$> ( Manifest.RecordArrayOps rank (nameToText elemtype)
+            <$> recordArrayProjectFunctions types desc fs vds
+            <*> recordArrayZipFunctions types desc fs vds rank
+            <*> recordArrayIndexFunctions space types desc rank elemtype vds
+            <*> recordArrayShapeFunctions desc
+        )
+opaqueExtraOps space types desc (OpaqueArray rank elemtype _) vds =
+  Just . Manifest.OpaqueArray
+    <$> ( Manifest.OpaqueArrayOps rank (nameToText elemtype)
+            <$> opaqueArrayIndexFunctions space types desc rank elemtype vds
+            <*> opaqueArrayShapeFunctions desc
+        )
 
 opaqueLibraryFunctions ::
+  Space ->
   OpaqueTypes ->
   Name ->
   OpaqueType ->
-  CompilerM op s (Manifest.OpaqueOps, Maybe Manifest.RecordOps, Maybe Manifest.SumOps)
-opaqueLibraryFunctions types desc ot = do
+  CompilerM op s (Manifest.OpaqueOps, Maybe Manifest.OpaqueExtraOps)
+opaqueLibraryFunctions space types desc ot = do
   name <- publicName $ opaqueName desc
   free_opaque <- publicName $ "free_" <> opaqueName desc
   store_opaque <- publicName $ "store_" <> opaqueName desc
@@ -637,7 +891,7 @@
     (OpaqueDecl desc)
     [C.cedecl|$ty:opaque_type* $id:restore_opaque($ty:ctx_ty *ctx, const void *p);|]
 
-  (record, sumops) <- processOpaqueRecord types desc ot vds
+  extra_ops <- opaqueExtraOps space types desc ot vds
 
   -- We do not need to enclose most bodies in a critical section,
   -- because when we operate on the components of the opaque, we are
@@ -685,8 +939,7 @@
           Manifest.opaqueStore = store_opaque,
           Manifest.opaqueRestore = restore_opaque
         },
-      record,
-      sumops
+      extra_ops
     )
 
 generateArray ::
@@ -712,16 +965,20 @@
       pure Nothing
 
 generateOpaque ::
+  Space ->
   OpaqueTypes ->
   (Name, OpaqueType) ->
   CompilerM op s (T.Text, Manifest.Type)
-generateOpaque types (desc, ot) = do
+generateOpaque space types (desc, ot) = do
   name <- publicName $ opaqueName desc
   members <- zipWithM field (opaquePayload types ot) [(0 :: Int) ..]
   libDecl [C.cedecl|struct $id:name { $sdecls:members };|]
-  (ops, record, sumops) <- opaqueLibraryFunctions types desc ot
+  (ops, extra_ops) <- opaqueLibraryFunctions space types desc ot
   let opaque_type = [C.cty|struct $id:name*|]
-  pure (nameToText desc, Manifest.TypeOpaque (typeText opaque_type) ops record sumops)
+  pure
+    ( nameToText desc,
+      Manifest.TypeOpaque (typeText opaque_type) ops extra_ops
+    )
   where
     field vt@(ValueType _ (Rank r) _) i = do
       ct <- valueTypeToCType Private vt
@@ -734,7 +991,7 @@
 generateAPITypes arr_space types@(OpaqueTypes opaques) = do
   mapM_ (findNecessaryArrays . snd) opaques
   array_ts <- mapM (generateArray arr_space) . M.toList =<< gets compArrayTypes
-  opaque_ts <- mapM (generateOpaque types) opaques
+  opaque_ts <- mapM (generateOpaque arr_space types) opaques
   pure $ M.fromList $ catMaybes array_ts <> opaque_ts
   where
     -- Ensure that array types will be generated before the opaque
@@ -743,6 +1000,10 @@
     -- the innards to increment reference counts.
     findNecessaryArrays (OpaqueType _) =
       pure ()
+    findNecessaryArrays (OpaqueArray {}) =
+      pure ()
+    findNecessaryArrays (OpaqueRecordArray _ _ fs) =
+      mapM_ (entryPointTypeToCType Public . snd) fs
     findNecessaryArrays (OpaqueSum _ variants) =
       mapM_ (mapM_ (entryPointTypeToCType Public . fst) . snd) variants
     findNecessaryArrays (OpaqueRecord fs) =
diff --git a/src/Futhark/CodeGen/Backends/SimpleRep.hs b/src/Futhark/CodeGen/Backends/SimpleRep.hs
--- a/src/Futhark/CodeGen/Backends/SimpleRep.hs
+++ b/src/Futhark/CodeGen/Backends/SimpleRep.hs
@@ -24,6 +24,8 @@
     fromStorage,
     cproduct,
     csum,
+    allEqual,
+    allTrue,
     scalarToPrim,
 
     -- * Primitive value operations
@@ -36,13 +38,17 @@
   )
 where
 
+import Control.Monad (void)
 import Data.Char (isAlpha, isAlphaNum, isDigit)
 import Data.Text qualified as T
+import Data.Void (Void)
 import Futhark.CodeGen.ImpCode
 import Futhark.CodeGen.RTS.C (scalarF16H, scalarH)
 import Futhark.Util (hashText, showText, zEncodeText)
 import Language.C.Quote.C qualified as C
 import Language.C.Syntax qualified as C
+import Text.Megaparsec
+import Text.Megaparsec.Char (space)
 
 -- | The C type corresponding to a signed integer type.
 intTypeToCType :: IntType -> C.Type
@@ -132,28 +138,37 @@
   where
     ok c = isAlphaNum c || c == '_'
 
-isArrayName :: T.Text -> (Int, T.Text)
-isArrayName s =
-  if "[]" `T.isPrefixOf` s
-    then
-      let (k, s') = isArrayName (T.drop 2 s)
-       in (k + 1, s')
-    else (0, s)
+-- | Find a nice C type name name for the Futhark type. This solely
+-- serves to make the generated header file easy to read, and we can
+-- always fall back on an ugly hash.
+findPrettyName :: T.Text -> Either String T.Text
+findPrettyName =
+  either (Left . errorBundlePretty) Right . parse (p <* eof) "type name"
+  where
+    p :: Parsec Void T.Text T.Text
+    p = choice [pArr, pTup, pAtom]
+    pArr = do
+      dims <- some "[]"
+      (("arr" <> showText (length dims) <> "d_") <>) <$> p
+    pTup = between "(" ")" $ do
+      ts <- p `sepBy` pComma
+      pure $ "tup" <> showText (length ts) <> "_" <> T.intercalate "_" ts
+    pAtom = T.pack <$> some (satisfy (`notElem` ("[]{}()," :: String)))
+    pComma = void $ "," <* space
 
 -- | The name of exposed opaque types.
 opaqueName :: Name -> T.Text
 opaqueName "()" = "opaque_unit" -- Hopefully this ad-hoc convenience won't bite us.
 opaqueName s
-  | (k, s'') <- isArrayName s',
-    k > 0,
-    valid s'' =
-      "opaque_arr_" <> s'' <> "_" <> showText k <> "d"
+  | Right v <- findPrettyName s',
+    valid v =
+      "opaque_" <> v
   | valid s' = "opaque_" <> s'
   where
     s' = nameToText s
 opaqueName s = "opaque_" <> hashText (nameToText s)
 
--- | The 'PrimType' (and sign) correspond to a human-readable scalar
+-- | The 'PrimType' (and sign) corresponding to a human-readable scalar
 -- type name (e.g. @f64@).  Beware: partial!
 scalarToPrim :: T.Text -> (Signedness, PrimType)
 scalarToPrim "bool" = (Signed, Bool)
@@ -185,6 +200,19 @@
 csum (e : es) = foldl mult e es
   where
     mult x y = [C.cexp|$exp:x + $exp:y|]
+
+-- | An expression that is true if these are also all true.
+allTrue :: [C.Exp] -> C.Exp
+allTrue [] = [C.cexp|true|]
+allTrue [x] = x
+allTrue (x : xs) = [C.cexp|$exp:x && $exp:(allTrue xs)|]
+
+-- | An expression that is true if these expressions are all equal by
+-- @==@.
+allEqual :: [C.Exp] -> C.Exp
+allEqual [x, y] = [C.cexp|$exp:x == $exp:y|]
+allEqual (x : y : xs) = [C.cexp|$exp:x == $exp:y && $exp:(allEqual(y:xs))|]
+allEqual _ = [C.cexp|true|]
 
 instance C.ToIdent Name where
   toIdent = C.toIdent . zEncodeText . nameToText
diff --git a/src/Futhark/CodeGen/ImpGen.hs b/src/Futhark/CodeGen/ImpGen.hs
--- a/src/Futhark/CodeGen/ImpGen.hs
+++ b/src/Futhark/CodeGen/ImpGen.hs
@@ -487,6 +487,8 @@
 entryPointSignedness types (TypeOpaque desc) =
   case lookupOpaqueType desc types of
     OpaqueType vts -> map valueTypeSign vts
+    OpaqueArray _ _ vts -> map valueTypeSign vts
+    OpaqueRecordArray _ _ fs -> foldMap (entryPointSignedness types . snd) fs
     OpaqueRecord fs -> foldMap (entryPointSignedness types . snd) fs
     OpaqueSum vts _ -> map valueTypeSign vts
 
@@ -499,6 +501,8 @@
 entryPointSize types (TypeOpaque desc) =
   case lookupOpaqueType desc types of
     OpaqueType vts -> length vts
+    OpaqueArray _ _ vts -> length vts
+    OpaqueRecordArray _ _ fs -> sum $ map (entryPointSize types . snd) fs
     OpaqueRecord fs -> sum $ map (entryPointSize types . snd) fs
     OpaqueSum vts _ -> length vts
 
diff --git a/src/Futhark/IR/Parse.hs b/src/Futhark/IR/Parse.hs
--- a/src/Futhark/IR/Parse.hs
+++ b/src/Futhark/IR/Parse.hs
@@ -654,7 +654,7 @@
 pOpaqueType =
   (,)
     <$> (keyword "type" *> (nameFromText <$> pStringLiteral) <* pEqual)
-    <*> choice [pRecord, pSum, pOpaque]
+    <*> choice [pRecord, pSum, pOpaque, pRecordArray, pOpaqueArray]
   where
     pFieldName = choice [pName, nameFromString . show <$> pInt]
     pField = (,) <$> pFieldName <* pColon <*> pEntryPointType
@@ -676,6 +676,20 @@
           )
 
     pOpaque = keyword "opaque" $> OpaqueType <*> braces (many pValueType)
+
+    pRecordArray =
+      keyword "record_array"
+        $> OpaqueRecordArray
+        <*> (pInt <* lexeme "d")
+        <*> (nameFromText <$> pStringLiteral)
+        <*> braces (many pField)
+
+    pOpaqueArray =
+      keyword "array"
+        $> OpaqueArray
+        <*> (pInt <* lexeme "d")
+        <*> (nameFromText <$> pStringLiteral)
+        <*> braces (many pValueType)
 
 pOpaqueTypes :: Parser OpaqueTypes
 pOpaqueTypes = keyword "types" $> OpaqueTypes <*> braces (many pOpaqueType)
diff --git a/src/Futhark/IR/Pretty.hs b/src/Futhark/IR/Pretty.hs
--- a/src/Futhark/IR/Pretty.hs
+++ b/src/Futhark/IR/Pretty.hs
@@ -417,6 +417,15 @@
     "sum" <+> nestedBlock "{" "}" (stack $ pretty ts : map p cs)
     where
       p (c, ets) = hsep $ "#" <> pretty c : map pretty ets
+  pretty (OpaqueArray r v ts) =
+    "array" <+> pretty r
+      <> "d"
+        <+> dquotes (pretty v)
+        <+> nestedBlock "{" "}" (stack $ map pretty ts)
+  pretty (OpaqueRecordArray r v fs) =
+    "record_array" <+> pretty r <> "d" <+> dquotes (pretty v) <+> nestedBlock "{" "}" (stack $ map p fs)
+    where
+      p (f, et) = pretty f <> ":" <+> pretty et
 
 instance Pretty OpaqueTypes where
   pretty (OpaqueTypes ts) = "types" <+> nestedBlock "{" "}" (stack $ map p ts)
diff --git a/src/Futhark/IR/Syntax/Core.hs b/src/Futhark/IR/Syntax/Core.hs
--- a/src/Futhark/IR/Syntax/Core.hs
+++ b/src/Futhark/IR/Syntax/Core.hs
@@ -608,6 +608,11 @@
     -- represent that constructor payload. This is necessary because
     -- we deduplicate payloads across constructors.
     OpaqueSum [ValueType] [(Name, [(EntryPointType, [Int])])]
+  | -- | An array with this rank and named opaque element type.
+    OpaqueArray Int Name [ValueType]
+  | -- | An array with known rank and where the elements are this
+    -- record type.
+    OpaqueRecordArray Int Name [(Name, EntryPointType)]
   deriving (Eq, Ord, Show)
 
 -- | Names of opaque types and their representation.
diff --git a/src/Futhark/IR/TypeCheck.hs b/src/Futhark/IR/TypeCheck.hs
--- a/src/Futhark/IR/TypeCheck.hs
+++ b/src/Futhark/IR/TypeCheck.hs
@@ -560,6 +560,11 @@
       mapM_ (checkEntryPointType known . snd) fs
     check known (OpaqueSum _ cs) =
       mapM_ (mapM_ (checkEntryPointType known . fst) . snd) cs
+    check known (OpaqueArray _ v _) =
+      checkEntryPointType known (TypeOpaque v)
+    check known (OpaqueRecordArray _ v fs) = do
+      checkEntryPointType known (TypeOpaque v)
+      mapM_ (checkEntryPointType known . snd) fs
     check _ (OpaqueType _) =
       pure ()
     checkEntryPointType known (TypeOpaque s) =
diff --git a/src/Futhark/Internalise/Entry.hs b/src/Futhark/Internalise/Entry.hs
--- a/src/Futhark/Internalise/Entry.hs
+++ b/src/Futhark/Internalise/Entry.hs
@@ -7,7 +7,7 @@
 where
 
 import Control.Monad
-import Control.Monad.State
+import Control.Monad.State.Strict
 import Data.List (find, intersperse)
 import Data.Map qualified as M
 import Futhark.IR qualified as I
@@ -85,7 +85,11 @@
   case find ((== name) . fst) ts of
     Just (_, t')
       | t /= t' ->
-          error $ "Duplicate definition of entry point type " <> E.prettyString name
+          error . unlines $
+            [ "Duplicate definition of entry point type " <> E.prettyString name,
+              show t,
+              show t'
+            ]
     _ -> I.OpaqueTypes ts <> I.OpaqueTypes [(name, t)]
 
 isRecord ::
@@ -124,6 +128,23 @@
   where
     opaqueField e_t i_ts = snd <$> entryPointType types e_t i_ts
 
+opaqueRecordArray ::
+  VisibleTypes ->
+  Int ->
+  [(Name, E.EntryType)] ->
+  [I.TypeBase I.Rank Uniqueness] ->
+  GenOpaque [(Name, I.EntryPointType)]
+opaqueRecordArray _ _ [] _ = pure []
+opaqueRecordArray types rank ((f, t) : fs) ts = do
+  let (f_ts, ts') = splitAt (internalisedTypeSize $ E.entryType t) ts
+  f' <- opaqueField t f_ts
+  ((f, f') :) <$> opaqueRecordArray types rank fs ts'
+  where
+    opaqueField (E.EntryType e_t _) i_ts =
+      snd <$> entryPointType types (E.EntryType e_t' Nothing) i_ts
+      where
+        e_t' = E.arrayOf (E.Shape (replicate rank E.anySize)) e_t
+
 isSum ::
   VisibleTypes ->
   E.TypeExp E.Exp VName ->
@@ -159,6 +180,10 @@
       ets' <- map snd <$> zipWithM (entryPointType types) ets (map (map (ts !!)) is')
       pure $ zip ets' $ map (map (+ 1)) is' -- Adjust for tag.
 
+entryPointTypeName :: I.EntryPointType -> Name
+entryPointTypeName (I.TypeOpaque v) = v
+entryPointTypeName (I.TypeTransparent {}) = error "entryPointTypeName: TypeTransparent"
+
 entryPointType ::
   VisibleTypes ->
   E.EntryType ->
@@ -180,16 +205,31 @@
   | otherwise = do
       case E.entryType t of
         E.Scalar (E.Record fs)
-          | not $ null fs ->
+          | not $ null fs -> do
               let fs' = recordFields types fs $ E.entryAscribed t
-               in addType desc . I.OpaqueRecord =<< opaqueRecord types fs' ts
+              addType desc . I.OpaqueRecord =<< opaqueRecord types fs' ts
         E.Scalar (E.Sum cs) -> do
           let (_, places) = internaliseSumTypeRep cs
               cs' = sumConstrs types cs $ E.entryAscribed t
               cs'' = zip (map fst cs') (zip (map snd cs') (map snd places))
           addType desc . I.OpaqueSum (map valueType ts)
             =<< opaqueSum types cs'' (drop 1 ts)
+        E.Array _ shape (E.Record fs)
+          | not $ null fs -> do
+              let fs' = recordFields types fs $ E.entryAscribed t
+                  rank = E.shapeRank shape
+                  ts' = map (strip rank) ts
+                  record_t = E.Scalar (E.Record fs)
+              ept <- snd <$> entryPointType types (E.EntryType record_t Nothing) ts'
+              addType desc . I.OpaqueRecordArray rank (entryPointTypeName ept)
+                =<< opaqueRecordArray types rank fs' ts
+        E.Array _ shape et -> do
+          let ts' = map (strip (E.shapeRank shape)) ts
+          ept <- snd <$> entryPointType types (E.EntryType (E.Scalar et) Nothing) ts'
+          addType desc . I.OpaqueArray (E.shapeRank shape) (entryPointTypeName ept) $
+            map valueType ts
         _ -> addType desc $ I.OpaqueType $ map valueType ts
+
       pure (u, I.TypeOpaque desc)
   where
     u = foldl max Nonunique $ map I.uniqueness ts
@@ -197,6 +237,9 @@
       maybe (nameFromText $ prettyTextOneLine t') typeExpOpaqueName $
         E.entryAscribed t
     t' = E.noSizes (E.entryType t) `E.setUniqueness` Nonunique
+    strip k (I.Array pt (I.Rank r) t_u) =
+      I.arrayOf (I.Prim pt) (I.Rank (r - k)) t_u
+    strip _ ts_t = ts_t
 
 entryPoint ::
   VisibleTypes ->
diff --git a/src/Futhark/Optimise/InliningDeadFun.hs b/src/Futhark/Optimise/InliningDeadFun.hs
--- a/src/Futhark/Optimise/InliningDeadFun.hs
+++ b/src/Futhark/Optimise/InliningDeadFun.hs
@@ -11,8 +11,10 @@
 import Control.Monad.Identity
 import Control.Monad.State
 import Control.Parallel.Strategies
+import Data.Functor (($>))
 import Data.List (partition)
 import Data.Map.Strict qualified as M
+import Data.Maybe
 import Data.Set qualified as S
 import Futhark.Analysis.CallGraph
 import Futhark.Analysis.SymbolTable qualified as ST
@@ -106,17 +108,85 @@
             to_inline_later
 
 calledOnce :: CallGraph -> S.Set Name
-calledOnce = S.fromList . map fst . filter ((== 1) . snd) . M.toList . numOccurences
+calledOnce =
+  S.fromList . map fst . filter ((== 1) . snd) . M.toList . numOccurences
 
 inlineBecauseTiny :: Prog SOACS -> S.Set Name
 inlineBecauseTiny = foldMap onFunDef . progFuns
   where
     onFunDef fd
-      | (length (bodyStms (funDefBody fd)) < 2)
+      | (length (bodyStms (funDefBody fd)) <= k)
           || ("inline" `inAttrs` funDefAttrs fd) =
           S.singleton (funDefName fd)
       | otherwise = mempty
+      where
+        k = length (funDefRetType fd) + length (funDefParams fd)
 
+progStms :: Prog SOACS -> Stms SOACS
+progStms prog =
+  progConsts prog <> foldMap (bodyStms . funDefBody) (progFuns prog)
+
+data Used = InSOAC | InAD deriving (Eq, Ord, Show)
+
+directlyCalledInSOACs :: Prog SOACS -> M.Map Name Used
+directlyCalledInSOACs = flip execState mempty . mapM_ (onStm Nothing) . progStms
+  where
+    onBody :: Maybe Used -> Body SOACS -> State (M.Map Name Used) ()
+    onBody u = mapM_ (onStm u) . bodyStms
+    onStm u stm = onExp u (stmExp stm) $> stm
+    onExp (Just u) (Apply fname _ _ _) = modify $ M.insert fname u
+    onExp Nothing Apply {} = pure ()
+    onExp u e = walkExpM (walker u) e
+    onSOAC u soac = void $ traverseSOACStms (const (traverse (onStm u'))) soac
+      where
+        u' = max u $ Just $ usage soac
+    usage JVP {} = InAD
+    usage VJP {} = InAD
+    usage _ = InSOAC
+    walker u =
+      (identityWalker :: Walker SOACS (State (M.Map Name Used)))
+        { walkOnBody = const (onBody u),
+          walkOnOp = onSOAC u
+        }
+
+-- Expand set of function names with all reachable functions.
+withTransitiveCalls :: CallGraph -> M.Map Name Used -> M.Map Name Used
+withTransitiveCalls cg fs
+  | fs == fs' = fs
+  | otherwise = withTransitiveCalls cg fs'
+  where
+    look :: (Name, Used) -> M.Map Name Used
+    look (f, u) = M.fromList $ map (,u) (S.toList (allCalledBy f cg))
+    fs' = foldr (M.unionWith max . look) fs $ M.toList fs
+
+calledInSOACs :: CallGraph -> Prog SOACS -> M.Map Name Used
+calledInSOACs cg prog = withTransitiveCalls cg $ directlyCalledInSOACs prog
+
+-- Inline those functions that are used in SOACs, and which involve
+-- arrays of any kind, as well as any functions used in AD.
+inlineBecauseSOACs :: CallGraph -> Prog SOACS -> S.Set Name
+inlineBecauseSOACs cg prog =
+  S.fromList $ mapMaybe onFunDef (progFuns prog)
+  where
+    called = calledInSOACs cg prog
+    isArray = not . primType
+    inline _ InAD = True
+    inline fd InSOAC =
+      any (isArray . paramType) (funDefParams fd)
+        || any (isArray . fst) (funDefRetType fd)
+        || arrayInBody (funDefBody fd)
+    onFunDef fd = do
+      guard $ maybe False (inline fd) $ M.lookup (funDefName fd) called
+      Just $ funDefName fd
+    arrayInBody = any arrayInStm . bodyStms
+    arrayInStm stm =
+      any isArray (patTypes (stmPat stm)) || arrayInExp (stmExp stm)
+    arrayInExp (Match _ cases defbody _) =
+      any arrayInBody $ defbody : map caseBody cases
+    arrayInExp (Loop _ _ body) =
+      arrayInBody body
+    arrayInExp _ = False
+
 -- Conservative inlining of functions that are called just once, or
 -- have #[inline] on them.
 consInlineFunctions :: (MonadFreshNames m) => Prog SOACS -> m (Prog SOACS)
@@ -125,10 +195,10 @@
   where
     cg = buildCallGraph prog
 
--- Inline everything that is not #[noinline].
+-- Inline aggressively; in particular most things called from a SOAC.
 aggInlineFunctions :: (MonadFreshNames m) => Prog SOACS -> m (Prog SOACS)
 aggInlineFunctions prog =
-  inlineFunctions 3 cg (S.fromList $ map funDefName $ progFuns prog) prog
+  inlineFunctions 3 cg (inlineBecauseTiny prog <> inlineBecauseSOACs cg prog) prog
   where
     cg = buildCallGraph prog
 
diff --git a/src/Futhark/Pass/ExtractKernels.hs b/src/Futhark/Pass/ExtractKernels.hs
--- a/src/Futhark/Pass/ExtractKernels.hs
+++ b/src/Futhark/Pass/ExtractKernels.hs
@@ -289,8 +289,7 @@
            )
     unbalancedStm _ (BasicOp _) =
       False
-    unbalancedStm _ (Apply fname _ _ _) =
-      not $ isBuiltInFunction fname
+    unbalancedStm _ Apply {} = False
 
 sequentialisedUnbalancedStm :: Stm SOACS -> DistribM (Maybe (Stms SOACS))
 sequentialisedUnbalancedStm (Let pat _ (Op soac@(Screma _ _ form)))