你忘记了一个维度，以及分配内存的开销。所示的代码在第三个维度上分配内存的效率非常低，导致了太多的开销。

float*** a = new float**[N];

这将分配大约22000 * sizeof(float **)，大约176kb，可以忽略不计。

a[m] = new float*[M - 1];

这里的一次分配是给44099 * sizeof(float *)的，但是你会得到22000个这样的. 22000 * 44099 * sizeof(float *)，或者大约7.7gb的额外内存，这就是你停止计数的地方，但是你的代码还没有完成，还有很长的路要走。

a[m][n] = new float[2];

这是一个8字节的单次分配，但是这个分配将被完成22000 * 44099次。这是另一个*7.7gb被冲走。你现在大约需要分配超过15gig的应用程序所需的内存。
但是每个分配都不是免费的，并且new float[2]需要 * 超过 * 8个字节。每个单独分配的块必须由C++库在内部跟踪，以便delete可以回收它。最简单的基于链表的堆分配实现需要一个前向指针、一个后向指针、以及分配的块中有多少字节的计数。假设不需要为对齐目的填充任何内容，则在64位平台上，每次分配至少会有24字节的开销。
现在，由于第三个维分配了22000 * 44099个分配，第二个维分配了22000个分配，第一个维分配了一个分配：如果我用手指数，这将需要（22000 * 44099 + 22000 + 1） 24，或者另外22GB的内存，仅仅是为了消耗最简单、基本的内存分配方案的开销。
如果我没算错的话，我们现在使用最简单、最可能的堆分配跟踪所需的RAM高达38 GB，您的C++实现可能会使用稍微复杂一点的堆分配逻辑，但开销较大。
去掉new float[2]，计算矩阵的大小，new是一个7.7gb的块，然后计算其余指针应该指向哪里，另外，为矩阵的第二维分配一块内存，计算第一维的指针。
你的分配代码应该只执行三条new语句，一条用于第一个维度指针，一条用于第二个维度指针，还有一条用于组成第三个维度的大块数据。

下面的示例只是对已经给出的一个答案的补充，它基本上是对这里给出的关于如何创建连续2D数组的答案的扩展，并且只说明了对new[]的3次调用的用法。
优点是保留了通常用于三重指针的[][][]语法（尽管我强烈建议不要使用这样的“三星”来编写代码，但我们已经拥有了），缺点是为指针分配了更多的内存，增加了用于数据的单个内存池。

#include <iostream>
#include <exception>

template <typename T>
T*** create3DArray(unsigned pages, unsigned nrows, unsigned ncols, const T& val = T())
{
    T*** ptr = nullptr;  // allocate pointers to pages
    T** ptrMem = nullptr;
    T* pool = nullptr;
    try 
    {
        ptr = new T**[pages];  // allocate pointers to pages
        ptrMem = new T*[pages * nrows]; // allocate pointers to pool
        pool = new T[nrows*ncols*pages]{ val };  // allocate pool

        // Assign page pointers to point to the pages memory,
        // and pool pointers to point to each row the data pool
        for (unsigned i = 0; i < pages; ++i, ptrMem += nrows)
        {
            ptr[i] = ptrMem;
            for (unsigned j = 0; j < nrows; ++j, pool += ncols)
                ptr[i][j] = pool;
        }
        return ptr;
     }
     catch(std::bad_alloc& ex)
     {
         // rollback the previous allocations
        delete [] ptrMem;
        delete [] ptr;
        throw ex; 
    }
}

template <typename T>
void delete3DArray(T*** arr)
{
    delete[] arr[0][0]; // remove pool
    delete[] arr[0];  // remove the pointers
    delete[] arr;     // remove the pages
}

int main()
{
    double ***dPtr = nullptr;
    try 
    {
        dPtr = create3DArray<double>(4100, 5000, 2);
    }
    catch(std::bad_alloc& )
    {
        std::cout << "Could not allocate memory";
        return -1;
    }
    dPtr[0][0][0] = 10;  // for example
    std::cout << dPtr[0][0][0] << "\n";
    delete3DArray(dPtr);  // free the memory
}

