Schätzung - Universität der Bundeswehr München

Ein Verfahren zur Schätzung von Formfaktoren,
Magnus-Koeffizienten usw. aus Flugbahndaten
Prof. Dr. Andreas Rudolph
Universität der Bundeswehr
München
Institut für Mathematik und Informatik
FB BW
Werner-Heisenberg-Weg 39
85577 Neubiberg
Email: [email protected]
8. September 2015
Inhaltsverzeichnis
1 Überblick
1.1 Einführung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Mathematische Grundlagen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2 Diverse Algorithmen
2.1 Ein einführendes Beispiel - das Kaliber .308 . . . . . . . . . . . . . . .
2.2 Das Kaliber .300 Win Mag . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 Das Kaliber .338 Lapua Magnum . . . . . . . . . . . . . . . . . . . . .
2.4 Schätzung des Formfaktors mithilfe der Geschwindigkeit . . . . . . . .
2.5 Schätzung des Formfaktors mithilfe der Energie . . . . . . . . . . . . .
2.6 Schätzung des Formfaktors mithilfe der Zeiten . . . . . . . . . . . . . .
2.7 Verallgemeinerungen . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.7.1 Schätzung der Parameter aus Flughöhe und Seitenabweichung
2.7.2 Schätzung der Parameter aus der Geschwindigkeit . . . . . . .
2.7.3 Schätzung der Parameter aus der Geschossenergie . . . . . . .
2.8 Bereitstellung der Daten . . . . . . . . . . . . . . . . . . . . . . . . . .
2.9 Auswahl eines adäquaten Luftwiderstandsmodells . . . . . . . . . . . .
2.10 Erweiterungen des Verfahrens . . . . . . . . . . . . . . . . . . . . . . .
v
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
1
2
11
11
15
17
18
21
21
23
26
27
28
29
29
30
vi
INHALTSVERZEICHNIS
Kapitel 1
Überblick
1.1 Einführung
Bei allen Modellen zur Beschreibung der Flugbahn eines Geschosses liegt die Newton’sche Mechanik zugrunde, d. h. die Gesamtheit aller auf das Geschoss einwirkenden Kräfte F ist gleich
deren Masse m multipliziert mit seiner Beschleunigung a:
F = m·a = m
d2
x
dt 2
wenn x = x(t) die Ortskoordinate des Geschosses zum Zeitpunkt t darstellt. Hierbei unterstellt
man, dass das Projektil als Massenpunkt dargestellt werden kann - man beschreibt strenggenommen die Bahn seines Schwerpunktes1 .
In erster Näherung vernachlässigt man zuerst die Ausdehnung des Geschosses - Rotationseffekte
wie die Präzession, die Nutation, den Magnus-Effekt und den Spin-Dämpfungseffekt betrachtet
man dann in einem zweiten Schritt.
Formuliert man nun die Kräfte, die auf das Geschoss einwirken, und löst das obige Differentialgleichungssystem zweiter Ordnung
F = m·a = m
d2
x
dt 2
dann erhält man damit die Flugbahn des Geschosses in sogenannter parametrisierter Form damit kann man in Abhängigkeit von der Zeit alle Größen bestimmen, die für das sportliche,
aber auch für das militärische Schießen von Interesse sind.
1 Anscheinend
wird der Lagrange-Formalismus und die Hamilton-Jacobi Theorie bisher nicht verwendet - die
Möglichkeit allgemeinerer Koordinatensysteme bleibt unbenommen.
1
KAPITEL 1. ÜBERBLICK
2
1.2 Mathematische Grundlagen
Betrachtet man das vorliegende Differentialgleichungssystem zweiter Ordnung
F = m·a = m
d2
x
dt 2
dann stellt sich die Frage, ob und unter welchen Bedingungen ein derartiges System lösbar ist,
und wenn es lösbar ist, wann die Lösung eindeutig ist.
Dazu schreibt man das System um in ein System erster Ordnung - man erhält dadurch prinzipiell
ein System der Gestalt
y′ = f(t, y)
dem noch entsprechende Anfangsbedingungen y(t0 ) = y0 hinzugefügt werden.
Hierbei ist nun




y1 (t)
f1 (t, y1 , . . . , yn )




..
y(t) =  ...  und f(t, y) = 

.
yn (t)
fn (t, y1 , . . . , yn )
Gesucht ist damit ein Vektor y(t) von Funktionen y1 (t), . . . , yn (t), der das obige Differentialgleichungssystem erster Ordnung samt Anfangsbedingungen erfüllt.
Beispiel 1
Wir betrachten das sogenannte Vakuum-Modell, d. h. der Luftwiderstand wird vernachlässigt.
Damit haben wir als Differentialgleichungssystem
mẍ =
0
mÿ = −mg
wenn nur die Flugweite x = x(t) und die Flughöhe y = y(t) betrachtet werden sollen. Der Punkt
soll dabei die Ableitung nach der Zeit bedeuten.
Zusätzlich stellt man Anfangsbedingungen
x(0)
y(0)
vx (0)
vy (0)
=
=
=
=
x0
y0
v0 cos(θ )
v0 sin(θ )
wenn v0 die Anfangsgeschwindigkeit des Geschosses ist und g die Erdbeschleunigung.
Dieses Differentialgleichungssystem zweiter Ordnung überführen wir durch Division mit m in
1.2. MATHEMATISCHE GRUNDLAGEN
3
ẍ =
0
ÿ = −g
Wir definieren die Vektoren


x
 y 

y=
 ẋ  ,
ẏ


ẋ
 ẏ 

f(t, y) = 
 0 
−g
und erhalten damit

 
x
ẋ
  ẏ
d 
y
 =
dt  ẋ   0
ẏ
−g




oder
ẏ = f(t, y)
ein Differentialgleichungssystem erster Ordnung.
Damit stellt sich die Frage nach der Existenz und Eindeutigkeit einer Lösung eines derartigen
Differentialgleichungssystems erster Ordnung.
Ein Satz, der die Existenz und Eindeutigkeit von Lösungen von Differentialgleichungssystemen
erster Ordnung generell sicherstellt, ist der bekannte Satz von Picard und Lindelöf - dieser besagt:
Satz 1
Wir betrachten das Differentialgleichungssystem
y′ = f(x, y),
y(u) = v
Die Abbildung f sei in dem (n + 1)-dimensionalen Zylinder
| x − u |< a,
k y − v k< b
stetig und beschränkt, sie soll weiter in diesem eine Lipschitz-Bedingung erfüllen, d. h. es gelte
k f(x, y) k≤ A,
k f(x, y1 ) − f(x, y2 ) k≤ M k y1 − y2 k
Dann gibt es für das Differentialgleichungssystem
y′ = f(x, y)
KAPITEL 1. ÜBERBLICK
4
genau eine Lösung y = y(x), welche die Anfangsbedingung y(u) = v erfüllt. und für | x − u |< α
existiert, wobei
b
α = min{a, }
A
ist.
B EWEIS
Siehe zum Beispiel Kamke (1969), S. 74, Satz 1 und S. 104, Satz 3.
Beispiel 2
Betrachten wir nochmals das obige Beispiel des Vakuum-Modells, dann besagt der Satz von
Picard-Lindelöf, dass genau eine Lösung existiert.
Dies bedeutet, dass ein Vektor y(t) existiert mit


x(t)
 y(t) 

y(t) = 
 vx (t) 
