С Новым годом! Форум программистов, компьютерный форум, киберфорум
Python для начинающих
Войти
Регистрация
Восстановить пароль
Блоги Сообщество Поиск Заказать работу  
 
Рейтинг 4.86/222: Рейтинг темы: голосов - 222, средняя оценка - 4.86
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577

Кубический сплайн

04.05.2019, 22:49. Показов 48136. Ответов 11

Студворк — интернет-сервис помощи студентам
Вечер добрый. Поступило такое задание. Необходимо реализовать при использовании языка Python интерполяцию эмпирических данных кубическим сплайном с естественными граничными условиями. Есть у кого-то мысли от чего нужно отталкиваться? (Сухих теоретических данных маловато для непосредственной реализации). Спасибо.
P.S. График рекомендуют строить с помощью mathplotlib. Насколько это будет удобно и рационально?
0
Лучшие ответы (1)
Programming
Эксперт
39485 / 9562 / 3019
Регистрация: 12.04.2006
Сообщений: 41,671
Блог
04.05.2019, 22:49
Ответы с готовыми решениями:

Кубический корень от х
Подскажите, пожалуйста, как рассчитать кубический корень от х. pow не работает с 1/3. а math.exp(math.log(x)/3) не работает с числами...

График функций в виде буквы G. Сплайн Эрмита,Кубический сплайн Эрмита
Нужно построить график в виде буквы G.Гладкая кривая + 2 прямых.+ Написать к этим функция сплайны Эрмита и многочлен Лагранжа.В написание...

Кубический сплайн
В лекции на сайте нашел вот такой сплайн для коэф.Добеши Не получается построить график полученных кривых. Не понимаю второй аргумент:...

11
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
05.05.2019, 13:10  [ТС]
Нашел подобную реализацию на С++. Может кто-то помочь переделать под Python?
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
#include <cstdlib>
#include <cmath>
#include <limits>
 
class cubic_spline
{
private:
        // Структура, описывающая сплайн на каждом сегменте сетки
        struct spline_tuple
        {
                double a, b, c, d, x;
        };
 
        spline_tuple *splines; // Сплайн
        std::size_t n; // Количество узлов сетки
 
        void free_mem(); // Освобождение памяти
 
public:
        cubic_spline(); //конструктор
        ~cubic_spline(); //деструктор
 
        // Построение сплайна
        // x - узлы сетки, должны быть упорядочены по возрастанию, кратные узлы запрещены
        // y - значения функции в узлах сетки
        // n - количество узлов сетки
        void build_spline(const double *x, const double *y, std::size_t n);
 
        // Вычисление значения интерполированной функции в произвольной точке
        double f(double x) const;
};
 
cubic_spline::cubic_spline() : splines(NULL)
{
 
}
 
cubic_spline::~cubic_spline()
{
        free_mem();
}
 
void cubic_spline::build_spline(const double *x, const double *y, std::size_t n)
{
        free_mem();
 
        this->n = n;
 
        // Инициализация массива сплайнов
        splines = new spline_tuple[n];
        for (std::size_t i = 0; i < n; ++i)
        {
                splines[i].x = x[i];
                splines[i].a = y[i];
        }
        splines[0].c = 0.;
 
        // Решение СЛАУ относительно коэффициентов сплайнов c[i] методом прогонки для трехдиагональных матриц
        // Вычисление прогоночных коэффициентов - прямой ход метода прогонки
        double *alpha = new double[n - 1];
        double *beta = new double[n - 1];
        double A, B, C, F, h_i, z;
        alpha[0] = beta[0] = 0.;
        for (std::size_t i = 1; i < n - 1; ++i)
        {
                h_i = x[i] - x[i - 1], h_i1 = x[i + 1] - x[i];
                A = h_i;
                C = 2. * (h_i + h_i1);
                B = h_i1;
                F = 6. * ((y[i + 1] - y[i]) / h_i1 - (y[i] - y[i - 1]) / h_i);
                z = (A * alpha[i - 1] + C);
                alpha[i] = -B / z;
                beta[i] = (F - A * beta[i - 1]) / z;
        }
 
        splines[n - 1].c = (F - A * beta[n - 2]) / (C + A * alpha[n - 2]);
 
        // Нахождение решения - обратный ход метода прогонки
        for (std::size_t i = n - 2; i > 0; --i)
                splines[i].c = alpha[i] * splines[i + 1].c + beta[i];
 
        // Освобождение памяти, занимаемой прогоночными коэффициентами
        delete[] beta;
        delete[] alpha;
 
        // По известным коэффициентам c[i] находим значения b[i] и d[i]
        for (std::size_t i = n - 1; i > 0; --i)
        {
                double h_i = x[i] - x[i - 1];
                splines[i].d = (splines[i].c - splines[i - 1].c) / h_i;
                splines[i].b = h_i * (2. * splines[i].c + splines[i - 1].c) / 6. + (y[i] - y[i - 1]) / h_i;
        }
}
 