Live Example
这是一个基本的RAII类。它还没有经过彻底的测试，但是应该可以正常工作。

#include <exception>
#include <algorithm>
#include <iostream>

template <typename T>
class Array3D
{
    T*** data_ptr;
    unsigned m_rows;
    unsigned m_cols;
    unsigned m_pages;

    T*** create3DArray(unsigned pages, unsigned nrows, unsigned ncols, const T& val = T()) 
    {
        T*** ptr = nullptr;  // allocate pointers to pages
        T** ptrMem = nullptr;
        T* pool = nullptr;
        try
        {
            ptr = new T * *[pages];  // allocate pointers to pages
            ptrMem = new T * [pages * nrows]; // allocate pointers to pool
            pool = new T[nrows * ncols * pages]{ val };  // allocate pool

            // Assign page pointers to point to the pages memory,
            // and pool pointers to point to each row the data pool
            for (unsigned i = 0; i < pages; ++i, ptrMem += nrows)
            {
                ptr[i] = ptrMem;
                for (unsigned j = 0; j < nrows; ++j, pool += ncols)
                    ptr[i][j] = pool;
            }
            return ptr;
        }
        catch (std::bad_alloc& ex)
        {
            // rollback the previous allocations
            delete[] ptrMem;
            delete[] ptr;
            throw ex;
        }
    }

    public:
        typedef T value_type;
        T*** data() { return data_ptr;  }

        unsigned get_num_pages() const { return m_pages; }
        unsigned get_num_rows() const { return m_rows; }
        unsigned get_num_cols() const { return m_cols; }

        Array3D() : data_ptr{}, m_rows{}, m_cols{}, m_pages{} {}

        Array3D(unsigned pages, unsigned rows, unsigned cols, const T& val = T())
        {
            if ( pages == 0 )
                throw std::invalid_argument("number of pages is 0");
            if (rows == 0)
                throw std::invalid_argument("number of rows is 0");
            if (cols == 0)
                throw std::invalid_argument("number of columns is 0");
            data_ptr = create3DArray(pages, rows, cols, val);
            m_pages = pages;
            m_rows = rows;
            m_cols = cols;
        }

        ~Array3D()
        {
            if (data_ptr)
            {
                delete[] data_ptr[0][0]; // remove pool
                delete[] data_ptr[0];  // remove the pointers
                delete[] data_ptr;     // remove the pages
            }
        }

        Array3D(const Array3D& rhs) : m_rows(rhs.m_rows), m_cols(rhs.m_cols), m_pages(rhs.m_pages)
        {
            data_ptr = create3DArray(m_pages, m_rows, m_cols);
            std::copy(&rhs.data_ptr[0][0][0], &rhs.data_ptr[m_pages - 1][m_rows - 1][m_cols], &data_ptr[0][0][0]);
        }

        Array3D& operator=(const Array3D& rhs)
        {
            if (&rhs != this)
            {
                Array3D temp(rhs);
                swap(*this, temp);
            }
            return *this;
        }

        Array3D(Array3D&& rhs) noexcept : data_ptr(rhs.data_ptr), m_pages(rhs.m_pages), m_rows(rhs.m_rows), m_cols(rhs.m_cols)
        {
            rhs.data_ptr = nullptr;
        }

        Array3D& operator =(Array3D&& rhs) noexcept
        {
            if (&rhs != this)
            {
                swap(rhs, *this);
            }
            return *this;
        }

        void swap(Array3D& left, Array3D& right)
        {
            std::swap(left.data_ptr, right.data_ptr);
            std::swap(left.m_cols, right.m_cols);
            std::swap(left.m_rows, right.m_rows);
            std::swap(left.m_pages, right.m_pages);
        }

        T** operator[](unsigned page)
        {
            return data_ptr[page];
        }

        const T** operator[](unsigned page) const
        {
            return data_ptr[page];
        }
    };

