DE102019117558B4 - Verfahren und Computerprogramm zum Bestimmen einer Qualität einer Maske einer Fotolithographie-Vorrichtung - Google Patents

Abstract

Verfahren (1000) zum Bestimmen einer Qualität einer Maske einer Fotolithographie-Vorrichtung, umfassend:
Empfangen eines Designs (1010) der Maske,
Empfangen von optischen Eigenschaften (1020) der Fotolithographie-Vorrichtung,
paralleles Berechnen (1030), mit einer Vielzahl von Recheneinrichtungen, eines Referenz-Luftbildes basierend auf dem Design und den optischen Eigenschaften, wobei das parallele Berechnen (1030) des Referenz-Luftbildes basierend auf dem Pupil-Shift-Verfahren erfolgt, wobei das Pupil-Shift-Verfahren für eine parallele Berechnung mit der Vielzahl von Recheneinrichtungen implementiert ist,
Empfangen (1040) eines gemessen Luftbildes (400),
Vergleichen (1050) des Referenz-Luftbildes (200) mit dem gemessenen Luftbild (400),
Bestimmen (1060) einer Qualität der Maske (500) basierend auf dem Vergleichen.

Description

• TECHNISCHES GEBIET
• Verschiedene Ausführungsbeispiele betreffen das Ermitteln der Qualität von fotolithographischen Einrichtungen oder Komponenten dieser Systeme mittels parallelen Berechnens.
• HINTERGRUND
• Zur Bewertung von fotolithographischen Vorrichtungen, beispielsweise von Masken, die in solchen Vorrichtungen Verwendung finden, werden Bildgebungssimulationen von sogenannten Luftbildern durchgeführt. Solche Masken werden manchmal als Fotomasken bezeichnet.
• Bildgebungssimulationen solcher Luftbilder sind beim Herstellungsprozess moderner Fotomasken unabdingbar, da sie eine zuverlässige Vergleichsbasis für das reale Fabrikat bildet und es somit ermöglicht, Herstellungsfehler zu detektieren.
• Als Luftbild wird allgemein die Intensitätsverteilung in der Bildebene eines optischen Systems verstanden.
• Luftbilder sind beispielsweise in A. Garetto, Th. Scherübl and J. H. Peters, „Aerial imaging technology for photomask qualification: from a microscope to a metrology tool“, Adv. Opt. Techn., 2012; 1(4): 289-298, beschrieben.
• Luftbilder können in optischen Systemen, die eine Vielzahl von Bildebenen aufweisen, entsprechend definiert sein. Hierbei kann ein Luftbild Informationen für die Vielzahl von Bildebenen aufweisen und/oder die Intensitätsverteilung in der Vielzahl der Bildebenen kann durch eine Vielzahl von Luftbildern beschrieben werden.
• Um die Rechenanforderungen solcher Berechnungen zu verringern, sind Verfahren zur partiell kohärenten Bildgebungssimulation (Englisch „partially coherent image simulation“) bekannt und finden Anwendung für die Simulation von fotolithographischen Vorrichtungen, die eine Köhlersche Beleuchtung aufweisen. Solche Vorgehensweisen erlauben die Simulation von optischen Bildern einer virtuellen idealisierten Maske.
• Verschiedene Vereinfachungen, um die Durchführung solcher Simulationen zeiteffizienter zu gestalten und den benötigten Rechenaufwand zu verringern, sind für sich genommen bekannt.
• Aus Kenji Yamazoe, „Computation theory of partially coherent imaging by stacked pupil shift matrix“, J. Opt. Soc. Am. A, 25(12): 3111-3119, Dezember 2008, und Kenji Yamazoe, „Two models for partially coherent imaging“, J. Opt. Soc. Am. A, 29(12): 2591-2597, Dezember 2012, sind mathematische Methoden zur Berechnung von Luftbildern von Fotolithographie-Vorrichtungen mit Köhlerscher Beleuchtung mittels partiell kohärenter Bildgebungssimulation und mittels des Pupil-Shift-Verfahrens bekannt.
• Aus H. Kirchauer, Photolithography Simulation, Technische Universität Wien: Dissertationen der Technischen Universität Wien. Österreichischer Kunst- und Kulturverl., 2000, sind mathematische Verfahren bekannt, um idealisierte Masken mit Hilfe von Dreiecken für partiell kohärente Bildgebungssimulationen zu beschreiben.
• Die DE 10 2013 101 445 A1 offenbart ein Verfahren zum Bestimmen einer Maske, wobei ein Referenzluftbild basierend auf einem Maskendesign und optischen Eigenschaften einer Fotolithographie-Vorrichtung berechnet wird. Das Referenzluftbild wird mit einem gemessenen Luftbild verglichen, und basierend auf dem Vergleich wird eine Qualität einer Maske bewertet.
• Ein derartiges Verfahren wird auch in der nachveröffentlichten DE 10 2018 207 880 A1 als bekannt beschrieben.
• Die US 7 027 143 B1 offenbart eine Vorrichtung zur Maskencharakterisierung mit einer Messvorrichtung zum Erfassen eines Luftbildes einer Maske und eine Rechenvorrichtung, die parallele Prozessoren aufweisen kann.
• Trotz dieser für sich genommen bekannten Vorgehensweisen und Vereinfachungen sind solche optischen Bildgebungssimulationen sehr rechenintensiv und in bekannten Implementierungen auch mit schnellen Rechnern daher zeitaufwändig bei der Berechnung. In vielen Anwendungsfällen wird eine Vielzahl von Simulationen für verschiedene Bedingungen, beispielsweise verschiedene Einstellungen der fotolithographischen Einrichtung, durchgeführt. Dies kann mit erheblich mehr Rechenzeitbedarf einhergehen.
• Ein weiterer Nachteil der bekannten theoretischen Verfahren ist, dass basierend auf theoretischen Annahmen simulierte Bilder von idealisierten theoretischen Masken erhalten werden, eine Interpretation der Simulationsergebnisse im Hinblick auf eine real existierende Fotolithographie-Vorrichtung mithin schwierig sein kann.
• Deshalb besteht ein Bedarf für Verfahren und Vorrichtungen, die ein Bestimmen einer Qualität von Masken von Fotolithographie-Vorrichtungen verbessern.
• KURZBESCHREIBUNG
• Es wird ein Verfahren wie in dem unabhängigen Patentanspruch definiert bereitgestellt. Die Merkmale der abhängigen Patentansprüche definieren weitere Ausführungsbeispiele.
• Ein Verfahren zum Bestimmen einer Qualität einer Maske einer Fotolithographie-Vorrichtung umfasst:
• Empfangen eines Designs der Maske und Empfangen von optischen Eigenschaften der Fotolithographie-Vorrichtung.
• Das Verfahren umfasst ferner paralleles Berechnen, mit einer Vielzahl von Recheneinrichtungen, eines Referenz-Luftbildes basierend auf dem Design und den optischen Eigenschaften auf einer Vielzahl von Recheneinrichtungen sowie Empfangen eines gemessenen Luftbildes.
• Ferner umfasst das Verfahren: Vergleichen des Referenz-Luftbildes mit dem gemessenen Luftbild und Bestimmen einer Qualität der Maske basierend auf dem Vergleichen.
• Die Fotolithographie-Vorrichtung kann hierbei ein System mit einer Köhlerschen Beleuchtung sein.
• Die Messvorrichtung kann einen Sensor umfassen, der in der Bildebene der Fotolithographie-Vorrichtung platziert ist und eingerichtet ist, eine räumliche Intensitätsverteilung zu messen.
• Die Vielzahl von Recheneinrichtungen können eine Vielzahl von Prozessoren und/oder eine Vielzahl von Parallelprozessoreinheiten und/oder eine programmierbare Hardwarekomponente umfassen. Hierbei können die Vielzahl von Recheneinrichtungen verschiedene Speicher umfassen. Beispielsweise können Parallelprozessoreinheiten gemeinsam Zugriff auf einen gemeinsamen Speicher (englisch: shared memory) haben und ferner Zugriff auf einen der jeweiligen Parallelprozessoreinheit zugehörigen Speicher haben. Hierbei kann der Zugriff auf den zugehörigen Speicher schneller sein als auf den gemeinsamen Speicher. Parallelprozessoreinheiten können in manchen Ausführungsbeispielen als Grafikprozessoreinheiten (englisch: Graphics Processing Unit, GPU) realisiert sein. Die Grafikprozessoreinheiten können in einer Rechenvorrichtung in Kommunikation mit einem oder mehreren Zentralprozessoreinheiten (englisch: Central Processing Unit, CPU) stehen.
• Die Vielzahl von Recheneinrichtungen kann hierbei eine Parallelberechnungskarte umfassen oder sein.
• Die Parallelberechnungskarte kann beispielsweise CUDA fähig sein, wie auf der Seite „CUDA“. In: Wikipedia, Die freie Enzyklopädie. Bearbeitungsstand: 7. Mai 2019, 19:50 UTC. URL: https://de.wikipedia.org/w/index.php?title=CUDA&oldid=188338305 (Abgerufen: 13. Mai 2019, 12:32 UTC) beschrieben. Beispielsweise kann die Version CUDA 9.1 verwendet werden.
• Unter Referenz-Luftbild wird ein simuliertes Luftbild (englisch. aerial image) verstanden. Als Luftbild wird allgemein die Intensitätsverteilung in der Bildebene eines optischen Systems verstanden. Solche Referenz-Luftbilder simulieren das Verhalten von fotolithographischen Vorrichtungen und sind beispielsweise in A. Garetto, Th. Scherübl and J. H. Peters, „Aerial imaging technology for photomask qualification: from a microscope to a metrology tool“, Adv. Opt. Techn., 2012; 1(4): 289-298, bekannt.
• Das gemessene Luftbild umfasst eine räumliche Intensitätsverteilung in einer Bildebene der Fotolithographie-Vorrichtung. Das gemessene Luftbild kann beispielsweise mit einem Sensor und optional weiteren optischen Komponenten direkt in der Fotolithographie-Vorrichtung bestimmt werden. Es kann alternativ oder ergänzend in einer dafür speziell eingerichteten Vorrichtung, beispielsweise in einer Fotomaskenqualifikationsvorrichtung wie dem Aerial Image Measurement System, AIMS, von Zeiss, bestimmt werden.
• Mit anderen Worten kann das gemessene Luftbild mit einer Messvorrichtung in der Fotolithographie-Vorrichtung und/oder in einer externen Vorrichtung wie einer Fotomaskenqualifikationsvorrichtung gemessen werden.
• Hierbei kann in den Fällen, in denen die Messung in der Fotolithographie-Vorrichtung und in der Fotomaskenqualifikationsvorrichtung erfolgt, entweder mit einer Messvorrichtung gearbeitet werden, oder es können Informationen verschiedener Messvorrichtungen kombiniert werden. In solchen Fällen kann die Messvorrichtung eine erste Messvorrichtung und eine zweite Messvorrichtung umfassen.
• Das Verfahren kann den Vorteil haben, dass durch die Vielzahl von Recheneinrichtungen die Bearbeitung beschleunigt wird.
• Das parallele Berechnen des Referenz-Luftbildes erfolgt basierend auf dem Pupil-Shift-Verfahren, wobei das Pupil-Shift-Verfahren für eine parallele Berechnung mit der Vielzahl von Recheneinrichtungen implementiert sein kann.
• Das Pupil-Shift-Verfahren wird manchmal auch als Verfahren nach Hopkins bezeichnet und ist für sich genommen beispielsweise aus Kenji Yamazoe: „Computation theory of partially coherent imaging by stacked pupil shift matrix“, J. Opt. Soc. Am. A, 25(12): 3111-3119, Dezember 2008, Kenji Yamazoe: „Two models for partially coherent imaging“, J. Opt. Soc. Am. A, 29(12): 2591-2597, Dezember 2012, sowie H. Kirchauer, „Photolithography Simulation“, Technische Universität Wien: Dissertationen der Technischen Universität Wien, Österreichischer Kunst- und Kulturverl., 2000, bekannt. Das Pupil-Shift-Verfahren bzw. die eigentliche Verschiebung der Koordinaten im Doppelintegral wird manchmal auch als zirkuläre Verschiebung (englisch: circular shift) bezeichnet.
• Das parallele Berechnen kann hierbei umfassen, dass das zu berechnende Referenz-Luftbild in Blöcke unterteilt wird, wobei die einzelnen Blöcke jeweils parallel berechnet werden.
• Das Pupil-Shift-Verfahren kann hierbei eine gemeinsame Intensitätsfunktion, eine Maskentransmissionsfunktion, und eine diskrete Pupillenfunktion verwenden oder äquivalente oder reduzierte Repräsentationen dieser Funktionen verwenden wie nachfolgend am Beispiel der Transmissions-Kreuzkoeffizientenmatrix erläutert.
• Die gemeinsame Intensitätsfunktion beschreibt die Beleuchtung des Objekts.
• Die gemeinsame Intensitätsfunktion wird manchmal mit „J“ bezeichnet. Die gemeinsame Intensitätsfunktion wird manchmal auch als Intensitätsfunktion oder als gemeinsame Intensitätsverteilungsfunktion bezeichnet.
• Die diskrete Pupillenfunktion beschreibt eine komplexe skalare Amplitudenverteilungsfunktion der Bildgebung der Fotolithographie-Vorrichtung. Die diskrete Pupillenfunktion wird manchmal mit „K“ bezeichnet. Manchmal wird die diskrete Pupillenfunktion auch als Pupillenfunktion bezeichnet.
• Die Maskentransmissionsfunktion ist ein zweidimensionales komplexes Feld aus welchem mittels Fourier-Transformation eine ortsabhängige Transmission der Maske bestimmt werden kann, wobei die Transmission der Maske für eine zweidimensionale Maske angegeben ist. Die lokale Transmission der Maske kann als reale Transmissionscharakteristik und/oder als komplexe Transmissionscharakteristik bei phasenverändernden Masken angegeben werden. Die Transmission der Maske wird manchmal mit „t“, die Maskentransmissionsfunktion wird manchmal mit „T“ bezeichnet. Die Maskentransmissionsfunktion kann durch Maskentransmissionskoeffizienten t_k beschrieben werden, wie nachfolgend näher erläutert.
• Hierbei kann die diskrete Pupillenfunktion basierend auf mindestens einer Messung der Messvorrichtung bestimmt werden.
• Beispielsweise kann die diskrete Pupillenfunktion in einer oder mehreren Bildebenen über eine Kalibriermessung bestimmt werden. Die diskrete Pupillenfunktion kann hierbei unabhängig von einer Maske und/oder einer Maskentransmissionsfunktion bestimmt werden.
• Die Maskentransmissionskoeffizienten t_k können ebenfalls über eine unabhängige Kalibriermessung für die betreffende Maske bestimmt werden.
• Hierdurch kann die Interpretation des Referenz-Luftbildes vereinfacht werden, da weniger Annahmen über die Eingangsparameter der Simulation getroffen werden müssen.
• Die diskrete Pupillenfunktion und die gemeinsame Intensitätsfunktion können als eine Transmissions-Kreuzkoeffizientenmatrix, TCC-Matrix, dargestellt werden.
• Die TCC-Matrix kann hierbei nach Gleichung (2.7) bestimmt werden, wie nachfolgend beschrieben. J kann in manchen Ausführungsbeispielen als konstant und/oder bekannt angenommen werden.
• Dies kann den Vorteil haben, dass der Speicherbedarf der Vielzahl von Recheneinrichtungen insgesamt oder der zugehörigen Speicher verringert werden kann oder die Anzahl an Zugriffen auf den Speicher oder den zugehörigen Speicher verringert werden kann.
• Das Verfahren kann ferner umfassen:
• Bestimmen einer Vielzahl von Eigenwerten der Transmissions-Kreuzkoeffizientenmatrix,
• Auswählen einer Menge von Eigenwerten aus der Vielzahl von Eigenwerten,
• wobei die Anzahl der Menge von Eigenwerten kleiner ist als die Anzahl der Vielzahl von Eigenwerten,
• Durchführen des Verfahrens mit der Menge der Eigenwerte.
• Die Anzahl der Menge von Eigenwerten kann um einen Faktor a kleiner sein als eine Anzahl von verwendeten Lichtquellen-Punkten.
• Lichtquellen-Punkte sind hierbei punktförmige Lichtquellen in der zweidimensionalen Beleuchtungsebene. Diese können mit einer bestimmten Intensität und unter einem bestimmten Winkel auf die Maske leuchten. Der Winkel kann abhängig von der Position des Punktes auf der Pupillenebene sein. Lichtquellen-Punkte können auch als eine Punktlichtquelle ausgeführt sein, also eine von Winkeln unabhängige Intensitätsverteilung aufweisen.
• Lichtquellen-Punkte werden manchmal mit dem Symbol „L“ bezeichnet. Hierbei kann L die Koordinaten und/oder die Intensitätsverteilung von Lichtquellen-Punkten, beispielsweise dargestellt als Einstrahlwinkel, umfassen. Das Licht mehrerer solcher Lichtquellen-Punkte kann durch die gemeinsame Intensitätsfunktion, wie zuvor beschrieben, manchmal mit dem Symbol „J“ bezeichnet, beschrieben werden.
• Hierbei kann a > 5, beispielsweise a > 10, beispielsweise a > 15, beispielsweise a > 20, beispielsweise a > 25, a > 30, sein.
• Alternativ oder ergänzend können zur Beschreibung einer Maske zwischen 50 und 1000 Eigenwerte, insbesondere zwischen 100 und 500 Eigenwerte, verwendet werden.
• Solche Vorgehensweisen, um die Anzahl der verwendeten Eigenwerte zu verringern, können den Vorteil haben, dass der Speicherbedarf der Vielzahl von Recheneinrichtungen insgesamt oder der zugehörigen Speicher verringert werden kann und/oder die Anzahl an Zugriffen auf den Speicher oder den zugehörigen Speicher verringert werden kann und/oder die Berechnungsdauer verkürzt werden kann.
• Das Verfahren kann umfassen:
• Darstellen des Designs der Maske als eine Menge von k Polygonen,
• Durchführen einer analytischen Spektrumsberechnung für die Menge von k Polygonen,
• Bestimmen der Maskentransmissionsfunktion basierend auf der analytischen Spektrumsberechnung für Polygone mittels Fourier-Transformation.
• Hierbei können die Polygone beliebig gewählt werden.
• Das Design der Maske kann hierbei eine Maskentransmission im Ortsraum sein, manchmal als t(x, y) bezeichnet.
• Unter analytischer Spektrumsberechnung wird hierbei verstanden, basierend auf der Maskentransmission im Ortsraum durch Fourier-Transformation die Maskentransmissionsfunktion T zu bestimmen. Die Maskentransmissionsfunktion ist hierbei im Frequenzraum definiert und wird daher manchmal auch als Maskenspektrum oder Maskenspektralfunktion bezeichnet. Wird die Maskentransmission durch Dreiecke beschrieben, kann die Fourier-Transformation für einzelne Dreiecke durchgeführt werden. Dies wird später näher erläutert.
• Die Polygone können Dreiecke sein. In solchen Fällen kann das Verfahren umfassen:
• Darstellen des Designs der Maske als eine Menge von k Dreiecken,
• Durchführen einer analytischen Spektrumsberechnung für die Menge von k Dreiecken,
• Bestimmen der Maskentransmissionsfunktion basierend auf der analytischen Spektrumsberechnung mittels Fourier-Transformation.
• Eine mögliche Vorgehensweise zur analytischen Spektrumsberechnung von Dreiecken ist für sich genommen aus H. Kirchauer, Photolithography Simulation, Technische Universität Wien: Dissertationen der Technischen Universität Wien. Österreichischer Kunst- und Kulturverl., 2000, bekannt, vgl. B.1. Die Gleichungen weisen hierbei Fallunterscheidungen auf. Dies kann in Implementierungen dieser Berechnungen sogenannte Verzweigungen (englisch: branching) erforderlich machen.
• Üblicherweise werden solche Verzweigungen in Gleichungen für parallelisierte Berechnungsverfahren als nachteilhaft angesehen. Eine solche Vorgehensweise zu verwenden, kann jedoch den Vorteil haben, dass bei manchen Implementierungen, beispielsweise mittels CUDA-Verfahren, ein solches Verzweigen kein flaschenhalsartiges Verhalten (englisch: bottle neck) aufweist. Dies wird später näher erläutert.
• Nachfolgend wird ein allgemeines Verfahren für Polygone beschrieben. Dieses Verfahren kann ergänzend oder alternativ zu dem zuvor beschriebenen Verfahren für Dreiecke angewandt werden. Es kann auch für Polygone mit einer höheren Eckenanzahl eingesetzt werden. Das betrachtete Polygon kann beispielsweise ein Tetragon, Pentagon, Hexagon, Heptagon oder ein Oktogon sein. Aber auch höhere Eckenzahlen sind möglich.
• Für das Verfahren werden die folgenden Bezeichnungen gewählt:
• Der Index k = 0, 1, 2..., N-1 (zyklisch mod N) bezeichnet, sofern nichts anderes angegeben ist, den Index der Ecken des Polygons.
• Zum Beispiel läuft k von 0 bis 2 Im Falle eines Dreiecks und von 0 bis 4 im Falle eines Pentagons.
• „Zyklisch mod N“ bedeutet hier beispielsweise p_N=p0, p_N+1=p1 usw.
• Die jeweiligen Eckpunkte werden mit dem Vektor p_k bezeichnet, wobei k der zuvor beschriebene Index ist.
• Die Kanten können gemäß [p_k, p_k+1] = {(1-t) p_k + t p_k+1: 0 ≤ t ≤ 1} parametrisiert werden.
• Mit $m_k=\frac{p_{k+1}+p_k}{2}$

  

  können die Vektoren der Kantenmittelpunkte bezeichnet werden. $t_k=\frac{p_{k+1}-p_k}{2}$

  

  bezeichnet die halben Kantenvektoren.
  n̂_k bezeichnet die nach außen gerichteten Einheitsnormalenvektoren der Kanten.