double cubic_spline::f(double x) const
{
        if (!splines)
                return std::numeric_limits<double>::quiet_NaN(); // Если сплайны ещё не построены - возвращаем NaN
 
        spline_tuple *s;
        if (x <= splines[0].x) // Если x меньше точки сетки x[0] - пользуемся первым эл-тов массива
                s = splines + 1;
        else if (x >= splines[n - 1].x) // Если x больше точки сетки x[n - 1] - пользуемся последним эл-том массива
                s = splines + n - 1;
        else // Иначе x лежит между граничными точками сетки - производим бинарный поиск нужного эл-та массива
        {
                std::size_t i = 0, j = n - 1;
                while (i + 1 < j)
                {
                        std::size_t k = i + (j - i) / 2;
                        if (x <= splines[k].x)
                                j = k;
                        else
                                i = k;
                }
                s = splines + j;
        }
 
        double dx = (x - s->x);
        return s->a + (s->b + (s->c / 2. + s->d * dx / 6.) * dx) * dx; // Вычисляем значение сплайна в заданной точке по схеме Горнера (в принципе, "умный" компилятор применил бы схему Горнера сам, но ведь не все так умны, как кажутся)
}
 
void cubic_spline::free_mem()
{
        delete[] splines;
        splines = NULL;
}
0
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
06.05.2019, 21:30  [ТС]
Вот еще на Python нашел. Но не совсем понятно является ли этот кубический сплайн именно с естественными граничными условиями.

Python
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
from collections import defaultdict
from matplotlib import mlab
import bisect, pylab, math
 
def drange(start, stop, step):
    while start < stop: yield start; start += step;
 
class Dot:
    def __init__(self, x, y): self.x, self.y = [x, y]
 
class Tuple: a, b, c, d, x = [0., 0., 0., 0., 0.]
 
def buildSpline(dots):
    for i in range(len(dots)): splines[i].x, splines[i].a = dots[i].x, dots[i].y
 
    alpha, beta = [defaultdict(lambda: 0.), defaultdict(lambda: 0.)]
 
    for i in range(1, len(dots)-1):
        C = 4. * in_step
        F = 6. * ((dots[i + 1].y - dots[i].y) / in_step - (dots[i].y - dots[i - 1].y) / in_step)
        z = (in_step* alpha[i - 1] + C)
        alpha[i] = -in_step / z
        beta[i] = (F - in_step* beta[i - 1]) / z
 
    for i in reversed(range(1, len(dots) - 1)): splines[i].c = alpha[i] * splines[i+1].c + beta[i]
 
    for i in reversed(range(1, len(dots))):
        hi = dots[i].x - dots[i-1].x
        splines[i].d = (splines[i].c - splines[i-1].c) / hi
        splines[i].b = hi * (2.0 * splines[i].c + splines[i - 1].c) / 6.0 + (dots[i].y - dots[i-1].y) / hi
 
def calc(x):
    distribution = sorted([t[1].x for t in splines.items()])
    indx = bisect.bisect_left(distribution, x)
    if indx == len(distribution): return 0
    dx = x - splines[indx].x
    return splines[indx].a + splines[indx].b * dx + splines[indx].c * dx**2 / 2. + splines[indx].d * dx**3 / 6.
#============================================
 
in_func =  lambda x:  math.sin(x)
in_min_x = 0
in_max_x = 25 
in_step = 2.5
 
out_min_x = 0
out_max_x = 25
out_step = 0.1
 
#============================================
 
