Ayuda: Graficar splines cúbicos
Publicado por Fernando Montiel (1 intervención) el 25/04/2019 02:47:52
Mi programa genera la matriz de los datos y obtiene los coeficientes más no he logrado tener un gráfica de los polinomios generados. Les comparto el código.
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
%Datos
x = [-1; 1; 2; 5];
y = [-1; 1; 5; 3];
%%%%%%%%% 3.- Splines cúbicos %%%%%%%%%%%%
N = length(y); %Número de datos
Num_eq = N - 1; %Número de ecuaciones
coef = 4 * Num_eq; %Coeficientes
Num_filas = 2 * (N - 2) + 2;
BY = zeros(2*(N-2)+2 + (N-2) + (N-2) + 2,1); %Vector B
H0 = zeros(Num_filas, coef); % H0 -> Matriz de las ec.
%Orden por grado x^3, x^2, x ...
for idx = 1:(Num_filas/2)
col = (idx-1)*4+1;
fila = (idx-1)*2+1;
for jdx = col:(col+3)
H0(fila,jdx) = x(idx)^(3 - jdx + col);
H0(fila+1,jdx) = x(idx+1)^(3 - jdx + col);
end
BY(fila,1) = y(idx);
BY(fila+1,1) = y(idx+1);
end
H1 = zeros(N - 2, coef); %H1 -> Matriz primera derivada
for idx = 2:(N - 1)
col = 1 + (idx - 2)*4;
H1(idx - 1, col) = 3*x(idx)^2;
H1(idx - 1, col + 1) = 2*x(idx);
H1(idx - 1, col + 2) = 1;
H1(idx - 1, col + 3) = 0;
H1(idx - 1, col + 4) = -3*x(idx)^2;
H1(idx - 1, col + 5) = -2*x(idx);
H1(idx - 1, col + 6) = -1;
end
H2 = 0*H1; %H2 -> Matriz segunda derivada
for idx = 2:(N - 1)
col = 1 + (idx - 2)*4;
H2(idx - 1, col) = 6*x(idx);
H2(idx - 1, col + 1) = 2;
H2(idx - 1, col + 2) = 0;
H2(idx - 1, col + 3) = 0;
H2(idx - 1, col + 4) = -6*x(idx);
H2(idx - 1, col + 5) = -2;
H2(idx - 1, col + 6) = 0;
end
H3 = zeros(2, coef); %H3 -> Matriz de extremos
H3(1,1) = 6*x(1,1);
H3(1,2) = 2;
H3(2,end-2) = 2;
H3(2,end-3) = 6*x(end);
%Revisar H3
H = [H0;H1;H2;H3];%Matriz
%Solución a coeficientes
X = H\BY;
%Gráfica
L = 1;U = 1;kx = size(x,1);%Tamaño datos x
KX = size(X,1);%Tamaño vector coef
while L < KXsyms s;
d = X(L,1) * s.^3;
L = L +1;
c = X(L,1) * s.^2;
L = L +1;
b = X(L,1) * s;
L = L+1;
a = X(L,1);
L = L+1;
funcion = d + c + b + a;%Polinomio
%Gráfica del polinomio
while U < kxdi = x(U,1);
din = x(U+1,1);
dom = [di din];
fplot(funcion,dom);
U = U+1;
grid on
end
end
Valora esta pregunta
0
