C++ Error: Access violation reading location 0x004475c4
<p>I am struggling to debug my program and keep getting the following error: </p>

<p>Unhandled exception at 0x004160b5 in 3Dmesher.exe: 0xC0000005: Access violation reading location 0x004475c4.</p>

<p>From my very limited understanding of C++, this is probably caused by a pointer trying to access an area in memory it doesn't have access to?? I will include the code, but be warned...it is awful!! (still <em>very</em> new to programming). It's supposed to build a 3D mesh from a 2D slice stored in a csv, and then perform FDTD update equations on each node in the mesh, depending on the nodetype (algorithm selected via switch statement based on index associated with each node).</p>

<p>Whilst trying to get this to run, I <em>am</em> simultaneously reading a beginner's book to learning c++, so please bear with me!</p>

<p>Thanks in advance :-)</p>

<p>Code here:</p>

<pre><code>// SonicFDTDmodeller.cpp : Defines the entry point for the console application.
//

//#include "stdafx.h"
#include &lt;fstream&gt;
#include &lt;iostream&gt;
#include &lt;string.h&gt;
#include &lt;sstream&gt;

using namespace std;
typedef struct {
int x;
int y;
int z;
int i;
} pos;

//int _tmain(int argc, _TCHAR* argv[])
int main(int argc, char* argv[])
{
    // Declare Variables
    //fstream csv;
    string filename = "";
    string root = "C:\\Documents and Settings\\shnh500\\My Documents\\Matlab\\Alex's Code\\";   
    string num;
    stringstream ss;
    int L = 0;
    int W = 0;
    int H;
    int *M;
    int *MESH;
    int cnt = 0;
    int nodeType;

    // Ask for Desired Height-------------
    while (true)
    {
        cout &lt;&lt; endl &lt;&lt; "Enter Desired Height: ";
        getline(cin,num);
        ss.clear();
        ss.str("");
        ss &lt;&lt; num;      //give input to stringstream
        if (ss &gt;&gt; H)    //try to convert stringstream to int
            break;      //conversion sucessful exit while loop 
        cout &lt;&lt; "Invalid number, please try again" &lt;&lt; endl;
    }
    // End Ask for desired height --------

    // Host parameters (air, constant temperature) (Most not used but included for posterity)
    //float Rho = 1.21;                 // Medium density (kg/m^3) of host
    //float K = 1.42 *pow (10.0,5.0);   // Bulk modulus (resistance to compression) of medium (Pa) of host
    //float c = sqrt(K/Rho);            // Ambient speed of sound in host (air), (m/s)
    float dx = 0.005;                   // Spatial resolution (m) for 2 dimensions
    //float dt = dx/(sqrt (3.0) *c);    // The timestep (to satisfy courrant limit)
    //int Fs = ceil (1/dt);             // Sampling frequency
    float courant = 1 / sqrt (3.0);     // Courant limit

    cout &lt;&lt; "\n Read CSV File" &lt;&lt; endl &lt;&lt; endl;
    cout &lt;&lt; "Enter File Name: ";
    getline(cin,filename);

    // Read CSV file determine L and W---------
    std::ifstream csv;
    std::string line;
    std::string cell;
    std::stringstream nodeStream;
    std::stringstream lineStream;

    csv.open(root + filename.c_str()); // Remember to input double slashes in the file path input i.e. D:\\Folder\\test.csv

    if (csv.good())
    {   
        while(std::getline(csv,line))
        {
            stringstream lineStream(line);
            string cell;
            while(getline(lineStream,cell,','))
            {               
                // Count Node W
                if (L == 0){W++;}
            }
            // Count Node L
            L++;
        }
    }
    csv.close();
    csv.clear();
    // End Read  CSV ----------

    // Now Define 2D node grid using L and W
    // We will store the 2D grid as a 1D array of L*W entries

    M = new int [L*W];

    // open File again and read the values

    csv.open(root + filename.c_str()); // Remember to input double slashes in the file path input i.e. D:\\Folder\\test.csv

    if (csv.good())
    {   
        while(std::getline(csv,line))
        {
            lineStream.clear();
            // Pop the string into the linestream so it's easier to deal with.
            lineStream &lt;&lt; line;
            // Grab the bit up to the next ',' delimiter and stick it in a string
            while(getline(lineStream,cell,','))
            {
              nodeStream.clear(); // reset the stream like before
              nodeStream &lt;&lt;  cell; // pop that string into the node stream
              nodeStream &gt;&gt;  nodeType; // and then dump the node into the int container
              M[cnt] = nodeType;
              cnt++;              
            }
        }
    }
    csv.close();
    csv.clear();

    // End Read  CSV ----------

    cout &lt;&lt; "L :" &lt;&lt; L &lt;&lt; endl;
    cout &lt;&lt; "W :" &lt;&lt; W &lt;&lt; endl;

    // Now extend the 2D slice to 3D
    MESH = new int [L*H*W];
    int st;
    int en;

    // Do all slice apart from the top and bottom
    st = (L*W); // start
    en = (L*H*W)-(L*W); // end
    cnt = 0;
    for (int i = st;i &lt; en;i++)
    {
        MESH[i] = M[cnt];        
        cnt++;
        if (cnt == L*W)
        {  // reset counter
            cnt = 0;
        }
    }

    // Can delete M now
    delete []M;

    // Do top and bottom slices
    // First make ALL nodes in top and bottom slices surface boundaries
    for (int i = 0; i &lt; L*W; i++)
    {
        MESH[i] = 1; // bottom slice
        MESH[en + i] = 1; //top slice
    }

    // now do the edge nodes
    for (int i = 0; i &lt; L;i++)
    {
        //bottom slice
        MESH[i*W] = 15; //left side edge
        MESH[((i+1) *W)-1] = 15; //right side edge
        //top slice
        MESH[en + (i*W)] = 15; //left side edge
        MESH[en + ((i+1) *W) -1] = 15; //right side edge
    }
    for (int i = 0; i &lt; W;i++)
    {
        //bottom slice
        MESH[i] = 15; // front side edge
        MESH[(L*W)-W + i] = 15; //back side edge
        //top slice
        MESH[en + i] = 15;  // front side edge
        MESH[en + (L*W) - W + i] = 15; //right side edge

    }
    // now do the eight corner nodes
    MESH[0] = 6;
    MESH[W-1] = 6;
    MESH[(L*W)-W] = 6;
    MESH[(L*W)-1] = 6;
    MESH[en] = 6;
    MESH[en+W-1] = 6;
    MESH[en+(L*W)-W] = 6;
    MESH[en+(L*W)-1] = 6;

    cout &lt;&lt; "3D Meshing Complete." &lt;&lt; endl;

    // Perform Model    

    // Number of timesteps
    int T;
    T = 1;  

    // Create an array to hold the correction signal...
    float *h;
    h = new float [T];

    // Initiate excitation signal array...
    float *excitation;
    excitation = new float [T];
    for (int i = 0; i &lt; T; i++)
    {
        if (i==0)
        {
            excitation[i] = 1.0;
        } else {
            excitation[i] = 0.0;
        }
    }

    // Absorption coefficient (1 for max, 0 for free-field)
    int alpha = 1;  

    // Assign pressure memory into 1D array (store past and present states)
    float *p; // Declare pointer to hold address of float
    p = new float [3*L*W*H];

    // Set receiver position in nodes, index from 0 to L-1 , 0 to W-1
    pos rec;
    rec.y = W/2;
    rec.x = ceil (2.16/dx);
    rec.z = (H/2)-1;
    // Calculate Position in the 1D array
    rec.i = (rec.z*L*W)+(rec.x*W)+rec.y;

    for (int i = 0; i &lt; 3*L*W*H; i++)
    {   // initialize to zero 
        p[i] = 0.0;
    }

    // Set input position in nodes, index from 0 to L-1 , 0 to W-1
    pos src;
    src.y = (W/2)-1;
    src.x = (floor (0.5/dx));
    src.z = (H/2)-1;
    // Calculate Position in the 1D array
    src.i = (src.z*L*W)+(src.x*W)+src.y;

    int n,m,l,val;

    // Prepare to write results to csv file
    std::ofstream my_output("C:\\Documents and Settings\\shnh500\\My Documents\\Matlab\\Alex's Code\\3D_P.csv");
    std::ofstream my_ir("C:\\Documents and Settings\\shnh500\\My Documents\\Matlab\\Alex's Code\\3D_IR.csv");
    std::ofstream my_corr("C:\\Documents and Settings\\shnh500\\My Documents\\Matlab\\Alex's Code\\3D_H.csv");

    // At each time step perform the following updates...
    for (int t = 0; t &lt; T; t++)
    {   
        // Perform AIR Update (update and store pressure at past, present and future timesteps)
        for (l=0;l&lt;H;l++)
        {
            std::cout &lt;&lt; '\n' &lt;&lt; endl;
        for (n=0;n&lt;L;n++)
        {
            std::cout &lt;&lt; '\n' &lt;&lt; endl;
            for (m=0;m&lt;W;m++)
            {    
                std::cout &lt;&lt; MESH[(l*W*L)+(n*W+m)];
                switch(MESH[(l*W*L)+(n*W+m)])
                {
                // For all INCLUSION EDGE BOUNDARY nodes do this...             
                case 12:            
                    p[(l*W*L)+(n*W+m)+(2*L*W*H)] = 0;
                    break;

                // For all INCLUSION CORNER nodes do this...
                case 11:
                    p[(l*W*L)+(n*W+m)+(2*L*W*H)] = 0;
                    break;

                // For all AIR nodes INSIDE INCLUSION do this...
                case 13:
                    //p[(l*W*L)+(n*W+m)] = 0;
                    break;

                // For all other AIR nodes do this...
                case 0:                     
                    // Capture/apply the correction signal
                    if ((l*W*L)+(n*W+m)==(src.z*W*L)+(src.x*W)+src.y)
                    {
                        h[t] = pow (courant,2) * (p[((l+1) *W*L)+(n*W+m)] + p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))] ) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] - p[(l*W*L)+(n*W+m)+(L*W*H)] + excitation[t];
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = pow (courant,2) * (p[((l+1) *W*L)+(n*W+m)] + p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))] ) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] - p[(l*W*L)+(n*W+m)+(L*W*H)] + excitation[t];//- src_compensation;
                    } else {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = pow (courant,2) * (p[((l+1) *W*L)+(n*W+m)] + p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))] ) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] - p[(l*W*L)+(n*W+m)+(L*W*H)];
                    }
                    break;

                // For all SURFACE BOUNDARY nodes do this...
                case 1: 
                    // Right
                    if (l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    // Left
                    } else if (l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    // Front
                    } else if (n==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    // Back
                    } else if (n==L-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+(n*W+(m-1))])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    // Top
                    } else if (m==W-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    // Bottom
                    } else if (m==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[((l+1) *W*L)+(n*W+m)])) + 2*(1-(3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+alpha*courant);
                    }
                    break;

                // For all CORNER nodes do this...
                case 6:
                    // Top front left 
                    if (m==0 &amp;&amp; n==0 &amp;&amp; l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n+1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top front right 
                    } else if (m==0 &amp;&amp; n==L-1 &amp;&amp; l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom front left
                    } else if (m==0 &amp;&amp; n==0 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n+1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom front right
                    } else if (m==0 &amp;&amp; n==L-1 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom back left
                    } else if (m==W-1 &amp;&amp; n==0 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top back left
                    } else if (m==W-1 &amp;&amp; n==0 &amp;&amp; l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom back right
                    } else if (m==W-1 &amp;&amp; n==L-1 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top back right
                    } else if (m==W-1 &amp;&amp; n==L-1 &amp;&amp; l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) * (2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1 - (3*pow (courant,2))) * p[(l*W*L)+(n*W+m)] + (3*alpha*courant-1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    }
                    break;

                    // For all EDGE BOUNDARY nodes do this...
                case 15: 
                    // Top front
                    if (l==H-1 &amp;&amp; m==0 &amp;&amp; n!=0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom front
                    } else if (l==0 &amp;&amp; m==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom back
                    } else if (l==0 &amp;&amp; m==W-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top back
                    } else if (l==H-1 &amp;&amp; m==W-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((n-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top left
                    } else if (l==H-1 &amp;&amp; n==0)
                    {
            ! IT BREAKS HERE !          p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((m+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Top right
                    } else if (n==L-1 &amp;&amp; l==H-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((m-1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom left
                    } else if (n==0 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((l+1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Bottom right
                    } else if (n==L-1 &amp;&amp; l==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+((m-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Front left
                    } else if (n==0 &amp;&amp; m==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Front right
                    } else if (n==L-1 &amp;&amp; m==0)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m+1))] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Back left
                    } else if (n==0 &amp;&amp; m==W-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n+1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    // Back right
                    } else if (n==L-1 &amp;&amp; m==W-1)
                    {
                        p[(l*W*L)+(n*W+m)+(2*L*W*H)] = (pow (courant,2) *(2*(p[((l+1) *W*L)+(n*W+m)] + p[(l*W*L)+(n*W+(m-1))] + p[(l*W*L)+((n-1) *W+m)] + p[(l*W*L)+((l-1) *W+m)])) + 2*(1-3*pow (courant,2)) * p[(l*W*L)+(n*W+m)] + (2*alpha*courant - 1) * p[(l*W*L)+(n*W+m)+(L*W*H)])/(1+3*alpha*courant);
                    }
                    break;
                }                   
            }
        }
        }

        // Record pressure at receiver location
        my_corr &lt;&lt; h[t] &lt;&lt; ",";
        my_ir &lt;&lt; p[rec.i] &lt;&lt; ",";

        // Past &lt; Present :: Present &lt; Future
        for (l=0;l&lt;H;l++)
        {
            for (m=0;m&lt;W;m++)
            {
                for (n=0;n&lt;L;n++)
                {
                    p[(l*W*L)+(n*W+m)+(L*W*H)] = p[(l*W*L)+(n*W+m)];
                    p[(l*W*L)+(n*W+m)] = p[(l*W*L)+(n*W+m)+(2*L*W*H)];
                }
            }
        }
    }       

    delete []h;
    delete []my_ir;
    delete []p; 
    delete []MESH;
    system("pause");
    return 0;
}
</code></pre>
