using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.IO;

namespace HeatEquation {
    class Function {
        double h;
        double tao;

        public int K;
        public int J;

        public double lStepH;
        public double lStepT;

        double[] mu;
        double[] mu1;
        double[] mu2;

        public string[,] graphAcc;

        double e = 2.71828182846;
        public double a = 0.018;

        public Function(double hi, double taoi) {
            h = hi;
            tao = taoi;

            J = Convert.ToInt32(1 / h) + 1;
            K = Convert.ToInt32(1 / tao) + 1;

            lStepH = h - 1 - h * J;
            lStepT = 1 - tao * K;

            double t = 0;
            double x = 0;

            mu = new double[J + 1];
            mu1 = new double[K + 1];
            mu2 = new double[K + 1];

            graphAcc = new string[11, 11];
            double taog = 0.1;
            double hg = 0.1;
            for (int i = 0; i < 11; ++i)
                for (int j = 0; j < 11; ++j)
                    graphAcc[i, j] = Convert.ToString(Value(i * hg, j * taog));

            for (int i = 0; i < J; ++i) {
                mu[i] = Value(x, 0);
                x += h;
            }
            //mu[J] = Value(1, 0);

            for(int i = 0; i < K; ++i) {
                mu1[i] = Value(0, t);
                mu2[i] = Value(1, t);
                t += tao;
            }

            mu1[K] = Value(0, t);
            mu1[K] = Value(0, t);
        }

        public double Value(double x, double t) {
            return 2 * Math.Pow(x, 4) - 3 * Math.Pow(t, 3) + 3 * Math.Pow(t, 2) * x - 2 * Math.Pow(e, x);
        }

        public double f(double x, double t) {
            return -9 * Math.Pow(t, 2) + 6 * t * x - a * (24 * Math.Pow(x, 2) - 2 * Math.Pow(e, x));
        }

        public double dTwoFx(double x, double t) {
            return -a * (48 - 2 * Math.Pow(e, x));
        }

        public double F(double x, double t) {
            return f(x, t) + h * h * dTwoFx(x, t) / 12;
        }


        public double[] Mu() {
            return mu;
        }

        public double[] Mu1() {
            return mu1;
        }

        public double[] Mu2() {
            return mu2;
        }
    }

    class Program {
        static double tao = 0.01;
        static double h = 0.02;

        static void Main(string[] args) {

            Function u = new Function(h, tao);

            double a = u.a;
            double sigma = 0.5 - h * h / (12 * a * tao);

            int J = u.J;
            int K = u.K;

            int xg = Convert.ToInt32(J / 2);

            string[,] graphApp = new string[J, K];

            double[] uNext = new double[J];
            double[] uCur = new double[J];

            double[] error = new double[J];
            double[] rightPart = new double[J];

            double x = h;
            double t = 0;

            double mainDiag = 1 / tao + 2 * a * sigma / (h * h);
            double secnDiag = -a * sigma / (h * h);

            double[] al = new double[J];
            double[] bt = new double[J];

            for (int i = 0; i < J; ++i) {
                uCur[i] = u.Mu()[i];
                graphApp[i, 0] = Convert.ToString(u.Mu()[i]);
            }

            double maxError = 0;

            double k = a * (1 - sigma) / (h * h);

            for (int l = 1; l < K; ++l) {
                rightPart[1] = u.F(x, t + tao / 2) + k * uCur[2] + k * uCur[0] + (1 / tao - 2 * k) * uCur[1] + a * sigma / (h * h) * u.Mu1()[l];
                for (int i = 2; i < J - 2; ++i) {
                    x += h;
                    rightPart[i] = u.F(x, t + tao / 2) + k * uCur[i + 1] + k * uCur[i - 1] + (1 / tao - 2 * k) * uCur[i];
                }
                x += h;
                rightPart[J - 2] = u.F(x, t + tao / 2) + k * uCur[J - 1] + k * uCur[J - 3] + (1 / tao - 2 * k) * uCur[J - 2] + a * sigma / (h * h) * u.Mu2()[l];

                x = h;

                al[1] = secnDiag / mainDiag;
                bt[1] = rightPart[1] / mainDiag;
                for (int i = 2; i < J - 1; ++i) {
                    al[i] = secnDiag / (mainDiag - al[i - 1] * secnDiag);
                    bt[i] = (rightPart[i] - bt[i - 1] * secnDiag) / (mainDiag - al[i - 1] * secnDiag);
                }

                uNext[J - 1] = u.Mu2()[l];
                uNext[J - 2] = bt[J - 2];
                for (int i = J - 3; i > 0; --i) {
                    uNext[i] = -al[i] * uNext[i + 1] + bt[i];
                }
                uNext[0] = u.Mu1()[l];

                t += tao;

                for (int i = 0; i < J; ++i) {
                    error[i] = u.Value(i * h, t) - uNext[i];
                    uCur[i] = uNext[i];
                    graphApp[i, l] = Convert.ToString(uNext[i]);
                }

                double curError = 0;
                for (int i = 0; i < J; ++i) {
                    error[i] = Math.Abs(error[i]);

                    if (curError < error[i]) {
                        curError = error[i];
                    }

                    //Console.WriteLine(error[i]);
                    // Console.ReadKey();
                }

                if (maxError < curError)
                    maxError = curError;

                Console.WriteLine(maxError);
                //Console.ReadKey();
            }
            Console.WriteLine(maxError);

            using (StreamWriter app = new StreamWriter(@"C:\UserDocs\HeatEquationApproximation.txt")) {
                app.WriteLine("Approximation:");
                for (int i = 0; i < J; ++i) {
                    app.WriteLine(i + ".   ");
                    for (int j = 0; j < K; ++j) {
                        app.WriteLine(graphApp[i, j]);
                    }
                    app.WriteLine();
                }

            }

            using (StreamWriter acc = new StreamWriter(@"C:\UserDocs\HeatEquationAccuraci.txt")) {
                acc.WriteLine("Accuraci:");
                for (int i = 0; i < 11; ++i) {
                    acc.WriteLine(i + ".   ");
                    for (int j = 0; j < 11; ++j) {
                        acc.WriteLine(u.graphAcc[i, j]);
                    }
                    acc.WriteLine();
                    acc.WriteLine();
                }
            }

          //  System.IO.File.WriteAllLines(@"C:\UserDocs\HeatEquation.txt", u.graphAcc);
            //System.IO.File.WriteAllLines(@"C:\UserDocs\HeatEquation.txt", graphApp);

            // Console.ReadKey()  ;
            Console.WriteLine();
            //Console.WriteLine(maxError);
            Console.ReadKey();
            Console.ReadKey();
        }
    }
}