• Die (skalierten) Normalenvektoren der Kanten können gegeben sein durch: $${\text{n}}_k=\left|{\text{p}}_{k+1}-{\text{p}}_k\right|{\widehat{n}}_k$$

  
• Mit den zuvor genannten Konventionen kann das Spektrum des Polygons wie folgt analytisch berechnet werden:
• Der Gradient der Indikatorfunktion $1p$
  
  von P ist gegeben durch
$$\nabla 1p={\displaystyle \sum _k{\widehat{n}}_k}{\delta }_{\left[p_k,p_{k+1}\right]}$$

und hat die Fourier-Transformation $$i\omega \widehat{1}p\left(\omega \right)={\displaystyle \sum _k{\text{n}}_k}{\displaystyle \underset{\left[p_k,p_{k+1}\right]}{\int }{e}^{-i\omega r}\text{dr}.}$$

• Dann ergibt sich für das Integral über die Kante: $$\begin{array}{c}\frac{1}{\left|p_{k+1}-p_k\right|}{\displaystyle \underset{\left[p_k,p_{k+1}\right]}{\int }{e}^{-i\omega r}\text{dr}}={\displaystyle \underset{0}{\overset{1}{\int }}{e}^{-i\omega \left(\left(1-t\right)p_k+t{p}_{k+1}\right)}\text{d}t=}{\left[\frac{e^{-i\omega \left(\left(1-t\right)p_k+i{p}_{k+1}\right)}}{-i\omega \left(p_{k+1}-p_k\right)}\right]}_0^1\\ =\frac{e^{-i\omega p_{k+1}}-e^{-i\omega p_k}}{-i\omega \left(p_{k+1}-p_k\right)}\end{array}$$

  
• Somit gilt: $${\omega }^2\widehat{1}p\left(\omega \right)={\sum _k}_{\omega t_k}^{\omega n_k}\left(e^{-\omega p_{k+1}}-e^{-i\omega p_k}\right).$$

  
• Bezeichnet man den Winkel zwischen der k-ten Kante und dem Frequenzvektor mittels θ_k,so dass gilt: $$\text{sin}{\theta }_k=\widehat{\omega }{\widehat{n}}_k$$

  

  $$\text{cos}{\theta }_k=\widehat{\omega }{\widehat{t}}_k,$$

  

  dann gilt: $$\widehat{1}p\left(\omega \right)=\frac{1}{{\omega }^2}{\displaystyle \sum _k\text{tan}{\theta }_k}\left(e^{-\omega p_{k+1}}-e^{-i\omega p_k}\right)$$

  

  $$=\frac{1}{{\omega }^2}{\displaystyle \sum _k\left(\text{tan}{\theta }_{k-1}-\text{tan}{\theta }_k\right)}{e}^{-i\omega p_k}.$$

  
• Das Kantenintegral kann auch geschrieben werden als: $$\frac{1}{\left|p_{k+1}-p_k\right|}{\displaystyle \underset{\left[p_k,p_{k+1}\right]}{\int }{e}^{-i\omega r}\text{dr}}=-e^{-i\omega \frac{p_{k+1}+p_k}{2}}\frac{e^{+iw\frac{p_{k+1}-p_k}{2}}-e^{-i\omega \frac{p_{k+1}-p_k}{2}}}{2i\omega \frac{p_{k+1}-p_k}{2}},$$

  
• Unter Berücksichtigung der Konvention $$\text{sinc}\left(x\right):=\frac{\text{sin}x}{x}$$

  

  gilt: $$\widehat{1}p\left(\omega \right)=-\frac{1}{{\omega }^2}{{\displaystyle \sum _k\omega n_{k}e}}^{-i\omega m_k}\text{sinc}\omega {\text{t}}_k$$

  