int main()
{
    try
    {
        Array3D<double> dPtr(10, 10, 10);
        dPtr[0][0][0] = 20;
        dPtr[0][0][3] = -23;
        std::cout << dPtr[0][0][0] << " " << dPtr[0][0][3] << "\n";
        Array3D<double> test = dPtr;
        std::cout << test[0][0][0] << " " << test[0][0][3] << "\n";
        Array3D<double> test2;
        test2 = test;
        std::cout << test2[0][0][0] << " " << test2[0][0][3] << "\n";
    }
    catch (std::exception& ex)
    {
        std::cout << ex.what();
    }
}

Live Example

这可能是问题的简化版本，但您使用的数据结构（“三星”数组）几乎从来都不是您想要的。如果您创建一个像这样的密集矩阵，并为每个元素分配空间，那么进行数百万次微小分配根本没有任何优势。如果您想要一个稀疏矩阵，通常需要压缩稀疏行这样的格式。
如果数组是“矩形”的（或者我认为3D数组是“四四方方”的），并且所有的行和列都是相同大小的，那么与分配单个内存块相比，这种数据结构纯粹是浪费，因为要执行数百万次微小的分配，为数百万个指针分配空间，并且丢失内存的局部性。
这个样板文件为动态3-D数组创建了一个零成本的抽象。同时存储底层一维std::vector的长度和各个维度是多余的。）API使用a(i, j, k)作为a[i][j][k]的等价物，使用a.at(i,j,k)作为边界检查的变量。
这个API还有一个选项，可以用索引函数f(i,j,k)填充数组。如果你调用a.generate(f)，它会设置每个a(i,j,k) = f(i,j,k)。理论上，这个强度减少了内部循环中的偏移量计算，使它更快。API还可以将生成函数作为array3d<float>(M, N, P, f)传递给构造函数。你可以随意扩展它。

#include <cassert>
#include <cstddef>
#include <cstdlib>
#include <functional>
#include <iomanip>
#include <iostream>
#include <vector>

using std::cout;
using std::endl;
using std::ptrdiff_t;
using std::size_t;

/* In a real-world implementation, this class would be split into a
 * header file and a definitions file.
 */
