Interpolacja funkcja sklejana

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <math.h>

using namespace std;

const int n = 5;
const int dn = 2;
const int rm = n + dn;
const double xt[ n ] = { - 4.0, - 2.0, 0.0, 2.0, 4.0 };
const double yt[ n ] = { 1155.0, 79.0, 3.0, 63.0, 931.0 };
const double dxt[ 2 ] = { - 4.0, 4.0 };
const double dyt[ 2 ] = { - 1124.0, 940.0 };
const double xv = 1.0;
const int dpos[] = { 0, 4 };

void rozwiaz_uklad_rownan( double macierz[ rm ][ rm + 1 ], double wynik[ rm ] )
{
    for( int i = 0; i < rm - 1; i++ )
    {
        for( int j = i + 1; j < rm; j++ )
        {
            double m = - macierz[ j ][ i ] / macierz[ i ][ i ];
            for( int k = i + 1; k <= rm; k++ )
            {
                macierz[ j ][ k ] += m * macierz[ i ][ k ];
            }
        }
    }
   
    for( int i = rm - 1; i >= 0; i-- )
    {
        double sp = macierz[ i ][ rm ];
        for( int j = rm - 1; j >= i + 1; j-- )
        {
            sp -= macierz[ i ][ j ] * wynik[ j ];
        }
        wynik[ i ] = sp / macierz[ i ][ i ];
    }
}

int main()
{
    double macierz[ rm ][ rm + 1 ];
   
    for( int i = 0; i < n; i++ )
    {
        for( int j = 0; j < n; j++ )
        {
            macierz[ i ][ j ] = powf( xt[ i ], j );
        }
       
        for( int j = 1;( n + j <= rm ) &&( j <= i ); j++ )
        {
            macierz[ i ][ n - 1 + j ] = pow( xt[ i ] - xt[ j ], 3.0 );
        }
    }
   
    for( int i = n; i < dn + n; i++ )
    {
        for( int j = 0; j < n; j++ )
        {
            macierz[ i ][ j ] = j * powf( xt[ dpos[ dpos[ i - n ] ], j - 1 );
        }
        for( int j = 1;( i + j <= rm + 1 ) &&( j < dpos[ i - n ] ); j++ )
        {
            macierz[ i ][ n + j - 1 ] = 3.0 * powf( xt[ dpos[ dpos[ i - n ] ] - xt[ j ], 2.0 );
        }
    }
   
    for( unsigned i = 0; i < n; i++ )
    {
        macierz[ i ][ rm ] = yt[ i ];
    }
   
    for( unsigned i = n; i < rm; i++ )
    {
        macierz[ i ][ rm ] = dyt[ i - n ];
    }
   
    for( unsigned y = 0; y < 7; y++ )
    {
        cout << "Macierz:" << endl;
        for( int x = 0; x < 8; x++ )
        {
            cout << macierz[ y ][ x ] << "\t";
        }
        cout << endl;
    }
   
    double wyniki[ 7 ];
   
    rozwiaz_uklad_rownan( macierz, wyniki );
   
    double w0 = wyniki[ 0 ];
    double w1 = wyniki[ 1 ];
    double w2 = wyniki[ 2 ];
    double w3 = wyniki[ 3 ];
    double w4 = wyniki[ 4 ];
    double w5 = wyniki[ 5 ];
    double w6 = wyniki[ 6 ];
   
    cout << "Wartosci zmiennych:" << endl;
    for( unsigned i = 0; i < 7; i++ )
    {
        cout << wyniki[ i ] << endl;
    }
   
    double wynik_koncowy = w0 + w1 * xv + w2 * powf( xv, 2 ) + w3 * powf( xv, 3 ) +
    w4 * powf( xv, 4 );
   
    cout << "Wynik f(" << xv << ") = "
    << wynik_koncowy << endl;
   
    while( 1 ) { }
   
    return 0;
}