• Das Verfahren kann ferner umfassen:
• Berechnen einer Vielzahl von Referenz-Luftbildern für eine zugehörige Vielzahl von Fokussiereinstellungen, FE, der Lithographie-Vorrichtung, wobei die Transmissions-Kreuzkoeffizientenmatrix jeweils konstant gehalten wird. Das Berechnen der Vielzahl von Referenz-Luftbildern kann ferner umfassen: Bestimmen einer Vielzahl von diskreten Pupillenfunktionen, K(FE).
• Hierbei ist jede diskrete Pupillenfunktion K(FE) eine diskrete Pupillenfunktion K des Pupil-Shift-Verfahrens, wie zuvor beschrieben, wobei die Unterschiede der jeweiligen diskreten Pupillenfunktion K(FE) durch die jeweilige hierdurch beschriebene Fokussiereinstellung FE der Lithographie-Vorrichtung bestimmt sind. Alternativ kann K(FE) auch mit K_z bezeichnet werden. Die Intensitätsverteilung des zu der jeweiligen Fokussiereinstellung zugehörigen Luftbildes kann als I_z bezeichnet werden. Einer Fokussiereinstellung kann in manchen Ausführungsbeispielen eine zugehörige Fokusebene im abbildenden System, beispielsweise in der Fotolithographie-Vorrichtung, zugeordnet werden.
• Hierbei kann jede diskrete Pupillenfunktion der Vielzahl von diskreten Pupillenfunktionen jeweils einer Fokussiereinstellung der Vielzahl von Fokussiereinstellungen zugeordnet sein. Ferner kann eines der vorstehend genannten Verfahren die Vielzahl der diskreten Pupillenfunktionen jeweilig durchgeführt werden.
• Hierbei kann das Verfahren ferner umfassen: Darstellen jeder der diskreten Pupillenfunktionen der Vielzahl von diskreten Pupillenfunktionen jeweils als Matrix; Laden, in einen zugehörigen Speicher, der Vielzahl von diskreten Pupillenfunktionen; Berechnen der Vielzahl von Referenz-Luftbildern ohne Nachladen von Matrizen in den zugehörigen Speicher.
• Solche Vorgehensweisen können den Vorteil haben, dass keine Zeitverluste durch das Nachladen von Informationen von dem gemeinsamen Speicher in den zugehörigen Speicher auftreten. Hierdurch können die kürzeren Zugriffszeiten der zugehörigen Speicher in solchen Vorrichtungen effektiv genutzt werden.
• Das Pupil-Shift-Verfahren kann hierbei eine oder mehrere Multiplikationen von Matrizen umfassen. Hierbei kann mindestens eine der Multiplikationen von Matrizen mittels Skalarmultiplikation durchgeführt werden.
• Dies kann den Vorteil haben, dass Skalarmultiplikationen auf manchen Hardwarekomponenten für parallele Berechnungen besonders zeiteffizient durchgeführt werden können und somit das Verfahren schneller durchgeführt werden kann.
• Das Verfahren kann mit Werten einfacher Genauigkeit in der Vielzahl von Recheneinrichtungen durchgeführt werden. Die Werte können hierbei als Gleitkommazahlen dargestellt werden, das Verfahren kann die Verwendung von Gleitkommaarithmetik umfassen.
• Dies kann den Vorteil haben, dass die Berechnungen, beispielsweise gegenüber von Berechnungen mit Werten doppelter Genauigkeit, beschleunigt durchgeführt werden können und dennoch ein hinreichend genaues Ergebnis erzielt werden kann.
• Hierbei kann ausgenutzt werden, dass manche Recheneinrichtungen, beispielsweise GPUs, effizienter mit Gleitkommaarithmetik von einfacher Genauigkeit umgehen können. Bei Bildgebungssimulationen, wie hier beschrieben, kann die einfache Genauigkeit den Vorteil haben, dass der erforderliche Speicherplatz im Vergleich zu doppelter Genauigkeit halbiert werden kann.
• Einfache Genauigkeit kann alternativ oder ergänzend bei der analytischen Spektrumsberechnung verwendet werden. Dies kann den Vorteil haben, dass Recheneinrichtungen der Vielzahl von Recheneinrichtungen, beispielsweise die GPUs, beispielsweise MSPs eines oder mehrerer GPUs, gleichmäßiger ausgelastet sein können, da die Berechnung mit einfacher Genauigkeit auf den Recheneinrichtungen deutlich weniger Rechenzeit benötigt als doppelte Genauigkeit. Dies kann einen zusätzlichen Geschwindigkeitsvorteil bringen, da die Spektrumsberechnung auf eine Vielzahl trigonometrischer Berechnungen, beispielsweise Sinus- und Cosinus-Funktionen, wie nachfolgend und zuvor beschrieben, angewiesen ist.
• Es wird ein Computerprogramm mit einem Programmcode bereitgestellt. Der Programmcode bewirkt ein Durchführen eines der Verfahren, wenn der Programmcode auf einer Vielzahl von Prozessoren und/oder einer Vielzahl von Parallelprozessoreinheiten und/oder einer programmierbaren Hardwarekomponente parallelisiert abläuft.
• Es wird eine für sich genommen nicht beanspruchte Vorrichtung zur Maskencharakterisierung bereitgestellt. Diese umfasst mindestens eine Rechenvorrichtung. Die Rechenvorrichtung umfasst eine Vielzahl von Recheneinrichtungen, wobei die Vielzahl von Recheneinrichtungen zum Durchführen eines der Verfahren eingerichtet ist.
• Die Vorrichtung umfasst ferner eine Messvorrichtung. Die Messvorrichtung ist eingerichtet, das gemessene Luftbild zu bestimmen und an die Rechenvorrichtung bereitzustellen.
• Die Messvorrichtung kann in einer Fotolithographie-Vorrichtung und/oder in einer Fotomaskenqualifikationsvorrichtung angeordnet sein.
• Hierbei kann in den Fällen, in denen die Messung in der Fotolithographie-Vorrichtung und in der Fotomaskenqualifikationsvorrichtung erfolgt, entweder mit einer Messvorrichtung gearbeitet werden, oder es können Informationen verschiedener Messvorrichtungen kombiniert werden. In solchen Fällen kann die Messvorrichtung eine erste Messvorrichtung und eine zweite Messvorrichtung umfassen.
• Das gemessene Luftbild kann mit einem Sensor, der eine räumliche Intensitätsverteilung in der Bildebene misst, bestimmt werden.
• Die Messvorrichtung kann hierbei einen Sensor umfassen. Der Sensor kann in eine Abbildungsebene eingebracht werden. Ergänzend oder alternativ kann die Vorrichtung zur Maskencharakterisierung einen variablen Strahlengang umfassen, wobei der variable Strahlengang in einer ersten Konfiguration das Licht in eine Bildebene abbildet und in einer zweiten Konfiguration das Licht auf einen Sensor, der sich nicht in der Bildebene der ersten Konfiguration befindet, abbildet.
• Figurenliste
• 1 zeigt schematisch eine Fotolithographie-Vorrichtung.
• 2 zeigt einen schematischen Aufbau einer Fotolithographie-Vorrichtung und ein mögliches Modell zur Beschreibung der Fotolithographie-Vorrichtung.
• 3 zeigt ein Ausführungsbeispiel eines Verfahrens.
• 4 zeigt ein Ausführungsbeispiel eines Verfahrens.
• 5 zeigt ein Diagramm zur Veranschaulichung eines Teils eines Verfahrens gemäß einem Ausführungsbeispiel.
• 6 zeigt anhand eines Ausführungsbeispiels, wie eine parallelisierte Berechnung vorgenommen werden kann.
• 7 zeigt ein Ausführungsbeispiel für eine Unterteilung des Referenz-Luftbildes in Blöcke.
• 8A zeigt einen Pseudocode für eine mögliche Vorgehensweise zum parallelen oder sequentiellen Berechnen eines Luftbildes.
• 8B zeigt eine mögliche Ablaufsteuerung eines Algorithmus gemäß verschiedener Ausführungsbeispiele.
• 8C zeigt einen Pseudocode für eine mögliche Vorgehensweise zum parallelen Berechnen eines Luftbildes.
• 9 zeigt einen Kernel, der eine skalare Matrixmultiplikation gemäß Gleichung (2.6) implementiert.
• 10 zeigt eine mögliche Implementierung des Kernels, welcher die Ausgabe der schnellen Fourier-Transformation zu dem Referenz-Bild akkumuliert.
• 11A und 11B zeigen benötigte Berechnungszeiten verschiedener Implementierungen gemäß verschiedener Ausführungsbeispiele als Funktion der Anzahl von Lichtquellen-Punkten.
• 12A und 12B zeigen eine Güte der parallelisierten Berechnungen gemäß verschiedener Ausführungsbeispiele.
• 13A und 13B zeigen Geschwindigkeitsvorteile der parallelen Berechnung mit einfacher Genauigkeit.
• AUSFÜHRLICHE BESCHREIBUNG
• Die oben beschriebenen Eigenschaften, Merkmale und Vorteile dieser Erfindung sowie die Art und Weise, wie diese erreicht werden, werden klarer und deutlicher verständlich im Zusammenhang mit der folgenden Beschreibung der Ausführungsbeispiele, die im Zusammenhang mit den Zeichnungen näher erläutert werden.
• Nachfolgend werden Ausführungsbeispiele unter Bezugnahme auf die Zeichnungen näher erläutert. In den Figuren bezeichnen gleiche Bezugszeichen gleiche oder ähnliche Elemente. Die Figuren sind schematische Repräsentationen verschiedener Ausführungsbeispiele. In den Figuren dargestellte Elemente sind nicht notwendigerweise maßstabsgetreu dargestellt. Vielmehr sind die verschiedenen in den Figuren dargestellten Elemente derart wiedergegeben, dass ihre Funktion und genereller Zweck dem Fachmann verständlich werden. In den Figuren dargestellte Verbindungen und Kopplungen zwischen funktionellen Einheiten und Elementen können auch als indirekte Verbindung oder Kopplung implementiert werden. Eine Verbindung oder Kopplung kann drahtgebunden oder drahtlos implementiert sein. Funktionale Einheiten können als Hardware, Software oder eine Kombination aus Hardware und Software implementiert werden.
• Nachfolgend werden im Zusammenhang mit 1 und 2 Vorrichtungen und Verfahren im Zusammenhang mit Fotolithographie erläutert. Hierbei werden die spezifischen Anwendungserfordernisse näher beschrieben.
• Im Zusammenhang mit den 3 und 4 werden verschiedene Ausführungsbeispiele für Verfahren erläutert.
• Im Zusammenhang mit den 5 bis 7 wird anschließend auf verschiedene Ausführungsbeispiele zur Implementierung von parallelen Berechnungen, unter anderem im Zusammenhang mit parallelen Berechnungen mittels Grafikprozessoreinheiten, eingegangen.
• Anschließend werden mögliche Vorteile von parallelisierten Vorgehensweisen im Zusammenhang mit den 11A und 11B bis 13A und 13B diskutiert.
• 1 zeigt schematisch eine Fotolithographie-Vorrichtung 90. Die Fotolithographie-Vorrichtung 90 umfasst die Komponenten: Beleuchtung, beispielsweise umgesetzt durch eine Beleuchtungskomponente 410, Objektebene 415, Bildgebung 420 und Bildebene 440. In der Bildebene 440 kann beispielsweise eine Abbildung auf einen Wafer 401 zum Zwecke der Fotolithographie erfolgen.
• 2 zeigt einen detaillierten schematischen Aufbau der Fotolithographie-Vorrichtung 90 mit Köhlerscher Beleuchtung und die mathematischen Annahmen zur Beschreibung der Fotolithographie-Vorrichtung 90. Die Beleuchtungskomponente 410 umfasst hierbei eine Lichtquelle 411, die beispielsweise als eine Menge von Punktlichtquellen angenommen werden kann, eine kollimierende Optik 412, eine Beleuchtungsebene 413 sowie eine Fokussierlinse 414. Die Beleuchtungsebene 413 wird hierbei manchmal auch als Beleuchtungspupillenebene oder Pupillenebene bezeichnet. Die Fokussierlinse 414 fokussiert das Licht in die Objektebene 415. Der abbildende Teil 420 setzt sich ebenfalls aus zwei Linsen 421, 423 zusammen. Hierbei kann der abbildende Teil 420 eine Bildblende 422 umfassen. Der abbildende Teil 420 bildet die Objektebene 415 auf die Bildebene 440 ab. In der Objektebene 415 kann eine Maske platziert sein um ein Bild in der Bildebene 440 zu erzeugen. In der Bildebene 440 kann in einer ersten Betriebsart ein Sensor 404 angeordnet werden, um ein gemessenes Luftbild 400 der Fotolithographie-Vorrichtung 90 basierend auf der Maske in der Bildebene 415 zu erhalten. In einer zweiten Betriebsart kann in der Bildebene 440 ein Wafer, beispielsweise der Wafer 401 aus 1, angeordnet sein. Hierbei kann die erste Betriebsart zu Diagnosezwecken der Fotolithographie-Vorrichtung verwendet werden und die zweite Betriebsart verwendet werden, um auf den Wafer fotolithographisch zu schreiben. In anderen Ausführungsbeispielen kann auch eine unabhängige Vorrichtung bereitgestellt werden, die lediglich die zweite Betriebsart bereitstellt. In solchen Fällen werden zwei unabhängige Vorrichtungen verwendet. Hierbei kann der Sensor 404 in nur einer oder in beiden Vorrichtungen verwendet werden. Alternativ oder ergänzend können verschiedene Sensoren verwendet werden, beispielsweise jeweils ein Sensor pro unabhängiger Vorrichtung.
• Hierbei können die folgenden Bezeichnungen verwendet werden:
• • S kann eine reale skalare Funktion sein, welche die Intensität der Lichtquelle, beispielsweise als Funktion der Koordinaten in der Beleuchtungsebene 413, beispielsweise in der Form von Lichtquellen-Punkten beschreibt.
• • F kann eine komplexe skalare Amplitudenverteilungsfunktion der Beleuchtungskomponente 410 sein.
• • J kann die gemeinsame Intensitätsfunktion bezeichnen, die das Objekt beleuchtet.
• • t kann eine komplexe skalare Funktion sein, die die Amplitude der Transmittanz, beispielsweise die Transmissionscharakteristik des Objekts, beispielsweise der Maske, beschreibt.
• • T kann das Spektrum des Objekts im Frequenzraum beschreiben, beispielsweise die Fourier-Transformation von t sein. T kann auch als Maskentransmissionsfunktion bezeichnet werden.
• • K kann die komplexe skalare Amplitudenverteilungsfunktion der Bildgebung 420 bezeichnen. K kann als diskrete Pupillenfunktion bezeichnet werden.
• • I kann die Intensitätsverteilung in der Bildebene 440 bezeichnen. I kann das Referenz-Luftbild, beispielsweise dargestellt als I(u,v), bezeichnen.
• Mögliche Implementierungen von S, F, J und K sind beispielsweise in:
• Konstantinos Adam: Domain Decomposition Methods for the Electromagnetic Simulation of Scattering from Three-Dimensional Structures with Applications in Lithography. PhD thesis, EECS Department, University of California, Berkeley, 2001,
• H. Kirchauer, Photolithography Simulation, Technische Universität Wien: Dissertationen der Technischen Universität Wien. Österreichischer Kunst- und Kulturverl., 2000,
• Bahaa E. A. Saleh and Malvin Carl Teich: Fundamentals of Photonics, Wiley, New York, 1991, und
• Herbert Gross, Wolfgang Singer, and Totzeck Michael: Handbook of Optical Systems - Volume 2: Physical Image Formation. Wiley-VCH Verlag GmbH & Co. KGaA, New York, 2005, beschrieben.
• Unter zahlreichen vereinfachenden Annahmen, insbesondere, dass die Lichtquelle 411 eine Punktlichtquelle oder eine Vielzahl von Punktlichtquellen ist und dass die Maske in der Objektebene 415 als unendlich dünn angenähert werden kann, kann das Lichtverhalten lokal zweidimensional in den Ortsvariablen x, y bzw. auf der Abbildungsseite mit den Ortsvariablen u, v sowie in Abhängigkeit von den Ortskoordinaten α, β der Lichtquelle 411 beschrieben werden.
• Dies ermöglicht eine deutliche Vereinfachung der Beschreibung gegenüber einer vollständigen Beschreibung mittels elektromagnetischer- und/oder Wellenoptik. Als charakteristisches Merkmal der Köhlerschen Beleuchtung trifft das Licht in der Objektebene 415 als ebene Wellen ein. Um das Verhalten des Lichts in der Beleuchtungskomponente 410 hinreichend genau zu beschreiben, kann somit die Funktion F(x, y, α, β) herangezogen werden, wobei die Funktion S(α, β) die Intensität der Lichtquelle, beispielsweise als Funktion der Koordinaten in der Beleuchtungsebene 413, beispielsweise in der Form von Lichtquellen-Punkten beschreiben kann. Hierbei kann S(α, β) insbesondere durch eine Superposition einer Vielzahl von Lichtquellen-Punkten, die in der Beleuchtungsebene 413 liegen, beschrieben werden.
• Die Abbildung einer Maske mit Transmissionscharakteristik t(x, y, α, β) als Funktion der Ortskoordinaten x und y und der lokalen Koordinaten α, β durch die Bildgebung 420 kann entsprechend durch eine diskrete Pupillenfunktion K(u, v, x, y) beschrieben werden, wobei x, y, u und v wie oben stehend definiert sind. Hierdurch kann eine Intensitätsverteilung I(u, v) in der Bildebene 440 berechnet werden. Diese Intensitätsverteilung I(u, v) kann das Referenz-Luftbild beschreiben.
• Die Messvorrichtung kann einen Sensor umfassen, der in der Bildebene der Fotolithographie-Vorrichtung platziert ist und eingerichtet ist, eine räumliche Intensitätsverteilung zu messen.
• Dies ist beispielsweise in den zuvor genannten Quellen von Yamazoe und Kirchauer näher beschrieben. Hierdurch kann der Rechenaufwand bereits gegenüber vollen elektromagnetischen Simulationen und/oder vollständigen wellenoptischen Simulationen des Lichtverlaufs stark verringert werden.
• 3 zeigt ein Ausführungsbeispiel eines Verfahrens 1000. Das Verfahren 1000 wird im Zusammenhang mit den 3 und 4 erläutert.
• Das Verfahren 1000 umfasst bei 1010 ein Empfangen eines Designs einer Maske. Das Design der Maske kann das Design der Maske 110 der 4 sein.
• Das Verfahren umfasst bei 1020 ein Empfangen von optischen Eigenschaften. Die optischen Eigenschaften können die optischen Eigenschaften 120 der 4 sein und können optische Eigenschaften der Fotolithographie-Vorrichtung 90 der 1 und 2 umfassen. Die optischen Eigenschaften können beispielsweise gemäß dem Pupil-Shift-Verfahren, wie zuvor und nachfolgend erläutert, beschrieben werden.
• Bei 1030 erfolgt ein paralleles Berechnen eines Referenz-Luftbildes 300 mit einer Vielzahl von Recheneinrichtungen basierend auf dem Design der Maske 110 und den optischen Eigenschaften 120. Das Design der Maske 110 kann hierbei basierend auf Dreiecken, wie zuvor und nachfolgend erläutert, beschrieben werden.
• Bei 1040 erfolgt ein Empfangen eines gemessenen Luftbildes 400. Das gemessene Luftbild kann mit einer Messvorrichtung in der Fotolithographie-Vorrichtung gemessen worden sein. Beispielsweise kann das gemessene Luftbild 400 durch Platzierung eines Sensors 404 in der Bildebene 440 bestimmt werden.
• Bei 1050 erfolgt ein Vergleichen des Referenz-Luftbilds 200 mit dem gemessenen Luftbild 400.
• Das Vergleichen kann ein Vergleichen der Intensitätswerte I(u, v) des Referenz-Luftbildes mit den gemessenen räumlichen Intensitätswerten des gemessenen Luftbildes umfassen. Das Vergleichen der Intensitätswerte kann pixelweise erfolgen, wenn die jeweiligen Koordinaten des Referenz-Luftbildes und des gemessenen Luftbildes einander entsprechen, beispielsweise weil die Koordinaten des Referenz-Luftbildes vor der Berechnung entsprechend der Koordinaten des gemessenen Luftbildes gewählt wurden. In anderen Ausführungsbeispielen kann auch eine oder mehrere Koordinatentransformation durchgeführt werden, beispielsweise von dem Referenz-Luftbild und/oder dem gemessenen Luftbild, um einen genauen Vergleich zu ermöglichen.
• Zur Quantifizierung des Vergleichens zwischen dem Referenz-Luftbild und dem gemessenen Luftbild stehen zahlreiche Möglichkeiten zur Verfügung, die für sich genommen von Bildvergleichen bekannt sind. Beispielsweise kann eine Pixel-zu-Pixel Differenz betrachtet werden oder Merkmal-basierte (englisch: feature-based) Differenzen betrachtet werden. Hierbei kann die Quantifizierung in einer Zahl (beispielweise entsprechend einem Abweichungsquadrat) oder in eine Vielzahl von Zahlen, beispielsweise verschiedene Zahlen für verschiedene Merkmale oder verschiedene Bereiche von Referenz-Luftbild und gemessenem Luftbild, erfolgen.
• Bei 1060 erfolgt ein Bestimmen einer Qualität der Maske 500 basierend auf dem Vergleichen bei 1050.
• Hierbei kann die Qualität der Maske basierend auf der Quantifizierung des Vergleichens erfolgen. Die Qualität der Maske kann hierbei beispielsweise durch einen vorgewählten Schwellenwert bestimmt werden und eine binäre Entscheidung sein, beispielsweise „Maske in Ordnung“ und „Maske nicht in Ordnung“. In anderen Ausführungsbeispielen kann die Qualität der Maske als eine Zahl quantifiziert werden, beispielsweise als eine Gütezahl zwischen 100% für eine perfekte Übereinstimmung und 0% für eine vollständige nicht-Übereinstimmung. In wieder anderen Ausführungsbeispielen kann die Qualität der Maske als eine Vielzahl von Zahlen und/oder Parametern beschrieben werden, beispielsweise als „Intensitätsentsprechung 99,5%, Qualität der Leiterbahnen-Strukturen - sehr gut, Qualität von Transistor-Strukturen: nicht akzeptabel“.
• Die Qualität kann dann auf Grundlage von Mindesterfordernissen beurteilt werden, beispielsweise kann im vorherigen Beispiel die Qualität der Maske mit „nicht in Ordnung“ beurteilt werden, weil die Qualität der Transistor-Strukturen nicht akzeptabel ist. Hierbei kann dann insbesondere die Intensitätsentsprechung von 99,5% unberücksichtigt bleiben. Mit anderen Worten kann ein einzelner lokaler Fehler zur Beurteilung der Qualität herangezogen werden.
• Das parallele Berechnen, wie in 3 bei Schritt 1030 erläutert, kann basierend auf dem Pupil-Shift-Verfahren erfolgen. Hierbei kann das Pupil-Shift-Verfahren für eine parallele Berechnung mit der Vielzahl von Recheneinrichtungen implementiert sein.
• Dabei kann das Referenz-Luftbild 300 in der Bildebene 400 berechnet werden, indem die vereinfachten Annahmen des Pupil-Shift-Verfahrens vorteilhaft angewendet werden. In manchen Ausführungsbeispielen kann das Pupil-Shift-Verfahren basierend auf einer gemeinsamen Intensitätsfunktion J, einer Maskentransmissionsfunktion T und einer diskreten Pupillenfunktion K durchgeführt werden. In weiteren Ausführungsbeispielen kann die gemeinsame Intensitätsfunktion und die Maskentransmissionsfunktion als eine Transmissions-Kreuzkoeffizientenmatrix, TCC-Matrix, dargestellt werden.
• Die Intensität I(u, v) des Referenz-Luftbildes kann dann durch $$I\left(u,v\right)={\left|{\displaystyle \underset{-\infty }{\overset{\infty }{\int }}{\displaystyle \int K\left(u,v,x,y\right)}}t\left(x,y\right)J\left(x,y\right)dxdy\right|}^2$$

  

  beschrieben werden. Hierbei wird als Approximation verwendet, dass die Amplitude der Maskentransmission t(x, y) unabhängig von der Lichtquelle S(α, β) ist. Aufgrund dieser Annahme ist es möglich, die Maske als die Maskentransmissionsfunktion T als ein zweidimensionales komplexes Feld zu beschreiben, was durch Fourier-Transformation der Maskentransmission t(x, y) wie folgt bestimmt werden kann: $$T\left(n,m\right)=\frac{M^2}{ab}{\displaystyle {\int }_0^a{\displaystyle {\int }_0^{b}t}}\left(x,y\right)e^{-j2\pi \left(nxM/a+myM/b\right)}dxdy$$

  
