Die Lösung der elastischen Bewegungsgleichung als ... · Kraftstoß in Achsenrichtung eines Stabes...

Die Lösung der elastischenBewegungsgleichung als Modellfall der

Schallemission

vorgelegt vonDiplom – Physiker

Thomas – Rex Heinrich

zurErlangung des akademischen Grades

doctor rerum naturalium (Dr. rer. nat.)

an der Fakultät V für Verkehrs – und Maschinensystemeder

Technischen Universität Berlin

Promotionsausschuss:

Vorsitzender: Prof. Dr. U. von Wagner1. Gutachter: Prof. Dr. W. Müller2. Gutachter: Dr. Chr. Marotzke

Tag der wissenschaftlichen Aussprache: 31. August 2009

Berlin 2009D 83

1

Inhalt

Inhalt.................................................................................................................1

Zusammenfassung..........................................................................................4

1. Einleitung

1.1. Wellen in Physik und Technik………….. …………………………………............6

1.2. Zerstörungsfreie Materialprüfung…………………………………………………..7

1.3. Die Schallemission als Methode der zerstörungsfreien Prüfung…...…………..8

1.4. Probleme der Schallemission……………………………………...……………...11

1.5. Geschichtliches……………………………………………………………………..15

1.6. Das Ziel der Arbeit...........................................................................................17

1.7. Grundlagenliteratur und weiterführende Artikel...............................................18

2. Die Bewegungsgleichung elastischer Wellen auf einem endlichenGebiet mit inhomogenen Randbedingungen

2.1. Die Bewegungsgleichung…………………………………………………….…....23

2.2. Das Anfangs-Randwertproblem……………………………….………….….…...30

3. Numerische Lösung des 3d-Problems

3.1. Warum Numerik?.............................................................................................34

3.2. Normierung der Bewegungsgleichung....………………………………………...35

3.3. Das finite Differenzenverfahren…………………………………………………...37

3.4. Diskretisierung des Gebietes und der Bewegungsgleichung....…………….…41

3.5. Die Randwerte

a) Die Achse……………………………………………………………………46

b) Der Rand ohne Ecken……………………………..……………………….46

c) Die Ecken…………………………………………………………………....47

3.6. Die Anfangswerte…………………………………………………………………...49

3.7. Stabilität……………………….……………………………………………………..50

2

3.8. Numerische Lösung des 3d-Problems (2.2.12)………...……………………….51

3.9. Interpretation der numerischen Rechnungen……………………………………63

3.10. Folgerungen für Schallemissionssignale…………………………………………65

4. Wellenausbreitung auf einem zweidimensionalen Gebiet

4.1. Der Ebene Spannungszustand einer dünnen Scheibe………………..…….....68

4.2. Numerik der vollständigen Gleichungen des EbenenSpannungszustandes…………...……………………………...…………..……...70

a) Der Rand ohne Ecken ………………………………………………..……72

b) Die Ecken…………………………………………………………………....74

4.3. Numerische Lösung der 2d-Bewegungsgleichungen.……………………...…..77

4.4. Interpretation der numerischen Rechnungen………………………………..…..88

4.5. Folgerungen für Schallemissionssignale………………………………….……..90

5. Ein Experiment zur Prüfung der 1d-Reihenlösungen(6.3.33), (6.3.40) für einen zeitbegrenzten Spannungsimpuls………….....93

6. Anhang

6.1. Longitudinal- und Transversalwellen eines Diracimpulses

a) Die Quellenfunktion q………………………………………………….….112

b) Die Greensche Funktion G………………………………………..….….113

c) Aufspaltung des Verschiebungsfeldes ineinen divergenz- und rotationsfreien Anteil…………………………….117

d) Longitudinal- und Transversalwellen……………………………….…...118

6.2. Das Kontaktproblem..………………………………………………....................122

6.3. Analytische Lösung der 1d-Wellengleichung auf einemendlichen Intervall mit inhomogenen Randbedingungen

a) Zeit während des Stoßes……………………………………………..….127

b) Zeit nach dem Stoß……………………………………………………….135

c) Dispersionseigenschaften………………………………………..………136

3

6.4. Lösung des Differentialgleichungsproblems (6.3.13) mit HilfeGreenscher Funktionen………………………………………...............….……137

6.5. Beweis, dass (6.3.33) das Problem (6.3.1) löst ……………………….………139

6.6. Beweis, dass (6.3.40) das Problem (6.3.36) löst ……………………………..140

6.7. 2d-Bewegungsgleichungen für den EbenenSpannungszustand und für Axialsymmetrie……………...…………………….141

6.8. Analytische Lösung desAnfangs-Randwertproblems auf einer Scheibe………………….….…………145

6.9. Beweis, dass (6.8.28) das Problem (6.8.2) löst……………..…………………151

6.10. Fourierkoeffizienten……………………………………………………………….153

6.11. Die d’Alembertlösung des 1d-Problems........................................................154

6.12. Das Spektrum des Kraftimpulses..................................................................158

6.13. 3d-Programmcode.........................................................................................160

6.14. 2d-Programmcode……………………………………………………………...…166

7. Symbolverzeichnis……………………………………………………….……...172

8. Literatur…………………………………………………………………………....173

Danksagung.................................................................................................176

Erklärung......................................................................................................177

Lebenslauf....................................................................................................178

4

Zusammenfassung

Ausgehend von einer kritischen Analyse der Schallemission als einer Methode derzerstörungsfreien Prüfung kommt der Verfasser zu dem Ergebnis, dass die Abbil-dung elastischer Wellen durch Schallemissionssensoren nur dann beurteilt werdenkann, wenn das Eingangssignal durch Analyse der Wellenausbreitung bekannt ist,und das berechnete Ausgangssignal vom gemessenen Signal bestätigt wird.Dazu wird ein lineares Modell aufgestellt, das einen zeitbegrenzten punktförmigenKraftstoß in Achsenrichtung eines Stabes und radial auf den Rand einer Scheibe be-schreibt. Die Ränder sind frei. Gesucht werden das Verschiebungsfeld und das resul-tierende Sensorsignal. Die analytische Lösung des Anfangs-Randwertproblems mitden vollständigen Wellengleichungen gelingt nicht. Es wird darum eine numerischeRechnung durchgeführt, die mit der experimentell durch einen berührungslos mes-senden optischen Sensor ermittelten Verschiebung verglichen wird. Für die Scheibeerhält man folgende Verschiebungsverläufe:

Verschiebung Platte

1,05

1,07

1,09

1,11

1,13

1,15

1,17

1,19

1,21

1,23

1,25

0,003 0,0035 0,004 0,0045 0,005 0,0055

t [s]

U[V

]

Abbildung 1: Errechnete (oben) und gemessene Verschiebung der Scheibe (unten).Die Verschiebung bezieht sich auf den der Erregung gegenüber liegenden Punkt.

Die in Abbildung 1 gezeigten Verschiebungsverläufe stimmen sehr gut überein, undzeigen, dass der numerischen Rechnung vertraut werden darf.

Zweifache zeitliche Ableitung des numerisch ermittelten Signals und Faltung mit ei-ner idealisierten Impulsantwort des Messsystems ergeben das zu erwartende Signal(oben), das dem realen Sensorsignal (unten) gegenübergestellt wird:

5

Sensorsignal

1

1,05

1,1

1,15

1,2

1,25

1,3

1,35

1,4

1,45

1,5

0,01305 0,0131 0,01315 0,0132 0,01325 0,0133 0,01335

t [s]

U[V

]

Abbildung 2: Errechnetes (oben) und gemessenes Sensorsignal (unten) der Scheibe.

Es zeigt sich, dass Sensoren in der ersten Schwingung das theoretisch zu erwarten-de Signal produzieren, oder anders herum, dass ein theoretisch erwartetes Sensor-signal unter kontrollierten experimentellen Bedingungen auch gemessen wird.

Über die Wellenausbreitung wird festgestellt:

1. Das Signal ist bei einem Laufweg s allseitig nach der Zeit / Lt s c zu beobachten.

Es handelt sich darum um das Signal der Longitudinalwelle. In einem zeitlichen Ab-stand entsprechend der Transversalgeschwindigkeit folgt die Oberflächenwelle, wieman in Abbildung 3 sehen kann.

Abbildung 3: Errechnetes Signal der Longitudinal (1)- und Oberflächenwelle (2) aufder Kante einer Scheibe.

2. Das Signal der Oberflächenwelle ist auf derjenigen Kante, auf der die Erregungstattfindet, Hauptträger der Energie. D.h., dass eine kumulative Energiemessung anentsprechendem Ort hauptsächlich die Oberflächenwelle nachweißt. Abbildung 3

1

2

6

1. Einleitung

1.1. Wellen in Physik und Technik

Das immer wiederkehrende Grundthema in Physik und Technik ist die mathemati-sche Beschreibung von Diffusions-, Potential- und Wellenphänomenen. Gründe dafürgibt es genug: Bei dynamischer Belastung von Bauteilen kann der mit der Wellen-ausbreitung einhergehende Energietransport erhebliche, bisweilen zerstörerische,Auswirkungen haben. Die hinlänglich bekannten Bilder von laufenden Wellen indurch Windböen zum Schwingen erregten Brücken sind ein extremes Beispiel dafür.Weniger spektakulär ist die Optimierung von Wärmeleitungsvorgängen in den Ent-wicklungsabteilungen vieler Unternehmen, an die der technische Fortschritt in ande-ren Bereichen oft eng gekoppelt ist. Wellen sind also durch raum-zeitliche Ausbrei-tung von Schwingungszuständen der betrachteten Feldgrößen gekennzeichnet, derohne Massentransport, aber mit Energietransport verbunden ist. Es sei darauf hin-gewiesen, dass zur Ausbreitung von Wellen nicht unbedingt ein stoffliches Mediumexistieren muss, wie die Ausbreitung elektromagnetischer Wellen lehrt. Vielmehr isthier der Raum selbst Träger der Schwingungen. Wellen werden durch hyperbolischepartielle Differentialgleichungen zweiter Ordnung beschrieben. Das bekannteste Bei-spiel dafür ist sicher die klassische d’Alembertsche Wellengleichung:

22

2c

t

uu .

Es sei angemerkt, dass die Gleichung beim Übergang zwischen zueinander beweg-ten Bezugssystemen unter Galilei-Transformation nicht invariant ist, da sie die Wel-lengeschwindigkeit c enthält. Die Form der Naturgesetze darf sich aber beim Wech-sel des Bezugssystems nicht ändern. Das ist der Ansatzpunkt der Einsteinschen Re-lativitätstheorie. In dieser Arbeit werden die durch die Dynamik sich ausbreitenderRisse erzeugten elastischen Wellen in homogenen Körpern untersucht. Grundlagehierfür ist eine verallgemeinerte Wellengleichung, die Lamésche Differentialgleichung(hier mit verschwindenden Volumenkräften):

2

2 2 2

2grad div

uu uT L Tc c c

t.

Ihre Ableitung und ein konkretes Differentialgleichungsproblem, das den Riss durcheinen Spannungsstoß auf den Rand des Körpers ersetzt, werden in Kapitel 2 ange-geben. Die zentrale Gleichung der Quantenmechanik ist eine Wellengleichung ersterOrdnung in der Zeit. Es ist die zeitabhängige Schrödinger-Gleichung, welche die zeit-

liche Entwicklung der Wellenfunktion , t x determiniert:

2

2, ,2

i t V tt m

x x x

.

Andere Beispiele wären die Korteweg-de Vries Gleichung, für welche Solitonen exak-te Lösungen darstellen, die Klein-Gordon Gleichung, die Nichtlineare Schrödinger-Gleichung, die Ginzburg-Landau Gleichung usw. [1]

7

1.2. Zerstörungsfreie Materialprüfung

Jeder professionelle Einsatz von Material ist an dessen Prüfung gebunden, um Risi-ken, die durch Belastungen aller Art entstehen, begegnen zu können.

Zum einen muss man sich der physikalischen und chemischen Eigenschaften desMaterials versichern, um es sinnvoll einsetzen zu können. Es bedarf keiner Erläute-rung, dass mechanische Parameter wie der Elastizitätsmodul usw. bekannt seinmüssen, um für eine mechanische Anwendung entscheiden zu können, welches Ma-terial eingesetzt werden kann. Andererseits muss geprüft werden können, ob wäh-rend eines technologischen Prozesses oder während des Einsatzes eines BauteilsFehler an demselben entstanden sind.

Die Beantwortung der ersten Frage zielt auf die Qualitätssicherung in der Produktion,die Beantwortung der zweiten Frage auf ein Monitoring in sicherheitsrelevanten Fra-gen.

Eine Prüfung kann durch zerstörende oder zerstörungsfreie Methoden vorgenommenwerden. Im ersten Fall ist eine Weiterverwendung des Bauteils für andere Untersu-chungen, etwa Festigkeitsbestimmungen, nicht möglich, im zweiten Fall ist es mög-lich. Es ist klar, dass zerstörende Methoden immer auf Stichprobenkontrolle be-schränkt sein müssen, während zerstörungsfreie Methoden die Grundgesamtheiterfassen können.Natürlich ist es auch möglich, beide Methoden an ausgewählten Bauteilen zu kombi-nieren, was insbesondere bei sensiblen Fällen, etwa kleineren Druckbehältern, auchgetan wird.

Es können nun drei Gruppen von Prüfverfahren für den Fall zerstörungsfreier Prü-fung unterschieden werden [2]:

1. die Prüfung mit Strahlung, insbesondere mit Röntgen- und - Strahlen,

2. die Prüfung mit Ultraschall,3. die magnetischen und magnetinduktiven Prüfverfahren.

Besonders wichtige Beispiele sind für

1. Untersuchungen von Schweißnähten auf Poren, Schlackeeinschlüsse, Ris-se,…

2. Prüfungen von Blechen auf Dopplungen; Gussstücke auf Lunker, Gasblasen,Risse, Poren, Einschlüsse; Schmiedestücke auf Risse

3. Prüfung von Härterissen, Schmiederissen, Schleifrissen, Oberflächenrissen inSchweißnähten; Prüfung von ferromagnetischen und nicht ferromagnetischenmetallischen Werkstücken auf Risse, Abmessungen, Härte,…

Diese Arbeit widmet sich der Prüfung durch Schallemission, wobei, im Unterschiedzum Ultraschallverfahren, elastische Wellen nicht aktiv in das Bauteil eingeschallt,sondern die durch Mikrobrüche im Bauteil erzeugten akustischen Wellen passiv emp-fangen und analysiert werden. Es sind hier Frequenzen von wenigen Hertz bis heraufzu 2 MHz relevant. Auch hier steht die Prüfung auf Rissbildung im Vordergrund, wo-bei, rein theoretisch, den Prüfteilen keinerlei Beschränkungen in Größe und Materialauferlegt sind.

8

1.3. Die Schallemission als Methode der zerstörungsfreien Prüfung

Ausgangspunkt dieser Arbeit ist die Schallemission. Was ist das Prinzip der Mes-sung?

Ein Prüfkörper wird belastet. Auftretende Spannungen verursachen Materialbrüche,die im Anfangsstadium mikroskopisch sind. Im Bruchbereich wird das Material lokalentlastet. Der Spannungsimpuls erzeugt elastische Wellen, die sich allseitig ausbrei-ten und an den Grenzflächen des Körpers vielfältige Reflexionen mit Modenumwand-lungen erfahren. Überlagerung erzeugt ein unübersichtliches Verschiebungsfeld, dasnicht einfach zu berechnen ist. An der Oberfläche des Prüfkörpers sindPiezosensoren angebracht, welche die unter ihnen auftretenden Beschleunigungenin elektrische Spannungen umwandeln. Die Signale werden über Vorverstärker undFilter auf einen computergestützten Messplatz geleitet.

Abbildung 4: Prinzipieller Aufbau einer Schallemissionsmessung.

Insbesondere ist technisch von Interesse, wie die beobachtete elektrische Signal-

spannung U t mit den mechanischen Spannungen am Erregungsort zusammen-

hängen. Gedacht ist dabei an die einfache Situation, dass zeitlich veränderliche äu-ßere Kräfte lokale Randspannungen hervorrufen. Diese Störung des Gleichgewichtsist Ursache elastischer Wellen, welche durch Sensoren detektiert werden können.

,BU tr

0 ,σ r t

0

Abbildung 5: Eine zeitlich veränderliche Randspannung σ t verursacht in einem

piezoelektrischem Sensor eine zeitlich veränderliche elektrische Spannung U t .

VerstärkerFilter

VerstärkerFilter

VerstärkerFilter

9

Sehen wir uns die Modellbildung genauer an. Allgemein bezeichne 1M einkontinuumsmechanisches Modell, das als Abbildung der elektrischen Spannung

,BU tr am festen Beobachtungsort Br auf den mechanischen Spannungstensor

0 ,σ r t am ebenfalls festen Ort der Wellenerregung 0r definiert ist:

10: , , r σ rBM U t t .

Da das elektrische Spannungssignal des Sensors den Beschleunigungen und damitden zweifachen Zeitableitungen der Verschiebung u am Beobachtungsort folgt, geht

1M über in

10: , , u r σ rBM t t .

Das ist die von der Praxis gestellte Grundaufgabe der Schallemission, die aber theo-

retisch als Differentialgleichungsmodell sicher nicht zu bewältigen ist. Da in 1M Ur-sache und Wirkung vertauscht sind, spricht man auch von einem inversen Problem.

Anders sieht es beim umgekehrten, kausalen Fall aus. Wenn von der Ursache aufdie Wirkung geschlossen wird, ist aus der Spannung am Erregungsort die Verschie-bung am Beobachtungsort zu bestimmen:

0: , ,σ r u r BM t t .

Das ist die theoretische Grundaufgabe der Schallemission. Formal muss sie als Dif-ferentialgleichungsproblem formuliert werden. Inhaltlich liegt ein Wellenausbrei-tungsproblem vor, wie sich herausstellen wird.

In einem ersten Schritt ist der zu beschreibende Vorgang in einem Anfangs-Randwertproblem zu formulieren. Dabei ist zu beachten, dass das Problem korrektgestellt wird. [3]

Das bedeutet:

1. Existenz: Es gibt wenigstens eine Lösung , tu r , welche alle Bedingungen

erfüllt.

2. Eindeutigkeit: Es gibt höchstens eine Lösung , tu r , welche alle Bedingun-

gen erfüllt.

3. Stabilität: Die Lösung , tu r hängt stetig von den Anfangs- und Randbedin-

gungen ab. Das bedeutet: Ändert man die Daten geringfügig, so ändert sichdie zugehörige Lösung nur wenig.

4. Zusätzlich muss gefordert werden, dass die Anfangs- und Randbedingungenverträglich sind, dass also die Anfangsbedingungen am Rand gleich denRandbedingungen zum Anfangszeitpunkt sind.

Im Folgenden ist also ein kontinuumsmechanisches Differentialgleichungsproblemals Anfangs-Randwertproblem aufzustellen und zu lösen.

10

Bei der Aufstellung und analytischen Lösung derartiger Probleme werden Vereinfa-chungen vorgenommen, welche das Modell in seinem Gültigkeitsbereich begren-zen. Solange man die Grenzen der Lösungen beachtet, entstehen aber keine Prob-leme.

Die wenigsten Differentialgleichungsprobleme können analytisch gelöst werden.(Auch wenn eine Lösung existiert, wie es der reale Vorgang in der Natur beweist, istes meistens nicht möglich, die mathematischen Schwierigkeiten ohne zu große Ver-einfachungen zu überwinden.)

In diesem Falle bleibt der Weg über numerische Methoden, die nur (oder immerhin)eine zahlenmäßige Lösung des vorliegenden Problems zulassen.Numerische Lösungen haben ihre eigenen Schwierigkeiten. Insbesondere kann mansich der Konvergenz und Stabilität der Lösung von vornherein nicht sicher sein. Essind Maßnahmen notwendig, um diese Eigenschaften zu sichern.

Die Möglichkeiten der Rechentechnik setzen den Rechnungen in der Regel schnelleine Grenze. Hier werden Fortschritte in der Hardwareentwicklung von erweitertenProblemstellungen eingeholt.

Die Numerik ist an möglichst einfachen Fragestellungen zu testen. „Möglichst ein-fach“ heißt, dass das Ergebnis am besten a priori bekannt ist oder gesicherte funkti-onale Zusammenhänge existieren.

In der Praxis geht man dann so vor: Es wird im Modell ein Spannungszustand überdie Randbedingungen vorgegeben und das berechnete Verschiebungsfeld mit denMessungen verglichen. Wenn keine gute Übereinstimmung besteht, werden dieRandbedingungen neu justiert und es erfolgt ein weiterer Versuch. Durch Versuchund Irrtum kann man sich an die gemessenen Werte annähern. Gelingt es, die Be-obachtungen mit dem Modell zu reproduzieren, sind die Ursachen des Messaus-gangs verstanden.

Alle qualitativen und quantitativen Zusammenhänge liefert nun das Modell.

11

1.4. Probleme der Schallemission

Die Schallemission ist ein passives Detektionsverfahren mit einer Vielzahl von Prob-lemen, derer man sich bewusst werden muss. Die folgende Abbildung vermittelt einBild von der komplexen Struktur des detektierten Signals, die einer einfachen Inter-pretation entgegensteht.

reales Schallemissionssignal

-1

-0,8

-0,6

-0,4

-0,2

0

0,2

0,4

0,6

0,8

0,00015 0,00017 0,00019 0,00021 0,00023 0,00025 0,00027 0,00029 0,00031 0,00033 0,00035

t [s]

U[V

]

Abbildung 6: Reales Signal eines Schallemissionssensors im Zeitbereich.

Schon die Anzahl registrierter Verschiebungsimpulse ist unsicher. Zwei schnell auf-einander folgende Impulse können in geringer Entfernung nicht sicher getrennt wer-den, da der Nachhall des Ersten den zweiten überlagert. In der Schallemissionspricht man darum nur von registrierten „Ereignissen“. Dagegen können Ereignissegezählt und ihr Ort durch Lokalisierungsalgorithmen auf der Basis von Laufzeitunter-schieden bestimmt werden.

Die Energiemessung eines detektierten Ereignisses bezieht sich auf die gemesseneelektrische Energie an einem Punkt 0r :

2el

0

1

E U t dtR

.

( R bezeichnet den elektrischen Widerstand des Messkreises, U die Spannung)

Sie stimmt natürlich nicht mit der Energie des räumlich ausgedehnten elastischenWellenfeldes überein. Auch wird ein solcher Energiewert keinem elementaren Span-nungsimpuls zuzuordnen sein, weil sich dieser zeitlich nicht auflösen lässt und nurnoch „Ereignisse“ beobachtet werden. Hinzu kommt, dass der Einfluss der material-spezifischen Dämpfung nicht sicher bekannt ist.

12

Also kann der Messung nur sicher entnommen werden, dass „etwas“ am Ort 1r

stattfand, das relativ zu einem anderen Etwas am Ort 2r ein mehr oder weniger star-

kes elektrisches Signal erzeugte.

Beim Durchlaufen der Messkette ist das Signal an einigen Stellen Umwandlungenunterworfen, die das Signal verändern:

Das erste Glied der Messkette, der Sensor, stellt ein System gekoppelterSchwinger dar, dessen Eigenschwingungsverhalten das Signal stark beein-flusst. Der Sensor wirkt als Bandpassfilter.

Im folgenden Vorverstärker wird das Signal hochpassgefiltert. Vor der Digitalisierung wird das Signal tiefpassgefiltert (Antialiasing-Filter) und

rechentechnisch verarbeitet.

Ein erstes Beispiel widmet sich dem Sensor.

Es wird gezeigt, wie ein Schallemissionssensor auf eine direkt auftreffende Stahlku-gel reagiert. Die Reaktion des Sensors auf den Stoß prägt sich als Impulsantwortdurch Faltung dem Signal auf. Das Nachschwingen des Sensors, also sein „Eigen-leben“, ist gut zu beobachten.

Impulsantwort eines Sensors

-0,5

-0,4

-0,3

-0,2

-0,1

0

0,1

0,0005 0,0007 0,0009 0,0011 0,0013 0,0015 0,0017 0,0019

t [s]

U[V

]

Abbildung 7: Transiente Reaktion eines Schallemissionssensors auf eine direkt auf-prallende Stahlkugel. Messung ohne Vorverstärker und Filter.

Die Abbildung 7 ist ohne Vorverstärker und ohne Filter eines modernen Messplatzes,also oszilloskopisch aufgenommen worden.

Die Auswirkung eines Vorverstärkers auf ein ideales Modellsignal sei an den Abbil-dungen 8 und 9 dargestellt. Dabei wurde ein Generatorsignal anstelle des Sensor-signals auf die Messkette gegeben. Es kann Invertierung und Differenzierung beo-bachtet werden:

13

Rechteckimpuls

-0,5

0

0,5

1

1,5

2

2,5

1,50E-05 1,90E-05 2,30E-05 2,70E-05 3,10E-05 3,50E-05

t [s]

U[V

]

Abbildung 8: Generatorrechteck von 2,0V und einer Breite von 5T µ s .

verformtes Rechteck

-3

-2

-1

0

1

2

3

1,50E-05 1,70E-05 1,90E-05 2,10E-05 2,30E-05 2,50E-05 2,70E-05 2,90E-05 3,10E-05 3,30E-05 3,50E-05

t [s]

U[V

]

Abbildung 9: Rechteck von Abbildung 8 nach Durchgang durch die Messkette.

Die Interpretation der Kurve in Abbildung 9 ist:

1. Das Rechteck wird durch eine entsprechende Verstärkerschaltung invertiert.Man kann es sich durch die Überlagerung von zwei Sprungfunktionen zu-

sammengesetzt denken: 1 2H t t H t t .

14

U

1t 2t t

2. Das invertierte Signal wird differenziert. Entsprechend seiner Ableitung

1 2 1 2 2 1

dH t t H t t t t t t t t

dt

entstehen an den Flanken bei 1 2,t t Diracimpulse. Intervalle konstanter Span-

nung verschwinden.3. Den Diracimpulsen folgen Ausschaltkurven, die exponentiell abklingen.4. Die Elektronik wirkt als RC -Differenzierglied, und ist somit ein Hochpass.

Demzufolge ist in Abbildung 9 zu sehen, was ein Hochpass mit der Grenzfre-quenz 1gf MHz von einem Rechteck übrig lässt.

Derartiges Verhalten verstellt den Charakter der aufgenommenen Signale und er-schwert deren Interpretation sehr. Um das Verhalten des Sensors überhaupt beurtei-len zu können, ist es natürlich wünschenswert, über ein mathematisches Modell desVerschiebungsfeldes zu verfügen. Daraus lässt sich das Beschleunigungssignal ab-leiten, das den Sensor anregt.

Was wären die Vorteile eines solchen Modells?

1. Durch die Berechnung des Verschiebungsfeldes können absolute Werte fürVerschiebung, Beschleunigung usw. für alle Punkte des Körpers und für alleZeiten angegeben werden.

2. Das Beschleunigungssignal am Sensor kann bei Kenntnis des Spannungs-zustandes auf dem Rand vorhergesagt werden. (Lösung der Grundaufgabe)Durch diskrete Faltung dieses Beschleunigungssignals mit der Impulsantwortdes Sensors kann das Sensorsignal konstruiert werden. Damit können die Pa-rameter der Messkette auf das zu erwartende Signal im Vorfeld der Messungabgestimmt werden.

Die Grenzen eines Modells liegen nach dem Gesagten auf der Hand:

1. Die Kenntnis der Anfangs- und Randbedingungen ist oft lückenhaft oder hypo-thetisch.

2. Randbedingungen sind numerisch oft schwer zu behandeln.3. Die Lösung der Bewegungsgleichungen ist mit numerischen Unsicherheiten

behaftet und stößt schnell an die Kapazitätsgrenze der Rechentechnik.

