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
|
Function IntegralP(ByVal X As Range, ByVal Y As Range)
'Integral d'une fonction donnée en des points discrets non équidistants
'interpolation par un polynôme du second degré
'identique à la méthode de Simpson si point équidistants
'les X doivent être croissant
'Tailles matrices
Dim L As Long, L2 As Long, C As Long, C2 As Long
L = X.Rows.Count
L2 = Y.Rows.Count
C = X.Columns.Count
C2 = Y.Columns.Count
'Erreur de taille
If C > 1 Or C2 > 1 Then IntegralP = "#COLONNE!": Exit Function
If L <> L2 Then IntegralP = "#LIGNE!": Exit Function
If L < 3 Then IntegralP = "#POINT!": Exit Function
'Calcul initial
Dim a0 As Double, b0 As Double, k As Double
a0 = ((Y(1) - Y(3)) * (X(2) - X(3)) - (X(1) - X(3)) * (Y(2) - Y(3))) / ((X(1) ^ 2 - X(3) ^ 2) * (X(2) - X(3)) - (X(1) - X(3)) * (X(2) ^ 2 - X(3) ^ 2))
b0 = (Y(2) - Y(3) - a0 * (X(2) ^ 2 - X(3) ^ 2)) / (X(2) - X(3))
k = 2 * a0 * X(1) + b0
'Algo
Dim i As Long, ai As Double, bi As Double, ci As Double
For i = 2 To L
'coef de l'eq polynome second degré : ax^2+bx+c
ai = (Y(i) - Y(i - 1) - k * (X(i) - X(i - 1))) / (X(i) - X(i - 1)) ^ 2
bi = k - 2 * ai * X(i - 1)
ci = Y(i) - ai * X(i) ^ 2 - bi * X(i)
'integration : [ax^3/3+bx^2/2+cx]
IntegralP = IntegralP + ai * (X(i) ^ 3 - X(i - 1) ^ 3) / 3 + bi * (X(i) ^ 2 - X(i - 1) ^ 2) / 2 + ci * (X(i) - X(i - 1))
'dérivé au point de départ pour continuité au segment suivant
k = 2 * ai * X(i) + bi
Next i
End Function |