#build model
splines = defaultdict(lambda: Tuple())
buildSpline(map(lambda x: Dot(x, in_func(x)), [x for x in drange(in_min_x, in_max_x+1, in_step)]))
 
#print result
for x in drange(out_min_x, out_max_x, out_step):
    print str(x) + ';' + str(calc(x))
 
#build graphics
xlist = mlab.frange (out_min_x, out_max_x, out_step)
pylab.plot(xlist, [calc(x) for x in xlist])
pylab.plot(xlist, [in_func(x) for x in xlist])
pylab.show()
0
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
09.05.2019, 00:34  [ТС]
Может с С# кто-то сможет?
C#
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
// Интерполирование функций естественными кубическими сплайнами
 
using System;
 
class CubicSpline
{
    SplineTuple[] splines; // Сплайн
    
    // Структура, описывающая сплайн на каждом сегменте сетки
    private struct SplineTuple
    {
        public double a, b, c, d, x;
    }
 
    // Построение сплайна
    // x - узлы сетки, должны быть упорядочены по возрастанию, кратные узлы запрещены
    // y - значения функции в узлах сетки
    // n - количество узлов сетки
    public void BuildSpline(double[] x, double[] y, int n)
    {
        // Инициализация массива сплайнов
        splines = new SplineTuple[n];
        for (int i = 0; i < n; ++i)
        {
            splines[i].x = x[i];
            splines[i].a = y[i];
        }
        splines[0].c = splines[n - 1].c = 0.0;
        
        // Решение СЛАУ относительно коэффициентов сплайнов c[i] методом прогонки для трехдиагональных матриц
        // Вычисление прогоночных коэффициентов - прямой ход метода прогонки
        double[] alpha = new double[n - 1];
        double[] beta  = new double[n - 1];
        alpha[0] = beta[0] = 0.0;
        for (int i = 1; i < n - 1; ++i)
        {
            double hi  = x[i] - x[i - 1];
            double hi1 = x[i + 1] - x[i];
            double A = hi;
            double C = 2.0 * (hi + hi1);
            double B = hi1;
            double F = 6.0 * ((y[i + 1] - y[i]) / hi1 - (y[i] - y[i - 1]) / hi);
            double z = (A * alpha[i - 1] + C);
            alpha[i] = -B / z;
            beta[i] = (F - A * beta[i - 1]) / z;
        }
 
        // Нахождение решения - обратный ход метода прогонки
        for (int i = n - 2; i > 0; --i)
        {
            splines[i].c = alpha[i] * splines[i + 1].c + beta[i];
        }
        
        // По известным коэффициентам c[i] находим значения b[i] и d[i]
        for (int i = n - 1; i > 0; --i)
        {
            double hi = x[i] - x[i - 1];
            splines[i].d = (splines[i].c - splines[i - 1].c) / hi;
            splines[i].b = hi * (2.0 * splines[i].c + splines[i - 1].c) / 6.0 + (y[i] - y[i - 1]) / hi;
        }
    }
 
    // Вычисление значения интерполированной функции в произвольной точке
    public double Interpolate(double x)
    {
        if (splines == null)
        {
            return double.NaN; // Если сплайны ещё не построены - возвращаем NaN
        }
        
        int n = splines.Length;
        SplineTuple s;
        
        if (x <= splines[0].x) // Если x меньше точки сетки x[0] - пользуемся первым эл-тов массива
        {
            s = splines[0];
        }
        else if (x >= splines[n - 1].x) // Если x больше точки сетки x[n - 1] - пользуемся последним эл-том массива
        {
            s = splines[n - 1];
        }
        else // Иначе x лежит между граничными точками сетки - производим бинарный поиск нужного эл-та массива
        {
            int i = 0;
            int j = n - 1;
            while (i + 1 < j)
            {
                int k = i + (j - i) / 2;
                if (x <= splines[k].x)
                {
                    j = k;
                }
                else
                {
                    i = k;
                }
            }
            s = splines[j];
        }
        
        double dx = x - s.x;
        // Вычисляем значение сплайна в заданной точке по схеме Горнера (в принципе, "умный" компилятор применил бы схему Горнера сам, но ведь не все так умны, как кажутся)
        return s.a + (s.b + (s.c / 2.0 + s.d * dx / 6.0) * dx) * dx; 
    }
}
0
1293 / 677 / 367
Регистрация: 07.01.2019
Сообщений: 2,300
09.05.2019, 04:42
За 5 дней можно было уже самому разобраться