15

1.5. Geschichtliches

In dieser Arbeit werden elastische Wellen behandelt. Die Kontinuumstheorie und dieUntersuchung elastischer Wellen hat eine relativ lange Geschichte. Einer der erstenfundamentalen Beiträge stammt von Robert Hooke (1635-1703) mit der Angabe desZusammenhangs zwischen Spannung und Verzerrung (1678), dem Materialgesetz.Weitere, auch technisch wichtige Resultate verdanken wir Jakob Bernoulli (1654-1705). Er befasste sich mit Studien zur Balkenbiegung und stellte 1694 die Differen-tialgleichung der Biegelinie

3

2 21

w xC M x

w x

,

auf, wobei M das Biegemoment ist. Überhaupt führten im 18. Jahrhundert vorwie-gend Schwingungs- und Biegungsprobleme die Mathematiker auf Differentialglei-chungen höherer Ordnung. So auch auf das Problem der schwingenden Saite, dasfür diese Arbeit eine zentrale Bedeutung hat, und dessen Entwicklung darum etwasausführlicher beleuchtet werden soll. Einen ersten Beitrag hierzu leistete Brook Tay-lor 1713, der fand, dass „Momentaufnahmen“ der Schwingung immer gegeben wür-den durch

sin 0x

y A x LL

.

Eine Zeitabhängigkeit berücksichtigte er noch nicht. Ein großer Fortschritt in der Be-handlung des Problems der schwingenden Seite war durch die Arbeit von Jean Bap-tiste le Rond d’Alembert (1717-1783) zu verzeichnen. Ihm haben wir die Gleichungder schwingenden Saite

2 22

2 2

y yc

t z

(1747) und deren allgemeine Lösung in der Form fortschreitender Wellen

,y z t f ct z g ct z

zu verdanken, die darum auch d’Alembertlösungen genannt werden. Ein Jahr späterpublizierte Leonhard Euler (1707-1783) seine Einsichten in die Bedingungen, denen

eine Lösung ,y z t der Wellengleichung genügen muss. Die noch unbestimmten

Funktionen ,f g werden, neben den Einspannbedingungen, durch die Anfangswerte

festgelegt. Mit den Worten von Bernard Riemann (1826-1866): “Er bemerkte, dassder Natur des Problems nach die Bewegung der Saite vollständig bestimmt sei, wennfür irgendeinen Zeitpunkt die Form der Saite und die Geschwindigkeit jedes Punktes(also y und /y t ) gegeben seien, und zeigte, dass sich, wenn man diese beiden

Functionen sich durch willkürlich gezogene Curven bestimmt denkt, daraus stets

durch eine einfache geometrische Construktion die d’Alembertsche Function f z

finden lässt.“ Mit Daniel Bernoulli (1700-1782) wurde 1753 eine andere Seite des

16

Problems beleuchtet. Er betonte, mehr physikalisch interessiert, die Existenz vonOberschwingungen, welche schon aus der Musik bekannt waren. Als erster weist ermathematisch Eigenschwingungen am Objekt der hängenden Kette (!) nach. DenAufbau komplizierter Bewegungen durch Überlagerung von Eigenschwingungen be-schreibt dann wieder Leonhard Euler. Er gibt in einer Arbeit aus dem Jahre 1748 fürden Fall einer sinusförmigen Anfangsauslenkung

sin 0nn

n xg x c x L

L

(ob die Reihe unendlich ist oder nicht, ist hier noch nicht klar) und verschwindenderAnfangsgeschwindigkeit die Bewegung der Saite

, sin cosnn

n x n ty x t c

L L

als Superposition der Eigenschwingungen sin / cos /nc n x L n t L an. Nun war

die Frage zu beantworten, ob Anfangsauslenkungen der obigen Form allgemein sind,d.h. ob jede beliebige Funktion sich in eine trigonometrische Reihe entwickeln lässt,oder ob es sich um eine spezielle Anfangsauslenkung handelt. Damit hatte das Prob-lem der schwingenden Saite die Frage nach der Entwicklung einer willkürlichenFunktion in eine trigonometrische Reihe hervorgebracht. Dieses Problem ist mit JeanBaptiste Fourier (1768-1830) auf das engste verbunden. Fourier verwandte in seinerersten Arbeit über die Wärme aus dem Jahre 1807 das Konzept, eine willkürlich ge-gebene Funktion durch eine trigonometrische Reihe auszudrücken, was eine neueEpoche in der mathematischen Physik einläutete. Allerdings war schon 1798 Eulerauf die Formeln zur Berechnung der nc gestoßen, was sich in deren Namen („Euler-

Fouriersche Formeln“) widerspiegelt. Das Konzept des Fourieransatzes hat in dertheoretischen Physik eine überragende Bedeutung und wird auch in dieser Arbeitverwendet. ( 6.3, 6.8)) Auf Euler gehen auch die ersten Untersuchungen zumKnickverhalten von elastischen dünnen Stäben zurück. Die Elastizitätstheorie in ihrerheutigen Form begründete Augustin Louis Cauchy (1789-1857). Er führte die Begriffedes Spannungszustandes und den des Verzerrungszustandes ein. In diesem Zu-sammenhang ist die Einführung des nach ihm benannten Spannungstensors zu nen-nen. Simeon Denis Poisson (1781-1840) stellte eine nach ihm benannte Differential-gleichung auf, die es gestattet, Potentiale aus ihren Quellen und Randbedingungenzu bestimmen, G.B. Airy (1801-1892) führte die nach ihm benannte Spannungsfunk-tion der Bipotentialgleichung ein, und so wären noch viele Namen zu nennen. Mei-lensteine waren die Fundamentallösungen der Elastizitätstheorie:

Lösung des Problems einer Punktlast auf der Oberfläche eines elastischenHalbraums durch Joseph Boussinesq 1878.

Lösung des Problems einer Linienlast auf der Oberfläche eines elastischenHalbraums durch Alfred Flamant, einem Studienfreund von Boussinesq.

Lösung des Problems einer Punktlast, die im Inneren eines unendlichen elas-tischen Raumes wirkt, durch William Thompson (Lord Kelvin) 1848.

Lösung des Problems einer tangential auf der Oberfläche eines elastischenHalbraums wirkenden Punktlast durch den Italiener V. Cerrutti 1882.

Lösung des Problems einer (normalen oder tangentialen) Punktlast, die im In-neren eines elastischen Halbraumes wirkt, durch Raymond Mindlin 1936.

17

1.6. Das Ziel der Arbeit

Es soll herausgefunden werden, ob ein piezoelektrischer Sensor die Anregung durcheine elastische Welle korrekt abbildet. Das ist nicht von vornherein klar, denn Schall-emissionssensoren sind zusammengesetzte mechanische Systeme und haben einfür gekoppelte mechanische Oszillatoren typisches kompliziertes Eigenschwingungs-verhalten. Nur durch diese Brille sehen wir auf die objektiven Bewegungen des Kon-tinuums und können darum nicht sicher sein, zutreffende Informationen über dessenDeformationszustand zu erhalten.Ein Kopf-an-Kopf-Betrieb zweier Sensoren kann die Frage leider nicht entscheiden,denn das elektrische Generatorsignal wird bei der Wandlung im Aktuator in unbe-kannter Weise verändert, sodass nicht genau gesagt werden kann, was für eine me-chanische Schwingung als Eingangssignal dient. Eine Beantwortung der Frage kanndaher nur durch die Berechnung des Verschiebungsfeldes eines Probekörpers, d.h.in der Auffindung von Lösungen der elastischen Bewegungsgleichung unter be-stimmten Anfangs- und Randbedingungen, erreicht werden. Wenn Sensoren im Mo-dellexperiment die zweifachen zeitlichen Ableitungen des berechneten Verschie-bungsfeldes reproduzieren, muss davon ausgegangen werden, dass sie Beschleuni-gungen, und damit natürlich auch Verschiebungen, richtig abbilden und es sich imWesentlichen nicht um Eigenschwingungen der Sensoren handelt.Damit wird die Beantwortung der Frage, ob Schallemissionssensoren den Bewe-gungszustand eines Probekörpers richtig abbilden, vom Glauben an die Versiche-rungen der Hersteller unabhängig gemacht und auf eine objektive Ebene gestellt.

Weiterhin werden aus der Dynamik des Verschiebungsfeldes Schallemissionssignalefür verschiedene Messpunkte von Prüfkörpern abgeleitet und verglichen. Es sollauch etwas über die Ausbreitung elastischer Wellen in den Prüfkörpern herausge-funden werden.

Ein methodisches Ziel der Arbeit besteht in der möglichst nachvollziehbaren Entwick-lung mathematischer Grundlagen, ohne den Bereich der komplexen Zahlen zu ver-wenden oder Integraltransformationsmethoden zu nutzen, die oft schwer interpretier-bare Lösungen produzieren. Auch der Übergang zu Potentialen wurde auf ein Mini-mum beschränkt.

Die Arbeit ist in der Form so gestaltet, dass alle Zwischenrechnungen ausführlichangegeben und Beweise nicht unterschlagen werden. Der umfangreiche Anhangenthält Ergänzungen, die benutzte Formeln transparenter machen. Die Ableitungs-schritte sind konsequent durchnummeriert, um das Nachrechnen zu erleichtern. Ein-zelne Abschnitte können darum als „Kochrezepte“ dienen, die einfach auf ein ande-res Problem übertragen werden müssen, um es zu lösen.

18

1.7. Grundlagenliteratur und weiterführende Artikel

Es folgt eine Zusammenstellung der Bücher, die der Verfasser zum Studium von par-tiellen Differentialgleichungen und elastodynamischer Probleme für unverzichtbarhält, gefolgt von einer kurzen Liste von Artikeln, die zu erweiterten Fragestellungenvon Anisotropie, Inhomogenität, Dämpfung, Streuung und Materialbruch hinführen.

Bücher:

Achenbach, J.D.: Wave Propagation In Elastic Solids,Verlag North-Holland bzw. ©Elsevier Science Publishers 1984

Duffy, Dean: Green’s Functions with Applications,Chapman & Hall/CRC 2001

Landau, Lev; Lifschitz, Evgenij: Lehrbuch der Theoretischen Physik,ElastizitätstheorieAkademie - Verlag 1989

Schmutzer, Ernst: Grundlagen der Theoretischen Physik, Teil I,Bibliographisches Institut - Wissenschaftsverlag 1989

Trim, Donald W.: Applied Partial Differential Equatoins,PWS Publishing Company Boston 1989

Tychonoff, A.N., Samarski, A.A.: Differentialgleichungen der mathematischen Physik,VEB Deutscher Verlag der Wissenschaften 1959

Artikel:

Albuquerque EL, Mauriz PW:“Elastic waves propagation in damped media”COMPUTATIONAL SCIENCE AND ITS APPLICATIONS - ICCSA 2003, PT 2,PROCEEDINGS LECTURE NOTES IN COMPUTER SCIENCE 2668: 774-781 2003

Alem L, Chorfi L:“Uniqueness theorem for an elastic waves problem”COMPTES RENDUS MATHEMATIQUE 336 (6): 525-530 MAR 15 2003

Auriault JL:“Body wave propagation in rotating elastic media”MECHANICS RESEARCH COMMUNICATIONS 31 (1): 21-27 JAN-FEB 2004

Bhaskar A:“Waveguide modes in elastic rods”PROCEEDINGS OF THE ROYAL SOCIETY OF LONDON SERIES A-MATHEMATICAL PHYSICAL AND ENGINEERING SCIENCES 459 (2029): 175-194JAN 8 2003

19

Blanchard JP:“Elastic waves induced by surface heating in a half-space”JOURNAL OF APPLIED MECHANICS-TRANSACTIONS OF THE ASME 70 (4): 613-615 JUL 2003

Cao Y, Yang KQ:“Hamiltonian system approach for simulation of acoustic and elastic wave propaga-tion”ACTA PHYSICA SINICA 52 (8): 1984-1992 AUG 2003

Cattani C.:“Solitary waves in dispersive media “NUOVO CIMENTO DELLA SOCIETA ITALIANA DI FISICA B-GENERAL PHYSICSRELATIVITY ASTRONOMY AND MATHEMATICAL PHYSICS AND METHODS 119(7-9): 875-888 JUL-SEP 2004

Cortes-Vega LA:“Resonant frequencies for a system of time-harmonic elastic wave”JOURNAL OF MATHEMATICAL ANALYSIS AND APPLICATIONS 279 (1): 43-55MAR 1 2003

Danilov VN:“Scattering of elastic pulses by a crack”RUSSIAN JOURNAL OF NONDESTRUCTIVE TESTING 39 (1): 48-53 JAN 2003

El-Raheb M, Wagner P:“Transient flexural waves in a disk and square plate from off-center impact”JOURNAL OF THE ACOUSTICAL SOCIETY OF AMERICA 110 (6): 2991-3002 DEC2001

Escalona JL, Mayo JM, Dominguez J:“Influence of reference conditions on the analysis of impact-induced elastic waves”MULTIBODY SYSTEM DYNAMICS 7 (2): 209-228 MAR 2002

Gachter GK, Grote MJ:“Dirichlet-to-Neumann map for three-dimensional elastic waves”WAVE MOTION 37 (3): 293-311 MAR 2003

Grimal Q, Naili S, Watzky A:“A method for calculating the axisymmetric response of a two-layered half-space un-der dynamic concentrated loading”JOURNAL OF SOUND AND VIBRATION 276 (3-5): 755-780 SEP 22 2004

Gsell D, Dual J:“Simulation of elastic wave propagation in cylindrical structures including excitationby piezoelectric transducers”ULTRASONICS 40 (1-8): 181-186 MAY 2002

20

Guy ED, Radzevicius SJ, Conroy JP:“Computer programs for application of equations describing elastic and electromag-netic wave scattering from planar interfaces”COMPUTERS & GEOSCIENCES 29 (5): 569-575 JUN 2003

Guz AN:“Elastic waves in bodies with initial (residual) stresses”INTERNATIONAL APPLIED MECHANICS 38 (1): 23-59 JAN 2002

Hussain W, Ogden RW:“The effect of pre-strain on the reflection and transmission of plane waves at an elas-tic interface”INTERNATIONAL JOURNAL OF ENGINEERING SCIENCE 39 (8): 929-950 MAY2001

Jin B:“Dynamic displacements of an infinite beam on a poroelastic half space due to amoving oscillating load”ARCHIVE OF APPLIED MECHANICS 74 (3-4): 277-287 DEC 2004

Kawashita M, Kawashita W:“Analyticity of the resolvent for elastic waves in a perturbed isotropic half space”,MATHEMATISCHE NACHRICHTEN 278 (10): 1163-1179 2005

Kononov AV, De Borst R:“Radiation emitted by a constant load in a circular motion on an elastically supportedMindlin plate”JOURNAL OF SOUND AND VIBRATION 245 (1): 45-61 AUG 2 2001

Krawczuk M, Zak A, Ostachowicz W, Cartmell M:“Propagation of elastic waves in beams - Including damping effects”MODERN PRACTICE IN STRESS AND VIBRATION ANALYSIS MATERIALSSCIENCE FORUM 440-4: 179-186 2003

Langenberg KJ, Marklein R:“Transient elastic waves applied to nondestructive testing of transversely isotropiclossless materials: a coordinate-free approach”WAVE MOTION 41 (3): 247-261 MAR 2005

Lavrov NA, Pavlovskaya EE:“Transient behavior of a few dies on an elastic half-space”ACTA MECHANICA 144 (3-4): 185-195 2000

Lavrov AV:“Fracture localization in space as a cause of variations in the amplitude distribution ofacoustic emission signals”ACOUSTICAL PHYSICS 51 (3): 321-327 MAY-JUN 2005

Lee DS:“Scattering of torsional waves by a spherical cavity in a long circular elastic cylinder”ACTA MECHANICA 164 (1-2): 47-59 2003

21

Le Rousseau JH, de Hoop MV:“Generalized-screen approximation and algorithm for the scattering of elastic waves”QUARTERLY JOURNAL OF MECHANICS AND APPLIED MATHEMATICS 56: 1-33Part 1 FEB 2003

Liu SW, Huang JH:“Transient dynamic responses of a cracked solid subjected to in-plane loadings”INTERNATIONAL JOURNAL OF SOLIDS AND STRUCTURES 40 (18): 4925-4940SEP 2003

Morales A, Mendez-Sanchez RA, Flores J:“Topological defects in 1D elastic waves”PHYSICA E-LOW-DIMENSIONAL SYSTEMS & NANOSTRUCTURES 19 (3): 289-291 AUG 2003

Pan Y, Rossignol C, Audoin B:“Identification of laser generated acoustic waves in the two-dimensional transient re-sponse of cylinders”JOURNAL OF THE ACOUSTICAL SOCIETY OF AMERICA 117 (6): 3600-3608 JUN2005

Pan Y, Rossignol C, Audoin B:“Acoustic waves generated by a laser line pulse in a transversely isotropic cylinder”APPLIED PHYSICS LETTERS 82 (24): 4379-4381 JUN 16 2003

Pierce AD:“Vibrations of spherical inclusions in elastic solids”ACOUSTICAL PHYSICS 51 (1): 5-10 JAN-FEB 2005

Sun WM, Yang GS, Li DX:“Exact analysis of wave propagation in an infinite rectangular beam”APPLIED MATHEMATICS AND MECHANICS-ENGLISH EDITION 25 (7): 768-778JUL 2004

Tyutekin VV:“Simulation and synthesis of elastic wave absorbers at the boundary of a rigid body:Incidence of longitudinal waves”ACOUSTICAL PHYSICS 47 (4): 461-467 JUL-AUG 2001

Widehammar S:“Estimation of 3D field quantities and energy flux associated with elastic waves in abar”JOURNAL OF SOUND AND VIBRATION 259 (4): 893-915 JAN 23 2003

Widehammar S, Gradin PA, Lundberg B:“Approximate determination of dispersion relations and displacement fields asso-ciated with elastic waves in bars - Method based on matrix formulation of Hamilton'sprinciple”JOURNAL OF SOUND AND VIBRATION 246 (5): 853-876 OCT 4 2001

22

Wu CY, Chang JS, Wu KC:“Analysis of wave propagation in infinite piezoelectric plates”JOURNAL OF MECHANICS 21 (2): 103-108 JUN 2005

Yang LY, Turner JA:“Scattering of elastic waves in damaged media”JOURNAL OF THE ACOUSTICAL SOCIETY OF AMERICA 113 (6): 2992-3000 JUN2003

Ziv M:“Source signature and elastic waves in a half-space under a momentary shear lineimpulse”INTERNATIONAL JOURNAL FOR NUMERICAL AND ANALYTICAL METHODS INGEOMECHANICS 27 (3): 233-258 MAR 2003

Ziv M:“Source signature and elastic waves in a half-space under a sustainable line-concentrated impulsive normal force”INTERNATIONAL JOURNAL FOR NUMERICAL AND ANALYTICAL METHODS INGEOMECHANICS 26 (4): 373-406 APR 10 2002

23

2. Die Bewegungsgleichung elastischer Wellen auf einemendlichen Gebiet mit inhomogenen Randbedingungen

2.1. Die Bewegungsgleichung

In diesem Kapitel wird die Bewegungsgleichung für das Verschiebungsfeld abgeleitetund das Anfangs - Randwertproblem formuliert, wobei die Darstellung ErnstSchmutzer [4] und Landau/Lifschitz [5] folgt.

Man stelle sich Folgendes vor: Ein beliebiger Bewegungsvorgang des Kontinuums

versetzte ein Massenelement, das sich zur Zeit 0t an der Stelle P mit dem Ortsvek-

tor x r e befinde, zu einer Stelle P mit dem Ortsvektor , ,t t r r r u r bzw.

, ,x t x u t r r . (2.1.1)

Dabei ist der Verschiebungsvektor , tu r :

,u t u e r r r . (2.1.2)

, tu r Momentankonfiguration

Referenzkonfiguration

r , t r r u r

0

Abbildung 10: Der Vektor r bezeichnet ein bestimmtes Massenelement. Dessen

Verbleib, die Deformation, beschreibt der Verschiebungsvektor , tu r .

Bilde nun das Differential der Verschiebungen:

udu dx H dx

x

. (2.1.3)

Die Größe H nennt man den Verschiebungsgradient. Aus der Verschiebung

können die materiellen Koordinaten x eliminiert werden:

, ,u u x x t t . (2.1.4)

Dann gilt unter der Voraussetzung

P

P

24

1u

x

, (2.1.5)

die geometrische Linearisierung heißt [6], für den Verschiebungsgradienten:

u u x u x uH

x x x x x x

u u

x x

,

denn das Kroneckerdelta ist 1 für , und 0 sonst. Es ist also

u uH

x x

, (2.1.6)

was bedeutet, dass zwischen den Ableitungen nach x und nach x nicht mehr un-

terschieden werden muss. Stelle nun den Verschiebungsgradient in der folgendenWeise dar:

1 1

2 2

D

u u uu uH

x x x x x

. (2.1.7)

Es ist

1

2

u uD

x x

(2.1.8)

der antisymmetrische Rotationstensor, und

1

2

u u

x x

(2.1.9)

der symmetrische Deformationstensor. Dessen Spur hat die Form

1Tr div

2

u u u

x x x

u . (2.1.10)

Den Zusammenhang von Spannung und Dehnung liefert das Hooksche Gesetz,auch Materialgesetz genannt, für den linearen elastischen Körper:

el C . (2.1.11)

Der Tensor C ist der Elastizitätstensor, der als absolutes Glied aus der Taylor-

25

Entwicklung des allgemeinen Hookschen Materialtensors hervorgeht. Der Elastizi-tätstensor für isotrope elastische Festkörper ist gegeben durch

C . (2.1.12)

2, /N m sind die Laméschen Konstanten. Einsetzen in (2.1.11) ergibt

el

el

el

Da ist und für die Spur (2.1.10) gilt, kann man

2el u

x

schreiben. Der Deformationstensor wird nun mit Hilfe von (2.1.9) eliminiert:

el u u u

x x x

. (2.1.13)

Betrachtet man nun die Dynamik eines N - Massenpunktsystems, so gilt in einembeschleunigungsfreien Bezugssystem die Bewegungsgleichung

1

Ne

k kk

m

r F . (2.1.14)

Die N Massepunkte mögen sich in einem kleinen Volumen V befinden. eF isteine äußere Kraft, die auf die Masse im Inneren des Volumens wirkt. Sie lässt sich

als Integral über die Volumenkraftdichte , ,t f t f r e r schreiben:

e

V

dV

F f . (2.1.15)

Eine weitere Kraft, die Flächenkraft AF , resultiert aus der Wechselwirkung des um-gebenden Kontinuums mit dem Volumenelement und greift an dessen Oberfläche

an. Sie ergibt sich als Integral einer Flächenkraftdichte , tP r , auch Spannungsvek-

tor (Spannung) genannt, über die Oberfläche A des Volumenelements V :

A

V

dA

F P . (2.1.16)

Die Flächenkraftdichte kann immer in Normal- und Tangentialkomponente (Normal-und Tangentialspannung) zerlegt werden:

26

n tP P P n t . (2.1.17)

Um ein Kräftegleichgewicht zu beschreiben, müssen wir (2.1.14) um die von derUmgebung des Volumenelements herrührende Flächenkraft ergänzen:

1

Ne A

k kk

m

r F F . (2.1.18)

Es wird nun eine lokale Massendichte 0

, lim /V

t m V

r eingeführt, wobei der

Grenzübergang so zu verstehen ist, dass jedes noch so kleine Volumen V immernoch so groß ist, dass es genügend viele Teilchen enthält. Dann kann die Gesamt-masse des Volumenelementes V auch integral ausgedrückt werden:

1

,N

kk V

m t dV

r .

Diese Masse des betrachteten Volumenelementes wird mit dessen Beschleunigungr multipliziert, womit (2.1.18) übergeht in

e A

V

dV

r F F . (2.1.19)

Benutzt man die Kraftdarstellungen (2.1.15), (2.1.16), kann auch die rechte Seiteintegral formuliert werden:

V V V

dV dV dA

r f P . (2.1.20)

Subtraktion vonV

dV f ergibt

V V

dV dA

r f P , (2.1.21)

oder in Komponenten

V V

x f dV P dA

. (2.1.22)

Die rechte Seite muss nun auch in ein Volumenintegral umgeformt werden. Es wurdefestgestellt, dass P im allgemeinen nicht die Richtung des Flächennormalenvektorsn besitzt. Es soll P mit dem Spannungstensor durch n ausgedrückt werden:

P n . (2.1.23)

Folglich ist mit dA n dA

27

V V V

P dA n dA dA

. (2.1.24)

Allgemein kann ein Oberflächenintegral mit dem Gaußschen Satz

divV V V V

ad dV a dA dV

x

a A a

in ein Volumenintegral umgewandelt werden. Daraus folgt für (2.1.24)

V V V

P dA dA dVx

. (2.1.25)

Ersetzen der rechten Seite in (2.1.22) ergibt

V V

x f dV dVx

.

Durch Subtraktion der rechten Seite entsteht ein Integral, das für jede beliebige Stel-le dV des Volumens V gilt, und darum identisch verschwindet:

0V

x f dVx

.

Daraus kann auf das Verschwinden des Integranden geschlossen werden, und esresultiert die lokale Form der Bewegungsgleichung der Kontinuumsmechanik

x fx

. (2.1.26)

Wegen v x kann das als

vv

v f

t x(2.1.27)

geschrieben werden. Da für Festkörper das in der Geschwindigkeit quadratischekonvektive Glied vernachlässigbar ist [4], gilt

v

ft x

. (2.1.28)

Aus (2.1.1) und der Zeitunabhängigkeit von r folgt

v , r t x 2 2

2 2, , ,

r r r

u uu t t t

t t,

28

da bei kleinen Verschiebungsbeträgen näherungsweise nicht mehr zwischen r undr unterschieden werden muss. [19] Also folgt für linke Seite von (2.1.28):

2

2

uf

t x

. (2.1.29)

Des Weiteren wird eine konstante Massendichte vorausgesetzt:

0, r t .

Eliminiere mit (2.1.13) die Spannungen:

2 2 2 2

0 2 2

u u uu u u uf f

t x x x x x x x x x

.

Da nur für von Null verschieden ist, folgt

2 2 2 2

0 2 2

u u u uf

t x x x x x

.

Umbenennung der Indizes in der Form lässt die Gleichung

2 2

0 2

u uf u

t x x

(2.1.30)

entstehen. In Komponenten:

2 2 2 231 1 1 1 1 2

0 12 2 2 21 2 3 1 1 2 3

2 2 2 232 2 2 2 1 2

0 22 2 2 21 2 3 2 1 2 3

2 2 2 23 3 3

0 32 2 21 2

uu u u u u uf

t x x x x x x x

uu u u u u uf

t x x x x x x x

u u uf

t x x

3 31 223 3 1 2 3

u uu u

x x x x x

(2.1.31)

Dabei müssen die Volumenkraftdichten f vorgegeben werden. Die Gleichungen

können in Vektorschreibweise ausgedrückt werden:

2

0 1 2 3 1 2 3 1 2 3 1 2 321 2 3 1 2 3

grad div

, , , , , , , , , , , ,u u u f f f u u u u u ut x x x x x x

u f u

u

Division durch 0 ergibt die endgültige Form

29

2

20 0 0

grad divt

u fu u . (2.1.32)

Die Lösung ist nur unter Angabe von Anfangs- und Randbedingungen eindeutig:

,0 und ,0

Rand, oder Rand,el

t

t P n t

uu r r

u

(2.1.33)

Die Bewegungsgleichung (2.1.32) lässt sich noch in einer anderen Form schreiben.Dazu benötigt man die Identität grad div rot rot u u u .

Beweis:

Voraussetzung: a b c b a c b a c

Behauptung: grad div rot rot u u u

