next up previous contents
Next: Feldberechnung mit FDTD auf Up: Zentrale Systeme Previous: Feldberechnung mit CONCEPT auf

Simulation polykristalliner Festkörper auf der Hydra

Erläuterung zum wissenschaftlichen Thema:

Bei dem zu parallelisierenden Programm handelt es sich um ein Simulationsprogramm, welches die Spannungen, Dehnungen und kristallographischen Rotationen in einem polykristallinen Festkörper beschreibt, der einer plastischen Verformung unterliegt.

Die einfachste Theorie zu diesem Thema ist die Theorie nach Taylor, Bishop und Hill (vereinfachte Darstellung):

  1. Gegeben ist ein makroskopischer Deformationstensor, der für alle Körner in dem Polykristall als gleich angenommen wird (damit ist die Kontinuität des Körpers gewährleistet).

    displaymath1556

  2. Weiterhin sei die Fließfläche für das betrachtete Material in einkristalliner Form bekannt. Eine Fließfläche ist der geometrische Ort der Punkte im Spannungsraum, die die an der Fließgrenze liegenden Spannungszustände (Fließspannung S) repräsentieren.
  3. Es gilt das Prinzip des ``maximum plastic work'', d.h. zu einem vorgegebenen Deformationstensor D stellt sich die Fließspannung S stets so ein, daß gilt:

    displaymath1564

Mit Hilfe dieser Theorie ist es möglich, die Spannungsverteilung und (sofern die mikroskopischen Verschiebungstensoren bekannt sind) die plastischen Gitterrotationen in einem Polykristall zu berechnen. Für die Orientierungsverteilungen, die sich aus einer solchen Simulation ergeben, hat man in einigen Fällen gute qualitative Übereinstimmung zwischen Theorie und Experiment gefunden. Für viele technisch wichtige Fälle, wie z.B. die Verformung eines zweiphasigen Polykristalls, muß diese Theorie jedoch versagen, da sie auf einer etwas gewaltsamen Annahme (D=const.) basiert, was in einem Material, bestehend aus einer harten und einer weichen Phase, sicher nicht erfüllt ist. Daher gibt es schon seit längerem andere Modelle, die eine Berechnung von Spannungen und Dehnungen mit realistischeren Annahmen gestatten. In diesen Modellen ergibt sich der makroskopische Deformationstensor durch eine Mittelung über die Körner des Polykristalls:

  equation322

Entsprechendes gilt auch für den Spannungstensor:

  equation327

Die Deformationstensoren tex2html_wrap_inline1566 können also verschieden sein. In dem sogenannten selbstkonsistenten Modell nach Molinari, Canova et al. sind dann die Spannungen und Dehnungen durch folgende Gleichung gegeben:

  equation332

Der Tensor M (Eshelby-Tensor) beschreibt die Wechselwirkung eines Korns mit einer als gemittelt gedachten Umgebung (homogenes, äqivalentes Medium) und ergibt sich in einer hier nicht näher erläuterten Weise aus dem plastischen Sekantenmodul L:

equation338

equation343

L ergibt sich aus eingegebenen, materialspezifischen Daten (Gleitsysteme, Verfestigung, Geschwindigkeitsempfindlichkeit). Gleichung (3) läßt sich durch eine iterative Prozedur lösen, d.h. die Größen in (3) müssen so gewählt werden, daß S und D aus (1), (2) erfüllt sind und konsistent mit dem Tensor M zusammenpassen.

Dabei kann der Rechenaufwand für den Tensor M, der in einer fortgeschrittenen Version des Modells ein komplizierter Tensor 4. Stufe ist, sehr groß werden.

Bei der Berechnung von M lassen sich Teile des Programms parallelisieren, d. h. es werden Prozeduren für jedes einzelne Korn durchgeführt, die von den Ergebnissen der übrigen Körner nicht abhängen. Die Parallelisierung lohnt sich dabei erst, wenn der Rechenaufwand für ein einzelnes Korn sehr groß wird, was in einfacheren Versionen des beschriebenen Modells noch nicht der Fall ist.

Erläuterung zum Programm:

Das Programm ist in Fortran77 geschrieben und wurde von Ricardo Lebensohn am Institut für Physik (Rosario, Argentinien) entwickelt. Eine Analyse des Laufverhaltens mit cxpa sowie der Programmstrukur ergab, daß fast die gesamte Rechenzeit in einer einzigen DO-Schleife verbraucht wird, die für alle Körner (zur Zeit sind es 800) wiederholt wird. Dabei empfiehlt es sich, daß die Schleife in einer separaten Subroutine abgearbeitet wird. Die benötigten Variablen (Felder, Konstanten usw.) müssen dabei in der Parameterliste (und nicht in einem COMMON-Block) übergeben werden. Die Schleife selbst (DO-Label 2096) wird mit einer Serie von einfachen Compilerdirektiven (C$DIR ...) parallelisiert:

  
SUBROUTINE LOOP2096(......)
          .
          .
          .
C$DIR loop_parallel
C$DIR loop_private(DUMMY,c4gr,r0gr,ax1,ax2,dummye,dummyr,egr,rgr,e)
C$DIR loop_private(e11,e22,e12,e21,xmcbc1,xmcbc2,bc1,bc2)
C$DIR loop_private(i,j,m,n,i1,j1,m1,n1,ll,k1,iq,ip)

      DO 2096 K=1,NGR    !DO-Schleife über alle Körner (NGR maximal 800)
          .
          .
      DO 7729 I=1,3
      DO 7729 J=1,3
      DO 7729 M=1,3
      DO 7729 N=1,3
          .
          .
7729  CONTINUE
          .
          .
      CALL INTER(AX1,AX2,XJAC,R0GR,C4GR,EGR,RGR,0)
          .
          .
2096  CONTINUE
          .
          .
      RETURN
END

Dabei muß beachtet werden, daß hydra mit einem shared memory arbeitet, d. h. Zahlen, die für alle Körner in einem seriellen Programm gemeinsam addressiert sind, müssen durch eine loop_private-Direktive separaten Adressen zugewiesen werden. So braucht z.B. ein Feld der Form ARRAY(I,J,K) nicht in diese Liste aufgenommen werden, da Zahlen, die zu verschiedenen Indizes K gehören, auch getrennt adressiert sind. Dagegen müssen aber z.B. alle inneren DO-Schleifen-Parameter unter loop_private aufgelistet sein.

Die parallelisierte Routine wurde mit folgenden Optionen compiliert:

OPTIONS= -O3 -cxpa -noautopar -mrl

Die übrigen Routinen wurden nur mit der Option -cxpa compiliert.

Nach Beachtung dieser Hinweise lief das Programm einwandfrei mit einer erheblichen Zeitersparnis gegenüber der seriellen Version. Zwei Testläufe mit cxpa mit jeweils sechs CPUs ergaben folgendes Resultat:

Zahl der Körner CPU-Time Benutzerzeit Effizienz der
Verformungsschritte Parallelisierung
24 /2 481 sec 111 sec 4.3
100 /2 1720 sec 286 sec 6.0
800 /20 ca. 32 h 5:20 h ca. 6

Das bedeutet, bei entsprechend vielen Körnern ergab sich die maximal mögliche Zeitersparnis von einem Faktor 6 gegenüber der seriellen Version.

Hierbei sollte auch beachtet werden, daß die Eingabedaten nur für einen relativ einfachen Fall (nur kubische Kristallstruktur, kleine Verfestigung) gelten. Kompliziertere Fälle würden sich mit vertretbarem Zeitaufwand nur noch auf einem Parallelrechner rechnen lassen.

Ergebnis:

Für den simulierten Verbundwerkstoff ergeben sich die in Fig. 5,6 dargestellten Ergebnisse:

  figure365
Abbildung 4:   Spannungs-Dehnungskurven eines zweiphasigen Faser-Verbundwerkstoffs.

     figure369
Abbildung 5: (111)-Polfiguren der Faser
                                  
Abbildung 6: 111-Polfiguren der Matrix

Wie zu erwarten, übernimmt die harte Phase größere Spannungen und kleinere Dehnungen als der Verbund. Die weiche Phase verhält sich umgekehrt. Eine genaue Analyse des Spannungs-Dehnungsverhaltens ergibt außerdem für die weichere Phase einen durch die Texturentwicklung bedingten nichtlinearen Anstieg der Verfestigung. Die Texturen der beiden Phasen lassen ebenfalls (trotz gleicher Kristallstruktur) eine unterschiedliche Entwicklung der Vorzugslagen erkennen.

Christian Hartig (Arbeitsbereich Werkstoffphysik), Tel.: 2276


next up previous contents
Next: Feldberechnung mit FDTD auf Up: Zentrale Systeme Previous: Feldberechnung mit CONCEPT auf

Marco Budde
Fri Apr 12 18:53:17 MESZ 1996