Python
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
# Структура, описывающая сплайн на каждом сегменте сетки
class SplineTuple:
    def __init__(self, a, b, c, d, x):
        self.a = a
        self.b = b
        self.c = c
        self.d = d
        self.x = x
 
# Построение сплайна
# x - узлы сетки, должны быть упорядочены по возрастанию, кратные узлы запрещены
# y - значения функции в узлах сетки
# n - количество узлов сетки
def BuildSpline(x, y, n):
    # Инициализация массива сплайнов
    splines = [SplineTuple(0, 0, 0, 0, 0) for _ in range(0, n)]
    for i in range(0, n):
        splines[i].x = x[i]
        splines[i].a = y[i]
    
    splines[0].c = splines[n - 1].c = 0.0
    
    # Решение СЛАУ относительно коэффициентов сплайнов c[i] методом прогонки для трехдиагональных матриц
    # Вычисление прогоночных коэффициентов - прямой ход метода прогонки
    alpha = [0.0 for _ in range(0, n - 1)]
    beta  = [0.0 for _ in range(0, n - 1)]
 
    for i in range(1, n - 1):
        hi  = x[i] - x[i - 1]
        hi1 = x[i + 1] - x[i]
        A = hi
        C = 2.0 * (hi + hi1)
        B = hi1
        F = 6.0 * ((y[i + 1] - y[i]) / hi1 - (y[i] - y[i - 1]) / hi)
        z = (A * alpha[i - 1] + C)
        alpha[i] = -B / z
        beta[i] = (F - A * beta[i - 1]) / z
  
 
    # Нахождение решения - обратный ход метода прогонки
    for i in range(n - 2, 0, -1):
        splines[i].c = alpha[i] * splines[i + 1].c + beta[i]
    
    # По известным коэффициентам c[i] находим значения b[i] и d[i]
    for i in range(n - 1, 0, -1):
        hi = x[i] - x[i - 1]
        splines[i].d = (splines[i].c - splines[i - 1].c) / hi
        splines[i].b = hi * (2.0 * splines[i].c + splines[i - 1].c) / 6.0 + (y[i] - y[i - 1]) / hi
    return splines
 
 
# Вычисление значения интерполированной функции в произвольной точке
def Interpolate(splines, x):
    if not splines:
        return None # Если сплайны ещё не построены - возвращаем NaN
    
    n = len(splines)
    s = SplineTuple(0, 0, 0, 0, 0)
    
    if x <= splines[0].x: # Если x меньше точки сетки x[0] - пользуемся первым эл-тов массива
        s = splines[0]
    elif x >= splines[n - 1].x: # Если x больше точки сетки x[n - 1] - пользуемся последним эл-том массива
        s = splines[n - 1]
    else: # Иначе x лежит между граничными точками сетки - производим бинарный поиск нужного эл-та массива
        i = 0
        j = n - 1
        while i + 1 < j:
            k = i + (j - i) // 2
            if x <= splines[k].x:
                j = k
            else:
                i = k
        s = splines[j]
    
    dx = x - s.x
    # Вычисляем значение сплайна в заданной точке по схеме Горнера (в принципе, "умный" компилятор применил бы схему Горнера сам, но ведь не все так умны, как кажутся)
    return s.a + (s.b + (s.c / 2.0 + s.d * dx / 6.0) * dx) * dx;
    
spline = BuildSpline([1, 3, 7, 9], [5, 6, 7, 8], 4)
print(Interpolate(spline, 5))
1
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
09.05.2019, 13:30  [ТС]
tooru, подскажите кто-нибудь как можно это отобразить графически? Как вариант, в mathplotlib.
0
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
09.05.2019, 23:17  [ТС]
tooru, Нашел код, который реализует сплайн и отображается графически. Подскажите, пожалуйста кто-нибудь, является ли он сплайном с естественными граничными условиями.
Python
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
import numpy as np
from math import sqrt
 