Beweis: In rot rot u u ist dann , , a b c u . Also folgt

rot rot grad div u c c u u , womit grad div rot rot u u u bewiesen

ist. ■

Mit den Quadraten

2 2

0 0

2L Tc c

(2.1.34)

der longitudinalen und transversalen Phasengeschwindigkeit ( Helmholtzsche Zer-legung (6.1.24)) kann (2.1.32) in der modifizierten Form

22 2

20

grad div rot rotL Tc ct

u fu u . (2.1.35)

geschrieben werden [4].

30

2.2. Das Anfangs-Randwertproblem

Das gesamte Gebiet eines Stabes ist in Zylinderkoordinaten wie folgt definiert:

: , , 0 0 2 0G r z r R z L . (2.2.1)

( bedeutet „und“, R ist der Radius des Stabes, L dessen Länge)

I

L r

z

R

F t

Abbildung 11: Die Bezeichnung der Koordinaten und ein rz -Schnitt.

Der Stoß erfolgt mit der zeitabhängigen Kraft F t allein im Punkt 0z der Symmet-

rieachse I und wird damit idealisiert als punktförmig betrachtet. Daraus folgt, dassfür das Problem Axialsymmetrie vorausgesetzt werden muss:

0 0u

. (2.2.2)

Äußere Volumenkraftdichten sollen nicht vorhanden sein:

, 0t f r . (2.2.3)

Aus (2.1.35) folgt unter den genannten Voraussetzungen:

31

2

2 2

2grad div rot rotr r z z L r r z z T r r z zu u c u u c u u

t

e e e e e e (2.2.4)

Das ist natürlich wegen der Linearität der Differentialoperatoren

2 22 2 2 2

2 2grad div grad div rot rot rot rotr z

r z L r r L z z T r r T z z

u uc u c u c u c u

t t

e e e e e e (2.2.5)

Bilde nun die Terme der rechten Seite unter der Voraussetzung (2.2.2):

2 2

2 2

1 1grad div

1 1 1 1 1

grad div

rr r r z r r z r

r r r r r rr r z r r r z

z z r z

uu ru u

r z r r r z r r

u u u u u uu u u

r r r z r r r r r r r z r z

ur

e e e e e

e e e e

e e e2 2

2

2 2 2

2 2

2 2

1 1rot rot rot

1 1rot rot rot

z z zr z

r r r r r rr r r r r z z r

z z z zz z z z r z r z

u u u

z z r z z

u u u u u uu u r

z z r r z r z r z z

u u u uu u r

r r z r r r r z r

e e

e e e e e e e

e e e e e e e2

2z zu u

r r

Einsetzen in die Bewegungsgleichung (2.1.35) ergibt

2 2 2 22

2 2 2 2

2 2 2 2 2 22 2 2

2 2 2

1 1 1

1 1

r z r r r rr z L r r z

z z r r r z z zL r z T z r T r z

u u u u u uc u

t t r r r r r z r z

u u u u u u u uc c c

r z z r z r z z r z r r r

e e e e

e e e e e e

Ordnen nach Einheitsvektoren:

2 2 2 2 2 22 2

2 2 2 2 2

2 2 2 22 2

2 2

1 1

1 1 1

r z r r z r zr z r L r T

r r z z r r zz L T

u u u u u u uc u c

t t r r r r r z z r z

u u u u u u uc c

r z r z z r r r z r z r

e e e

e

(2.2.6)

Vergleich linke Seite-rechte Seite ergibt die endgültigen Bewegungsgleichungen:

32

2 2 2 2 22 2

2 2 2 2

2 2 2 2 22 2

2 2 2

1 1

1 1 1

r r r z r zL r T

z z r r z z r rL T

u u u u u uc u c

t r r r r r z z r z

u u u u u u u uc c

t z r z r z r r r r z r z

(2.2.7)

Des Weiteren wird ein rz -Schnitt betrachtet, der in Abbildung 12 dargestellt ist.

r

rK

zK *G zK

z0 I

Abbildung 12: rz -Schnitt mit den Kanten und der Symmetrieachse.

Das gesamte rechteckige Gebiet G , die Kanten K und die Symmetrieachse I sind:

: , 0 0

: , 0 0

: , 0

: , 0

: , 0 0

z

r

z

G r z r R z L

K r z r R z

K r z r R z L K

K r z r R z L

I r z r z L

(2.2.8)

Der Rand ist dann die Vereinigung

z r zK K K K . (2.2.9)

Durch den Aufprall einer Kugel wirke auf den Stab die zeitlich begrenzte Kraft

0 0sin sin 0

0

tF F t t T

F t T

t T

(2.2.10)

( 0F ist die Amplitude, T die Stoßzeit, die zugehörige Kreisfrequenz.)

Es werde der zeitliche Druckverlauf

02

sin sin

z

Ff t H T t t H T t p t

a(2.2.11)

33

über der Kontaktzone angenommen, wobei a der Kontaktradius des Kontaktprob-lems Kugel-Stab ist. (Spektrum Abb. 64, 6.12.) zz ist in Flächennormalenrichtung

gerichtet, nach Konvention positiv, sowie der äußeren Kraft entgegengerichtet:

I

f t zz

Abbildung 13: Druck und Spannung im Angriffspunkt der äußeren Kraft für t T .

Das Anfangs-Randwertproblem ist mit (2.2.7):

2 2 2 2 22 2

2 2 2 2

2 2 2 2 22 2

2 2 2

1 1

1 1 1

RB: : 0,

: 0, 0

:

r r r z r zL r T

z z r r z z r rL T

z rz zz

z rz zz

r rz

u u u u u uc u c

t r r r r r z z r z

u u u u u u u uc c


K r f t

K

K

0, 0

AB: , ,0 0 , ,0 0

rr

r z r zt

uu

(2.2.12)

Die Lösung für das Verschiebungsfeld ist bei Angabe von Rand- und Anfangsbedin-gungen eindeutig. (Beweis in [7]) Da die Spannungen an den Rändern für konkreteRechnungen in den Verschiebungen formuliert werden müssen, sind die folgendenBeziehungen in Zylinderkoordinaten wichtig [8]:

1 12

1 12

1 12

1 1

1

r zrr r

r zr

z rzz r

rr r

r zrz zr

z z

uu uu

r r r z

u u uu

r r r z

uu uu

z r r r

uuu

r r r

u u

z r

u

z r

zu

(2.2.13)

34

3. Numerische Lösung des 3d-Problems

3.1. Warum Numerik?

Um sich einen Überblick über den dreidimensionalen Verschiebungsverlauf zu ver-schaffen, müssen die vollständigen Gleichungen (2.2.12) gelöst werden, was abermit analytischen Mitteln aus zwei Gründen nicht gelingt.

Das Ausgangsproblem kann zum einen nicht durch eine Transformation auf eineForm gebracht werden, bei der homogene Randbedingungen und eine inhomogeneDifferentialgleichung (PDGL mit Störterm) vorliegen.

Der oft gesuchte Ausweg, einen auf dem Rand wirkenden Spannungsstoß direkt inden Kraftdichteterm für die äußeren Volumenkräfte hineinzuschreiben [9], ist nichtzulässig. [10]

Auch wenn eine Transformation des Problems, analog der eindimensionalen Lösungdes Anfangs-Randwertproblems im Anhang, gelänge, hat man es mit gekoppeltenGleichungen zu tun, die so nicht lösbar sind.

Es kann jedoch versucht werden, eine zahlenmäßige Lösung zu erhalten, wie es nungeschehen soll. Die Störung wird hierbei ganz zwanglos als Randbedingung berück-sichtigt.

Das weitere Vorgehen sei in einigen Punkten dargelegt:

Eingangs wird die Bewegungsgleichung dimensionslos gemacht. Dann sind die Bewegungsgleichung, das Gebiet und die Zeit zu diskretisieren. Es muss ein Algorithmus angegeben werden, nach dem die Lösung in der Zeit

entwickelt wird. Die Randbedingungen einschließlich der Störung sind umzusetzen. Die Anfangswerte sind zu berücksichtigen. Der numerische Algorithmus ist an der eindimensionalen Theorie zu kontrollie-

ren. Es ist sicherzustellen, dass der Algorithmus stabil arbeitet.

Bei der Umsetzung dieser Punkte wird sich herausstellen, dass die umfangreichenDifferenzensterne auf dem Rand und der Achse besondere Schwierigkeiten bereiten,was jedoch durch geeignete Handhabung der Randbedingungen aufgelöst werdenkann.

35

3.2. Normierung der Bewegungsgleichung

Die dynamischen Gleichungen (2.2.12) sind:

2 2 2 2 2

2 2 2 22 2

2 2 2 2 2

2 2 2

1 1

1 1 1

r r r z r zr

L T

z z r r z z r r

u u u u u uu

t r r r r r z z r zc c

u u u u u u u u


Des Weiteren werden alle Längen auf den Radius R als charakteristische Länge desSystems, und die Zeit auf die Ausbreitungszeit / LR c normiert:

Lcr zr z t t

R R R . (3.2.1)

Damit gehen die dimensionsbehafteten Verschiebungen über in

, , , , , , , ,r r z zu r z t R u r z t u r z t R u r z t . (3.2.2)

Die Differentialoperatoren sind nun

2 2

2 2 2

2 2

2 2 2

22 2

2 2 2

1 1

1 1

L L

r

r r r R r r R r

z

z z z R z z R z

c ct

t t t R t t R t

(3.2.3)

Einsetzen in die Bewegungsgleichungen ergibt

2Lc R 2

2 2

2

r

L

u

R t

c R

2

2

2 2

L

z

R

cu

R t

2

2 2ru R

R r

2

1 ru R

R r r

2 2

1r

Ru

R r

2

2zu

R r z

R

2

2 2zu R

R z

2

1 ru R

R r z

2

2

2

T

r

R

cu

R r z

2

2 2ru R

R z

2

2zu

R r z

R

2

2 2zu R

R r

2

1 zu R

R r r

2

1 ru R

R r z

2

2ru

R r z

Nach Multiplikation mit 2 2/ LR c und Weglassen der Striche erhält man folgende di-

mensionslose Differentialgleichungen:

2 2 2 2 2

22 2 2 2

22 2 2 2 2

2 2 2

1 1

1 1 1

r r r z r zr

T

Lz z r r z z r r

u u u u u uu

ct r r r r r z z r z

cu u u u u u u u


(3.2.4)

36

Die Koeffizienten 2 1Lc und 2 2 2/T L Tc c c übernehmen die Funktion von dimensionslo-

sen Wellengeschwindigkeiten im normierten Gitter. Diese Gleichungen bilden denAusgangspunkt der numerischen Rechnungen.

Die Randspannungen sind invariant gegenüber der Normierung, da die Normie-rungsgröße R in Nenner und Zähler steht:

2zz

R

Rzu R

z

R

1r

Ru

r

Rru

r

(3.2.5)

Da in der Störung im Argument der Sinusfunktion sich der transformierende Faktorkürzt, und die Sprungstelle der Sprungfunktion sich durch Multiplikation mit einer po-sitiven Zahl nicht verändert, behält die Störung ihre Form:

sinL

z

L

R

cRf t H T t p

c

L

t

R

c

sinz

tH T t p

TT

.

Auch hier werden die Striche wieder weggelassen:

sinz

tf t H T t p

T

. (3.2.6)

37

3.3. Das finite Differenzenverfahren

Es soll das Finite Differenzenverfahren zur Anwendung gebracht werden, das eini-ge Vorteile hat. Differenzenverfahren sind sehr nahe an der Theorie und gestattendie unmittelbare Umsetzung des Formalismus. Zum anderen ist das Verfahren leichtverständlich und es kann gut auf erweiterte Probleme ausgeweitet werden. [11] Au-ßerdem lässt sich die Methode leicht auf einem Rechner umsetzen. Mit dem Verfah-ren können sowohl gewöhnliche als auch partielle Differentialgleichungen gelöstwerden. Leider lässt sich das Verfahren nur sinnvoll auf rechteckige Gebiete anwen-den, was im vorliegenden Problem der Fall ist. Des Weiteren soll das explizite Diffe-renzenverfahren benutzt werden. [12]

Die Verfahrensweise sei im eindimensionalen Fall an folgendem Schema erläutert.

Stelle die Bewegungsgleichung und die Rand - und Anfangsbedingungen auf.

2 2

2

2 2RB: 0, Störung , 0

u u u uc t t L t

t z z z

Ersetze in der Bewegungsgleichung und den Randbedingungen alle Differential-operatoren durch Differenzenquotienten und diskretisiere das Gebiet.

2

2 2

, 1 2 , , 1 1, 2 , 1,

1, 1, 2 Störung 1, 1,

u j t u j t u j t u j t u j t u j tc

t h

u t u t h t u n t u n t

Stelle die diskretisierte Bewegungsgleichung nach dem Zeitschritt 1t um.

2

2

2, 1 1, 2 , 1, 2 , , 1

tu j t c u j t u j t u j t u j t u j t

h

Die Anfangsverschiebungen legen das Feld zum Zeitpunkt 0 fest. Die Anfangsge-schwindigkeiten legen das Feld durch einen einfachen Rückwärtsschritt in der Zeitzum Zeitpunkt 1 fest.

,0 , 1

,0 ,0 , 1u j u j

u j f j v j g j u j f j t g jt

Berechnung des Zeitschritts 1, indem die diskretisierte Bewegungsgleichung für uexplizit ausgerechnet wird. Außerhalb des Gebiets liegende Werte werden durch dieRandbedingungen aus den Differenzensternen eliminiert. So wird auch die Störungberücksichtigt.

38

2

2

2,1 1,0 2 ,0 1,0 2 ,0 , 1

tu j c u j u j u j u j u j

h

Berechnung des Zeitschritts 2 aus den Zeitschritten 0,1 usw. ....

Die konkrete Rechnung stellt sich so dar:

,r zu u sind für die Zeitschritte 1,t t bekannt.

Setze für 1t die Störung am Rand

Berechne für 1t die Verschiebung aus den Zeitschritten 1,t t

Ist das zeitliche Ende erreicht?

Ja Nein: 1t t

Ende

Grundlage des Diskretisierungsverfahrens ist die Approximation der Lösungsfunktionan den Stützstellen eines Gitters. [13] Als Approximation der Lösung wird eine Git-

terfunktion u eingeführt, welche zunächst nur an den Gitterpunkten ,i tx t definiert

ist. Zwischen den Gitterpunkten kann interpoliert werden. Das Gitter sei in einer Di-mension der Länge 1 durch folgendes Vorgehen definiert:

Definition 1

Wähle eine natürliche Zahl 1n . Wähle in x - Richtung die Schrittweite 1/x n .

Die x - Koordinaten der Gitterpunkte lauten 0,1,2, ,ix i x i n .

Der Zeitschritt 0t werde vorgegeben. Aus ihm errechnen sich die Zeitstu-fen pt t t .

Um Näherungsausdrücke für die Differentialoperatoren der Differentialgleichung undder Randbedingungen zu erhalten, betrachtet man die Taylorentwicklung einerviermal stetig differenzierbaren Funktion.

2 3 43 4

1

2 3 43 4

2

2 6 24

2 6 24

x x xu x x u x x u x u x u x u x x

x x xu x x u x x u x u x u x u x x

(3.3.1)

39

wobei 10 x x und 20 x x ist. Addition und Subtraktion der Gleichungen

(3.3.1) ergibt

44 42

1 2

3 43 4 4

1 2

224

23 24

xu x x u x x u x x u x u x x u x x

x xu x x u x x x u x u x u x x u x x

(3.3.2)

Nach Umstellung erhält man die Approximation der 1. und 2. Ableitung:

2 33 4 4

1 2

24 4

1 22

2 6 48

2

24

h

h

E x

E x

u x x u x x x xu x u x u x x u x x

x

u x x u x u x x xu x u x x u x x

x

(3.3.3)

Im Fall der zweiten Ableitung genügt der Fehlerterm hE x der Ungleichung

2

4sup

12

g

h gx

M xE x M u x

. (3.3.4)

Für jede vorgegebene Funktion u konvergiert der Fehlerterm gegen Null, wenn x

gegen Null geht. Bei Polynomen 3. Grades ist 4 0u , was bedeutet, dass der Feh-ler für alle x verschwindet. Die zentralen Differenzenquotienten sind also

2

2 2

2

2

du x u x x u x x d u x u x x u x u x x

dx x dx x

(3.3.5)

Da bei partieller Ableitung eine Variable konstant bleibt, können die Approximationen(3.3.5) für den Fall mehrerer Veränderlicher übernommen werden.

2

2 2

, , , , , 2 , ,

2

u x y u x x y u x x y u x y u x x y u x y u x x y

x x x x

(3.3.6)

Wendet man den Differenzenquotient der 1. Ableitung auf den einer anderen Variab-len an, erhält man die Approximation der gemischten Ableitung [14]:

2

2

2

, , , ,, 2 2

2

, , , ,,

4

, , , , ,

4

u x x y y u x x y y u x x y y u x x y yu x y x x

y x y

u x x y y u x x y y u x x y y u x x y yu x y

y x x y

u x y u x x y y u x x y y u x x y y u x x y y

y x x y

(3.3.7)

40

Die partiellen Ableitungen bis zur 2. Ordnung werden auf einem äquidistanten Gitter

x y h

durch die folgenden zentralen Differenzenquotienten approximiert, wobei die durch

den Abbruch entstehende Ungenauigkeit immer 2O h bzw. 2O t ist:

2

2 2

2

2 2

2

2

, , , , , ,

2

, , , , 2 , , , ,

, , , , 2 , , , ,

, , , , , , , , , ,

4

u x y t u x h y t u x h y t

x h

u x y t u x y t t u x y t u x y t t

t t

u x y t u x h y t u x y t u x h y t

x h

u x y t u x h y h t u x h y h t u x h y h t u x h y h t

x y h

(3.3.8)

41

3.4. Diskretisierung des Gebietes und der Bewegungsgleichung

Ein äquidistantes Raum-Zeitgitter wird nach Definition 1 aufgebaut, wobei t auseinem Stabilitätskriterium folgt, das weiter unten ( 3.7.) begründet wird. Die nor-mierten Maße des Stabes und die dimensionslose Gesamtzeit sind als gestricheneKonstanten gegeben mit

1 L LE E

c cR LR L T T T T

R R R R . (3.4.1)

Der Raum- und Zeitschritt ist mit 1/16 ( 3.7.), max 1c (maximale Wellenge-

schwindigkeit im normierten Gitter) und Nn :

1 1 1

16 2 16 2L

R h hr z h t

n n c

. (3.4.2)

Die Striche werden nun weggelassen. Die Gitterpunkte sind mit , , Ni j t

0 1

0 1

0 1

: 0 1 1

: 0 1

: 0 1

i

j

t E

r i h r h r h n h

z j h z h z h l h L

t t t t t t t q t T

(3.4.3)

Definition 2

Des Weiteren seien die Punkte des Raum-Zeitgitters durch die ihnen zugeordnetennatürlichen Zahlen bezeichnet. Der Zeitschritt t wird zur Abkürzung weggelassen,

und stattdessen , ,, , : ,r z r zu i j t u i j geschrieben.

Ersetzen der Differentialoperatoren in (3.2.4) durch (3.3.8) ergibt für die erste Glei-chung:

22

2

2

22

2 2

, , 1 2 , , , , 1

1, 2 , 1, 1, 1,1 1,

2

1, 1 1, 1 1, 1 1, 1

4

, 1 2 , , 1 1, 1

r r r

r r r r rr

z z z z

r r r z zT

L

u i j t u i j t u i j t

u i j u i j u i j u i j u i ju i j

h i h h i ht

u i j u i j u i j u i j

h

u i j u i j u i j u i j u ict

c h

2

1, 1 1, 1 1, 1

4z zj u i j u i j

h

(3.4.4)

(Man beachte, dass die ganze rechte Seite zum Zeitschritt t gehört.)

Und ebenso wird mit der anderen Gleichung von (3.2.4) verfahren:

42

22

2

2 22

2

, , 1 2 , , , , 1

, 1 2 , , 1 , 1 , 11

2

1, 1 1, 1 1, 1 1, 1

4

1, 2 , 1, 1, 1,1 1

2

z z z

z z z r r

r r r r

z z z z z

T

L


u i j u i j u i j u i j u i j

h i h ht


h


c h i h h it

c

2

, 1 , 1

2

1, 1 1, 1 1, 1 1, 1

4

r r

r r r r

u i j u i j

h h


h

(3.4.5)

Mit den Abkürzungen

22 2

2 2 24 4T

L

ct ta b

h h c

erhält man die folgenden Differenzengleichungen:

2

, , 1 2 , , , , 1

2 44 1, 2 , 1, 1, 1, ,

1, 1 1, 1 1, 1 1, 1

4 , 1 2 , , 1 1, 1 1, 1 1, 1 1, 1

,

r r r

r r r r r r

z z z z

r r r z z z z

z


u i j u i j u i j u i j u i j u i ji ia


b u i j u i j u i j u i j u i j u i j u i j

u i

, 1 2 , , , , 1

24 , 1 2 , , 1 , 1 , 1

1, 1 1, 1 1, 1 1, 1

24 1, 2 , 1, 1, 1,

2, 1 , 1 1, 1 1, 1

z z

z z z r r

r r r r

z z z z z

r r r r

j t u i j t u i j t

u i j u i j u i j u i j u i jia


u i j u i j u i j u i j u i ji

b

u i j u i j u i j u i ji

1, 1 1, 1r ru i j u i j

Weiterhin wird in den 1/ i -Termen zu i eine kleine Zahl mit 1h addiert, unddamit die Achse in eine - Umgebung eingebettet:

1/ 1/i i . (3.4.6)

Mit welchem Ausdruck die Einbettung realisiert wird, ist willkürlich. Wichtig ist allein,dass sie folgende Eigenschaften hat:

00

lim1/ 1/ 1/ 1/

i

i i i .

43

Das etwas größer als h sein muss hängt mit dem Abbruchfehler zusammen, der

durch die Behandlung der Ecken von der Ordnung O h ist, wie man 3.5.c) ent-

nehmen kann.

Die für ganz *G zu lösenden Differenzengleichungen sind:

2

, , 1 2 , , , , 1

2 44 1, 2 , 1, 1, 1, ,

1, 1 1, 1 1, 1 1, 1

4 , 1 2 , , 1 1, 1 1, 1 1, 1 1, 1

r r r

r r r r r r

z z z z

r r r z z z z


u i j u i j u i j u i j u i j u i ji ia



, , 1 2 , , , , 1

24 , 1 2 , , 1 , 1 , 1

1, 1 1, 1 1, 1 1, 1

24 1, 2 , 1, 1, 1,

2, 1 , 1 1, 1

z z z

z z z r r

r r r r

z z z z z

r r r


u i j u i j u i j u i j u i jia


u i j u i j u i j u i j u i ji

b

u i j u i j u i ji

1, 1 1, 1 1, 1r r ru i j u i j u i j

(3.4.7)

Die rechten Seiten der Gleichungen können auch anschaulicher in Sternform darge-stellt werden:

za b u 4 rbu za b u

12 2 ra u

i

2

82

r

aa b u

i

12 2 ra u

i

za b u 4 rbu za b u

Abbildung 14: Erste Differenzengleichung (3.4.7) in Sternform.

44

ra b u 2

4 z rau a b ui

ra b u

12 2 zb u

i

8 za b u

12 2 zb u

i

ra b u 2

4 z rau a b ui

ra b u

Abbildung 15: Zweite Differenzengleichung (3.4.7) in Sternform.

Man stelle sich das Gitter nun folgendermaßen durchnummeriert vor:

j

0,2 1,2 2,2

0,1 1,1 2,1

i

0,0 1,0 2,0

Abbildung 16: Durchnummerierung des Gitters.

Nun wird das Gleichungssystem der Differenzengleichungen (3.4.7) in Vektor-Matrixschreibweise formuliert.

Dazu wird das Variablenfeld auf eine lineare Struktur, den Vektor u , abgebildet. Dassind gleichzeitig die Spalten der Matrix M :

: 0,0 , 1,0 , , , , 0,0 , 1,0 , , ,r r r z z zu u u n l u u u n lu . (3.4.8)

Jede Zeile der Matrix M entspricht nun einem Gitterpunkt. Dies entspricht der Aus-wertung des Sterns aus Abbildung 14 auf dem Gitter.Dann werden der Matrix M noch einmal soviel Zeilen angefügt. Dies entspricht derAuswertung des Sterns aus Abbildung 15 auf dem Gitter.

45

Es entsteht eine quadratische Matrix mit 2

2 1 1n l Elementen:

0,0 , 0,0 ,

0,0

: ,

0,0

,

r r z zu u n l u u n l

n l

n l

M

Abbildung 17: Schema der Matrix.

Die Randstörung tritt in Form von Zahlen auf, die nicht als Koeffizienten von Ver-schiebungen fungieren, und durch die Elimination hereinkommen. Sie bilden die In-homogenität I .

Mit dem Vektor (3.4.8), der Matrix M und der Inhomogenität I können die Gleichun-gen (3.4.7) kompakt geschrieben werden:

21 2 1t t t t t u u u I M u . (3.4.9)

ru

zu

46

3.5. Die Randwerte

Auf dem gesamten Rand K können die Sterne nicht unmittelbar berechnet werden,

da dort Verschiebungen außerhalb des Gebietes benötigt werden. Das Problem wirdgelöst, indem aus den Rand- und Symmetriebedingungen Beziehungen, „Eliminati-onsbeziehungen“, gewonnen werden, mit deren Hilfe die überhängenden Verschie-bungen aus den Differenzensternen eliminiert werden können. [14]

a) Die Achse 0i :

Die sich aus den Symmetriebedingungen ergebenden Werte 0, 0ru z sollen nicht

direkt verwendet, sondern numerisch errechnet werden. Dabei werden zur Eliminati-on die Symmetriebedingungen

1, 1, 1, 1,r r z zu j u j u j u j (3.5.1)

genutzt. Diese bedeuten ja gerade, dass die Funktionen bezüglich der Symmetrie-achse gerade fortgesetzt werden können. [14] Der Parameter steuert die Güte der

Näherung 0

lim , 0ru z

.

b) Der Rand ohne Ecken:

Auf dem Rand existieren keine Normalspannungen, außer im Punkt 0,0 :

1

2 z rr zz

u uu f t f t

z r r

,

was unter Beachtung der Ableitungsdifferenzen (3.3.8)

22

, 1 ,1 ,0 1,0 1,02 2

z z r r r

h f tu i u i u i u i u i

i

(3.5.2)

ergibt. Für alle anderen Randpunkte verschwindet die Inhomogenität. Für ru muss

die Tatsache der Schubspannungsfreiheit verwendet werden:

0 0r zrz

u u

z r

. (3.5.3)

Auch hier kann nach Diskretisierung nach den überhängenden Verschiebungen um-gestellt werden:

, 1 ,1 1,0 1,0r r z zu i u i u i u i . (3.5.4)

Entsprechend wird auf den anderen Rändern vorgegangen.

47

c) Die Ecken:

Da in den Differenzengleichungen gemischte Ableitungen enthalten sind, erhält manan den Ecken Gitterpunkte, die in r - und z -Richtung gleichzeitig über den Rand hin-ausgehen, d.h. diagonal sind. Für deren Elimination stelle man nach dem Schema(3.3.7) mit Hilfe einfacher Vorwärts- und Rückwärtsschritte Differenzenapproximatio-

nen für die gemischten Ableitungen auf. Für die Ecke , 0n wird für die Ableitung in

r -Richtung natürlich eine Rückwärtsdifferenz, und für die Ableitung in z -Richtungeine Vorwärtsdifferenz verwendet:

2

2

2

2

2

,0

, , , ,, ,

, , , , ,

,1 1,1 ,0 1,0

n

u r z h u r h z h u r z u r h zu r z u r z h hr z z r h

