Решение уравнения теплопроводности с помощью неявной разностной схемы

#include "math.h"
#include "iostream"
#include"fstream"
using namespace std;
 
 
 
double func( double x , double t )
{ return exp ( - 3.14 * 3.14 * t ) * sin ( 3.14 * x ) ; }
 
 
int main()
{
int n1 = 11 ;
int n,k;
double z;
double A[n1][n1];
double B[11];
double eps[n1];
double X[n1];
double et[n1];
double cou, tau=0.005,h=1.0/(n1-1);
 
cou = 2.0 * tau / ( h * h );
 
for ( int l = 0 ; l < n1 ; l++ )
for( int m = 0 ; m < n1 ; m++ )
A[m][l]=0;
 
for ( int j = 0 ; j < n1 + 1 ; j++ )
{
A[j][j-1] = cou / 2 ;
A[j][j] = - 1 - cou ;
A[j][j+1] = cou / 2 ;
}
 
for ( int i = 0; i < n1; i++ )
{
B[i] = - sin ( i * h * 3.14 );
}
 
n = n1 - 1 ;
eps[0] = - A[0][1] / A[0][0] ;
et[0] = B[0] / A[0][0] ;
 
for(int i=1;i<n;i++)
{
z=A[i][i]+A[i][i-1]*eps[i-1];
eps[i]=-A[i][i+1]/z;
et[i]=(B[i]-A[i][i-1]*et[i-1])/z;
}
 
X[n]=(B[n]-A[n][n-1]*et[n-1])/(A[n][n]+A[n][n-1]*eps[n-1]);
 
for(int i=n-1;i>=0;i--)
X[i]=eps[i]*X[i+1]+et[i];
 
 double max = 0;
for(int i=0; i<n1-1;i++)
if(fabs(func(i*h, tau)-X[i])>max)
max=fabs(func(i*h, tau)-X[i]);
printf("%f", max); 
 
return 0;
}

const int n1 = 11 ;