template <typename T>
  class array3d {
    public:
    using value_type = T;
    using size_type = size_t;
    using difference_type = ptrdiff_t;
    using reference = T&;
    using const_reference = const T&;
    using pointer = T*;
    using const_pointer = const T*;
    using iterator = typename std::vector<T>::iterator;
    using const_iterator = typename std::vector<T>::const_iterator;
    using reverse_iterator = typename std::vector<T>::reverse_iterator;
    using const_reverse_iterator = typename
      std::vector<T>::const_reverse_iterator;

/* For this trivial example, I don’t define a default constructor or an API
 * to resize a 3D array.
 */
    array3d( const ptrdiff_t rows,
             const ptrdiff_t cols,
             const ptrdiff_t layers )
    {
      const ptrdiff_t nelements = rows*cols*layers;

      assert(rows > 0);
      assert(cols > 0);
      assert(layers > 0);
      assert(nelements > 0);

      nrows = rows;
      ncols = cols;
      nlayers = layers;
      storage.resize(static_cast<size_t>(nelements));
    }

/* Variant that initializes an array with bounds and then fills each element
 * (i,j,k) with a provided function f(i,j,k).
 */
    array3d( const ptrdiff_t rows,
             const ptrdiff_t cols,
             const ptrdiff_t layers,
             const std::function<T(ptrdiff_t, ptrdiff_t, ptrdiff_t)> f )
    {
      const ptrdiff_t nelements = rows*cols*layers;

      assert(rows > 0);
      assert(cols > 0);
      assert(layers > 0);
      assert(nelements > 0);

      nrows = rows;
      ncols = cols;
      nlayers = layers;
      storage.reserve(static_cast<size_t>(nelements));

      for ( ptrdiff_t i = 0; i < nrows; ++i )
        for ( ptrdiff_t j = 0; j < ncols; ++j )
          for ( ptrdiff_t k = 0; k < nlayers; ++k )
            storage.emplace_back(f(i,j,k));

      assert( storage.size() == static_cast<size_t>(nelements) );
    }

    // Rule of 5:
    array3d( const array3d& ) = default;
    array3d& operator= ( const array3d& ) = default;
    array3d( array3d&& ) = default;
    array3d& operator= (array3d&&) = default;

    /* a(i,j,k) is the equivalent of a[i][j][k], except that the indices are
     * signed rather than unsigned.  WARNING: It does not check bounds!
     */
    T& operator() ( const ptrdiff_t i,
                    const ptrdiff_t j,
                    const ptrdiff_t k ) noexcept
    {
      return storage[make_index(i,j,k)];
    }

    const T& operator() ( const ptrdiff_t i,
                          const ptrdiff_t j,
                          const ptrdiff_t k ) const noexcept
    {
      return const_cast<array3d&>(*this)(i,j,k);
    }

    /* a.at(i,j,k) checks bounds.  Error-checking is by assertion, rather than
     * by exception, and the indices are signed.
     */
    T& at( const ptrdiff_t i, const ptrdiff_t j, const ptrdiff_t k )
    {
      bounds_check(i,j,k);
      return (*this)(i,j,k);
    }

    const T& at( const ptrdiff_t i,
                 const ptrdiff_t j,
                 const ptrdiff_t k ) const
    {
      return const_cast<array3d&>(*this).at(i,j,k);
    }

/* Given a function or function object f(i,j,k), fills each element of the
 * container with a(i,j,k) = f(i,j,k).
 */
    void generate( const std::function<T(ptrdiff_t,
                                         ptrdiff_t,
                                         ptrdiff_t)> f )
    {
      iterator it = storage.begin();

      for ( ptrdiff_t i = 0; i < nrows; ++i )
        for ( ptrdiff_t j = 0; j < ncols; ++j )
          for ( ptrdiff_t k = 0; k < nlayers; ++k )
            *it++ = f(i,j,k);

      assert(it == storage.end());
    }

/* Could define a larger API, e.g. begin(), end(), rbegin() and rend() from the STL.
 * Whatever you need.
 */

    private:
    ptrdiff_t nrows, ncols, nlayers;
    std::vector<T> storage;

    constexpr size_t make_index( const ptrdiff_t i,
                                 const ptrdiff_t j,
                                 const ptrdiff_t k ) const noexcept
    {
      return static_cast<size_t>((i*ncols + j)*nlayers + k);
    }

    // This could instead throw std::out_of_range, like STL containers.
    constexpr void bounds_check( const ptrdiff_t i,
                                 const ptrdiff_t j,
                                 const ptrdiff_t k ) const
    {
      assert( i >=0 && i < nrows );
      assert( j >= 0 && j < ncols );
      assert( k >= 0 && k < nlayers );
    }
};

// In a real-world scenario, this test driver would be in another source file:

constexpr float f( const ptrdiff_t i, const ptrdiff_t j, const ptrdiff_t k )
{
  return static_cast<float>( k==0 ? 1.0 : -1.0 *
                             ((double)i + (double)j*1E-4));
}

int main(void)
{
  constexpr ptrdiff_t N = 2200, M = 4410, P = 2;
  const array3d<float> a(N, M, P, f);

  // Should be: -1234.4321
  cout << std::setprecision(8) << a.at(1234,4321,1) << endl;

  return EXIT_SUCCESS;
}

值得注意的是，这段代码在技术上包含未定义的行为：它假定带符号整数乘法溢出产生负数，但实际上，如果程序在运行时请求一些不合理的存储器量，则编译器有权生成完全损坏的代码。
当然，如果数组边界是常量，只需声明它们为constexpr，并使用一个具有固定边界的数组。
不幸的是，每一个新的C++程序员首先学习char** argv，因为这使人们认为“二维”数组是指向行的指针的“不规则”数组。
在真实的世界中，这几乎从来都不是最适合这项工作的数据结构。