• Die Rücktransformation ist entsprechend gegeben durch: $$t\left(x,y\right)={\displaystyle {\int }_0^{a\text{/}M}{\displaystyle {\int }_0^{b\text{/}M}T}}\left(n,m\right)e^{j2\pi \left(nxM\text{/}a+myM\text{/}b\right)}dndm$$

  
• Hierbei sind jeweils a/M und b/M die Objektperioden, so dass a und b die Bildperioden für ein Projektionssystem mit einer Vergrößerung M=1 in x- und y-Richtung sind.
• Hierbei wird ferner angenommen, dass die Maske periodisch sowohl in x- als auch in y-Richtung ist, so dass gilt: $$t\left(x+a\text{/}M,y+b\text{/}M\right)=t\left(x,y\right)$$

  
• Durch Fourier-Transformation von Gleichung (2.2) kann die Faltung als eine Multiplikation in dem Frequenzraum geschrieben werden, welche entsprechend der Beleuchtungsrichtung von J verschoben (englisch: shifted) ist. Somit kann Gleichung 2.2 geschrieben werden als: $$I\left(u,v\right)=\stackrel{\infty }{{\displaystyle \underset{-\infty }{\iint }J}}\left(p,q\right){\left|{\displaystyle {\int }_0^{a/M}{\displaystyle {\int }_0^{b/M}T}}\left(n,m\right)K\left(n,m\right)e^{j2\pi \left(\left(n-p\right)u/a+\left(m-q\right)v/b\right)}dndm\right|}^{2}dpdq$$

  
• Führt man nun die Verschiebung in der diskreten Pupillenfunktion K durch, erhält man die als Pupil-Shift-Verfahren bezeichnete Gleichung für die Intensitätsverteilung I(u, v) des Luftbildes: $$I\left(u,v\right)=\stackrel{\infty }{{\displaystyle \underset{-\infty }{\iint }J}}\left(p,q\right){\left|{\displaystyle {\int }_0^{a/M}{\displaystyle {\int }_0^{b/M}T}}\left(n,m\right)K\left(n-p,m-q\right)e^{j2\pi \left(nu\text{/}a+mv\text{/}b\right)}dndm\right|}^{2}dpdq$$

  
• Alternativ ist es möglich, T anstelle von K zu verschieben. Dies wird als Spektral-Shift-Verfahren (engl.: spectral shift method) bezeichnet: $$I\left(u,v\right)=\stackrel{\infty }{{\displaystyle \underset{-\infty }{\iint }J}}\left(p,q\right){\left|{\displaystyle {\int }_0^{a/M}{\displaystyle {\int }_0^{b/M}T}}\left(n-p,m-q\right)K\left(n,m\right)e^{j2\pi \left(nu\text{/}a+mv\text{/}b\right)}dndm\right|}^{2}dpdq$$

  
• Zur effektiven Berechnung kann, wie in Kenji Yamazoe, „Two models for partially coherent imaging“, J. Opt. Soc. Am. A, 29(12): 2591-2597, Dezember 2012, beschrieben, eine Eigenfunktionsanalyse durchgeführt werden. Dies wird als Transmissions-Kreuzkoeffizientenmatrix-,TCC-, Ansatz bezeichnet. Ein Umformen von Gleichung 2.5 liefert hierbei: $$\begin{array}{l}I\left(u,v\right)={\displaystyle \iiint {\displaystyle \int TCC}}\left(x_1,y_1,x_2,y_2\right)T\left(x_1,y_1\right)T*\left(x_2,y_2\right)e^{j2\pi \left(\left(x_1-x_2\right)u+\left(y_1-y_2\right)v\right)}d{x}_{1}d{x}_{2}d{y}_{1}d{y}_2\\ TCC\left(x_1,y_1,x_2,y_2\right)={\displaystyle \iint J\left(f,g\right)}K\left(f+x_1,g+y_1\right)K*\left(f+x_2,g+y_2\right)dfdg\end{array}$$

  
• Die Transmissions-Kreuzkoeffizientenmatrix TCC enthält hierbei die optischen Einstellungen der Lithographie-Vorrichtung und kann für statische optische Einstellung vorberechnet werden. Die Transmissions-Kreuzkoeffizientenmatrix kann mittels ihrer Eigenwerte µ_j und Eigenfunktionen Φ_j diagonalisiert werden.
• Nimmt man TCC als eine hermitesche Matrix an, kann das Eigenwert- und Eigenfunktionsproblem folgendermaßen formuliert werden: $$I\left(u,v\right)={\displaystyle \sum _{j=1}^E{u}_j}{\left|{\displaystyle \iint T\left(n,m\right)}{\varphi }_j\left(n,m\right)e^{j2\pi \left(nu+mv\right)}dndm\right|}^2$$

  
• Durch die Darstellung als Summe von Eigenwerten kann - bei entsprechender Sortierung der Darstellung nach betraglicher Größe der Eigenwerte - ausgenutzt werden, dass der Beitrag der einzelnen Terme der Summe zu dem Referenz-Luftbild abnimmt. Hierdurch kann die Anzahl an Fourier-Transformationen, die benötigt werden, um das Referenz-Luftbild I(u, v) hinreichend genau zu berechnen, verringert werden. Insbesondere kann E ≤ N für eine Anzahl von N Lichtquellen-Punkten angenommen werden.
• Da die Fourier-Transformationen einen großen Beitrag zu der benötigten Laufzeit beitragen, kann die Reduktion der benötigten Fourier-Transformationen zu einer starken Beschleunigung des Verfahrens führen. Hierbei kann ausgenutzt werden, dass die TCC-Matrix mittels Eigenfunktionen und Eigenwerten wie zuvor beschrieben als eine positiv definite hermetische Matrix diagonalisiert werden kann. Dies kann, beispielsweise einmalig, in externen Programmen durchgeführt werden, beispielsweise in MATLAB. Eine solche Vorgehensweise kann dann vorteilhaft sein, wenn die optischen Einstellungen, welche durch die TCC-Matrix beschrieben werden, konstant bleiben.
• Hierbei kann die TCC-Matrix in diagonalisiert approximiert werden. Beispielsweise können nur die betraglich größten Eigenwerte verwendet werden und die restlichen Eigenwerte auf null gesetzt werden. Hierdurch kann die Berechnung weiter beschleunigt werden. Dies wird nachfolgend näher erläutert. Ergebnisse zum Einfluss der Eigenwertauswahl auf Genauigkeit und Laufzeit wird in den 12A und 12B gezeigt.
• In manchen Ausführungsbeispielen wird die Maskentransmissionsfunktion T aus Designdaten bestimmt. Wie in Gleichung (2.2b) beschrieben ist die Maskentransmissionsfunktion T die Fourier-Transformation der Maskentransmission t(x, y). Die Maskentransmission an einem beliebigen Punkt kann entsprechend des in H. Kirchauer, Photolithography Simulation, Technische Universität Wien: Dissertationen der Technischen Universität Wien. Österreichischer Kunst- und Kulturverl., 2000 beschriebenen Verfahrens folgendermaßen dargestellt werden:
• Der nachfolgend in der Beschreibung, im Zusammenhang mit den Formeln (2.9) bis (2.19), verwendete Index k ist dabei ein Index für Muster und unterscheidet sich vom zuvor verwendeten Index für Ecken des Polygons.
$$t\left(x,y\right)={\displaystyle \sum _k{t}_k}{\text{PAT}}_k\left(x,y\right)$$

