API tools faq
paste
Advertisement
osipyonok

Untitled

osipyonok
Nov 14th, 2016
143
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C# 13.78 KB | None | 0 0
raw download clone embed print report
  1. using System;
  2. using System.Collections.Generic;
  3. using System.Linq;
  4. using System.Text;
  5. using MaximaSharp;
  6. using System.Collections;
  7.  
  8. namespace _4m_lab3
  9. {
  10.     class Program
  11.     {
  12.         static SplineTuple[] splines;
  13.  
  14.         private struct SplineTuple
  15.         {
  16.             public double a, b, c, d, x;
  17.         }
  18.         static void Main(string[] args)
  19.         {
  20.             //Maxima.GnuPlot("plot x/cos(x)");
  21.             //Maxima.GnuPlot("plot 1.0341-0.21178*(x-2)-0.52526*(x-2)**2+0.3029*(x-2)**3");
  22.             //Maxima.GnuPlot("plot 1.0341+-0.21178*(x-2)+-0.52526*(x-2)**2++0.3029*(x-2)**3");
  23.             /*Maxima.GnuPlot(@"
  24.                 f(x)=  x>=1 & x<2 ? 1.0002+0*(x-1)+0.3134*(x-1)**2+-0.2795*(x-1)**3 : x>=2 && x<3 ? 1.0341+-0.21178*(x-2)+-0.5252*(x-2)**2+0.3029*(x-2)**3 : 0
  25.                 plot [0:10] f(x)
  26.             ");*/
  27. //           Maxima.GnuPlot(@"
  28. //                f(x)=0<x<1 ? (x-1)**2 : x>7 ? x : (x-8)**2
  29. //                plot [0:10] f(x)
  30. //            ");          
  31.          /*   Maxima.GnuPlot(@"
  32.                 set parametric
  33.                 set pm3d depthorder hidden3d
  34.                 set isosamples 15, 15
  35.                 splot [-pi:pi][-pi:pi] cos(u)*(cos(v)+3), sin(u)*(cos(v)+3), sin(v) w pm
  36.             ");*/
  37.             string fx;
  38.             int a, b, n;
  39.             //fx = Console.ReadLine();
  40.             //n = Int32.Parse(Console.ReadLine());
  41.             //a = Int32.Parse(Console.ReadLine());
  42.             //b = Int32.Parse(Console.ReadLine());
  43.             n = 6;
  44.             a = 1;
  45.             b = 6;
  46.  
  47.             List<double> points = new List<double>();
  48.             //points.Add(double.MinValue);
  49.             Random r = new Random( );
  50.  
  51.             /*for (int i = 0; i < n - 2; ++i)
  52.             {
  53.                 points.Add(r.NextDouble() * (b - a) + a);
  54.             }*/
  55.             for (double i = 1.0; i <= 6.0; ++i)
  56.             {
  57.                 points.Add(i);
  58.             }
  59.             //points.Add(a);
  60.             //points.Add(b);
  61.             points.Sort();
  62.             List<double> y = new List<double>();
  63.             //for (int i = 0; i < n; ++i)
  64.             //{
  65.             //    string func = fx.Replace("x", points[i].ToString().Replace(',','.'));
  66.             //    y.Add(double.Parse(Maxima.Eval(func).Replace('.',',')));
  67.             //    //Console.WriteLine(Maxima.Eval(func).Replace('.',','));
  68.             //}
  69.             //y.Add(0);
  70.             y.Add(1.0002);
  71.             y.Add(1.0341);
  72.             y.Add(0.6);
  73.             y.Add(0.40105);
  74.             y.Add(0.1);
  75.             y.Add(0.23975);
  76.  
  77.             foreach(double u in y){
  78.                 Console.WriteLine(u);
  79.             }
  80.             Spline(points, y, n);
  81.             //BuildSpline(points, y, n);
  82.             //Cubiq(points, y, n);
  83.             string all = "f(x)=";
  84.             for(int i = 0 ; i < n - 1 ; ++i){
  85.                 SplineTuple sp = splines[i];
  86.                 //Console.WriteLine(sp.a + " " + sp.b + " " + sp.c + " " + sp.d);
  87.                 string cur = "";
  88.                 cur += sp.a.ToString();
  89.                 cur += "+" + sp.b + "*(x - " + sp.x + ")";
  90.                 cur += "+" + sp.c + "*(x - " + sp.x + ")**2";
  91.                 cur += "+" + sp.d + "*(x - " + sp.x + ")**3";
  92.                 /*cur += sp.a.ToString();
  93.                 cur += Sgn(sp.b, false);
  94.                 cur += sp.b + "*(x" + Sgn(-sp.x) + -sp.x + ")";
  95.                 cur += Sgn(sp.c, false);
  96.                 cur += "(" + sp.c + ")" + "*(x" + Sgn(-sp.x) + -sp.x + ")**2";
  97.                 cur += Sgn(sp.d, false);
  98.                 cur += "(" + sp.d + ")" + "*(x" + Sgn(-sp.x) + -sp.x + ")**3";*/
  99.            
  100.                 cur = cur.Replace(',', '.');
  101.                 if(i != 0 && i != n - 2)all += " x>=" + points[i] + " & x<" + points[i + 1] + " ?";
  102.                 if (i == 0) all += "x<" + points[i + 1] + " ?";
  103.                 //if(i == n-2) all += "x
  104.                 all += cur;
  105.                 if(i != n - 2)all += ":";
  106.                 //if (i + 1 == n - 1) all += "0";
  107.                 //Console.WriteLine(cur);
  108.                 //Maxima.GnuPlot("plot [0:10] " + cur);
  109.                 //Console.ReadKey();
  110.             }
  111.             Console.WriteLine(all);
  112.             all += "\n";
  113.             all += "plot [1:6] f(x)";
  114.             Maxima.GnuPlot(all);
  115.  
  116.             //string fplot = "plot [" + a.ToString() + ":" + b.ToString() + "] " + fx + " x+1";
  117.             //Maxima.GnuPlot(fplot);
  118.             //Maxima.Eval("plot2d([[parametric, cos(t), sin(t), [t, 0, 2*%pi]], -abs(x)], [x, -sqrt(2) , sqrt(2)] , same_xy)");
  119.  
  120.            //for (int i = 0; i < n; ++i)Console.WriteLine(points[i]);
  121.            
  122.             //while (true)
  123.             //{
  124.             //    s = Console.ReadLine();
  125.             //    Console.Write(Maxima.Eval(s) + "\n");
  126.             //}
  127.             Console.ReadKey();
  128.          //   System.
  129.            
  130.         }
  131.  
  132.         public static string Sgn(double d , bool flag = true)
  133.         {
  134.             if (d >= 0) return flag ? "-" : "+";
  135.             return "";
  136.         }
  137.  
  138.         //public static double[] Tridiagonal(double[,] A , double[] B , int n)
  139.         //{
  140.         //    double[] x = new double[n];
  141.         //    double[] alpha = new double[n - 1];
  142.         //    double[] beta = new double[n - 1];
  143.  
  144.         //    alpha[0] = -A[2,0] / A[1,0];
  145.         //    beta[0] = B[0] / A[1, 0];
  146.         //    for (int i = 1; i < n - 1; ++i)
  147.         //    {
  148.         //        alpha[i] = -A[2, i] / (A[1, i] + A[0, i] * alpha[i - 1]);
  149.         //        beta[i] = (B[i] - A[0, i] * beta[i - 1]) / (A[1, i] + A[0, i] * alpha[i - 1]);
  150.         //    }
  151.         //    x[n - 1] = (B[n - 1] - A[0, n - 1] * beta[n - 2]) / (A[1, n - 1] + A[0, n - 1] * alpha[n - 2]);
  152.         //    for (int i = n - 2; i > - 1; --i)
  153.         //    {
  154.         //        x[i] = x[i + 1] * alpha[i] + beta[i];
  155.         //        //Console.WriteLine((i + 1).ToString() + " " + i.ToString());
  156.         //    }
  157.         //    return x;
  158.         //}
  159.  
  160.         public static void Spline(List<double> x, List<double> y, int n)
  161.         {
  162.             splines = new SplineTuple[n-1];
  163.             for (int i = 0; i < n - 1; ++i)
  164.             {
  165.                 splines[i].x = x[i];
  166.                 splines[i].a = y[i];
  167.             }
  168.  
  169.             double[] alpha = new double[n-2];
  170.             double[] beta = new double[n-2];
  171.             double[] gamma = new double[n-2];
  172.             double[] delta = new double[n-2];
  173.  
  174.             for (int i = 1; i <= n - 3; ++i)
  175.             {
  176.                 double hi = x[i + 1] - x[i];
  177.                 alpha[i] = 1.0 / hi;
  178.             }
  179.             for (int i = 0; i <= n - 3; ++i)
  180.             {
  181.                 double hi = x[i + 1] - x[i];
  182.                 double hii = x[i + 2] - x[i + 1];
  183.                 beta[i] = 2.0 * ((1.0 / hi) + (1.0 / hii));
  184.             }
  185.             for (int i = 0; i <= n - 4; ++i)
  186.             {
  187.                 double hii = x[i + 2] - x[i + 1];
  188.                 gamma[i] = 1.0 / hii;
  189.             }
  190.             for (int i = 0; i <= n - 3; ++i)
  191.             {
  192.                 double hi = x[i + 1] - x[i];
  193.                 double hii = x[i + 2] - x[i + 1];
  194.                 double fi = (y[i + 1] - y[i]) / (x[i + 1] - x[i]);
  195.                 double fii = (y[i + 2] - y[i + 1]) / (x[i + 2] - x[i + 1]);
  196.                 delta[i] = 3.0 * ((fii / hii) + (fi / hi));
  197.             }
  198.             double[,] A = new double[n-2, n-2];
  199.             A[0,0] = beta[0];
  200.             A[0,1] = gamma[0];
  201.             A[n - 3, n - 3] = beta[n - 3];
  202.             A[n - 3, n - 4] = alpha[n - 3];
  203.             for (int i = 1; i <= n - 4; ++i)
  204.             {
  205.                 A[i, i - 1] = alpha[i];
  206.                 A[i, i] = beta[i];
  207.                 A[i, i + 1] = gamma[i];
  208.             }
  209.             /*for (int i = 0; i < alpha.Length; ++i)
  210.             {
  211.                 Console.Write(alpha[i] + " ");
  212.             }*/
  213.             for (int i = 0; i <= n - 3; ++i)
  214.             {
  215.                 for (int j = 0; j <= n - 3 ; ++j)
  216.                 {
  217.                     Console.Write(A[i, j] + " ");
  218.                 }
  219.                 Console.WriteLine(delta[i]);
  220.             }
  221.  
  222.             double[] tmp = Tridiagonal(A, delta , n - 2);
  223.             double[] b = new double[n];
  224.             b[0] = 0;
  225.             b[n - 1] = 0;
  226.             for (int i = 1; i <= n - 2; ++i)
  227.             {
  228.                 b[i] = tmp[i - 1];
  229.             }
  230.            
  231.             double[] c = new double[n - 1];
  232.             for(int i = 0 ; i < n - 1 ; ++i){
  233.                 double fi = (y[i + 1] - y[i]) / (x[i + 1] - x[i]);
  234.                 double hi = x[i + 1] - x[i];
  235.                 c[i] = (3 * fi - b[i+1] - 2 * b[i])/hi;
  236.             }
  237.             double[] d = new double[n-1];
  238.             for (int i = 0; i < n - 1; ++i)
  239.             {
  240.                 double fi = (y[i + 1] - y[i]) / (x[i + 1] - x[i]);
  241.                 double hi = x[i + 1] - x[i];
  242.                 d[i] = (b[i] + b[i + 1] - 2 * fi) / Math.Pow(hi, 2);
  243.             }
  244.             for (int i = 0; i < n - 1; ++i)
  245.             {
  246.                 splines[i].a = y[i];
  247.                 splines[i].b = b[i];
  248.                 splines[i].c = c[i];
  249.                 splines[i].d = d[i];
  250.                 //Console.WriteLine(y[i] + " " + b[i] + " " + c[i] + " " + d[i] + " ");
  251.             }
  252.            
  253.         }
  254.  
  255.         public static double[] Tridiagonal(double[,] A, double[] b, int n)
  256.         {
  257.             double[] x = new double[n];
  258.             double[] alpha = new double[n+10];
  259.             double[] beta = new double[n+10];
  260.             alpha[1] = -A[0,1] / A[0,0];
  261.             beta[1] = b[0] / A[0,0];
  262.             for (int i = 1; i <= n - 2; ++i)
  263.             {
  264.                 alpha[i + 1] = -A[i, i + 1] / (A[i, i - 1] * alpha[i] + A[i, i]);
  265.                 beta[i + 1] = (b[i] - (A[i, i - 1] * beta[i])) / (A[i, i - 1] * alpha[i] + A[i, i]);
  266.             }
  267.             x[n - 1] = (b[n - 1] - (A[n - 1,n - 2] * beta[n - 1])) / (A[n - 1,n - 1] + (A[n - 1,n - 2] * alpha[n - 1]));
  268.             for (int i = n - 2; i >= 0; --i)
  269.             {
  270.                 x[i] = (alpha[i + 1] * x[i + 1]) + beta[i + 1];
  271.             }
  272.             return x;
  273.         }
  274.  
  275.         //public static void BuildSpline(List<double> x, List<double> y, int n)
  276.         //{
  277.         //    splines = new SplineTuple[n];
  278.         //    for (int i = 0; i < n ; ++i)
  279.         //    {
  280.         //        splines[i].x = x[i];
  281.         //        splines[i].a = y[i];
  282.         //    }
  283.         //    double[] a = new double[n - 1];
  284.         //    double[] b = new double[n];
  285.         //    double[] c = new double[n - 1];
  286.         //    double[] d = new double[n - 1];
  287.         //    double[] delta = new double[n - 1];
  288.         //    double[] h = new double[n - 1];
  289.         //    double[,] A = new double[3, n];
  290.         //    double[] f = new double[n];
  291.         //    double x3, xn;
  292.         //    if (n < 3) return;
  293.         //    x3 = splines[2].x - splines[0].x;
  294.         //    xn = splines[n - 1].x - splines[n - 3].x;
  295.  
  296.         //    for (int i = 0; i < n - 1; ++i)
  297.         //    {
  298.         //        a[i] = splines[i].a;
  299.         //        h[i] = splines[i + 1].x - splines[i].x;
  300.         //        delta[i] = (splines[i + 1].a - splines[i].a) / h[i];
  301.         //        A[0, i] = i > 0 ? h[i] : x3;
  302.         //        f[i] = i > 0 ? 3 * (h[i] * delta[i - 1] + h[i - 1] * delta[i]) : 0;
  303.         //    }
  304.         //    A[1, 0] = h[0];
  305.         //    A[2, 0] = h[0];
  306.         //    int k = n - 1;
  307.         //    for (int i = 1; i < n - 1; ++i)
  308.         //    {
  309.         //        A[1, i] = 2 * (h[i] + h[i - 1]);
  310.         //        A[2, i] = h[i];
  311.         //    }
  312.         //    A[1, n - 1] = h[n - 2];
  313.         //    A[2, n - 1] = xn;
  314.         //    A[0, n - 1] = h[n - 2];
  315.  
  316.         //    f[0] = ((h[0] + 2 * x3) * h[1] * delta[0] + Math.Pow(h[0], 2) * delta[1]) / x3;
  317.         //    f[n - 1] = (Math.Pow(h[n - 2], 2) * delta[n - 3] + (2 * xn + h[n - 2]) * h[n - 3] * delta[n - 2]) / xn;
  318.  
  319.         //    b = Tridiagonal(A, f,n);
  320.         //    for (int i = 0; i < n - 1; ++i)
  321.         //    {
  322.         //        d[i] = (b[i + 1] + b[i] - 2 * delta[i]) / (h[i] * h[i]);
  323.         //        c[i] = 2 * (delta[i] - b[i]) / h[i] - (b[i + 1] - delta[i]) / h[i];
  324.  
  325.         //        splines[i].a = a[i];
  326.         //        splines[i].b = b[i];
  327.         //        splines[i].c = c[i];
  328.         //        splines[i].d = d[i];
  329.         //    }
  330.         //}
  331.         public static void Cubiq(List<double> x, List<double> y, int n)
  332.         {
  333.             splines = new SplineTuple[n];
  334.             for (int i = 0; i < n; ++i)
  335.             {
  336.                 splines[i].x = x[i];
  337.                 splines[i].a = y[i];
  338.             }
  339.             splines[0].c = 0.0;
  340.             splines[n - 1].c = 0.0;
  341.  
  342.             double[] alpha = new double[n - 1];
  343.             double[] beta = new double[n - 1];
  344.             alpha[0] = 0.0;
  345.             beta[0] = 0.0;
  346.             for (int i = 1; i < n - 1; ++i)
  347.             {
  348.                 double hi = x[i] - x[i - 1];
  349.                 double hi1 = x[i + 1] - x[i];
  350.                 double A = hi;
  351.                 double C = 2.0 * (hi + hi1);
  352.                 double B = hi1;
  353.                 double F = 6.0 * ((y[i + 1] - y[i]) / hi1 - (y[i] - y[i - 1]) / hi);
  354.                 double z = (A * alpha[i - 1] + C);
  355.                 alpha[i] = -B / z;
  356.                 beta[i] = (F - A * beta[i - 1]) / z;
  357.             }
  358.             for (int i = n - 2; i > 0; --i)
  359.             {
  360.                 splines[i].c = alpha[i] * splines[i + 1].c + beta[i];
  361.             }
  362.             for (int i = n - 1; i > 0; --i)
  363.             {
  364.                 double hi = x[i] - x[i - 1];
  365.                 splines[i].d = (splines[i].c - splines[i - 1].c) / hi;
  366.                 splines[i].b = hi * (2.0 * splines[i].c + splines[i - 1].c) / 6.0 + (y[i] - y[i - 1]) / hi;
  367.             }
  368.  
  369.         }
  370.     }
  371. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!