Решение уравнения теплопроводности численными методами

...
#define X 20
...
{
        double  dt = 0.01,
                   dx = 0.1,
                   dy = 0.02,
                   a  = 0.0001;
        int       tsc; //время с сек.
       double dT[20][17],
                  T[20][17];
 //значение функции f(x) в точках по x
   int f[X] = {8,9,8,9,11,13,15,15,15,14,14,14,13,13,12,12,11,11,10,10};
//----------------------------Обнуление массивов---------------------------
 for (int i=0; i<19; ++i)
   {
       for (int j=0; j<16; ++j)
        {
                dT[i][j] = T[i][j] = 0;
        }
   }
//-----------------------------------------------------------------------------
for (int tsc=1; tsc<5000; ++tsc) //цикл по времени
{
         for (int i=1; i<18; ++i)
         {
               for (int f_iter=1; f_iter<X; ++f_iter)
                 {
                        int j = f[f_iter];
                        //уравнение теплопроводности в дескретной форме
                        double x = ( T[i-1][j] - 2 * T[i][j] + T[i+1][j] );
                        double y = ( T[i][j-1] - 2 * T[i][j] + T[i][j+1] );
                        dT[i][j] = a * dt * ( x/(dx*dx) + y/(dy*dy) + T[i][j]);
                 }
        }
//------------------граничные условия---------------------------------------
                       //сверху+воздействие
                       for (int x=0; x<19; ++x)
                       {
                                T[x][f[x]]  = 100;
                       }
                       //слева
                       for (int y_iter=0; y_iter<f[0]; ++y_iter)
                       {
                                int y = f[y_iter];
                                T[0][y]  = 0;
                       }
                       //снизу
                       for (int x=0; x<19; ++x)
                       {
                                T[x][0] = T[x][1];
                       }
                       //справа
                       for (int y_iter=0; y_iter<f[X-1]; ++y_iter)
                       {
                                int y = f[y_iter];
                                T[19][y]  = 0;
                       }
                       
//--------------------приращение-------------------------------------------
                       for (int x=0; x<19; ++x)
                       {
                                for (int y_iter=0; y_iter<X; ++y_iter)
                                {
                                        int y = f[y_iter];
                                        T[x][y] = T[x][y]+dT[x][y];
                                }
                       }
         Series1 -> AddXY(tsc, T[2][2]);
         Series2 -> AddXY(tsc,f[X-2]);
 }
}