wobei $${\text{PAT}}_k\left(x,y\right)=\{\begin{array}{ll}1\hfill & \text{falls}\left(\text{x},\text{y}\right)\text{innerhalb des Musters k}\hfill \\ 0\hfill & \text{sonst}\hfill \end{array}$$

gilt. PAT_k(x, y) beschreibt hierbei eine beliebiges elementares Muster der Maske. Der Transmissionskoeffizient t_k in Gleichung (2.9) ist gegeben durch: $$t_k\{\begin{array}{ll}\equiv 1\hfill & \text{f}\ddot{\text{u}}\text{r bin}\ddot{\text{a}}\text{re Masken}\hfill \\ \in \text{C}\hfill & \text{f}\ddot{\text{u}}\text{r phasenverschiebende Masken}\hfill \end{array}$$

• Wenn a und b die Breite und Höhe der Maske bezeichnen, kann T für eine 1:1 Abbildung beschrieben werden durch: $$T\left(n,m\right)=\frac{1}{ab}{\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}t\left(x,y\right)}}{e}^{-j2\pi \left(nx\text{/}a+my\text{/}b\right)}dydx$$

  
• Durch Einsetzen von Gleichungen (2.9), (2.10) und (2.11) in Gleichung (2.12) erhält man: $$T\left(n,m\right)={\displaystyle \sum _k{t}_k}\left[\frac{1}{ab}{\displaystyle \underset{0}{\overset{a}{\int }}{\displaystyle \underset{0}{\overset{b}{\int }}{\text{PAT}}_k}}\left(x,y\right)e^{-j2\pi \left(nx\text{/}a+my\text{/}b\right)}dydx\right]$$

  
• Hierbei bezeichnet der Ausdruck in eckigen Klammern das elementare Musterspektrum P_k(n, m): $$T\left(n,m\right)={\displaystyle \sum _k{t}_k}{P}_k\left(n,m\right)$$

  
• Dieses kann analytisch für jedes gegebene Dreieck k, welches aus den zugehörigen drei Punkten $$\left\{\left(x_{k1},y_{k1}\right),\left(x_{k2},y_{k2}\right),\left(x_{k3},y_{k3}\right)\right\}\in {ℝ}^2$$

  

  besteht, berechnet werden gemäß: $$P_k\left(n,m\right)=\frac{\left|\text{det A}\right|}{ab}{e}^{-j2\pi c_{nm}}W\left(a_{nm},b_{nm}\right)$$

  

  $$W\left(a_{nm},b_{nm}\right)={\displaystyle \underset{0}{\overset{1}{\int }}{\displaystyle \underset{0}{\overset{1-\xi }{\int }}{e}^{-j2\pi \left(a_{nm}\xi +b_{nm}\eta \right)}}}d\eta d\xi$$

  

  wobei $$A=\left(\begin{array}{cc}{x}_{k2}-x_{k1}& x_{k3}-x_{k1}\\ y_{k2}-y_{k1}& y_{k3}-y_{k1}\end{array}\right)$$

  

  die Transformationsmatrix bezeichnet, die beliebige Dreiecke in das Einheitsdreieck {(0,0), (1,0), (0,1)} überführt.
