Theoretische und numerische Untersuchungen zu den geodätischen Hauptaufgaben

Armin Retterath

Inhalt

• Inhalt
• 1. Theoretische Untersuchungen zu den geodätischen Hauptaufgaben
  • 1.1 Die geometrische Bedeutung der geodätischen Linie
    • 1.1.1 Grundlegende Begriffe
    • 1.1.2 Das Problem der kürzesten Verbindung zweier Raumpunkte
    • 1.1.3 Möglichkeiten zur Lösung der Differentialgleichung der geodätischen Linie
  • 1.2 Das Bewegungsproblem in der Mechanik
    • 1.2.1 Begriffsdefinitionen
    • 1.2.2 Erweiterte Methoden der Mechanik
  • 1.3 Herausstellen der Analogie zwischen dem geometrischen und mechanischen Problem
    • 1.3.1 Ableitung der Differentialgleichung der geodätischen Linie über die Eulersche Gleichung
    • 1.3.2 Herleitung der expliziten Form der Lagrange-Gleichungen zweiter Art
  • 1.4 Eine Möglichkeit der Lösung der geodätischen Hauptaufgaben über die kanonischen Gleichungen der geodätischen Linie
• 2. Numerische Untersuchungen zu den geodätischen Hauptaufgaben
  • 2.1 Numerische Verfahren zur näherungsweisen Bestimmung der Lösung von Differentialgleichungen erster Ordnung als Anfangswertproblem
    • 2.1.1 Beschreibung des Problems
    • 2.1.2 Lösung über Reihenentwicklung
    • 2.1.3 Runge-Kutta-Verfahren
    • 2.1.4 Verallgemeinerung des Runge-Kutta-Verfahrens auf höhere Ordnungen und Systeme von Differentialgleichungen
  • 2.2 Realisierung eines Programms zur Berechnung algebraisch vorliegender Differentialgleichungssysteme nach dem Runge-Kutta-Verfahren 7. Ordnung
  • 2.3 Lösung der 1. geodätischen Hauptaufgabe mit Hilfe der Funktion RK7multi.m
    • 2.3.1 Beschreibung der theoretischen Grundlagen
    • 2.3.2 Aufruf des Lösers über die graphische Benutzeroberfläche guiawp_graphik.m
    • 2.3.3 Ableitung einer empirischen Funktion zur Bestimmung der Zahl der nötigen Iterationen
    • 2.3.4 Diskussion der mit dem Programm guiawp_graphik.m zu erzielenden Ergebnisse
  • 2.4 Lösung der 2. geodätischen Hauptaufgabe nach dem Schießverfahren
    • 2.4.1 Prinzip des Schießverfahrens
    • 2.4.2 Berechnung der Näherungswerte
    • 2.4.3 Praktische Realisierung eines Programms zur Lösung der 2. geodätischen Hauptaufgabe nach dem Schießverfahren
    • 2.4.4 Diskussion der mit dem Programm guirwp_graphik.m zu erzielenden Ergebnisse
  • 2.5 Beschreibung des Programmsystems
    • 2.5.1 Die Darstellung und Abspeicherung der Ephemeriden
    • 2.5.2 Auslesen und Plotten der Datensätze
    • 2.5.3 Berechnung eines Parallelkoordinatensystems in Soldnerscher Anordnung
  • 2.6 Beispiel zur Anwendung des Programmsystems
  • 2.7 Anmerkungen zum Anhang
• Abbildungsverzeichnis
• Literatur

1. Theoretische Untersuchungen zu den geodätischen Hauptaufgaben

1.1 Die geometrische Bedeutung der geodätischen Linie

1.1.1 Grundlegende Begriffe

In diesem Abschnitt werden zunächst die Begriffe Mannigfaltigkeit und Riemannscher Raum mit ihrer differentialgeometrischen Bedeutung kurz erläutert. Außerdem soll die Interpretation der Christoffelsymbole als Elemente des affinen Zusammenhangs in der von einer Einbettung unabhängigen Riemannschen Geometrie herausgestellt werden.

1.1.1.1 Punktraum-Vektorraum

Wenn der Mensch den Begriff des Raumes verwendet, so verbindet er mit diesem einen dreidimensionalen ausgedehnten Raum, welchen er mit seinen Sinnen erfassen kann. Er kann den, seinen Sensoren zugänglichen Punkten des Anschauungsraumes, durch Messen von Richtungen und Abständen eine gegenseitige Lage zuordnen. Es sind nun, um die Eindeutigkeit der Zuordnungen zu gewährleisten, Vereinbarungen über den Maßstab und den Bezugspunkt der Messungen nötig. Dies führt zu einer Einführung von gewissen Berechnungsvorschriften, die dem Anschauungsraum eine einfache Ordnung geben. Um Berechnungen zwischen den Raumpunkten durchführen zu können, müssen die Punkte des Raumes einer mathematischen Beschreibung zugänglich gemacht werden.

Betrachtet man eine Menge von Punkten in einem n-dimensionalen Raum , welche im folgenden mit Mannigfaltigkeit bezeichnet wird und einen Vektorraum , so wird durch die Abbildung von in jeder Punkt einem Vektor zugeordnet. Die Zuordnung setzt eine Ordnung des Raumes durch eine Vereinbarung über ein System voraus, das als Koordinatensystem bezeichnet wird. Ein Koordinatensystem sei durch ein Paar , welches aus einem festen Punkt und einer Basis besteht, definiert. Von diesem ausgehend wird durch die bijektive Abbildung

$$f:\left( O,X\left( x_{1},x_{2},...,x_{n}\right) \right) \; \... ...rray}{c} x_{1}\\ x_{2}\\ \vdots \\ x_{n} \end{array}\right]$$

jedem Punkt bezüglich der Ortsvektor zugeordnet. (REINHARDT/SOEDER 1984) Zuerst wird jedoch vom Begriff des Vektors, als einfachster Art einer Beschreibung der gegenseitigen Lage zweier Punkte ausgegangen. Ein Vektor bezeichnet demnach eine Äquivalenzklasse parallelgleicher Pfeile. Die Pfeile lassen sich als geordnete Punktepaare interpretieren, wobei Angriffspunkte und Zielpunkte bezeichnen. Wichtig für den so eingeführten Begriff des Vektors ist dabei die Unabhängigkeit dieses geometrischen Objektes vom gewählten Koordinatensystem^1.1, d.h. seine Invarianz gegenüber bestimmten Basistransformationen. Um die Ordnung des Raumes aufrechtzuerhalten reicht die Vereinbarung über ein gewähltes Koordinatensystem nicht aus. Es müssen im weiteren Rechenregeln^1.2 für die Elemente des Raumes aufgestellt werden. Dabei gelten die folgenden Gesetze:

1. Summe
   1. 
   2. 
   3. , so daß
   4. , zu jedem so daß
2. gemischtes Produkt
   1. 
   2. 
   3. 
   4. 
3. Skalarprodukt
   1. 
   2. 
   3. 
   4. Wenn und beliebig, dann gilt

Sind nun die Gesetze 1, 2 und 3 erfüllt, so wird der Vektorraum als Euklidscher Vektorraum bezeichnet, gelten nur 1 und 2 so bezeichnet man ihn als Affinen Vektorraum.

1.1.1.2 Tensorielle Darstellung von Funktionen

Da in den folgenden Abschnitten häufig die tensorielle Darstellung Verwendung findet, soll hier zunächst ein kurzer Überblick gegeben werden. Dabei wird zwischen skalaren und vektorwertigen Funktionen unterschieden und deren partielle Ableitungen werden diskutiert. Die Einführung der Christoffelsymbole an dieser Stelle, dient dem Herausstellen der Analogie zur kovarianten Behandlung von Bewegungsgleichungen in der Form der Lagrangegleichungen zweiter Art. Die Darstellung der Funktionen hält sich dabei weitestgehend an die von KLINGBEIL (1984) verwandte Notation. Ein Index welcher nur der Unterscheidung von Objekten dient, wird als kovarianter Index mit einem Punkt als Trennsymbol zu den Summationsindizes angefügt.

1. Ein geometrischer Vektor läßt sich nach der Summenkonvention bezüglich der geradlinigen Basis sowie der krummlinigen Basis durch seine kontravarianten Komponenten wie folgt darstellen
   
   $$\mathbf{x}=x^{i}\mathbf{e}_{i}=u^{i}\mathbf{b}_{i}$$
   
   

2. skalarwertige Funktionen
   1. Für eine Skalarfunktion gilt bezüglich einer festen Basis
      
      $$s=s\left( x^{i}\right)$$
      
      
      Nach der Vereinbarten Notation ergibt sich dann für die Ableitungen nach den ko- bzw. kontravarianten Komponenten
      
      $$s_{,i}=\frac{\partial s}{\partial x^{i}}$$
      
      
      
      
      $$s_{,}^{i}=\frac{\partial s}{\partial x_{i}}$$
      
      
      Hierbei ist zu beachten, daß sich bei Ableitungen eines Skalars nach den im allgemeinen kontravarianten Koordinaten die kovarianten Komponenten eines Vektors (Tensor 1. Stufe) ergeben. Für das totale Differential gilt dann
      
      $$ds=\frac{\partial s}{\partial x^{i}}dx^{i}=s_{,i}dx^{i}$$
      
      

   2. Wird nun eine ortsabhängige Basis eingeführt, so kann am Transformationsverhalten gezeigt werden, daß es sich bei ebenfalls um (kovariante) Komponenten eines Tensors 1. Stufe handelt. Aus diesem Grund wird eine s.g. kovariante Ableitung
      

      $$s_{,i}=s\left\vert _{i}\right.$$ (1.1)

      
      
      eingeführt, welche durch ihre Eigenschaft der Invarianz gegenüber erlaubten Parametertransformationen definiert ist.
3. vektorwertige Funktionen
   1. Eine vektorwertige Funktion wird bezüglich einer festen kartesischen Basis beschrieben durch
      
      $$\mathbf{r}=\mathbf{r}\left( x^{i}\right) =r^{k}\left( x^{i}\right) \mathbf{e}_{k}$$
      
      
      Für ihre Ableitungen nach den kartesischen Koordinaten ergibt sich dann
      
      $$\frac{\partial \mathbf{r}}{\partial x^{i}}=\frac{\partial r^... ...l x^{i}}\mathbf{e}_{k}=r^{k}_{,i}\mathbf{e}_{k}=\mathbf{r}_{,i}$$
      
      
      Und das totale Differential lautet in diesem Fall
      

      $$d\mathbf{r}=\frac{\partial \mathbf{r}}{\partial x^{i}}dx^{i}=\mathbf{r}_{,i}dx^{i}$$ (1.2)

      
      
      mit , wobei die die gemischt ko/kontravarianten Komponenten eines Tensors zweiter Stufe bilden.
   2. Es sei nun eine vektorwertige Funktion der krummlinigen Koordinaten, so folgt analog zu () für das totale Differential
      
      $$d\mathbf{r}=\mathbf{r}_{,k}du^{k}$$
      
      
      Für die Ableitung nach den kartesischen Koordinaten gilt aber aufgrund der Kettenregel
      
      $$\mathbf{r}_{,i}=\mathbf{e}_{i}=\frac{\partial \mathbf{r}}{\partial u^{k}}\frac{\partial u^{k}}{\partial x^{i}}\; ,$$
      
      
      so daß sich die ortsabhängige Basis über
      
      $$\mathbf{b}_{k}=\mathbf{e}_{i}\frac{\partial x^{i}}{\partial u^{k}}$$
      
      
      bzw. die kartesische Basis nach
      
      $$\mathbf{e}_{i}=\mathbf{b}_{k}\frac{\partial u^{k}}{\partial x^{i}}$$
      
      
      bestimmen läßt. Im Gegensatz zu () kann aber gezeigt werden, daß im Falle eines Tensors erster oder höherer Stufe, dessen partielle Ableitung kein Tensor mehr ist. Am Beispiel eines Tensors 1. Stufe gilt
      
      $$r^{i}_{,k}\neq r^{i}\left\vert _{k}\right.$$
      
      
      Dies kommt daher, daß, wenn man bspw. die Ableitung eines Vektors in einem Punkt des krummlinigen Raumes bestimmen will, man Informationen über den Raum selbst, d.h. speziell über die Nachbarschaft des Punktes benötigt. Diese Informationen, die einen Punkt im Raum mit den weiteren, den Raum bildenden Punkten verknüpft, müssen sich aus der Metrik des Raumes ableiten lassen. Es ist nun möglich, aus der geforderten Invarianz der Ableitungen, einen Korrekturterm für die partiellen Ableitungen herzuleiten, so daß eine kovariante Ableitung im krummlinigen Raum gefunden werden kann. Für eine vektorwertige Funktion ergibt sich diese zu
      

      $$r^{i}\left\vert j\right. =r^{i}_{,j}+\Gamma ^{i}_{jk}r^{k}$$ (1.3)

      
      
      wobei die ^1.3 selbst keine Tensoreigenschaft besitzen. Mithilfe dieser Ableitung wird ein Absolutes Differential über
      

      $$Dr^{i}=r^{i}\left\vert _{j}\right. du^{j}$$ (1.4)

      
      
      definiert. Man kann zeigen, daß die Komponenten des vollständigen Differentials bezüglich der krummlinigen Basis bildet, d.h. es gilt dann
      
      $$d\mathbf{r}=Dr^{i}\mathbf{b}_{i}$$
      
      

1.1.1.3 Die Verallgemeinerung des Raumbegriffs

Nach dem theorema egregium von Gauß ist bekannt, daß die Krümmung einer Fläche allein vom Metriktensor abhängt. Das bedeutet, am Beispiel einer in einem dreidimensionalen kartesischen Raum eingebetteten Fläche, daß das Krümmungsmaß dieser Fläche ohne Informationen bezüglich einer Einbettung, d.h. allein aus Messungen auf der Fläche bestimmbar ist. Dies führte zu einer vom Einbettungsraum unabhängigen Flächentheorie, deren Grundlagen in den folgenden Abschnitten kurz aufgeführt sind.

1.1.1.3.1 Mannigfaltigkeit und Affiner Tangentialraum

Der von Riemann eingeführte verallgemeinerte Raumbegriff ist die Mannigfaltigkeit. Die Definition einer solchen erfolgt über die Transformationsgleichungen zweier Koordinatensysteme

$$u^{i} = u^{i}(\overline{u}^{1},...,\overline{u}^{n})$$ (1.5)

$$\overline{u}^{i} = \overline{u}^{i}(u^{1},...,u^{n})$$ (1.6)

Sind die Transformationsgleichungen stetig differenzierbar, so handelt es sich bei der durch (1.5) und () gegebenen Mannigfaltigkeit um eine n-dimensionale differenzierbare Mannigfaltigkeit. In jedem ihrer Punkt, können die totalen Differentiale der Koordinaten über und umgekehrt bestimmt werden. Mit diesem Schritt erhält man die linearen Transformationen der Koordinatendifferentiale über die für einen Tensor erster Stufe geltenden Transformationskoeffizienten.

Abbildung: Geometrische Interpretation des Affinen Tangentialraums

Mit Hilfe des Tensorcharakters der Koordinatendifferentiale läßt sich diesen in einem Punkt ein Vektor zuordnen.

$$du^{i}\longleftrightarrow d\mathbf{u}=du^{i}\mathbf{b}_{i}$$ (1.7)

Ordnet man einem zweiten Tensor (bzw. Koordinatendifferential) ebenfalls einen Vektor zu, so kann man durch Erklärung einer Vektorsumme

$$ds^{i}=du^{i}+dv^{i}\longleftrightarrow d\mathbf{s}=d\mathbf{u}+d\mathbf{v}=\left( du^{i}+dv^{i}\right) \mathbf{b}_{_{i}}$$

und eines Produktes mit einem Skalar

$$p^{i}=\alpha t^{i}\longleftrightarrow \mathbf{p}=\alpha \mathbf{t}=\alpha t\mathbf{b}_{i}$$

der Mannigfaltigkeit im Punkt einen Vektorraum anhängen. Es wird vereinbart, daß für die Elemente dieses Raumes, die Gesetze (1) und (2) aus Abschnitt 1.1.1.1 erfüllt sein sollen. Da kein Skalarprodukt erklärt ist, läßt sich nur ein affiner Vektorraum anhängen. Abb. stellt eine geometrische Interpretation dieses Vektorraumes dar. Es sei eine Raumkurve gegeben, außerdem existiert ein Punkt auf der Kurve, für den diese den Parameterwert annimmt. Ein benachbarter Punkt auf der Kurve habe nun den Parameterwert und die Koordinaten . Dem Koordinatendifferential wird gemäß () der Vektor zugeordnet der Endpunkt des Vektors liegt nicht in . Außerdem entspricht dieser Vektor nur im Differentiellen dem durch definierten Tangentenvektor an die Kurve.

Dies führt zu dem Problem, daß einem Tensor in Form eines Koordinatendifferentials der Mannigfaltigkeit, unterschiedliche Vektoren bezüglich der affinen Tangentialräume zugeordnet werden. Die Tangentialräume in zwei verschiedenen Punkten sind in diesem Fall nicht zusammenhängend.

Um zu einem geeigneten Zusammenhang zwischen den Räumen zu gelangen, ist es notwendig die Parallelität von Vektoren im verallgemeinerten Raum auf andere Weise zu definieren. Für eine gekrümmte Fläche, als Beispiel für eine zweidimensionale Mannigfaltigkeit, läßt sich das Prinzip einfach verdeutlichen. Der affine Tangentialraum in einem Punkt kann in einem solchen Fall mit der Tangentialebene an die Fläche in diesem Punkt identifiziert werden. Versucht man einen Vektor in der Tangentialebene in , parallel in einen benachbarten Punkt zu verschieben, so ist dies nach der euklidschen Definition von Parallelität, aufgrund der Ortsabhängigkeit der Basisvektoren, nicht möglich. Eine 'parallelste' Verschiebung kann jedoch auf diese Weise definiert werden, daß man räumlich parallel (im euklidschen Sinn) nach versetzt und diesen im Anschluß daran in die Tangentialebene in projeziert. Diese Parallelverschiebung gilt natürlich nur für eine differentiell kleine Entfernung zwischen und . Es kann gezeigt werden, daß diese Art der Verschiebung, der s.g. Parallelverschiebung eines Tensors 1. Stufe nach Levi Civita entspricht.

Ein Tensor 1. Stufe einer Mannigfaltigkeit heißt parallelverschoben nach Levi Civita, entlang einer Kurve mit dem Parameter , wenn folgende Differentialgleichungen gelten

$$\frac{Dv^{i}}{dt}=0$$

(d.h. sich sein absolutes Differential entlang der Kurve nicht ändert)

Mit () gilt die hierzu äquivalente Darstellung

$$v^{i}\left\vert _{j}\right. \dot{u}^{j}=0$$ (1.8)

Über () sind die Tangentialräume in den Punkten und miteinander verknüpft. Wie leicht zu erkennen ist, entspricht die Einführung der kovarianten Ableitung in () der Definition eines verallgemeinerten Paralleltransportes. Beide Operationen verknüpfen die affinen Tangentialräume und ermöglichen auf diese Weise die mathematische Beschreibung von Objekten in der Mannigfaltigkeit. Da diese Verknüpfung, nach Einführung einer Metrik, wie in () zu sehen, über die Christoffelsymbole geschieht, werden diese auch als Objekte des affinen Zusammenhangs im Riemannschen Raum bezeichnet.

1.1.1.3.2 Der Riemannsche Raum

Der Riemannsche Raum entsteht durch die Ausstattung einer n-dimensionalen Mannigfaltigkeit mit einem symmetrischen Tensorfeld , so daß die Länge einer Kurve zwischen den Parameterwerten und gegeben ist durch

$$\displaystyle s=\int _{t_{0}}^{t_{1}}\sqrt{g_{ij}\dot{v}^{i}\dot{v}^{j}}$$

Durch die Einführung dieses Feldes, kann den Tangentialräumen über die Definition eines Skalarproduktes

$$g_{ij}u^{i}v^{j}=\mathbf{u}\cdot \mathbf{v}$$

eine Euklidsche Metrik gegeben werden. Aus den partiellen Ableitungen des Metriktensors lassen sich dann, wie im folgenden noch gezeigt wird, die Christoffelsymbole der zweiten Art bestimmen. Sie sind die einzigen symmetrischen Objekte des Zusammenhangs, die die Länge eines parallelverschobenen Vektors unverändert lassen.