u r z u r z h u r h z h u r z u r h z

r z h

u n u n u n u nu

r z h

(3.5.5)

So auch für ,n l , nur dass jetzt für die Ableitung in z -Richtung ebenfalls eine Rück-

wärtsdifferenz genommen wird:

2

2

2

2

2

,

, , , ,,

, , , , ,

, 1, , 1 1, 1

n l

u r z u r h z u r z h u r h z hu r z h hr z h

u r z u r z u r h z u r z h u r h z h

r z h

u n l u n l u n l u n lu

r z h

(3.5.6)

Diese Approximationen werden auf den Ecken in den Differenzengleichungen (3.4.7)verwendet.

Es werden nun Eliminationsbeziehungen für die in einer Richtung überhängendenGitterpunkte benötigt. In z -Richtung nimmt man die in zentralen Differenzen formu-lierte Randbedingung des zK - und zK -Randes. Dies produziert zusätzlich Punkte,

die in r -Richtung überhängen. Diese werden mit einseitig formulierten Differenzendes rK -Randes eliminiert.

Die Richtigkeit dieses Vorgehens bestätigt ein Test, in dem das Verschiebungsfeldeines Spannungsstoßes über den gesamten Querschnitt bei fehlender Querdehnungberechnet wird. Da der Verschiebungsverlauf a priori klar ist, besteht hier einestrenge Kontrollmöglichkeit für die Numerik.

48

Tabelle 1: Eliminationsbeziehungen des Randes.

Mit den Eliminationsregeln der Tabelle 1, sowie mit den modifizierten Differenzen fürdie gemischten Ableitungen kann das vollständige Gebiet berechnet werden.

zK

0

22, 1 ,1 ,0 1,0 1,0

2 2z z r r r

i

h f tu i u i u i u i u i

i

, 1 ,1 1,0 1,0r r z zu i u i u i u i

rK 1, 1, , 1 , 1z z r ru n j u n j u n j u n j

2

1, 1, , , 1 , 12

r r r z zu n j u n j u n j u n j u n jn

zK 2

, 1 , 1 , 1, 1,2

z z r r ru i l u i l u i l u i l u i li

, 1 , 1 1, 1,r r z zu i l u i l u i l u i l

,0n 2

, 1 ,1 ,0 1,0 1,02

z z r r ru n u n u n u n u nn

, 1 ,1 1,0 1,0r r z zu n u n u n u n

1,0 ,0 ,1 ,0z z r ru n u n u n u n

1

1,0 ,0 ,0 ,1 ,02

r r r z zu n u n u n u n u nn

,n l 2

, 1 , 1 , 1, 1,2

z z r r ru n l u n l u n l u n l u n ln

, 1 , 1 1, 1,r r z zu n l u n l u n l u n l

1, , , , 1z z r ru n l u n l u n l u n l

1

1, , , , , 12

r r r z zu n l u n l u n l u n l u n ln

Achse 1, 1,r ru j u j

Achse 1, 1,z zu j u j

49

3.6. Die Anfangswerte

Da das vorliegende Problem ein Anfangs-Randwertproblem ist, muss die Lösungnicht nur die Differentialgleichungen und die Randwerte, sondern auch einen vorge-gebenen Anfangszustand erfüllen.

Und da es sich um Differentialgleichungen zweiter Ordnung handelt, sind zwei Zeit-schritte vorzugeben.

Das Verschiebungsfeld des Zeitschritts 0t ist auf *G vorgegeben:

, ,0 , ,0

, ,0 , ,0

r

z

u i j e i j

u i j f i j

(3.6.1)

Das Geschwindigkeitsfeld des Zeitschritts 0t ist ebenfalls vorgegeben, und seidurch einen Rückwärtsschritt in der Zeit approximiert:

, ,0 , , 1, ,0

, ,0 , , 1, ,0

r r

z z

u i j u i jg i j

t

u i j u i jh i j

t

(3.6.2)

Das Verschiebungsfeld des Zeitschritts 1t ist dann mit (3.6.1):

, , 1 , ,0 , ,0

, , 1 , ,0 , ,0

r

z

u i j e i j t g i j

u i j f i j t h i j

(3.6.3)

Da in dem zu lösenden Problem (2.2.12) die Anfangsverschiebungen- und Ge-schwindigkeiten verschwinden, also

, , , 0e f g h

sind, gilt dies nach (3.6.3) auch für die Komponenten des Verschiebungsvektors desZeitschritts 1t :

, , 1 0

, , 1 0

r

z

u i j

u i j

(3.6.4)

50

3.7. Stabilität

Die numerische Lösung darf kein chaotisches oder oszillatorisches Verhalten zeigen.Um dies zu gewährleisten, muss der Zeitschritt t zum Raumschritt h in einembestimmten Verhältnis stehen. Sie sind über die (maximale) Wellengeschwindigkeitgekoppelt. Ihren quantitativen Zusammenhang gibt die Courant-Friedrichs-Lewy-Bedingung an. Sie ist eine notwendige Bedingung für die Stabilität des Verfahrens.Arbeitet ein Algorithmus unter Einhaltung dieser Bedingung stabil, spricht man vonoptimaler Stabilität. [15]

Theorem (Courant-Friedrichs-Lewy-Bedingung, CFL-Bedingung):

Der Abhängigkeitsbereich der kontinuierlichen Gleichungen muss im Abhängigkeits-bereich der diskreten Gleichungen enthalten sein. So gilt bei Verwendung einesäquidistanten Gitters für ein nichtdissipatives homogenes isotropes Medium

max

c t h

d, (3.7.1)

wobei eine Konstante mit 0 1 ist, d der Dimension des Gitters entspricht und

maxc die maximale Wellengeschwindigkeit im normierten Gitter ist. Bei einem zweidi-

mensionalen Gitter wird also 2d gewählt. [16],[17]Anschaulich gesprochen bedeutet dies, dass sich die Welle im Gitter während einesZeitschritts nicht weiter ausbreiten darf, als der -te Teil der halben 2d-Gitterzellendiagonale oder als der -te Teil eines Drittels der 3d-Gitterzellendiagonale. Wird dies verletzt, kommt es zu Instabilität.

Weiterhin ist eine obere Schranke des Zeitschritts t durch das Abtasttheorem vonKotelnikov (1933) und Shannon (1948) gegeben [16]. Dieses besagt, dass die Ab-tastfrequenz sf größer sein muss, als das Doppelte der höchsten im Signal enthal-

tenden Frequenz f (Nyquist-Frequenz), um ein zeitkontinuierliches Signal eindeutig

aus seinen Abtastwerten rekonstruieren zu können:

2 sf f . (3.7.2)

Die Impulsdauer für die Numerik beim Stab ist 0,2T µs . Aus der Spektralformel

(6.12.3) ist zu entnehmen, dass dann mit einer Nyquist-Frequenz von 7 MHz zu

rechnen ist. Es folgt aus dem Abtasttheorem eine Ungleichung für den Zeitschritt t :

61 12 7 10 0,0714

t s

t s. (3.7.3)

Die maximale Wellengeschwindigkeit im Gitter ist auf 1 normiert, und alle Zeiten sindmit dem Normierungsfaktor /Lc R dimensionslos gemacht. Mit 8R mm und 1/ 2

ist dann durch (3.7.1) eine obere (dimensionslose) Schranke von 0,15 h gegeben.

In der Rechnung ist aber ein um den Faktor 1/10 besseres h realisiert: 0,02 h .

51

3.8. Numerische Lösung des 3d-Problems (2.2.12)

Das Programm ist in Mathematica-Code der Version 4.0 geschrieben, der auf einemC-Kern basiert. Der vollständige Programmcode ist im Anhang 6.13. angefügt, umdie Rechnungen bei Interesse nachvollziehen zu können.

Es wird auf einem Gitter von 49 61 Punkten gerechnet, da die Rechenzeiten beieiner größeren Anzahl von Punkten zu lang werden. Bei einem Radius von 8mm und

einer Länge von 10 mm wird ein h von 0,166 mm in beiden Richtungen erreicht, was

Reserven für die Erfassung eines schnellen Prozesses mit seinen hohen Frequenzenlässt, sowie eine Darstellung in guter Auflösung ermöglicht. Die Dauer der Störung istauf 0, 2T µs gesetzt, um ringförmige Wellen sichtbar zu machen.

Input:

Die Geometrie des Körpers wird durch die Konstanten

,L R : Länge des Stabes, Radius des Stabes charakterisiert.

Das Material ist durch die Parameter

0, , : Lamésche Konstanten, konstante Dichte festgelegt.

Der Stoß wird berücksichtigt durch

0, 2 ,T µs f t , die Stoßzeit und die Störungsfunktion.

Output:

In den folgenden Graphiken sind Dichte- und Schnittplots dargestellt.

Die Graphik stellt in der obersten Reihe links die Verteilung der zu -

Komponente des Verschiebungsfeldes auf einem rz -Schnitt räumlich undrechts flächenhaft dar.

In der Reihe darunter ist die ru -Komponente in der gleichen Weise und für

den gleichen Zeitschritt dargestellt.

Auf die Dichteplots folgen die rz -Schnitte. Es handelt sich dabei für beideVerschiebungskomponenten um Schnitte in z -Richtung bei 0, 4,8r mm und

in r -Richtung bei 0,5,10z mm .

Als Erstes wird ein Testfall gerechnet, welcher der Belastung eines Stabes über dengesamten Querschnitt, einem Hammerschlag, entspricht. Das Ergebnis wird durchRechnung und Messung bestätigt, die im Anhang 6.3. und Kapitel 5 zu finden sind.

52

Abbildung 18: Testrechnung, die den Fall der Belastung eines Stabes über den ge-samten Querschnitt darstellt. Reflexion bei 10z mm am freien Ende, 0ru .

53

Abbildung 19: Dichte- und Schnittplots für die Zeit von 0,07 µs.

54


55


56


57


58


59


60


61


62


63

3.9. Interpretation der numerischen Rechnungen

Das Maximum der ru -Komponente wird während des Spannungsstoßes in

45°-Richtung in den Halbraum abgestrahlt. Nachdem sich das ru -Maximum von der Zylindergrundfläche abgelöst hat,

wechselt auf dieser noch während des Spannungsstoßes die ru -Verschiebung

hinter der Wellenfront das Vorzeichen von Plus zu Minus. Das Maximum der zu -Komponente breitet sich während des Spannungssto-

ßes in Form eines Ellipsoids um den Koordinatenursprung aus, dessen großeHalbachse mit der z -Achse zusammenfällt.

Man kann in dem folgenden Dichteplot entsprechend 0 1,35t µs sehen, dass

das Einflussgebiet der Verschiebung einem Halbkreis um den Koordinatenur-sprung mit dem Radius 0Lc t entspricht. D.h. die Wellenfront breitet sich allsei-

tig mit Longitudinalgeschwindigkeit Lc in Form einer Kugelwelle aus.

Abbildung 29: Lage der Wellenfront für 0 1,35t µs ist 0Lc t .

Der Verschiebungsvektor muss daher an der Wellenfront immer senkrecht aufdieser stehen:

L z zuu e

L r r z zu u u e e

ze

re L r ruu e

Abbildung 30: Longitudinale Komponente an der Wellenfront.

Die Rechnung bestätigt, dass auf der Achse die ru -Komponente verschwin-

det, und die zu -Komponente sich mit Lc fortpflanzt. Ebenso verschwindet auf

der Grundfläche des Zylinders an der Wellenfront die zu -Komponente:

64

Abbildung 31: ru hat die Ecke 8, 0r z nach 1, 2 µs als Longitudinalwelle (1) er-

reicht, zu beeinflusst diese aber erst nach 1, 4 µs (2), und ist bis zum Einsatz der

Oberflächenwelle (3) unbedeutend.

Longitudinal- und Transversalwelle des freien Volumens überlagern sich ent-sprechend den Randbedingungen eines Halbraums zu einer Oberflächenwel-le, auch Rayleighwelle genannt. Sie übernimmt die Rolle der Transversalwelleauf der Grundfläche des Zylinders. Die Massenelemente bewegen sich hier inder Zeit auf elliptischen Bahnen. Im isotropen Material ist die Geschwindigkeitdieser Oberflächenwelle ungefähr die der Transversalwelle: O Tc c .

Die Wellenfront der zu -Komponente erfährt am der Erregung gegenüber lie-

genden Punkt 0, 10r z eine heftige Reflexion ohne Phasensprung wie an

einem freien Ende, und wandert als Kugelwelle in die umgekehrte Richtung.

Die Wellenfront der ru -Komponente erfährt bei 0, 10r z natürlich keine

Reflexion, denn sie verschwindet ja auf der Achse, sondern sie wandert als re-flektierte Kugelwelle, von der Zylinderdeckfläche geführt, zur oberen Ecke,und interferiert dort mit jenem Teil der Welle, der den Mantel heraufgewandertist.

Die Wellenzüge laufen nicht auseinander, sondern bleiben klar erkennbar zu-sammen. Das sie also keine Dispersion erfahren ist eine Folge der Vorausset-zung: Die zugrunde liegende Theorie ist linear.

1

32

65

3.10. Folgerungen für Schallemissionssignale

Des Weiteren soll die Frage untersucht werden, welche Zeitsignale ein Schallemissi-onssensor an bestimmten Punkten des zylindrischen Prüfkörpers registrieren würde.

Das weitere Vorgehen ist klar: Da die Signale von Piezosensoren durch Verformung,d.h. Krafteinwirkung, hervorgerufen werden, müssen die auftretenden dynamischenBeschleunigungen berechnet werden. An den interessierenden Punkten wird alsodas numerisch errechnete Verschiebungsfeld zweimal nach der Zeit differenziert. Esresultiert das Beschleunigungssignal, je nach differenzierter Verschiebungskompo-nente ra oder za .

Hierbei ist anzumerken, dass eine numerische Differentiation einen vorhandenenFehler verdoppelt. Demzufolge vervierfacht zweifache numerische Differentiation denvorhandenen numerischen Fehler.

Das Sensorsignal wird durch diskrete Faltung der Beschleunigung mit der idealisier-ten Impulsantwort des Messsystems gewonnen. (Zur Technik der diskreten Faltungsiehe [18].)Die Faltung bedeutet eine Glättung, was die starke Betonung der Fehler durch dienumerische zweifache Zeitableitung verringert.

Wie die idealisierte Impulsantwort aussehen muss, kann Abbildung 7 entnommenwerden. Der Verlauf wird gut durch eine Gaußfunktion approximiert. Vor dem Einset-zen der Erregung sei die Impulsantwort Null.

Abbildung 32: Angenommene Impulsantwort des Messsystems, deren Standardab-

weichung 71,5 10 s entspricht. Vor Einsetzen der Erregung ist sie Null.

Die Beschleunigungen und Sensorsignale sind an den Punkten 0, 0r z und

0, 10r z aus der zu -Verschiebung, an den Punkten 8, 0r z und

8, 5r z aus der ru -Verschiebung gebildet:

66

Abbildung 33: Verschiebung (blau) und das Sensorsignal (rot) für die Punkte

0, 0r z , 0, 10r z , 8, 0r z und 8, 5r z .

Entsprechend den bisher getroffenen Aussagen ist folgendes über zu erwartendeSchallemissionssignale zu sagen:

Das Longitudinalsignal wird bei einem Laufweg s allseitig nach der Zeit

/ Lt s c zu beobachten sein.

Die Signalform des Longitudinalsignals (erste vollständige Schwingung) ist füralle Messpunkte ähnlich, und folgt dem Signal im Krafteinleitungspunkt. (Abbildung 33)

Das Signal am Rand der Zylindergrundfläche bei 8, 0r z ist dem Signal-

verlauf am Krafteinleitungspunkt am unähnlichsten, und nimmt eine Sonder-stellung ein.

Der Einsatz der Transversalwelle geht mit signifikanter Amplitudenvergröße-

rung einher. Beispiel: Das in der Achse des Zylinders, am Punkt 0, 10r z ,

berechnete Signal zeigt nach Durchgang des Verschiebungsringes, der sichmit Transversalgeschwindigkeit ausbreitet , starke Amplitudenzunahme:

67

Abbildung 34: Signal der zu -Verschiebung am Punkt 0, 10r z . Einsatz der

Longitudinalwelle und der Transversalwelle ist mit 1 und 2 gekennzeichnet.

Das theoretisch zu erwatende Signal kann nicht einem realen Signal gegenüber ge-stellt werden, da Messungen an Volumen nicht vorliegen. Ein Vergleich errechnetes-gemessenes Signal findet sich im Abschnitt 4.5, in dem Messungen an einer Scheibemit numerischen Rechnungen an einem Scheibenmodell verglichen werden.

1

2

1

2

68

4. Wellenausbreitung auf einem zweidimensionalen Gebiet

4.1. Der Ebene Spannungszustand einer dünnen Scheibe

Die allgemeine Bewegungsgleichung ohne Volumenkräfte lautet

0 ,u . (4.1.1)

Ausgeschrieben:

0 1 1 , 11,1 12,2 13,3

0 2 2 , 21,1 22,2 23,3

0 3 3 , 31,1 32,2 33,3

u

u

u

(4.1.2)

Im Ebenen Spannungszustand (ESZ) ist 13 23 33 0 , d.h. alle Spannungen

mit z -Komponente verschwinden. [19] Weiterhin gilt für den Spannungstensor

13 31 23 32 :

0 1 11,1 12,2 13,3u

0 2 21,1 22,2 23,3u

0 3 31,1u 32,2 33,3

(4.1.3)

Also folgt

0 1 11,1 12,2

0 1,1 2,2 ,

0 2 21,1 22,2

3

, 1, 2

0

uu

u

u

(4.1.4)

Die Elimination der Spannungen aus (4.1.4) erfolgt mit

2 . (4.1.5)

Nun ist aber 11 22 33 . Das Hooksche Gesetz liefert für den ESZ

33 11 22 33 33 33 11 220 22

. (4.1.6)

Einsetzen in (4.1.5) ergibt

11 22 11 22 11 22 ,

22 2 1 2

2 2 2u

(4.1.7)

Um in (4.1.4) die Spannungen durch die Verschiebungen auszudrücken, müssen die

69

Elemente des Spannungstensors und die Dehnungen nach x abgeleitet werden:

/

, , , ,

/

, , , , ,

2 22 2

2 2

1 1

2 2

x

x

u u

u u u u

(4.1.8)

Eingesetzt in (4.1.4) erhält man

0 , , , , , , , ,

2 1 2 22 2

2 2 2 2u u u u u u u u

Ausgeschrieben:

0 1 1,11 1,11 1,22 2,12 1,11 2,21

0 2 2,11 1,21 2,22 2,22 1,12 2,22

2

2

2

2

u u u u u u u

u u u u u u u

oder

22 2 2

2 2 20 0 0

2 2 2 2

2 2 20 0 0

4 2 3

2 2

4 2 3

2 2

yx x x

y y y x

uu u u

t x y x y

u u u u

t y x x y

(4.1.9)

Es werden weiterhin die Abkürzungen

2 2 2 2

20 0 0 0 0 0

4 1 1 2 3 1

2 2 1 2 2 11L T L T

E E Ec c c c

(4.1.10)

verwendet. Damit lauten die Bewegungsgleichungen (4.1.9)

22 2 22 2 2 2

2 2 2

2 2 2 22 2 2 2

2 2 2

yx x xL T L T

y y y xL T L T

uu u uc c c c

t x y x y

u u u uc c c c

t y x x y

(4.1.11)

Die Gleichung der zu -Komponente wurde weggelassen, da sie entkoppelt ist.

Eine analytische Lösung der Gleichungen (4.1.11) unter der Annahme von Radial-symmetrie findet man im Anhang 6.8..

70

4.2. Numerik der vollständigen Gleichungen des Ebenen Spanungs-zustandes

In diesem Abschnitt sollen die Gleichungen des zweidimensionalen Ebenen Span-nungszustandes (4.1.11)

2 22 2 22 2

2 2 2

2 2 22 22 2

2 2 2

y yx x xL T

y y yx xL T

u uu u uc c

t x x y y x y

u u uu uc c

t y x y x x y

numerisch ausgewertet werden. Des Weiteren werden alle Längen auf die Breite Bals charakteristische Länge des Systems, und die Zeit auf die Ausbreitungszeit

/ LT B c normiert:

Lcx yx y t t

B B B (4.2.1)

Damit gehen die dimensionsbehafteten Verschiebungen über in

, , , , , , , ,x x y yu x y t B u x y t u x y t B u x y t .

Die Differentialoperatoren sind nun

22 2 2 2 2 2 2 2

2 2 2 2 2 2 2 2 2 2 2 2 2 2

1 1 1 Lc

x B x y B y x y B x y t B t

Einsetzen in die Bewegungsgleichungen ergibt nach Fortlassen der Striche

2 22 2 22 2 2

2 2 2 2 2 2

2 2 22 22 2 2

2 2 2 2 2 2

y yx x xL L T

y y yx xL L T

u uu u uc c c

B t B x x y B y x y

u u uu uc c c

B t B y x y B x x y

Nach Multiplikation mit 2 2/ LB c erhält man folgende dimensionslose Differentialglei-

chungen:

2 22 22

2 222

22 2 22 2

2 2 2

y yx xx

T

Ly y yx x

u uu uu

x x y y x yct

cu u uu u

t y x y x x y

(4.2.2)

71

Diese Gleichungen können sofort mit den Differenzenapproximationen (3.3.8)

diskretisiert werden. Mit 2 2 2 2 2 2/ 4 , / 4 /T La t h b t h c c erhält man die folgen-

den Differenzengleichungen:

, , 1 2 , , , 1

4 1, 2 , 1, 1, 1 1, 1 1, 1 1, 1

4 , 1 2 , , 1 1, 1 1, 1 1, 1 1, 1

, , 1 2 , , , 1

4 , 1 2 ,

x x x

x x x y y y y

x x x y y y y

y y y

y y

u i j t u i j u i j t

a u i j u i j u i j u i j u i j u i j u i j


u i j t u i j u i j t

a u i j u i j

, 1 1, 1 1, 1 1, 1 1, 1

4 1, 2 , 1, 1, 1 1, 1 1, 1 1, 1

y x x x x

y y y x x x x



(4.2.3)

Diese Gleichungen können anschaulicher in Sternform dargestellt werden:

ya b u 4 xbu ya b u

4 xau 8 xa b u 4 xau

ya b u 4 xbu ya b u

Abbildung 35: Erste Differenzengleichung (4.2.3) in Sternform.

72

xa b u 4 yau xa b u

4 ybu 8 ya b u 4 ybu

xa b u 4 yau xa b u

Abbildung 36: Zweite Differenzengleichung (4.2.3) in Sternform.

Die Scheibe der Breite B und der Länge L liege in einem kartesischen Koordinaten-system. Die Ränder werden wieder durch zwei Indizes bezeichnet, von denen dasErste die Richtung der Normale bezeichnet, und das Zweite den Richtungssinn der

Normale: , , ,x x y yK K K K . Beispiel: x

0 yK B

a) Der Rand ohne Ecken

Zur Aufstellung von Eliminationsregeln werden der Spannungstensor und der Verzer-rungstensor für den ESZ aus (4.1.8) herangezogen:

, , ,

2 12 , 1, 2

2 2u u u

.

Beispielsweise muss am Rand yK in y -Richtung eliminiert werden. Betrachte dazu

die Normalspannung yy mit , , / 2 /yy y y y y yu u u y :

, ,

2 22 2

2 2

y yxyy yy x x y y

u uuu u

y x y

. (4.2.4)

Am unteren Rand yK greift bei / 2x B ein Spannungsstoss des zeitlichen Verlaufs

sinyf t H T t p t (4.2.5)

73

mit Dauer T in positiver y -Richtung an. yy ist in Flächennormalenrichtung gerich-

tet, also positiv:

2

2 4

y xyy

u uf t f t

y x

Ersetze die Differentialoperatoren durch Differenzen nach (3.3.8):

,1 , 1 1,0 1,0 2

2 2 2 4

y y x xu i u i u i u i

f th h

.

Stelle nach , 1yu i um:

2

, 1 ,1 1,0 1,0 22 4

y y x xu i u i u i u i h f t

.

Schubspannungsfreiheit liefert:

2 0y yx xxy xy

u uu u

y x y x

.

Ersetzen der Differentialoperatoren durch Differenzen ergibt:

1,0 1,0,1 , 1

2 2

y yx xu i u iu i u i

h h

Stelle nach , 1xu i um:

, 1 ,1 1,0 1,0x x y yu i u i u i u i .

Am Rand xK muss in negative x -Richtung eliminiert werden. Es gilt für die Normal-

spannung

0

2

yxxx

uu

x y

,

was diskretisiert

1, 1, 0, 1 0, 12

x x y yu j u j u j u j

ergibt. Für die yu -Komponente erhält man aus der Schubspannungsbedingung

74

1, 1, 0, 1 0, 1

2 2

y y x xu j u j u j u j

h h

,

was nach 1,yu j umgestellt wird:

1, 1, 0, 1 0, 1y y x xu j u j u j u j .

b) Die Ecken

Nun müssen mit einfachen Vorwärts-Rückwärtsschritten neue Approximationen fürdie gemischten Ableitungen auf den Ecken gefunden werden. Diese sind vom drei-dimensionalen Fall her für die Ecken ,0 ,

,n n l

E E bekannt:

2

,0 2

2

, 2

, , 1 1, 1 , 1,:

, , 1, , 1 1, 1:

n

n l

u x y u i j u i j u i j u i jE

x y h


x y h

(4.2.6)

Für die Ecken 0, 0,0,

lE E findet man:

2

0, 2

2

0,0 2

, 1, , 1, 1 , 1:

, 1, 1 , 1 1, ,:

l


x y h


x y h

(4.2.7)

Auf den Ecken werden die Differentialgleichungen (4.2.2) mit diesen gemischten Ab-leitungen diskretisiert. Damit sind die diagonal überhängenden Gitterpunkte beseitigt.

Es werden nun Eliminationsbeziehungen für die in einer Richtung überhängendenGitterpunkte benötigt.

In y -Richtung nimmt man die in zentralen Differenzen formulierte Randbedingung

des yK - und yK -Randes. Dies produziert zusätzlich Punkte, die in x -Richtung

überhängen. Diese werden mit einseitig formulierten Differenzen der xK -Ränder

eliminiert.

Zusammenfassend sind in der Tabelle 2 alle Eliminationsbeziehungen aufgeführt:

75

yK

/ 2

2 2, 1 ,1 1,0 1,0

2 4y y x x

i n

h f tu i u i u i u i

, 1 ,1 1,0 1,0x x y yu i u i u i u i

yK

, 1 , 1 1, 1,2

y y x xu i l u i l u i l u i l

, 1 , 1 1, 1,x x y yu i l u i l u i l u i l

xK 1, 1, 0, 1 0, 1y y x xu j u j u j u j

1, 1, 0, 1 0, 12

x x y yu j u j u j u j

xK 1, 1, , 1 , 1y y x xu n j u n j u n j u n j

1, 1, , 1 , 12

x x y yu n j u n j u n j u n j

0,0E

0, 1 0,1 1,0 1,02

y y x xu u u u

0, 1 0,1 1,0 1,0x x y yu u u u

1,0 0,0 0,1 0,02

x x y yu u u u

1,0 0,0 0,1 0,0y y x xu u u u

,0nE

, 1 ,1 1,0 1,02

y y x xu n u n u n u n

, 1 ,1 1,0 1,0x x y yu n u n u n u n

1,0 ,0 ,1 ,02

x x y yu n u n u n u n

1,0 ,0 ,1 ,0y y x xu n u n u n u n

,n lE

, 1 , 1 1, 1,2

y y x xu n l u n l u n l u n l

, 1 , 1 1, 1,x x y yu n l u n l u n l u n l

1, , , , 12

x x y yu n l u n l u n l u n l

