Äquidistant - Warum so beliebt?
Warum sind eigentlich "äquidistante Stützstellen" so beliebt? Sicherlich ist es intuitiv eine "gerechte" Einteilung eines Intervalls [a,b]. Gestellt hat sich mir die Frage bei der Lektüre eines Buches über numerische Verfahren. Gestartet wird mit Polynominterpolation und festgestellt, dass die beste Stützstellenwahl die Tschebyscheff-Knoten (*) wären.
(*)
Mit diesem Wissen geht es zur numerischen Integration. Dort wird die Idee "integriere nicht f, sondern das Interpolationspolynom" verwendet. Als konkrete Beispiele werden dann die Newton-Cotes-Formeln verwendet, die wieder äquidistante Stützstellen verwenden. Warum, wenn diese Wahl doch nicht optimal ist? Sind die Tschebyscheff-Knoten "zu teurer/zu fehleranfällig" in der Berechnung?
Es soll hier nicht darum gehen, dass erst Intervallunterteilungen (Spline, zusammengesetzte Integrationsformeln) wirklichen "Erfolg" bringen
Gruß,
tigerbine 
|
RE: Äquidistant - Warum so beliebt?
Push.
Mal etwas losgelöster gefragt, wie aufwendig ist denn die Berechnung von Cosinus (und damit wohl auch sinus) Werten? Wird dazu eine Reihendarstellung benutzt? |
Wenn du die Geduld aufbringst, kannst du ja die Quellen der GNU Scientific Library studieren. Die offizielle Bibliotheksdokumentation allein ist nämlich leider nicht sehr hilfreich: Da werden nur die Schnittstellen, aber nicht die intern verwendeten Algorithmen beschrieben - schade.
Auf jeden Fall eine reiche Quelle zum Studium dessen, wie diverse transzendente Funktionen vergleichsweise effizient numerisch umgesetzt werden.
Bei Sinus, zu fínden in /specfunc/trig.c , steht da was von "Chebyshev expansion" - hier ein Auszug aus diesem File:
| code: |
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:
|
/* Chebyshev expansion for f(t) = g((t+1)Pi/8), -1<t<1
* g(x) = (sin(x)/x - 1)/(x*x)
*/
static double sin_data[12] = {
-0.3295190160663511504173,
0.0025374284671667991990,
0.0006261928782647355874,
-4.6495547521854042157541e-06,
-5.6917531549379706526677e-07,
3.7283335140973803627866e-09,
3.0267376484747473727186e-10,
-1.7400875016436622322022e-12,
-1.0554678305790849834462e-13,
5.3701981409132410797062e-16,
2.5984137983099020336115e-17,
-1.1821555255364833468288e-19
};
static cheb_series sin_cs = {
sin_data,
11,
-1, 1,
11
};
/* Chebyshev expansion for f(t) = g((t+1)Pi/8), -1<t<1
* g(x) = (2(cos(x) - 1)/(x^2) + 1) / x^2
*/
static double cos_data[11] = {
0.165391825637921473505668118136,
-0.00084852883845000173671196530195,
-0.000210086507222940730213625768083,
1.16582269619760204299639757584e-6,
1.43319375856259870334412701165e-7,
-7.4770883429007141617951330184e-10,
-6.0969994944584252706997438007e-11,
2.90748249201909353949854872638e-13,
1.77126739876261435667156490461e-14,
-7.6896421502815579078577263149e-17,
-3.7363121133079412079201377318e-18
};
static cheb_series cos_cs = {
cos_data,
10,
-1, 1,
10
};
/*-*-*-*-*-*-*-*-*-*-*-* Functions with Error Codes *-*-*-*-*-*-*-*-*-*-*-*/
/* I would have prefered just using the library sin() function.
* But after some experimentation I decided that there was
* no good way to understand the error; library sin() is just a black box.
* So we have to roll our own.
*/
int
gsl_sf_sin_e(double x, gsl_sf_result * result)
{
/* CHECK_POINTER(result) */
{
const double P1 = 7.85398125648498535156e-1;
const double P2 = 3.77489470793079817668e-8;
const double P3 = 2.69515142907905952645e-15;
const double sgn_x = GSL_SIGN(x);
const double abs_x = fabs(x);
if(abs_x < GSL_ROOT4_DBL_EPSILON) {
const double x2 = x*x;
result->val = x * (1.0 - x2/6.0);
result->err = fabs(x*x2*x2 / 100.0);
return GSL_SUCCESS;
}
else {
double sgn_result = sgn_x;
double y = floor(abs_x/(0.25*M_PI));
int octant = y - ldexp(floor(ldexp(y,-3)),3);
int stat_cs;
double z;
if(GSL_IS_ODD(octant)) {
octant += 1;
octant &= 07;
y += 1.0;
}
if(octant > 3) {
octant -= 4;
sgn_result = -sgn_result;
}
z = ((abs_x - y * P1) - y * P2) - y * P3;
if(octant == 0) {
gsl_sf_result sin_cs_result;
const double t = 8.0*fabs(z)/M_PI - 1.0;
stat_cs = cheb_eval_e(&sin_cs, t, &sin_cs_result);
result->val = z * (1.0 + z*z * sin_cs_result.val);
}
else { /* octant == 2 */
gsl_sf_result cos_cs_result;
const double t = 8.0*fabs(z)/M_PI - 1.0;
stat_cs = cheb_eval_e(&cos_cs, t, &cos_cs_result);
result->val = 1.0 - 0.5*z*z * (1.0 - z*z * cos_cs_result.val);
}
result->val *= sgn_result;
if(abs_x > 1.0/GSL_DBL_EPSILON) {
result->err = fabs(result->val);
}
else if(abs_x > 100.0/GSL_SQRT_DBL_EPSILON) {
result->err = 2.0 * abs_x * GSL_DBL_EPSILON * fabs(result->val);
}
else if(abs_x > 0.1/GSL_SQRT_DBL_EPSILON) {
result->err = 2.0 * GSL_SQRT_DBL_EPSILON * fabs(result->val);
}
else {
result->err = 2.0 * GSL_DBL_EPSILON * fabs(result->val);
}
return stat_cs;
}
}
} |
|
Man erkennt "Reihenbruchstücke" - nicht genau Taylor - und jede Menge Fallunterscheidungen. 
|
Oha. Dann leg ich mir schon einmal ein Snickers raus, wenn's mal wieder länger dauert. 
Dann vielleicht doch lieber die "äquidistanten Stützstellen". Aber vielleicht gibt es das das Pendant zur Tschebycheffknoten-Wahl ja doch, und ich "sehe" es gerade nicht (sicher muss sich eine Theorie nicht auf den roten Faden des Autors meines Buches gründen). Dort wird nach Newton Cotes gleich auf die zusammengesetzten Formeln eingegangen (entspricht ja dem Spline Ansatz) und erst dann die Frage gestellt, welche maximale Ordnung Quadratur-Formeln haben können.
Da kommt der Name Tschebyscheff ja nun wieder vor (die Tschebycheff-Polynome sind orthogonal) und solche braucht man ja bei der Gauß-Quadratur. Vielleicht steckt da das Pendant? |
Historie
Recherche bezgl. eines anderen Begriffes brachte:
| Zitat: |
Für den Fall äquidistanter Stützstellenverteilungen hat Cotes 1722 in seiner Harmonia mensurarum79 bis n = 10 die entsprechenden Gewichte ohne die zugehörige Rechnung angegeben, wozu Gauß bemerkt:
Bei dieser Methode liegt durchaus die Voraussetzung gleicher Abstände zwischen den Ordinaten zum Grunde. Allerdings scheint beim ersten Anblick diese Voraussetzung am einfachsten und natürlichsten zu sein und es war noch nicht in Frage gekommen, ob es nicht dem ungeachtet noch vorteilhafter sein könne, Ordinaten in ungleichen Abständen zu Grunde zu legen. Um diese Frage zu entscheiden, musste zuerst die Theorie der Quadraturcoefficienten in unbeschränkter Allgemeinheit entwickelt und der Grad der Genauigkeit bestimmt werden. |
|