• Die Koeffizienten a_nm, b_nm und c_nm sind hierbei gegeben durch: $$\begin{array}{r}\hfill a_{nm}=\frac{n\left(x_{k2}-x_{k1}\right)}{a}+\frac{m\left(y_{k2}-y_{k1}\right)}{b}\\ \hfill b_{nm}=\frac{n\left(x_{k3}-x_{k1}\right)}{a}+\frac{m\left(y_{k3}-y_{k1}\right)}{b}\\ \hfill c_{nm}=\frac{n{x}_{k1}}{a}+\frac{m{y}_{k1}}{b}\end{array}$$

  
• Durch Betrachtung von Grenzwerten und durch Transformation zu Einheitsdreiecken kann das Integral von Gleichung 2.16 umgeformt werden zu: $$W\left(a_{nm},b_{nm}\right)=\{\begin{array}{ll}\frac{1}{2}\hfill & \text{falls }{a}_{nm}=b_{nm}=0\hfill \\ \frac{j}{2\pi a_{nm}}\left(e^{-f2\pi a_{nm}}-sinc\left(\pi a_{nm}\right)e^{-j\pi a_{nm}}\right)\hfill & \text{falls }{a}_{nm}=b_{nm}\ne 0\hfill \\ \frac{j}{2\pi }\left(\frac{sinc\left(\pi a_{nm}\right)}{a_{nm}-b_{nm}}{e}^{-j\pi a_{nm}}+\frac{sinc\left(\pi b_{nm}\right)}{b_{nm}-a_{nm}}{e}^{-j\pi b_{nm}}\right)\hfill & \text{sonst}\hfill \end{array}$$

  
• Die Punkte von Dreiecken, die die Maske darstellen, können beispielweise als ein Array von dreidimensionalen realen Gleitkommazahlen mit doppelter Genauigkeit (real double precission float array) bereitgestellt werden.
• Die Transmittanz jeweiliger Dreiecke kann in einem eindimensionalen komplexen Array von realen Gleitkommazahlen mit doppelter Genauigkeit bereitgestellt sein.
• Gitter spektraler Koordinaten in x- und y-Richtung können jeweils als eindimensionale Arrays von realen Gleitkommazahlen mit doppelter Genauigkeit bereitgestellt werden.
• Die Größe des Referenz-Luftbildes in Pixeln in x- und y-Richtung können jeweils in der Form von realen Gleitkommazahlen mit doppelter Genauigkeit bereitgestellt werden.
• Aber auch andere Arten, die Werte bereitzustellen, sind möglich.
• Durch solche Vorgehensweisen kann das parallele Berechnen 1030 mit einer Vielzahl von Recheneinrichtungen des Referenz-Luftbildes 300 beschleunigt werden. Dies wird im Folgenden näher beschrieben.
• In manchen Ausführungsbeispielen können einige der benötigten Koeffizienten der Transformationsmatrix A aus Gleichung (2.17) konstant sein. Ferner kann der Transmittanzkoeffizient des Dreiecks, beispielsweise der Transmissionskoeffizient t_k, konstant sein.
• Um diese konstanten Koeffizienten effektiv bei der parallelen Berechnung teilen zu können, beispielsweise zwischen einzelnen Threads, kann ein Gitter so gewählt werden, dass jeder Block ein Dreieck berechnet. In solchen Ausführungsbeispielen können alle Threads dieses Blocks den Beitrag des jeweiligen Dreiecks für alle spektralen Koordinaten berechnen. Schließlich kann dieses so erhaltene Zwischenergebnis einem globalen Spektrum, beispielsweise für alle Dreiecke zusammen, hinzugefügt werden. Hierdurch kann dann T bestimmt werden.
• Hierbei kann es nachteilhaft sein, falls verschiedene Blöcke zeitgleich auf die gleiche spektrale Koordinate zugreifen. Daher kann es notwendig sein, diese Additionen atomar (englisch: atomic addition), beispielsweise mittels einer atomaren Bibliothek, beispielsweise mittels der atomaren CUDA-Funktionalitäten, durchzuführen.
• Alternativ oder ergänzend kann die Berechnung von trigonometrischen Funktionen im Ausdruck W(a_nm, b_nm) vereinfacht werden wie nachfolgend beschrieben.
• Der Ausdruck $$\text{sin}\left(\pi x\right)e^{-j\pi x}$$

  

  kann vereinfacht werden zu: $$sinc\left(\pi x\right)e^{-j\pi x}=\{\begin{array}{ll}1\hfill & \text{falls }\left|\pi x\right|\le {10}^{-e}\hfill \\ \frac{cos\left(\pi x\right)sin\left(\pi x\right)-jsi{n}^2\left(\pi x\right)}{\pi x}\hfill & \text{sons}\hfill \end{array}$$

  
• Hierbei sind die Argumente, die an Sinus und Cosinus übergeben, werden positiv. Der Term -πx, der ab e^-j übergeben wird, kann unter Verwendung der trigonometrischen Identität $$cos\left(-\pi x\right)+jsin\left(-\pi x\right)=cos\left(\pi x\right)-jsin\left(\pi x\right)$$

  

  eliminiert werden.
• Ebenfalls kann der Term e^-j2πx, der benötigt wird, um Fall 2 von W(a_nm, b_nm) zu berechnen, folgendermaßen umgeformt werden: $$e^{-j2\pi x}=cos\left(2\pi x\right)-jsin\left(2\pi x\right)=co{s}^2\left(\pi x\right)-si{n}^2\left(\pi x\right)-2jsin\left(\pi x\right)cos\left(\pi x\right)$$

  
• Dies kann gestatten, die Koeffizienten von Sinus und Cosinus nur einmal zu berechnen und für alle drei Fälle von W(a_nm, b_nm) zu verwenden. Hierdurch kann das Verfahren weiter beschleunigt werden.
• Ein Diagramm zur Veranschaulichung eines Teils eines Verfahrens gemäß einem Ausführungsbeispiels ist in 5 gezeigt.
• Das im Zusammenhang mit 5 beschriebene Verfahren kann bei Schritt 1030 der 3 verwendet werden.
• Paralleles Berechnen kann mittels einzelnen Threads parallelisiert durchgeführt werden. Hierbei bezeichnen Threads einzelne Prozessschritte, die parallelisiert oder sequentiell abgearbeitet werden können. Blöcke können in Abhängigkeit der vorhandenen Hardware gebündelte Threads bezeichnen. Das Gesamtproblem, welches berechnet werden soll, kann dann auf einem Gitter (engl.: grid) definiert sein und eine Bündelung von Blöcken bezeichnen. Ein konkretes Problem kann mithin auf einem Gitter definiert werden und dann vorteilhaft parallelisiert berechnet werden.
• Diese Vorgehensweise, das Gitter abzuarbeiten, kann als Gitter-Stride (von englisch stride - durchlaufen, abschreiten) bezeichnet werden. Dies kann beispielsweise mittels einer Gitter-Stride-Schleife realisiert werden.
• Solche Vorgehensweisen können die effiziente Durchführung eines Problems auch auf verschiedenen Hardwarekomponenten erleichtern, weil die Verteilung des Gitters auf die Blöcke in Abhängigkeit der vorhandenen Hardware gewählt werden kann und somit weitestgehend oder vollständig unabhängig von dem konkreten Problem vorgenommen werden kann.
• Vorliegend kann das zu simulierende Referenz-Luftbild 300 auf einem Gitter G beschrieben werden und hierzu in eine Anzahl von Blöcken B0 bis B8 aufgeteilt werden. Das Referenz-Luftbild kann als Intensität von zwei räumlichen Koordinaten, beispielsweise u und v, berechnet werden, beispielsweise als Intensitätsverteilung I(u, v) wie zuvor beschrieben. Die Intensitätsverteilung kann auf dem Gitter definiert werden.
• Hierbei kann es möglich sein, die einzelnen Pixel der Intensitätsverteilung unabhängig voneinander zu berechnen. In solchen Fällen kann das parallele Berechnen pixelweise erfolgen. Hierzu können die Pixel gruppiert und jeweiligen Blöcken, beispielsweise den Blöcken B0 bis B8, zugeordnet werden.
• Exemplarisch ist der Block B8 zusätzlich vergrößert in 5 dargestellt. Im gezeigten Ausführungsbeispiel ist das simulierte Referenz-Luftbild 300 in neun Blöcke B0 bis B8 unterteilt. Jeder der Blöcke B0 bis B8 umfasst wiederum neun Pixel. Es wird ein Ausführungsbeispiel für den Block B8 gezeigt, der die Pixel 0,0 bis Pixel 2,2 umfasst. Dieses Ausführungsbeispiel ist hierbei nur schematisch und die Anzahl der Pixel kann unterschiedlich sein, insbesondere müssen weder das simulierte Referenz-Luftbild 300 noch die Blöcke quadratisch aufgebaut sein.
• 6 zeigt am Ausführungsbeispiel eines Grafikprozessors, wie eine parallelisierte Berechnung der Blöcke der 5 mit Hilfe von verschiedenen Recheneinrichtungen vorgenommen werden kann.
• Im gezeigten Ausführungsbeispiel der 6 setzt sich ein Programm P0 zur parallelisierten Simulation eines Referenz-Luftbilds aus acht Blöcken B0 bis B7 zusammen. Mit anderen Worten verwendet das in 6 gezeigte Ausführungsbeispiel einen Block weniger als das in 5 gezeigte Ausführungsbeispiel.
• In einem ersten Ausführungsbeispiel A der 6 steht als die Vielzahl von Recheneinrichtungen ein Grafikprozessor GPU-1 mit zwei Streaming Multiprozessoren (englisch: „Streaming Multiprocessors“) SMs SM0, SM1 zur Verfügung. In zeitlicher Abfolge T1 werden nun die jeweiligen Blöcke berechnet. Da zwei Streaming Multiprozessoren zur Verfügung stehen, werden zunächst die Blöcke 0 und 1, anschließend die Blöcke 2 und 3, gefolgt von den Blöcken 4 und 5, abschließend gefolgt von den Blöcken 6 und 7 berechnet. Die Reihenfolge der Blöcke ist hier nur zu Illustrationszwecken angeführt.
• In einem zweiten Ausführungsbeispiel B steht eine andere Hardwarekonfiguration zur Verfügung. Die Vielzahl der Recheneinrichtungen ist im Ausführungsbeispiel B ein Grafikprozessor GPU-2 mit vier Streaming Multiprozessoren SM0 bis SM3. Entsprechend können parallel vier Blöcke des Referenz-Luftbildes gleichzeitig berechnet werden. Im gezeigten Ausführungsbeispiel der 6, Ausführungsbeispiel B, werden zunächst die Blöcke 0 bis 3 und anschließend die Blöcke 4 bis 7 berechnet. Hierbei ist die gezeigte Reihenfolge und Zuordnung zu den Streaming Multiprozessoren, ebenfalls wie im Ausführungsbeispiel A, als Ausführungsbeispiel angegeben und beliebig modifizierbar.
• 7 zeigt ein Ausführungsbeispiel einer Architektur eines Computers mit einem Grafikprozessor.
• Die gezeigte Architektur kann die Vielzahl von Recheneinrichtungen bereitstellen, wie nachfolgend erläutert wird.
• Der Computer umfasst einen Speicher MEM-HOST sowie eine Grafikkarte GPU-1. Die Grafikkarte GPU-1 der 7 kann beispielsweise die Grafikkarte GPU-1 der 6 sein. Aber auch eine andere Grafikkarte ist möglich. Im gezeigten Ausführungsbeispiel der 7 umfasst die Grafikkarte GPU-1 Streaming-Multiprozessoren SM0, SM1. Die Streaming-Multiprozessoren SM0, SM1 umfassen jeweils eine Vielzahl von Skalarprozessoren SKP sowie zugehörige geteilte Speicher SM-MEMO bzw. SM-MEM1. Ferner umfasst die Grafikkarte GPU-1 einen Grafikkartenspeicher MEM-GPU.
• In manchen Ausführungsbeispielen kann hierbei der geteilte Speicher SM-MEMO, SM-MEM1 schneller sein als der Grafikkartenspeicher MEM-GPU und/oder der Arbeitsspeicher MEM-HOST. Bei solchen Parallelprozessoreinheiten können die hier beschriebenen Verfahren, die ohne Nachladen von Matrizen in den Speicher auskommen, den Vorteil haben, dass die Vielzahl von diskreten Pupillenfunktionen in einen schnellen Speicher, beispielsweise den geteilten Speicher SM-MEMO und/oder SM-MEM1, geladen werden können. Hierdurch können die Berechnungen von Referenz-Luftbildern, insbesondere bei Veränderungen der optischen Eigenschaften der Fotolithographie-Vorrichtung, besonders schnell durchgeführt werden.
• Die parallele Berechnung kann auf verschiedene Arten implementiert werden. Sie kann beispielsweise mittels der MEX-Schnittstelle von MATLAB in einer Hochsprache, beispielsweise C/C++ implementiert werden. Wie zuvor beschrieben, kann die partiell kohärente Bildgebungssimulation mittels verschiedener Funktionen durchgeführt werden. Diese können als ein-, zwei- und/oder dreidimensionale Arrays, beispielsweise als komplexe Gleitkommazahlen mit einfacher oder doppelter Genauigkeit implementiert werden.
• Beispielsweise kann T als ein zweidimensionales Array komplexer Zahlen mit doppelter Genauigkeit, K als ein dreidimensionales Array komplexer Zahlen mit doppelter Genauigkeit, L als ein eindimensionales Array komplexer Zahlen mit einfacher Genauigkeit und J als ein eindimensionales Array realer Zahlen mit doppelter Genauigkeit angegeben werden.
• Die MEX-Funktion kann ein oder mehrere Referenz-Luftbilder parallel berechnen. Das Referenz-Luftbild kann als ein zweidimensionales Array realer Zahlen mit einfacher oder doppelter Genauigkeit dargestellt werden.
• 8A zeigt einen Pseudocode für eine mögliche Vorgehensweise zum parallelen oder sequentiellen Berechnen eines Luftbildes.
• Der gezeigte Pseudocode ist in der Lage, gemäß $$I_z\left(u,v\right)={\displaystyle \sum _{i}J\left({\alpha }_i,{\beta }_i\right)}{\left|\text{IFT}\left[T\left(f-{\alpha }_i,g-{\beta }_i\right)K_z\left(f,g\right)\right]\right|}^2$$

  

  das Luftbild I_z für eine jeweilige Fokussiereinstellung, bezeichnet mit K_z, zu berechnen. Hierbei wird der Index „z“ zur Bezeichnung der Fokussiereinstellung verwendet. Der Index z kann der Notation „FE“ für die Fokussiereinstellung, wie zuvor erläutert, entsprechen.
• Es können beispielsweise Z=20 verschiedene Fokussiereinstellungen, beispielsweise bezeichnet mit z=0..19, berechnet werden.
• Hierbei bezeichnet IFT eine inverse zweidimensionale Fourier-Transformation.
• (α_i, β_i) bezeichnen die Koordinaten des i-ten Lichtquellen-Punktes und (f, g) bezeichnet die Koordinatennotation der IFT.
• Die zwei äußeren Integrale über die Lichtquelle wurden hierbei in eine einzige Summation, die über eine diskretisierte zweidimensionale Lichtquelle, beispielsweise durch Superposition von einer Vielzahl von Lichtquellen-Punkten dargestellt, läuft, zusammengefasst.
• Der in 8A gezeigte Pseudocode benötigt vier Eingabeargumente und liefert eine Ausgabe.
• Die Ausgabe I ∈ ℝ^N×M×Z ist ein Tensor dritter Stufe, welcher für jede Fokusebene z das Referenz-Luftbild beschreibt.
• Die vier Eingabeargumente sind:
• ein Vektor, J ∈ ℝ^S, der die Lichtintensitäten für jeden Lichtquellen-Punkt beschreibt,
• ein Vektor L ∈ ℂ^S, der die Koordinaten jedes Lichtquellen-Punktes beschreibt,
• eine Matrix T ∈ ℂ^N×M, die das Fourier-Spektrum der Maske, beispielsweise eine Maskentransmissionsfunktion, beschreibt,
• ein Tensor dritter Stufe, K ∈ ℂ^N×M×Z der die zweidimensionale Pupillenfunktion, beispielsweise die diskrete Pupillenfunktion, des abbildenden Systems für Z Fokusebenen beschreibt.
• Der Pseudocode durchläuft mit einer Schleife alle Lichtquellen-Punkte und führt folgende Schritte durch:
1. 1. Verschieben (englisch: shifting) des Objektspektrums T entsprechend der Beleuchtungsbedingungen,
2. 2. Anwenden der diskreten Pupillenfunktion K des abbildenden Systems auf das verschobene Array und Durchführen von Z inversen Fourier-Transformationen,
3. 3. Akkumulieren des Gains der diskreten Fourier-Transformation, gewichtet mit der gemeinsamen Intensitätsfunktion J zu dem Referenz-Bild I.
• Das Pupil-Shift-Verfahren, wie hier und im Ausführungsbeispiel der 8A implementiert, durchläuft Schleifen über eine zweidimensionale gemeinsame Intensitätsfunktion J und ein zweidimensionales Objekt T und ein abbildendes System K.
• Die innersten zweidimensionalen Fourier-Transformationen haben, beispielsweise in ihrer reinsten Form, eine Komplexität von O((M × N)²). Um ein einzelnes Referenz-Bild, beispielsweise für eine Fokussiereinstellung z, zu simulieren, muss die Fourier-Transformation S mal berechnet werden, wobei S die Dimensionalität der zweidimensionalen Lichtquelle ist, beispielsweise S = S_x × S_y.
• Im Pseudocode der 8A wurde das Pupil-Shift-Verfahren dahingehend erweitert, Referenz-Bilder für Z verschiedene Fokusebenen zu berechnen. Daher muss das Verfahren Z mal durchgeführt werden. Vom Rechenaufwand her ist es - unter der Annahme, dass die Kosten für Speicherzugriff ignoriert werden - irrelevant, ob die Bildgebungssimulation in allen Fokusebenen, wie in 8A gezeigt, oder jeweils separat für die einzelnen Fokusebenen durchgeführt wird.
• Somit ergibt sich eine Gesamtkomplexität von: $$\mathcal{O}\left(S_x\times S_y\times \left(Z\times {\left(M\times N\right)}^2\right)\right)=\mathcal{O}\left(Z\times \left(S_x\times S_y\times {\left(M\times N\right)}^2\right)\right)\approx \mathcal{O}\left(Z\times n^6\right)=\mathcal{O}\left(n^6\right)$$

  
• Hierbei kann Z als konstanter Faktor, der nicht von der Komplexität des konkreten Simulationsproblems abhängt, ausgelassen werden, da die verschiedenen Fokusebenen nur erfordern, dass die Simulation Z mal durchgeführt wird.
• Die benötigten Rechenoperationen sind skalare Multiplikationen und skalare Additionen sowie die bereits erwähnten Fourier-Transformationen, welche den rechenzeitintensivsten Teil des Problems darstellen.
• Eine mögliche Ablaufsteuerung (englisch: control flow) dieses Algorithmus ist in 8B gezeigt.
• Der in 8B gezeigte Algorithmus weist keine Verzweigungen in Abhängigkeit von Bedingungen bei den Berechnungen auf.
• In 8B wird die Abhängigkeit einer Implementierung des Pupil-Shift-Verfahrens, beispielsweise gemäß der 8A, für einen Lichtquellen-Punkt gezeigt. Die Berechnung von voneinander nicht abhängigen Lichtquellen-Punkten kann gleichzeitig, beispielsweise parallel, durchgeführt werden, da die Lichtquellen-Punkte voneinander unabhängig sind.
• Das gezeigte Verfahren umfasst folgende Schritte:
• Bei 810 wird die Verschiebung bestimmt.
• Bei 820 erfolgt die Multiplikation von T mit K.
• Bei 830 wird eine inverse Fourier-Transformation durchgeführt.
• Bei 840 wird das Bild durch Akkumulation der Ergebnisse der Fourier-Transformationen durchgeführt.
• Die Schritte bei 810, 820 und 840 beinhalten nur skalare Matrixoperationen. Somit können die jeweiligen Schritte gleichzeitig, beispielsweise parallelisiert, durchgeführt werden.
• Ebenfalls ist festzustellen, dass die Eingabedaten für den Algorithmus nicht neu geladen werden müssen, sondern für alle Lichtquellen-Punkte verwendet werden können. Die Speichertransferkosten (englisch: memory transfer costs) der Eingangsdaten können hergeleitet werden, indem die Dimensionen eines typischen Datensatzes betrachtet werden. Diese sind gegeben durch: $$\text{N}=\text{M}=\mathrm{300,}\text{ Z}=7\text{ und S}=17000.$$

  
• Geht man von Gleitkommazahlen mit doppelter Präzision aus, kann hieraus der insgesamt benötigte Speicher wie folgt angegeben werden: $$\text{sizeof}\left(J\right)=\text{sizeof}\left(double\right)\times S=\mathrm{136,000}bytes$$

  

  $$\text{sizeof}\left(L\right)=2\times \text{sizeof}\left(double\right)\times S=\mathrm{272,000}bytes$$

  

  $$\text{sizeof}\left(T\right)=2\times \text{sizeof}\left(double\right)\times N\times M=\mathrm{1,440.000}bytes$$

  

  $$\text{sizeof}\left(K\right)=2\times \text{sizeof}\left(double\right)\times N\times M\times Z=\mathrm{10,080.000}bytes$$

  
• Somit ergeben sich circa 11,4 MiB benötigte Eingabedatengröße und sizeof(K)/2 Bytes für das Referenz-Bild I.
• Wie im Zusammenhang mit den 8A und 8B erläutert, ist das Problem der Referenz-Bild-Simulation aufgrund folgender Eigenschaften gut geeignet, um parallelisiert berechnet zu werden, beispielsweise auf einer GPU oder mehreren GPUs:
• • vernachlässigbare Speichertransferkosten,
• • zahlreiche (größtenteils skalare) Operationen, die auf die Eingabedaten angewendet werden und
• • eine einfache, zugleich hochgradig parallelisierbare Ablaufsteuerung.
• Eine mögliche Vorgehensweise zur parallelisierten Berechnung wird nachfolgend näher beschrieben.
• Die 8C zeigt einen Pseudocode für eine mögliche Vorgehensweise zum parallelen Berechnen eines Luftbildes.
• Der Pseudocode der 8C kann auf einem Grafikprozessor durchgeführt werden, aber auch andere parallele Implementierungen sind möglich.
• Es werden die folgenden Schritte ausgeführt:
1. 1. Eine Schleife durchläuft alle Lichtquellen-Punkte S, die auf einer Rechenvorrichtung gespeichert sind.
2. 2. Zwischenergebnisse werden in einem Speicher, auf den die Variable fft_temp verweist, zwischengespeichert. Der Speicher kann beispielsweise der globale Grafikkartenspeicher MEM-GPU der 7 sein.
3. 3. Es werden Z zweidimensionale inverse Fourier-Transformationen durchgeführt. Die Berechnungen können beispielsweise parallelisiert auf einer Grafikprozessoreinheit durchgeführt werden, beispielsweise auf SM0 und SM1 der 7 und/oder auf SM0 und SM1 von GPU-1 oder auf SM0-SM3 der GPU-2 der 6.
4. 4. Das rücktransformierte Zwischenergebnis in fft_temp wird von der Funktion accumulate_image akkumuliert, beispielsweise aufaddiert, um das Referenz-Luftbild zu erhalten, beispielsweise nachdem die Schleife vollständig durchlaufen wurde.
• Die hier beschriebene Vorgehensweise kann den Vorteil haben, dass der Speicherbedarf gering ausfallen kann. Beispielsweise kann der Speicherbedarf der Variable fft_temp auf N * M * Z* 1 Element begrenzt werden, wie zuvor im Zusammenhang mit 8B beschrieben.
• Bei einer Implementierung mit einer Grafikkarte kann es notwendig sein, drei Kontext-Wechsel (englisch: context switches) zwischen der CPU und dem GPU durchzuführen. Für große Werte der Anzahl von Lichtquellen-Punkten S kann der Overhead dieser Kontext-Wechsel dominant werden und daher diese Vorgehensweise ineffizient machen.
• In solchen Fällen können andere Vorgehensweisen, wie nachfolgend beschrieben, schnellere Berechnungen ermöglichen. Dies wird unten als die weiter verbesserte Vorgehensweise beschrieben.
• 9 zeigt einen Kernel, der eine skalare Matrixmultiplikation gemäß Gleichung (2.6) implementiert.
• Ein Kernel ist ein Programm, welches auf einer Rechenvorrichtung, beispielsweise einer Grafikkarte, parallele Berechnungen durchführt, um ein bestimmtes Problem zu lösen.
• So kann ein Kernel beispielsweise eine Luftbildgebungssimulation nach Hopkins durchführen, ein anderer Kernel kann beispielsweise analytisch ein Spektrum der Fotomaske berechnen.
• Die äußeren beiden Schleifen, die verwendet werden, um das Gitter zu durchlaufen, werden als Grid-Stride-Schleifen (englisch: grid-stride-loops) bezeichnet. Diese können ein gegebenes Thread-Gitter durchlaufen, wie im Zusammenhang mit 5 erläutert. Ist das Gitter groß genug, um einen gesamten Datensatz, welcher zur Berechnung benötigt wird, zu speichern, werden diese Schleifen nur einmal durchlaufen. In anderen Fällen gewährleisten die Schleifen, dass bei anderen Hardwarekonfigurationen die entsprechende Anzahl von Berechnungsschritten durchgeführt wird, wie beispielsweise oben im Zusammenhang mit 6 für die Fälle A und B erläutert.
• Die innere Schleife führt den Shift von T entsprechend den Beleuchtungsbedingungen durch. Hierbei erhält jeder Thread das geshiftete Element von T durch Shiften des Indexes mit dem auf T zugegriffen wird. Um implementierungsbedingtes Verhalten zu vermeiden, falls die Operatoren der verbleibenden Operationen im Vorzeichen in den Zeilen 10 und 11 von 9 abweichen, wird der Dividend durch Addition der Variablen size.x und size.y dazu gezwungen, jeweils positiv zu sein.
• Das geshiftete Element wird mit dem jeweiligen Element jeder Fokusebene der diskreten Pupillenfunktion K multipliziert. Das komplexe Produkt wird in dem Puffer fft_temp abgelegt.
• Hierbei bleibt der Speicherzugriff auf T zusammenhängend (engl.: coalesced), weil angrenzende Threads auf angrenzende Elemente von T zugreifen.
• Um die inverse Fourier-Transformation zu berechnen, kann beispielsweise die cuFFT Bibliothek, die Bestandteil des CUDA Softwareentwicklungs-Sets (engl.: Software Development Kit, SDK) ist, verwendet werden. Aber auch andere Bibliotheken zur Fourier-Transformation, beispielsweise zur schnellen Fourier-Transformation (englisch: fast fourier transform, FFT), können verwendet werden.
• 10 zeigt eine mögliche Implementierung des Kernels, welcher die Ausgabe der schnellen Fourier-Transformation zu dem Referenz-Bild akkumuliert.
• Die zwei äußeren Schleifen von 10 sind die bereits beschriebenen Grid-Stride-Schleifen.
• Die zwei inneren Schleifen berechnen das Quadrat der Betragswerte des Ergebnisses der komplexen schnellen Fourier-Transformation für jede Fokusebene und addiert das jeweilige Ergebnis zu dem Referenz-Luftbild für die jeweilige Fokusebene. In jeder Fokusebene verwaltet ein Thread ein Element von fft_input.
• Die Kernel, die im Zusammenhang mit den 9 und 10 beschrieben wurden, können Schleifen über Puffer laufen lassen, welche die Matrizen enthalten, die in eindimensionalem Speicher abgelegt werden.
• In den gezeigten Beispielen erfolgt bei beiden Kerneln der 9 und 10 die Iteration über die Matrix Zeile für Zeile mittels zweier Schleifen.
• Wird eine solche Implementierung beispielsweise mit CUDA vorgenommen, können beispielsweise 32 Threads gleichzeitig auf einer GPU ablaufen. In Fällen, in denen die Matrizen kein Vielfaches von 32 sind, kann dies dazu führen, dass die letzten Threads einer oder beider Schleifen untätig sind, beispielsweise wenn die Bedingung x < size.x nicht erfüllt ist. Würden die Threads auf die Elemente zugreifen, wären die Zugriffe jedoch wohldefiniert, da die Matrizen in linearen Puffern abgelegt werden. Die zweite Schleife kann daher als eine unnötiger weise divergierender Zweig (englisch: diverging branch) betrachtet werden. Die Vorgehensweise kann mithin weiter dadurch verbessert werden, dass die beiden äußeren Schleifen in eine einzige Schleife umgeformt werden.
• Hierdurch kann die benötigte Rechenzeit weiter verringert werden.
• Der Code in 10 führt die Summierung des Referenz-Luftbildes I im globalen Speicher durch. Der „+=“ Operator kann ein Read-Modify-Write Operator sein. In den globalen Speicher schreiben kann zeitintensiv (engl.: „expensive“) sein und kann dazu führen, dass ein Cache, der I enthält, ungültig wird. In Fällen von ausgerichtetem (engl.: unaligned) oder unzusammenhängendem Zugriff kann der nachfolgende Lesezugriff unter dem Zugriff auf nicht-zwischengespeicherten (englisch: uncached) globalen Speicher leiden. Hierdurch kann die Rechenzeit verlängert werden.
• Weiter verbesserte Vorgehensweise: Parallelisierte Iteration über die Lichtquellen-Punkte
• Zur weiteren Verbesserung kann die Iteration über die Lichtquellen-Punkte ebenfalls parallelisiert durchgeführt werden, beispielsweise auf einer GPU. Dies kann Vorteile wie das Vermeiden unnötiger Kernel-Aufrufe und/oder verringertem Overhead haben und somit die benötigte Rechenzeit weiter reduzieren.
• Die vorherstehend beschriebenen Verfahren können beispielsweise sieben Kernel benötigen, um einen Lichtquellen-Punkt abzubilden. So kann mul_circshift3D_kernel und accumulate_image_kernel jeweils einen Kernel und cuFFT bis zu fünf Kernel ausführen. Somit würden 17000 Lichtquellen-Punkte 119000 Kernel-Starts herbeiführen. Somit kann es vorteilhaft sein, jede Schleife für jeden Lichtquellen-Punkt auf einem GPU oder einen anderen Vielzahl von Recheneinrichtungen durchführen zu können, um die Anzahl von Kernel-Aufrufen und hierdurch beispielsweise den benötigten Overhead zu reduzieren.
• Hierbei kann dieser Ansatz jedoch das Speichern von fft_temp in den zugehörigen Speicher der GPU, beispielsweise MEM-GPU der 7, erforderlich machen.
• Um alle Transformationen innerhalb eines Durchlaufs durchführen zu können, kann es erforderlich sein, dass der Puffer die Größe des Speicherbedarfs der Elemente N*M*Z*S, wie zuvor im Zusammenhang mit der 8B beschrieben, hat. Dies kann dazu führen, dass die Speichergröße von fft_temp in der Größenordnung mehrerer Gigabyte, beispielsweise bei 80 GB liegen kann.
• Hierbei ist zu beachten, dass diese Speicheranforderungen höher sind als die Hardwarefähigkeiten gegenwärtig verfügbarer Grafikkarten. Auch kann eine Simulation einer größeren Anzahl von Lichtquellen-Punkten notwendig und/oder wünschenswert sein, so dass auch andere Beispiele einer Vielzahl von Recheneinrichtungen, die einen größeren Speicher, beispielsweise mit mehr als 100 GB aufweisen, an ihre Grenzen stoßen.
• Daher kann es vorteilhaft sein, das Problem zu unterteilen (englisch: slice), um mit einem gegebenen globalen GPU Speicher, beispielsweise MEM-GPU der 7, zu arbeiten.
• Hierzu können die Kernel mul_circshift3d_kernel und accumulate_image_kernel in einen eindimensionalen Thread-Block umgeformt werden, indem die äußeren beiden der Gitter-Stride-Schleifen in eine einzelne Schleife umgeformt werden. Wie zuvor erläutert, können hierdurch unnötige divergierende Zweige vermieden werden.
• Da in dieser Variante beide Kernel über alle Lichtquellen-Punkte auf der GPU iterieren, wird bei jeder Iteration nur jeweils ein Element L(i) der Lichtquellen-Punkte L und J von allen Threads benötigt. Für jedes Element wird eine vollständige Cache Line der GPU geladen, jedoch nur ein einzelnes Element verwendet.
• Greift das Programm auf ein einzelnes Element im Speicher oder auch auf nur eine einzelne Speicherzelle zu, so kann in manchen Ausführungsbeispielen nicht nur dieses einzelne Element oder diese einzelne Zelle geladen werden, sondern ein bestimmter Anteil an benachbartem Speicher gleich mit. Ziel kann es hierbei sein, diesen Speicher, der quasi ohne zusätzliche Kosten mitgeladen wurde, in unmittelbar folgenden Rechenschritten mit zu benutzen.
• Hierdurch kann in manchen Ausführungsbeispielen vermieden werden, dass im nächsten Rechenschritt an eine ganz andere Stelle im Speicher gesprungen wird und somit der „kostenlos“ geladenen Speicher, beispielsweise die Cache Line, verworfen wird. Die Größe einer Cache Line ist hardwareabhängig und kann eine Potenz von 2 sein. Typische Größen bei gegenwärtigen Rechenvorrichtungen sind 64 Byte, manchmal auch 128 Byte.
• Um das Verwerfen des Rests der Cache Line zu vermeiden, kann geteilter Speicher, beispielsweise SM-MEMO und SM-MEM1 der 7, verwendet werden, um die Daten zwischenzuspeichern. Um eine bessere Nutzung des geteilten Speichers zu gewährleisten, können die verwendeten Datentypen, beispielsweise der Datentyp von L, optimiert werden. Beispielsweise kann L als short2 abgelegt werden. Dies kann möglich sein, wenn für alle Koordinaten der Lichtquellen-Punkte gilt: $$\le \sqrt{N\times M}$$

  
• Im Falle von vorzeichenbehafteten 16 bit Zweierkomplement-Ganzzahlen, können die Koordinaten im Interval [-32768; 32767] liegen.
• Das Iterieren über alle Lichtquellen-Punkte auf der GPU kann den weiteren Vorteil haben, dass der accumulate_image_kernel von 10 die Summierung von I in Thread-Registern anstelle von globalem Speicher durchführt. Hierdurch kann der Beitrag von allen Lichtquellen-Punkten in Registern akkumuliert, beispielsweise aufaddiert, werden und zu dem Referenz-Bild nur einmalig im globalen Speicher akkumuliert und/oder kopiert werden. Hierdurch kann das Verfahren weiter beschleunigt werden, da zu einem größeren Anteil der schnellere Speicher SM-MEMO und SM-MEM1 der 7 verwendet wird.
• Die 11A und 11B zeigen benötigte Berechnungszeiten verschiedener Implementierungen gemäß verschiedener Ausführungsbeispiele als Funktion der Anzahl von Lichtquellen-Punkten.
• In 11B sind die entsprechenden Zahlenwerte zu den in 11A grafisch dargestellten Zusammenhängen angegeben.
• Eine Kurve 601 zeigt eine nicht-parallelisierte Berechnung, wie sie aus dem Stand der Technik bekannt ist. Im gezeigten Ausführungsbeispiel erfolgt die Berechnung mittels eines CPUs ohne Parallelisierung. Eine Kurve 602 zeigt die zuvor beschriebene parallelisierte Vorgehensweise, und eine Kurve 603 zeigt die zuvor beschriebene weiter verbesserte Vorgehensweise und eine Kurve 604 die zuvor beschriebene weiter verbesserte Vorgehensweise, ergänzt um die Verwendung von geteiltem Speicher, wie ebenfalls zuvor beschrieben, gemäß verschiedener Ausführungsbeispiele. Wie aus 11A zu ersehen ist, ist durch die hier beschriebene Parallelisierung eine erhebliche Beschleunigung um einen Faktor >10 bis hin zu einem Faktor 52 für die Referenz-Luftbild-Berechnung erzielt.
• Wie aus den Kurven der 11A weiter zu ersehen ist, skalieren die Implementierungen proportional mit der Größe des Problems. Dies wird dadurch erreicht, dass bei den beschriebenen Vorgehensweisen kein Kommunikations- oder Synchronisations-Overhead benötigt wird und alle Berechnungen gleichzeitig durchgeführt werden können. Die weiter verbesserte Implementierung ist gegenüber der verbesserten Implementierung um einen Faktor 3 schneller.
• Die Verwendung des geteilten Speichers in der weiter verbesserten Vorgehensweise zeigt sich in diesem Ausführungsbeispiel jedoch als nicht-vorteilhaft. Im vorliegenden Beispiel führt das zuvor beschriebene Slicen dazu, dass Synchronisationen auf der Ebene des Thread-Blöcke notwendig werden. Dieser benötigte Synchronisations-Overhead kann dazu führen, dass diese Vorgehensweise geringfügig langsamer sein kann, als die weiter verbesserte Vorgehensweise an sich. Somit kann es in manchen Situationen vorteilhaft sein, auf die Verwendung des geteilten Speichers, wie zuvor beschrieben, zu verzichten.
• Die 12A und 12B zeigen eine Güte der parallelisierten Berechnungen gemäß verschiedener Ausführungsbeispiele.
• In den 12A und 12B werden hierbei Vorgehensweisen, welche die TCC-Matrix, wie zuvor beschrieben, verwenden, untersucht.
• In 12A ist auf einer ersten y-Achse die Berechnungszeit als Funktion der Zeit und auf einer zweiten y-Achse eine normierte Intensitätsdifferenz zwischen dem berechneten Bild und dem parallelisiert berechneten Bild als eine Funktion von verwendeter Anzahl von Eigenwerten N der TCC-Matrix aufgetragen.
• Die Kurve 701 zeigt hierbei, dass die Berechnungszeit mit der Anzahl von Eigenwerten die für das parallelisierte Berechnungsverfahren zur Ermittlung der TCC-Matrix herangezogen werden, ansteigt. Die Kurve 702 zeigt die mittlere Abweichung zwischen einem mit maximaler Präzision gerechneten Referenz-Luftbild und einem mit reduzierter Anzahl von Eigenwerten N berechneten Referenz-Luftbild. Die Kurve 703 zeigt die maximale Abweichung. Wie aus der 12A hervorgeht, kann N entsprechend den Anforderungen der Anwendung und der benötigten Genauigkeit gewählt werden. Die Zahlenwerte der Kurven der 12A sind in 12B aufgetragen.
• Als ein Ausführungsbeispiel kann N = 60 gewählt werden, wenn eine Genauigkeit von 10^-1 als Maximalabweichung 703 für die Anwendung akzeptabel ist. Dies kann zu einer erheblichen Steigerung der Berechnungsgeschwindigkeit führen.
• Die Anzahl der Eigenwerte kann je nach Anwendungsfall gewählt werden. Es kann vorteilhaft sein, die Anzahl der Menge von Eigenwerten in Abhängigkeit der Anzahl von verwendeten Lichtquellen-Punkte zu wählen. Hierbei kann a > 5 sein. Hierbei kann a > 5, beispielsweise a > 10, beispielsweise a > 15, beispielsweise a > 20, beispielsweise a > 25, a > 30 sein.
• Alternativ oder ergänzend können zur Beschreibung einer Maske zwischen 50 und 1000, insbesondere zwischen 100 und 500 Eigenwerte verwendet werden.
• Die 13A und 13B zeigen Geschwindigkeitsvorteile der parallelen Berechnung mit einfacher Genauigkeit.
• In manchen Ausführungsbeispielen, bei denen gemeinsamen Speicher verwendet wird, können zehn Gleitkommazahlen im Thread-Register eingespart werden. Hierdurch kann die Anforderung an diese Register, beispielsweise Register der SKPs in 7, verringert werden. Hierdurch können beispielsweise nur 63 Register pro Thread benötigt werden. Hierdurch kann es möglich sein, ein Maximum von 1024 threads pro Block zu verwenden.
• In Ausführungsbeispielen, die den gemeinsamen Speicher nicht oder nicht zu diesem Zweck verwenden, kann es vorteilhaft sein, die Anzahl von Threads pro Block zu reduzieren. Sonst könnte sog. register spilling auftreten und die Laufzeit negativ, beispielsweise um einen Faktor 100 gesteigert, beeinflussen.
• Um unnötige Multiplikationen während der Laufzeit zu vermeiden, kann es vorteilhaft sein, Klammern von Berechnungen geschickt zu wählen. Gleitkommaberechnungen sind nicht assoziativ, beispielsweise gilt: $$\left(2\cdot \pi \right)\cdot x\ne 2\cdot \left(\pi \cdot x\right)$$

  
• Um einem optimierendem Kompiler eine Konstantenpropagation (englisch constant folding) zu ermöglichen, kann es vorteilhaft sein, Klammern um Ausdrücke wie $$2\cdot \pi$$

  

  zu setzen.
• In einigen Ausführungsbeispielen verschlechtert die atomare Addition der Spektren einzelner Dreiecke zu dem globalen Spektrum die Performance nicht.
• In einigen Ausführungsbeispielen ist die Berechnung der Sinus- und Cosinus-Terme am zeitaufwändigsten. Daher kann es vorteilhaft sein, die Sinus- und Cosinus-Terme nicht separat zu berechnen, sondern, beispielsweise durch die CUDA spezifische sincosf()-Funktion, beide Terme gemeinsam in einem einzelnen Rechenschritt zu berechnen.
• In manchen Ausführungsbeispielen wird die Berechnung von W(a_nm, b_nm) gemäß Gleichung (2.19) mittels IF-Blöcken implementiert. Dies kann zu Verzweigungen führen. In Ausführungsbeispielen, bei denen die Vielzahl von Recheneinrichtungen eine sogenannte Single Instruction, Multiple Data-, SIMD-, Struktur aufweisen, wie dies beispielsweise bei GPUs der Fall sein kann, können solche Verzweigungen die Performance negativ beeinflussen, wenn die Ausführungen verschiedener Threads verzweigt.
• Einige der Ausführungsbeispiele zeigen jedoch keine negativen Effekte, da die Fälle 1 und 2 von W(a_nm, b_nm) nur sehr unwahrscheinlich wahr sind. Daher führen die Verzweigungen in einigen Ausführungsbeispielen nicht zu flaschenhalsartigem Verhalten.
• in 13A ist die Berechnungszeit als Funktion der Anzahl von Dreiecken aufgetragen. Die Kurve 801 zeigt die Abhängigkeit für ein nicht-paralleles Rechenverfahren aus dem Stand der Technik. Die Kurve 802 zeigt den Zusammenhang für Berechnungen gemäß einem Ausführungsbeispiel mit Werten mit einfacher Genauigkeit, die Kurve 803 ein anderes Ausführungsbeispiel für Werte mit doppelter Genauigkeit. Wie aus der 13A und den tabellierten Zahlenwerten der 13B zu ersehen ist, kann durch eine parallelisierte Berechnung mit doppelter Genauigkeit eine Verkürzung der Berechnungszeit um zwei Größenordnungen erreicht werden. Durch Verwendung einer Berechnung mit einfacher Genauigkeit kann nochmals eine Steigerung um mehr als eine Größenordnung erreicht werden.
• In einem Ausführungsbeispiel mit einfacher Genauigkeit, vgl. Kurve 802, wird für 62.656 Dreiecke ein Intensitätsunterschied zu einer Implementierung aus dem Stand der Technik von -1,892837 × 10^-6 + i1,046527 × 10^-5 erreicht.
• Hierbei tritt das Maximum des Intensitätsunterschiedes in der Nähe der Frequenz von 0 auf.
• In einem Ausführungsbeispiel mit doppelter Genauigkeit wird für gleiche Bedingungen eine Intensitätsunterschied zu einer Implementierung aus dem Stand der Technik von 4,67629401697351 × 10^-10 -i4,43115737363925 × 10^-12 erreicht, vgl. Kurve 803.
• In manchen Ausführungsbeispielen kann die Genauigkeit des Ausführungsbeispiels mit einfacher Genauigkeit hierbei ausreichend sein. In manchen dieser Ausführungsbeispiele propagiert der Fehler nicht zu dem Ergebnis Referenz-Luftbild durch.
• Selbstverständlich können die Merkmale der vorab beschriebenen Ausführungsbeispiele und Aspekte der Ausführungsbeispiele miteinander kombiniert werden. Insbesondere können die Merkmale nicht nur in den beschriebenen Kombinationen, sondern auch in anderen Kombinationen oder für sich genommen verwendet werden, ohne das Gebiet der Ausführungsbeispiele zu verlassen.
• Die hier beschriebenen Vorgehensweisen können verwendet werden, um akkurater Bildgebungssimulation jeder beliebigen unendlich dünnen Struktur durchzuführen, um Referenz-Luftbilder zu bestimmen.
• Hierbei wurden Vorrichtungen, Verfahren und Computerprogramme bereitgestellt, die die zugrundeliegenden partiell kohärenten Bildgebungssimulation hochparallelisiert implementiert und auf aktuell verfügbaren Grafikkarten lauf- und skalierfähig bereitgestellt wird. Hierdurch können Simulation um Faktor 54 gegenüber herkömmlichen Implementierungen beschleunigt werden.
• In Fällen, beispielsweise bei Seriensimulationen, bei denen das zugrundeliegende optische System nicht verändert wird, erlaubt die beschriebene Optimierung mittels der TCC-Matrix eine weitere deutliche Beschleunigung der Simulation um Faktor 20.
• Die beschriebenen Vorgehensweisen ermöglichen es, die Hardwarespezifikationen zu variieren, ohne die Funktion der vorgestellten Lösung zu beeinträchtigen.
• Die hier beschriebenen Vorrichtungen, Computerprogramme und Verfahren lassen sich jedoch auf jedes beliebige zweidimensionale (d.h. unendlich dünne) Objekt anwenden.

Claims (14)

1. Verfahren (1000) zum Bestimmen einer Qualität einer Maske einer Fotolithographie-Vorrichtung, umfassend: Empfangen eines Designs (1010) der Maske, Empfangen von optischen Eigenschaften (1020) der Fotolithographie-Vorrichtung, paralleles Berechnen (1030), mit einer Vielzahl von Recheneinrichtungen, eines Referenz-Luftbildes basierend auf dem Design und den optischen Eigenschaften, wobei das parallele Berechnen (1030) des Referenz-Luftbildes basierend auf dem Pupil-Shift-Verfahren erfolgt, wobei das Pupil-Shift-Verfahren für eine parallele Berechnung mit der Vielzahl von Recheneinrichtungen implementiert ist, Empfangen (1040) eines gemessen Luftbildes (400), Vergleichen (1050) des Referenz-Luftbildes (200) mit dem gemessenen Luftbild (400), Bestimmen (1060) einer Qualität der Maske (500) basierend auf dem Vergleichen.
2. Verfahren nach Anspruch 1, wobei das Pupil-Shift-Verfahren umfasst: eine gemeinsame Intensitätsfunktion und eine Maskentransmissionsfunktion, und eine diskrete Pupillenfunktion.
3. Verfahren nach Anspruch 2, wobei das Bestimmen der gemeinsamen Intensitätsfunktion und/oder der Maskentransmissionsfunktion basierend auf mindestens einer Messung der Messvorrichtung erfolgt.
4. Verfahren nach Anspruch 2 oder 3, wobei die diskrete Pupillenfunktion und die gemeinsame Intensitätsfunktion als eine Transmissions-Kreuzkoeffizientenmatrix, dargestellt werden.
5. Verfahren nach Anspruch 4, wobei das Verfahren umfasst: Bestimmen einer Vielzahl von Eigenwerten der Transmissions-Kreuzkoeffizientenmatrix, Auswählen einer Menge von Eigenwerten aus der Vielzahl von Eigenwerten, wobei die Anzahl der Menge von Eigenwerten kleiner ist als die Anzahl der Vielzahl von Eigenwerten, Durchführen des Verfahrens mit der Menge der Eigenwerte.
6. Verfahren nach Anspruch 5, wobei die Anzahl der Menge von Eigenwerten um einen Faktor a kleiner ist als eine Anzahl von verwendeten Lichtquellen-Punkten, wobei a > 5 ist.
7. Verfahren nach einem der Ansprüche 4-6, wobei das Verfahren umfasst: Berechnen einer Vielzahl von Referenz-Luftbildern für eine zugehörige Vielzahl von Fokussiereinstellungen, FE, der Lithographie-Vorrichtung, wobei die Transmissions-Kreuzkoeffizientenmatrix jeweils konstant gehalten wird und das Berechnen der Vielzahl von Referenz-Luftbildern ferner umfasst: Bestimmen einer Vielzahl von diskreten Pupillenfunktionen, K(FE), wobei jede diskrete Pupillenfunktion der Vielzahl von diskreten Pupillenfunktionen jeweils einer Fokussiereinstellung der Vielzahl von Fokussiereinstellungen zugeordnet ist, und jeweiliges Durchführen des Verfahrens gemäß einem der Ansprüche 4-6 für die Vielzahl der diskreten Pupillenfunktionen.
8. Verfahren nach Anspruch 7, ferner umfassend: Darstellen jeder der diskreten Pupillenfunktionen der Vielzahl von diskreten Pupillenfunktionen jeweils als Matrix, Laden, in einen zugehörigen Speicher, der Vielzahl von diskreten Pupillenfunktionen und Berechnen der Vielzahl von Referenz-Luftbildern ohne Nachladen von Matrizen in den zugehörigen Speicher.
9. Verfahren nach einem der Ansprüche 1-8, wobei das Verfahren umfasst: Darstellen des Designs der Maske als eine Menge von k Polygonen, Durchführen einer analytischen Spektrumsberechnung für die Menge von k Polygonen, Bestimmen der Maskentransmissionsfunktion basierend auf der analytischen Spektrumsberechnung für Polygone mittels Fourier-Transformation.
10. Verfahren nach Anspruch 9, wobei die Polygone Dreiecke sind.
11. Verfahren nach einem der Ansprüche 1-10, wobei das Pupil-Shift-Verfahren ein oder mehrere Multiplikationen von Matrizen umfasst, wobei mindestens eine der Multiplikationen von Matrizen mittels Skalarmultiplikation durchgeführt wird
12. Verfahren nach einem der Ansprüche 1-11, wobei das Verfahren mit Werten einfacher Genauigkeit in der Vielzahl von Recheneinrichtungen durchgeführt wird.
13. Verfahren nach einem der Ansprüche 1-12, wobei das gemessen Luftbild: in der Fotolithographie-Vorrichtung und/oder in einer Fotomaskenqualifikationsvorrichtung gemessen wurde.
14. Computerprogramm mit einem Programmcode, der ein Durchführen eines der Verfahren der Ansprüche 1-13 bewirkt, wenn der Programmcode auf einer Vielzahl von Recheneinrichtungen abläuft.

Images