1, , , , 1y y x xu n l u n l u n l u n l

76

Tabelle 2: Eliminationsbeziehungen des Randes.

Mit den Eliminationsregeln der Tabelle 2, sowie mit den modifizierten Differenzen fürdie gemischten Ableitungen kann das vollständige Gebiet berechnet werden.

0,lE

0, 1 0, 1 1, 1,2

y y x xu l u l u l u l

0, 1 0, 1 1, 1,x x y yu l u l u l u l

1, 0, 0, 0, 12

x x y yu l u l u l u l

1, 0, 0, 0, 1y y x xu l u l u l u l

77

4.3. Numerische Lösung der 2d-Bewegungsgleichungen

Über den Programmiercode gelten die gleichen Bemerkungen, wie sie unter 3.8. fürdie dreidimensionale Rechnung gegeben wurden. Der vollständige Programmcodeist im Anhang 6.14. angefügt, um die Rechnungen bei Interesse nachvollziehen zukönnen.

Die Störung dauert hier 0,6T µs , wie sie im Experiment durch einen Kugelaufprall

realisiert wurde.

Input:

Die Geometrie der Scheibe wird durch die Konstanten

,B L : Breite der Scheibe, Länge der Scheibe charakterisiert.

Das Material ist durch die Parameter

0, , : Lamésche Konstanten, konstante Dichte festgelegt.

Der Stoß wird berücksichtigt durch

0,6 ,T µs f t , die Stoßzeit und die Störungsfunktion.

Output:

Die Graphik stellt in der obersten Reihe links die Verteilung der yu -

Komponente des Verschiebungsfeldes auf einem xy -Schnitt räumlich und

rechts flächenhaft dar.

In der Reihe darunter ist die xu -Komponente in der gleichen Weise und für

den gleichen Zeitschritt dargestellt.

Auf die Dichteplots folgen die xy -Schnitte. Es handelt sich dabei für beide

Verschiebungskomponenten um Schnitte in y -Richtung bei 4,6,8x mm und

in x -Richtung bei 0,5,10y mm .

Auch in diesem Fall kann zum Test ein Spannungsstoß über die gesamte yK -Kante

gerechnet werden, was hier aber nicht, des Platzes wegen, angegeben werden soll.

78


79


80


81


82


83


84


85


86


87


88

4.4. Interpretation der numerischen Rechnungen

In einem ersten Schritt soll die Zuverlässigkeit der Rechnung beurteilt werden. Dazuwird ein einfaches Experiment vorgestellt, dass wie folgt aufgebaut ist:

1191L mm

Sensor1300B mm

0,012d m

Abbildung 47: Aufnahme der Verschiebung an einer PMMA-Platte, die in der Ebeneeinen Stoß erfährt.

Die Stahlkugel trifft die PMMA-Platte genau auf der Kante bei 650B mm . Der Kraft-

stoß erfolgt also in der Ebene. In der folgenden Abbildung ist das Verschiebungssig-nal des Experimentes dem der numerischen Rechnung gegenüber gestellt:

Verschiebung Platte

1,05

1,07

1,09

1,11

1,13

1,15

1,17

1,19

1,21

1,23

1,25

0,003 0,0035 0,004 0,0045 0,005 0,0055

t [s]

U[V

]

Abbildung 48: Vergleich der errechneten yu -Verschiebung am Punkt 4, 10x y

(unten) mit der yu -Verschiebung aus dem Experiment nach Abbildung 47 (oben).

89

Der Verschiebungsverlauf der numerischen Rechnung ist wegen der Empfindlich-keitscharakteristik des optischen Sensors mit 1 multipliziert. Die Übereinstimmungvon Rechnung und Experiment ist sehr gut. Deshalb werden den numerischen Re-sultaten vertraut, und diese zur Interpretation herangezogen.

Die in Abschnitt 3.9. zur Wellenausbreitung gemachten Aussagen gelten auch hier:

Maximum und Minimum der xu -Komponente werden während des Span-

nungsstoßes in 45°-Richtung in den Halbraum abgestrahlt. Man kann in dem folgenden Dichteplot entsprechend 0 0,68t µs sehen, dass

das Einflussgebiet der Verschiebung einem Halbkreis um den Koordinatenur-sprung mit dem Radius 0Lc t entspricht. D.h. die Wellenfront breitet sich allsei-

tig mit Longitudinalgeschwindigkeit Lc in Form einer Kreiswelle aus:

Abbildung 49: Lage der Wellenfront für 0 0,68t µs ist 0Lc t .

Das Maximum der yu -Komponente breitet sich während des Spannungssto-

ßes in Form eines Ellipsoids um den Koordinatenursprung aus, dessen großeHalbachse mit der y -Achse zusammenfällt.

Auf der Achse existiert nur die yu -Komponente, die sich dort als reine Longi-

tudinalwelle ausbreitet. Auf der Kante 0y existieren beide Verschiebungskomponenten. Die xu -

Komponente breitet sich hier, wie erwartet, mit Longitudinalgeschwindigkeitaus. Die yu -Komponente ist kurz hinter der Wellenfront extrem geschwächt.

In einer Entfernung entsprechend der Transversalwellengeschwindigkeit folgtdann mit stark erhöhter Amplitude die Oberflächenwelle, die aus der xu - und

yu -Komponente gebildet wird.

Die Wellenfront der yu -Komponente erfährt am der Erregung gegenüber lie-

genden Punkt eine heftige Reflexion ohne Phasensprung, wie an einem freienEnde, und wandert als Kreiswelle in die umgekehrte Richtung.

Die Wellenzüge unterliegen, entsprechend der linearen Theorie, keiner Dis-persion, und laufen nicht auseinander.

90

4.5. Folgerungen für Schallemissionssignale

In der folgenden Abbildung sind die errechneten Verschiebungs- und Signalverläufefür ausgewählte Punkte wiedergegeben.

Abbildung 50: Verschiebung und das Sensorsignal für die Punkte

4, 0 , 4, 10 , 8, 0x y x y x y und 8, 5x y .

Alle in Abschnitt 3.10. gemachten Bemerkungen über die Signale gelten auch hier:

Das Longitudinalsignal ist bei einem Laufweg s nach der Zeit / Lt s c allseitig

zu beobachten. Die Signalform des Longitudinalsignals (erste vollständige Schwingung) ist für

alle Messpunkte ähnlich, und folgt dem Signal im Krafteinleitungspunkt. (Abbildung 50)

Das an der unteren Ecke 8, 0x y auftretende Signal wird durch die Ober-

flächenwelle bestimmt und nimmt eine Sonderstellung ein, da es dem Signal-verlauf am Krafteinleitungspunkt am unähnlichsten ist:

91

Abbildung 51: Durch die Oberflächenwelle geprägtes Signal am Punkt 8, 0x y .

Der Einsatz der Transversalwelle geht auch im Fall der Scheibe mit signifikan-

ter Amplitudenvergrößerung einher. Beispiel: Das am Punkt 4, 10x y be-

rechnete Signal zeigt nach der Überlagerung zweier von links und rechtskommender yu -Maxima starke Amplitudenzunahme:

Abbildung 52: Das Signal der yu -Verschiebung am Punkt 4, 10x y wurde

wegen der Vergleichbarkeit mit dem Experiment mit 1 multipliziert. Der Einsatzvon Longitudinal- und Transversalwelle ist mit 1 und 2 gekennzeichnet.

Das errechnete Sensorsignal vom Punkt 4, 10x y aus Abbildung 50 soll nun mit

dem Sensorsignal des oben beschriebenen Experimentes verglichen werden:

1

2

1

2

92

Sensorsignal

1

1,05

1,1

1,15

1,2

1,25

1,3

1,35

1,4

1,45

1,5

0,01305 0,0131 0,01315 0,0132 0,01325 0,0133 0,01335

t [s]

U[V

]

Abbildung 53: Vergleich des errechneten Sensorsignals am Punkt 4, 10x y

(oben) mit dem Signal aus dem Experiment nach Abbildung 47 (unten).

Die Rechnung gibt den zeitlichen Verlauf, so auch die Asymmetrie zwischen ersterund zweiter Halbwelle, gut wieder.

Somit kann behauptet werden, dass der verwendete Piezosensor das Wellenfeldkorrekt abbildet.

93

5. Ein Experiment zur Prüfung der 1d-Reihenlösungen(6.3.33), (6.3.40) für einen zeitbegrenzten Spannungsimpuls

Im Anhang 6.3. werden Reihenlösungen für das 1d-Problem eines Spannungsimpul-ses der Zeitdauer T auf einen freien Stab abgeleitet. Ein Experiment soll klären, obdie Reihen (6.3.33) und (6.3.40) richtige Voraussagen über die Verschiebung amfreien Ende eines Stabes machen. Damit ist die Motivation des Experimentes erklärt:Es handelt sich um eine Prüfung theoretischer Überlegungen und nicht um ein„experimentum crucis“, wie es manchmal in der Physik vorkommt. [20] Dort soll einstrategisches Experiment zwischen zwei konkurrierenden Theorien entscheiden, hiergeht es allein um die Bestätigung oder Widerlegung von Rechnungen, die auf gesi-cherten, wenn auch unübersichtlichen mathematischen Grundlagen beruhen.

Die Schwierigkeit des auszuführenden Experimentes besteht im Wesentlichen auszwei Punkten:

Die zu erwartenden Verschiebungen liegen im Submikrometerbereich. Dahermuss die Messkette eine Auflösung in diesem Bereich garantieren.

Erfahrungen aus dem Schallemissionsbereich haben gezeigt, dass piezoelekt-rische Sensoren ihre Eigenschwingungen jeder Messung aufprägen. Daherwird von der Verwendung von Piezosensoren abgegangen und ein berüh-rungslos messender optischer Sensor verwendet.

Der prinzipielle Versuchsaufbau ist in der folgenden Abbildung erläutert.

3

Stahlstab Ø 0,0156m Sensor

k k k2L m

Abbildung 54: Versuchsanordnung zur Messung eines Verschiebungsimpulses aneinem freien Ende.

Ein Stahlstab von 2m Länge sei auf drei schmalen Elastomerlagern der Feder-konstante k gelagert. Ein am linken Rand aufgestelltes Fadenpendel der Länge

0, 465l m wird um 3 ausgelenkt. Die Kugel der Masse 0,0083Km kg und mit

dem Radius 0,0063kR m ruft einen longitudinalen Verschiebungsimpuls hervor.

Dessen Zeitsignal wird berührungslos von einem optischen Sensor am anderen, frei-en Ende aufgenommen. Das Signal wird auf die Messkarte eines PCs (Taktfrequenz

94

2MHz ) geleitet, die im Modus eines Speicheroszilloskops arbeitet. Der Sensor ist alsGlasfaserkabel ausgeführt, das fest mit einem nicht invertierenden Verstärker ver-bunden ist.

Abbildung 55: Faseroptischer Sensor D 20 von PhilTec mit Verstärkereinheit.

Durch einen Teil der Glasfasern wird Licht auf die zu messende Fläche geführt. Derandere Teil der Fasern absorbiert einen Teil des reflektierten Lichts und führt es demVerstärker zu. Eine typische Empfindlichkeitskurve ist in den Abbildungen 56 a und bdargestellt, die der Internetdarstellung der Firma PhilTec entnommen ist:

Abbildung 56 a: Typischer Verlauf der Empfindlichkeitskurve des optischen Ver-schiebungssensors Typ D 20 von PhilTec. Die Nahflanke ist die empfindlichere Seite

und wird in ihrem linearen Teil genutzt.

95

Abbildung 56 b: Typischer Verlauf der Fernflanke des optischen Verschiebungssen-sors Typ D 20.

Der Ausgleich unterschiedlicher Reflexionseigenschaften wird über die Verstärkungvorgenommen: Es wird das Maximum auf der Empfindlichkeitskurve aufgesucht unddann das Signal am Verstärker auf 5V geregelt. Anschließend wählt man einen Ar-beitspunkt auf dem linearen Teil der empfindlichen Flanke.

Das eigentliche Nutzsignal ist eine Wechselspannung, die der eingestellten Gleich-spannung am Arbeitspunkt aufmoduliert ist. Es kann mit dem verwendeten Modell D20 auf der Nahflanke eine Empfindlichkeit von

912, 215 10m

mV (5.1)

erreicht werden. Auf der unempfindlicheren Fernflanke ist es immerhin noch

9130 10m

mV . (5.2)

Lenkt man die Kugel, wie gesagt, um 3 aus ihrer Gleichgewichtslage aus und trifftdiese die Querschnittsfläche etwa zentrisch, wird folgendes Zeitsignal registriert:

96

0.05 0.10 0.15 0.20 0.25

-0.3

-0.2

-0.1

0.0

0.1

0.2

0.3

Kugel 3° AuslenkungU

[V]

t [s]

Abbildung 57: Gesamtes Zeitintervall. Man beachte die exponentielle Dämpfung unddie Modulation der langperiodischen Schwingung.

0.050 0.075 0.100 0.125 0.150

-0.3

-0.2

-0.1

0.0

0.1

0.2

0.3

Kugel 3° Auslenkung

U[V

]

t [s]

Abbildung 58: Vergrößert dargestellter Anfangsbereich: 0,0500 0,1500s t s , Gitter-

schritt 0,01t s .

97

0.050 0.055 0.060 0.065 0.070

-0.3

-0.2

-0.1

0.0

0.1

0.2

0.3

Kugel 3° AuslenkungU

[V]

t [s]

Abbildung 59: Wiederum vergrößerter Anfangsbereich: 0,0500 0,0725s t s , Gitter-

schritt 0,0025t s .

0.050 0.052 0.054 0.056

-0.3

-0.2

-0.1

0.0

0.1

0.2

0.3

Kugel 3° Auslenkung

U[V

]

t [s]

Abbildung 60: Nochmals vergrößerter Anfangsbereich: 0,0500 0,0575s t s , Gitter-

schritt 0,001t s .

98

Beginnt man mit der Betrachtung der größten Zeitauflösung, kann man unschwererkennen, dass die Verschiebung in Stufen erfolgt. Steile Verschiebungssprünge lö-sen sich mit Ruhephasen von ungefähr 0,7 ms ab. Das ist aber ungefähr die Zeit, die

zwischen zwei Reflexionen liegt, wenn man eine longitudinale Wellengeschwindigkeitvon 6000 /m s für Stahl zugrunde legt:

2 40,67

6000L

L mt ms

mcs

. (5.3)

Bei den Stufen handelt es sich also um die Reflexionen des Impulses am rechtenStabende. Dass die Spannung bei kleiner werdender Lücke zwischen Stab und Sen-sor sinkt, ist nicht verwunderlich, da auf der Nahflanke der Empfindlichkeitskurve einkleiner werdender Spalt mit einer sinkenden Spannung beantwortet wird. ( Abbil-dung 56) Um die Formeln (6.3.33) und (6.3.40) für den Verschiebungsverlauf aus-zuwerten, muss geklärt werden, welche dynamische Kraft durch den Kugelaufprallerzeugt wird. Betrachten wir also zunächst den Kraftstoß. Die Kugel wird um 3

ausgelenkt und losgelassen.

0, 465l m

24mmz

kR

kR

Abbildung 61: Geometrische Verhältnisse vor und nach dem Kugelaufprall.

Nun ist der Winkel in dem rechtwinkligen Dreieck gegeben durch

°

sin rad=180

z

l

. (5.4)

Die Horizontalentfernung des Kugelmittelpunkts von der Gleichgewichtslage ist

30,465 0,024

180z l m m

.

Im Rahmen idealelastischer Eigenschaften erreicht die Kugel wieder die Anfangsaus-lenkung, allerdings ist der Winkel entgegengesetzt zu orientiert:

99

. (5.5)

Durch die Auslenkung gewinnt die Kugel potentielle Energie, die beim Durchgang

durch die Gleichgewichtslage 0 vollständig in kinetische umgewandelt ist. Da-

raus kann die Endgeschwindigkeit vK der Kugel beim Aufprall errechnet werden.

lx

h

Abbildung 62: Hebung h der Kugel am Fadenpendel.

Also ist

Pot KinE E m1

2 g h m 2v v 2 K K g h . (5.6)

Die durch die Auslenkung gewonnene Höhe der Kugel ist

cos 1 cosh l x l l l . (5.7)

Damit ist die Kugelgeschwindigkeit vor dem Stoß bekannt:

v 2 1 cos K g l .

Dieser Ausdruck kann umgeformt und für kleine vereinfacht werden.

1 1

v 1 cos sin v 0,1122 2 2 2

K K

mg l g l g l g l

s.(5.8)

Aus (5.5) folgt für die Geschwindigkeit nach dem Stoß:

v

v v

K

K Kg l g l . (5.9)

Ausgangspunkt für das Weitere ist das 2. Newtonsche Axiom:

d

dt

pF . (5.10)

100

Die zwischen Kugel und Stab wirkende Kraft ist die Ursache der Impulsänderung.Gleichung (5.10) wird bei bekannter Impulsänderung über die Stoßzeit integriert:

0

0 0 0

T T TT

dp Fdt dp Fdt p Fdt . (5.11)

Der Kraftimpuls auf den Stab hat nach (2.2.10) die Form

0 sin F t F t .

Als Ursache für die Impulsänderung der Kugel ist er entgegen gerichtet zu nehmen:

0 00

0 0

20 sin v v cos

TT

K K K

F Ftp T p F t dt m

T. (5.12)

vK kann mit (5.9) eliminiert werden:

0v sin 02

v v 2 v0

K K

K K K K K

m t t TFm m F t

t T(5.13)

Die Stoßzeit erhält man als Lösung des Hertzschen Kontaktproblems. Die im Anhang6.2. abgeleitete Formel für die Kontaktzeit ist

1

22 5

2

3 142,94 ,

v 5 2

kK

K

RmT k D

k D E.

Mit den bekannten Konstanten

11

20,0083 v 0,112 0,0063 2,1 10 0,3 K K k

m Nm kg R m E

s m(5.14)

folgt für die Stoßzeit als oberer Integrationsgrenze

1

2 52 2

2

225 12,94 69

64 v

K

k K

mT µs

E R. (5.15)

In der folgenden Abbildung ist der Stoß vergrößert dargestellt, so dass das gemes-sene Verschiebungs- und Zeitintervall abgelesen werden kann.

101

Verschiebungssprung

-50

-45

-40

-35

-30

-25

-20

-15

-10

-5

0

0,0005 0,00052 0,00054 0,00056 0,00058 0,0006 0,00062 0,00064 0,00066 0,00068 0,0007

t [s]

U[m

V]

Abbildung 63: Stoßamplitude mit 41 7,5 33,5U mV mV mV

entsprechend 409 nm und der Stoßzeit 640 570 70t µs µs µs .

Im Anhang 6.12. wird die spektrale Darstellung des Kraftstoßes abgeleitet:

Abbildung 64: Kraftstoß im Zeit - und Frequenzraum.

Für die erste Nullstelle bei 24gf kHz und die Wellengeschwindigkeit 6000 /m s re-

sultiert eine minimale Wellenlänge von 0, 25m , was sehr viel größer als der Stabradi-

us von 0,0078m ist. Daher verhält sich der Stab wie ein eindimensionales System.

Mit den Laméschen Konstanten , , der Dichte von Stahl, sowie dem Stabradius,

11 2 10 2 301, 2 10 / 8,1 10 / 7800 / 0,0078N m N m kg m R m [19], (5.16)

können die Lösungsfunktionen (6.3.33), (6.3.40) mit

7

2 10 2 11 2 2 2

43,177,85 10

2 2 8,1 10 / 1, 2 10 / 0,0078

F NA

R N m N m m

(5.17)

102

ausgewertet werden. Zwischen den Abbildungen liegt ein Zeitintervall von 100t µs

in dem sich der Impuls um 0,6 m fortbewegt. Der Amplitudensprung nach der Refle-

xion von 400 nm stimmt gut mit der gemessenen ersten Reflexion überein, die in

Abbildung 63 dargestellt ist.

Abbildung 65: Auswertung der Reihenlösungen (6.3.33), (6.3.40). Ausbreitung desVerschiebungsimpulses bei freien Enden mit Reflexion bei 2z m .

103

Die Messung bezieht sich auf die Verschiebung des Stabendes. Darum wird dieRechnung für z L wiederholt.

Abbildung 66: Mit der Lösungsfunktion (6.3.40) errechnete Reflexionen am rechtenEnde des Stabes.

Es ist offensichtlich, dass das Ergebnis mit Abbildung 60, bis auf die durch die Ver-stärkercharakteristik hervorgerufene Invertierung, übereinstimmt. Der gemesseneVerschiebungssprung von 409 nm aus Abbildung 63 wird richtig beschrieben, und die

in (5.3) errechnete Laufzeit des Impulses mit 0, 67 ms wird ebenfalls reproduziert.

Wenn man die zeitliche Entwicklung der Verschiebung am rechten Stabende be-trachtet, wie sie in Abbildung 66 dargestellt ist, gelangt man zu dem Schluss, dassder Stab eine diskontinuierliche Bewegung, ähnlich der eines Inch Worm, vollführenmuss. Warum? An den Enden, die spannungsfrei sind, wächst die Verschiebung aufden doppelten Wert, während hinter dem durch den Stab wandernden Impuls dieVerschiebung nicht mehr auf ihren Ausgangswert zurück geht. So stößt einmal dasrechte Ende um den doppelten Betrag der verursachten Verschiebung vor. Dannzieht der Stab das linke Ende um den gleichen Betrag nach:

0 zu L

u 2u

3u 2u

3u 4u

Abbildung 67: Diskontinuierliche Bewegung eines Stabes der Länge L mit freien En-den durch einen Verschiebungsimpuls der Amplitude u .

104

Die in dem zugrunde gelegten Differentialgleichungsproblem angenommene Lage-rung des Stabes ist zweifellos pathologisch: Es gibt keine! Wie wird sich eine Lage-rung auf schmalen Elastomerstreifen, wie im Experiment, auswirken? Der Stab kannnicht endlos nach vorne rücken. Solange die Haftung zwischen den Elastomer-oberflächen und dem Stahl nicht überschritten wird, werden die Elastomere durchden Stab wie Federn gespannt. Es baut sich eine rückstellende Kraft auf, die denStab in seiner Vorwärtsbewegung zunehmend behindert. Wenn die Kugel vom Stababgehoben hat, kann man von einer freien Schwingung des Stabes auf den Lagernausgehen. Und genau das kann man in den Abbildungen 57 bis 60 sehen. Abhängigvon der Zeitskalierung ist entweder der Impuls oder die langperiodische Schwingungdes Stabes dominierend. In Abbildung 60 mit der feinsten Zeitauflösung ist ab derdritten Ruhephase ein leicht zunehmendes Rückstellen zwischen den Reflexionen zusehen. In Abbildung 57, die das vollständige Zeitintervall darstellt, können mehrerePerioden der freien Schwingung gesehen werden, auf welche die Reflexionen desLongitudinalimpulses aufmoduliert sind.Es soll eine Beschreibung für die langperiodische Schwingung des Stabes auf denElastomeren im Rahmen eines Starrkörpermodells gefunden werden. Das Systemkann durch einen Federschwinger mit der Federkonstante k ersetzt werden, auf den

die Kraft F t wirkt:

F t k

z z0

Abbildung 68: Modell für die Stabschwingung auf den Elastomerlagern.

Wie groß ist die den drei Elastomerlagern zugeordnete Federkonstante k ?

Dazu wurde (in diesem Fall) ein bekannter Kraftgeber in Form einer Spiralfeder mit629 /Fk N m benutzt. Mit diesem Kraftgeber wird nun untersucht, um welchen Weg

u sich der Stab auf den Elastomeren verschiebt. Dabei ist auf der Fernflanke desSensors zu messen und darauf zu achten, dass der Stab nicht rutscht.

opt. Sensor

Fk

l u u

Abbildung 69: Bestimmung der den Elastomerlagern zugeordneten Federkonstante.

Dann ist die den Elastomerstücken zugeordnete Federkonstante:

5

4

5

629 / 0,0045 6,76 104,1242 10

6,76 10F

N m m mk l u Nk

u m m

. (5.18)

105

Es gilt für t T die Bewegungsgleichung

22

20S

d zz

dt (5.19)

mit /S Sk m ( Sm Masse des Stabes), und für t T

2

2

2sinS

d zz B t

dt (5.20)

mit der Abkürzung v / K K SB m m . Die inhomogene Differentialgleichung hat nach

[21] die Lösung

1 2 2 2

1 2

sin cos sin

sin cos cos2

S S S

S

S S S S

S

BC t C t t

z tB

C t C t t t

(5.21)

Da k bekannt ist, kann ein Wert für S angegeben werden:

20

1118 53

S

S

k kT ms

m R L s. (5.22)

Da 45530 / s beträgt, ist nur der Fall S zu verfolgen. Die Konstanten

1 2,C C werden durch die Anfangsbedingungen festgelegt:

2 1 12 2 2 20 0 0 0S

S S S

B Bz C z C C

Also resultiert für die Zeit während des Stoßes

2 2

vsin sin

K KS

S S S

mz t T t t

m. (5.23)

Für die Zeit t T wird eine Zeitverschiebung t t T vorgenommen, und der Strichanschließend wieder weggelassen. Die Lösung von (5.19) ist hinlänglich bekannt:

1 2sin cosS Sz t C t C t .

Wieder erfolgt die Anpassung an die Anfangsbedingungen:

2 2

2 1 2 22 2

v v0 sin 0 1 cos

K K K K

S S S

S S SS S

m mz C T z C T

m m

106

Daraus folgt nach Rücksubstitution der Zeitvariablen die Lösung nach dem Stoß:

2

2 2

vsin sin

K K

S S

S S S

mz t T t T t

m. (5.24)

Um die Dämpfung zu simulieren, wird der Schwingung eine Exponentialfunktion te

überlagert. Dabei wird der Abklingkoeffizient über das logarithmische Dekrement der freien Schwingung bestimmt:

1 1 3ln ln

2

z t

T T z t T T

. (5.25)

Die fallende Charakteristik der Empfindlichkeitskurve auf der Nahflanke wird durchMultiplikation der Lösungsfunktion mit 1 berücksichtigt. Um ein Spannungssignal in

Volt zu erhalten, wird durch die Empfindlichkeit des Sensors von 912, 215 10 /m mV

und dann noch einmal durch 1000 dividiert.

6

1

12,215 10tV

U t e z tm

(5.26)

Abbildung 70: Freie Schwingung des Stabes auf drei Elastomerlagern.

Nun soll der freien Schwingung der Verschiebungsimpuls überlagert werden. Dazumüssen aus den Fourier-Reihen (6.3.33) und (6.3.40) die in z absoluten Glieder ge-strichen werden:

22 2

2 2 2 21

0

2 2

2 2 21

0

sin sin cos2, sin

cos2 2, sin sin

n n n n

n n n

nnn n

n n n

t t k zAc Ac Au z t T t t

L L L k

k zAc T Ac Au z t T t T t t T

L L L k

(5.27)

107

Warum? Sie entsprechen der über den Stab gemittelten Verschiebung, denn in einerFourier-Reihe

01

cosn nn

f z a a k z

mit 0

0

1L

a f z dzL

wird 0a wie das Mittel der Funktion f auf dem Intervall L berechnet. Diese mittlere

Verschiebung entspricht der Bewegung des Schwerpunktes. Letzteres kann mandurch die folgende Überlegung nachweisen:

Das System Kugel-Stab kann man natürlich unter dem Blickwinkel von Impuls- undEnergieerhaltung betrachten. Nach dem Ende des Kontaktes gilt / 0d dt p , der Im-

puls bleibt erhalten. Unter diesen Voraussetzungen und aus den Geschwindigkeiten

der Kugel und des Stabes vor dem Stoss v 0,112 / , v 0K Sm s können die Ge-