vy (t)
der die Anfangsbedingungen x(0) = x0 , y(0) = y0 , vx (0) = v0 cos(θ ) und vy (0) = v0 sin(θ ) erfüllt,
wenn v0 die Anfangsgeschwindigkeit darstellt.
Diese Lösung kann man durch entsprechende Integrationen oder durch einen Potenzreihenansatz erhalten2 - es ergibt sich:
x(t) = v0 cos(θ )t
y(t) = y0 + v0 sin(θ )t − 21 gt 2
und durch Differentiation
vx (t) = v0 cos(θ )
vy (t) = v0 sin(θ ) − gt
Damit hängt die Lösung von einem Parameter ab - dem Erhebungswinkel θ .
Fordert man nun, dass in einer vorgegebenen Entfernung x(tZiel ) gelten soll
y(tZiel ) = 0
dann stellt sich damit die Frage, wie der Winkel θ gewählt werden muss, um dies zu erreichen.
Insbesondere stellt sich die Frage, ob in jedem Fall ein entsprechender Winkel existiert.
2 Es
gibt auch die Möglichkeit, die Picard-Lindelöf-Iterationen durchzuführen - eine Frage des Geschmacks
1.2. MATHEMATISCHE GRUNDLAGEN
5
Generell gilt hierzu folgender
Satz 2
Sei G(x, y) ein Gebiet, c = (c1 , . . . , ck ) ein k-dimensionaler Vektor, C sei eine offene Menge in dem
c-Raum und G × C das Kartesische Produkt von G und C - d. h. der Bereich der Punkte x, y, c eines
(n + k + 1)-dimensionalen Raumes mit x, y ∈ G und c ∈ C.
Die Abbildung f(x, y, c) sei in G × C stetig, so dass für jedes c ∈ C durch jeden Punkt x, y ∈ G mindestens eine Integralkurve der Differentialgleichung
y′ = f(x, y, c)
geht. Es wird weiter vorausgesetzt, dass es auch nur eine solche von Rand zu Rand laufende Integralkurve gibt, die mit
y = Y(x, u, v, c)
bezeichnet wird.
Dann gilt folgendes:
(a) Wenn die Funktion y = Y(x, u, v, c) an einer Stelle x0 , u0 , v0 , c0 existiert und wenn α ≤ x ≤ β ein
abgeschlossenes Intervall ist, dass x0 , u0 als innere Punkte enthält und in dem Y(x, u0 , v0 , c0 )
existiert, dann existieren die Funktionen y = Y(x, u, v, c) in diesem Intervall auch für alle hinreichend nahe an u0 , y0 , c0 gelegenen u, y, c.
(b) Die Funktionen y = Y(x, u, v, c) sind an jeder Stelle ihres Existenzbereiches stetige Funktionen
ihrer Argumente x, u, v, c.
B EWEIS
Siehe Kamke (1969), S. 116, Satz 2.
Beispiel 3
Wir betrachten wiederum das Vakuum-Modell mit der Lösung
x(t) = v0 cos(θ )t
y(t) = y0 + v0 sin(θ )t − 21 gt 2
Offensichtlich hängen die Lösungen in stetiger Weise von t und θ ab.
Gibt man nun ein x(tZiel ) = x0 vor und fordert, dass y(tZiel ) = 0, dann liegt es nahe, das folgende
Funktional zu minimieren:
1 2
))2
h(θ ) = (x(tZiel ) − v0 cos(θ )tZiel )2 + (y(tZiel ) − (y0 + v0 sin(θ )tZiel − gtZiel
2
Da es eine maximale Flugzeit gibt und ebenso einen maximalen Winkel für θ , wird die Funktion
h(.) auf einem kompakten Intervall zu betrachten sein.
KAPITEL 1. ÜBERBLICK
6
Aus Stetigkeitsgründen wird somit ein entsprechender Winkel existieren müssen.
Um den notwendigen Erhebungswinkel θ zu berechnen - will man in einer vorgegebenen Entfernung x0 treffen - eliminieren wir zuerst den Parameter t, indem wir die erste Gleichung nach
t auf lösen
t=
x(t)
v0 cos(θ )
und dann t in der zweiten Gleichung eliminieren:
1
x2
y(x) = y0 + tan(θ )x − g 2 2
2 v0 cos (θ )
Wir beachten nun die Beziehung
1
cos2 (θ )
= 1 + tan2 (θ )
und erhalten
0 = y(x0 ) = y0 + tan(θ )x0 −
gx20 gx20 2
−
tan (θ )
2v20 2v20
Durch Umstellen ergibt sich hieraus eine quadratische Gleichung in tan(θ )
gx20
gx2
tan2 (θ ) − tan(θ )x0 = y0 − 20
2
2v0
2v0
oder
tan2 (θ ) − 2
v20 x0
2v20 y0
θ
)
=
tan(
−1
gx20
gx20
Mittels einer quadratischen Ergänzung ergibt sich hieraus
2 2 2 2
v0
v20
v0
2v2 y0
tan (θ ) − 2
=
+ 02 − 1
tan(θ ) +
gx0
gx0
gx0
gx0
2
und damit
tan(θ ) =
v20
gx0
±
s
v20
gx0
2
+
2v20 y0
−1
gx20
1.2. MATHEMATISCHE GRUNDLAGEN
7
woraus sich mithilfe der Umkehrfunktion des Tangens der gesuchte Winkel θ ergibt - sofern man
fordert, dass der Winkel zwischen 0 und 45◦ liegen soll. Man nimmt hierzu die negative Wurzel.
In dem obigen Beispiel konnten wir den notwendigen Winkel θ analytisch bestimmen - dies ist
sicherlich die Ausnahme.
Man kann aber das obige Beispiel verallgemeinern.
Sei hierzu ein Vektor yZiel gegeben, dann ist ein xziel und cZiel gesucht mit
h(xZiel , cZiel ) =k yZiel − Y(xZiel , u, v, cZiel ) k2 = min
Das Funktional
h(x, c) =k yZiel − Y(x, u, v, c) k2
soll somit minimiert werden.
Unterstellen wir nun, dass
- das ballistische Problem auf einer kompakten Teilmenge D von G × C betrachtet werden
kann
- yZiel aus dem Bildraum Y(x, u, v, c)(D) entstammt
dann muss es nach dem Satz von Weierstrass vom Minimum und Maximum ein xZiel und ein
cZiel geben, welches das obige Funktional minimiert - das obige Funktional h(c) ist ja nach dem
vorhergehenden Satz stetig.
Das wesentliche Problem besteht damit darin, wie die Größe cZiel berechnet werden kann.
Dieses Optimierungsproblem kann aber mit Hilfe eines Verfahrens wie dem BFGS-Verfahren oder
dem DFP-Verfahren gelöst werden - siehe hierzu die Funktion fminunc aus der Programmersprache Octave.
Für die weiteren Betrachtungen wollen wir das vorliegende Problem etwas verallgemeinern hierbei übernimmt der Rohrerhebungswinkel (und gegebenenfalls der Seitenwinkel) die Rolle
eines Hilfsparameter (im folgenden Satz der Parametervektor b) und der Formfaktor (und gegebenenfalls weitere Parameter wie der Magnus-Koeffizient, Lift-Koeffizient etc.) die Rolle des
eigentlichen Parameters (im folgenden Satz der Parametervektor c).
Der folgende Satz stellt dabei auch noch die Differenzierbarkeit der Integralkurven
y = Y(x, u, v, c)
sicher.
Dabei sei c = (c1 , . . . , ck ) ein k-dimensionaler Vektor, b = (b1 , . . . , bl ) ein l-dimensionaler Vektor, C
eine offene Menge im c-Raum und B eine offene Menge im b-Raum.
KAPITEL 1. ÜBERBLICK
8
Ist G(x, y) ein Gebiet, dann sei G × C × B die Menge aller Punkte x, y, c, b mit x, y ∈ G, c ∈ C, b ∈ B.
Es gilt dann
Satz 3
Die Abbildung f(x, y; c, b) sei in G × C × B stetig und besitze stetige partielle Ableitungen bzgl y1 , . . . ,
yn , c1 , . . . , ck .
Für jeden Punkt c, b von C × B gibt es dann durch den Punkt u, v ∈ G eine in G von Rand zu Rand
verlaufende Integralkurve
y = Y(x, u, v, c, b)
der Differentialgleichung
y′ = f(x, y; c, b)
Die Integralkurve y = Y(x, u, v, c, b) ist an jeder Stelle x, u, v, c, b ihres Existenzbereiches nach x, u, v1 ,
. . . , vn , c1 , . . . , ck partiell differenzierbar, und diese partiellen Ableitungen sind stetige Funktionen
von x, u, v, c, b.
Ferner ist die partielle Ableitung
wp =
∂Y
,
∂ cp
1≤ p≤k
die Lösung der linearen Differentialgleichung
w′ = fc p (x, Y(x, u, v, c, b)) + F(x, Y(x, u, v, c, b), c, b)w
mit w p (u) = 0
Hierbei ist F die Matrix
F(x, y, c, b) =
∂ fp
,
∂ yq
p, q = 1, . . . , n
B EWEIS
Siehe Kamke (1969), S. 122, Satz 2.
Anmerkung
Ist die Abbildung f(x, y, c, b) entsprechend oft partiell differenzierbar, dann überträgt sich dies
auf die Integralkurven
y = Y(x, u, v, c, b)
siehe hierzu Kamke (1969), S. 124, Satz 3, und S. 125, Satz 4.
1.2. MATHEMATISCHE GRUNDLAGEN
Die angeführten Sätze bilden nun die Arbeitspferde für alle weiteren Betrachtungen.
9
10
KAPITEL 1. ÜBERBLICK
Kapitel 2
Diverse Algorithmen
2.1 Ein einführendes Beispiel - das Kaliber .308
Wir betrachten das Kaliber .308 und hier das Produkt D46 vom Hersteller Lapua mit den folgenden ballistischen Daten (zuerst die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800
760
704
650
597
452 369
3466 2974 2533 2141 1228 815
Das Geschossgewicht beträgt hier 185 Grains oder 12 Gramm.
Für dieses Produkt ergeben sich folgende Flughöhen (nach Herstellerangaben):
100m 200m 300m
600m
800m
100m
0 −157 −547 −3588 −7908
300m
182
207
0 −2495 −6450
600m
598 1039 1247
0 −3123
Als Luftwiderstandsgesetz unterstellen wir G7 - siehe zum Beispiel Peelen (2014) - und gesucht
sei der zugehörige Formfaktor.
Das zugehörige Differentialgleichungssysstem lautet wie folgt mit der Abkürzung
Ĉd∗ =
ρπ Cd
1
· ρ · S ·Cd =
,
2m
8 C
in Komponentenschreibweise
v˙x = −Ĉd∗ vvx
v˙y = −Ĉd∗ vvy − g
v˙z = −Ĉd∗ vvz
11
C=
m
d2
KAPITEL 2. DIVERSE ALGORITHMEN
12
Abstrahieren wir von der z-Komponente und separieren wir Ĉd∗ in ein Produkt aus Formfaktor
und Luftwiderstandsbeiwert Cd , dann ergeben sich mindestens folgende Möglichkeiten:
- es werden Zeitpunkte t1 , . . . , tn vorgegeben (diese Zeitpunkte müssen dann sinnvoll gewählt werden) und zu diesen Zeitpunkte ti werden die zugehörigen Flugweiten x(ti ) und
die zugehörigen Flughöhen y(ti ) gemessen.
Der Formfaktor wird dann so bestimmt, dass an den vorgegebenen Zeitpunkten ti die Flugweiten x(ti ) und die Flughöhen y(ti ) möglichst gut übereinstimmen.
- es werden Flugweiten x1 , . . . , xn vorgegeben und zu diesen Flugweiten werden Flughöhen
y1 , . . . , yn gemessen.
Der Formfaktor wird dann so bestimmt, dass die Lösungskomponente für die Flughöhe
möglichst gut an die gemessenen Flughöhen angepasst wird.
Da die Komponente für die Flugweite eine streng monotone und differenzierbare Funktion
ist, muss diese Komponente umkehrbar sein, d. h. zu vorgegebenen Flugweiten xi müssen
sich die zugehörigen Zeitpunkte ti berechnen lassen.
Demzufolge können die Flughöhen als Funktion der Flugweite prinzipiell dargestellt werden.
- wieder werden sinnvoll gewählte Zeitpunkte t1 , . . . , tn vorgegeben und die zugehörigen
Geschwindigkeiten vi gemessen.
Der Formfaktor wird dann so bestimmt, dass die gemessenen Geschwindigkeiten möglichst
gut an die durch die Lösung des zugehörigen Differentialgleichungssystems bestimmten
Beträge des Geschwindigkeitsvektors angepasst sind.
- es werden Flugweiten x1 , . . . , xn vorgegeben und die zugehörigen Geschwindigkeiten v1 ,
. . . , vn gemessen.
Der Formfaktor wird dann so bestimmt, dass die gemessenen Geschwindigkeiten möglichst
gut an die durch die Lösung des zugehörigen Differentialgleichungssystems bestimmten
Beträge des Geschwindigkeitsvektors angepasst sind.
Auch hier gilt:
Da die Komponente für die Flugweite eine streng monotone und differenzierbare Funktion
ist, muss diese Komponente umkehrbar sein, d. h. zu vorgegebenen Flugweiten xi müssen
sich die zugehörigen Zeitpunkte ti berechnen lassen.
Demzufolge können die Geschwindigkeiten als Funktion der Flugweite prinzipiell dargestellt werden.
- Da die kinetischen Energien Funktionen der Geschwindigkeit sind, ist auch der Weg über
die kinetischen Energien denkbar.
Wir erhalten damit folgendes Differentialgleichungssystem
2.1. EIN EINFÜHRENDES BEISPIEL - DAS KALIBER .308
13
1
· ρ · S · i ·Cd vvx
2m
1
· ρ · S · i ·Cd vvy − g
= −
2m
v˙x = −
v˙y
wobei i der zu bestimmende Formfaktor ist (siehe auch das Buch von Peelen (2014)), ρ die
Luftdichte und S die Referenzfläche.
v
Cd ist der von der Mach-Zahl abhängige Luftwiderstandsbeiwert, demzufolge Cd = Cd ( ), wenn
a
a die Schallgeschwindigkeit ist.
Wir überführen nun das Differentialgleichungssystem zweiter Ordnung in ein System erster Ordnung, indem wir definieren:



x
 y 

′


y=
 vx  ⇒ y = 
vy

vx
vy 

v′x 
v′y
und damit

′

vx
vy
x

 y  



 =
1
y′ = 
 vx   −
· ρ · S · i ·Cd vvx

 
2m

1
vy
−
· ρ · S · i ·Cd vvy − g
2m
wobei v =
q








v2x + v2y , m die Geschossmasse und ρ die Luftdichte ist.
Nehmen wir zuerst an, wir geben eine Fleckschussentfernung von zum Beispiel x(t100 ) = 100
Metern vor, und die Flughöhen sind gegeben.
Ebenso haben wir Anfangsbedingungen x(0) = 0, y(0) = y0 , vx (0) = v0 cos(θ ) und vy (0) = v0 sin(θ ),
wobei der Rohrerhebungswinkel θ als ein Parameter fungiert.
Die Lösung des Differentialgleichungssystems besteht dann aus zwei Funktionen x(t) und y(t)
wie auch Geschwindigkeiten vx (t) und vy (t).
Diese Funktionen hängen damit von dem Formfaktor i ab, d. h.
x(t)
y(t)
vx (t)
vy (t)
=
=
=
=
Gegeben seien nun folgende Daten vom Anfang:
x(t; i)
y(t; i)
vx (t; i)
vy (t; i)
KAPITEL 2. DIVERSE ALGORITHMEN
14
100m 200m 300m
600m
800m
100m
0 −157 −547 −3588 −7908
300m
182
207
0 −2495 −6450
600m
598 1039 1247
0 −3123
Wollen wir somit die Flughöhen möglichst genau vorhersagen, müssen wir den Formfaktor so
bestimmen, dass das folgende Funktional minimal wird:
h(i) =
(0 − y(t100 ; i))2 + (−0, 157 − y(t200 ; i))2 + (0, 547 − y(t300 ; i))2
+(−3, 588 − y(t600 ); i))2 + (−7, 908 − y(t800 ); i))2
Dabei muss zuerst der Rohrerhebungswinkel für die Fleckschussweite bestimmt und danach das
angegebene Funktional minimiert werden.
Eine geeignete Funktion in Octave ist die Funktion fminbnd - hier liegt ja eine eindimensionale
Funktion vor, sie basiert auf dem Verfahren des goldenen Schnitts.
Aufgrund der vorhergehenden Sätze muss eine optimale Lösung für den Formfaktor existieren es ist dann sicher interessant, einen Vergleich anzustellen, welche Ergebnisse mit der Angabe für
den Formfaktor aus der Literatur erhalten werden und welche Ergebnisse mit dem optimierten
Formfaktor erzielt werden.
Zu den oben angegebenen Daten ergibt sich durch Rechnung mit dem Formfaktor 1, 081 aus dem
Buch von Litz (2011), S. 530, und dem Luftwiderstandsgesetz G7:
100m 200m 300m
600m
800m
100m
0 −158 −549 −3584 −7863
300m
183
208
0 −2486 −6400
600m
597 1068 1243
0 −3085
Die erste Spalte beinhaltet die Fleckschussweiten.
Berechnen wir gemäß der vorhergehenden mathematischen Sätze den optimalen Formfaktor zu
1, 08944, so erhalten wir die folgende Tabelle:
100m 200m 300m
600m
800m
100m
0 −158 −550 −3597 −7905
300m
183
208
0 −2498 −6439
600m
600 1041 1249
0 −3109
Hier nochmals zum Vergleich die Herstellerdaten:
2.2. DAS KALIBER .300 WIN MAG
15
100m 200m 300m
600m
800m
100m
0 −157 −547 −3588 −7908
300m
182
207
0 −2495 −6450
600m
598 1039 1247
0 −3123
Wir betrachten nun die Daten mit dem Produkt-Code GB491 von Lapua (Kaliber .308). Hier
haben wir folgende ballistische Daten (zuerst die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800 1000m
860
794
730
669
503 408
331
3698 3151 2663 2235 1264 833
548
Das Geschossgewicht beträgt hier 155 Grains oder 10 Gramm.
Für dieses Produkt ergeben sich folgende Flughöhen (nach Herstellerangaben, die erste Spalte
beinhaltet wieder die Fleckschussweiten):
100m 200m 300m
600m
800m
1000m
100m
0 −115 −414 −2813 −6265 −12106
300m
138
161
0 −1984 −5160 −10726
600m
469
822
992
0 −2515 −7418
und berechnet mit dem Formfaktor 0, 988 aus dem Buch von Litz (2011), S. 527:
100m 200m 300m
600m
800m
1000m
100m
0 −116 −414 −2795 −6206 −11984
300m
138
161
0 −1967 −5101 −10603
600m
466
816
983
0 −2479 −7325
Optimieren wir jetzt gemäß unseren mathematischen Sätzen den Formfaktor und rechnen damit
das Differentialgleichungssystem durch, dann erhalten wir folgende Tabelle:
100m
300m
600m
100m 200m 300m
600m
800m
1000m
0 −116 −416 −2811 −6253 −12111
139
161
0 −1979 −5145 −10725
468
821
990
0 −2506 −7427
Der aus den Daten geschätzte Formfaktor ist hierbei 0, 99852.
Die Berechnungen wurden mit dem Laufdurchmesser 7, 83 mm durchgeführt (siehe hierzu das
Buch von Peelen (2014) auf S. 19, wo eine Tabelle mit den Laufdurchmessern angegeben ist).
2.2 Das Kaliber .300 Win Mag
Als weiteres Beispiel betrachten wir das Kaliber .300 Win Mag.
KAPITEL 2. DIVERSE ALGORITHMEN
16
RWS gibt in der Serie Target Elite Plus für dieses Kaliber folgende Daten an (Geschossgewicht
13 Gramm, ballistischer Koeffizient 0.555):
0m 200m 400m 600m 800 1000m
870
762
661
568 484
407
4904 3758 2833 2091 1521
1073
Die zweite Zeile umfasst die Geschwindigkeiten, die dritte Zeile die kinetische Energie.
Für die Trajektorien gibt RWS folgende Daten an (in Millimetern):
100m
300m
500m
200m 400m
600m
800m 1000m
−98 −827 −2446 −5345 −9924
144 −343 −1721 −4377 −8715
507
382 −633 −2927 −6903
Leider macht RWS keine Angaben (?) über die Höhe des Zielfernrohrs über der Laufachse, es
steht aber zu vermuten, dass es 4 oder 5 Zentimeter sein könnten (dies ist die übliche Höhe).
Bei 5 Zentimetern ergibt sich mit dem vorgestellten Verfahren und Luftwiderstandsgesetz G7 ein
Formfaktor von 1, 07755 - rechnet man damit das Differentialgleichungssystem durch, so erhält
man folgende Tabelle:
100m
300m
500m
200m 400m
600m
800m 1000m
−98 −827 −2463 −5332 −9927
146 −340 −1732 −4358 −8710
508
384 −646 −2909 −6899
Bei einer Höhe des Zielfernrohrs von 4 Zentimetern werden die Ergebnisse schlechter - von
daher wurden die Ergebnisse von RWS wohl mit einer Höhe von 5 Zentimetern gewonnen.
Der angegebenen ballistische Koeffizient scheint für das Luftwiderstandsgesetz G1 zu gelten
- rechnet man das Differentialgleichungssystem mit G1 durch, so erhält man Ergebnisse, die
leidlich zu den angegebenen Daten passen.
Allerdings - so Litz (2011), S. 13 ff. - ist das Luftwiderstandsgesetz G7 geeigneter für Büchsenpatronen mit entsprechenden Geschossformen.
Das Luftwiderstandsgesetz G1 ist eher geeignet für Revolver- und Pistolengeschosse aufgrund
ihrer eher zylindrischen Form - dafür sind die ballistischen Koeffizienten größer (?) .
Anmerkung
Die analoge Rechnung wurde bei gleicher Zielfernrohrhöhe auch für das Kaliber .308 durchgeführt (Geschossgewicht 10,69 Gramm, Mündungsgeschwindigkeit 805 Meter pro Sekunde,
ebenfalls Target Elite plus von RWS).
Hier der Vergleich - zuerst die Herstellerangaben:
2.3. DAS KALIBER .338 LAPUA MAGNUM
17
0m 200m 400m 600m 800 1000m
805
676
555
447 362
311
3529 2487 1675 1088 715
528
für Geschwindigkeit und Energie, nun die Flughöhen
200m
400m
600m
800m
1000m
100m −131 −1078 −3336 −7644 −15035
300m
186 −443 −2385 −6375 −13449
500m
677
538 −913 −4413 −10996
Der optimale Formfaktor ergibt sich zu 1, 12667 - mit diesem Formfaktor das Differentialgleichungssystem mit dem Luftwiderstandsgesetz G7 durchgerechnet ergibt die Trajektorien:
200m
400m
600m
800m
1000m
100m −131 −1087 −3342 −7608 −15048
300m
188 −450 −2386 −6334 −13456
500m
676
526 −922 −4382 −11015
2.3 Das Kaliber .338 Lapua Magnum
Die Munition mit dem Produkt-Code GB528 von Lapua hat folgende ballistische Daten (zuerst
die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800 1000m
830
789
750
712
605 539
478
6696 6054 5465 4925 3558 2827
2223
Das Geschossgewicht beträgt hier 300 Grains oder 19,4 Gramm.
Für dieses Produkt ergeben sich folgende Flughöhen (nach Herstellerangaben):
100m
300m
600m
100m 200m 300m
600m
800m 1000m
0 −117 −413 −2583 −5231 −9289
136
155
0 −1710 −4146 −7932
421
725
855
0 −1865 −5082
100m
300m
600m
100m 200m 300m
600m
800m 1000m
0 −117 −406 −2506 −5178 −9169
135
154
0 −1694 −4094 −7814
418
719
847
0 −1836 −4992
und berechnet:
KAPITEL 2. DIVERSE ALGORITHMEN
18
Auch hier beinhaltet die erste Spalte die Fleckschussweiten.
Der Formfaktor wurde dem Buch von Litz (2011), Seite 557, entnommen zu 0, 956.
Optimieren wir nun den Formfaktor, so erhalten wir hierfür den Wert 0, 98351, und es ergibt sich
folgende Tabelle nach dem Lösen des Differentialgleichungssystems:
100m 200m 300m
600m
800m 1000m
100m
0 −117 −408 −2526 −5233 −9297
300m
136
155
0 −1710 −4145 −7936
600m
421
725
855
0 −1865 −5086
Hier wurde mit dem Laufdurchmesser 8, 59 mm gearbeitet.
2.4 Schätzung des Formfaktors mithilfe der Geschwindigkeit
Wir betrachten nochmals das Kaliber .308 und hier das Produkt D46 vom Hersteller Lapua mit
den folgenden ballistischen Daten (zuerst die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800
760
704
650
597
452 369
3466 2974 2533 2141 1228 815
Das Geschossgewicht beträgt hier 185 Grains oder 12 Gramm.
Falls es schwierig, wenn nicht unmöglich ist, die Flughöhen zu messen, liegt es nahe, stattdessen
die Geschwindigkeiten für die Schätzung des gesuchten Formfaktors zu verwenden.
Die Waffe sei hierzu auf 100 Meter eingeschossen.
Wieder betrachten wir das Differentialgleichungssystem

′
vx
x

vy
 y  



 
1
y′ = 
 vx  =  −
· ρ · S · i ·Cd vvx
 

2m

1
vy
−
· ρ · S · i ·Cd vvy − g
2m

wobei v =
q
v2x + v2y , m die Geschossmasse und ρ die Luftdichte ist.
Die Lösungen seien wie vorhin
x(t)
y(t)
vx (t)
vy (t)
=
=
=
=
x(t; i)
y(t; i)
vx (t; i)
vy (t; i)








2.4. SCHÄTZUNG DES FORMFAKTORS MITHILFE DER GESCHWINDIGKEIT
19
sie hängen somit von dem Formfaktor i ab.
Wir bilden nun die Geschwindigkeit v = v(t) =
tional
p
vx (t; i)2 + vy (t; i)2 und damit das folgende Funk-
h(i) =
(704 − v(t100 ; i))2 + (650 − v(t200 ; i))2 + (597 − v(t300 ; i))2
+(452 − v(t600 ); i))2 + (369 − v(t800 ); i))2
Dabei werden die Zeitpunkte t j , j = 100, j = 200, j = 300, j = 600 und j = 800, so bestimmt, dass
x(t j ; i) = j.
Dieses Funktional muss nun bzgl. des Formfaktors i minimiert werden - es liegt nahe, wie in
den vorhergehenden Beispielen vorzugehen und die Funktion fminbnd aus dem Programmpaket
Octave zu verwenden.
Wir erhalten dann wiederum den optimalen Formfaktor, mit dem dann die Flugbahn berechnet
werden kann - hier ergibt sich für den Formfaktor der Wert 1, 09748.
Rechnet man hierfür das Differentialgleichungssystem durch, dann erhält man folgende Tabelle:
100m 200m 300m
600m
800m
100m
0 −159 −551 −3610 −7945
300m
184
209
0 −2508 −6476
600m
602 1045 1254
0 −3132
Hier nochmals zum Vergleich die Herstellerdaten:
100m 200m 300m
600m
800m
100m
0 −157 −547 −3588 −7908
300m
182
207
0 −2495 −6450
600m
598 1039 1247
0 −3123
Betrachten wir nun noch die zweite Laborierung des Kalibers .308 - hier hatten wir folgende
Angaben für die Geschwindigkeit und die kinetische Energie:
0m 100m 200m 300m 600m 800 1000m
860
794
730
669
503 408
331
3698 3151 2663 2235 1264 833
548
Für dieses Produkt ergeben sich folgende Flughöhen (nach Herstellerangaben, die erste Spalte
beinhaltet wieder die Fleckschussweiten):
KAPITEL 2. DIVERSE ALGORITHMEN
20
100m 200m 300m
600m
800m
1000m
100m
0 −115 −414 −2813 −6265 −12106
300m
138
161
0 −1984 −5160 −10726
600m
469
822
992
0 −2515 −7418
und berechnet mit dem Formfaktor 0, 988 aus dem Buch von Litz (2011), S. 527:
100m 200m 300m
600m
800m
1000m
100m
0 −116 −414 −2795 −6206 −11984
300m
138
161
0 −1967 −5101 −10603
600m
466
816
983
0 −2479 −7325
Optimieren wir jetzt gemäß unseren mathematischen Sätzen den Formfaktor und rechnen damit
das Differentialgleichungssystem durch, dann erhalten wir folgende Tabelle:
100m
300m
600m
100m 200m 300m
600m
800m
1000m
0 −116 −416 −2811 −6253 −12111
139
161
0 −1979 −5145 −10725
468
821
990
0 −2506 −7427
Optimieren wir nun den Formfaktor mittels der Geschwindigkeit, so erhalten wir für den Formfaktor den Wert 0, 99847 und damit folgende Tabelle
100m
300m
600m
100m 200m 300m
600m
800m
1000m
0 −116 −416 −2810 −6253 −12110
139
161
0 −1979 −5145 −10725
468
821
990
0 −2506 −7426
Es verbleibt jetzt noch, die Daten für das Kaliber .338 Lapua Magnum durchzurechnen.
Hier hatten wir folgende Herstellerdaten (zuerst die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800 1000m
830
789
750
712
605 539
478
6696 6054 5465 4925 3558 2827
2223
Das Geschossgewicht beträgt hier 300 Grains oder 19,4 Gramm.
Für dieses Produkt ergeben sich folgende Flughöhen (nach Herstellerangaben):
100m
300m
600m
100m 200m 300m
600m
800m 1000m
0 −117 −413 −2583 −5231 −9289
136
155
0 −1710 −4146 −7932
421
725
855
0 −1865 −5082
2.5. SCHÄTZUNG DES FORMFAKTORS MITHILFE DER ENERGIE
21
Nach dem Optimieren erhielten wir seinerzeit einen Formfaktor mit dem Wert 0, 98351, und
folgende Tabelle nach dem Lösen des Differentialgleichungssystems:
100m 200m 300m
600m
800m 1000m
100m
0 −117 −408 −2526 −5233 −9297
300m
136
155
0 −1710 −4145 −7936
600m
421
725
855
0 −1865 −5086
Optimieren wir den Formfaktor über die Geschwindigkeit, so erhalten wir für den Formfaktor
den Wert 0, 98521, und folgende Tabelle nach dem Lösen des Differentialgleichungssystems:
100m
300m
600m
100m 200m 300m
600m
800m 1000m
0 −117 −408 −2528 −5237 −9305
136
155
0 −1711 −4148 −7944
421
725
856
0 −1867 −5092
2.5 Schätzung des Formfaktors mithilfe der Energie
Die kinetische Energie ist gegeben durch
1
· m · v2
2
von daher besteht die Möglichkeit, die kinetische Energie
- zu vorgegebenen Zeitpunkten t1 , . . . , tn zu messen und ein entsprechendes Funktional bzgl.
des Formfaktors zu minimieren
- zu vorgegebenen Flugweiten x1 , . . . , xn zu messen und ebenfalls ein entsprechendes Funktional bzgl. des Formfaktors zu minimieren
Da die kinetische Energie eine Funktion der Geschwindigkeit ist, steht allerdings zu erwarten,
dass sich bei dieser Vorgehensweise keine neuen Erkenntnisse gewinnen lassen.
Aus praktischer Sicht ist es eine Frage der Bequemlichkeit, ob man lieber mit der Geschwindigkeit oder der kinetischen Energie arbeiten möchte.
2.6 Schätzung des Formfaktors mithilfe der Zeiten
Wir betrachten wieder das Kaliber .308 und hier das Produkt D46 vom Hersteller Lapua mit den
folgenden ballistischen Daten (zuerst die Geschwindigkeit, dann die kinetische Energie):
0m 100m 200m 300m 600m 800
760
704
650
597
452 369
3466 2974 2533 2141 1228 815
KAPITEL 2. DIVERSE ALGORITHMEN
22
Das Geschossgewicht beträgt hier 185 Grains oder 12 Gramm.
Nun gehen wir davon aus, dass die Zeiten gemessen werden, die zum Erreichen der Entfernungen 100, 200, 300, 600 und 800 Meter benötigt werden - diese seien mit t˜100 , t˜200 , t˜300 , t˜600 und t˜800
bezeichnet.
Wiederum soll der Formfaktor i optimal geschätzt werden.
Die Waffe sei hierzu auf 100 Meter eingeschossen.
Wieder betrachten wir das Differentialgleichungssystem

′

vx
vy
x

 y  



 =
1
y′ = 
 vx   −
· ρ · S · i ·Cd vvx

 
2m

1
vy
−
· ρ · S · i ·Cd vvy − g
2m
wobei v =
q








v2x + v2y , m die Geschossmasse und ρ die Luftdichte ist.
Die Lösungen seien wie vorhin
x(t)
y(t)
vx (t)
vy (t)
=
=
=
=
x(t; i)
y(t; i)
vx (t; i)
vy (t; i)
sie hängen somit von dem Formfaktor i ab.
Wir berechnen nun zu den Entfernungen 100, 200, 300, 600 und 800 die Zeiten aus den Gleichungen 100 = x(t100 ; i), 200 = x(t200 ; i), 300 = x(t300 ; i), 600 = x(t600 ; i) und 800 = x(t800 ; i) - die so
gewonnenen Zeitpunkte t j hängen dann natürlich noch von dem Formfaktor i ab.
Wir betrachten dann folgendes Funktional
h(i) =
(t˜100 − t100 )2 + (t˜200 − t200 )2 + (t˜300 − t300 )2
+(t˜600 − t600 )2 + (t˜800 − t800 )2
Dieses Funktional muss dann bzgl. i minimiert werden - dies ergibt dann die optimale Schätzung
für den gesuchten Formfaktor.
Allerdings sind in dieser Situation keine Daten verfügbar - Auswertungen sind somit nicht möglich.
2.7. VERALLGEMEINERUNGEN
23
2.7 Verallgemeinerungen
Der Ausgangspunkt des einführenden Beispiels war das Differentialgleichungssystem
1
· ρ · S · i ·Cd vvx
2m
1
= −
· ρ · S · i ·Cd vvy − g
2m
v˙x = −
v˙y
wobei i der zu bestimmende Formfaktor ist (siehe auch das Buch von Peelen (2014)).
Die erste Erweiterung besteht darin, anstatt des zweidimensionalen Systems ein dreidimensionales System zu betrachten.
1
· ρ · S · i ·Cd vvx
2m
1
· ρ · S · i ·Cd vvy − g
= −
2m
1
= −
· ρ · S · i ·Cd vvz
2m
v˙x = −
v˙y
v˙z
Dieses System zweiter Ordnung kann man natürlich analog in ein System erster Ordnung überführen, indem wir definieren:




y=



x
y
z
vx
vy
vz








 ⇒ y′ = 






vx
vy
vz
v′x
v′y
v′z








und damit






′
y =





wobei v =
q
x
y
z
vx
vy
vz
′

vx
vy

 
 
 
vz
 
 
1
 =
· ρ · S · i ·Cd vvx
−
 

2m
 
 
  − 1 · ρ · S · i ·Cd vvy − g
  2m

1
−
· ρ · S · i ·Cd vvz
2m















v2x + v2y , m die Geschossmasse und ρ die Luftdichte ist.
In dieser Formulierung ist es jetzt möglich, auch den Einfluss des Seitenwindes zu modellieren.
KAPITEL 2. DIVERSE ALGORITHMEN
24
Dieses Differentialgleichungssystem kann man nun um Rotationseffekte erweitern - siehe hierzu
McCoy (1999) oder Carlucci und Jacobson (2008).
Diese Erweiterung beruht auf dem Drehimpulssatz - hierzu müssen wir folgende Größen definieren - siehe McCoy (1999), S. 189:
Cd
CL
CNp
CMp
CMq
Cl p
CMα
=
=
=
=
=
=
=
Luftwiderstandskoeffizient
Auftriebskoeffizient
Koeffizient für die Magnus-Kraft
Koeffizient für das Magnus-Moment
Koeffizient für die Kraft der Pitch-Dämpfung
Koeffizient der Spin-Dämpfung
Koeffizient des Überschlagsmoments
Um in den Differentialgleichungen nicht ständig die Masse wie auch die Spin-Rate mitschleppen
zu müssen, bildet man folgende abgekürzte Größen:
C̃d =
C̃L =
C̃Np
=
C̃Mp
=
C̃Mq
=
C̃l p
=
C̃Mα
=
ρ vSCd
2m
ρ vSCL
2m
ρ vSCNp p
2m
ρ Sd 2CMp p
2Iy
ρ vSdCMq
2m
ρ vSd 2Cl p p
2Iy
ρ vSdCMα
2Iy
wobei
p=
Iy
(h · x)
Ix
die axiale Spin-Rate des Projektils
v1 = V1 −W1 , v2 = V2 −W2 , v3 = V3 −W3 , v =
cos(α ) =
q
v21 + v22 + v23
v1 x1 + v2 x2 + v3 x3
v
der Kosinus des Anstellwinkels des Projektils ist, Wi , i = 1, 2, 3 die Komponenten der Windgeschwindigkeit und ẋi = vi .
2.7. VERALLGEMEINERUNGEN
25
Vi sind die Komponenten des Geschwindigkeitsvektors des Projektils bzgl. eines erdfesten Koordinatensystem, vi die Komponenten des Geschwindigkeitsvektors des Projektils bzgl der Luft.
Ix und Iy sind die Trägheitsmomente des Projektils in x- bzw. y-Richtung des Projektils.
Hiermit können wir dann folgende sechs gewöhnliche Differentialgleichungen aufstellen:
v̇1 = −C̃d v1 + C̃L v2 x1 − vv1 cos(α ) + C̃Np (x2 v3 − x3 v2 ) + C̃Nq (h2 x3 − h3 x2 ) + g1
v̇2 = −C̃d v2 + C̃L v2 x2 − vv2 cos(α ) + C̃Np (x3 v1 − x1 v3 ) + C̃Nq (h3 x1 − h1 x3 ) + g2
v̇3 = −C̃d v3 + C̃L v2 x3 − vv3 cos(α ) + C̃Np (x1 v2 − x2 v1 ) + C̃Nq (h1 x2 − h2 x1 ) + g3
ḣ1 = C̃l p x1 + C̃Mα (v2 x3 − v3 x2 ) + C̃Np (v1 − vx1 cos(α )) + C̃Mq
ḣ2 = C̃l p x2 + C̃Mα (v3 x1 − v1 x3 ) + C̃Mp (v2 − vx2 cos(α )) + C̃Mq
Ix p
x2
h2 −
Iy
ḣ3 = C̃l p x3 + C̃Mα (v1 x2 − v2 x1 ) + C̃Mp (v3 − vx3 cos(α )) + C̃Mq
Ix p
x3
h3 −
Iy
Ix p
h1 −
x1
Iy
wobei g1 , g2 und g3 die Komponenten des Gravitationsvektors sind.
Dabei ist h = (h1 , h2 , h3 )T und H der Drehimpulsvektor - h ist der Drehimpulsvektor H durch Iy
devidiert.
Beachtet man noch, dass v̇i = V̇i , dann nehmen die sechs Differentialgleichungen die folgende
Form an:
V̇1 = −C̃d v1 + C̃L v2 x1 − vv1 cos(α ) + C̃Np (x2 v3 − x3 v2 ) + C̃Nq (h2 x3 − h3 x2 ) + g1
V̇2 = −C̃d v2 + C̃L v2 x2 − vv2 cos(α ) + C̃Np (x3 v1 − x1 v3 ) + C̃Nq (h3 x1 − h1 x3 ) + g2
V̇3 = −C̃d v3 + C̃L v2 x3 − vv3 cos(α ) + C̃Np (x1 v2 − x2 v1 ) + C̃Nq (h1 x2 − h2 x1 ) + g3
ḣ1 = C̃l p x1 + C̃Mα (v2 x3 − v3 x2 ) + C̃Np (v1 − vx1 cos(α )) + C̃Mq
ḣ2 = C̃l p x2 + C̃Mα (v3 x1 − v1 x3 ) + C̃Mp (v2 − vx2 cos(α )) + C̃Mq
Ix p
x1
h1 −
Iy
Ix p
x2
h2 −
Iy
KAPITEL 2. DIVERSE ALGORITHMEN
26
ḣ3 = C̃l p x3 + C̃Mα (v1 x2 − v2 x1 ) + C̃Mp (v3 − vx3 cos(α )) + C̃Mq
Ix p
x3
h3 −
Iy
Dieses System zweiter Ordnung kann man aber leicht umschreiben in ein System erster Ordnung.
Dazu definieren wir (plus geeignete Anfangsbedingungen)




y=



X1
X2
X3
h1
h2
h3








 ⇒ ẏ = 






Ẋ1
Ẋ2
Ẋ3
ḣ1
ḣ2
ḣ3








 = ẏ = 






V1
V2
V3
ḣ1
ḣ2
ḣ3








Fasst man nun die rechten Seiten der sechs gewöhnlichen Differentialgleichungen zu einer Funktion f(y) zusammen, so erhält man offensichtlich ein Differentialgleichungssystem
ẏ = f(y)
dem noch geeignete Anfangsbedingungen hinzuzufügen sind.
Dies bedeutet, dass eine Fleckschussentfernung gewählt worden sein muss, auf die sich die
Anfangsbedingungen beziehen.
Bei dem Differentialgleichungssystem ist nun
f(y) = f(y; C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα ),
wobei C̃d eine Funktion des Formfaktors i ist.
Diese sieben Parameter C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα müssen somit aus geeigneten Daten geschätzt werden.
Setzen wir nun voraus, dass geeignete Daten vorliegen und dass
y(t; C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα )
die Lösung des Differentialgleichungssystems darstellt.
Um die Parameter zu schätzen, muss ein geeignetes Funktional definiert werden - analog zu den
vorhergehenden Beispielen.
2.7.1 Schätzung der Parameter aus Flughöhe und Seitenabweichung
Nehmen wir an, dass die zweite Komponente im Vektor y die Flughöhe und die dritte Komponente im Vektor y die seitliche Abweichung darstellen.
2.7. VERALLGEMEINERUNGEN
27
Sind nun y(ti ) für i = 1, . . . , n die Daten für die Flughöhen und z(ti ) für i = 1, . . . , n die Daten
für die seitlichen Abweichungen bei bestimmten Flugweiten x(ti ), dann liegt es nahe, folgendes
Funktional zu minimieren:
n
n
i=1
i=1
h(C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα )) = ∑ (y(ti ) − X2 (ti ))2 + ∑ (z(ti ) − X3(ti ))2
wobei X2 (t) und X3 (t) jeweils von C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα abhängen.
Dabei sind die ti aus den zugehörigen Flugweiten zu berechnen.
Die Minimierung dieses Funktionals entspricht prinzipiell der Lösung eines nichtlinearen Ausgleichsproblems - die Funktion fminunc aus dem Programm Octave oder Matlab löst diese Aufgabe.
2.7.2 Schätzung der Parameter aus der Geschwindigkeit
Wir kürzen den Parametervektor durch
C̃ = (C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα )T
ab.
Analog zu der Vorgehensweise in den vorausgegangenen Beispielen sei
v(t j ; C̃) =
q
vx (t j ; C̃)2 + vy (t j ; C̃)2 + vz (t j ; C̃)2
der Betrag des Geschwindigkeitsvektors zur Zeit t j , j = 1, . . . , n.
Da die Geschwindigkeit in gewissen Entfernungen j = x(t j ; i) gemessen wird, ergeben sich
hieraus die zugehörigen Beträge der Geschwindigkeitsvektoren.
Seien nun in bestimmten Entfernungen o. B. d. A. z. B. j = 100, j = 200, j = 300, j = 400, j = 500,
j = 600, j = 800 und j = 1000 Meter die zugehörigen Geschwindigkeiten v100 , v200 , v300 , v400 , v500 ,
v600 , v800 und v1000 gemessen.
Wir betrachten dann das Funktional
h(C̃) = (v100 − v(t100 ; C̃))2 + (v200 − v(t200 ; C̃))2 + (v300 − v(t300 ; C̃))2 + (v400 − v(t400 ; C̃))2
+(v500 − v(t500 ; C̃))2 + (v600 − v(t600 ; C̃))2 + (v800 − v(t800 ; C̃))2 + (v1000 − v(t1000 ; C̃))2
KAPITEL 2. DIVERSE ALGORITHMEN
28
oder ganz allgemein
n
h(C̃) = ∑ (v j − v(t j ; C̃))2
i=1
Dieses Funktional muss nun in Abhängigkeit von C̃ minimiert werden - eine derartige Aufgabe
löst die Funktion fminunc aus dem Programm Octave.
2.7.3 Schätzung der Parameter aus der Geschossenergie
Analog zu der Vorgehensweise in den vorausgegangenen Beispielen sei
v(t j ; C̃) =
q
vx (t j ; C̃)2 + vy (t j ; C̃)2 + vz (t j ; C̃)2
der Betrag des Geschwindigkeitsvektors zur Zeit t j , j = 1, . . . , n.
Hieraus kann die kinetische Energie
1
e(t j ; C̃) = mv(t j ; C̃)2
2
berechnet werden.
Seien nun in bestimmten Entfernungen o. B. d. A. z. B. j = 100, j = 200, j = 300, j = 400, j = 500,
j = 600, j = 800 und j = 1000 Meter die zugehörigen kinetischen Energien e100 , e200 , e300 , e400 ,
e500 , e600 , e800 und e1000 gemessen.
Wir betrachten dann das Funktional
h(C̃) = (e100 − e(t100 ; C̃))2 + (e200 − e(t200 ; C̃))2 + (e300 − e(t300 ; C̃))2 + (e400 − e(t400 ; C̃))2
+(e500 − e(t500 ; C̃))2 + (e600 − e(t600 ; C̃))2 + (e800 − e(t800 ; C̃))2 + (e1000 − e(t1000 ; C̃))2
oder ganz allgemein
n
h(C̃) = ∑ (e j − e(t j ; C̃))2
i=1
Dieses Funktional muss nun in Abhängigkeit von C̃ minimiert werden - eine derartige Aufgabe
löst die Funktion fminunc aus dem Programm Octave.
2.8. BEREITSTELLUNG DER DATEN
29
2.8 Bereitstellung der Daten
Bisher wurde vorausgesetzt, dass die Daten zur Schätzung des Formfaktors etc. bereits vorliegen.
Normalerweise müssen diese aber erst erzeugt werden - es muss zum Beispiel auf bestimmten
Entfernungen erst eine bestimmte Anzahl von Schüssen abgegeben werden.
Damit läuft man aber sofort in eine Stichprobenproblematik hinein, d. h. es stellt sich die Frage,
wieviele Treffer muss man aufnehmen in welcher Entfernung.
Eine Möglichkeit besteht darin, ein Konfidenzintervall vorgegebener Länge bzw. ein Konfidenzellipsoid vorgegebener Länge zu betrachten und hierüber den notwendigen Stichprobenumfang
zu berechnen.
Diese Problematik ist aber in der Literatur bekannt, siehe zum Beispiel das Buch von Morrison
(1976) oder das Buch von Seber (1984).
2.9 Auswahl eines adäquaten Luftwiderstandsmodells
Die vorgestellten Verfahren ermöglichen offensichtlich, unter gegebenen Modellen ein optimales
Luftwiderstandsmodell auszuwählen.
Berechnet man hierzu für optimalen Formfaktor die Standardabweichungen des jeweiligen Modells, so liegt es nahe, dasjenige Modell zu wählen, bei dem die Standardabweichung minimal
ist.
Sind beispielsweise G1 bis G8 die zur Auswahl stehenden Luftwiderstandsmodelle mit zugehörigen Standardabweichungen s1 bis s8 , dann würde das Modell gewählt werden, deren zugehörige
Standardabweichung minimal ist.
Eine weitere Möglichkeit besteht darin, das Funktional selbst zur Entscheidung heranzuziehen
- je besser die Anpassung durch das Luftwiderstandsmodell, umso kleiner muss das Funktional
zum Beispiel
h(i) =
(0 − y(t100 ; i))2 + (−0, 157 − y(t200 ; i))2 + (0, 547 − y(t300 ; i))2
+(−3, 588 − y(t600 ); i))2 + (−7, 908 − y(t800 ); i))2
(siehe das Beispiel vom Anfang) in Abhängigkeit von dem Luftwiderstandsgesetz sein. Ein derartiges Funktional beschreibt ja gerade eines Teils der Lösung des Differentialgleichungssystems
an die Daten - im vorliegenden Fall an die Flughöhen.
Betrachten wir hierzu nochmals die Daten des Kalibers .308 der Target Elite Plus von RWS:
KAPITEL 2. DIVERSE ALGORITHMEN
30
200m
400m
600m
800m
1000m
100m −131 −1078 −3336 −7644 −15035
300m
186 −443 −2385 −6375 −13449
500m
677
538 −913 −4413 −10996
Mit der Fleckschussweite von 100 Metern erhalten wir für das Luftwiderstandsgesetz G7 einen
Formfaktor von 0, 12667 und für das Funktional einen Wert von 0, 0015967 - bei einer Zielfernrohrhöhe von 5 Zentimetern.
Bei gleicher Fleckschussweite und gleicher Ziefernrohrhöhe ergibt das Funktional für das Luftwiderstandsgesetz G1 einen Wert von 0, 0038125 und für den Formfaktor einen Wert von 1, 13100.
Der Wert des Funktionals ist für das Luftwiderstandsgesetz G1 mehr als doppelt so groß, damit
passt das Luftwiderstandsgesetz G7 besser zu den Daten als das Luftwiderstandsgesetz G1, was
die Behauptung von Bryan Litz stützt.
2.10 Erweiterungen des Verfahrens
Das vorgestellte Verfahren wurde bisher nur für Munition für Handfeuerwaffen konzipiert.
Erweiterungen sind vorstellbar:
- Wählt man ein geeignetes Luftwiderstandsgesetz (ein Luftwiderstandsgesetz, welches Artilleriegeschosse adäquat beschreibt), dann eröffnet sich die Möglichkeit, die entsprechenden Parameter
C̃ = (C̃d , C̃L , C̃Np , C̃Mp , C̃Mq , C̃l p , C̃Mα )T
optimal zu schätzen. Damit wäre auch eine potentielle Anwendung für die Artillerie aufgezeigt.
- Ebenso liegt es nahe, mit weiteren Luftwiderstandsgesetzen zu arbeiten - das vorgeschlagene Verfahren ist ja unabhängig von den Gesetzen G1 bis G7.
Literaturverzeichnis
[Carlucci und Jacobson 2008] CARLUCCI, D. E. ; JACOBSON, S. S.: Ballistics Theory and Design
of Guns and Ammunition. CRC Press, 2008
[Kamke 1969]
[Litz 2011]
KAMKE, E.: Differentialgleichungen. Akademische Verlagsgesellschaft, 1969
L ITZ, B.: Applied Ballistics for Long-Range Shooting. Applied Ballistics, LLC, 2011
[McCoy 1999]
M C C OY, R. L.: Modern Exterior Ballistics. Schiffer Military History, 1999
[Morrison 1976]
pany, 1976
M ORRISON, D. F.: Multivariate Statistical Methods. McGraw-Hill Book Com-
[Peelen 2014] P EELEN, J.: Ballistik kleiner Kaliber: der Luftwiderstand. Verlagshaus Monsenstein und Vannerdat OHG Münster, 2014
[Seber 1984]
S EBER, G. A. F.: Multivariate Observations. Wiley, 1984
31
32
LITERATURVERZEICHNIS