Writing Several 2D Arrays Into Text

[AliMoradi]

I've written a 2D Jacobi solver and I'm trying to write all my 2D arrays into a text file. I'm trying to write all the 2D arrays into a text file. I tried putting the entire MxN array into a vector of size M*N. This method seems a bit archaic. Is there a better way to do this?

#include <iostream>
#include <math.h>
#include <stdio.h>
#include <fstream>
#include <cmath>

using namespace std;

const int N = 21;
const int M = 21;
const int MAX = 101;
double ax = 0;
double bx = 1;
double ay = 0;
double by = 1;
double hx = (bx - ax)/(N - 1);
double hy = (by - ay)/(M - 1);

//define PI
double const pi = 3.141592653589793;

int main()
{
    //Create X and Y Vectors
    double x [N];
    double y [N];
    for (int i = 0; i < N; i++)
    {
        x[i] = ax + hx*i;
        y[i] = ay + hy*i;
    }

    //Define Vectors
    double PHI [N][M]; //Previous Iteration
    double PHIn [N][M]; //New Iteration
    double PHI50 [N][M]; //At 50 Iterations
    double PHI100 [N][M]; //At 100 Iterations

    //Set Boundary Conditions
    for (int i = 0; i < N; i++)
    {
        for (int j = 0; j < M; j++)
        {
            PHI[i][j] = 0;
            PHIn[i][j] = 0;
            PHI50[i][j] = 0;
            PHI100[i][j] = 0;
        }
    }

    //Function For b vector
    //double func(double x,double y)
    //{
        //return 2*(2 - x*x - y*y);
    //}

    //Perform Jacobi Method
    int step = 1;
    while (step < MAX)
    {
        for (int j = 1; j < M-1; j++)
        {
            for (int i = 1; i < N-1; i++)
            {
               PHIn[i][j] = 0.25*(PHI[i-1][j] + PHI[i+1][j] + PHI[i][j-1] + PHI[i][j+1]) - 0.25*(2*(2 - x[i]*x[i] - y[j]*y[j]));

                if (step == 50)
                {
                    PHI50[i][j] = PHIn[i][j];
                }
                if (step == 100)
                {
                    PHI100[i][j] = PHIn[i][j];
                }
            }
        }

        for (int i = 1; i < N-1; i++)
        {
            for (int j = 1; j < M-1; j++)
            {
                PHI[i][j] = PHIn[i][j];
            }
        }

        step++;
    }

    //Output Data Into Text File
    ofstream outputdata("jacobidim.txt");
    for (int i = 0; i<N; i++)
    {
            outputdata << x[i] << "    " << y[i] << endl;
    }
    outputdata.close();

    ofstream outputdata2("jacobi.txt");
    for (int i = 0; i<N; i++)
    {
        for (int j = 0; j<M; j++)
        {
            outputdata2 << PHIn[i][j] << "    " << PHI50[i][j] << "    " << PHI100[i][j] << endl;
        }
    }
    outputdata2.close();

    return 0;
}