1.1.2 Das Problem der kürzesten Verbindung zweier Raumpunkte

Eine Aufgabe der Geodäsie besteht in der Anlage von Lagefestpunktfeldern. Die Bezugsfläche ist ein festgelegtes Rotationsellipsoid und die Lagekoordinaten eines geodätischen Punktes werden, nach einer geeigneten Projektion desselben auf das Bezugsellipsoid, in einem System Gaußscher Flächenparameter definiert. Die numerischen Werte dieser Parameter entsprechen somit den Lagekoordinaten der Punkte. Die am häufigsten Verwendung findenden Parametersätze können Tabelle entnommen werden.

Tabelle: Flächenparametersysteme auf dem Rotationsellipsoid

Bezeichnung
geographische Koordinaten	(ellipsoidische Breite)	(geographische Länge)
geodätische Polarkoordinaten	(Länge der geodätischen Linie)	(Azimut der geodätischen Linie)
geodätische Parallelkoordinaten:
nach Soldner	(Abzissenwert)	(Ordinatenwert)
Gauß'sche isotherme Flächenkoordinaten	(Gauß'sche Abzisse)	(Gauß'sche Ordinate)

Da meist mehrere solcher Systeme nebeneinander Verwendung finden, ist es notwendig, Flächenparametertransformationen zwischen den einzelnen Systemen durchzuführen. Eine Möglichkeit dieser Transformation beruht auf der Differentialgleichung der geodätischen Linie, d.h der i.A. eindeutigen, bezüglich eines Linienelementes kürzesten, Verbindung zweier Punkte auf einer Fläche. Intuitiv ergeben sich die Bedingungen für das Vorliegen einer kürzesten Verbindung zweier Punkte entsprechend Abb..

Abbildung: Verlauf der geodätischen Linie auf einer gekrümmten Fläche

Im kartesischen Koordinatenraum wird die kürzeste Verbindung zweier Punkte mit den Ortsvektoren und durch einen Vektor realisiert. Sind diese beiden Punkte nun Elemente einer Fläche, so kann man direkt folgern, daß die kürzeste Verbindung der Punkte über eine Flächenkurve durch die Projektion der geradlinigen Verbindung auf die Fläche erhalten wird. Diese Projektion erfolgt so, daß jedem einzelnen Punkt von ein Lotfußpunkt auf der Fläche zugeordnet wird. Die Lotfußpunkte bilden dann die geodätische Linie zwischen den beiden Punkten. Dies gilt streng genommen nur für stetige Flächen. Außerdem gibt es auf beliebigen Flächen keine eindeutige Zuordnung, da es mehrere geodätische Linien zwischen zwei Flächenpunkten geben kann. Weiterhin kann eine Flächenkurve nur dann die kürzeste Verbindung zwischen 2 Punkten darstellen, wenn sie bezüglich der Fläche selbst keine Krümmung aufweist, d.h. ihre Krümmung in jedem Punkt mit der Krümmumg der Fläche identisch ist. Damit folgt, daß der Hauptnormalenvektor der geodätischen Linie in jedem Punkt dem Normalenvektor der Fläche entspricht. In den beiden folgenden Abschnitten wird nun die Bestimmungsleichung für den Verlauf der geodätischen Linie auf verschiedene Weise hergeleitet. In der ersten Herleitung wird der Begriff der Krümmung einer Verbindunglinie zweier Punkte auf gekrümmte Räume verallgemeinert und damit die Differentialgleichung der geodätischen Linie als Bedingung für eine bezüglich des gekrümmten Raumes krümmungslose Verbindung abgeleitet. In der zweiten Herleitung sei die Analogie zwischen einer Geraden im kartesischen Raum und der geodätischen Linie auf eine andere Weise herausgestellt. Die gesuchte Differentialgleichung ergibt sich als Bedingung für den Paralleltransport eines Tangenteneinheitsvektors. Anders ausgedrückt, ist die geodätische Linie die Linie, analog zur Geraden im kartesischen Raum, entlang welcher der Tangenteneinheitsvektor eine konstante Richtung besitzt.

1.1.2.1 Herleitung der Differentialgleichung der geodätischen Linie über die geradlinige Verbindung zweier Punkte

Für die Bestimmung der kürzesten Verbindung zweier Örter im dreidimensionalen Raum gilt allgemein folgender Zusammenhang

$$\frac{d^{2}\mathbf{x}}{ds^{2}}=0$$ (1.9)

Diese Gleichung beinhaltet, daß die Krümmung der durch definierten Kurve Null ist. Das Linienelement

$$ds^{2}=d\left( x^{1}\right) ^{2}+d\left( x^{1}\right) ^{2}+d\left( x^{1}\right) ^{2}$$

beschreibt in einem kartesischen Koordinatensystem die Rechenvorschrift zur Bestimmung des Abstandes. Die Lösung von () ergibt sich nach zweifacher Integration zu

$$\mathbf{x}\left( s\right) =\mathbf{a}\cdot s+\mathbf{b}$$

der Punktrichtungsform der Geradendarstellung im mit dem Parameter . Um zu einer allgemein gültigen Darstellung der Differentialgleichung () für beliebige n-dimensionale Koordinatensysteme zu gelangen, muß zweimal integriert werden. Für gilt nach dem totalen Differential

$$\frac{d\mathbf{x}}{ds}=\frac{\partial \mathbf{x}}{\partial u... ...^{n}\frac{\partial \mathbf{x}}{\partial u^{i}}\frac{du^{i}}{ds}$$

Nach der Summenkonvention ergibt sich

$$\frac{d\mathbf{x}}{ds}=\frac{\partial \mathbf{x}}{\partial u^{i}}\frac{du^{i}}{ds}$$

und nach der zweiten Differentiation

$$\frac{d^{2}\mathbf{x}}{ds^{2}}=\frac{d}{ds}\left( \frac{\par... ...\partial \mathbf{x}}{\partial u^{i}}\frac{d^{2}u^{i}}{ds^{2}}=0$$

Die Multiplikation mit führt zu

$$\frac{\partial ^{2}\mathbf{x}}{\partial u^{i}\partial u^{j}}... ...\partial \mathbf{x}}{\partial u^{k}}\frac{d^{2}u^{i}}{ds^{2}}=0$$

Hierin kann nun durch den Metriktensor ersetzt werden. Um die totalen Differentiale des Ortsvektors und der Basisvektoren des Dreibeins in die Basis zurückzuschreiben, d.h. vom einbettenden Raum unabhängig zu werden, muß der Term

$$\frac{\partial ^{2}\mathbf{x}}{\partial u^{i}\partial u^{j}}\frac{\partial \mathbf{x}}{\partial u^{k}}$$ (1.10)

durch eine nur von den krummlinigen Koordinaten abhängigen Größe ersetzt werden. Mit der Gleichung

$$\frac{\partial g_{ik}}{\partial u^{j}}=\frac{\partial \left(... ...i}}\frac{\partial ^{2}\mathbf{x}}{\partial u^{j}\partial u^{k}}$$

und zwei weiteren, welche durch zyklische Vertauschung aus der obigen abzuleiten sind, erhält man nach einfacher Rechnung den folgenden Zusammenhang für ()

$$\frac{1}{2}\left( \frac{\partial g_{ik}}{\partial u^{j}}+\fr... ...\frac{\partial \mathbf{x}}{\partial u^{k}}\equiv \Gamma _{kij}$$ (1.11)

Die mit bezeichneten dreidimensionalen Felder werden als Christoffelsymbole erster Art bezeichnet. Mit ihnen ergibt sich für die Gleichung der kürzesten Verbindung zweier Punkte

$$g_{ik}\frac{d^{2}u^{i}}{ds^{2}}+\Gamma _{kij}\frac{du^{i}}{ds}\frac{du^{j}}{ds}=0$$

Über die Multiplikation mit den kontravarianten Metrikkoeffizienten und die anschließende Einführung der Christoffelsymbole der zweiten Art

$$\Gamma _{ij}^{k}\equiv g^{kl}\Gamma _{kij}=\frac{1}{2}g^{kl}... ...\partial u^{i}}-\frac{\partial g_{ij}}{\partial u^{k}}\right)$$ (1.12)

läßt sich diese weiter vereinfachen zu

$$\frac{d^{2}u^{l}}{ds^{2}}+\Gamma _{jk}^{l}\frac{du^{j}}{ds}\frac{du^{k}}{ds}=0$$ (1.13)

1.1.2.2 Herleitung der Differentialgleichung der geodätischen Linie über Parallelverschiebung von Vektoren

Allgemein bedeutet eine parallele Verschiebung eines Vektors, daß Richtung und Länge desselben von der Verschiebung unberührt bleiben. Im euklidschen Raum bleibt ein Vektor und somit dessen kartesische Komponenten durch Paralleltransport unverändert. Um diese Unveränderlichkeit mathematisch zu beschreiben, setzt man das Differential des Vektors zu 0:

$$d\mathbf{v}=0$$

in allgemeinen Koordinaten gilt für diesen Vektor .

Für die Basisvektoren gilt:

$$\mathbf{b}_{i}=\frac{\partial \mathbf{x}}{\partial u^{i}}$$ (1.14)

Das totale Differential ergibt sich damit nach der Produktregel zu

$$d\mathbf{v}=\mathbf{b}_{i}dv^{i}+v^{i}d\mathbf{b}_{i}=0$$

Die Multiplikation mit den Basisvektoren und die Einführung des Metriktensors führt dann zu

$$g_{ji}dv^{i}+\mathbf{b}_{j}d\mathbf{b}_{i}v^{i}=0$$

Mit dem totalen Differential der Basisvektoren und () folgt

$$\mathbf{b}_{j}d\mathbf{b}_{i}=\frac{\partial \mathbf{x}}{\pa... ... u^{i}\partial u^{k}}\frac{\partial \mathbf{x}}{\partial u^{j}}$$

Die rechte Seite der Gleichung entspricht () und kann dementsprechend durch die Christoffelsymbole der ersten Art ersetzt werden. Dies führt analog zum vorigen Kapitel, nach Einführung der Christoffelsymbolen zweiter Art, zur allgemeinen Form des Paralleltransportes eines kontravarianten Vektors:

$$dv^{l}+\Gamma _{jk}^{l}v^{j}du^{k}=0$$ (1.15)

Ist nun ein Tangentenvektor an eine Raumkurve durch gegeben, dann folgt für dessen Parallelverschiebung, wenn

$$\frac{dv^{l}}{ds}+\Gamma _{jk}^{l}v^{j}\frac{du^{k}}{ds}=0$$

das Äquivalent zu () bezüglich des Paralleltransportes entlang einer Kurve beschreibt direkt die Differentialgleichung der geodätischen Linie:

$$\frac{d^{2}u^{l}}{ds^{2}}+\Gamma _{jk}^{l}\frac{du^{j}}{ds}\frac{du^{k}}{ds}=0$$

1.1.3 Möglichkeiten zur Lösung der Differentialgleichung der geodätischen Linie

1.1.3.1 Lösung über Reihenentwicklung in Legendre Reihen

Für beliebige orthogonale Flächenparametersysteme ergibt sich die Gleichung der geodätischen Linie vereinfacht aus folgenden 3 Differentialgleichungen 1. Ordnung:

$$\frac{du}{ds} = \frac{1}{\sqrt{g_{11}}}\cos \left( \alpha \right)$$

$$\frac{dv}{ds} = \frac{1}{\sqrt{g_{22}}}\cos \left( \alpha \right)$$ (1.16)

$$\frac{d\alpha }{ds} = \frac{g_{11,2}}{2g_{11}\cdot \sqrt{g_{22}}}\cos \left( \alpha \right) -\frac{g_{22,1}}{2g_{22}\sqrt{g_{22}}}\sin \left( \alpha \right)$$

Damit lassen sich die die Koordinaten eines beliebigen Punktes auf der geodätischen Linie, die im Pol des lokalen geodätischen Polarkoordinatensystems den Richtungswinkel bezüglich der Parameterlinie besitzt, in Abhängigkeit der Länge der geodätischen Linie aus der Lösung von () über eine Taylorentwicklung nach

$$\begin{eqnarray*} u^{n}=u\left( s\right) & = & u\left( 0\right) +\sum _{k=1}^{\i... ...}s^{k}\cdot \frac{d^{k}\alpha }{ds^{k}}\left\vert _{s=0}\right. \end{eqnarray*}$$

bestimmen. Da die Ableitungen von und nach der Bogenlänge Funktionen sind, ergeben sich die Ableitungen höherer Ordnung durch schrittweises Anwenden der Kettenregel am Beispiel für nach:

$$\frac{d^{2}u}{ds^{2}}=\frac{d}{ds}\left( f\left( u\left( s\r... ...artial f}{\partial \alpha }\frac{\partial \alpha }{\partial s}$$ (1.17)

Für beliebige Gauß'sche Flächenparametersysteme hängen die oben genannten Differentialquotienten linear von und ab. Somit kann die Lösung der 1. geodätischen Hauptaufgabe als Reihenentwicklung dargestellt werden:

$$\begin{eqnarray*} \triangle u=u^{n}-u^{0} & = & A_{10}s\cos \left( \alpha ^{0}\r... ...) \right) ^{i}\left( \sin \left( \alpha ^{0}\right) \right) ^{k} \end{eqnarray*}$$

Nach Einführung der Riemannschen Normalkoordinaten und , ergeben sich die s.g. Legendreschen Reihen

$$\begin{eqnarray*} \triangle u & = & A_{10}\cdot \chi _{c}+\sum _{i+k=2}^{\infty ... ...k}\cdot \left( \chi _{c}\right) ^{i}\left( \chi _{s}\right) ^{k} \end{eqnarray*}$$

Der Wert für ist ebenfalls auf diese Weise bestimmbar, jedoch existiert beim Beispiel der geographischen Koordinaten auf dem Rotationsellipsoid eine strenge Rechenvorschrift, welche sich aus dem Satz von Clairaut ableiten läßt und der Reihenentwicklung aus numerischen Gesichtspunkten vorzuziehen ist:

$$\alpha \left( u^{0},u^{n},\alpha ^{0}\right) =\arcsin \left(... ... ^{2}u^{0}}{1+e'+\tan ^{2}u^{n}}\cdot \sin \alpha ^{0}}\right)$$

(wobei der geographischen Breite entspricht) Der Vollständigkeit halber sei hier noch die Lösung der 2. geodätischen Hauptaufgabe, welche sich aus der Umkehrung der Reihen mittels Reiheninversion ergibt, angegeben:

$$\begin{eqnarray*} \chi _{c} & = & \overline{A}_{10}\triangle u+\sum ^{\infty }_{... ...{ik}\left( \triangle u\right) ^{i}\left( \triangle v\right) ^{k} \end{eqnarray*}$$

Die benötigten Koeffizienten ergeben sich nach () sukzessive aus den Ableitungen. Die ersten Koeffizienten^1.4 seien in den Tabellen () und () kurz dargestellt.

Tabelle: Koeffizienten der Legendreschen Reihe für 1. geod. Hauptaufgabe

	Berechnungsvorschrift

Tabelle: Koeffizienten der Legendreschen Reihe für 2. geod. Hauptaufgabe

	Berechnungsvorschrift

1.1.3.2 Lösung über numerische Integration des Differentialgleichungssystems

Die Entwicklung immer schnellerer Rechenanlagen in den letzten 20 Jahren führte dazu, daß die klassische Lösung der geodätischen Hauptaufgaben über die im vorigen Abschnitt beschriebene Reihenentwicklung mehr und mehr in den Hintergrund trat. Zwar weist sie bezüglich der Rechengeschwindigkeit Vorteile gegenüber numerischen Verfahren auf, doch in Anbetracht der heutigen Prozessorleistungen von bis zu 1GHz, relativiert sich dies. Die numerische Integration kann nach Aufstellen der Differentialgleichung der geodätischen Linie für jede beliebig gekrümmten Fläche erfolgen, während die Lösung über eine Reihenentwicklung zunächst einen großen Aufwandt an algebraischen Umformungen erfordert, der mit steigender Anforderung an die Genauigkeit auch von heute erhältlichen Computer Algebra Systemen (CAS) noch nicht in befriedigender Form durchgeführt werden kann.

Es steht eine Vielzahl verschiedener Integratoren^1.5 zur numerischen Quadratur von Differentialgleichungen zur Verfügung. Diese unterscheiden sich in der Zahl ihrer Funktionsaufrufe, sowie in der Zuverlässigkeit. Einige arbeiten mit einer Fehlerschätzung, so daß sie, bei Vorgabe einer bestimmten Genauigkeit, ohne großen Zeitverlust, die Zahl der nötigen Iterationen zur Einhaltung der Vorgabe selbst abschätzen. An dieser Stelle sei auf das zweite Kapitel verwiesen, das sich intensiver mit dieser Problematik auseinandersetzt.

1.2 Das Bewegungsproblem in der Mechanik

In diesem zweiten Abschnitt des ersten Teils soll nun auf die Möglichkeiten der Formulierung des Bewegungsproblems in der Mechanik eingegangen werden. Hierzu seien zunächst kurz einige Begriffsdefinitionen aufgeführt. Im Anschluß daran wird auf die alternativen Darstellungsformen der Bewegungsgleichungen nach Lagrange und Hamilton eingegangen.

1.2.1 Begriffsdefinitionen

• Geschwindigkeit

Ein Teilchen bewege sich entlang einer Bahn , dann ist seine Position zur Zeit im Punkt durch den Ortsvektor gekennzeichnet. Als Momentangeschwindigkeit des Teilchens im Punkt definiert man . Ist der Ortsvektor gegeben durch , so folgt für die Geschwindigkeit Der Betrag ergibt sich folglich zu .

• Beschleunigung

Als Beschleunigung eines Teilchens im Punkt wird definiert .

• Impuls

$$\mathbf{p}=m\cdot \mathbf{v}=m\cdot \frac{dx^{i}}{dt}\mathbf{e}_{i}$$

• Kraft

Greift an einem Körper der Masse eine äußere Kraft an und bewegt sich dieser infolge der Kraftwirkung mit der Geschwindigkeit , dann gilt , bei zeitunabhängiger Masse gilt dann das 2. Newtonsche Axiom .

• Arbeit

Bewirkt eine Kraft eine Verschiebung des Massepunktes um , dann leistet sie die Arbeit . Die gesamte Arbeit die entlang eines Weges von einem Vektorfeld bei der Bewegung eines Massepunktes verrichtet wird ergibt sich aus dem Kurvenintegral

$$\displaystyle W=\int ^{P_{2}}_{P_{1}}\mathbf{F}\cdot d\mathb... ...t ^{\mathbf{x}_{2}}_{\mathbf{x}_{1}}\mathbf{F}\cdot d\mathbf{x}$$

• Kinetische Energie

Die gesamte Arbeit bei der Bewegung eines Masseteilchens von nach , wobei die Geschindigkeit zum Zeitpunkt und die zum Zeitpunkt ist, ergibt sich nach

$$W=\int ^{\mathbf{x}_{2}}_{\mathbf{x}_{1}}\mathbf{F}\cdot d\mathbf{x}=\frac{1}{2}mv_{2}^{2}-\frac{1}{2}mv_{1}^{2}=T_{2}-T_{1}$$

bezeichnet hier die kinetische Energie des Teilchens zu den jeweiligen Zeitpunkten .

• Konservatives Kraftfeld

Unter der Annahme, es existiere eine Skalarfunktion , so daß gilt , beträgt die bei der Bewegung eines Massepunktes geleistete Arbeit . In diesem Fall ist die geleistete Arbeit unabhängig vom zurückgelegten Weg, daher wird ein solches Kraftfeld, das sich aus einem Potential ableiten läßt als konservatives Kraftfeld bezeichnet.

Folgende Bedingungen müssen erfüllt sein, damit ein konservatives Kraftfeld vorliegt:

1. 
2. 

• Potentielle Energie oder Potential

Der Skalar wird als potentielle Energie oder als Potential eines Teilchens in einem konservativem Kraftfeld bezeichnet. Das Potential ist nur bis auf eine additive Konstante festgelegt. Unter der Annahme für , gilt für dieses

$$V=-\int _{\mathbf{x}_{0}}^{\mathbf{x}}\mathbf{F}\cdot d\mathbf{x}$$

• Energieerhaltung

Für ein konservatives Kraftfeld gilt . Die Summe aus kinetischer Energie und potentieller Energie wird als Gesamtenergie bezeichnet. In einem konservativem Kraftfeld ist die Gesamtenergie .

1.2.2 Erweiterte Methoden der Mechanik

Im Allgemeinen wird bei der Formulierung von Bewegungsproblemen hauptsächlich auf die Newton/Eulersche Bewegungsgleichung

$$m\ddot{\mathbf{x}}=\mathbf{F}\left( \alpha ,\mathbf{x}\left( t\right) ,\dot{\mathbf{x}}\left( t\right) ,t\right)$$ (1.18)

zurückgegriffen. Wird hier das Kraftgesetz als quantitativer Einfluß der Umgebung auf einen Massepunkt vorgegeben, so liefert () eine vektorielle Differentialgleichung 2. Ordnung zur Bestimmung der Bahn dieses Masseteilchens. Die Herausstellung der Analogie zu den geodätischen Hauptaufgaben erfordert aber eine andere Betrachtungsweise des Bewegungsproblems, so daß in den folgenden Abschnitten zunächst drei weitere Methoden zur Lösung von Bewegungsproblemen dargestellt werden.

1.2.2.1 Lösung von Bewegungsproblemen über die Lagrangegleichungen erster Art

Innere und äußere Kräfte, die auf ein System von Massepunkten wirken, schränken die Bewegungsmöglichkeiten der einzelnen Teilchen meist nicht ein. Je nach Wahl der 6 Anfangswerte sind die Massepunkte in einem solchen Fall in der Lage, jeden Punkt des Raumes zu erreichen. Die meisten Bewegungen verlaufen jedoch nicht frei, d.h. die Massepunkte sind in ihrer Bewegung durch gewisse Nebenbedingungen eingeschränkt. Diese Nebenbedingungen sind zunächst geometrischer Art, bspw. kann als Nebenbedingung für die Bewegung eines starren Pendels, zum einen die feste Schwingebene, zum anderen der kugelförmige mögliche Aufenthaltsraum des Pendelkörpers festgelegt werden. Die Bewegung erfolgt dann entlang des Schnittes der beiden geometrischen Körper. Man kann nun über die Formulierung der s.g. Lagrangegleichung erster Art die Bewegungsgleichung lösen, indem man virtuelle Zwangskräfte , d.h. Kräfte einführt, deren Wirkungslinien senkrecht zur einschränkenden Geometrie stehen und somit keine Arbeit verrichten.

Dies ist äquivalent zum s.g D'Alembertschen Prinzip^1.6:

$$\sum _{k}\left( m_{k}\ddot{\mathbf{x}}_{.k}-\mathbf{F}_{.k}\right) \cdot \delta \mathbf{x}_{.k}=0$$ (1.19)

Mit den Zwangskräften .

ist dabei so zu interpretieren, daß sich die eingeprägte Kraft mit der Trägheitskraft im Gleichgewicht befindet.

Die Zwangskräfte dienen nur der geometrischen Festlegung der Bewegung und tragen, da sie keine Arbeit verrichten, nicht zur Bewegung selbst bei^1.7. Allgemein lassen sich die Lagrangegleichungen erster Art in folgender Form darstellen

$$m_{k}\ddot{\mathbf{x}}_{.k}=\mathbf{F}_{.k}+\mathbf{Z}_{.k}$$

$$\Leftrightarrow \: m_{k}\ddot{\mathbf{x}}_{.k}=\mathbf{F}_{.k}+\sum _{i=1}^{r}\lambda _{i}\nabla g_{.i}$$ (1.20)

Die Zwangskräfte werden, wie bereits erwähnt, derart eingeführt, daß sie orthogonal zur bewegungseinschränkenden Geometrie stehen, und somit in () dargestellt werden können als Vektoren

$$\mathbf{Z}_{i}=\lambda _{i}\nabla g_{.i}$$

in Richtung der Flächennormalen der restriktiven geometrischen Figuren, welche in der impliziten Darstellung

$$g_{.i}\left( \mathbf{x}_{.k},t\right) =0$$ (1.21)

gegeben sind.

Abbildung: Körper gleitet unter Einfluß der Schwerkraft eine Ebene herab

Beispiel:Ein Körper gleitet unter Einfluß der Schwerkraft eine schiefe Ebene hinunter:

Wie in Abb. zu sehen, bildet die Schnittgerade der schiefen Ebene mit der -Ebene den Winkel mit der -Achse. Zum einen gilt daher zu jedem Zeitpunkt

$$g_{1}\left( \mathbf{x},t\right) =x(t)\cdot \cos \left( \alpha \right) -z\left( t\right) \cdot \sin \left( \alpha \right) =0$$

zum anderen

$$g_{2}\left( \mathbf{x},t\right) =y=0$$

Für () ergibt sich in diesem Fall

$$m\cdot \ddot{\mathbf{x}}=\mathbf{F}+\lambda _{1}\nabla g_{.1}+\lambda _{2}\nabla g_{.2}$$ (1.22)

Um die Multiplikatoren zu bestimmen, kann man, da die Nebenbedingungen nicht explizit von der Zeit abhängen, d.h. gilt, man in diesem Fall die 2. Zeitableitungen der Nebenbedingungen über ermitteln. Es folgt dann

$$\frac{d^{2}g_{1}}{dt^{2}}= \ddot{x}\sin \left( \alpha \right) -\ddot{z}\cos \left( \alpha \right) =0$$ (1.23)

$$\frac{d^{2}g_{i}}{dt^{2}}= \ddot{y} =0$$ (1.24)

Aus () ergibt sich damit

$$m\cdot \left( \begin{array}{c} \ddot{x}\\ \ddot{y}\\ \ddot... ...mbda _{1}\cdot \cos \left( \alpha \right) \end{array}\right)$$ (1.25)

Nach Einsetzen der zweiten Zeitableitungen der Koordinaten aus () in () entstehen die folgenden 3 Bestimmungsgleichungen für die Lagrangemultiplikatoren .

$$\begin{eqnarray*} m\cdot \frac{\ddot{z}\cos \left( \alpha \right) }{\sin \left( ... ...lpha \right) } & = & -mg+\lambda _{1}\cos \left( \alpha \right) \end{eqnarray*}$$

$$\Longrightarrow \begin{array}{ccc} \lambda _{1} & = & -mg\cdot \cos \left( \alpha \right) \\ \lambda _{2} & = & 0 \end{array}$$

Die Ergebnisse werden nun in die allgemeine Lagrangegleichung 1. Art () eingesetzt und es ergibt sich

$$m\cdot \left( \begin{array}{c} \ddot{x}\\ \ddot{y}\\ \ddot... ...) \\ 0\\ -\sin ^{2}\left( \alpha \right) \end{array}\right)$$

Diese Gleichung läßt sich direkt integrieren. Nach der geeigneten Wahl der Konstanten folgt die Bewegungsgleichung:

$$\mathbf{x}\left( t\right) =\left( \begin{array}{c} s\left( t... ...\left( t\right) \sin \left( \alpha \right) \end{array}\right)$$

mit

$$s\left( t\right) =-\frac{1}{2}g\cdot t^{2}\cdot \sin \left( \alpha \right) +v_{0}\cdot t+s_{0}$$

Die Weiterentwicklung dieses Lösungsprinzips führt dann zu den Lagrangegleichungen zweiter Art, sowie zum Hamiltonformalismus, die in den folgenden Abschnitten eingehender behandelt werden.

1.2.2.2 Lösung des Bewegungsproblems über die Lagrangegleichung zweiter Art

Ein Ansatz, der ebenfalls auf Lagrange zurück geht, versucht die Bewegungsgleichung zu lösen, indem der Lösungsraum durch die Einführung einer anderen Geometrie eingeschränkt wird. Die Einführung der Geometrie erfolgt über die Definition eines neuen Koordinatensystems, dessen Zahl unabhängiger Koordinaten eben der Zahl der Freiheitsgrade des Systems von Massepunkten entspricht, für welche die Bewegungsgleichungen zu lösen sind. Die Zahl der Freiheitsgrade eines Systems entspricht dabei der Zahl der Koordinaten, die nötig sind, um die Ortsfunktionen der Massepunkte des Systems zu beschreiben. Über das zu wählende Koordinatensystem kann frei verfügt werden, jedoch bestimmt die Art des Systems die Komplexität der Lösung des Bewegungsproblems. Die Zahl der Freiheitsgrade verringert sich nun mit der Zahl der gegebenen Nebenbedingungen. Die minimal zur Beschreibung des Problems notwendigen Koordinaten werden als generalisierte Koordinaten^1.8 bezeichnet. Die Bewegung muß in diesen Koordinaten zwangsfrei verlaufen. Ein System von Massepunkten in einem 3-D Raum wird durch kartesische Koordinaten für und beschrieben. Sei die Zahl der Randbedingungen, so ergibt sich die Zahl der Freiheitsgrade des Systems nach . Hieraus folgt, daß im Falle die Bewegung des Punktsystems durch die Bedingungsgleichungen eindeutig bestimmt ist.

Der Zusammenhang zwischen den rechtwinkligen kartesischen und den generalisierten Koordinaten eines Massepunktes läßt sich wie folgt darstellen

$$\begin{eqnarray*} x_{1.i} & = & x_{1.i}\left( q^{1},...,q^{f},t\right) \\ x_{2.... ...right) \\ x_{3.i} & = & x_{3.i}\left( q^{1},...,q^{f},t\right) \end{eqnarray*}$$

Bei diesen Gleichungen ist zu beachten, daß die Transformationsformeln stetig sind und ebenfalls stetige Ableitungen besitzen^1.9. Man klassifiziert^1.10 nun verschiedene mechanische Systeme folgendermaßen:

1. skleronome und rheonome Systeme
   1. skleronom: explizite Zeitabhängigkeit in den Transformationsgleichungen
   2. rheonom: keine explizite Zeitabhängigkeit
2. holonom und nichtholonome Systeme
   1. holonom: Nebenbedingungen lassen sich in der Form für darstellen
   2. nichtholonom andernfalls
3. konservative und nicht konservative Systeme
   1. konservativ: alle auf ein System wirkenden Kräfte lassen sich aus einer Potentialfunktion herleiten
   2. andernfalls

In den folgenden beiden Abschnitten werden die Lagrangegleichungen zweiter Art auf verschiedene Arten hergeleitet. Die erste Herleitung geht von den Lagrangegleichungen erster Art aus. Die zweite verfolgt einen Ansatz über das Variationsproblem. Bei dieser wird der Beweis zunächst nur für Integrale der Form

$$\displaystyle \int _{a}^{b}F\left( x,y,y'\right) dx$$

angegeben, das Ergebnis ist jedoch leicht auf die in diesem Zusammenhang auftretenden Integrale der Form

$$\displaystyle \int _{a}^{b}F\left( x,y_{1},y_{2},...,y_{f},y_{1}',y_{2}',...,y_{f}'\right) dx$$

zu übertragen.

1.2.2.2.1 Herleitung der Lagrange Gleichungen zweiter Art über die Newton/Eulersche Bewegungsgleichung

Der Ortsvektor eines Massepunktes bezüglich des rechtwinklig kartesischen Sytems ist gegeben durch

$$\mathbf{x}_{.i}=\left[ \begin{array}{c} x^{1}_{.i}\\ x^{2}_... ...{.i} \end{array}\right] =x^{j}_{.i}\mathbf{e}_{j}\; ,\; j=1,2,3$$

Unter der Annahme, daß sich die generalisierten Koordinaten eines Massepunktes sich aufgrund einer geleisteten Arbeit um verändern, folgt für die daraus entstehende Verschiebung

$$d\mathbf{x}_{.i}=\frac{\partial \mathbf{x}_{.i}}{\partial q^{k}}dq^{k}_{.i}=\mathbf{x}_{,k}dq^{k}_{.i}$$

Für die geleistete Arbeit in einem System von Massepunkten ergibt sich nach

$$dW=\sum _{i=1}^{n}\mathbf{F}_{.i}d\mathbf{x}_{.i}=\sum _{i=1... ..._{.i}g_{lk.i}dq^{k}_{.i}=\sum _{i=1}^{n}\Phi _{k.i}dq^{k}_{.i}$$ (1.26)

mit als der mit der ortsabhängigen Koordinate assozierten generalisierten Kraft, welche bestimmt werden kann über

$$\Phi _{k}=\sum _{i=1}^{n}\overline{F}^{l}_{.i}g_{lk.i}=\sum ... ...}\mathbf{F}_{.i}\frac{\partial \mathbf{x}_{.i}}{\partial q^{k}}$$

Für einen Massepunkt folgt aus () und mit dem totalen Differential

$$dW=\frac{\partial W}{\partial q^{k}}dq^{k}=\sum _{i=1}^{n}\P... ...ft( \Phi _{k}-\frac{\partial W}{\partial q^{k}}\right) dq^{k}=0$$

Übertragen auf ein System von Massepunkten ergibt sich damit ein Gleichungssystem in dem, da die unabhängig voneinander sind, die Koeffizienten alle Null sein müssen.

$$\Rightarrow \; \; \Phi _{k}=\frac{\partial W}{\partial q^{k}}$$

Die totale kinetische Energie eines Systems ergibt sich nach

$$T=\frac{1}{2}\sum _{i=1}^{n}m_{i}\dot{\mathbf{x}}_{i}^{2}$$ (1.27)

oder als quadratische Form der generalisierten Geschwindigkeiten für den skleronomen Fall zu

$$T=\frac{1}{2}\sum _{i=1}^{n}m_{i}\frac{\partial \mathbf{x}_{... ...partial \mathbf{x}_{.i}}{\partial q^{l}}\dot{q}^{k}\dot{q}^{l}$$ (1.28)

Im Folgenden werden zwei zur Herleitung der Lagrange-Gleichungen notwendige Rechenregeln abgeleitet:

1. Kürzen von Produkten

   In verallgemeinerten Koordinaten gilt für die Geschwindigkeit aufgrund der Kettenregel
   

   $$\dot{\mathbf{x}}_{.i}=\sum _{k=1}^{f}\left( \frac{\partial \... ...dot{q}^{k}\right) +\frac{\partial \mathbf{x}_{.i}}{\partial t}$$ (1.29)

   
   
   Die Ableitung nach den generalisierten Geschwindigkeiten führt zum gesuchten Zusammenhang
   

   $$\frac{\partial \dot{\mathbf{x}}_{.i}}{\partial \dot{q}^{k}}=... ...ac{\partial \mathbf{x}_{.i}}{\partial q^{k}}=\mathbf{x}_{,k.i}$$ (1.30)

   
   

2. Vertauschen von Operatoren

   Nach () folgt

   $$\begin{eqnarray*} \frac{\partial \dot{\mathbf{x}}_{.i}}{\partial q^{l}} & = & \s... ...) +\frac{\partial ^{2}\mathbf{x}_{.i}}{\partial t\partial q^{l}} \end{eqnarray*}$$
   
   
   
   
   
   Es ergibt sich dann unter der Voraussetzung, daß die zweiten partiellen Ableitungen von stetig sind:
   

   $$\frac{d}{dt}\left( \frac{\partial \mathbf{x}_{.i}}{\partial ... ...\right) =\frac{\partial \dot{\mathbf{x}}_{.i}}{\partial q^{l}}$$ (1.31)

   
   

Der Sinn der Einführung der neuen Koordinaten liegt nun darin, daß man die Nebenbedingungen () für die Lagrangegleichungen der ersten Art auch interpretieren kann als Gleichungen von -dimensionalen Hyperflächen in einem, durch die Punkte des mechanischen Systems definierten -dimensionalen Konfigurationsraum. Der Schnitt der Hyperflächen bildet einen -dimensionalen Riemannschen Teilraum mit den Parametern . Dieser ist der Lösungsraum des Bewegungsproblems, mit anderen Worten, es ist der Raum, der noch von den Bahnkurven der Massepunkte erreichbar ist.

Wenn man die Lagrangegleichungen der ersten Art in den generalisierten Koordinaten ausdrückt, so fallen die Zwangskräfte weg. Dies liegt daran, daß durch Anwenden der Kettenregel bei der Ableitung von

$$g_{.i}\left( \mathbf{x}_{.k}\left( q^{1}\left( t\right) ,...,q^{f}\left( t\right) \right) ,t\right) =0$$

nach den generalisierten Koordinaten zwei zueinander senkrecht stehende vektorielle Funktion entstehen, die miteinander skalar multipliziert werden. Zum einen ergeben sich Vektoren in Richtung der Flächennormalen der durch

$$g_{.i}\left( \mathbf{x}_{.k},t\right) =0$$

definierten Nebenbedingungen, zum anderen ergeben sich die Tangentenvektoren an den, durch die Parameter definierten Lösungsraum der Bewegung.

Die Darstellung von () in generalisierten Koordinaten lautet

$$m_{i}\ddot{\mathbf{x}}_{.i}\frac{\partial \mathbf{x}_{.i}}{\... ...f{x}_{.i}}\cdot \frac{\partial \mathbf{x}_{.i}}{\partial q^{l}}$$

Hierin entsprechen die den senkrecht zu der durch die Nebenbedingung definierten Fläche stehenden Vektoren, während die die Tangentenvektoren an den Lösungsraum sind. Damit gilt

$$\sum _{l=1}^{r}\lambda _{l}\frac{\partial g_{.i}\left( \math... ...x}_{.i}}\cdot \frac{\partial \mathbf{x}_{.i}}{\partial q^{l}}=0$$

und die Darstellung der Bewegungsgleichung in den generalisierten Koordinaten reduziert sich auf die Erweiterung der Newton/Eulerschen Bewegungsgleichung mit . Die Herleitung der Lagrangegleichungen der zweiten Art folgt dann über

$$m_{i}\ddot{\mathbf{x}}_{.i}\cdot \mathbf{x}_{,k.i}=\mathbf{F}_{.i}\cdot \mathbf{x}_{,k.i}$$ (1.32)

Wendet man die Produktregel und () an so gilt allgemein

$$\frac{d}{dt}\left( \dot{\mathbf{x}}_{.i}\frac{\partial \math... ...{.i}\cdot \frac{\partial \dot{\mathbf{x}}_{.i}}{\partial q^{k}}$$

Nach Umstellung und Multiplikation mit folgt mit () unter Annahme der Konstanz der Masse

$$\frac{d}{dt}\left( m_{i}\cdot \dot{\mathbf{x}}_{.i}\frac{\pa... ...f{x}}_{.i}\cdot \frac{\partial \mathbf{x}_{.i}}{\partial q^{k}}$$

Die rechte Seite der obigen Gleichung entspricht dabei der generalisierten Kraft (vgl. ()):

$$m_{i}\cdot \ddot{\mathbf{x}}_{.i}\cdot \frac{\partial \mathb... ...{.i}\frac{\partial \mathbf{x}_{.i}}{\partial q^{k}}=\Phi _{k.i}$$

Nach der Summation über alle Massepunkte folgt:

$$\frac{d}{dt}\left( \sum _{i}\left\{ m_{i}\cdot \dot{\mathbf{... ...mathbf{x}}_{.i}}{\partial q^{k}}\right\} =\sum _{i}\Phi _{k.i}$$ (1.33)

mit () und () ergibt sich andererseits

$$\begin{eqnarray*} \frac{\partial T}{\partial \dot{q}^{k}} & = & \sum _{i}\left\{... ...t \frac{\partial \dot{\mathbf{x}}_{.i}}{\partial q^{k}}\right\} \end{eqnarray*}$$

was eingesetzt in () zu den Lagrangeschen Bewegungsgleichungen zweiter Art führt :

$$\frac{d}{dt}\left( \frac{\partial T}{\partial \dot{q}^{k}}\right) -\frac{\partial T}{\partial q^{k}}=\Phi _{k}$$ (1.34)

Falls sich die generalisierten Kräfte aus einer Potentialfunktion nach

$$\Phi _{k}=\frac{\partial W}{\partial q^{k}}=-\frac{\partial V}{\partial q^{k}}$$ (1.35)

herleiten lassen, wobei das Potential eine reine Ortsfunktion ist, kann man, mit Einführung der Lagrangefunktion

$$L=T-V$$ (1.36)

für () schreiben

$$\frac{d}{dt}\left( \frac{\partial L}{\partial \dot{q}^{k}}\right) -\frac{\partial L}{\partial q^{k}}=0$$ (1.37)

1.2.2.2.2 Herleitung der Euler/Lagrange-Gleichung über das Variationsproblem

Ein dem differentiellen D'Alembertschen Prinzip gleichwertiges, integrales Prinzip ergibt sich aus der Integration der s.g. Lagrangeschen Zentralgleichung^1.11. Es handelt sich um das Hamiltonsche Prinzip^1.12:

Ein konservatives mechanisches System bewegt sich zwischen den Zeitpunkten und so, daß das s.g. Wirkungsintegral

$$\int _{t_{1}}^{t_{2}}Ldt$$

ein Extremum^1.13 annimmt. Anders ausgedrückt gilt

$$\delta ^{1}\int _{t_{1}}^{t_{2}}Ldt=0$$

mit als Zeichen für die erste Variation.

Wie schon zuvor für das D'Alembertsche Prinzip, können auch hier die Lagrangegleichungen zweiter Art abgeleitet werden. Der Einfachheit halber sei dies im folgenden, wie bereits erwähnt, für Funktionen des Typs dargelegt.

Gesucht ist die Bedingung die eingehalten werden muß, damit ein Integral der Art

$$\int _{a}^{b}F\left( x,y,y'\right) dx$$ (1.38)

Extremwerte animmt. Dazu sei zunächst eine Kurve durch im Intervall vorgegeben welche durch die Punkte und verläuft. Eine benachbarte Kurve, ebenfalls durch und kann dann durch

$$y=Y\left( x\right) +\varepsilon \eta \left( x\right) =Y+\varepsilon \eta$$

unter der Nebenbedingung gegeben sein. sei in diesem Fall ein von unabhängiger Parameter. Das zu () äquivalente Integral der benachbarten Kurve lautet damit

$$I\left( \varepsilon \right) =\int _{a}^{b}F\left( x,Y+\varepsilon \eta ,Y'+\varepsilon \eta '\right) dx$$

Die Bedingung für ein Extremum ist dann . Nach der Durchführung der Differentiation unter dem Integral, ergibt sich die notwendige Bedingung zu

$$\frac{\partial I\left( \varepsilon \right) }{\partial \varep... ...artial y}\eta +\frac{\partial F}{\partial y'}\eta '\right) dx=0$$

Partielle Integration führt dann zu

$$\int _{a}^{b}\eta \left( \frac{\partial F}{\partial y}-\frac{d}{dx}\left( \frac{\partial F}{\partial y'}\right) \right) dx=0$$

Da beliebig gewählt werden darf, gilt entsprechend auch

$$\frac{d}{dx}\left( \frac{\partial F}{\partial y'}\right) +\frac{\partial F}{\partial y}=0$$ (1.39)

(s.g Eulersche Gleichung)

1.2.2.2.3 Beispiel zur Lösung von Bewegungsproblemen über die Lagrangegleichungen zweiter Art

Das Prinzip der Lösung von Bewegungsproblemen über die Lagrangegleichungen zweiter Art sei am, schon bei den Lagrangegleichungen erster Art besprochenen, einfachen Beispiel eines, eine schiefe Ebene reibungslos herabgleitenden Körpers, dargestellt. Wie in Abb. ersichtlich, kann im vorliegenden kartesischen Koordinatensystem eine generalisierte Koordinate eingeführt werden, die die Bewegung des Körpers ausreichend beschreibt. In diesem Fall ist leicht erkennbar, daß das System nur einen Freiheitsgrad besitzt und somit eine generalisierte Koordinate zur Beschreibung der Bewegung ausreicht. Die Koordinate unterliegt keinerlei Zwängen, d.h. sie kann je nach Wahl der Anfangsbedingungen jeden Wert annehmen. Die Transformationsgleichungen lauten in diesem Fall

$$\begin{eqnarray*} x & = & s\cdot \cos \left( \alpha \right) \\ y & = & 0\\ z & = & s\cdot \sin \left( \alpha \right) \end{eqnarray*}$$

Nun wird zur Lösung des Problems die Lagrangefunktion nach () aufgestellt. Für die kinetische Energie ergibt sich

$$T=\frac{1}{2}m\left( \dot{x}^{2}+\dot{y}^{2}+\dot{z}^{2}\right) =\frac{1}{2}m\dot{s}^{2}$$

Die potentielle Energie ist gegeben durch

$$V=m\cdot g\cdot z=m\cdot g\cdot s\cdot \sin \left( \alpha \right)$$

Die Lagrangefunktion lautet somit

$$L=\frac{1}{2}m\dot{s}^{2}-m\cdot g\cdot s\cdot \sin \left( \alpha \right)$$

und mit und ergibt sich für ()

$$\frac{d}{dt}\left( \frac{\partial L}{\partial \dot{q}^{k}}\r... ...k}}=m\cdot \ddot{s}-m\cdot g\cdot \sin \left( \alpha \right) =0$$

Daraus ergibt sich die folgende Differentialgleichung zweiter Ordnung

$$\ddot{s}=g\cdot \sin \left( \alpha \right)$$

Die Integration führt schließlich zur Bewegungsgleichung:

$$s\left( t\right) =-\frac{1}{2}gt^{2}\sin \left( \alpha \right) +\dot{s}_{0}t+s_{0}$$

Mit Hilfe des Lagrangeformalismus kann man die Zahl, der für das Bewegungsproblem zu lösenden Differentialgleichungen von auf reduzieren. Jedoch bleibt hier das Problem Konstanten der Bewegung zu bestimmen, die in der Form einer Differentialgleichung erster Ordnung vom Typ die Gesamtordnung des Systems jeweils um verringern. Man kann zwar durch Betrachten von , zyklische Koordinaten bestimmen, d.h. Koordinaten, von der die Lagrangefunktion nicht explizit abhängt^1.14, aber dieses Verfahren ist recht umständlich. Um das Problem zu vereinfachen, führt man zusätzlich weitere Koordinaten so ein, daß diese konstant werden, wenn die ihnen zugehörige generalisierte Koordinate zyklisch wird. Nach () ist dies der Fall bei den generalisierten Impulsen

$$p^{k}=\frac{\partial L}{\partial \dot{q}^{k}}=0$$

Die Erweiterung des Lösungsraumes auf einen -dimensionalen Raum^1.15 mit den Koordinaten für führt zur s.g. Hamilton Theorie mit der eine Möglichkeit gegeben wird, zyklische Koordinaten und damit Konstanten der Bewegung auf einfache Weise zu bestimmen.

1.2.2.3 Lösung mechanischer Probleme über den Hamiltonformalismus

Die Definition der Hamiltonfunktion mit Hilfe der Lagrangefunktion lautet

$$H=\sum _{k=1}^{f}p^{k}\dot{q}^{k}-L$$ (1.40)

Diese Funktion entspricht im Falle eines konservativen Systems, der gesamten Energie

$$H=E_{ges}=E_{kin}+E_{pot}=T+V$$

dieses Systems. Damit ist ersichtlich, daß es sich bei (), unter Voraussetzung daß die potentielle Energie unabhängig von der Zeit ist, nur um eine Form des Energieerhaltungsatzes handelt. Allgemein sollte die Hamiltonfunktion sich mit Hilfe der Lagrangegleichungen als eine Funktion der generalisierten Koordinaten und der generalisierten Impulse darstellen lassen. Man kann nun die Bewegungsgleichungen auf folgende Weise in symmetrischer Form angeben

$$\begin{eqnarray*} \dot{p}^{k} & = & -\frac{\partial H}{\partial q^{k}}\\ \dot{q}^{k} & = & \frac{\partial H}{\partial p^{k}} \end{eqnarray*}$$

Diese Gleichungen werden als kanonische oder Hamiltonsche Bewegungsgleichungen bezeichnet. Die Symmetrie der Impuls- und Ortskoordinaten ist leicht zu erkennen. Jede zyklische Koordinate für die gilt verringert nun die Ordnung der kanonischen Gleichungen um 2. Diese zyklische Koordinate läßt sich durch direkte Integration bestimmen. Ist es nun möglich, durch eine Koordinatentransformation alle zyklisch zu machen, so ist das Bewegungsproblem auf die Integration von Bestimmungsgleichungen der zurückgeführt. Der Vorteil des Hamiltonformalismus besteht nun darin, daß es möglich ist aufgrund der Symmetrie, durch Einführung neuer Koordinaten bzw. über Koordinatentransformationen, zyklische Koordinaten auf einfache Weise zu ermitteln. Jedoch ist die kanonische Darstellung nicht invariant bezüglich beliebiger Transformationen ,wie es bei den Lagrangeschen Gleichungen zweiter Art der Fall ist. Eine Transformation, für die die kanonische Darstellung invariant ist, wird als kanonische Transformation bezeichnet. Seien die neuen Koordinaten durch

$$\begin{eqnarray*} \overline{p}^{k} & = & \overline{p}^{k}\left( p^{k},q^{k},t\ri... ...erline{q}^{k} & = & \overline{q}^{k}\left( p^{k},q^{k},t\right) \end{eqnarray*}$$

gegeben, und gibt es eine neue Funktion , so daß

$$\begin{eqnarray*} \dot{\bar{p}}^{k} & = & \frac{\partial \overline{H}}{\partial ... ...k} & = & \frac{\partial \overline{H}}{\partial \overline{p}^{k}} \end{eqnarray*}$$

gilt, dann existiert eine Lagrangefunktion , für die dann folgt

$$\displaystyle \overline{H}=\sum \overline{p}^{k}\overline{\dot{q}}^{k}-\overline{L}$$

Aus dem Hamiltonschen Prinzip

$$\begin{eqnarray*} \displaystyle \int _{t_{1}}^{t_{2}}Ldt & \equiv & min\\ \int _{t_{1}}^{t_{2}}\overline{L}dt & \equiv & min \end{eqnarray*}$$

läßt sich eine notwendige Existenz einer s.g. Erzeugenden Funktion mit ableiten. Hängt diese von , und ab, so wird sie mit bezeichnet. Hier gilt dann

$$\begin{eqnarray*} \overline{H}=\frac{\partial \widetilde{S}}{\partial t}+H, & \o... ...{p}}^{k}=\frac{\partial \overline{H}}{\partial \overline{q}^{k}} \end{eqnarray*}$$

Wenn man nun erreicht, daß und zyklische Variablen werden, so läßt sich, wie zuvor bereits erwähnt, das mechanische Problem einfach lösen. Es existiert in Form der Hamilton-Jacobischen Gleichung eine Bedingung für die s.g. Erzeugende Funktion einer solchen kanonischen Transformation. Um das Bewegungsproblem zu lösen, wird eine Funktion derart bestimmt, daß sie die Bedingung oder

$$\frac{\partial \widetilde{S}}{\partial t}+H\left( \frac{\partial \widetilde{S}}{\partial t},q^{k},t\right) =0$$ (1.41)

erfüllt. Die Lösung kann in der Form erhalten werden, dabei bezeichnen die die übrigen freien Parameter^1.16, welche bis auf eine gemeinsame additive Konstante festgelegt sind. Die Ortskoordinaten können über und bestimmt werden. Falls die Hamiltonfunktion nicht explizit von der Zeit abhängt, läßt sich der von der Zeit unabhängige Teil der Erzeugenden, nach Separation der Variablen, schreiben als

$$\widehat{S}=\widehat{S}_{1}\left( q^{1}\right) +\widehat{S}_{2}\left( q^{2}\right) +...+\widehat{S}_{f}\left( q^{f}\right)$$

Damit vereinfacht sich die Hamilton-Jacobische Gleichung zu

$$H\left( \frac{\partial \widehat{S}}{\partial q^{k}},q^{k}\right) =E=const.$$

worin die Gesamtenergie des Systems bezeichnet.

1.3 Herausstellen der Analogie zwischen dem geometrischen und mechanischen Problem

1.3.1 Ableitung der Differentialgleichung der geodätischen Linie über die Eulersche Gleichung

Die Differentialgleichung der geodätischen Linie kann nun, nach den in den vorigen Abschnitten beschriebenen Verfahren, abgeleitet werden. Die Eigenschaft der geodätischen Linie, die kürzeste Verbindung in beliebigen, durch ihre Metrik ausgezeichneten Riemannschen Räumen darzustellen, führt zu einer Möglichkeit der Ableitung ihrer Bestimmungsgleichung über die Minimierung des durch

$$s_{AB}=\int _{A}^{B}ds=\int _{A}^{B}\sqrt{g_{ij}du^{i}du^{j}... ...s_{1}}^{s_{2}}\sqrt{g_{ij}\frac{du^{i}}{ds}\frac{du^{j}}{ds}}ds$$

gegebenen Abstandsmaßes. Eingesetzt^1.17 in ergibt sich, wenn zunächst angenommen wird und somit die Wurzel vernachlässigt werden kann,

$$\frac{d}{ds}\left( \frac{\partial \left( g_{ij}\frac{du^{i}}... ...}\frac{du^{i}}{ds}\frac{du^{j}}{ds}\right) }{\partial u^{k}}=0$$ (1.42)

$$\Rightarrow \frac{d}{ds}\left( 2\cdot g_{ij}\frac{du^{i}}{ds... ...u^{k}}\frac{du^{i}}{ds}\frac{du^{j}}{ds}=0\left\vert /2\right.$$

Nach Umordnung folgt dann

$$g_{ij}\ddot{u}^{i}+\frac{1}{2}\left( g_{ik,j}+g_{kj,i}-g_{ij,k}\right) \dot{u}^{i}\dot{u}^{j}=0$$

Jetzt werden wieder die Christoffelsymbole der ersten Art eingeführt und analog zu den in den vorigen Kapiteln beschriebenen Verfahren ergibt sich die folgende Differentialgleichung als Bestimmungsgleichung^1.18 für den Verlauf der geodätischen Linie im durch definierten Riemannschen Raum.

$$\frac{d^{2}u^{l}}{ds^{2}}+\Gamma _{jk}^{l}\frac{du^{j}}{ds}\frac{du^{k}}{ds}=0$$

1.3.2 Herleitung der expliziten Form der Lagrange-Gleichungen zweiter Art

In den vorhergegangenen Abschnitten wurde zum einen die Herleitung der Lagrangeschen Gleichungen aus der Erweiterung der Newton/Eulerschen Bewegungsgleichung in generalisierten Koordinaten aufgezeigt, zum anderen wurde gezeigt, daß die Bedingungsgleichungen des Variationsproblems

$$\displaystyle \delta ^{1}\int ^{b}_{a}F\left( t,q_{1},q_{2},...,q_{f},\dot{q}_{1},\dot{q}_{2},...,\dot{q}_{f}\right) dt=0$$

den Differentialgleichungen der geodätischen Linie ebenso äquivalent, sind wie auch den Lagrangegleichungen der zweiten Art. Es handelt sich also bei der Lösung der geodätischen Hauptaufgaben, wie auch bei der Lösung der Bewegungsgleichung über die Lagrangegleichungen der zweiten Art, um eine Lösung desselben Problems. Um die Analogie etwas besser herausstellen zu können, sei an dieser Stelle die Herleitung der expliziten Form der Lagrangegleichungen über die Einführung der generalisierten Koordinaten angeführt. Dieses schon zuvor beschriebene Verfahren, das Bewegungsproblem über eine Einschränkung des Lösungsraumes, durch die Wahl von generalisierten Koordinaten, zu vereinfachen, formuliert das mechanische Problem in einer -dimensionalen differenzierbaren Mannigfaltigkeit. Über die Einführung einer Metrik wird diese Mannigfaltigkeit zu einem Riemannschen Raum. Die Bedingung welche die Metrik erfüllen muß, ist an die Existenz eines symmetrischen Tensorfeldes gebunden, so daß sich die Länge einer Kurve zwischen 2 Parameterwerten und berechnet läßt nach

$$s=\int _{t_{0}}^{t_{1}}\sqrt{g_{ij}\dot{q}^{i}\dot{q}^{j}}dt$$

Wie in () gezeigt, erfüllt die kinetische Energie diese Forderung, so daß man einen Metriktensor der Form

$$g_{ij}\left( q^{k}\right) =\sum _{l}m_{l}\frac{\partial \mat... ...\partial q^{i}}\frac{\partial \mathbf{x}_{.l}}{\partial q^{j}}$$ (1.43)

einführen kann. Hieraus ergibt sich die Transformation wobei das auftretende Linienelement als kinematisches Linienelement bezeichnet wird. Das Hamiltonsche Prinzip besagt, daß falls ein konservatives System vorliegt und sich daher die Massenkräfte nach () berechnen lassen, die erste Variation von () verschwindet. Die Lagrangefunktion ergibt sich dann mit () und zu

$$L=\frac{1}{2}g_{ij}\left( q^{k}\right) \dot{q}^{i}\dot{q}^{j}-V\left( q^{k}\right)$$ (1.44)

Setzt man diese in die Eulerschen Gleichungen ein, so folgt analog zu ()

$$\frac{d}{dt}\left( g_{ij}\dot{q}^{i}\right) -\frac{1}{2}\fra... ...rtial g_{ij}}{\partial q^{k}}\dot{q}^{i}\dot{q}^{j}-\Phi _{k}=0$$

Wie leicht ersichtlich, ergibt sich dann

$$\frac{d^{2}q^{l}}{dt^{2}}+\Gamma _{jk}^{l}\frac{dq^{j}}{dt}\frac{dq^{k}}{dt}=\Phi _{k}$$

Diese Gleichung entspricht (), damit gilt für den kräftefreien Fall, d.h. für den Fall mit Potential die Differentialgleichung der geodätischen Linie:

$$\frac{d^{2}q^{l}}{dt^{2}}+\Gamma _{jk}^{l}\frac{dq^{j}}{dt}\frac{dq^{k}}{dt}=0$$

Somit ist die Lösung der Differentialgleichung der geodätischen Linie vollkommen analog zur Lösung der Bewegungsgleichung eines Massepunktes auf einer Fläche mit dem Potential Null und der Bogenlänge als unabhängigem Parameter. Dies führt dazu, daß Lösungsstrategien zwischen den beiden Problemen übertragbar sind. Vor allem läßt sich damit die Differentialgleichung der geodätischen Linie in kanonischer Form angeben.

1.4 Eine Möglichkeit der Lösung der geodätischen Hauptaufgaben über die kanonischen Gleichungen der geodätischen Linie

Aufgrund der zuvor beschriebenen Analogie, ist es möglich, den Hamiltonformalismus zur Lösung der geodätischen Hauptaufgaben zu verwenden. Angewandt auf die Problematik der geodätischen Hauptaufgaben auf einem Rotationsellipsoid, folgt aus () und wenn das Potential den Wert erhält, daß die Gesamtenergie, die in einem solchen Fall durch die Hamiltonfunktion gegeben ist, nur aus dem Anteil der kinetischen Energie besteht. Wichtig ist hierbei die Betrachtung des kinematischen Linienelementes , welches beim Lagrangeformalismus das geometrische Linienelement ersetzt. Falls nun der Hamiltonformalismus angewandt werden soll, so muß die quadratische Form für die kinetische Energie durch die Metrik des Rotationsellipsoides ersetzt werden. Demzufolge ergibt sich analog zu () für die Hamiltonfunktion

$$H=T+V=T=\frac{1}{2}g_{ij}\dot{q}^{i}\dot{q}^{j}$$ (1.45)

Im Vergleich zum allgemeinen Linienelement sieht man, daß in der kinetischen Energie und damit auch im kinematischen Linienelement die Geschwindigkeiten auftreten. Bei der Allgemeinen Darstellung sind jedoch nur die Differentiale vorhanden. Aus diesem Grund ergibt sich die Hamiltonfunktion zur Bestimmung der kürzesten Verbindung zweier Punkte in beliebigen Riemannschen Räumen zu

$$H=\frac{1}{2}\widetilde{g}_{ij}\frac{du^{i}}{ds}\frac{du^{j}}{ds}$$

Um das Problem lösen zu können, müssen die generalisierten Geschwindigkeiten in () eliminiert werden, damit sich die Hamiltonfunktion in der benötigten Form ergibt. Dabei kann man die Bestimmungsgleichung für die generalisierten Impulse

$$p^{k}=\frac{\partial T}{\partial \dot{q}^{k}}$$ (1.46)

verwenden.

Für den Fall der Lösung der geodätischen Hauptaufgaben folgt mit () und

$$\left[ \widetilde{g}_{ij}\left( u^{\alpha }\right) \right] =... ... 0 & N^{2}\cos ^{2}\left( \varphi \right) \end{array}\right]$$

die Hamiltonfunktion

$$H=\frac{1}{2}\left( M^{2}\cdot \left( \frac{d\varphi }{ds}\r... ...i \right) \cdot \left( \frac{d\lambda }{ds}\right) ^{2}\right)$$

oder für ein allgemeines orthogonales Flächenparametersystem

$$H=\frac{1}{2}\left( E\cdot \left( \frac{du^{1}}{ds}\right) ^{2}+G\cdot \left( \frac{du^{2}}{ds}\right) ^{2}\right)$$ (1.47)

Unter Verwendung von () für den Fall eines Potentials , d.h. es gilt , ergeben sich die Bestimmungsgleichungen für die 'generalisierten Impulse'

$$\begin{eqnarray*} p_{u^{1}} & = & E\frac{du^{1}}{ds}\\ p_{u^{2}} & = & G\frac{du^{2}}{ds} \end{eqnarray*}$$

Einsetzen der quadrierten Impulse in () führt zur Hamiltonfunktion des Systems, die dann im Falle der geodätischen Hauptaufgaben auf einem Rotationsellipsoid folgende Form annimmt

$$H=\frac{1}{2}\left( \frac{p_{\varphi }^{2}}{M^{2}}+\frac{p_{... ...right) \right) \cos ^{2}\left( \varphi \right) }{a^{4}}\right)$$

Die kanonischen Gleichungen lauten dann

$$\frac{d\varphi }{ds}=\frac{\partial H}{\partial p_{\varphi }}$$

$$\frac{d\lambda }{ds}=\frac{\partial H}{\partial p_{\lambda }}$$

sowie

$$\frac{dp_{\varphi }}{ds}=-\frac{\partial H}{\partial \varphi }$$

$$\frac{dp_{\lambda }}{ds}=-\frac{\partial H}{\partial \lambda }$$

Auf einfache Weise lösen läßt sich das Problem, wenn man die Erzeugende einer Transformation findet, für die die neuen Koordinaten zyklisch werden. Diese muß dann () erfüllen, wobei die Separation der Variablen zur gesuchten Lösung für die Erzeugende führt (es ergeben sich elliptische Funktionen). Aus den konstanten neuen 'Impulsen' und Koordinaten lassen sich über die Transfomationsformeln sowie und in Abhängigkeit von berechnen, so daß damit das Problem gelöst ist.

An dieser Stelle sei noch auf zwei weitere Ansätze zur Lösung der geodätischen Hauptaufgaben über die kanonische Darstellung der geodätischen Linie verwiesen. Der Ansatz von SCHÜRER (1975) geht im Unterschied zum zuvor dargelegten, von der reduzierten Breite mit aus, was an () nichts ändert. Der Metriktensor hat dann die etwas einfachere Form

$$ds^{2}=a^{2}\left( 1-e^{2}\cos ^{2}\left( \beta \right) \rig... ...eta ^{2}+a^{2}\cos ^{2}\left( \beta \right) \cdot d\lambda ^{2}$$

und werden entsprechend dem obigen Ansatz eingeführt, jedoch umgeht er die Lösung der elliptischen Funktionen, indem er die Hamiltonfunktion in zwei Anteile zerlegt. Die Unterteilung erfolgt derart, daß die erste Teilfunktion ohne Probleme mit dem Hamiltonformalismus gelöst und die Restfunktion als Störfunktion in Analogie zur Störungsrechnung in der Himmelsmechanik betrachtet werden kann. Die Lösung von () mit führt zur 'Bewegungsgleichung' der geodätischen Linie auf einer Kugelfläche. Die hier abgeleiteten Formeln eignen sich gut zur Untersuchung des Verlaufs der geodätischen Linien über mehrere Umläufe hinweg.

Einen weitere Lösung kann KLOTZ (1991) entnommen werden. Er gelangt über den Hamiltonformalismus zu einer integralen Darstellung der ellipsoidischen Länge in Abhängigkeit der von reduzierten Breite .

2. Numerische Untersuchungen zu den geodätischen Hauptaufgaben

Im Rahmen dieses zweiten Teils der Arbeit wurde mit MATLAB unter WINDOWS ein Programmsystem erstellt, das dem Benutzer die Möglichkeit bietet numerische und geometrische Eigenschaften von geodätischen Linien auf einem Rotationsellipsoid zu untersuchen. Hierbei wurde für die Berechnung der geodätischen Hauptaufgaben auf die numerische Integration eines Differentialgleichungssystems der geodätischen Linie in kartesischen Koordinaten zurückgegriffen. Um die geometrischen Eigenschaften von geodätischen Linien zu untersuchen, ist es notwendig ihren Verlauf auf der Bezugsfläche graphisch zu visualisieren. Das erstellte Programmsystem bietet hierzu die Möglichkeit, nach der Berechnung der geodätischen Hauptaufgaben, den Verlauf der zugehörigen Geodäte auf der Bezugsfläche, stellvertretend durch eine frei zu wählende Zahl von Stützstellen, dreidimensional darzustellen. Der Blickpunkt kann beliebig gewählt werden und es besteht die Möglichkeit einer Vergrösserung bzw. Verkleinerung. Für die Aufgaben der Landesvermessung ist es nötig, die geodätischen Hauptaufgaben unter Vorgabe einer bestimmten Fehlerschranke zu lösen. In den meisten numerischen Ansätzen geschieht dies durch eine Fehlerschätzung. In dieser Arbeit wird jedoch ein anderer Weg beschritten. Es wird eine, speziell auf die Lösung der geodätischen Linien auf einem Rotationsellipsoid bezogene Funktion, mit der die Zahl von Iterationen des Differentialgleichungslösers nach Vorgabe der geforderten Genauigkeit berechnet werden kann, empirisch abgeleitet.

Als ein Beispiel für die Untersuchung numerischer Eigenschaften der Geodäten auf einem Rotationsellipsoid bietet das Programmsystem weiterhin die Möglichkeit ein geodätisches Parallelkoordinatensystem in Soldnerscher Anordnung zu berechnen und darzustellen.

2.1 Numerische Verfahren zur näherungsweisen Bestimmung der Lösung von Differentialgleichungen erster Ordnung als Anfangswertproblem

2.1.1 Beschreibung des Problems

Allgemein ist eine Differentialgleichung erster Ordnung gegeben durch , die Werte sowie seien als bekannt vorausgesetzt. Das Problem besteht nun darin, den Funktionswert der Funktion an der Stelle zu bestimmen. Im weiteren sei dieser mit bezeichnet. Die zugehörige Integralgleichung lautet damit

$$\displaystyle \overline{y}=y_{0}+\int _{x_{0}}^{x_{0}+h}f\left( x,y\right) dx$$

.

2.1.2 Lösung über Reihenentwicklung

Eine Möglichkeit der Lösung besteht darin, die Funktion an der Stelle in eine Taylorreihe

$$y=g\left( x_{0}\right) +\left( x-x_{0}\right) g'\left( x_{0}... ...rac{1}{6}\left( x-x_{0}\right) ^{2}g'''\left( x_{0}\right) +...$$

zu entwickeln. Mit der Ausgangsgleichung gilt und für die höheren Ableitungen folgt

$$y''\left( x\right) =\frac{\partial f}{\partial x}+f\frac{\partial f}{\partial y}$$

$$y'''\left( x\right) =\frac{d}{dx}\left( \frac{\partial f}{\p... ...rtial y}\right) ^{2}+f^{2}\frac{\partial ^{2}f}{\partial y^{2}}$$

$$...$$

Wenn man nun die partiellen Ableitungen durch , , , , ersetzt, so erhält man den folgenden übersichtlichen Ausdruck für die Taylorentwicklung an der Stelle :

$$\overline{y}=y_{0}+h\cdot f_{0}+\frac{1}{2}\left( p_{0}+f_{0... ...dot s_{0}+f_{0}\cdot q_{0}^{2}+f_{0}^{2}\cdot t_{0}\right) +...$$

Für Approximationszwecke ist es ausreichend die ersten beiden Glieder der Taylorentwicklung zu betrachten. Wenn man das Integrationsintervall h in infinitisimal kleine Teilintervalle der Breite , unterteilt, so kann man auf diese Weise den Funktionswert der Funktion an der Stelle iterativ bestimmen. Der Rechenaufwand kann durch Berücksichtigung höherer Glieder der Taylorentwicklung zwar verringert werden, jedoch erhöht sich dabei die Komplexität der Funktionen.

2.1.3 Runge-Kutta-Verfahren

Eine andere Alternative ergibt sich, wenn angenommen wird, daß bspw. die Werte , , der Funktion an den Stellen , und vorliegen. Mit diesen Werten folgt mit der Prismaregel^2.1 :

$$k=\overline{y}-y_{0}\int ^{x_{0}+h}_{x_{0}}f\left( x,y\right... ...rac{1}{2}h,y_{1}\right) +f\left( x_{0}+h,y_{2}\right) \right)$$ (2.1)

Die in Wirklichkeit, bis auf , unbekannten Funktionswerte von , werden mittels einer nach Runge bezeichneten Methode wie folgt approximiert:

$$y_{1}\approx \widetilde{y_{1}}=y_{0}+\frac{1}{2}hf_{0}$$

$$y_{2}\approx \widetilde{y_{2}}=y_{0}+hf\left( x_{0}+h,y_{0}+hf_{0}\right)$$

Eingesetzt in (), die genäherte Gleichung für , ergibt sich dann

$$\displaystyle k\approx \frac{h}{6}\left( f_{0}+4f\left( x_{0... ...}+h,y_{0}+hf\left( x_{0}+h,y_{0}+hf_{0}\right) \right) \right)$$

Hierbei wird iterativ über

$$k\approx \frac{1}{6}\left( k_{1}+4k_{4}+k_{3}\right)$$

mit

, , und

bestimmt. Je nach Zahl und Wahl der Approximationen für , kann das Verfahren abgeändert werden. Die Abbildung gibt das Prinzip des Verfahrens nach Runge wieder.

Abbildung: Prinzip des Verfahrens nach Runge

2.1.4 Verallgemeinerung des Runge-Kutta-Verfahrens auf höhere Ordnungen und Systeme von Differentialgleichungen

Wie leicht zu erkennen ist, läßt sich dieses Verfahren allgemeiner darstellen mit , wobei sich als eine Funktion der Integrationsweite und einer diskreten Näherungsfunktion, die im weiteren mit bezeichnet wird, darstellen läßt . Für ^2.2 gilt:

$$\displaystyle \Phi \left( x,y_{0},h,f\right) =\sum ^{m}_{j=1}c_{j}\cdot k_{j}\left( x,y_{0},h\right)$$

Es werden dabei von und unabhängige Zahlen , , , mit und unter der Voraussetzung der Bedingungen und so bestimmt, daß die Beziehungen

$$k_{j}\left( x,y_{0}\right) =f\left( x+a_{j}h,y_{0}+h\cdot \left( \sum ^{j-1}_{l=1}b_{ij}k_{j}\right) \right)$$

für und

$$\displaystyle \overline{y}=y_{0}+h\cdot \sum ^{m}_{j=1}c_{j}k_{j}$$

eine möglichst hohe Konsistenzordnung ergeben. Da diese Näherung über die Approximation der Funktionswerte nicht zu einer exakten Lösung führen kann, wird das Integrationsintervall, wie bereits bei der Integration über Reihenentwicklung erwähnt, in Teilintervalle eingeteilt. Auf diesem Weg erhält man das Ergebnis über

$$y_{i}=y_{i-1}+h_{i-1}\cdot \sum ^{m}_{j=1}c_{j}k_{j}\left( y_{i-1}\right)$$

für .

Das Verfahren wird auf Systeme von Differentialgleichungen übertragen, indem man für die Koeffizienten eine Matrix ansetzt, in der für jede Gleichung des Systems eine Spalte sukzessive in den Iterationsschritten bestimmt wird.

2.2 Realisierung eines Programms zur Berechnung algebraisch vorliegender Differentialgleichungssysteme nach dem Runge-Kutta-Verfahren 7. Ordnung

Für die Realisierung unter MATLAB 5.0 wird, wegen des benötigten Zugriffs auf algebraische Objekte die symbolische Toolbox von MATLAB verwandt. Diese beruht auf einem MAPLE Kern und ermöglicht es, das Programm als Funktion zu realisieren, welche symbolische Vektoren als Eingabewerte verarbeiten kann. Die Berechnung erfolgt nach dem im vorigen Abschnitt beschriebenen Verfahren, wobei die Zahl der Iterationen vorgegeben werden muß, da keine automatische Schrittweitensteuerung implementiert wurde. Der Aufruf der Funktion (Abb. zeigt ein Flußdiagramm) von der MATLAB Eingabeaufforderung geschieht über:

z=RK7multi(DGLSYS,i,h,Naeherung,x0);

Abbildung: Flußdiagramm der Funktion RK7multi.m

Dabei bezeichnet DGLSYS einen symbolischen Spaltenvektor der Form sym('[y1;y2;...;yn]'), der die rechten Seiten der Differentialgleichungen enthält. Die Variable i steht für die Zahl der Iterationen, h bezeichnet die Integrationsweite, Y ist der symbolische Vektor der Variablennamen in der Form sym('[Y1;Y2;...;Yn]'), Naeherung ist ein Zeilenvektor mit den bekannten Werten der Variablen an den Näherungstellen und x0 bezeichnet den Startwert der unabhängigen Variabeln .

Beispiel: Es sei folgendes Differentialgleichungssytem erster Ordnung gegeben:

$$\begin{array}{c} \frac{\partial y}{\partial x}=x+\sqrt{z}\\ \frac{\partial z}{\partial x}=y-\sqrt{z} \end{array}$$

Für dieses sind und an der Stelle , für diejenige partikuläre Lösung zu prädizieren, für die gilt

$$y=0.5\mid _{x=0.2}$$

und

$$z=0\mid _{x=0.2}$$

Der Funktionsaufruf lautet dann:  
z=RK7multi(sym('[x+sqrt(z);y-sqrt(z)]'),1,0.1,sym('[y;z]'),[0.5,0],0.2);
(für i kann hier 1 angesetzt werden) Es ergibt sich

$$\begin{eqnarray*} y & = & 0.5384\\ z & = & 0.0376 \end{eqnarray*}$$

2.3 Lösung der 1. geodätischen Hauptaufgabe mit Hilfe der Funktion RK7multi.m

2.3.1 Beschreibung der theoretischen Grundlagen

ARINGER (1994) beschreibt die Möglichkeit der Lösung der geodätischen Hauptaufgaben über eine Darstellung der Differentialgleichungen der geodätischen Linie in kartesischen Koordinaten. Der Grund für deren Einführung liegt in den Problemen, die sich aufgrund von Singularitäten bei der Parameterdarstellung regulärer Flächen ergeben. Werden beispielsweise geographische Breite und geographische Länge als Parameter zur Darstellung einer Einheitskugel verwendet, so ergeben sich diese in den Polen. Eine andere Möglichkeit bietet die Wahl von Mongeschen Flächenstücken, wobei hier Fallunterscheidungen zu treffen sind, da sonst Singularitäten an den begrenzenden Großkreisen auftreten können. Diese Problematiken lassen sich ebenfalls auf das Ellipsoid übertragen. Dort kommt es dann in der Nähe der Unstetigkeitsstellen zu numerischen Unsicherheiten bei der Berechnung der geodätischen Hauptaufgaben über Quadraturverfahren. Als Lösung schlägt ARINGER die Wahl der impliziten Form

$$Q\left( x_{1},x_{2},x_{3}\right) =0$$

der Flächendarstellung in kartesischen Koordinaten vor. Über die Bestimmungsgleichung der geodätischen Linie

$$\frac{\partial ^{2}\mathbf{x}}{\partial s^{2}}=\kappa \cdot \mathbf{n}=\kappa _{n}\cdot \mathbf{n}$$

mit der Normalkrümmung

$$k_{n}=\frac{-1}{\sqrt{\frac{x_{1}^{2}}{a^{4}}+\frac{x_{2}^{2... ...2}}\left( \frac{\partial x_{3}}{\partial s}\right) ^{2}\right)$$

und dem zugehörigen Normalenvektor , wobei , leitet er das folgende Differentialgleichungssystem her:

$$\left[ \begin{array}{c} DY_{1}=Y_{2}\\ DY_{2}=\frac{-Y_{1}}... ...}^{2}+\frac{a^{2}}{b^{2}}Y_{6}^{2}\right) \end{array}\right]$$ (2.2)

Die Differentiale sind mit , die Variablen mit bezeichnet. stehen hierbei für die kartesischen Koordinaten des Anfangspunktes, die drei übrigen Variablen bezeichnen die Komponenten des Tangenteneinheitsvektors im Anfangspunkt, die zunächst über die gegebenen Größen zu berechnen sind. Wie bereits erwähnt, besitzt die implizite Art der Flächendarstellung und damit das aus dieser abgeleitete Differentialgleichungssystem, den Vorteil, daß es frei von Singularitäten ist. Den meisten Verfahren, die bei der numerischen Integration der ersten Hauptaufgabe Verwendung finden, liegt das, in der Literatur meist angegebene, Gleichungssystem

$$\left[ \begin{array}{c} \frac{\partial u}{\partial s}=\frac{... ...}\sqrt{g_{11}}}\sin \left( \alpha \right) \end{array}\right]$$ (2.3)

zugrunde. Für das Linienelement in ellipsoidischen Koordinaten gilt

$$ds^{2}=M^{2}\cdot d\varphi ^{2}+N^{2}\cos ^{2}\left( \varphi \right) \cdot d\lambda ^{2}$$ (2.4)

mit dem Meridiankrümmungshalbmesser und dem Querkrümmungshalbmesser . bezeichnet den Polkrümmungshalbmesser und es gilt , mit der zweiten numerischen Exzentrizität . Für den kovarianten Metriktensor folgt dementsprechend

$$\left[ g_{ij}\right] =\left[ \begin{array}{cc} M^{2} & 0\\ 0 & N^{2}\cos ^{2}\left( \varphi \right) \end{array}\right]$$

Eingesetzt in () ergibt sich folgendes Differentialgleichungssystem erster Ordnung für die geodätische Linie auf dem Rotationsellipsoid:

$$\left[ \begin{array}{c} \frac{\partial \varphi }{\partial s}... ...varphi \right) \sin \left( \alpha \right) \end{array}\right]$$ (2.5)

Dieses Gleichungssystem kann zur Lösung der 1. geodätischen Hauptaufgabe direkt integriert werden, wobei hier jedoch auf die schon erwähnten Singularitäten in den Polen, hervorgerufen durch die Wahl der ellipsoidischen Koordinaten, zu achten ist.

2.3.2 Aufruf des Lösers über die graphische Benutzeroberfläche guiawp_graphik.m

Mit MATLAB wurde eine graphische Benutzerschnittstelle entwickelt, die die Berechnung der 1. geodätischen Hauptaufgabe in kartesischen Koordinaten auf einem vorgegebenen Rotationsellipsoid ermöglicht. Dabei wird zur Lösung des Gleichungssystems () über eine Anfangswertaufgabe auf die zuvor eingeführte Funktion RK7multi zurückgegriffen. Der Funktionsaufruf lautet in diesem Fall:

z=RK7multi(DGLSYS,schritte,s,Y,[x1 t(1) x2 t(2) x3 t(3)],0);

mit den Variablen

• DGLSYS: Symbolische Zeilenmatrix mit den rechten Seiten der Differentialgleichungen (). (kann dem Quelltext entnommen werden)
• Y=sym('[y1;y2;y3;y4;y5;y6]');
• schritte bezeichnet die vorzugebende Zahl der Iterationen
• s entspricht dem Integrationsintervall
• [x1 t(1) x2 t(2) x3 t(3)] gibt die Komponenten des Anfangsvektors und des Tangentenvektors im Anfangspunkt an

Es sei hier zu bemerken, daß in DGLSYS die Konstantenbezeichner vor Aufruf der Funktion noch durch ihre numerischen Werte zu ersetzen sind. Die in Abb. dargestellte Benutzerschnittstelle (G(raphical)U(ser)I(nterface)) kann durch Aufruf der Funktion guiawp_graphik.m an der MATLAB Eingabeaufforderung gestartet werden.

Abbildung: GUI für Lösung des Anfangswertproblems

Die Felder für die Anfangswerte und die vorzugebenden Ellipsoidparameter lassen sich editieren, und die Funktion geo_haupt_1.m, welche die Übergabe der Daten zwischen dem GUI und dem Script RKF7multi.m übernimmt, wird über den Button 'Berechnung Starten' aufgerufen. Das Programm selbst bietet zum einen die Möglichkeit eine konstante Zahl von Iterationen vorzugeben, zum anderen kann aus verschiedenen Genauigkeiten ausgewählt werden. Mit der Ableitung, des für diese Aufgabe empirisch zu bestimmenden Zusammenhangs zwischen Genauigkeit und vorzugebender Zahl der Iterationen, beschäftigt sich der folgende Abschnitt. Hier sei zu erwähnen, daß man durch Wahl des Schalters 'Genauigkeit auswählen' unter 9 verschiedenen Genauigkeiten der Lösung der Anfangswertaufgabe wählen kann, und mit Übernahme dieser Wahl, die Zahl der benötigten Iterationen automatisch, nach Auswerten der Iterationsfunktion, in die Oberfläche übernommen werden. Auf die Bedeutung der rechts im Fenster angeordneten Schalter 'Ephemeridenberechnung' und 'Zurück zur Hauptauswahl' wird in den folgenden Abschnitte noch eingegangen. Die Rechenzeit^2.3 des Programms (0.1''-2''/AWP) liegt in der Regel beim 2- bis 9- fachen eines kompilierten FORTRAN-Programms.

Abbildung: Rechenzeit in Abhängigkeit von der Entfernung

Verkürzen läßt sich diese, indem man entweder vorhandene FORTRAN- bzw. C-Funktionen über ein API (Application Program Interface) in die MATLAB-Routinen einbindet oder einen MATLAB-Compiler wie MATCOM oder MIDEVA^2.4 benutzt, der die m-Scripte in C++ übersetzt und dann compiliert.

2.3.3 Ableitung einer empirischen Funktion zur Bestimmung der Zahl der nötigen Iterationen

Da bei der Realisierung des Programms auf eine automatische Schrittweitensteuerung^2.5 verzichtet wurde, ist es wichtig, die für eine vorgegebene Genauigkeit der Endpunktkoordinaten nötige Schrittweite, bzw. die Zahl der Iterationen in Abhängigkeit der Integrationsweite des Runge-Kutta-Verfahrens, empirisch abzuleiten. Über ein Hilfsprogramms wird dabei die Zahl der Iterationen solange erhöht, bis sich die Endpunktkoordinaten der Lösung mit Iterationen von denen der Lösung mit Iterationen nur noch um eine vorgegebene Schranke unterscheiden. Setzt man beispielsweise als Schranke für die Genauigkeit der kartesischen Koordinaten der Endpunkte einen Wert von an, so ergibt der in Abb. dargestellte, wie zu erwarten lineare, Zusammenhang zwischen der Anzahl der benötigten Iterationen und der Länge der geodätischen Linie.

Abbildung: Iterationen in Abhängigkeit von der Integrationsweite s

Im weiteren bestehen Abhängigkeiten der nötigen Iterationen von der Abplattung, dem Azimut sowie von der geforderten Genauigkeit . In den oberen beiden Graphen der Abb.^2.6 ist gut zu erkennen, daß der Einfluß der Abplattung auf die Iterationen bei Integration in Polrichtung am stärksten ist, wohingegen bei Integration in Äquatorrichtung diese Abhängigkeit nicht mehr auftritt. Dies war von vornherein klar, da mit geringerem Azimut die Abweichung der geodätischen Linie von einem ebenen Kreisabschnitt immer größer wird. In Äquatorrichtung entspricht die geodätische Linie dann dem Äquator. Aufgrund der Tatsache, daß die maximale Zahl der benötigten Iterationen in Polrichtung auftritt, wird bei der weiteren Ableitung einer Näherungsfunktion immer von dieser Integrationsrichtung ausgegangen und somit der Einfluß des Azimuts unter Inkaufnahme längerer Rechenzeiten unberücksichtigt gelassen. Die beiden unteren Graphen zeigen zum einen den Einfluß der Abplattung auf die Funktion , die im weiteren die Abhängigkeit der nötigen Iterationen von der Integrationsweite beschreibt, sowie den Einfluß einer variierenden Integrationsweite auf die Funktion , durch die der Zusammenhang zwischen den Iterationen und dem Azimut erfaßt wird.

Abbildung: Abhängigkeiten von Azimut und Abplattung

Der nächste Schritt besteht darin die Einflüsse der übrigen 3 Parameter, Integrationsweite , Abplattung und geforderte Genauigkeit in einer Funktion zusammenzufassen. Es bietet sich dabei an, die Einflüsse von und auf den in Abb. dargestellten linearen Zusammenhang zu untersuchen. Im ersten Graphen von Abb. sind zunächst verschiedene Kurven für unterschiedliche Abplattungswerte in dargestellt. Darunter folgt derselbe Graph für variierende ( ). Aus diesen beiden Abbildungen folgt, daß die Einflüsse der anderen beiden Parameter auf nichts an deren linearer Eigenschaft ändern. Im weiteren sind daher die Veränderungen der Parameter und des Polynoms erster Ordnung bei Änderung von und zu untersuchen. Für die Abhängigkeit der Steigung der Geraden von der geforderten Genauigkeit (vgl. Abb. ), ergibt sich, nach Ermittlung der über lineare Regression, der folgende logarithmische Zusammenhang :

$$a^{I\left( s\right) }_{1}\left( \varepsilon \right) =4.44\cdot 10^{-4}\cdot \varepsilon ^{-0.126}$$

Die Abhängigkeit des y-Achsenabschnittes der Kurven kann dem 2. Graphen in Abb. entnommen werden. Im Folgenden sei die Änderung von bei Variation von , unter Festhalten von , mit bezeichnet.

Abbildung: Lineare Zusammenhänge zwischen und sowie

Der Koeffizient läßt sich durch ein Polynom dritten Grades in approximieren:

$$a^{I\left( s\right) }_{1}\left( f\right) =10^{-3}\left( 0.7431+0.0447f+0.001f^{2}+0.0002f^{3}\right)$$ (2.6)

Abb. zeigt nun, analog zu Abb. , die sich für diesen Fall ergebenden Graphen.

Abbildung: Änderung der Parameter von in Abhängigkeit von

Abbildung: Änderung der Parameter von in Abhängigkeit von

Man kann zur Ableitung der Funktion 2 Wege beschreiten, zum einen kann unter Variation von untersucht werden, zum anderen ist der umgekehrte Weg möglich. Aufgrund der Einfachheit der Polynome 3. Grades wird der erste Weg beschritten, und 4 Funktionen, welche mit für bezeichnet werden, lassen sich aus den Variationen der in () eingeführten Koeffizienten des Polynoms unter Änderung von ableiten.

Abbildung: Steigungs- und y-Achsenabschnittsfunktion in Abhängigkeit von der zu fordernden Genauigkeit

In Abb. sind wiederum analog zu Abb. die sich ergebenden Steigungsfunktionen sowie die y-Achsenabschnittsverläufe für verschiedene dargestellt. Die daraus folgenden Koeffizienten der angesetzten Polynome können der Tabelle entnommen werden.

Tabelle: Polynomkoeffizienten

	0.000024	0.000017	0.000013	0.000006	0.000013	0.000011	0.000006	0.000004	0.000003
	-0.00112	-0.00080	-0.00064	-0.00026	-0.00082	-0.00059	-0.00023	-0.00013	-0.00014
	0.04785	0.03450	0.02811	0.01695	0.02392	0.01656	0.00921	0.00845	0.00855
	0.71968	0.55341	0.41451	0.32873	0.21614	0.16594	0.12830	0.08221	0.04088

Abb. zeigt eine graphische Darstellung der aus Tabelle . Wie man dieser entnehmen kann, lohnt es aufgrund der Sprünge nicht, die sich ergebenden Zusammenhänge funktional zu approximieren.

Abbildung: Darstellung der Approximationspolynome

Die im Programm implementierte Iterationsfunktion wird folglich für die diskreten Werte , wie sie der Tabelle zu entnehmen sind, bestimmt. Die Möglichkeit der Wahl der Genauigkeit zwischen und reicht für die Anwendung des Programms vollkommen aus.

Nach diesen Ergebnissen der Untersuchung der Iterationsfunktion erfolgt die automatische Berechnung der, für eine vorgegebene Genauigkeit benötigten, Zahl der Iterationen des Differentialgleichungslösers, indem zunächst mit dem zu erreichenden Genauigkeitswert aus Tabelle die zugehörigen Polynomkoeffizienten bestimmt, dann der Funktionswert dieses Polynoms an der Stelle der gegebenen Abplattung ermittelt und mit diesem die Funktion ausgewertet wird. Wie aus Abb. unten ersichtlich, ist es angebracht den y-Achsenabschnitt, der scheinbar keinem funktionalen Zusammenhang folgt ,als Konstante mit dem Wert einzuführen. Damit ergibt sich folgende Iterationsfunktion:

$$I\left( s,f,\varepsilon \right) =ceil\left( \left( a_{0}\lef... ...a_{3}\left( \varepsilon \right) f^{3}\right) \cdot s+2\right)$$ (2.7)

(Diese Funktion wurde ausreichend verprobt.^2.7)

Abbildung: Graphische Darstellung der Polynomkoeffizienten

2.3.4 Diskussion der mit dem Programm guiawp_graphik.m zu erzielenden Ergebnisse

Um die Leistungsfähigkeit des Programms zu testen, werden die in ARINGER (1994) aufgeführten Beispiele berechnet. In Tabelle 8 seiner Dissertation vergleicht er die Ergebnisse verschiedener Integratoren. Da das hier benutzte Programm nach dem Runge-Kutta-Verfahren 7. Ordnung arbeitet, werden die Ergebnisse mit denjenigen verglichen, die ARINGER nach dem, im weiteren mit RKF7 abgekürzten, Runge-Kutta-Fehlberg-Verfahren 7. Ordnung ermittelte. Der Unterschied der beiden Verfahren liegt in der, beim RKF7 Verwendung findenden Fehlerschätzung^2.8 nach FEHLBERG. Bei geod_haupt_1.m ist diese Fehlerschätzung durch die im vorigen Abschnitt abgeleitete empirische Iterationsfunktion ersetzt worden, so daß das Ergebnis ebenfalls nach Vorgabe einer gewünschten Genauigkeit ermittelt werden kann^2.9. in Tabelle sind nun die Ergebnisse des RKF7^2.10 zum Vergleich in der letzten Spalte angeordnet. Bei der Berechnung über geod_haupt_1.m wurde zunächst die Zahl der Iterationen auf 35 festgelegt. Die Werte in den Klammern enthalten im Vergleich hierzu die Änderungen, die sich nach Vorgabe einer Genauigkeit von über die Auswertung der empirischen Iterationsfunktion zeigten. Die Zahl der sich für diese Genauigkeit nach () ergebenden nötigen Iterationen ist zusätzlich angegeben. Der Tabelle kann entnommen werden, daß die Ergebnisse für die Koordinaten sich bei Längen bis zu nicht unterscheiden. Werden die Entfernungen jedoch größer treten hier Abweichungen im Bereich einiger zehntel Millimeter auf. Der Azimut im Endpunkt wird nur für identisch erhalten. Für treten Abweichungen im Bereich von Sekunden auf. Da Abweichungen in dieser Größenordnung auch bei Aringers Vergleich der verschiedenen Integratoren vorkommen^2.11, ist es schwierig auf deren Ursache zu schließen. Ein möglicher Grund könnte in der unterschiedlichen Größenordnung der Anfangswerte liegen, die jedoch durch eine Änderung der Skalierung umgangen werden kann.

Tabelle: Beispiele für die 1. geodätische Hauptaufgabe

	Anfangspunkt	Ergebnisse nach guiawp_graphik.m	Ergebnisse nach RKF7
Beispiel 1
	4045649.321	4110102.0794 (*.*)	4110102.079
	713357.130	759450.1051 (*.*)	759450.105
	4862882.427	4801881.8161 (*.*)	4801881.816
	160	160 21'24.925395'' (*.*)	160 21'24.925395''
	100000
	3
Beispiel 2
	4892928.819	4202463.8663 (*.*)	4202463.866
	0	172937.6598 (*.*)	172937.660
	4078053.805	4778979.7695 (*.*)	4778979.770
	10	11 39'15.778910'' (*.*)	11 39'15.778910''
	1000000
	3
Beispiel 3
	2662030.662	-2185590.6824 (*.*)	-2185590.682
	469387.829	13443.8838 (*.*)	13443.884
	5757839.704	5972056.2814 (*.*)	5972056.281
	5	173 48'43.328944'' (*.*)	173 48'43.328938''
	5000000
	6
Beispiel 4
	1894899.301	3786572.2845 (*.*)	3786572.285
	-1094020.622	-4701166.891 (*.8909)	-4701166.890
	5971179.947	-2053338.1932 (*.*)	-2053338.194
	200	187 07'17.004099'' (*.*)	187 07'17.004102''
	10000000
	10
Beispiel 5
	4045649.321	-757346.0989 (*.*)	-757346.098
	713357.130	2808023.9817 (*.9816)	2808023.982
	4862882.427	-5657616.2719 (*.*)	-5657616.272
	140	114 46'41.483865'' (*.483863)	114 46'41.483900''
	15000000
	14
Beispiel 6
	1405039.264	-1477765.6895 (*.6894)	-1477765.690
	-3860313.652	3782381.2759 (*.2758)	3782381.275
	-4862882.427	4902141.1762 (*.1760)	4902141.176
	310	230 48''04.1889954 (*.1889935)	230 48''04.188997
	19900000
	17

2.4 Lösung der 2. geodätischen Hauptaufgabe nach dem Schießverfahren

2.4.1 Prinzip des Schießverfahrens

Wie im vorigen Abschnitt gezeigt, lassen sich nichtlineare Anfangswertprobleme auf einfache Weise durch numerische Integration lösen. Bei Randwertaufgaben ist dies im allgemeinen nicht der Fall, da es sich meist um Probleme mit freiem Rand, d.h. um Probleme handelt, bei denen das Integrationsintervall nicht vorgegeben, sondern Teil der Lösung des Randwertproblems ist. Betrachtet man das Gleichungssystem () so trifft dies hier ebenfalls zu. Die gesuchte Entfernung tritt als Integrationsvariable auf. Eine Möglichkeit der Lösung dieses Problems besteht darin, ein anderes gleichwertiges System aufzustellen, in dem das Integrationsintervall einer neu einzuführenden unabhängigen Variablen gegeben ist. ARINGER stellt ein solches System für die Differentialgleichung der geodätischen Linie in kartesischen Koordinaten in seiner Dissertation^2.12 vor.
Meist werden solche Systeme jedoch nach einem als empirisch zu bezeichnenden Verfahren, dem s.g. Schießverfahren gelöst. Voraussetzung für dieses Verfahren sind gute Näherungswerte für die Lösung des Randwertproblems. Im ersten Schritt wird mit den Näherungswerten eine Anfangswertaufgabe gelöst. Im Fall der zweiten geodätischen Hauptaufgabe sind die Koordinaten des Anfangs-, sowie die des Endpunktes gegeben. Man benötigt dann sowohl Näherungswerte für die Tangentenrichtung an die geodätische Linie im Anfangspunkt, als auch für die zu bestimmende Entfernung. Nach Berechnung des AWP's können die Sollwerte der Endpunktkoordinaten mit der Lösung verglichen und eine Abweichung bestimmt werden.^2.13 Um das Prinzip mathematisch zu verwirklichen, werden homogene Randbedingungen der Form aufgestellt. Hierbei bezeichne die Lösung der Anfangswertaufgabe mit dem Näherungsvektor in der Iteration und die vorgegebenen Randwerte für das RWP (geg. Endpunktkoordinaten). Im Falle der zweiten geodätischen Hauptaufgabe ergeben sich die folgenden Randbedingungen:

$$\left[ \begin{array}{c} R_{1}\\ R_{2}\\ R_{3}\\ R_{4} \en... ... =\left[ \begin{array}{c} 0\\ 0\\ 0\\ 0 \end{array}\right]$$ (2.8)

mit und

Als 4. Randbedingung wird der Betrag des Tangenteneinheitsvektors eingeführt. Damit ist () eine nichtlineare Gleichung in , die mit Hilfe des allgemeinen Newtonverfahrens^2.14 iterativ gelöst werden kann. Die Bestimmung der Näherungswerte für die te Iteration erfolgt über:

$$\mathbf{L}_{0}^{i+1}=\mathbf{L}_{0}^{i}-\mathbf{d}^{i}$$ (2.9)

mit den Verbesserungen

$$\mathbf{d}^{i}=\left( \frac{\partial \mathbf{R}}{\partial \m... ...\right) ^{-1}\cdot \mathbf{R}\left( \mathbf{L}_{0}^{i}\right)$$ (2.10)

Da die Differentation in () nicht explizit durchgeführt werden kann, muß an dieser Stelle der Differentialquotient durch den Differenzenquotienten ersetzt werden.

Abbildung: Prinzip des Schießverfahrens

Im Startvektor , welcher durch eine geeignete Näherung zu ermitteln ist, wird die -te Komponente um einen noch festzulegenden Wert variiert und erneut eine Anfangswertaufgabe berechnet. Nach der Iteration über erhält man damit ein lineares Gleichungssystem für die Ermittlung des Vektors der Verbesserungen . Die Startwerte für den 'ten Schritt erhält man über (). Nach einer erneuten Berechnung des AWP's kann ein Abbruchkriterium () derart eingefügt werden, daß falls die euklidsche Distanz des berechneten Endpunktes vom Sollwert eine Schranke unterschreitet, das Schießverfahren beendet wird. Die Abb. zeigt das Prinzip des Schießverfahrens am Beispiel der 2. geodätischen Hauptaufgabe. Es sind im weiteren noch geeignete Näherungswerte für den Startvektor und die Wahl der Größenordnung der Variationen der Näherungsvektoren sowie deren Änderungen in der Iterationsphase zu bestimmen. Bei den Größenordnungen der Variationen handelt es sich um empirisch abzuleitende Werte, auf die im Hinblick auf den Umfang der Arbeit, hier nicht weiter eingegangen werden soll.

2.4.2 Berechnung der Näherungswerte

2.4.2.1 Berechnung der Näherungswerte für den Tangentenvektor im Anfangspunkt

Eine geeignete Näherung für den Tangentenvektor an die geodätische Linie stellt die Projektion des räumlichen Verbindungsvektors von Anfangspunkt und Endpunkt in die Tangentialebene der Fläche im Anfangspunkt dar. Hierzu muß zunächst der Normalenvektor der Fläche im Anfangspunkt A bestimmt werden. Mit der Normierung

$$\rho _{a}=\frac{1}{\sqrt{\left( \frac{x_{a}}{a^{2}}\right) ^... ...a}}{a^{2}}\right) ^{2}+\left( \frac{z_{a}}{b^{2}}\right) ^{2}}}$$

ergibt dich dieser zu

$$\mathbf{n}_{a}=\rho _{a}\left[ \begin{array}{c} \frac{x_{a}}... ... \frac{y_{a}}{a^{2}}\\ \frac{z_{a}}{b^{2}} \end{array}\right]$$

Der orthogonale Abstand des Endpunkts von der Tangentialebene im Anfangspunkt bestimmt sich aus dem Skalarprodukt

$$p_{\perp }=\left[ \begin{array}{c} x_{b}-x_{a}\\ y_{b}-y_{a}\\ z_{b}-z_{a} \end{array}\right] ^{T}\cdot \mathbf{n}_{a}$$

und der Tangentenvektor selbst über

$$\mathbf{t}_{a}=\frac{\left( \mathbf{x}_{b}-\mathbf{x}_{a}\ri... ...ft\vert p_{\perp }\right\vert \cdot \mathbf{n}_{a}\right\vert }$$

Die Komponenten dieses Vektors bilden die Näherungswerte für den ersten Schritt des Schießverfahrens.

Abbildung: Skizze zur Bestimmung des Tangentenvektors an direkte Verbindung von A und E

2.4.2.2 Berechnung der Näherungswerte für die Länge der geodätischen Linie zwischen Anfangs- und Endpunkt

Die einfachste Möglichkeit besteht in der Annäherung der geodätischen Linie über einen Kreisabschnitt. Dazu werden aus den kartesischen Koordinaten nach bekannten Verfahren die ellipsoidischen Koordinaten bestimmt und über den Kosinussatz der sphärischen Geometrie, unter Betrachtung der ellipsoidischen Koordinaten als Kugelkoordinaten, sowie der Annahme , ein Näherungswert für die Länge der geodätischen Linie ermittelt. Es zeigt sich, daß diese Näherung jedoch, aufgrund der lokalen Konvergenz des Schießverfahrens, nur ausreichend ist für Strecken bis 5000km, so daß im weiteren die Länge eines Vertikalschnittes als Näherung Verwendung findet. Die Vertikalschnittebene wird von und in A aufgespannt. Um die Berechnung zu vereinfachen reicht es aus, das Problem in einem 2-D Koordinatensystem, in dem sich der Vertikalschnitt als ebene Kurve eindeutig darstellen läßt, zu betrachten. Ohne näher auf die Transformation einzugehen, sei erwähnt, daß der Ursprung des Koordinatensystems in den Anfangspunkt A verschoben wurde, die neue x-Achse in Richtung und die neue z-Achse in Richtung weist. (vgl. hierzu Abb. ) Es ergibt sich die folgende Gleichung für die Vertikalschnittkurve aus der Gleichung des Ellipsoides im gedrehten ebenen System:

$$\widehat{x}^{2}b_{11}+\widehat{z}^{2}b_{33}+2\widehat{x}\widehat{z}b_{13}+2\widehat{x}c_{1}+2\widehat{z}c_{3}=0$$

Die Auflösung der quadratischen Gleichung führt zu

$$\widehat{z}\left( \widehat{x}\right) =\frac{-\left( \widehat... ...widehat{x}^{2}b_{11}+2\widehat{x}c_{1}\right) b_{33}}}{b_{33}}$$ (2.11)

und

$$\widehat{x}\left( \widehat{z}\right) =\frac{-\left( \widehat... ...widehat{z}^{2}b_{33}+2\widehat{z}c_{3}\right) b_{11}}}{b_{11}}$$ (2.12)

Abbildung: Intervallteilung zur Berechnung des Vertikalschnitts

Zur Berechnung der Länge des Verikalschnitts erfolgt eine Einteilung der z-Achse des ebenen Systems gem. Abb. in Teilintervalle . Über () werden die zugehörigen Koordinatendifferentiale ermittelt und mit den erhaltenen Werten dann die Streckendifferentiale . Die gesuchte Strecke ergibt sich demzufolge nach

$$\displaystyle s=\sum _{i=1}^{n}\sqrt{d\widehat{x}_{i}^{2}+d\widehat{z}_{i}^{2}}$$

Die Abbildung zeigt, daß es numerisch sinnvoll ist, das gesamte Integrationsintervall in 3 Teilintervalle zu ordnen. Im ersten, sowie im dritten wird , im zweiten als unabhängige Variable gewählt. Damit findet bei der Berechnung der Koordinatendifferentiale im zweiten Intervall Gleichung (), in den übrigen Gleichung () Verwendung. Im ersten Intervall erhält () das Vorzeichen ``+'', im dritten ``-``. Die Genauigkeit der Berechnung ist abhängig von der Zahl der Teilintervalle bzw. .

Mit dem, auf diese Weise ermittelten Näherungswert für die Länge der geodätischen Linie zwischen Anfangspunkt A und Endpunkt E auf einem durch seine Parameter und gegebenem Ellipsoid, kann die 2. geodätische Hauptaufgabe nach dem Schießverfahren gelöst werden^2.15.

2.4.3 Praktische Realisierung eines Programms zur Lösung der 2. geodätischen Hauptaufgabe nach dem Schießverfahren

2.4.3.1 Beschreibung der Benutzeroberfläche

Die Vorgehensweise des MATLAB-Programms guirwp_graphik.m zur Berechnung der 2. geodätischen Hauptaufgabe entspricht weitestgehend den in den Abschnitten zuvor erwähnten Verfahren. Die Anfangswertaufgaben werden mit RK7multi.m gelöst und zusätzlich wird zu Beginn, um die Intervallteilung bei der Berechnung des Vertikalschnitts durchzuführen, die Länge eines Kreisabschnitts ermittelt. Das zur Bedienung des Programms erstellte GUI (siehe Abb. ) enthält Eingabefelder für

• vorzugebende kartesische Koordinaten des Anfangs- und Endpunktes als Randwerte
• Maße des zugrundezulegenden Ellipsoides
• Schranke für den euklidschen Abstand der ermittelten Endpunktkoordinaten von den vorgegebenen Randwerten als Abbruchkriterium für das Schießverfahren
• die zu wählende Intervallbreite bei der Berechnung der Länge des Vertikalschnittes
• Größe der relativen Variationen von (in der Oberfläche mit dsigma bezeichnet) in der ersten Iteration des Schießverfahrens
• einen Änderungsfaktor für die relativen Variationen

Ausgegeben werden

• die Ergebnisse der Berechnung des Vertikalschnittes
• die Zahl der benötigten Iterationen
• der über die Kreisabschnittsberechnung abgeschätzte Entfernungsbereich
• die Länge, sowie die Azimute der geodätischen Linie zwischen A und E

2.4.3.2 Beschreibung des Berechnungsablaufs

Zunächst sei erwähnt, daß alle auftretenden Routinen als Funktionen realisiert sind, d.h., daß die in den Zwischenschritten berechneten Ergebnisse am Ende des Programmablaufs keiner weiteren Anwendung mehr zur Verfügung stehen^2.16. Um die Berechnung zu starten, kann entweder an der Eingabeaufforderung das Script guirwp_graphik.m aufgerufen werden oder man startet das komplette Paket durch den Aufruf von geod_haupt.m. Beide Scripte dienen dem Aufbau der graphischen Benutzeroberflächen. Mit der Berechnung der 2. geodätischen Hauptaufgabe, wird nach Eingabe der Randwerte und der Steuerparameter in die Oberfläche, durch Anklicken des Startbuttons mit automatischem Aufruf der Funktion geo_haupt_2_graphik.m begonnen. Diese Funktion besitzt keine Ausgabewerte sondern dient nur zur Berechnung der Randwertaufgabe sowie als Schnittstelle für das GUI. Nach Auslesen der eingegebenen Daten wird zunächst das Fenster initialisiert und eine erste Abschätzung der Entfernung wird über die Länge eines Kreisbogens vorgenommen. Diese Abschätzung dient der Näherungswertberechnung über den Vertikalschnitt zur Bestimmung der Intervallteilung^2.17. Die Berechnung erfolgt durch die Funktion kreisbogen.m welche die Länge eines Kreisbogens zwischen den Vektoren des Anfangs- und des Endpunktes bezüglich einer Kugel mit dem Radius ermittelt.

Abbildung: Flußdiagramm zur Berechnung der 2. geodätischen Hauptaufgabe

Nach dieser groben Ermittlung wird die Funktion vertikalschnitt.m, die die Länge des Vertikalschnittes und dessen Richtung am Anfangs- sowie am Endpunkt berechnen soll, aufgerufen. Übergeben werden hier nur die Koordinaten des Anfangs- und Endpunktes, die übrigen Steuerparameter wurden in der aufrufenden Funktion zuvor als global deklariert. Wie im Abschnitt auf Seite beschrieben, müssen, vor der Berechnung der Strecken, die Tangentenvektoren für die Transformation der Ellipse in das ebene System bestimmt werden. Dies geschieht durch die Funktion kart_lokal.m, die als Eingabewerte wiederum die Koordinaten der Randpunkte verarbeitet. Da die Tangentenvektoren weiterhin verwendet werden, besitzt kart_lokal.m ebenfalls keine Ausgabeparamter, sondern ändert nur die von vornherein als global deklarierten Variablen. Die Berechnung der Länge des Vertikalschnittes erfolgt nach Ablauf von kart_lokal.m wieder in vertikalschnitt.m. Je nach Lage des Endpunktes müssen unter Berücksichtigung der über die Funktion kreisbogen.m erhaltenen Abschätzung die Streckenanteile in den einzelnen Integrationsintervallen berechnet werden. Ohne weiter auf die Berechnung, die schon zuvor beschrieben wurde, einzugehen, sei erwähnt, daß diese von einzelnen kleineren Hilfsfunktionen, deren Quellcode der beigefügten Diskette entnommen werden kann, durchgeführt wird. Zum Schluß der Berechnung werden die zur Lösung des Anfangswertproblems nötigen Näherungswerte in einem Vektor an die Ursprungsfunktion geo_haupt_2_graphik.m übergeben. Vor dem Start der Iteration erfolgt dann dort, mit Initialisierung der Abbruchbedingung (siehe Seite ), die erste Berechnung des AWP's mit den o.g. Näherungswerten. Die Berechnung der Randbedingungen (), bzw. der Soll-Ist Abgleich, wird von der gleichnamigen Funktion randbedingungen.m durchgeführt. Danach startet das Schießverfahren mit den einzelnen Berechnungen der Anfangswertaufgaben, wobei in jedem Schritt auf 2 verschiedene Arten eine Überprüfung auf Divergenz, sowie eine Abfrage nach Zutreffen des Abbruchkriteriums erfolgt. Wird dieses erfüllt endet die Iteration, die Ergebnisse werden an die Benutzerschnittstelle übergeben und es erscheint die Meldung Berechnung beeendet.

Abbildung: Graphische Benutzeroberfläche zur Lösung der 2. geodätischen Hauptaufgabe

2.4.4 Diskussion der mit dem Programm guirwp_graphik.m zu erzielenden Ergebnisse

In Tabelle werden die Ergebnisse, wie auch schon bei der ersten Hauptaufgabe, Aringers Daten gegenübergestellt. Die der Berechnung zugrundeliegenden Ellipsoidparameter waren

$$a=6378388\; m$$

und

$$b=6356911.94613\; m$$

Die Randwerte für die Berechnung über guirwp_graphik.m, bildeten die von Aringer mit 3 Nachkommastellen angegebenen Koordinaten. Die vorzugebende Genauigkeit zur Lösung der Anfangswertaufgaben während des Schießverfahrens wurde mit angesetzt. Das Programm beendete die Iteration, wenn der berechnete Endpunkt um weniger als vom Sollwert entfernt war. Die für die Vertikalschnittberechnung nötige Intervallbreite betrug . Für die relative Größe der Variationen der Komponenten des Tangenteneinheitsvektors, hat sich ein Wert von , für Variation der der Strecke von als zweckmäßig erwiesen^2.18. Der Veränderungsfaktor wurde mit angesetzt. An dieser Stelle sei darauf hingewiesen, daß alle Veränderungen der mit den Variationen der Anfangswerte zusammenhängenden Steuerparameter, einen sehr großen Einfluß auf das Konvergenzverhalten des Schießverfahrens haben. Um hier jedoch einen funktionalen Zusammenhang abzuleiten, bedarf es weiterer Untersuchungen die den Umfang dieser Arbeit überschreiten.

Wie der Tabelle zu entnehmen ist, differieren die Ergebnisse für die Azimute im Bereich von Sekunden. Die Entfernung wird im Gegensatz dazu identisch erhalten. Es verwundert, daß Aringer bei den Beispielen 2 - 6 das Azimut auf die 7-te Nachkommastelle genau betimmen kann (vgl. Beispiele für die 1. geodätische Hauptaufgabe), obwohl er die Koordinaten des Endpunktes, die er aus der Berechnung der Anfangswertaufgaben ermittelt hat, nur mit 3 Nachkommastellen ansetzt. Das führt dazu, daß der Vergleich der Ergebnisse nicht besonders aussagekräftig ist.

Um die Integrität des Programmsystems zu prüfen, wurden einige Anfangswertaufgaben berechnet und die Ergebnisse ungerundet in das Programm zur Lösung der Randwertaufgabe übernommen. Dabei konnte keine Differenz, außerhalb des Rahmens der Rechengenauigkeit, festgestellt werden.

Tabelle: Beispiele für die Lösung der 2. geodätischen Hauptaufgabe

Randwerte	über geod.Linie		über Vertikalschnitt
Beispiel 1	guirwp_graphik.m	Aringer	guirwp_graphik.m	Aringer
	159 59'59.999985''	... ..'59.999984''	159 59'59.992231''	... ...'59.992362''
	160 21'24.925248''	... ..'24.925379''	160 21'24.940530''	... ..'24.940661''
	100000	100000	1000000	100000
	benötigte Iterationen: 2	benötigte Iterationen: 3
	Iterationen zur Lösung der AWP's: 3
Beispiel 2
	10 00'00.000026''	... ..'00.000000''	10 00'00.555657''	... ..'00.555635''
	11 39'15.778941''	... ..'15.778911''	11 39'14.690175''	... ..'14.690145''
	1000000	1000000	1000000	1000000
	benötigte Iterationen: 2	benötigte Iterationen: 3
	Iterationen zur Lösung der AWP's: 3
Beispiel 3
	5 0'00.0000047''	. ..'00.000000''	5 00''00.01311050''	. ..'01.311046''
	173 48'43.328933''	... ..'43.328938''	173 48'41.789155''	... ..'41789160''
	4999999.999	5000000	4999999.993	5000000
	benötigte Iterationen: 3	benötigte Iterationen: 3
	Iterationen zur Lösung der AWP's: 6
Beispiel 4
	199 59'59.999965''	200 00'00.000000''	200 00'59.157638''	... ..'59.157672''
	187 07'17.004081''	... ..'17.004101''	187 05'18.805158''	... ..'18.805171''
	10000000.001	10000000	10000000.2482	10000000.2606
	benötigte Iterationen: 4	benötigte Iterationen: 4
	Iterationen zur Lösung der AWP's:10
Beispiel 5
	140 00'00.000017''	... ..'00.000000''	139 42'53.422259''	... ..'53.422242''
	114 46'41.483912''	... ..'41.483903''	115 11'27.067537''	... ..'27.067529''
	14999999.999	14999999.999	15000052.7846	15000052.8051
	benötigte Iterationen: 5	benötigte Iterationen: 6
	Iterationen zur Lösung der AWP's: 14
Beispiel 6
	310 00'00.000633''	... ..'00.000000''	281 52'14,760975''	... ..'14.759261''
	230 48'04.188394''	... ..'04.189003''	259 09'25,780566''	... ..'25.782254''
	19899999,998	19899999.999	19911966,8419	19911966.8693
	benötigte Iterationen: 12	benötigte Iterationen: 11
	Iterationen zur Lösung der AWP's: 17

Abbildung: Aufbau des Programmsystems

2.5 Beschreibung des Programmsystems

Abbildung: Benutzeroberfläche des Programmsystems

Der Aufgabenstellung entsprechend wurde unter MATLAB ein Programmsystem^2.19 erstellt, das einem Nutzer erlaubt, geometrische, wie auch numerische Eigenschaften von geodätischen Linien auf einem Rotationsellipsoid zu untersuchen. Um auch im Umgang mit MATLAB unerfahrenen Benutzern die Anwendung zu ermöglichen, wurde das System mit einer graphischen Benutzeroberfläche (GUI) versehen. Diese erlaubt einen interaktiven Zugriff auf MATLABs interne Grafikfunktionen, ohne daß ein Benutzer sich vorher intensiv mit der Programmiersystem auseinandersetzen muß. Abb. zeigt das Hauptfenster des Systems, wie es sich, nach Eingabe von geod_haupt am MATLAB-Prompt, öffnet. Mit Hilfe dieser Oberfläche lassen sich nun folgende Aufgaben lösen:

1. Berechnung der 1. geodätischen Hauptaufgabe
2. Berechnung der 2. geodätischen Hauptaufgabe
3. Berechnung von Punkten auf einer geodätischen Linie und deren Abspeicherung in einer Datei
4. Eine dreidimensionale Darstellung des Verlaufs von geodätischen Linien^2.20 auf dem zugrundeliegenden Ellipsoid
5. Eine dreidimensionale Visualisierung eines Parallelkoordinatensystems in Soldnerscher Anordnung

Da die Berechnung der Hauptaufgaben schon in den vorangegangenen Abschnitten erläutert wurde, sei im folgenden kurz auf die Punkte 3-5 eingegangen.

2.5.1 Die Darstellung und Abspeicherung der Ephemeriden

Das Programm zur Darstellung der Ephemeriden wird jeweils nach der Berechnung der ersten bzw. zweiten geodätischen Hauptaufgabe aus den Oberflächen heraus gestartet. Nach der Wahl des Tastenfeldes Ephemeridenberechnung öffnet sich das Fenster in Abb. . Hier wird zum einen die gewünschte Zahl der Stützpunkte^2.21 der geodätischen Linie eingetragen, zum anderen besteht die Möglichkeit zwischen zwei Arten der Darstellung des Ellipsoides zu wählen. In der Gitterdarstellung werden nur die, über eine Hilfsfunktion berechneten Stützpunkte eines geographischen ellipsoidischen Koordinatensystems mit Linien verbunden, dargestellt. Über die Zahl der Stützstellen läßt sich frei verfügen. Bei der Flächendarstellung hingegen werden diese Gitterelemente farbig ausgefüllt, was dazu führt, daß die später zwischen die Stützpunkte gezeichneten Verbindungslinien nicht mehr komplett sichtbar sind.

Abbildung: Auswahlfenster zur Ephemeridenberechnung

Die Genauigkeit, mit der die Zwischenpunkte berechnet werden, richtet sich nach der Genauigkeit, die für die Berechnung der Hauptaufgaben angesetzt wurde. Zur Verkürzung der Rechenzeit wird die empirische Iterationsfunktion für jeden Stützpunkt einzeln ausgewertet. Die Endpunktkoordinaten werden außerdem zur Kontrolle ein zweites Mal ermittelt. Nach der Bestätigung der Eingabe in Fenster Abb. startet das Script ephemeridenberechnung.m, es öffnet sich das Fenster mit der dreidimensionalen Darstellung des Ellipsoides, sowie das Fenster in Abb. . Diesem kann während der Berechnung die Nummer des aktuell ermittelten Punktes entnommen werden. Parallel dazu werden im Grafikfenster die Stützpunkte als rote Kreuze geplottet.

Abbildung: GUI zum Steuern der Ephemeridendarstellung und Abspeichern der Stützpunktkoordinaten

Die Richtung des Blickpunktes wird automatisch in Richtung des Anfangspunktvektors gesetzt, so daß die Linie vollständig im Blickfeld des Betrachters erscheint. Anfangs- sowie Endpunkt werden durch einen, grün dargestellten, räumlichen Vektor hervorgehoben.

Am Ende der Berechnungen können die kartesischen Koordinaten der Zwischenpunkte in einer binären MATLAB-Workspace (*.mat) Datei abgespeichert werden^2.22. Hierbei ist zu beachten, daß falls mehrere Linien in einer Datei zusammengefaßt werden, diese sich auf dasselbe Ellipsoid beziehen sollten, da die Ellipsoidparameter ebenfalls dort abgelegt werden. Der Dateiname muß mit voller Pfadangabe, aber ohne Erweiterung angegeben werden.

Über die Grafikoptionen im unteren Bereich des Steuerungsfensters läßt sich der Blickpunkt sowie ein gewünschter Zoomfaktor für die Darstellung festlegen. Abbildung zeigt das graphische Ergebnis einer Ephemeridenberechnung mit 20 Stützstellen bei einer Länge der geodätischen Linie von 15000km.

Abbildung: Darstellung einer Ephemeride

2.5.2 Auslesen und Plotten der Datensätze

Im Hauptfenster gelangt man durch Wahl des Schalters geod_plot zu einer Möglichkeit mehrere, im Anschluß an die Berechnung der Hauptaufgaben abgespeicherte Linienverläufe,

Abbildung: Darstellungsschnittstelle

darzustellen. (Abb. zeigt die zugehörige Oberfläche) Durch Angabe des gewünschten Dateinamens im Rahmen Datensatz und anschließender Wahl des Schalters Laden wird der Inhalt der Datei im darunter liegenden Auswahlfenster angezeigt. Sollen mehrere Linien innerhalb einer Grafik dargestellt werden, so muß unter den Plotoptionen das Feld hinzufügen angewählt sein. Das im Anschluß daran vorliegende dreidimensionale Modell kann mit Hilfe der Maus beliebig um seine 3 Achsen gedreht werden. Die Art der Darstellung der Ephemeride ist frei wählbar.

Tabelle: Zuordnung der Linientypen

Symbol	Linientyp
-	durchgezogen
-	gestrichelt
:	gepunktet
-.	abwechselnd
none	keine Darstellung

Tabelle: Zuordnung der Punktsymbole

Symbol	Punktbeschreibung
+	Plus
o	Kreis
*	Stern
x	Kreuz
s	Quadrat
d	Raute
h	Hexagramm
.	Punkt

Tabelle: Zuordnung der Farben

Abbkürzung	Farbe
k	Schwarz
w	Weiß
r	Rot
g	Grün
b	Blau
y	Gelb
c	Cyan
m	Magenta

Da es programmiertechnisch einfacher war in den Auswahlfenstern die MATLAB spezifischen Abbkürzungen für Farbe und Form der graphischen Objekte anzeigen zu lassen, sind die entsprechenden Zuordnungen in den Tabellen - kurz aufgeführt.

Abbildung: Darstellung mehrerer Ephemeriden nach Auslesen der Datensätze

Ein Beispiel für die Darstellung mehrerer Ephemeriden zeigt Abb. .

2.5.3 Berechnung eines Parallelkoordinatensystems in Soldnerscher Anordnung

Die zum Programmsystem gehörende Funktion Soldner.m dient der Untersuchung des Verlaufs von geodätischen Linien in der Nähe des Querpols eines globalen Parallelkoordinatensystems in Soldnerscher Anordnung.

2.5.3.1 Grundlagen

Bei einem allgemeinen geodätischen Parallelkoordinatensystem auf einem Rotationsellipsoid findet als Abzissenlinie, der durch den Koordinatenanfangspunkt verlaufende Meridian Verwendung. Die Lage eines Punktes wird durch die Bogenlänge auf der Abzissenlinie, sowie die Länge einer geodätischen Linie durch diesen Punkt, die die Abzissenlinie senkrecht schneidet, bestimmt.

Abbildung: Prinzip eines Parallelkoordinatensystems in Soldnerscher Anordnung

Dabei entspricht der Bogenlänge auf der Abzissenlinie, gemessen vom Koordinatenanfangspunkt zum Fußpunkt und der Länge der geodätischen Linie durch und . (vgl. Abb. )

und werden in diesem Zusammenhang als die ellipsoidischen Parallelkoordinaten des Punktes bezeichnet. Sie lassen sich als Gaußsche Flächenparameter auf dem Ellipsoid interpretieren. Die als Ordinatenlinien bezeichneten geodätischen Linien verlaufen auf dem Ellipsoid weder eben noch sind es geschlossene Kurven. Ebenso bilden die Abzissenlinien oder geodätische Parallelen, zwar geschlossene, jedoch keine ebenen Kurven. An Stelle des Querpols, der sich bei Einführung eines Parallelkoordinatensystems als Schnittpunkt der Geodäten auf der Kugel ergibt, bildet sich im Falle eines Rotationsellipsoides ein, zum Äquator symmetrisches Gebiet, das durch eine asteroidenförmige Flächenkurve begrenzt wird. Mit Hilfe des im folgenden beschriebenen Programms läßt sich diese Figur anhand eines dreidimensionalen Gittermodells des Koordinatennetzes in Abhängigkeit gewisser Parameter graphisch darstellen.

2.5.3.2 Beschreibung des Programms

Abbildung: GUI zur Darstellung eines Soldner'schen Koordinatensystems

Das Programm kann entweder von der Eingabeaufforderung durch Aufruf des Skriptes soldner_graphik.m oder im Hauptfenster des Programmsystems, durch Wahl des Schalters soldner, gestartet werden. Über das Steuerfenster (siehe Abb. ) lassen sich die Ellipsoidparameter, die Zahl der Stützpunkte des Koordinatensystems, sowie die Genauigkeit der Berechnung festlegen. Die Berechnung wird durch Auswahl des Buttons Plotten, der das Script Soldner.m aufruft, begonnen. Zunächst wird der Bezugsmeridian entsprechend der Wahl der Breitenteilung in Intervalle aufgeteilt. Dabei sei anzumerken, daß durch eine solche Teilung zwar keine gleichen Längenteilung des Bezugmeridians erreicht wird, jedoch ist dies für die Darstellung der Figur auch nicht erforderlich. Da es sich bei dem zugrundegelegten Ellipsoid um einen rotationssymetrischen Körper handelt, reicht es aus, die Stützpunkte des Koordinatennetzes für nur einen Oktanten dieses Körpers zu berechnen und die anderen durch Spiegelungen abzuleiten. Im Programm selbst werden für die einzelnen Stützpunkte des Bezugsmeridians sukzessive Anfangswertaufgaben für ein Azimut von und einer Streckenteilung gemäß der Auswahl im Fenster für ein Intervall von gelöst. Dabei ist zu beachten daß, falls ein Rotationsellipsoid einer anderen Größenordnung als der eines Erdellipsoides gewählt wird, dieser Wert im Quellcode sinnvoll abzuändern ist. Die Auswahl der Genauigkeit dient in diesem Fall einer Steigerung der Rechengeschwindigkeit. Auch hier wird die Iterationsfunktion für jedes einzelne AWP erneut ausgewertet und die auf diese Weise ermittelte Zahl der notwendigen Iterationen wird der Berechnung zugrundegelegt. Tabelle zeigt jeweils 4 verschiedene Ansichten eines Modells, wie es sich nach Ablauf der Berechnung dem Betrachter darstellt.

Tabelle: Koordinatennetz eines Parallelkoordinatensystems in Soldnerscher Anordnung

2.6 Beispiel zur Anwendung des Programmsystems

Nach dem Satz von Clairaut läßt sich der Verlauf einer geodätischen Linie auf einem rotationsymmetrischen Körper beschreiben mit

$$p\cdot \sin \left( \alpha \right) =const$$ (2.13)

entspricht dabei dem Parallelkreishalbmesser, der sich ausgedrückt, durch die reduzierte Breite , darstellen läßt über

$$p=a\cdot \cos \left( \beta \right)$$

Gilt , und bezeichne das Azimut der geodätischen Linie am Äquator, so folgt und für () folgt

$$a\cdot \sin \left( \alpha _{\ddot{A}}\right) =const=p_{min}$$

Die Konstante entspricht dabei dem minimalen Parallelkreishalbmesser, der von der geodätischen Linie nicht unterschritten wird. Diesem entspricht die maximale reduzierte Breite ,die nicht überschritten wird. Damit schwingt die geodätische Linie, mit einer von der Abplattung und dem Azimut abhängigen Periode, symmetrisch zum Äquator zwischen den reduzierten Breiten und .

Tabelle:

Bezeichnung	Farbe
Linie 1	schwarz
Linie 2	rot
Linie 3	grün
Linie 4	blau

Mit Hilfe des Darstellungsprogramms für den Verlauf geodätischer Linien läßt sich dieses Verhalten graphisch veranschaulichen. Es wurden für ein Ellipsoid mit den Parametern

$$\begin{eqnarray*} a & = & 6378388\; m\\ b & = & 4500000\; m \end{eqnarray*}$$

4 Lösungen der ersten geodätischen Hauptaufgabe für die Azimute und bestimmt. Das Integrationsintervall betrug dabei . Der Ortsvektor des Startpunktes lautet

$$\mathbf{x}=\left[ \begin{array}{c} 6378388\; m\\ 0\; m\\ 0\; m \end{array}\right]$$

Abbildung bietet eine dreidimensionale Aufsicht der sich ergebenden Darstellung. Zur besseren Unterscheidung der Linien seien zusätzlich die 2 Projektionen in Form der Abbildungen und beigefügt. Die farbliche Zuordnung der Linien, sowie die zugehörigen Werte für die maximalen Breiten sind Tabelle ^2.23 zu entnehmen.

Abbildung: Veranschaulichung der geometrischen Bedeutung des Satzes von Clairaut für den Verlauf der geodätischen Linien

Abbildung: Projektion der Abb. in die -Ebene

Abbildung: Projektion der Abb. in die -Ebene

2.7 Anmerkungen zum Anhang

Der Arbeit sind als Anhang die wichtigsten MATLAB-Scripte, die schon im Flußdiagramm bezeichnet wurden, beigefügt. Zur Programmierung der Oberfläche kann leider kein Quellcode abgegeben werden, da diese interaktiv, per Drag & Drop, geschieht und hierbei zusätzlich entstehender Code in binären Dateien abgespeichert wird. Als Beispiel für ein, bei der Oberflächenprogrammierung automatisch von MATLAB erstelltes Script, kann der Ausdruck der Datei liniendarstellung.m, welche den Aufbau der Oberfläche aus Abb. regelt, ebenfalls dem Anhang entnommen werden.

Abbildungsverzeichnis

• . Geometrische Interpretation des Affinen Tangentialraums
• . Verlauf der geodätischen Linie auf einer gekrümmten Fläche
• . Körper gleitet unter Einfluß der Schwerkraft eine Ebene herab
• . Prinzip des Verfahrens nach Runge
• . Flußdiagramm der Funktion RK7multi.m
• . GUI für Lösung des Anfangswertproblems
• . Rechenzeit in Abhängigkeit von der Entfernung
• . Iterationen in Abhängigkeit von der Integrationsweite s
• . Abhängigkeiten von Azimut und Abplattung
• . Lineare Zusammenhänge zwischen und sowie
• . Änderung der Parameter von in Abhängigkeit von
• . Änderung der Parameter von in Abhängigkeit von
• . Steigungs- und y-Achsenabschnittsfunktion in Abhängigkeit von der zu fordernden Genauigkeit
• . Darstellung der Approximationspolynome
• . Graphische Darstellung der Polynomkoeffizienten
• . Prinzip des Schießverfahrens
• . Skizze zur Bestimmung des Tangentenvektors an direkte Verbindung von A und E
• . Intervallteilung zur Berechnung des Vertikalschnitts
• . Flußdiagramm zur Berechnung der 2. geodätischen Hauptaufgabe
• . Graphische Benutzeroberfläche zur Lösung der 2. geodätischen Hauptaufgabe
• . Aufbau des Programmsystems
• . Benutzeroberfläche des Programmsystems
• . Auswahlfenster zur Ephemeridenberechnung
• . GUI zum Steuern der Ephemeridendarstellung und Abspeichern der Stützpunktkoordinaten
• . Darstellung einer Ephemeride
• . Darstellungsschnittstelle
• . Darstellung mehrerer Ephemeriden nach Auslesen der Datensätze
• . Prinzip eines Parallelkoordinatensystems in Soldnerscher Anordnung
• . GUI zur Darstellung eines Soldner'schen Koordinatensystems
• . Veranschaulichung der geometrischen Bedeutung des Satzes von Clairaut für den Verlauf der geodätischen Linien
• . Projektion der Abb. in die -Ebene
• . Projektion der Abb. in die -Ebene

Literatur

1

Aringer, K. (1994): Geodätische Hauptaufgaben auf Flächen in kartesischen Koordinaten. DGK C 421.

1

Heck, B. (1993): Rechenverfahren und Auswertemodelle der Landesvermessung. H. Wichmann Verlag, Karlsruhe.

1

Heil M., Kitzka F. (1984): Grundkurs Theoretische Mechanik. B.G. Teubner, Stuttgart.

2

Heitz, S., Stöcker-Meier, E. (1998): Grundlagen der Physikalischen Geodäsie. Dümmler, Bonn.

3

Klingbeil, E. (1989): Tensorrechnung für Ingenieure. B.I.-Hochschultaschenbuch Band 197, Mannheim, Wien, Zürich.

4

Klotz, J. (1991): Eine analytische Lösung kanonischer Gleichungen der geodätischen Linie zur Transformation geodätischer Flächenkoordinaten. DGK C 385.

5

Moritz, H., Hofmann-Wellenhof, B. (1993): Geometry, Relativity, Geodesy. H. Wichmann Verlag, Karlsruhe.

6

Reinhardt F., Soeder H. (1984): dtv-Atlas zur Mathematik.Deutscher Taschenbuch Verlag, München.

6

Schürer, M. (1975): Die geodätischen Linien und die Mechanik. Vermessung, Photogrammetrie, Kulturtechnik; Fachblatt III/IV-75.

2

Spiegel, M. R. (1990): Allgemeine Mechanik. McGraw-Hill Book Company GmbH, Hamburg.

7

Script zur Vertiefungsvorlesung im Fach Astronomisch/Physikalische Geodäsie, WS 1997/98.

Über dieses Dokument ...

Theoretische und numerische Untersuchungen zu den geodätischen Hauptaufgaben

This document was generated using the LaTeX2HTML translator Version 99.1 release (March 30, 1999)

Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
Copyright © 1997, 1998, 1999, Ross Moore, Mathematics Department, Macquarie University, Sydney.

The command line arguments were:
latex2html -no_subdir -split 0 -show_section_numbers /home/hannah/Diplomarbeit/Text/Diplomarbeit.tex

The translation was initiated by root on 2000-10-19

Fußnoten

... Koordinatensystem^1.1

Ansonsten würde die Eindeutigkeit der Zuordnung der Punkte verloren gehen

... Rechenregeln^1.2

KLINGBEIL(1989)

... ^1.3

Es handelt sich hier um die Christoffelsymbole 2. Art, die Ableitung dieser 'Objekte des räumlichen Zusammenhanges' aus den Metrikkoeffizienten wird in einem folgenden Abschnitt eingehender behandelt.

... Koeffizienten^1.4

Aus MORITZ, H., HOFMANN-WELLENHOF, B. (1993).

... Integratoren^1.5

Ein Vergleich der Leistungsfähigkeit verschiedener Integratoren bei der Lösung der geodätischen Hauptaufgaben auf einem Rotationsellipsoid ist in ARINGER (1994) dargestellt.

... Prinzip^1.6

Das Prinzip () geht von den s.g. virtuellen Verschiebungen eines Massepunktes zu einem festen Zeitpunkt aus, und ist daher als differentielles Prinzip zu betrachten.

... bei^1.7

Eine gute Einführung in die unterschiedlichen Methoden der theoretischen Mechanik bietet HEIL/KITZKA (1984)

... Koordinaten^1.8

Bei diesen kann es sich auch um Winkel, Abstände oder daraus abgeleitete Größen handeln.

... besitzen^1.9

An dieser Stelle ist auf die Analogie zur Definition einer Mannigfaltigkeit hingewiesen.

... klassifiziert^1.10

Beispiel für eine Klassifikation eines mechanischen Systems: Ein Massepunkt bewege sich reibungslos auf einen langen Draht, welcher mit konstanter Winkelgeschwindigkeit um eine senkrechte Achse rotiert. In diesem Fall handelt es sich um ein rheonomes, holonomes und konservatives System.

... Zentralgleichung^1.11

Für die Herleitung dieser Gleichung sei auf HEITZ/STÖCKER-MEIER (1998) verwiesen.

... Prinzip^1.12

siehe SPIEGEL (1990)

... Extremum^1.13

Aufgrund der Tatsache, daß es sich bei diesem meist um ein Minimum handelt, wird das Prinzip oft als Hamiltonsches Prinzip der kleinsten Wirkung bezeichnet.

... abhängt^1.14

Für diese gilt daher

... Raum^1.15

Wird auch als pq - Phasenraum bezeichnet.

... Parameter^1.16

Diese werden mit den neuen Impulsen identifiziert.

... Eingesetzt^1.17

Hier sei anzumerken, daß zum besseren Vergleich gesetzt wurde

... Bestimmungsgleichung^1.18

Diese Gleichung beschreibt den Verlauf einer entfernungsmäßig kürzesten Verbindung, während in () durch die Transformation eine kürzeste Verbindung bezüglich des eingeführten Linienelementes beschrieben wird.

... Prismaregel^2.1

Für Runge-Kutta-Verfahren höherer Ordnung findet die Simpsonsche Regel als Verallgemeinerung der Prismaregel Verwendung.

... ^2.2

siehe hierzu ARINGER 1994

... Rechenzeit^2.3

siehe hierzu Abb.

... MIDEVA^2.4

Weitere Informationen zu den Compilern und freie Evaluations-Lizenzen sind erhältlich unter http://www.mathtools.com.

... Schrittweitensteuerung^2.5

Möglichkeiten zur Implementierung einer solchen sind in ARINGER 1994 aufgeführt.

...2d_abh_it_von_s_abpl^2.6

Bei der Berechnung wurde hier davon ausgegangen das der Anfangspunkt auf dem Äquator liegt. Das Azimut wird von Norden nach Osten als positiv angenommen. Der hier noch nicht variierende Parameter ist die Genauigkeit , welche für die Berechnung der Graphen mit angesetzt wurde.

... verprobt.^2.7

Es sei hier noch einmal bemerkt, daß diese Funktion nur für die diskreten Werte aus Tabelle , sowie für Abplattungen bis gilt. Ebenso wird der lineare Verlauf der Funktion auch für Entfernungen über angenommen.

... Fehlerschätzung^2.8

Die Schätzung erfolgt über 2 Näherungen der Lösung, welche jeweils über Diskretisierungsverfahren der globalen Fehlerordnungen und gewonnen werden. Praktisch wird beim RKF7 ein von Fehlberg angegebener Fehlerterm angesetzt, wobei und aus zwei weiteren Auswertungen der Differentialgleichung stammen. Damit ergibt sich der Fehlerterm aus dem Vergleich eines Verfahrens 8. mit dem eines 7. Ordnung. Über diesen wird bei ARINGER (1994) die automatische Schrittweitensteuerung implementiert.

... kann^2.9

Dies führt zu einer geringfügigen Verkürzung der Rechenzeit im Vergleich zu den Programmen mit Fehlerschätzung.

...RKF7^2.10

Aus ARINGER (1994), Tabelle 8, wobei hier zu bemerken sei, daß für diese Berechnung die Vorgabe der Genauigkeit der kartesischen Koordinaten mit angegeben wurde.

... vorkommen^2.11

Bei Vorgabe derselben geforderten Genauigkeit.

... Dissertation^2.12

siehe Gleichung (6.61)

... werden.^2.13

Um Unklarheiten die im Vergleich mit Aringers Arbeit auftreten können zu vermeiden, sei erwähnt, daß hier nur das normale Differentialgleichungssystem zur Lösung der Anfangswertaufgabe Verwendung findet. Der Schritt den Aringer beschreibt um, wie er es nennt, eine Randwertaufgabe mit freiem Rand auf eine gewöhnliche Randwertaufgabe zurückzuführen, scheint mir unnötig, außerdem bringt das auf diese Weise abgeleitete Differentialgleichungsystem, welches mit dem des AWP's bis auf eine Erweiterung mit der Entfernung identisch ist, keine brauchbaren Ergebnisse.

... Newtonverfahrens^2.14

Die Grundlage dieses Verfahrens folgt aus einer nach dem linearen Glied abgebrochenen Taylorentwicklung nach . Unter der Vorraussetzung, daß ergibt sich die gesuchte Größe aus . Die Gleichung ist dann die Verallgemeinerung für Funktionen von n-Veränderlichen.

... werden^2.15

ARINGER beschreibt in seiner Dissertation die Berechnung des Vertikalschnittbogens über eine Differentiation von () und () mit anschließender Integration (vgl. Gleichungen (6.55) und (6.56) seiner Arbeit). Der Sinn dieser Handlungsweise ist mir nicht klar, da ein konstanter Anteil dadurch wegfällt. Führt man die Integration nach (6.56) formal durch, so ergeben sich nicht die Ergebnisse die er in seinen Tabellen angibt, sondern man erhält eine zu kurze Strecke.

... stehen^2.16

Ausgenommen sind hier die über das GUI im Fenster sichtbar oder versteckt übergebenen Variablenwerte.

... Intervallteilung^2.17

Wenn das Programm für andere Ellipsoiddimensionen (6000km < a > 7000km) verwandt werden soll, so muß hier der Quelltext abgeändert werden, da die Intervallteilung durch Konstanten festgelegt wurde.

... erwiesen^2.18

Um in Beispiel 6 auf nur 12 Schritte zu gelangen, ist für die relative Variation der Strecke anzusetzen.

... Programmsystem^2.19

Ein Flußdiagramm das den Aufbau des Programmsystems erläutert, ist in Abb. dargestellt.

... Linien^2.20

Im folgenden wird der Verlauf der Linie in Analogie zur Satellitenbahn als Ephemeride bezeichnet.

... Stützpunkte^2.21

Die Koordinaten dieser Punkte werden über eine Anfangswertaufgabe ermittelt.

... werden^2.22

Der Vorteil der binären Dateien liegt darin, daß man in einer solchen, auf einfache Weise, verschiedene Linien unter eigenen Bezeichnern abspeichern kann. MATLAB bietet zusätzlich die Möglichkeit diese Datei im ASCII Format auszugeben.

...linienzuordnung_clairaut^2.23

Die Berechnung der geographischen Breiten aus den reduzierten erfolgt nach .