lab6

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

using namespace std;

#define N 1000

double Func(double x){
    return log(x * x) * sin(x * 0.5) * exp(pow(x, 1.0 / 7.0));
}

double *CoefA(double a, double b, int n){
    double *A = new double[n];
    double h = (b - a) / (double) n;
    for(int i = 0; i < n; i++)  {
        A[i] = Func(a + i * h);
    }
    return A;
}

double **Matrix(double a, double b, int n) {
    double h = (b - a) / (double)n;

    double **matrix = new double*[n];
    for (int i = 0; i < n; i++) {
        matrix[i] = new double[n + 1] {0};
    }

    matrix[0][0]= 1;
    matrix[n - 1][n - 1] = 1;

    for (int i = 1; i < n - 1; i++) {
        for (int j = 0; j <= n; j++) {
            matrix[i][j] = 0;
            if (j == (i - 1)) matrix[i][j] = h;
            if (j == i) matrix[i][j] = 4 * h;
            if (j == (i + 1)) matrix[i][j] = h;
            if (j == n) matrix[i][j] = 6 * (Func(a + h * (i - 1)) - 2 * Func(a + h * i) + Func(a + h * (i + 1)))/h;
        }
    }

    return matrix;
}

double *CoefC(double a, double b, int n) {
    double **matrix = Matrix(a, b, n);
    double *C = new double[n];

    for (int k = 0; k < n - 1; k++){
        double m = 1.0 / matrix[k][k];
        for (int j = 0; j <= n; j++) {
            matrix[k][j] = matrix[k][j] * m;
        }

        for (int i = k + 1; i < n; i++){
            double m = matrix[i][k];
            for (int j = 0; j <= n; j++) {
                matrix[i][j] = matrix[i][j] - matrix[k][j] * m;
            }
        }
    }

    for (int i = 0; i < n; i++) {
        C[n - i - 1] = matrix[n - i - 1][n];
        for (int j = 0; j < i; j++) {
            C[n - i - 1] -= matrix[n - i - 1][n - 1 - j] * C[n - i];
        }
    }
    return C;
}

double *CoefD(double a, double b, double *C, int n){
    double *D = new double[n];
    double h = (b - a) / (double) n;

    D[0] = C[0];
    for (int i = 1; i < n; i++) {
        D[i] = (C[i] - C[i - 1]) / h;
    }
    return D;
}

double *CoefB(double a, double b, double *C, double *D, int n){
    double *B = new double[n];
    double h = (b - a) / (double) n;
    B[0] = 0;

    for(int i = 1; i < n; i++) {
        B[i] = (h * C[i] / 2) - (h * h * D[i] / 2) + (Func(a + h * i) - Func(a + h * (i - 1))) / h;
    }
    return B;
}

double Spline(double x, double a, double b, double *A, double *B, double *C, double *D, int n){
    double h = (b - a) / (double) n;
    int i = 0;

    while (x > (a + i * h)) i++;
    --i;

    double xi = a + i * h;
    return A[i] + B[i] * (x - xi) + C[i] * (x - xi) * (x - xi) / 2 + D[i] * (x - xi) * (x - xi) * (x - xi) / 6;
}

int main(){
    double a = 1, b = 35;

    double *A = CoefA(a, b, N);
    double *C = CoefC(a, b, N);
    double *D = CoefD(a, b, C, N);
    double *B = CoefB(a, b, C, D, N);

    int k = (b - a) * 50;
    double step = (b - a) / (double)k;
    for (int i = 0; i <= k; i++) {
        double x = a + i * step;
        double y = Spline(x, a, b, A, B, C, D, N);

        printf("%4.2f %7f\n", x, y);
    }
}