def cubic_interp1d(x0, x, y):
    """
    Interpolate a 1-D function using cubic splines.
      x0 : a float or an 1d-array
      x : (N,) array_like
          A 1-D array of real/complex values.
      y : (N,) array_like
          A 1-D array of real values. The length of y along the
          interpolation axis must be equal to the length of x.
 
    Implement a trick to generate at first step the cholesky matrice L of
    the tridiagonal matrice A (thus L is a bidiagonal matrice that
    can be solved in two distinct loops).
 
    additional ref: www.math.uh.edu/~jingqiu/math4364/spline.pdf 
    """
    x = np.asfarray(x)
    y = np.asfarray(y)
 
    # remove non finite values
    # indexes = np.isfinite(x)
    # x = x[indexes]
    # y = y[indexes]
 
    # check if sorted
    if np.any(np.diff(x) < 0):
        indexes = np.argsort(x)
        x = x[indexes]
        y = y[indexes]
 
    size = len(x)
 
    xdiff = np.diff(x)
    ydiff = np.diff(y)
 
    # allocate buffer matrices
    Li = np.empty(size)
    Li_1 = np.empty(size-1)
    z = np.empty(size)
 
    # fill diagonals Li and Li-1 and solve [L][y] = [B]
    Li[0] = sqrt(2*xdiff[0])
    Li_1[0] = 0.0
    B0 = 0.0 # natural boundary
    z[0] = B0 / Li[0]
 
    for i in range(1, size-1, 1):
        Li_1[i] = xdiff[i-1] / Li[i-1]
        Li[i] = sqrt(2*(xdiff[i-1]+xdiff[i]) - Li_1[i-1] * Li_1[i-1])
        Bi = 6*(ydiff[i]/xdiff[i] - ydiff[i-1]/xdiff[i-1])
        z[i] = (Bi - Li_1[i-1]*z[i-1])/Li[i]
 
    i = size - 1
    Li_1[i-1] = xdiff[-1] / Li[i-1]
    Li[i] = sqrt(2*xdiff[-1] - Li_1[i-1] * Li_1[i-1])
    Bi = 0.0 # natural boundary
    z[i] = (Bi - Li_1[i-1]*z[i-1])/Li[i]
 
    # solve [L.T][x] = [y]
    i = size-1
    z[i] = z[i] / Li[i]
    for i in range(size-2, -1, -1):
        z[i] = (z[i] - Li_1[i-1]*z[i+1])/Li[i]
 
    # find index
    index = x.searchsorted(x0)
    np.clip(index, 1, size-1, index)
 
    xi1, xi0 = x[index], x[index-1]
    yi1, yi0 = y[index], y[index-1]
    zi1, zi0 = z[index], z[index-1]
    hi1 = xi1 - xi0
 
    # calculate cubic
    f0 = zi0/(6*hi1)*(xi1-x0)**3 + \
         zi1/(6*hi1)*(x0-xi0)**3 + \
         (yi1/hi1 - zi1*hi1/6)*(x0-xi0) + \
         (yi0/hi1 - zi0*hi1/6)*(xi1-x0)
    return f0
 
if __name__ == '__main__':
    import matplotlib.pyplot as plt
    x = np.linspace(0, 10, 11)
    y = np.sin(x)
    plt.scatter(x, y)
 
    x_new = np.linspace(0, 10, 201)
    plt.plot(x_new, cubic_interp1d(x_new, x, y))
 
    plt.show()
0
1293 / 677 / 367
Регистрация: 07.01.2019
Сообщений: 2,300
10.05.2019, 05:37
Лучший ответ Сообщение было отмечено Teylor как решение

Решение

Python
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
import matplotlib.pyplot as plt
 
# Структура, описывающая сплайн на каждом сегменте сетки
class SplineTuple:
    def __init__(self, a, b, c, d, x):
        self.a = a
        self.b = b
        self.c = c
        self.d = d
        self.x = x
 