schwindigkeiten nach dem Stoss, v , vK S , errechnet werden:

2 2v v v 0,111 v v v 0,000622

S K S S K K

K S K S S K

S K S K S K S K

m m m m m mm m

m m m m s m m m m s

Der den Stab repräsentierende Schwerpunkt bewegt sich nach dem Stoss, ohne denEinfluss einer äußeren Kraft, geradlinig gleichförmig mit der Geschwindigkeitv 0,00062 /S m s . Genau diesen Wert erhält man nach zeitlicher Differentiation der

Reihen (6.3.33), (6.3.40) aus den zu streichenden Termen:

2 22

v 1 cos v 0,00062Ac Ac m

t T t t TL L s

.

Diese über den Stab gemittelten Geschwindigkeiten sind in der folgenden Graphikdargestellt:

Abbildung 71: Beschleunigungsphase ( bis 0,069 ms ) und Translation des Stab-

schwerpunktes während und nach dem Stoss.

Die über den Stab gemittelte Verschiebung ist in der folgenden Abbildung zu erken-nen:

108

Abbildung 72: Die treppenförmigen Verschiebungssprünge bei z L , und die überden Stab gemittelte Verschiebung.

Zu den modifizierten Reihen (5.27) ist die freie Schwingung (5.26) zu addieren:

6

1,

12, 215 10tV

U t u L t z t em

(5.28)

Im Ergebnis erhält man eine überzeugende Annäherung an Abbildung 59, einerHalbwelle der freien Schwingung mit aufmodulierten Verschiebungssprüngen:

Abbildung 73: Halbwelle der freien Schwingung mit aufmodulierten Verschiebungs-sprüngen entsprechend der Messung aus Abbildung 59.

Damit ist die langperiodische Schwingung der Abbildung 57 als freie Schwingung desStabes auf den Elastomerlagern identifiziert, und die Gültigkeit der Fourier-Reihenexperimentell gesichert.

Nun wird mit ihrer Hilfe die Beschleunigung am Ende z L bestimmt. Dazu mussdie Reihe (6.3.40) für t T zweimal nach der Zeit abgeleitet werden:

32

2 2 2 21

cos2, sin sin nn

n nn n n

k Lu AL t t t T

t L k

. (5.29)

109

Dies wird für das Zeitintervall kurz vor bis kurz nach Ankunft des Impulses bei z Lfür eine wachsende Anzahl von Reihengliedern ausgewertet.

Abbildung 74: Beschleunigung für eine wachsende Zahl von Reihengliedern.

Offensichtlich konvergiert die Reihe gegen eine Grenzfunktion.

Diese ist im Anhang 6.11. mit der d’Alembertschen Methode ermittelt worden. Zwei-fache zeitliche Ableitung von (6.11.19) an der Stelle z L ergibt die Funktion

22

02

0 /

2, cos / /

2

0 /

t L c

puL t t kL L c t L c T

t k

t L c T

(5.30)

mit der graphischen Darstellung:

110

Abbildung 75: Grenzfunktion der unendlichen Reihe für die Beschleunigung beiz L , bestimmt mit Hilfe der d’Alembertschen Methode. ( Anhang 6.11.)

Des Weiteren wird das Beschleunigungssignal der Form (5.30) und der Dauer0,7T ms mit einem Gaußimpuls gefaltet, der die Impulsantwort der Messelektronik

auf einen Stoß simuliert. Die Breite der Impulsantwort ist / 4T , und die negativeAmplitude gibt die invertierende Charakteristik des Sensors wieder.

2

t

h t e

. (5.31)

Abbildung 76: Idealisierte Impulsantwort des Messsystems.

Die Faltung wird durch das Integral

2

/2

2

/

, costL c T

L c

us t L h t d kL e d

t

(5.32)

111

definiert [22], die graphisch folgendes ergibt:

Abbildung 77: Faltung eines Beschleunigungssignals der Form (5.30) mit der Impuls-antwort (5.31).

Die obige Abbildung soll mit dem Signal einer realen Schallemissionsmessung ver-glichen werden:

Sensorsignal

1

1,05

1,1

1,15

1,2

1,25

1,3

1,35

1,4

1,45

1,5

0,01305 0,0131 0,01315 0,0132 0,01325 0,0133 0,01335

t [s]

U[V

]

Abbildung 78: Signal einer realen Schallemissionsmessung.

Beide Abbildungen sind ähnlich. Die Unterschiede lassen sich durch eine abwei-chende reale Impulsantwort und damit erklären, dass das reale Wellenausbreitungs-problem höherdimensional ist.

112

6. Anhang

6.1. Longitudinal- und Transversalwellen eines Diracimpulses

a) Die Quellenfunktion q ( [3])

Den Lösungen der Wellengleichung entsprechen ebene Wellen, wenn zuzüglichzur Differentialgleichung keine weiteren Einschränkungen der Lösungsmenge be-stehen. Diese Wellen kommen aus dem Unendlichen, und pflanzen sich ins Un-endliche fort.

Die Lösungen der Wellengleichung können folgenden Einschränkungen unterliegen:

1. Die Lösungen müssen der Wellengleichung zusammen mit Anfangsbedin-gungen genügen. (Cauchyches Anfangswertproblem) Beispiel: Kugelwellen

2. Die Lösungen müssen der Wellengleichung zusammen mit Anfangsbedin-gungen und Randbedingungen genügen. (Anfangs-Randwertproblem) Bei-spiel: stehende Wellen

Wellen, die sich infolge einer dynamischen Punktlast unbeeinflusst von den Rändernausbreiten, werden durch Probleme vom Typ 1 beschrieben.

Die Quellenfunktion ,q tx ist durch das Anfangswertproblem auf der gesamten re-

ellen Zeitachse für einen deltaförmigen Geschwindigkeitsimpuls definiert:

2

2

2, , 0 ,0 0 ,0

q qt c q t q t

t t

x x x x x . (6.1.1)

Das Anfangswertproblem für eine stetige Anfangsbedingung ,

2

2

2, , 0 ,0 0 ,0

u ut c u t u

t t

x x x x x , (6.1.2)

wird durch Faltung der stetigen Anfangsbedingung mit q gelöst:

, , x x y y yu t q t d . (6.1.3)

Das eindimensionale Anfangswertproblem für q und 0t löst:

1

, 2

0

x ctq x t c

x ct

(6.1.4)

Beweis: Einerseits wird das Anfangswertproblem (6.1.2) für eine stetige Anfangsbe-dingung durch Formel (6.1.3) formal gelöst. „Formal“ bedeutet, dass die konkreteQuellenfunktion hier nicht bekannt ist. Andererseits wird das gleiche Anfangswert-

113

problem durch die d’Alembertsche Formel gelöst. Beide Lösungen können gleichge-setzt werden:

1

,2

x ct

x ct

q x y t y dy y dyc

.

,q x y t ist deshalb entweder gleich 1/ 2c oder gleich Null, abhängig davon, ob

x y im Intervall , ct ct liegt oder nicht. Wenn x y durch x ersetzt wird, folgt die

Behauptung. ■

Dieses Verhalten kann mit der Heavyside-Funktion und mit der Signum-Funktion fürdie gesamte reelle Zeitachse dargestellt werden:

2 2 21, sgn

2q x t t H c t x t

c . (6.1.5)

Das dreidimensionale Anfangswertproblem für q und 0t löst:

1

,4

xx

x

ctq t

c. (6.1.6)

(Die Ableitung gelingt mit Hilfe der Kirchhoffschen Formel.) Mit der Beziehung

2 2 10, 0

2a x a x a x

a

kann das für die gesamte reelle Zeitachse wie folgt geschrieben werden:

22 21,

2q t c t t

c

x x . (6.1.7)

Diese Funktion nennt man auch Riemann-Funktion. Ihr Zusammenhang mit derGreenschen Funktion G des Problems wird im nächsten Punkt untersucht.

b) Die Greensche Funktion G

Im Folgenden soll der Zusammenhang der Quellenfunktion q mit der Greenschen

Funktion G aufgeklärt werden. Dazu ist es notwendig, das grundlegende Konzeptder Lösungsmethode mit Greenschen Funktionen darzustellen. Nehmen wir eine in-homogene partielle DGL der Form

0 x x xL u f t , (6.1.8)

wobei x ein drei- oder mehrdimensionaler Ortsvektor, Lx ein linearer partieller Diffe-

rentialoperator in drei oder mehr Variablen mit konstanten Koeffizienten, und

114

,u fx x Funktionen von drei oder mehr unabhängigen Variablen sind. Die

Greensche Funktion dieses Problems erfüllt die definierende Gleichung

, x x x x xL G , (6.1.9)

und repräsentiert am Ort x die Wirkung einer deltaförmigen Störung am variablen

Ort x . Multiplikation von (6.1.9) mit der beliebigen Testfunktion f x , und Integrati-

on über den x -Raum zusammen mit der Regel f d f x x x x x ergibt:

, x x x x x x x x x xL G f d f d f . (6.1.10)

Nach Vertauschung von Differentialoperator Lx und Integral erhält man

, x x x x x xL G f d f . (6.1.11)

Ein einfacher Vergleich von (6.1.11) und (6.1.8) führt zu der formalen Lösung

, x x x x xu G f d . (6.1.12)

Bisher war der Differentialoperator noch unbestimmt geblieben. Es soll speziell derWellenoperator

22

2L c

t

x (6.1.13)

betrachtet werden. Es sei die Lösung der Wellengleichung mit Störterm gesucht,

2

2

2,

uc u f t

t

x , (6.1.14)

wobei ,f tx eine Quellendichte, etwa eine von außen wirkende Volumenkraftdichte,

darstellt. Die Gleichung ist von der Form (6.1.8). Die Lösung ist mit (6.1.12) gegeben,wenn dazu mit (6.1.9) die Greensche Funktion bestimmt wird:

2

2

2, , , , , ,

x x x x x xG t t c G t t t t

t. (6.1.15)

Die Lösung G wird auch Fundamentallösung des Wellenoperators genannt.([23][24][25][26]) Damit ist das Problem auf die Bestimmung der Greenschen Funkti-on verschoben, die nicht immer einfach ist. Für t t fordert die Kausalität:

, , , 0 , , , 0

x x x xG

G t t t tt

. (6.1.16)

115

Nun soll untersucht werden, wie die Greensche Funktion G mit der Quellenfunktionq zusammenhängt.

Satz 1: Die Greensche Funktion G ist mit der Quellenfunktion q für t t identisch,

d.h. die Definitionen (6.1.1) und (6.1.15) sind für t t äquivalent:

, , , , , , x x x xG t t q t t H t t (6.1.17)

Beweis: Die Definition von q ist für eine Störung zum Zeitpunkt t t :

2

2

20 , 0 ,

qc q q t q t

t

x x x x . (6.1.18)

Die definierende Gleichung von G ist für t t :

2

2

2

x x

Gc G t t

t. (6.1.19)

Beide Definitionen haben die gleiche Laplace-Transformierte, und damit die gleicheLösung.Dazu werden beide Gleichungen in die neue Zeitvariable 0t t transformiert.Anwendung der Laplace-Transformation auf Problem (6.1.18) ergibt

2 2

0

2 2

, ,0 ,0 , 0

, ,

x x

x x x x

x x x x

s q s s q q c q s

s q s c q s

Anwendung der Laplace-Transformation auf Gleichung (6.1.19) ergibt

2 2, , x x x x s G s c G s ,

wobei , q G die Laplace-Transformierten von ,q G sind. ■

Eindimensional ist mit (6.1.4) und , 0x t bereits die Greensche Funktion gegeben:

1,

2 G x t H ct x

c. (6.1.20)

Zweidimensional ist mit rx und , 0x t die Greensche Funktion:

2 2 2

1,

2

H ct rG x t

c c t r. (6.1.21)

Dreidimensional ist mit rx und , 0r t in (6.1.6) bereits die Greensche Funktion

116

gegeben:

1

,4

ct rG r t

c r. (6.1.22)

Beispiel: Finde die Lösung des folgenden Problems für einen an die Stelle a ver-schobenen Gaußschen Geschwindigkeitsimpuls der halben Breite :

2

22

2

20 ,0 0 ,0

x au u

c u u x x et t

.

Wegen der Beschränkung auf 0t kann die 1d-Greenfunktion (6.1.20) mit demGaußimpuls nach (6.1.3) gefaltet werden:

2 2

1 2 1, ,

2

x a x ax ct

x ct

u x t q x x t x dx H ct x x e dx e dxc c

.

Die Sprungfunktion legt die Integrationsgrenzen fest, denn sie stellt ist ein Rechteckin x dar, das im Intervall x ct x x ct gleich 1 ist, und sonst verschwindet.

Abbildung 79: Ausbreitung eines Gaußschen Geschwindigkeitsimpulses bei500a m mit der halben Breite 120 m und der Geschwindigkeit 600 /c m s nach

0,05 , 0, 25t s s und 1,5 s .

117

c) Aufspaltung des Verschiebungsfeldes in einen divergenz- und rotationsfrei-en Anteil

Die Lamésche Differentialgleichung in der Form (2.1.32) soll mit verschwindendenVolumenkräften in einem unbegrenzten Raum betrachtet werden:

2

2 2 2

2grad div

uu uT L Tc c c

t. (6.1.23)

Es sei an die Helmholtzsche Zerlegung des Verschiebungsfeldes in

L T u u u (6.1.24)

erinnert. Der tranversale Anteil genüge der Bedingung

div 0T u , (6.1.25)

der longitudinale Anteil erfülle

rot 0L u . (6.1.26)

Setze (6.1.24) in die Bewegungsgleichung ein, und beachte, dass div 0T u ist:

2 2 2 grad divL T T L T L T Lc c c u u u u u . (6.1.27)

Wende auf diese Gleichung den div -Operator an:

2 2 2div div div divgrad divL T T L T L T Lc c c

u u u u u

.

Laplace- und div -Operator vertauschen, und mit (6.1.25) erhält man

div divL Tu u 2 div divT L Tc u u 2 2

2

div

div div

L T L

L L L

c c

c

u

u u

Nun kann der div -Operator herausgezogen werden:

2div 0L L Lc u u . (6.1.28)

Laplace- und rot -Operator vertauschen, und mit (6.1.26) verschwindet rot -Bildungdes Klammerausdrucks (6.1.28) ebenfalls:

2rot 0L L Lc u u . (6.1.29)

Wenn Rotation und Divergenz eines Vektors im gesamten Raum verschwinden, soist der Vektor ein Nullvektor, und man kann folgern:

118

22

20L

L Lct

uu . (6.1.30)

Wende nun auf (6.1.27) den rot -Operator an:

rot Lu 2 2 2

0

rot rot rot grad divT T L T L T Lc c c u u u u

.

Laplace- und rot -Operator vertauschen, und mit (6.1.26) erhält man

2rot 0T T Tc u u . (6.1.31)

Offensichtlich verschwindet mit (6.1.25) auch die Divergenz des Klammerausdrucks:

2div 0T T Tc u u . (6.1.32)

Mit gleicher Argumentation, wie für die Wellengleichung (6.1.30) erhält man

22

20T

T Tct

uu . (6.1.33)

Die Gleichungen (6.1.30), (6.1.33) sind gewöhnliche dreidimensionale Wellenglei-chungen.

Tu hängt nicht mit Volumenänderungen zusammen div 0T u , dann aber mit

Scherungen.

Lu hängt nicht mit Scherungen zusammen rot 0L u , dann aber mit Kom-

pressionen und Dilatationen.

Die Wellengeschwindigkeiten ,L Tc c sind eindeutig den entsprechenden Ausbrei-

tungsmechanismen zugeordnet.

d) Longitudinal- und Transversalwellen

Mit den Gleichungen (6.1.30), (6.1.33) wurden also Bewegungsgleichungen für denrotationsfreien Anteil,

22

2rot 0 0L

L L Lct

uu u ,

und für den divergenzfreien Anteil,

22

2div 0 0T

T T Tct

uu u ,

abgeleitet. Beide Vektoren, Lu und Tu , besitzen drei skalare Komponenten:

119

, , , ,l l l t t tL i j k T i j ku u u u u u u u . (6.1.34)

Oben wurde das Problem einer deltaförmigen Störung,

2

2

2

x x

Gc G t t

t,

dreidimensional gelöst. Es resultierte die 3d-Greenfunktion (6.1.22):

1

,4

ct rG r t

c r.

Welches sind nun die Orte gleicher Phase, Orte, an denen die ,l ti iu u die gleiche Amp-

litude haben? Das ist genau dort der Fall, wo das Argument der Deltafunktion Null ist:

0ct r r ct . (6.1.35)

Da r unabhängig von der Richtung ist, sind das für festes t die Punkte einer Kugel-

schale. Jede der drei skalaren Komponenten von Lu und Tu breitet sich infolge einer

in Raum und Zeit hinlänglich konzentrierten Störung auf dünnen Kugelschalen mitdem Radius r ct aus.

Nun führen zwei allgemeine Sätze der Physik zur Lage der Wellenfront und Wellen-geschwindigkeit:

Satz 2: Orte gleicher Phase definieren eine Wellenfront. Im Raum periodische Wel-len haben unendlich viele Wellenfronten, Diracimpulse nur eine: die Kugelschale.

Satz 3: Die Geschwindigkeit von Orten gleicher Phase (Phasengeschwindigkeit)steht senkrecht auf der Wellenfront. Sie gibt die Ausbreitungsrichtung der Welle an.Legt man an in die Punkte der Kugelschalen Tangentialflächen, so stehen die Wel-lengeschwindigkeiten Lc und Tc senkrecht auf ihnen, und zeigen in re -Richtung.

,L Tc c

rre

Abbildung 80: Wellenfront und Lage der Wellengeschwindigkeiten einer Kugelwelle.

An dieser Stelle kann noch nicht gesagt werden, in welcher Richtung die Verschie-bung in Bezug auf die Bewegungsrichtung der Welle erfolgt. Der Verschiebungsvek-tor Lu wird wie folgt allgemein angesetzt:

0, lL Lt f c t u x x u . (6.1.36)

120

Über die Richtung von 0lu sei, wie gesagt, nichts angenommen. Wir haben in (6.1.26)

vorausgesetzt, dass Lu der Bedingung rot 0L u genügt:

0rot rot 0lL Lf c t u x u .

Mit den Regeln rot rot gradf f f a a a und a b b a erhält man

0 0 0

0

rot grad grad 0l l lL L Lf c t f c t f c t x u u x x u .

Nun ist grad xLf c t :

, , , , , ,L L

x y zf c t f c t f f f f x y z f

x y z

x x xx x x

.

Also erhält man

0 0lLf c t x x u . (6.1.37)

Weil Lf c t x im allgemeinen nicht verschwindet, muss das Vektorprodukt ver-

schwinden.

Untersuche nun den Verschiebungsvektor Tu :

0, tT Tt f c t u x x u . (6.1.38)

Wir haben in (6.1.25) vorausgesetzt, dass Tu der Bedingung div 0T u genügt.

Wende den div -Operator auf Tu an:

0div , div 0tT Tt f c t u x x u .

Mit der Regel div div gradf f f a a a erhält man

0 0

0

div grad 0t tT Tf c t f c t x u u x ,

und es bleibt mit Lf c t f x x :

0 0tTf c t x x u . (6.1.39)

Weil Tf c t x im allgemeinen nicht verschwindet, muss das Skalarprodukt ver-

schwinden.

121

Interpretiere nun die Gleichungen (6.1.37) und (6.1.39):

0 0 0 00 || 0l l t tr r r r e u u e e u u e . (6.1.40)

Der Einheitsverschiebungsvektor 0lu liegt also parallel zum Einheitsradiusvektor re .

Die Wellengeschwindigkeit Lc steht senkrecht auf der Wellenfront, ist also ebenfalls

parallel zu re . Da sowohl 0lu als auch Lc parallel zu re liegen, sind sie auch zueinan-

der parallel:

0lu

Lc

rre

Abbildung 81: Lage von Einheitsverschiebungsvektor 0lu und Wellengeschwindigkeit.

Der Amplitudenvektor dieser Wellen besitzt die Richtung der Wellenausbreitung. Mannennt solche Wellen Longitudinalwellen.

Der Amplitudenvektor 0tu steht senkrecht auf dem Einheitsradiusvektor re . Die Wel-

lengeschwindigkeit Tc steht senkrecht auf der Wellenfront, ist also parallel zu re . Al-

so steht 0tu senkrecht auf Tc , sie sind orthogonal:

0tu

Tc

rre

Abbildung 82: Lage von Einheitsverschiebungsvektor 0tu und Wellengeschwindigkeit.

Der Amplitudenvektor dieser Wellen steht senkrecht auf der Wellenausbreitungsrich-tung. Man nennt solche Wellen Transversalwellen. Natürlich kann der Amplituden-

vektor 0tu in der Blattebene oder senkrecht zur Blattebene senkrecht auf Tc stehen.

Im ersten Fall spricht man von hS -Wellen (horizontal polarisierte Sekundärwellen), im

zweiten Fall spricht man von vS -Wellen (vertikal polarisierte Sekundärwellen).

122

6.2. Das Kontaktproblem

Die Behandlung des Kontaktproblems folgt der Darstellung von Landau, Lifschitz [5].(Siehe auch [27].) Zwei Körper seien in einem Berührungspunkt O in Kontakt. Dortfallen die beiden Tangentialebenen in der xy -Ebene zusammen. Die z -Achse ver-

läuft positiv in das Innere der Körper und ist mit ,z z bezeichnet.

z

y

x

z

Abbildung 83: Lage der Kugeln im Koordinatensystem.

In der Nähe des Berührungspunktes kann die Oberfläche des oberen Körpers durchdie Gleichung

1 2, 1,2z g x x x x x y (6.2.1)

beschrieben werden. g ist ein metrischer Tensor, der die Krümmung der Fläche

charakterisiert. Analog gilt für die Oberfläche des unteren Körpers

1 2, 1,2z g x x x x x y . (6.2.2)

Beide Körper werden durch eine äußere Kraft F aufeinander gedrückt, sodass siesich um die Strecke h nähern. Für alle Punkte der Berührungsebene gilt

z zz u z u h . (6.2.3)

Einsetzen von (6.2.1) und (6.2.2) ergibt

z zg g x x u u h . (6.2.4)

123

Nun legt man die x - und y -Achse so, dass der Tensor g g Diagonalform be-

sitzt. Die Eigenwerte des Tensors heißen A und B und mit ihnen geht (6.2.4) über in

2 2z zAx By u u h . (6.2.5)

Die Größen ,A B hängen mit den Krümmungsradien 1 2,R R und 1 2,R R beider Flächen

über die Beziehungen

1 2 1 2

2 2

2

1 2 1 2 1 2 1 2

1 1 1 12

1 1 1 1 1 1 1 14 2cos 2

A BR R R R

A BR R R R R R R R

(6.2.6)

zusammen. ist der Winkel zwischen den Normalschnitten der Oberflächen, in de-

nen die Krümmungsradien 1 2,R R sind. Nun sei (6.2.6) auf den Fall eines Kugelauf-

pralls auf die ebene Fläche eines Zylinders angewendet.

z

y

x

z

Abbildung 84: Lage der Körper und Normalschnittebenen im Koordinatensystem.

Die Kugel habe die Krümmungsradien 1 2R R R , die Ebene 1 2R R R mit

R . Also geht (6.2.6) über in

22 2 2

2 lim 4 0R

A B A BR R R

. (6.2.7)

Also folgt für die Eigenwerte A und B des metrischen Tensors:

1 1 1 1

2 2A B

R R . (6.2.8)

124

Die Druckverteilung im Berührungsgebiet wird durch die Integralgleichung

2 22 2,1 1 1 zp x y

dx dy h Ax ByE E r

(6.2.9)

beschrieben, wobei ,v E die Poissonzahl und der Elastizitätsmodul sind. Für einen

Kontakt Stahl auf Stahl ist (6.2.9) natürlich

2

2 22 1 ,zp x y

dx dy h Ax ByE r

. (6.2.10)

Die Lösung findet man in [5]:

2 2

2 21z

x yp C

a a . (6.2.11)

a ist der Radius des Kontaktgebietes. Die Konstante C wird auf die vorgegebeneKraft F normiert:

2 22 2 2 3

2 2 2

Halbkugelvolumen

1 4 31

2 3 2z

x y C C Fp dxdy C dxdy a x y dxdy a F C

a a a a a

Also folgt für den Druck im Kontaktgebiet

2 2

2 2 2

31

2z

F x yp

a a a . (6.2.12)

Setzt man (6.2.12) in (6.2.10) ein

2 2

22 2

2 2

2

12 1 3

2

x yF a a dx dy h Ax By

E a r

und nutzt die Umformung [5]

2 22 2

2 2 22 2

20

11

2

x yx ya a aa a dx dy d

r a

,

ergibt sich die Integralgleichung

2 2

2 2 22 2

20

13 1

2

x yF a a

d h Ax ByE a

. (6.2.13)

125

Mit der Abkürzung

23 1

2D

E

(6.2.14)

kann die linke Seite von (6.2.13) als Summe dreier Integrale geschrieben werden:

2 2 2 2

2 22 2 20 0 0

1 1 1F D F D F Dd d x d y h Ax By

a a a

Die Gleichung muss im Kontaktgebiet 2 2 2x y a identisch erfüllt sein. Daher müs-

sen die sich entsprechenden Glieder der linken und rechten Seite gleich sein, wasfolgende Gleichungen ergibt:

2 22 2 2

0 0 0

1 1 1

F D F D F D

h d A d B da a a

.(6.2.15)

A und B sind bei einem kreisförmigen Kontaktgebiet gleich, was auch in (6.2.8) zumAusdruck kommt. Die letzten beiden Gleichungen drücken den Radius a des Kon-taktgebietes durch die vorgegebene Kraft F aus, denn A ist bekannt. (6.2.8)

1

2A

1 F D

R

23

3

1a F D R

a (6.2.16)

Anschließend kann die Annäherung der Körper h unter der Wirkung der Kraft F mitder ersten Gleichung von (6.2.15) berechnet werden:

F Dh

2 2

33

F D F Dh

a RF D R. (6.2.17)

Des Weiteren ist die potentielle Energie U der sich berührenden Körper

UF U Fdh

h

.

Stellt man (6.2.17) nach h um, ergibt einfache Integration

5

22

5

RU h

D . (6.2.18)

Wichtig für weitere Untersuchungen ist die Bestimmung der maximalen Annäherungder Körper 0h und der Kontaktzeit T . Ohne Energiedissipation ist die kinetische

Energie der Kugel kurz vor dem Aufprall gleich der potentiellen Energie der Verfor-mung beim Erreichen der maximalen Annäherung 0h der Körper. Mit (6.2.18) kann

man den Energiesatz formulieren:

126

52 2

0

2v

2 5K

K

m Rh

D. (6.2.19)

(Index K für die Kugelgeschwindigkeit) Einfaches Umstellen führt auf

2

52

0

5v

4

K K

Dh m

R. (6.2.20)

Die Kontaktzeit T folgt ebenfalls aus dem Energiesatz, wenn die Annäherungsge-schwindigkeit durch /dh dt und die potentielle Energie durch (6.2.18) ausgedrücktwird. Die Gesamtenergie, die durch die kinetische Energie der Kugel kurz vor demAufprall gegeben ist, ist in jedem Moment des Kontaktes gleich der Summe aus kine-tischer Energie und potentieller Energie:

2 2 25 5 52 2 22 2 2

52 2

5 52 22 2

2 4v v v

2 2 5 5

1 1v

v v

ges potkin

K KK K K K K

KkE EE

K

K

K K

K K

m m dh R dh R dh kh m m h h

dt D dt D dt m

dh k dth dt dh

dt m dh k kh h

m m

Die Integration erfolgt im Intervall 00 2h h . Ausnutzung der Symmetrie ergibt

0

50 2 2

12

v

h

K

K

T dhk