# Построение сплайна
# x - узлы сетки, должны быть упорядочены по возрастанию, кратные узлы запрещены
# y - значения функции в узлах сетки
# n - количество узлов сетки
def BuildSpline(x, y, n):
    # Инициализация массива сплайнов
    splines = [SplineTuple(0, 0, 0, 0, 0) for _ in range(0, n)]
    for i in range(0, n):
        splines[i].x = x[i]
        splines[i].a = y[i]
    
    splines[0].c = splines[n - 1].c = 0.0
    
    # Решение СЛАУ относительно коэффициентов сплайнов c[i] методом прогонки для трехдиагональных матриц
    # Вычисление прогоночных коэффициентов - прямой ход метода прогонки
    alpha = [0.0 for _ in range(0, n - 1)]
    beta  = [0.0 for _ in range(0, n - 1)]
 
    for i in range(1, n - 1):
        hi  = x[i] - x[i - 1]
        hi1 = x[i + 1] - x[i]
        A = hi
        C = 2.0 * (hi + hi1)
        B = hi1
        F = 6.0 * ((y[i + 1] - y[i]) / hi1 - (y[i] - y[i - 1]) / hi)
        z = (A * alpha[i - 1] + C)
        alpha[i] = -B / z
        beta[i] = (F - A * beta[i - 1]) / z
  
 
    # Нахождение решения - обратный ход метода прогонки
    for i in range(n - 2, 0, -1):
        splines[i].c = alpha[i] * splines[i + 1].c + beta[i]
    
    # По известным коэффициентам c[i] находим значения b[i] и d[i]
    for i in range(n - 1, 0, -1):
        hi = x[i] - x[i - 1]
        splines[i].d = (splines[i].c - splines[i - 1].c) / hi
        splines[i].b = hi * (2.0 * splines[i].c + splines[i - 1].c) / 6.0 + (y[i] - y[i - 1]) / hi
    return splines
 
 
# Вычисление значения интерполированной функции в произвольной точке
def Interpolate(splines, x):
    if not splines:
        return None # Если сплайны ещё не построены - возвращаем NaN
    
    n = len(splines)
    s = SplineTuple(0, 0, 0, 0, 0)
    
    if x <= splines[0].x: # Если x меньше точки сетки x[0] - пользуемся первым эл-тов массива
        s = splines[0]
    elif x >= splines[n - 1].x: # Если x больше точки сетки x[n - 1] - пользуемся последним эл-том массива
        s = splines[n - 1]
    else: # Иначе x лежит между граничными точками сетки - производим бинарный поиск нужного эл-та массива
        i = 0
        j = n - 1
        while i + 1 < j:
            k = i + (j - i) // 2
            if x <= splines[k].x:
                j = k
            else:
                i = k
        s = splines[j]
    
    dx = x - s.x
    # Вычисляем значение сплайна в заданной точке по схеме Горнера (в принципе, "умный" компилятор применил бы схему Горнера сам, но ведь не все так умны, как кажутся)
    return s.a + (s.b + (s.c / 2.0 + s.d * dx / 6.0) * dx) * dx;
    
 
x = [1, 3, 7, 9]
y = [5, 6, 7, 8]
new_x = 5
 
spline = BuildSpline(x, y, len(x))
 
plt.scatter(x, y)
plt.plot(x, y)
plt.scatter(new_x, Interpolate(spline, new_x))
plt.show()
0
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
10.05.2019, 16:06  [ТС]
tooru, благодарю!
0
14 / 14 / 0
Регистрация: 01.12.2017
Сообщений: 577
10.05.2019, 23:52  [ТС]
tooru, а в чем вообще суть сего задания (интерполяции), если график строится по заданным точкам? Понятно, что есть условия, когда точки не проходят проверку. Но в целом - просто по заданным точкам строится график.
0
1293 / 677 / 367
Регистрация: 07.01.2019
Сообщений: 2,300
11.05.2019, 08:43
Teylor, Допустим нам заданы точки [1, 3, 4, 5] и значения в этих точках [10, 30, 40, 50], какое значение имеет точка 2? Вот это мы и вычисляем, потом строим линию по известным точкам, а значение в требуемой точке наносим на график другим цветом, еще момент, линия это результат линейной интерполяции по точкам, а вычисленная точка результат расчета кубическим сплайном, видно, что центр точки немного не совпадает с линией
0
Эксперт Python
 Аватар для dondublon
4652 / 2072 / 366
Регистрация: 17.03.2012
Сообщений: 10,181
Записей в блоге: 6
11.05.2019, 10:05
Вообще кубическая интерполяция есть в scipy готовая.
0
Надоела реклама? Зарегистрируйтесь и она исчезнет полностью.
inter-admin
Эксперт
29715 / 6470 / 2152
Регистрация: 06.03.2009
Сообщений: 28,500
Блог
11.05.2019, 10:05
Помогаю со студенческими работами здесь