hm

.

Die Substitution 1

2 4 2 5v / K Kh m k x ergibt die neuen Integrationsgrenzen 0 1x .

Damit ist das Integral analytisch auszuwerten, auch wenn die Lösung aus nichtele-mentaren Funktionen aufgebaut ist:

1 112 25 5

2 250 2

7

1 52 2

9v v1 10

K K

K K

m mT dx

k kx

.

Gammafunktion, 4 / 5 x k R D Damit lautet der Ausdruck für die Kontaktzeit

12 5

22,94

v

K

K

mT

k. (6.2.21)

127

6.3. Analytische Lösung der 1d-Wellengleichung auf einem endli-chen Intervall mit inhomogenen Randbedingungen

Die eindimensionale Wellengleichung, über deren Ableitung man sich etwa in [28]informieren kann, wird zuerst für die Zeit während des Stoßes untersucht.

a) Zeit während des Stoßes: 0 t T

2 22

2 20 0 0

RB: 0, sin , 0

AB: ,0 0 ,0 0

u uc t T z L

t z

u ut A t L t

z z

uu z z

t

(6.3.1)

Mit der vorliegenden inhomogenen Randbedingung ist das Problem durch einen Se-parationsansatz nicht zu lösen. Einen Ausweg bietet die Transformation auf ein Prob-lem mit homogenen Randbedingungen und inhomogener Differentialgleichung. [29]

2 22

2 20

u uc

t z

RB 00,u

t tz

, L

uL t t

z

AB ,0u z f z , 0u

z g zt

, , ,u z t v z t z t

2 2

2

2 2,

v vc z t

t z

RB 0, 0tz

, 0

vL t

z

AB ,0 ,0v z f z z , 0 ,0v

z g z zt t

Das transformierte Problem kann durch Entwicklung nach Eigenfunktionen gelöstwerden. Es sei in allgemeiner Weise ein in den Randbedingungen inhomogenesProblem auf einem Intervall 0 z L gegeben. [30]

2 22

2 2

0

0 0

RB: 0, ,

AB: ,0 ,0

L

u uc t T

t z

u ut t L t t

z z

uu z f z z g z

t

(6.3.2)

128

Es wird eine Transformationsbeziehung, und damit eine neue Variable eingeführt:

, , , , , ,u z t v z t z t v z t u z t z t . (6.3.3)

Differenzieren nach z ergibt

, , ,v u

z t z t z tz z z

. (6.3.4)

Wegen der Forderung nach homogenen Randbedingungen kann man schreiben

!

0

!

0, 0, 0

, , 0L

vt t t

z z

vL t t L t

z z

(6.3.5)

Die Forderungen können nur erfüllt werden, wenn gilt:

00, , Lt t L t tz z

. (6.3.6)

Ein besonders einfacher Zusammenhang (er ist nicht der einzige), der dies reprodu-ziert, ist die Beziehung

0 0, L

zz t t t t

z L

.

Die transformierende Funktion ist die Stammfunktion zum obigen Ausdruck:

2

0 0,2

L

zz t z t t t C t

L .

Über die Funktion C t kann frei verfügt werden. Sie sei im Hinblick auf Vereinfa-

chungen des konkreten Problems in der folgenden Weise gewählt:

2 2

0 022

L Ac cC t t t t

L L

.

Dann kann als ein vollständiges Quadrat mit Zusatztermen geschrieben werden:

22 2 2

0 02,

2 2L

z Lz Ac cz t t t t t

L L L L

. (6.3.7)

Die Transformationsbeziehung lautet nun

22 2 2

0 02, ,

2 2L

z Lz Ac cu z t v z t t t t t

L L L L

. (6.3.8)

129

Einsetzen von (6.3.8) in (6.3.2) führt auf

2 2 2 22 2 2 2 2 22 0 0

02 2 2 2 2 2

,

2 2L

L

z t

z Lv v c c c zc

t z L t L t L L L t

. (6.3.9)

Das ist die transformierte PDGL mit einer Inhomogenität ,z t , die ausschließlich

die Randbedingungen des Ausgangsproblems enthält. Die neuen Anfangsbedingun-gen kann man aus (6.3.2) und (6.3.7) ohne Mühe bestimmen.

22 2

0 02

22 2 2

0 02

,0 ,0 ,0 ,0 0 0 02 2

,0 ,0 ,0 ,0 0 0 02 2

L

L

z Lz cv z u z z v z f z

L L L

z Lz Ac cv z u z z v z g z

L L L L

(6.3.10)

Damit ist die Transformation des Problems abgeschlossen:

2 2 2 22 2 2 22 2 20 0 0

2 2 2 2 2 2

22 2

0 02

22 2 2

0 02

2 2

RB: 0, 0 , 0

AB: ,0 0 0 02 2

,0 0 0 02 2

L L

L

L

z Lv v c zc c c

t z L t L t L L L t

v vt L t

z z

z Lz cv z f z

L L L

z Lv z Ac cz g z

t L L L L

(6.3.11)

Um die Funktionen 0 , , ( ), ( )Lt t f z g z zu bestimmen, vergleicht man das konkre-

te Problem (6.3.1) mit dem allgemeinen Problem (6.3.2).

0 sin 0 0 0Lt A t t f z g z . (6.3.12)

Geht man damit in das Problem (6.3.11) ein, ergibt sich

2 2 222

2 2

2

, sin2

RB: 0, 0 , 0

AB: ,0 0 ,02

v v Ac z t z L t

t z L

v vt L t

z z

v Av z z z L

t L

(6.3.13)

2 2

2

2

2, sin

2

Ac c Az t t z L t

L L

(6.3.14)

130

Wenn auf der rechten Seite der DGL (6.3.13) als Inhomogenität Deltafunktionen inOrt und Zeit stehen, nennt man ihre Lösung „Fundamentallösung“. Zur Lösung desProblems wird die Methode der Entwicklung nach Eigenfunktionen verwendet. [30],[31] Im ersten Schritt wird das Eigenfunktionssystem aufgefunden. Dazu wird die Se-

parationsmethode auf das Randwertproblem ohne Inhomogenität, , 0z t , ange-

wandt. Ansatz:

,v z t Z z T t (6.3.15)

Eingesetzt in die homogene PDGL (6.3.13) erhält man

2

Z T

Z c T

,

wobei eine beliebige Konstante ist. Jede Seite ist für sich gleich der Konstanten:

2I: 0 II: 0Z Z T c T . (6.3.16)

Die neuen Randbedingungen sind

0, 0 0 0 0

, ( ) ( ) 0 0

v t Z T t Z

v L t Z L T t Z L

(6.3.17)

Nun erfolgt eine Fallunterscheidung bezüglich der Separationskonstanten . Es sinddrei Fälle möglich:

0 0 0 . (6.3.18)

Die DGLen (6.3.16) haben in allen Fällen unterschiedliche Form. Der Fall 0 be-inhaltet keine Schwingungslösungen, weshalb er nicht weiter zu verfolgen ist:

2

2 2

2

0

0

I 0 I 0

II 0 II 0

k

T k c T T

Z k Z Z

(6.3.19)

Untersuche den Fall 2k : Die DGL II ist die Schwingungsgleichung, deren Lö-sung zum Standard der Theorie gewöhnlicher Differentialgleichungen gehört:

1 2cos sinZ z C kz C kz . (6.3.20)

Die Anwendung der Randbedingungen liefert:

2 2 10 0 0 sin 0Z C k C Z L C k kL . (6.3.21)

131

Wenn man der Annahme entgehen will, dass 1 0C ist (triviale Lösung), dann folgt

die interessante Tatsache, dass k folgender Diskretisierung genügen muss:

0,1,2,kL n n . (6.3.22)

Es folgen die Definitionen von Wellenzahl und Kreisfrequenz der n -ten Mode:

n n n

n n ck k c

L L

. (6.3.23)

Damit lauten die ortsabhängigen Eigenfunktionen

cos cosn n n n

n zZ z C C k z

L

. (6.3.24)

Untersuche den Fall 0 :

1 2 1 20 0Z z C z C Z C Z z C .

Insgesamt trägt der Fall 0 eine in z konstante Funktion zur Gesamtlösung bei,

die freilich im allgemeinen von der Zeit abhängen, und durch 0C t für den Index

0n berücksichtigt wird. Da die PDGL linear ist, setzt sich die allgemeine Lösungaus der Summe aller Einzellösungen zusammen. [3]

0 1 1cos cos k kZ z C t C t k z C t k z . (6.3.25)

Daraus resultiert der Fourieransatz

01

, cos

n nn

v z t b t b t k z . (6.3.26)

Der nach Ort und Zeit differenzierte Fourieransatz

20

1 1

, cos , cosn n n n nn n

v z t b t b t k z v z t b t k k z

wird in die PDGL (6.3.13) eingesetzt

2

22 20

1 1

cos cos sin2

n n n n nn n

Ab t b t k z c b t k k z z L t

L

.(6.3.27)

Um linke und rechte Seite vergleichen zu können, muss auch die rechte Seite in einecos - Fourierreihe entwickelt werden, was eine gerade Fortsetzung der Funktion überdas Intervall L z L hinaus bedeutet. Die Parabel ist um L nach rechts ver-schoben und gestreckt:

132

Abbildung 85: Gerade Fortsetzung der Parabel 22 / 2A L z L mit 2L .

Es ist

22

0 2 2 220

212

2

0

1

3 cos2

2 62 4cos

L

n

Ln n

n n

n

La z L dz

L k zA AL Az L

L L ka z L k z dz

L k

.

Eingesetzt in die Differentialgleichung (6.3.27) erhält man

2 2

2 20 2

1 1 1

2cos cos sin sin cos

6n n n n n n

n n n n

AL Ab t b t k z c b t k k z t t k z

Lk

Die Summen können zusammengefasst werden:

2 2

2 20 2

1 1

2cos sin sin cos

6n n n n n

n n n

AL Ab t b t c k b t k z t t k z

Lk

.

Daraus folgt durch Vergleich:

2

0

22 2

2

sin6

2sinn n n

n

ALb t t

Ab t c k b t t

Lk

(6.3.28)

Aus der ersten Gleichung findet man 0b . Sie kann unmittelbar integriert werden.

133

0 1

0 1 2

cos6

sin6

ALb t t C

ALb t t C t C

(6.3.29)

Die Integrationskonstanten sind durch die Anfangsbedingungen festgelegt.

0 01

2

01

,0 0 0 cos 0 0 0 0

,0 0 0 cos2

n n nn

n nn

v z b b k z b b

Av z b b k z z L

L

(6.3.30)

Nun wird die rechte Seite entwickelt, um die 0nb zu erhalten:

2

21

cos2

2 6n

n n

k zA AL Az L

L L k

.

Graphisch:

Abbildung 86: Gerade Fortsetzung der Parabel 2

/ 2A L z L .

Eingesetzt in die zweite Differentialgleichung (6.3.30) erhält man

0 21 1

cos20 0 cos

6n

n nn n n

k zAL Ab b k z

L k

.

Daraus folgt für die 0nb :

0 2

20 0

6n

n

AL Ab b

Lk

.

134

Nun können die Integrationskonstanten bestimmt werden:

0 2

0 1 1 0

0 0

0 0 sin6 6 6

b C

AL AL ALb C C b t t

(6.3.31)

Die zweite Gleichung von (6.3.28) ist eine inhomogene gewöhnliche DGL vom Typ

2

22 2 2

2

2sin sin

n

n

n n n n n

n

Ab t c k b t t y t y t t

Lk

.

Sie hat für n die Lösung

2

1 2 2 2 2

2sin cos sinn n n

n n

Ab t C t C t t

Lk

. [21]

Die Konstanten werden wieder unter Beachtung der Anfangsbedingungen festgelegt.

2

3 2

1 12 2 2 2 2 2 2

2 2

2 2 2 2 2 2

0 0

2 2 20 1

2 21 sin sin

n

n n

n n n n n n

n n

n n n n n

b C

A A Ab C C

Lk Lk L k

A Ab t t t

L k Lk

Einsetzen in den Fourieransatz (6.3.26) ergibt die Lösung des Teilproblems:

2 2 21

sin sin cos2, sin

6n n n

n n n

t t k zAL Av z t t

L k

. (6.3.32)

(Lösung mit Greenscher Funktion 6.4.) Die Lösung des Problems (6.3.1) setzt sich

aus , , ,u z t v z t z t zusammen:

2 2

2

2 2 2 21

sin sin cos2 2, sin sin

2 6n n n

n n n

v

t t k zAc c A AL Au z t t z L t t

L L L k

Nun wird 2

z L durch 2

21

cos4

3n

n n

k zL

k

ersetzt und zusammengefasst:

22 2

2 2 2 21

sin sin cos2, sin n n n n

n n n

t t k zAc Ac Au z t t t

L L L k

(6.3.33)

(Beweis, dass (6.3.33) Lösung von (6.3.1) ist: 6.5.)

135

b) Zeit nach dem Stoß: t T

Nach dem Ende des Impulses, beim Erreichen der Zeit T , ändert sich das Problem.Es liegt nunmehr ein in den Randbedingungen homogenes Problem vor. Darum wirdeine Verschiebung des Zeitnullpunktes an die Stelle t T vorgenommen:

t t T . (6.3.34)

Leite (6.3.33) nach der Zeit ab:

22 2

2 2 21

cos2, cos cos cos nn

nn n n

k zAc Ac Au z t t t t

L L L k

.(6.3.35)

Nun sind ,u z T und ,u z T die Anfangsbedingung für das neue Problem:

2 22

2 2

2

2 2 21

22

2 2 21

0 0

(0, ) 0 ( , ) 0

sin cos2( ,0)

1 cos cos2 2( ,0)

n n n

n n n

n n n

n n n

u uc t

t z

u ut L t

z z

T k zAc T Au z

L L k

T k zu Ac Az

t L L k

(6.3.36)

Da das Problem (6.3.36) in den Randbedingungen homogen ist, kann es durch einenSeparationsansatz gelöst werden. [22],[32]

Die Orts-Eigenfunktionen dieses Problems sind mit denen des transformierten Prob-lems identisch, und aus (6.3.25) zu entnehmen.

Die Zeitanteile ergeben sich aus der Untersuchung der DGL I (6.3.19). Das sind für2k die Lösung der Schwingungsgleichung, und für 0 ist es eine lineare

Funktion:

1 2 3 41

, cos sin cosn n n n nn

u z t C C t C t C t k z

. (6.3.37)

Die freien Konstanten sind durch die Anfangsbedingungen bestimmt:

2

1 3 2 2 21 1

sin2,0 cos cosn n

n n nn n n n

TAc T Au z C C k z k z

L L k

,

woraus sofort

136

2

1 3 2 2 2

sin2 n nn

n n

TAc T AC C

L L k

(6.3.38)

folgt. Vergleich der Anfangsgeschwindigkeiten

22

2 4 2 2 21 1

1 cos cos2 2,0 cos

n n nn n n

n n n n

T k zAc Au z C C k z

L L k

liefert

2

2 4 2 2 2

1 cos2 2 n n

n

n n

TAc AC C

L L k

. (6.3.39)

Damit ist die Lösungsfunktion für t T unmittelbar gegeben durch

2 2

2 2 21

cos2 2, sin sin nn

n nn n n

k zAc T Ac Au z t t T t t T

L L L k

(6.3.40)

(Beweis, dass (6.3.40) Lösung von (6.3.36) ist: 6.6.)

c) Dispersionseigenschaften

Es soll an den vorliegenden Wellenfunktionen untersucht werden, ob die Dispersi-onseigenschaften longitudinaler Wellen richtig wiedergegeben werden. Bekannt ist,dass dieser Wellentyp keine Dispersion zeigt. Dies ist immer dann der Fall, wenn dieGruppengeschwindigkeit gleich der Phasengeschwindigkeit ist:

gr phc c c . (6.3.41)

Die in (6.3.23) eingeführte Kreisfrequenz

n nk c (6.3.42)

stellt die Dispersionsbeziehung der n -ten Mode der Lösungsfunktion dar. Die Grup-pengeschwindigkeit [33] ist

ngr

n

dc c

dk

. (6.3.43)

Das Kriterium (6.3.41) trifft also zu. Daher bleibt die Form des Impulses in der Zeiterhalten. Ganz analoge Verhältnisse bestehen für die Lösungsfunktion, die für dieZeit nach dem Ende des Verschiebungsimpulses gilt.

137

6.4. Lösung des Differentialgleichungsproblems (6.3.13) mit HilfeGreenscher Funktionen

Es ist auch möglich, das inhomogene Anfangs - Randwertproblem (6.3.13)

2 2 222

2 2

2

, sin2

RB: 0, 0 , 0

AB: ,0 0 ,02

v v Ac z t z L t

t z L

v vt L t

z z

v Av z z z L

t L

mit Hilfe Greenschen Funktionen zu lösen.

Die allgemeine Lösungsformel für nicht verschwindende Anfangswerte lautet [21]

0 0 0 0

, , , , , , , ,L L t L

v z t G z t f d G z t g d G z t d dt

,

wobei G die Greenfunktion, und ,0 0, ,0v z f z v z g z die Anfangswerte

bezeichnen. Das Problem dieser Methode ist, dass die Greenfunktion für die gege-benen Randbedingungen konstruiert werden muss. Für das Problem (6.3.13) ist

1

2 1, , , cos cos sinn n n

n n

tG z t k k z t

L L

. [30]

Einsetzen der Greenfunktion, der Inhomogenität und der Anfangswerte in die Lö-sungsformel ergibt

22

10 0

2

10

2 1, cos cos sin sin

2

2 1cos cos sin

2

t L

n n nn n

L

n n nn n

t Av z t k k z t L d d

L L L

t Ak k z t L d

L L L

Nun werden die Integrale in

22

2

0 0

22

21 0 0

2 2

2 210 0

, sin2

1cos sin sin cos

1cos sin cos

2

t L

t L

n n nn n

L L

n n nn n

Av z t t L d d

L

Ak z t L k d d

L

A At L d k z t L k d

L L

gelöst. Die Integrale lassen sich analytisch auswerten.

138

2

0

2 2

0

32

0

2

2

0

sinsin

sin sinsin sin

3

2cos

t

tn n

n

n

L

L

n

n

t tt d

t tt d

LL d

LL k d

k

Einsetzen und Umformen ergibt:

2 2 21

sin sin cos2, sin

6n n n

n n n

t t k zAL Av z t t

L k

in Übereinstimmung mit (6.3.32). ■

139

6.5. Beweis, dass (6.3.33) das Problem (6.3.1) löst

Randwerte:

Bilde u in der Form v :

2 21

sin sin sin2sin n n n

n n n

t t k zA Au v z L t

L L k

. (6.5.1)

Dann ist 0, sinu t A t und , 0u L t .

Anfangswerte:

,0 0u z ist sofort zu sehen, und in (6.3.35) ist

22 2

2 2 21

00

cos2,0 1 1 0nn

n n n

k zAc Ac Au z

L L L k

.

Erfüllt die Funktion die Differentialgleichung?

Die zweifachen partiellen Ableitungen, von (6.5.1) nach z , und von (6.3.35) nach t :

22

2 2 21

2 21

cos2sin sin sin

sin sin2sin cos

nnn n

n n n

n nn

n n

k zAc Au t t t

L L k

t tA Au t k z

L L

Man kann unmittelbar2 2

2

2 20

u uc

t z

sehen, wenn 2 2 2/n nc k beachtet wird:

22

2 2 21

2 2

2 21

cos2sin sin sin

sin sin2sin cos 0

nnn n

n n n

n nn

n n

k zAc At t t

L L k

t tA Ac t c k z

L L

Damit ist die Richtigkeit der Lösung bewiesen. ■

140


Randwerte:

2 2

1

sin2sin sin nn

n nn n n

k zAu t t T

L k

Dann ist 0, 0u t und , 0u L t .

Anfangswerte:

2

2 2 21

22

2 2 21

sin cos2,

cos2 2, 1 cos

n n n

n n n

nnn

n n n

T k zAc T Au z T

L L k

k zAc Au z T T

L L k

was offensichtlich mit (6.3.36) übereinstimmt.


22

2 2 21

3

2 2 21

2 21

cos2 2cos cos

cos2sin sin

2sin sin cos

nnn n

n n n

nnn n

n n n

nn n n

n n

k zAc Au t t T

L L k

k zAu t t T

L k

Au t t T k z

L

Man kann unmittelbar2 2

2

2 20

u uc

t z

sehen, wenn 2 2 2/n nc k beachtet wird:

3

2 2 21

2

2 21

cos2sin sin

2sin sin cos 0

nnn n

n n n

nn n n

n n

k zAt t T

L k

Ac t t T k z

L

Damit ist die Richtigkeit der Lösung bewiesen. ■

141

6.7. 2d-Bewegungsgleichungen für den Ebenen Spannungszustandund für Axialsymmetrie

Die Bewegungsgleichung ist

u σ F , (6.7.1)

wobei σ der Spannungstensor ist. Dessen Divergenz ist

, 1, 1 2, 2 3, 3 σ e e e eij i j i i i i i i .

Komponentenweise ist (6.7.1):

,j ij i ju F . (6.7.2)

Nun wird in ,j i i j umindiziert. Wegen der Symmetrie des Spannungstensors

ij ji kann man schreiben:

,i ij j iu F . (6.7.3)

Die Bewegungsgleichung für beliebige Raummetrik ist

, j ij i jm j mi ji im miv F , (6.7.4)

mit den Christoffel-Symbolen erster und zweiter Art:

, , ,

1

2r rs

ijs si j sj i ij s ij ijsg g g g . (6.7.5)

Die Elemente des ko- und kontravarianten Metrischen Tensors sind

2 2

1 0 0 1 0 0

0 0 0 1/ 0

0 0 1 0 0 1

ijijg r g r

. (6.7.6)

Die Christoffel-Symbole erster Art sind also bestimmt durch die Ableitungen der Ele-mente des kovarianten Metrischen Tensors. Für Zylinderkoordinaten:

122 212

221

0 sonstijkr

r

(6.7.7)

Dann sind die Christoffel-Symbole zweiter Art für Zylinderkoordinaten:

1 11 1222 221 222

2 221 12

0 sonst

1/

ijkg g r

r

(6.7.8)

142

Schreibe die Komponenten der Bewegungsgleichung (6.7.4) hin:

1 1 1 1 1,

2 2 2 2 2,

3 3 3 3 3,

i i m mii im mi

i i m mii im mi

i i m mii im mi

v F

v F

v F

Summiere über i :

1 11 1 1 1 1 21 2 1 1 2 31 3 1 1 3 1,1 1 1 ,2 2 2 ,3 3 3

2 12 1 2 2 1 22 2 2 2 2 32 3 2 2 3 2,1 1 1 ,2 2 2 ,3 3 3

3 13 1 3 3 1 23 2 3 3 2,1 1 1 ,2 2 2 ,3

m m m m m mm m m m m m

m m m m m mm m m m m m

m m m mm m m m

v F

v F

v

33 3 3 3 3 33 3

m mm m F

Beachte bei der Summierung über m , dass alle rij außer 1 2 2

22 21 12, , verschwinden:

1 11 21 31 2 11 1 22 1,1 ,2 ,3 21 22

2 12 22 32 2 21 2 12 2,1 ,2 ,3 21 12

3 13 23 33 2 31 3,1 ,2 ,3 21

2

v F

v F

v F

Man erinnere sich an die Bedeutung der Christoffel-Symbole zweiter Art aus (6.7.8):

1 11 21 31 11 22 1,1 ,2 ,3

2 12 22 32 21 12 2,1 ,2 ,3

3 13 23 33 31 3,1 ,2 ,3

1

2 1

1

v r Fr

v Fr r

v Fr

Ersetze in diesen Gleichungen in der folgenden Weise:

11 22 33

2

12 21 13 31 23 32

1 2 3

1 2 3

1

1 1

rr zz

r rz z

r z

r z

r

r r

FF F F F F

r

uv u v v u

r

(6.7.9)

Dann sind die allgemeinen Bewegungsgleichungen in Zylinderkoordinaten:

143

1

1 2

1

rrr rz rrr r

r zr

zrz zz rzz z

u Fr r z r r

u Fr r z r

u Fr r z r

(6.7.10)

Des Weiteren sollen alle Volumenkräfte iF verschwinden.

In Sonderheit soll ein Ebener Spannungszustand (ESZ) betrachtet werden,in dem alle Spannungen in z -Richtung verschwinden: 0iz . Die Scheibe

kann sich in dieser Richtung frei bewegen:

1

1 2

0

rrr rrr

rr

z

ur r r r

ur r r

u

(6.7.11)

Weiterhin wird Axialsymmetrie vorausgesetzt. D.h. alle Größen sind, bei fes-tem ,r z , unabhängig vom Winkel, und es gilt / 0 .

2

rr rrr

rr

ur r r

ur r

(6.7.12)

Man beachte, dass die vorgegebenen Einschränkungen immer noch Torsion zulas-sen. Nun soll zu den Verschiebungen übergegangen werden, indem die Elementedes Spannungstensors mit (2.2.13) ersetzt werden:

1 12

1 12

1 1

r zrr r

r zr

rr

uu uu

r r r z

u u uu

r r r z

uuu

r r r

Eliminiere in ,rr die Dehnung /zu z durch

1 1 1

2 02

z r z rzz r r

uu u u uu u

z r r r z r r

.

Daraus folgt

144

4 42 1 1 2

2 2 2 2r r

rr r r

u uu u

r r r r

. (6.7.13)

Setze das und den Ausdruck für r in die Bewegungsgleichungen (6.7.12) ein:

2

2 2

2

2 2

4 1 1

2

1 1

r rr r

u uu u

r r r r

u uu u

r r r r

Definiere nun die Abkürzungen

2 24: :

2L Tc c

. (6.7.14)

Die Belastung sei so, dass keine Schubspannungen r existieren. Daraus folgt,

dass in (6.7.12) die Bewegungsgleichung für u entfällt. Die verbleibende Differenti-

algleichung

2 22

2 2 2

1 1r r rL r

u u uc u

t r r r r

(6.7.15)

ist keine Wellengleichung der Form 2r L ru c u , und darum auch nicht mit einem An-

satz in Form laufender Wellen zu lösen.

145

6.8. Analytische Lösung des Anfangs-Randwertproblems auf einerScheibe

Es wird nun versucht, das Problem als Anfangs-Randwertproblem auf einer Kreis-scheibe mit einem Radius R zu lösen. Die Anfangsbedingungen beinhalten einenradialen deltaförmigen Geschwindigkeitsimpuls bei r R :

,0 0 ,02

r

r

r Ruu r r

t r. (6.8.1)

Die Scheibe sei am Rand frei, was der Randbedingung , 0 rr R t entspricht:

2 22

2 2 2

1 10 0

RB: , 0

AB: ,0 0 ,02

r r rL r

rr

rr

u u uc u r

t r r r r

R t

r Ruu r r

t r

(6.8.2)

Es handelt sich um ein in den Anfangsbedingungen inhomogenes Problem. In einemersten Schritt wird das Eigenfunktionssystem des homogenen Randwertproblems

aufgefunden. Durch Einsetzen des Separationsansatzes ,ru r t X r T t in die

PDGL erhält man

2

2

1 1LX T c X T X T X T

r r

. (6.8.3)

Division durch X T und Gleichsetzen beider Seiten mit der Separationskonstanten führt auf

2

2

1 1L

T X Xc

T X r X r

.

Es resultieren die gewöhnlichen Differentialgleichungen

2 2

2

I 0

II 1 0L

T t T t

r X r rX r r X rc

(6.8.4)

Nun erfolgt eine Fallunterscheidung bezüglich der Separationskonstanten . Es sinddrei Fälle möglich:

0 0 0 . (6.8.5)

Die DGLen (6.8.4) haben in allen Fällen unterschiedliche Form. Der Fall 0 bein-haltet keine Schwingungslösungen, weshalb er nicht weiter zu verfolgen ist:

146

2

2

22 2 2

2

0

0

I 0 I 0

II 1 0 II 0L

T T T

r X rX r X r X rX Xc

(6.8.6)

Die Gesamtlösung ist die Summe beider Fälle.

Untersuche nun für 2 die DGL II. Durch Einführung der Wellenzahl

22

2:

L

kc

(6.8.7)

geht die Gleichung über in

2 2 2 1 0r X rX k r X . (6.8.8)

Durch die Substitutionen

2 22

2 2

d d dz d d dkr z k k

dr dz dr dz dr dz(6.8.9)

entsteht die Besselgleichung [34]:

2 2 1 0z X z zX z z X z . (6.8.10)

Sie wird von Besselfunktionen erster und zweiter Art und erster Ordnung gelöst:

1 1 2 1X z C J z C Y z .

Für 0r soll die Lösung beschränkt bleiben. Da 10

lim

r

Y kr gilt, ist 2 0C . Es

bleibt

1 1X r C J kr . (6.8.11)

Die Randbedingung ergibt nach (6.7.13) für den Ortsanteil:

4 2 1 1

, 0 02 2 2

rrr r

r R r R

uR t u X r X r

r r r.(6.8.12)

Nun werden für X die Eigenfunktionen (6.8.11) eingesetzt, wobei k die Variable ist:

147

1 1

10

2

J kR J kR

kR. (6.8.13)

Mit 1 0 1

1 J kr J kr J kr

krist das

0 1

2 10

2

J kR J kR

kR. (6.8.14)

Für diese Bestimmungsgleichung existieren positive Nullstellen ik k :

1 20 mk k k .

Die ersten zehn Nullstellen sind hier für 1R numerisch berechnet:

0 1 2 3 4

5 6 7 8 9

2,04792 5,38908 8,57168 11,7317 14,8837

18,0321 21,1785 24,3236 27,4679 30,6117

k k k k k

k k k k k

(6.8.15)

Für alle anderen Nullstellen können in sehr guter Näherung die Nullstellen von 0J k

verwendet werden, da in (6.8.14) der Term proportional 1/ k praktisch verschwindet.

Die vollständige Lösung ist dann die Summe aller Einzellösungen:

1 10

ii

X r C J k r . (6.8.16)

Untersuche nun für 0 die DGL II:

2 0r X rX X . (6.8.17)

Das ist die Eulergleichung

2 0x y x axy x by x

mit der Lösung [21]

11

2 221 2

11 4

2

a

y x x C x C x a b

. (6.8.18)

Offensichtlich sind in diesem Fall 1 1a b und die Betragsstriche in der Lösungs-

funktion überflüssig, da 0r :

1 2

1X r C r C

r .

148

Für 0r wächst der zweite Summand unbeschränkt. Daher muss 2 0C sein. Die

Randbedingung ergibt die Einschränkung

1 1 1

2 30

2 2

C C C .

Das kann nur gelten, wenn 1 0C ist, und die Lösung für 0 verschwindet.

Untersuche für 2 die DGL I von (6.8.6): 2 0T T . Das ist die Schwingungs-gleichung, deren Lösung hinlänglich bekannt ist:

1 2sin cosT t C t C t . (6.8.19)

Also besitzt ru die folgende Darstellung:

1 2 10

, sin cosr i L i L ii

u r t C k c t C k c t J k r

. (6.8.20)

Die freien Konstanten werden durch die Anfangsbedingungen bestimmt:

2 1 20

1 10

,0 0 0

,02

r ii

r i L ii

u r C J k r C

r Ru r C k c J k r

r

(6.8.21)

Die Deltafunktion wird ebenfalls im Eigenfunktionssystem 1 iJ k r entwickelt [23]:

10 2

i i

i

r Ra J k r

r. (6.8.22)

Multiplikation mit der Summe über alle 1 jJ k r ergibt

1 1 10 0 0 2

i i j j

i j j

r Ra J k r J k r J k r

r.

Das Produkt links kann zu einer Doppelsumme zusammengefasst werden:

1 1 10 0 0 2

i i j j

i j j

r Ra J k r J k r J k r

r

Nun wird Integriert:

1 1

0 0 0 2

R

i i ji j

r Ra J k r J k r rdr

r 1 jJ k r r

0 0

R

j

dr . (6.8.23)

149

Die Orthogonalitätsbeziehung für Besselfunktionen ist

21 1 2

0 1 0 2

0

2

R

i j

i i i

i j

J k r J k r rdr RJ k R J k R J k R i j

(6.8.24)

Das Integral ergibt nur etwas für i j , und das Integral über die Deltafunktion kann

problemlos ausgewertet werden:

2

21 0 2 1

0 0

1

2 2

i i i i ii i

Ra J k R J k R J k R J k R .

Zusammenfassen der Summen ergibt

2 21 0 2 1

0

10

i i i i ii

R a J k R J k R J k R J k R .

Da die Summe Gliedweise verschwinden muss, erhält man

1

2 21 0 2

1

i

i

i i i

J k Ra

R J k R J k R J k R, (6.8.25)

und dann für die Anfangsbedingung

1

12 20 1 0 2

1

2

i

ii i i i

J k R r RJ k r

R J k R J k R J k R r. (6.8.26)

Dann erhält man für die freie Konstante 1C aus (6.8.21)

11 1 12 2

0 0 1 0 2

1

1 2 21 0 2

1

1

ii L i i

i i i i i

i

L i i i i

J k RC k c J k r J k r

R J k R J k R J k R

J k RC

c R k J k R J k R J k R

(6.8.27)

Dann ergibt sich endlich für die Verschiebung:

1 1

2 20 1 0 2

sin1,

i i L i

riL i i i i

J k R k c t J k ru r t

c R k J k R J k R J k R. (6.8.28)

(Beweis, dass (6.8.28) das Problem (6.8.2) löst: Anhang 6.9.)

In Abbildung 87 ist die Ausbreitung eines -förmigen Geschwindigkeitsimpulses bei1r m zu sehen. Die 2d-Fundamentallösung wird am Koordinatenursprung wie an

einem festen Ende reflektiert. Es wurde hier mit 200 Termen gerechnet.

150

Abbildung 87: Ausbreitung einer -Störung bei 1r m als 2d-Fundamentallösung.

151


Randwerte:

1 1 1 1

2 2 2 20 01 0 2 1 0 2

1

2 21

4 2 1, ,

2 2

4 sin sin1 2 1

2 2

4 sin1

2

rrr r

r R

i i L i i i L i

i iL i i i L i i i i

i i L

L i

uR t u R t

r R

J k R k c t J k R J k R k c t J k R

c R J k R J k R J k R R c R k J k R J k R J k R

J k R k c t

c R J k R J

1 1

0 0 2

0

2

i ii i i i

J k R J k Rk R J k R k R

Anfangswerte:

112 2

0 1 0 2

6.8.261

12 20 1 0 2

sin 01,0 0

1,0

2

ir i

iL i i i i

ir i

i i i i

J k Ru r J k r

c R k J k R J k R J k R

J k R r Ru r J k r

R J k R J k R J k R r


21 1

2 2 20 1 0 2

1 1

2 20 1 0 2

21 1

2 2 20 1 0 2

sin1

sin1

sin1

i i L iri L

i i i i

i i L ir

iL i i i

i i L i ir

iL i i i

J k R k c t J k ruk c

t R J k R J k R J k R

J k R k c t J k ru

r c R J k R J k R J k R

J k R k c t k J k ru

r c R J k R J k R J k R

Eliminiere die Ableitungen der Besselfunktionen mit ik r z , 0 1J z J z und

1

1 0 1

1 0 1 0 1 1 1 1 0 12 2

1 02

1

1 1 1 1 1 1

2 11

J z

J z J z J zz

dJ z J z J z J z J z J z J z J z J z J z

dz z z z z z z

J z J zz z

152

21 1

2 2 20 1 0 2

10 12 2

0 1 0 2

21

2 2 2 2 20 1 0 2

sin1

sin1 1

sin1 21

i i L iri L

i i i i

i i Lri i

iL i i i i

i i Lri

iL i i i i

J k R k c t J k ruk c

t R J k R J k R J k R

J k R k c tuJ k r J k r

r c R J k R J k R J k R k r

J k R k c tuk

r c R J k R J k R J k R k r 1 0

1

i i

i

J k r J k rk r

Ist nun2 2

2

2 2 2

1 10r r r

L r

u u uc u

t r r r r

? Die linke Seite ist

11 02 2 2 2

0 1 0 2

1 1 120 12 2 2 2

0 01 0 2 1 0 2

sin1 2 11

sin sin1 1 1 1

i i Li i i

iL i i i i i

i i L i i i Li L L i

i ii i i L i i i i

J k R kc tk J k r J k r

c R J k R J k R J k R k r k r

J k R kc t J k r J k R kc tkc c J k r J k

R J k R J k R J k R r c R J k R J k R J k R k r

112 2 2

0 1 0 2

sin1 1

i

i i Li

iL i i i i

r

J k R kc tJ k r

r c R k J k R J k R J k R

Streiche die sich kompensierenden Terme:

1

2 2 2 20 1 0 2

1 1 2

2 20 1 0 2

sin1 2

sin1

i i Li

iL i i i i

i i L ii L L

i i i i

J k R kc tk

c R J k R J k R J k R k r

J k R kc t J krkc c

R J k R J k R J k R

1 0

11

i i

i

J kr J k rk r

102 2

0 1 0 2

sin1 1

i i L

iiL i i i

J k R kc tJ k r

r c R J k R J k R J k R 1

1 i

i

J k rkr

1

12 2 21 0 2

sin1 1

i i L

i

L i i i i

J k R kc tJ k r

r c R k J k R J k R J k R0

i

Es bleibt

1 1 1212 2 2 2

0 01 0 2 1 0 2

sin sin1 10

i i L i i i L

i L L i ii ii i i L i i i

J k R k c t J k r J k R k c tk c c k J k r

R J k R J k R J k R c R J k R J k R J k R

Damit ist bewiesen, dass die Reihe (6.8.28) das Differentialgleichungsproblem(6.8.2) erfüllt. ■

153

6.10. Fourierkoeffizienten

Angestrebt wird die Darstellung einer Funktion durch eine unendliche Reihe vonWinkelfunktionen.

Es wird im Folgenden vorausgesetzt, dass die Funktion f z gerade, und das Inter-

vall L z L definiert und 2L -periodisch ist. Es kann eine Kosinusreihe Ver-

wendung finden:

0 01

cos cosn kn

n z k zf z a a a a

L L

. (6.10.1)

Um die na zu berechnen wird (6.10.1) mit cos /n z L multipliziert und von

L z L unter Ausnutzung der Orthogonalitätseigenschaft

0cos cos

L

L

n kk x n xdx

L n kL L

integriert:

0cos cos cos cosL L L

k

L L L

n z n z k z n zf z dz a dz a dz

L L L L

Da die Funktion f z gerade ist, braucht nur über das Intervall 0 z L integriert

und mit 2 multipliziert zu werden, um die na zu erhalten:

0

2 cosL

n

n zf z dz L a

L

.

Damit können die na der Reihe (6.10.1) nach der Formel

0

2cos

L

n

n za f z dz

L L

, (6.10.2)

und das 0a mit

0

0

1L

a f z dzL

(6.10.3)

berechnet werden.

154

6.11. Die d’Alembertlösung des 1d-Problems

Es sei das Problem eines sinusförmigen Druckimpulses p t auf einen Stab ge-

stellt. [35],[36] Der Stab habe die Länge L und der Druckimpuls habe die Dauer T .

p t z

0 L

Abbildung 90: Ein Druckimpuls wirkt bei 0z auf einen Stab der Länge L .

Das Problem ist:

2 22

2 2

0

0

sin, 0 0,

2

,0 0 ,0 0

u uc

t z

H T t p tu uL t t

z z

uu z z

t

(6.11.1)

Nach d’Alembert ist die allgemeine Lösung der 1d-Wellengleichung:

, u z t f kz t g kz t . (6.11.2)

Die Anfangsbedingungen fordern

,0 0

,0 0

u z f kz g kz

u z f kz g kz(6.11.3)

(Die Striche bezeichnen Differentiation nach dem Argument.) Ableitung der erstenGleichung nach z und Addition beider Gleichungen ergibt:

00

0

f kz g kzg kz g kz A

f kz g kz. (6.11.4)

A ist eine Konstante. Es folgt für (6.11.2)

,0

f kz t A z ctu z t

z ct(6.11.5)

Für 0t gilt nach der ersten Gleichung von (6.11.3) und nach (6.11.4)

f kz g kz A . (6.11.6)

Also sind die Anfangsbedingungen erfüllt.

155

Die inhomogene Randbedingung ist im Zeitintervall 0 t T :

0 sin

0,2

p tut

z.

Andererseits ist mit (6.11.5) die linke Seite dieser Gleichung:

0,

ut kf t

z.

Die Ableitung der Verschiebung kann durch Gleichsetzen eliminiert werden:

0 sin

2

p tf t

k.

Jetzt wird eine Variablensubstitution vorgenommen:

s kz t . (6.11.7)

An der Stelle 0z entsprechend s t wird die Gleichung überführt in:

0 sin

2

p sf s

k. (6.11.8)

Die Stammfunktion ist

0 cos

2

p sf s B

k. (6.11.9)

Sie muss mit der Anfangsbedingung (6.11.6) verträglich sein. Für 0s erhält man:

002

pf B A

k.

Dann ist die Konstante B :

0

2

pB A

k.

Also ist die Funktion f s im Intervall 0 s entsprechend 0 t T :

01 1 cos

2

pf s s A

k. (6.11.10)

Für s t T ist das linke Ende spannungsfrei, woraus mit (6.11.8) folgt:

156

0 f s f s C .

Die Konstante C kann mit 1f s an der Stelle s bestimmt werden:

01

2

2

pf A C

k.

Also ist die Funktion f s für s :

02

2

2

pf s A

k. (6.11.11)

Die Verschiebung (6.11.5), formuliert in s , ist mit 1 2,f s f s :

2

1 0

0 0

f s A s

u s f s A s

s

Rücksubstitution von s ergibt die Lösung für den elastischen Halbraum:

0

0

20

2

, 1 cos2

0

pz c t T

k

pu z t kz t c t T z ct

k

z ct

(6.11.12)

In (6.11.1) wird an der Stelle z L die Randbedingung / 0 u z gefordert.

Diese Bedingung kann durch die Überlagerung des nach rechts laufenden Impulses

u kz t mit einem entgegengesetzt laufenden, virtuellen Impuls 2 u k z L t

eingehalten werden, denn es gilt:

0

0

0

2

1 cos 1 cos 22

sin sin 22

sin sin 02

z L z L

z L

z L

uu kz t u k z L t

z z

pkz t k z L t

k z

pk kz t k k z L t

k

pk kL t k kL t

k

157

Die Gesamtverschiebung der sich überlagernden Impulse ist bei z L :

0

0

, 1 cos 1 cos 22

2, 1 cos

2

z L

pu L t kz t k z L t

k

pu L t kL t

k

(6.11.13)

Die Amplitude verdoppelt sich also. Die Überlagerung ist im Folgenden graphischdurchgeführt:

Abbildung 91: Überlagerung von realem und virtuellem Verschiebungsimpuls. DieSpannungsbedingung bei 2m wird erfüllt, die Amplitude verdoppelt sich.

158

6.12. Das Spektrum des Kraftimpulses

Es soll das Spektrum eines zeitbegrenzten sinusförmigen Impulses gefunden wer-

den. Die Zeitbegrenzung wird durch ein Rechteck rect t realisiert:

1 / / 2 / 2 rect t T H t T H t T

t/ 2T 0 / 2T

Abbildung 90: Das Rechteck /rect t T der Breite T .

In diesem Zeitfenster g t kann ein sinusförmiger Impuls f t der festen Kreisfre-

quenz 0 /T und der Amplitude 0 1F formuliert werden:

0cos cos

t t tg t f t rect rect t

T T T. (6.12.1)

Das Spektrum S F g t f t wird nun unter Beachtung der Tatsache berech-

net, dass die Transformierte des Produktes zweier Funktionen im Zeitbereich derFaltung der zugehörigen Spektren im Frequenzbereich entspricht [37]:

F g t f t F g t F f t .

Die Einzelspektren können der Literatur entnommen werden [38]. Mit 0T ist:

0 0 0

2sin 2 2sinsin 2 2

2

1cos

2

T Tf

TftF rect T si Tf T

T Tf f

F t

(6.12.2)

Für die Faltung gilt das Distributivgesetz, durch dessen Anwendung man

0 0 0 0

sin sin sin2 2 2

T T T

S

erhält. Nun müssen die Faltungen berechnet werden, wobei ausgenutzt wird, dassdie Diracfunktion gerade ist:

159

x x .

Dann folgt mit der Abkürzung / 2 :T a :

0

0 0 0

0

0

0 0 0

0

sinsin sin sin

sinsin sin sin

aa a ad d

aa a ad d

Das Spektrum ist die Summe beider Faltungen:

0 0

0 0

sin sin

a aS .

Substitution der Abkürzung / 2a T ergibt die endgültige Form des Spektrums, denndie rechte Seite ist reell, und der Betrag muss nicht gebildet werden. Natürlich ist dasSpektrum noch mit der Kraftamplitude 0F aus (5.13) zu multiplizieren, die als reeller

Faktor in den obigen Rechnungen zu berücksichtigen ist:

0 0

0 0

0 0

sin sin2 2

T T

S FT

. (6.12.3)

Es ist zu beachten, dass

0 1si

ist, und dass das Spektrum eine Kraft pro Kreisfrequenz mit der Dimension N s

darstellt. [38] Integration des Spektrums über die Kreisfrequenz ergibt dann wiedereine Kraft.

160

6.13. 3d-Programmcode

166

6.14. 2d-Programmcode

172

7. Symbolverzeichnis (alphabetisch)

aussagenlogisches „und“, \ Mengenvereinigung, Mengendifferenz

u (fettes Zeichen) Vektor

A Abkürzung 2/ 2F R

,dA dA differentielles skalares, vektorielles Flächenelement

a KontaktradiusC Komponenten des Elastizitätstensors

D Komponenten des Rotationstensors

E Elastizitätsmodul

0F Kraftamplitude

f Volumenkraftdichte

,G G 3d-Gesamtgebiet, Gesamtgebiet des rz -Schnitts

H t Heaviside-Sprungfunktion

I Symmetrieachse des Körpers

, ,i j t Laufvariablen des Raumgittes ,i j und der diskreten Zeit t,zK Rand des rz -Schnitts, Flächennormalenvektor in j -Richtung

,yK Rand des xy -Gebietes, Normalenvektor in i -Richtung

L Länge des Stabesn Komponenten des Flächennormalenvektors

P Komponenten des Spannungsvektors

, ,r z Zylinderkoordinaten

, kR R Radius des Stabquerschnitts, Radius der Kugel

T Stoßzeit

ET Beobachtungszeit

1 2 3, ,u u u Komponenten der Verschiebung in indizierter Schreibweise

, ,r zu u u Komponenten der Verschiebung in Zylinderkoordinaten

kv Endgeschwindigkeit der Kugel vor dem Aufprall

1 2 3, , ;x x x t allg. Raumkoordinaten in indizierter Schreibweise; Zeit

, z Kroneckerdelta, Diracsche Deltadistribution

Komponenten des Verzerrungsstensors

, Lamésche Konstanten

0, lokale Dichte, konstante Dichte

σ Spannungstensor

Komponenten des Spannungstensors

nk Wellenzahl /n L der n -ten Fouriermode

Kreisfrequenz /T des Kraftimpulses

n Kreisfrequenz / nn c L k c der n -ten Fouriermode

173

8. Literatur

[1] Debnath, Lokenath: Nonlinear Partial Differential Equations,Birkhäuser, Boston, Basel, Berlin, 2. Aufl. 2005

[2] Schatt, Werner: Einführung in die Werkstoffwissenschaft,VEB Deutscher Verlag für Grundstoffindustrie, Leipzig 1972

[3] Strauss, Walter A.: Partielle Differentialgleichungen Vieweg 1995

[4] Schmutzer, Ernst: Grundlagen der Theoretischen Physik, Teil I,Bibliographisches Institut - Wissenschaftsverlag 1989

[5] Landau, Lev; Lifschitz, Evgenij: Lehrbuch der Theoretischen Physik, Elastizitäts-theorie,Akademie- Verlag Berlin 1989

[6] Gross, Dietmar; Hauger, Werner: Technische Mechanik 4,Springer - Verlag Berlin, Heidelberg, New York 2004

[7] Budo, Agoston: Theoretische Mechanik,Deutscher Verlag der Wissenschaften Berlin 1990

[8] Jung, Michael; Langer, Ulrich: Methode der finiten Elemente für Ingenieure,B.G.Teubner Stuttgart 2001

[9] Langenberg, K.J.: Modelle und Theorien für die Ultraschallprüfung 1990,DGZfP Berichtsband 23, Berlin, 5. und 6. November 1990

[10] mündliche Mitteilung Prof. Ch. Marotzke, Bundesanstalt für Materialforschung -und Prüfung Berlin

[11] Jonas, Piet: Instabilitäten und Strukturbildung in der positiven Säule von Neon-glimmentladungen, Diss. Univ. Greifswald 1998

[12] Fletcher, C.A.J.: Computational Techniques for Fluid Dynamiks, Volume 1,Springer - Verlag Berlin, Heidelberg, New York 1991

[13] Tveito, A.: Winther, Ragnar: Einführung in partielle Differentialgleichungen, Einnumerischer Zugang,Springer - Verlag Berlin, Heidelberg, New York 2002

[14] Schwarz, Hans Rudolf: Numerische Mathematik,B.G.Teubner Stuttgart 1997

[15] Bronstein; Semendjajew; Grosche; Ziegler: Teubner-Taschenbuch der Mathe-matik,B.G.Teubner Stuttgart 1996

174

[16] Marklein, René: Numerische Verfahren zur Modellierung von akustischen, elekt-romagnetischen, elastischen und piezoelektrischen Wellenausbreitungsproblemen imZeitbereich basierend auf der Finiten Integrationstechnik, D 34, (Diss. Univ. GH Kas-sel), Shaker Verlag, Aachen 1997

[17] Smith, G.D.: Numerical Solution of Partial Differential Equations: Finite Differ-ence Methods,Clarendon Press Oxford 1985

[18] Buttkus, Burkhard, Spektralanalyse und Filtertheorie in der angewandten Geo-physik,Springer-Verlag Berlin, Heidelberg, New York 1991

[19] Greve, Ralf: Kontinuumsmechanik,Springer - Verlag Berlin, Heidelberg, New York 2003

[20] Popper, Karl: Logik der Forschung,J.C.B. Mohr (Paul Siebeck) Tübingen 1989

[21] Polyanin, A.D.; Zaitsev, V. F.: Handbuch der linearen Differentialgleichungen,Spektrum Akademischer Verlag 1996

[22] Lang, Christian; Pucker, Norbert: Mathematische Methoden in der Physik,Spektrum Akademischer Verlag 1998

[23] Duffy, Dean: Green´s Functions with Applications,Chapmann & Hall/CRC 2001

[24] Preuss, W.; Kirchner, Hermann: Partielle Differentialgleichungen,VEB Fachbuchverlag Leipzig

[25] Tychonoff, A.; Samarski, A.: Differentialgleichungen der mathematischen Physik,VEB Deutscher Verlag der Wissenschaften Berlin 1959

[26] Wladimirow, W.S.: Gleichungen der mathematischen Physik,VEB Deutscher Verlag der Wissenschaften Berlin 1972

[27] Eringen, A. Cemal; Suhubi, Erdogan: Elastodynamics Volume II, Linear Theorie,Academic Press 1975

[28] Davis, Julian L.: Mathematics of Wave Propagation,Princeton University Press 2000

[29] Bronstein; Semendjajew; Musiol; Mühlig: Taschenbuch der Mathematik,Verlag Harri Deutsch 1995

[30] Trim, Donald W.: Applied Partial Differential Equations,PWS Publishing Company Boston 1989

[31] Egorov, Yu. V.; Shubin, M. A.: Partial Differential Equations I,Springer - Verlag Berlin, Heidelberg, New York 1991

175

[32] Lehner, Günther: Elektromagnetische Feldtheorie für Ingenieure und Physiker,Springer - Verlag Berlin, Heidelberg, New York 1994

[33] Greiner, Walter: Klassische Elektrodynamik,Verlag Harri Deutsch 1991

[34] Gradshteyn, L.S.; Ryznik, I.M.: Table of Integrals, Series and Products,Academic Press, 6. Aufl. 2000

[35] Achenbach, J.D.: Wave Propagation In Elastic Solids,Verlag North-Holland bzw. ©Elsevier Science Publishers 1984

[36] Brekhovskikh, Leonid; Goncharov, Valery: Mechanics of Continua and WaveDynamics,Springer-Verlag Berlin, Heidelberg, New York 1985

[37] Kories, Ralf; Schmidt-Walter, Heinz: Taschenbuch der Elektrotechnik,Verlag Harri Deutsch 2004

[38] Ohm, Jens-Rainer; Lüke, Hans Dieter: Signalübertragung,Springer-Verlag Berlin, Heidelberg, New York 2002

176

Danksagung

An erster Stelle möchte ich Prof. Dr. Chr. Marotzke, Bundesanstalt für Materialfor-

schung- und prüfung Berlin (BAM) und Prof. Dr. W. H. Müller, Technische Universität

Berlin, für ihren fachkundigen Rat, für die außerordentliche Gesprächsbereitschaft

und für die angenehme Betreuung der Arbeit danken.

Ebenso danke ich Dr. J. Bohse (BAM) für die Unterstützung in experimentellen Fra-

gen, und Prof. Dr. A. Hampe (BAM) für die Möglichkeit, die Arbeit ohne Zeitdruck zu

beenden. Herzlichen Dank gebührt Prof. Dr. G. Czichowski, Fachbereich Mathematik

der Universität Greifswald, für wertvolle Hinweise.

Ich danke meinen Eltern für die unbeschwerten Jahre meiner Jugend.

177

Erklärung

Hiermit erkläre ich, dass diese Arbeit bisher von mir weder der Mathematisch -

Naturwissenschaftlichen Fakultät der Technischen Universität Berlin noch einer

anderen wissenschaftlichen Einrichtung zum Zwecke der Promotion eingereicht

wurde.

Ferner erkläre ich, dass ich diese Arbeit selbständig verfasst und keine anderen als

die darin angegebenen Hilfsmittel benutzt habe.

Thomas-Rex Heinrich

178

Lebenslauf

Persönliche Angaben

Staatsangehörigkeit: Deutschland

Geburtsdatum: 14. März 1962

Geburtsort: Greifswald

Ausbildung und Lebensdaten

1968 - 1978Allgemeinbildende polytechnische Oberschule in Greifswald

Abschluss der 10. Klasse

1978 - 1986Technischer Assistent im Wetterdienst

1987 - 1989Militärdienst unter Verweigerung der Waffe

1989 - 1992Arbeit im sozialen Bereich; Besuch der Abiturstufe der Volkshochschule Greifs-wald

1992 Abschluss des Abiturs auf dem zweiten Bildungsweg

1992 - 1999Studium der Physik mit Spezialisierungsrichtung theoretische Physik

Abschluss als Diplomphysiker

2001 - 2005Doktorand an der Bundesanstalt für Materialforschung- und prüfung (BAM)

2006Wissenschaftlicher Mitarbeiter am Institut für Physik der Universität Greifswald

Die Lösung der elastischen Bewegungsgleichung als ... · Kraftstoß in Achsenrichtung eines Stabes...

Documents

Transcript of Die Lösung der elastischen Bewegungsgleichung als ... · Kraftstoß in Achsenrichtung eines Stabes...