Кубический сплайн
Помогите понять теорию. Сколько должно быть строк в трехдиагональной матрице? К примеру,у меня размер сплайна 6, тогда интервалов же 5, а...

Кубический сплайн
Есть простейшая прога для построения кубического сплайна #include &lt;iostream&gt; #include &lt;math.h&gt; using namespace std; int...

Кубический сплайн
Реализовал кубический сплайн. Кому интересно вот файл MathCAD:

Кубический сплайн
Всем доброго времени суток. Дали задание: написать програму в MS Visual C++ для построения кубического сплайна на основе данных в таблице....

Кубический сплайн
Для введённого значения аргумента х и заданной табличной зависимости у = f(x), решить задачу интерполяции заданным методом. Предусмотреть...


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
12
Ответ Создать тему
Новые блоги и статьи
Учёным и волонтёрам проекта «Einstein@home» удалось обнаружить четыре гамма-лучевых пульсара в джете Млечного Пути
Programma_Boinc 01.01.2026
Учёным и волонтёрам проекта «Einstein@home» удалось обнаружить четыре гамма-лучевых пульсара в джете Млечного Пути Сочетание глобально распределённой вычислительной мощности и инновационных. . .
Советы по крайней бережливости. Внимание, это ОЧЕНЬ длинный пост.
Programma_Boinc 28.12.2025
Советы по крайней бережливости. Внимание, это ОЧЕНЬ длинный пост. Налог на собак: https:/ / **********/ gallery/ V06K53e Финансовый отчет в Excel: https:/ / **********/ gallery/ bKBkQFf Пост отсюда. . .
Кто-нибудь знает, где можно бесплатно получить настольный компьютер или ноутбук? США.
Programma_Boinc 26.12.2025
Нашел на реддите интересную статью под названием Anyone know where to get a free Desktop or Laptop? Ниже её машинный перевод. После долгих разбирательств я наконец-то вернула себе. . .
Thinkpad X220 Tablet — это лучший бюджетный ноутбук для учёбы, точка.
Programma_Boinc 23.12.2025
Рецензия / Мнение/ Перевод Нашел на реддите интересную статью под названием The Thinkpad X220 Tablet is the best budget school laptop period . Ниже её машинный перевод. Thinkpad X220 Tablet —. . .
PhpStorm 2025.3: WSL Terminal всегда стартует в ~
and_y87 14.12.2025
PhpStorm 2025. 3: WSL Terminal всегда стартует в ~ (home), игнорируя директорию проекта Симптом: После обновления до PhpStorm 2025. 3 встроенный терминал WSL открывается в домашней директории. . .
Как объединить две одинаковые БД Access с разными данными
VikBal 11.12.2025
Помогите пожалуйста !! Как объединить 2 одинаковые БД Access с разными данными.
Новый ноутбук
volvo 07.12.2025
Всем привет. По скидке в "черную пятницу" взял себе новый ноутбук Lenovo ThinkBook 16 G7 на Амазоне: Ryzen 5 7533HS 64 Gb DDR5 1Tb NVMe 16" Full HD Display Win11 Pro
Музыка, написанная Искусственным Интеллектом
volvo 04.12.2025
Всем привет. Некоторое время назад меня заинтересовало, что уже умеет ИИ в плане написания музыки для песен, и, собственно, исполнения этих самых песен. Стихов у нас много, уже вышли 4 книги, еще 3. . .
От async/await к виртуальным потокам в Python
IndentationError 23.11.2025
Армин Ронахер поставил под сомнение async/ await. Создатель Flask заявляет: цветные функции - провал, виртуальные потоки - решение. Не threading-динозавры, а новое поколение лёгких потоков. Откат?. . .
Поиск "дружественных имён" СОМ портов
Argus19 22.11.2025
Поиск "дружественных имён" СОМ портов На странице: https:/ / norseev. ru/ 2018/ 01/ 04/ comportlist_windows/ нашёл схожую тему. Там приведён код на С++, который показывает только имена СОМ портов, типа,. . .
КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin
Copyright ©2000 - 2026, CyberForum.ru