Ein PDAE-Netzwerkmodell alsAbstraktes Differential-Algebraisches System

– Diplomarbeit –

Humboldt-Universitat zu Berlin

Mathematisch-Naturwissenschaftliche Fakultat II

Institut fur Mathematik

eingereicht von: Steffen Schulz

geboren am: 12.7.1976 in Berlin

Betreuerinnen: Prof. Dr. Roswitha Marz

Dr. Caren Tischendorf

Berlin, 30. September 2002

Inhaltsverzeichnis

1 Einleitung 1

2 Abstrakte Differential-Algebraische Systeme 4

2.1 Definitionen . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

2.2 Der Begriff des ADAS-Index . . . . . . . . . . . . . . . . . . . . 5

3 Die Netzwerkgleichung 8

3.1 Das elektrische Netzwerk als Graph . . . . . . . . . . . . . . . . 8

3.2 Die Modifizierte Knotenanalyse . . . . . . . . . . . . . . . . . . 10

3.3 Der Index der Netzwerkgleichung . . . . . . . . . . . . . . . . . 13

3.4 Zusammenfassung . . . . . . . . . . . . . . . . . . . . . . . . . . 15

4 Halbleitersimulation 17

4.1 Das Drift-Diffusionsmodell . . . . . . . . . . . . . . . . . . . . . 17

4.2 Das PDE-Modell . . . . . . . . . . . . . . . . . . . . . . . . . . 20

4.3 Bemerkungen zum PDE-Modell . . . . . . . . . . . . . . . . . . 23

5 PDAEs als Abstrakte Differential-Algebraische Systeme 25

5.1 Das gekoppelte System . . . . . . . . . . . . . . . . . . . . . . . 25

5.2 Das Abstrakte Differential-Algebraische System . . . . . . . . . 29

5.3 Aussagen zum ADAS-Index . . . . . . . . . . . . . . . . . . . . 31

5.3.1 ADAS ohne Driftkomponente . . . . . . . . . . . . . . . 32

5.3.2 ADAS mit Driftkomponente . . . . . . . . . . . . . . . . 35

5.3.3 Zusammenfassung . . . . . . . . . . . . . . . . . . . . . . 38

5.4 Beweise der Hauptresultate . . . . . . . . . . . . . . . . . . . . 40

5.5 Beispiele . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56

6 Thesen 61

A Funktionalanalytische Grundlagen 63

B Grundlagen aus der Linearen Algebra 68

C Verzeichnis der verwendeten Symbole 70

1 Einleitung

In der Industrie spielen integrierte Schaltkreise eine wesentliche Rolle. In prak-tisch allen elektronischen Produkten werden diese Schaltungen eingesetzt. ImRahmen von Entwicklung, Produktion und Test wird verstarkt auf die nume-rische Schaltungssimulation gesetzt, da durch die mathematische Modellierungnoch vor der eigentlichen Produktion Aussagen z.B. zum Zeitverhalten einerSchaltung gemacht werden konnen. Die Schaltungssimulation bietet eine ein-fache, effiziente Alternative, unterschiedliche Layouts der Schaltungen mitein-ander zu vergleichen und tragt so bei sehr geringen Kosten zu einer deutlichenOptimierung der Schaltungen und des Produktionsprozesses bei.

Bei der numerischen Schaltungssimulation haben sich vor allem Netzwerkan-satze bewahrt. Dabei gehen sowohl die Netzwerktopologie als auch die Bau-elemente der entsprechenden Schaltung in die mathematische Modellierung mitein. Die Netzwerktopologie wird durch den Netzwerkgraphen beschrieben. DasVerhalten der Bauelemente wird durch verschiedenste mathematische Modellewidergegeben.

Bauelemente, die durch algebraische Gleichungen oder Gewohnliche Differen-tialgleichungen beschrieben werden1, haben sich als relativ unproblematischerwiesen. Die Modellierung fuhrt auf Algebro-Differentialgleichungen2 (DAEs)

f(x(t), x′(t), t

)= 0

mit singularer Jacobi-Matrix f ′x. Die analytischen und numerischen Eigenschaf-ten solcher Gleichungen werden durch den Index charakterisiert (Kapitel 3.3).Bei Verwendung linearer Bauelemente hangt dieser Index ausschließlich vonder Netzwerktopologie ab (vgl. Satz 14).

Im Unterschied zu diesen Bauelementen werden Halbleiterbauelemente durchPartielle Differentialgleichungen3 (PDEs) beschrieben. Um auch Halbleiter mitAlgebro-Differentialgleichungen behandeln zu konnen, wird jedes Halbleiter-bauelement durch ein entsprechendes Kompaktmodell aus Widerstanden, Kon-densatoren, Spulen sowie Strom- und Spannungsquellen ersetzt (Abbildung 1).Dieses Kompaktmodell stellt dann praktisch eine Ortsdiskretisierung der PDEdar.

Mit fortschreitender Miniaturisierung der Bauelemente und immer geringerenEingangssignalen kann der Fehler, der bei der Modellierung mit Kompaktmo-dellen gemacht wird, in Zukunft nicht mehr vernachlassigt werden. Daruber

1Zu diesen Bauelementen zahlen Widerstande, Kondensatoren, Spulen aber auch un-abhangige Strom- und Spannungsquellen.

2engl.: Differential-Algebraic Equation (DAE); Meist werden Algebro-Differentialglei-chungen auch im deutschen Sprachgebrauch als DAEs bezeichnet.

3engl.: Partial Differential Equation (PDE)

1

elektrische Schaltungmit Halbleiterelement

Halbleiter-simulationmit PDEs

Halbleiter-simulationmit Ersatz-schaltbildern

Abbildung 1: Numerische Simulation elektrischer Schaltungen4

hinaus werden die Ersatzschaltbilder immer umfangreicher, so dass auch dieParameterextraktion fur die Kompaktmodelle schwieriger wird.

Als Ausweg bietet sich an, Halbleiterbauelemente direkt durch partielleDifferentialgleichungen zu beschreiben, also gekoppelte Systeme aus DAEs(Netzwerkgleichung) und PDEs (Bauelementsimulation) zu untersuchen. Sol-che gekoppelten Systeme werden als Partial Differential-Algebraic Equations(PDAEs) bezeichnet.

Die ersten Untersuchungen zu PDAEs in der Schaltungssimulation stammenbereits aus den 70er Jahren. Engl, Laur und Dirks entwickelten in Aachenden Simulator Medusa, der Teilsysteme unterschiedlicher Komplexitat lostund die Ergebnisse zur Behandlung umfangreicherer Schaltungen weiterver-wendet [3].Ende der 80er Jahre wurde mit Codecs von Mayaram an der Ber-keley University of California eine Simulatorkopplung vorgestellt [20]. Hierwird jeweils fur die Netzwerkgleichung und die Bauelementsimulation Spezi-alsoftware verwendet. Noch deutlicher wird dieser Ansatz bei der Kopplungdes Schaltungssimulators Spice mit dem Bauelementsimulator Pisces Mitteder 90er Jahre an der Stanford University [24]. Ein Uberblick uber bisherigeArbeiten wird auch in [23] gegeben.

Aktuell ist ein wachsendes Interesse an der mathematischen Untersuchung derPDAEs in der Schaltungssimulation zu erkennen [9, 13, 15]. Trotz der Existenzder o.g. Software ist die mathematische Struktur der gekoppelten Gleichungen

4Im Beispielnetzwerk aus Kondensator ( ), Halbleiterdiode ( ) und unabhangi-ger Spannungsquelle ( ) wird die Diode zunachst durch eine nichtlineare Stromquelle( ) und einen parallel geschalteten Kondensator ersetzt. Alternativ soll die Diodedurch ein ”Black-Box“-Element ( ) ersetzt werden, das durch PDEs beschrieben wird.

2

weitgehend ungeklart. So ist z.B. offen, in welcher Form der Indexbegriff derDAEs auf PDAEs erweitert werden kann. Auch Fragen der Existenz und Ein-deutigkeit von Losungen sind in Zukunft verstarkt zu untersuchen. In [1] konntebeispielsweise zunachst die eindeutige Losbarkeit des gekoppelten Systems fureinfache, eindimensionale Bauelemente im stationaren Fall gezeigt werden.

An der Humboldt-Universitat in Berlin wurde von Lamour, Marz und Ti-schendorf in [13] das Konzept derAbstrakten Differential-Algebraischen Syste-me (ADAS oder ADA-System) vorgestellt, das einen systematischen Zugangzum Studium von PDAEs bietet. Daruber hinaus kann im Rahmen der ADA-Systeme ein Index fur PDAEs, der ADAS-Index, definiert werden.

In der vorliegenden Arbeit werden PDAEs der Schaltungssimulation als Ab-strakte Differential-Algebraische Systeme untersucht. Es wird gezeigt, wie dasgekoppelte System aus Netzwerk- und Bauelementsimulation einheitlich alsADA-System geschrieben werden kann. Besonderes Augenmerk wird auf denBegriff des ADAS-Index gelegt. Ein Index fur das Gesamtsystem kann nur dannsinnvoll sein, wenn Eigenschaften der Schaltungen, wie sie aus der Modellie-rung mit Kompaktmodellen bekannt sind, korrekt widergespiegelt werden.

In Kapitel 2 wird das Konzept der Abstrakten Differential-AlgebraischenSysteme vorgestellt. Die Ausfuhrungen orientieren sich an [13]. Die Kapitel 3und 4 sind der Modellierung des Netzwerkes mit DAEs und der Bauelement-simulation mit PDEs gewidmet. Hier werden die Komponenten vorgestellt, diezusammen das PDAE-System bilden, das in Kapitel 5 als ADA-System unter-sucht wird. Im Abschnitt 5.3 sind die wesentlichen Aussagen zum ADAS-Indexzu finden. Die Beweise dieser Aussagen werden im Abschnitt 5.4 gefuhrt.

In der gesamten Arbeit werden Matrizen mit lateinischen Buchstaben A, Busw. bezeichnet. Fur Operatoren in unendlichdimensionalen Raumen werdendagegen Script-Buchstaben wie A, B usw. verwendet.

Funktionalanalytische Grundlagen, die zur Untersuchung dieser Operatorenverwendet werden, sind im Anhang A zusammengestellt. Anhang B enthaltAussagen der linearen Algebra, die haufig verwendet werden.

Das Gebiet der gekoppelten Netzwerk- und Bauelementsimulation ist noch sehrjung. Ich bin Frau Prof. R. Marz besonders dankbar, dass sie mich ermutigt hat,mich im Rahmen meiner Diplomarbeit mit einem so neuen und spannendenGebiet der Mathematik zu beschaftigen. Durch ihre engagierte Unterstutzungtrug sie sehr zum Gelingen der Arbeit bei. Gleichermaßen danke ich Frau Dr.C. Tischendorf fur die vielen Gesprache, die mir stets sehr geholfen haben.

3

2 AbstrakteDifferential-Algebraische Systeme

Partielle Differential-Algebraische Systeme (PDAEs) bestehen aus gekoppel-ten Systemen partieller Differentialgleichungen (PDEs) und Algebro-Differen-tialgleichungen (DAEs). Diese Systeme konnen im Rahmen der AbstraktenDifferential-Algebraische Systeme (ADAS) systematisch studiert werden.

In [13] werden ADA-Systeme als

A(t)(D(t)u(t)

)′+ B(t)u(t) = q(t), t ∈ I = [t0, T ] ⊂ IR, (1)

formuliert. Die Operatoren A(t) : Z → Y , D(t) : X → Z und B(t) : X → Ysind Abbildungen zwischen reellen Hilbert-Raumen. Losungen u : I → X sindPfade in dem (unendlichdimensionalen) Hilbert-Raum X.

Fur Systeme (1) kann ein Index definiert werden. Dabei dient die Definiti-on des Behandelbarkeitsindex fur Algebro-Differentialgleichungen als Orien-tierung [19].

2.1 Definitionen

Definition 1 H1 und H2 seien reelle Hilbert-Raume. Ein linearer OperatorL : H1 → H2 heißt

(i) normal losbar, wenn imL = imL, und

(ii) dicht losbar, wenn imL = H2 ist.

Mit M wird der Abschluss der Menge M ⊂ H im Hilbert-Raum H bezeichnet.

Das Symbol Lb(H1, H2) steht fur die Menge der linearen beschrankten Opera-toren H1 → H2. Abkurzend wird Lb(H) = Lb(H,H) geschrieben.

Abbildungen Q ∈ Lb(H), die Q2 = Q genugen, werden als Projektoren be-zeichnet. Projektoren sind also stets beschrankte lineare Abbildungen.

Q ist ein Projektor auf imQ langs kerQ. Mit Q ist P = I − Q ebenfalls einProjektor, der als zu Q komplementarer Projektor bezeichnet wird. Hier istI ∈ Lb(H) die Identitat.

X, Y und Z seien reelle Hilbert-Raume. Wir betrachten Abstrakte Differential-Algebraische Systeme (1) mit linearen Operatoren A(t) : Z→Y , D(t) : X→Z,B(t) : X → Y und rechter Seite q(t) ∈ Y . Die Operatoren A(t) und D(t)seien beschrankt und normal losbar. Der Operator B(t) sei auf einer dichtenTeilmenge DB ⊂ X definiert. Fur alle x ∈ DB sei B(·)x ∈ C(I, Y ) stetig.

Ein stetiger Pfad x ∈ C(I, X) mit x(t) ∈ DB ∀ t ∈ I heißt Losung von (1),wenn D x ∈ C1(I, Z) und (1) punktweise fur alle t ∈ I erfullt ist.

4

Definition 2 Der Hauptterm von (1) heißt proper formuliert, wenn die Ope-ratoren A und D in folgendem Sinne gut zusammenpassen:

(i) imA(t)D(t) = imA(t), kerA(t)D(t) = kerD(t), kerA(t)⊕ imD(t) = Z∀ t ∈ I.

(ii) Der Projektor R(t) mit imR(t) = imD(t) und kerR(t) = kerA(t)∀ t ∈ I hangt stetig differenzierbar von t ab.

Definition 3 Sei L ∈ Lb(H1, H2). Ein Operator L− ∈ Lb(H2, H1) heißt re-flexive verallgemeinerte Inverse von L, wenn

L−LL− = L− und LL−L = L.

Die Operatoren LL− ∈ Lb(H2) und L−L ∈ Lb(H1) sind Projektoren. In [13]wird gezeigt, dass fur normal losbare Operatoren L ∈ Lb(H1, H2) und beliebigeProjektoren R∈Lb(H2) und P∈Lb(H1) eine reflexive verallgemeinerte InverseL− ∈ Lb(H2, H1) konstruiert werden kann, fur die LL− = R und L−L = Pgilt. In diesem Fall ist L− eindeutig bestimmt.

2.2 Der Begriff des ADAS-Index

Es sei Q0(t) ∈ Lb(X) ein Projektor auf N0(t) = kerA(t)D(t), t ∈ I. Wie obenwird mit P0(t) = I −Q0(t) der zu Q0(t) komplementare Projektor bezeichnet.D(t)− : Z → X sei die reflexive verallgemeinerte Inverse mit

D(t)D(t)− = R(t) und D(t)−D(t) = P0(t).

Dabei ist R(t) der Projektor aus Definition 2. Mit G0(t) und B0(t) werden dieOperatoren

G0(t) = A(t)D(t) und B0(t) = B(t), t ∈ I,

bezeichnet. Fur i ≥ 0 definieren wir fur t ∈ I punktweise

Ni(t) = kerGi(t),

Qi(t) ∈ Lb(X), Qi(t)2 = Qi(t), imQi(t) = Ni(t), Pi(t) = I − Qi(t),

Wi(t) ∈ Lb(Y ), Wi(t)2 = Wi(t), kerWi(t) = imGi(t),

Si(t) = ker(Wi(t)Bi(t)

),

Gi+1(t) = Gi(t) + Bi(t)Qi(t),

Ci+1(t) = Gi+1(t)D(t)−(D(t)P1(t) · · · Pi+1(t)D(t)−

)′D(t),

Bi+1(t) = Bi(t)Pi(t)− Ci+1(t)P0(t)P1(t) · · · Pi(t).

5

Insbesondere Ci+1 kann nur formuliert werden, wenn die zur Definition desOperators verwendeten Ableitungen existieren. In den weiteren Uberlegungenwerden der Operator D und die Projektoren Pi konstant sein, so dass Ci+1 = 0fur alle i ≥ 0 folgt. Dann ist

G1(t) = G0(t) + B(t)Q0(t) und G2(t) = G1(t) + B(t)P0(t)Q1(t). (2)

Definition 4 Das Abstrakte Differential-Algebraische System (1) hat denIndex µ, wenn

(i) Wi(t) fur i = 0, . . . , µ − 1 stetig von t abhangt, dim(imWi(t)

)> 0 ist

und außerdem

(ii) der Operator Gµ(t) injektiv und dicht losbar ist.

Abschließend sollen die Begriffe am Beispiel der Warmeleitungsgleichung illus-triert werden [32].

Beispiel 5 Sei u(·, t) die Temperatur im Gebiet Ω ⊂ IR3 zur Zeit t.

j = −α grad u (3)

ist die Stromdichte des Warmeflusses in Ω. Werden Warmequellen mit f be-zeichnet, so gilt aus Grunden der Energieerhaltung die Beziehung

∂u

∂t+ div j = f in Ω× I. (4)

Wird (3) direkt in (4) eingesetzt, ergibt sich die klassische parabolische Warme-leitungsgleichung

∂u

∂t− α∆u = f.

Mit den Operatoren

A(t) =

(I0

), D(t) =

(I 0)

und B(t) =

(0 div

α grad I)

kann das gekoppelte System (3), (4) als Abstraktes Differential-AlgebraischesSystem

(I0

)((I 0) (

u(·, t)j(·, t)

))′+

(0 div

α grad I)(

u(·, t)j(·, t)

)=

(f0

)

geschrieben werden. Es sei

X = Y = L2(Ω)× [L2(Ω)]3, Z = L2(Ω).

6

Als Definitionsbereich des Operators B kann

DB = H1(Ω)× [H1(Ω)]3

gewahlt werden.

Es ist N0 = kerAD = 0 × [L2(Ω)]3, so dass Q0 =

0 00 I

ein Projektor5

auf N0 ist. Wegen G1 = AD + BQ0 = I div

0 I

mit G−11 =

I −div0 I

hat

das obige System den Index 1.

Soll die Temperatur u auf dem Rand ∂Ω vorgeschrieben werden, ist eineDirichlet-Randbedingung der Form

u = g auf ∂Ω× I (5)

zu fordern. Nach geeigneter Homogenisierung

u = u− u0

kann die Randbedingung im ADA-System durch die Wahl von

DB = H10 (Ω)× [H1(Ω)]3

realisiert werden. Dabei ist u0 ∈ C1(I, H1(Ω)

)mit u0 = g auf ∂Ω. Das ge-

koppelte System (3)-(5) wird durch das Abstrakte Differential-AlgebraischeSystem

(I0

)((I 0) (

u(·, t)j(·, t)

))′+

(0 div

α grad I)(

u(·, t)j(·, t)

)=

(f − ∂u0

∂t

−α grad u0

)

beschrieben.

5Die Identitat [H1(Ω)]k → [H1(Ω)]k wird fur k = 1 und k = 3 mit dem gleichen SymbolI bezeichnet.

7

3 Die Netzwerkgleichung

Bei der Entwicklung integrierter Schaltungen spielt die numerische Simulationeine wesentliche Rolle. Netzwerkmethoden haben sich als besonders geeig-net erwiesen, elektrische Schaltungen mathematisch zu beschreiben [5]. Da-bei werden auf Grundlage der Netzwerkstruktur der Schaltungen physikalischeGesetze wie Energie- und Ladungserhaltung mit den charakteristischen Glei-chungen der Bauelemente kombiniert. Die Topologie der Schaltung wird durchden Netzwerkgraphen beschrieben.

Diese Strategie wird zum Beispiel bei der Modifizierten Knotenanalyse6 (MNA)verfolgt, auf die im Abschnitt 3.2 kurz eingegangen wird.

Bei der MNA konnen die zu losenden Gleichungen automatisch aus den Netz-werkdaten generiert werden. Dabei muss aber im Allgemeinen mit einem Varia-blensatz gearbeitet werden, der fur das entsprechende System nicht minimalist. Die MNA fuhrt auf gekoppelte Systeme impliziter Differentialgleichun-gen und algebraischer Gleichungen, also auf Algebro-Differentialgleichungen(DAEs).

Die Theorie der Algebro-Differentialgleichungen ist aktueller Forschungsgegen-stand. Einblicke in die derzeitige Entwicklung geben Marz [18], Hairer undWanner [11], Brenan et al. [2], Strehmel/Weiner [27] und andere.

Bevor im Abschnitt 3.2 die Gleichungen der MNA formuliert werden, fasst derfolgende Abschnitt wichtige graphentheoretische Aussagen zusammen.

3.1 Das elektrische Netzwerk als Graph

Ein Graph besteht aus Knoten, die durch Zweige verbunden sind. Die Zwei-ge konnen dabei gerichtet oder auch ungerichtet sein. Pfade sind Tupel vonZweigen, die zwei Knoten verbinden [12].

Definition 6 Sei N 6= ∅ eine endliche Knotenmenge und B ⊂ N × N eineMenge von Zweigen7. Das Paar G = (N,B) wird als Graph bezeichnet. OhneEinschrankung wird N ⊂ IN angenommen.G heißt gerichteter Graph, wenn jeder Zweig b ∈ B mit einer Richtung ver-sehen ist8.

6engl.: Modified Nodal Analysis (MNA)7Im Allgemeinen sind Zweige ungeordnete Paare b ∈ N ×N .8G wird zum gerichteten Graph, wenn die Zweige b ∈ N ×N als geordnete Paare inter-

pretiert werden.

8

Im Folgenden werden ausschließlich Graphen betrachtet, bei denen alle Zweigeb = (n1, n2) jeweils zwei verschiedene Endpunkte n1 6= n2 besitzen.

Ein Tupel P = (b1, . . . , bk) ∈ Bk, k ∈ IN, heißt Pfad in G, wenn aufeinander-folgende Zweige bi, bi+1, i = 1, . . . , k− 1, durch einen gemeinsamen Knoten ni

verbunden sind und kein Knoten von G mehr als zwei Zweige aus P verbindet.

Ein Graph G = (N, B) heißt verbunden, wenn zwischen je zwei Knoten in Nein Pfad existiert.

Fur die folgenden topologischen Betrachtungen sind die Begriffe Schleife, Baumund Schnittmenge in einem Graphen von grundlegender Bedeutung.

Definition 7 Sei G = (N, B) ein verbundener Graph. Gt = (Nt, Bt) sei einTeilgraph von G, das heißt es sei Nt ⊂ N und Bt ⊂ Nt ×Nt ⊂ B.

(i) Gt heißt Schleife, wenn Gt verbunden ist und jeder Knoten von Gt genauzwei Zweige verbindet.

(ii) Gt heißt Baum, wenn Nt = N ist, Gt keine Schleifen enthalt und selbstein verbundener Graph ist.

(iii) Gt heißt Schnittmenge, wenn (N, B \ Bt) ein nicht-verbundener Graphist, aber fur beliebiges b ∈ Bt der Graph (N, (B \ Bt) ∪ b) verbundenist.

Wir bezeichnen mit n die Anzahl der Knoten von G. Um diese Knoten zueinem verbundenen Graphen zu verknupfen, sind mindestens n − 1 Zweigenotig. Jeder weitere Zweig wurde aber zu einer Schleife fuhren. Damit besitztjeder Baum in G genau n− 1 Zweige.

Die Lage der Zweige in G, also die Topologie des Netzwerkes, wird durch dieInzidenzmatrix beschrieben.

Definition 8 Sei G ein gerichteter, verbundener Graph mit n Knoten undk Zweigen. Die Matrix Aa ∈ MIR(n, k) mit Aa = (aij) heißt (vollstandige)Inzidenzmatrix von G, wenn

aij =

1, falls der Zweig j vom Knoten i wegfuhrt

−1, falls der Zweig j zum Knoten i hinfuhrt

0, sonst

.

Die Zeilen von Aa sind offenbar linear abhangig. Ohne Informationsverlustkann eine Zeile aus Aa gestrichen werden. Es ergibt sich so die reduzierte

9

Inzidenzmatrix A ∈ MIR(n−1, k). In [12] wird gezeigt, dass A vollen Zeilenrangbesitzt. In einem verbundenen Graphen ist damit stets

rang(A) = rang(Aa) = n− 1.

Der Knoten, der zu der gestrichenen Zeile gehort, wird Masseknoten genannt.

Im Folgenden ist mit dem Begriff der Inzidenzmatrix immer die reduzierteMatrix gemeint.

Das folgende Lemma gibt notwendige und hinreichende Bedingungen fur dieExistenz von Schleifen in einem Graphen an [12, 29].

Lemma 9 Sei A Inzidenzmatrix eines (nicht notwendig verbundenen) Gra-phen G. Der Graph enthalt genau dann Schleifen, wenn die Spalten von Alinear abhangig sind, d.h. wenn ker A 6= 0 ist.

Beweis: G bestehe aus r verbundenen Teilgraphen Gi. Die Zweige von

G seien o.E. so nummeriert, dass A =

(A1 0 ··· 00 A2 ··· 0...

......

...0 0 ··· Ar

)ist, wobei Ai Inzidenz-

matrix von Gi ist. Gi besitze ki Knoten. Dann ist rang(Ai) = ki − 1. (∗)Falls G eine Schleife enthalt, so gilt dies auch fur einen Teilgraphen Gi0 . Da Gi0

ein verbundener Graph ist, enthalt Gi0 folglich mindestens ki0 Zweige. Wegen(∗) sind die zugehorigen ki0 Spalten von Ai0 – und damit jene von A – linearabhangig.

Sind die Spalten von A linear abhangig, so sind auch die Spalten von Ai0 linearabhangig, i0 ∈ 1, . . . , r geeignet. Wurde Gi0 keine Schleife besitzen, so wareGi0 als verbundener Graph ein Baum, bestunde also aus genau ki0−1 Zweigen.Wegen der linearen Abhangigkeit der Spalten ware dann rang(Ai0) < ki0−1,was ein Widerspruch zu (∗) ist. Folglich enthalt Gi0 doch Schleifen. ¤

3.2 Die Modifizierte Knotenanalyse

Elektrische Netzwerke bestehen aus Bauelementen, die durch Knoten miteinan-der verbunden sind. Bei der mathematischen Beschreibung solcher Netzwerkewurden verschiedene Ansatze verfolgt [12, 31]. Das Netzwerk wird dabei stetsals verbundener Graph interpretiert. Die Netzwerkgleichungen ergeben sichaus der Netzwerktopologie und den mathematischen Modellen der Netzwerk-elemente. Als Vektor der Unbekannten werden Zweigstrome, Zweigspannungenund Knotenpotentiale9 gewahlt.

9Knotenpotentiale geben den Spannungsabfall gegenuber dem Masseknoten an.

10

Die Kirchhoffschen Gesetze geben Beziehungen zwischen den Zweigstromenund -spannungen an, die ausschließlich von der Netzwerktopologie abhangen.

(i) Kirchhoffsches Stromgesetz (KCL): Die Summe aller an einem Knotenzusammenfließenden Strome ist Null.

(ii) Kirchhoffsches Spannungsgesetz (KVL): Die Summe aller zu einer Schlei-fe gehorenden Spannungen ist Null.

Besteht das Netzwerk aus k Bauelementen, n Knoten und ist A ∈ MIR(n−1, k)die Inzidenzmatrix des Netzwerkgraphen, so folgen aus KCL und KVL dieBeziehungen

A i = 0 und v = AT e

(vgl. z.B. [12]). Hier ist i ∈ IRk der Vektor der Zweigstrome und v ∈ IRk bzw.e ∈ IRn−1 der Vektor der Zweigspannungen bzw. der Knotenpotentiale.

Diese Gleichungen werden durch die charakteristischen Gleichungen der Bau-elemente erganzt. Als Bauelemente wollen wir lineare zeitinvariante Wider-stande (R), Kapazitaten (C), Induktivitaten (L) und unabhangige Strom- undSpannungsquellen (I bzw. V ) zulassen. In diesem Fall spricht man auch vonRCL-Netzwerken.

Seien i =

(iRiCiLiIiV

)und v =

( vRvCvLvIvV

)die Strom- und Spannungsvektoren, wo-

bei iE der Vektor der Zweigstrome und vE der Vektor der Zweigspannungendurch bzw. an den Bauelementen E ∈ R, C, L, I, V ist. Fur die einzelnenBauelemente gelten dann die Beziehungen

iR = GvR, iC = Cd vC

dt, vL = L

d iLdt

, (6)

wenn G=diag(G1, . . . , GkG), C =diag(C1, . . . , CkC

) und L=diag(L1, . . . , LkL)

die Leitwert-, Kapazitats- und Induktivitatsmatrix ist. Mit kE bezeichnen wirdie Anzahl der Elemente E im Netzwerk.

Wenn fur die Inzidenzmatrix A = (ARACALAV AI) geschrieben wird, ist dasGleichungssystem

ARiR + ACiC + ALiL + AV iV + AIiI = 0, (KCL)

vR = ATRe, vC = AT

Ce, vL = ATLe, vV = AT

V e, vI = ATI e, (KVL)

iR = GvR, iC = Cd vC

dt, vL = L

d iLdt

(Bauelemente)

vV = vs, iI = is (unabhangige Quellen)

zu losen. is und vs bezeichnen unabhangige Strom- und Spannungsquellen.

11

Dieses System wird auch als Sparse-Tableau bezeichnet. Hier zahlen alle Zweig-strome und -spannungen sowie die Knotenpotentiale zu den Unbekannten.Schon bei kleinen Schaltungen fuhrt der Sparse-Tableau-Ansatz auf Systememit vielen Unbekannten.

Im Rahmen der modifizierten Knotenanalyse wird versucht, die Anzahl derUnbekannten zu reduzieren. Wegen KVL konnen alle Zweigspannungen durchKnotenpotentiale ausgedruckt werden. Daruber hinaus konnen Bauelementbe-ziehungen in Admittanzform direkt in KCL eingesetzt werden. Damit genugt

es, als Vektor der Unbekannten x =(

eiViL

)zu betrachten und die DAE

ACCATC 0 0

0 L 00 0 0

·

e(t)iL(t)iV (t)

′

+

ARGATR AL AV

−ATL 0 0

ATV 0 0

·

e(t)iL(t)iV (t)

=

−AI is(t)

0vs(t)

(7)

zu losen. Zur Bestimmung der Kerne und Bilder der beteiligten Matrizen istdas folgende Lemma hilfreich [28].

Lemma 10 Sei M ∈ MIR(k, k) positiv definit. Fur N ∈ MIR(l, k) gilt dann

(i) ker NMNT = ker NT (ii) im NMNT = im N.

Der Beweis des Lemmas ist im Anhang B angegeben.

Wie bei der Untersuchung der Abstrakten Differential-Algebraischen Systemeim Kapitel 2 werden Projektoren als zentrales Hilfsmittel zur Untersuchungvon (7) verwendet. Ist Q ∈ Lb(IR

k, IRk), k ∈ IN, ein Projektor, so wird mitP = Ik − Q wiederum der zu Q komplementare Projektor bezeichnet. Furjedes v ∈ IRk ist dann v = Pv + Qv, so dass im P = ker Q und ker P = im Qfolgt. Es ist PQ = QP = 0 und im P ⊕ ker Q = IRk.

Wir setzen voraus, dass die Kapazitatsmatrix C positiv definit ist. Sei dann QC

ein Projektor auf ker ATC = ker ACCAT

C (Lemma 10) und sei PC = In−1 −QC .Wir wahlen

A =

ACCATC 0

0 L0 0

, D =

(PC 0 00 IkL 0

), B =

ARGATR AL AV

−ATL 0 0

ATV 0 0

und q(t) =

−AI is(t)

0vs(t)

. Dann ist (7) eine lineare DAE der Form

A(Dx(t)

)′+ Bx(t) = q(t). (8)

Die Koeffizienten A ∈ C(I, Lb(IR

n−1+kL , IRm)), D ∈ C

(I, Lb(IR

m, IRn−1+kL))

und B ∈ C(I, Lb(IR

m, IRm))

sind konstante Matrixfunktionen. m=n−1+kL+kV

ist die Dimension des Systems (8).

12

Im Weiteren wird stets die elementbezogene Inzidenzmatrix (ARACALAV AI)verwendet. Die Bezeichnung der Koeffizientenmatrix mit dem Symbol ’A’ kanndaher nicht zu Verwechslungen mit der Inzidenzmatrix fuhren.

Wir suchen Losungen

x ∈ C1D

(I, IRm

)=

C

(I, IRm

) ∣∣ Dx ∈ C1(I, IRm

) ,

die fur jedes t ∈ I der Gleichung (8) genugen.

Statt (6) konnen die Bauelementbeziehungen auch ladungs- bzw. flussorientiertformuliert werden. In diesem Fall ergeben sich DAEs der Form

Ay′ + f(x) = s(t), y = q(x),

wobei y ein Vektor von Ladungen und Flussen ist [5].

3.3 Der Index der Netzwerkgleichung

Wir betrachten die DAE (8) auf dem Zeitintervall I = [t0, T ]. Besondere Auf-merksamkeit ist dem Hauptterm (A,D) zu widmen [19].

Definition 11 Seien A ∈ C(I, Lb(IR

n, IRm)), D ∈ C

(I, Lb(IR

m, IRn))

ste-tige Matrixfunktionen. Das geordnete Paar (A,D) heißt gut zusammenpassend,wenn

ker A(t)⊕ im D(t) = IRn, t ∈ I, (9)

und wenn diese Unterraume durch stetig differenzierbare Basisfunktionen auf-gespannt werden.Der Hauptterm von (8) heißt proper formuliert, wenn A und D in diesemSinne gut zusammenpassen.

Aufgrund der Bedingung (9) ist auch

ker AD = ker D und im AD = im A,

so dass sich Definition 11 fur die endlichdimensionenalen Raume X = Y = Rm

und Z = Rn als Spezialfall der Definition 2 ergibt.

Die Wahl der (konstanten) Matrixfunktionen A und D in Abschnitt 3.2 sichertdie propere Formulierung des Haupttermes von (8). Es ist

ker A = ker ACCATC × 0 = ker AT

C × 0 = im QC × 0

13

und

im D = im PC × IRkL .

Die Struktur von A, D und B bestimmt die Eigenschaften der Gleichung.Diese Matrixfunktionen werden daher zur Definition des Index µ der DAEherangezogen [19].

Sei dazu G0 = AD. Wegen im QC = ker ATC = ker ACCAT

C ist

Q0 =

(QC 0 00 0 00 0 IkV

)

ein Projektor auf ker G0. Es sei G1 = G0 + B Q0. Hier hangt G1 nicht von derZeit t ab. In diesem Fall wird G2 = G1 +B P0 Q1 definiert, wobei P0 = Im−Q0

ist. Außerdem ist Q1 ein Projektor auf ker G1.

Definition 12 Der Hauptterm von (8) sei proper formuliert. Dann besitzt(8) den (Behandelbarkeits-)Index µ ∈ 0, 1, 2, wenn Gi fur i = 0, . . . , µ − 1singular, aber Gµ regular ist.

Bemerkung 13 Die Definition 12 ist fur X = Y = Rm und Z = Rn ein Spe-zialfall der Definition 4 der Abstrakten Differential-Algebraischen Systeme.Die allgemeineren ADA-Systeme konnen DAEs als Komponenten enthalten,so dass in diesem Fall zum einen der Index des Gesamtsystems (Definition 4)und zum anderen der Index der involvierten DAE untersucht werden kann(Definition 12). In der vorliegenden Arbeit wird fur ADA-Systeme aus derSchaltungssimulation die Bedeutung des DAE-Index fur den Index des Ge-samtsystems, den ADAS-Index, untersucht.

Der Index µ der DAE (8) ist unabhangig von der speziellen Wahl von Q0 undQ1, und gibt, grob gesprochen, an, wie stark die DAE von einer gewohnlichenDifferentialgleichung abweicht. Ist µ = 0, so liegt eine gewohnliche Differen-tialgleichung vor. Ist aber µ = 1, so muss die Losung x algebraischen Zusatzbe-dingungen genugen. In diesem Fall konnen auch die Anfangswerte x(t0) = x0

nicht beliebig vorgegeben werden, sondern mussen Nebenbedingungen genugen(konsistente Anfangswerte).

Durch die Verwendung weiterer Projektoren kann die Matrixkette Gi fortge-setzt werden und so auch ein Index µ ≥ 3 definiert werden [25]. Im Folgendenwerden jedoch ausschließlich Index-1 Netzwerke betrachtet, so dass auf dieDefinition hoherer Indizes verzichtet werden kann.

Neben dem hier vorgestellten Behandelbarkeitsindex (tractability-index) [19]gibt es noch zahlreiche weitere Indexbegriffe fur DAEs wie den Differentia-tionsindex [2] oder den Storungsindex [9]. Ist (8) eine DAE mit konstanten

14

Koeffizientenmatrizen, so sind diese Indexbegriffe samtlich aquivalent. In die-sem Fall ergibt sich stets der Kronecker-Index des Matrixbuschels λG0 + B.

Fur DAEs (8), die sich aus der MNA ergeben, kann der Index µ aus der Struk-tur des Netzwerkgraphen bestimmt werden. In [22, 29] wird der folgende Satzbewiesen:

Satz 14 Fur ein RCL-Netzwerk seien die Matrizen G, C und L positiv defi-nit. Das Netzwerk enthalte keine Schnittmengen von Stromquellen und keineSchleifen von Spannungsquellen. Ist x0 ∈ IRm ein konsistenter Anfangswert,so ist die DAE (8) mit der Bedingung x(t0) = x0 eindeutig losbar und es gilt:

• Die DAE (8) besitzt genau dann den Index 0, wenn

(0.a) das Netzwerk keine Spannungsquelle enthalt und

(0.b) im Netzwerk ein Baum aus kapazitiven Zweigen existiert.

• Die DAE (8) besitzt genau dann den Index 1, wenn

(1.a) das Netzwerk eine Spannungsquelle enthalt odereine Schnittmenge nicht-kapazitiver Zweige existiert und

(1.b) das Netzwerk keine LI-Schnittmenge enthalt und

(1.c) das Netzwerk keine CV -Schleife mit mindestenseiner Spannungsquelle enthalt.

• In allen anderen Fallen besitzt die DAE (8) den Index 2. ¤

Eine LI-Schnittmenge besteht ausschließlich aus Zweigen, die zu Kapazitatenund Stromquellen gehoren. Analog wird eine CV -Schleife von kapazitiven Zwei-gen und Zweigen, die Spannungsquellen entsprechen, gebildet.

3.4 Zusammenfassung

Im Folgenden betrachten wir ein verbundenes RCL-Netzwerk mit n Knotenund der Inzidenzmatrix (ARACALAV AI). Die MNA fuhrt auf die Algebro-Differentialgleichung (8)

ACCATC 0

0 L0 0

((PC 0 00 IkL 0

)x(t)

)′

+

ARGATR AL AV

−ATL 0 0

ATV 0 0

x(t) =

−AI is(t)

0vs(t)

mit proper formuliertem Hauptterm. Dabei ist x(t) =

(e(t)iV (t)iL(t)

), t ∈ I, der

Vektor der Unbekannten. QC ist ein Projektor auf ker ATC und PC =In−1−QC .

15

Wir setzen voraus:

(V1) Die Leitwert-, Kapazitats- und Induktivitatsmatrizen G, C und L seienpositiv definit.

(V2) Das Netzwerk enthalte keine Schnittmengen von Stromquellen undkeine Schleifen von Spannungsquellen und habe den Index µ = 1.

(V3) Der Hauptterm der Netzwerk-DAE sei proper formuliert, d.h.ker A⊕ im D = IRn.

Mit dem Projektor

Q0 =

(QC 0 00 0 00 0 IkV

)

auf ker AD ist

G1 = AD + B Q0 =

(ACCAT

C+ARGATRQC 0 AV

−ATLQC L 0

ATV QC 0 0

).

Bemerkung 15 Die Voraussetzungen (V1)-(V3) haben vier, fur die weiterenUberlegungen besonders wichtige Konsequenzen:

(i) G1 ist regular,

(ii) ker AV = 0,(iii) ker QT

CAV = 0 und

(iv) ker(ACARAV )T = 0.

Beweis: (i) ist gerade das Index-1-Kriterium und (ii) folgt mit Lem-ma 9 aus (V2). (iii) druckt die Tatsache aus, dass das Netzwerk wegen (V2)und Satz 14 keine CV -Schleife mit mindestens einer Spannungsquelle enthalt[28]. Der Satz 14 sichert wegen (V2) ebenfalls, dass keine LI-Schnittmengeexistiert. Damit gibt es einen Baum, der keine induktiven Zweige oder Zwei-ge von Stromquellen enthalt. Dies zeigt rang(ACARAV ) = n − 1, so dassdim ker(ACARAV )T = n − 1 − rang(ACARAV ) = 0 und somit (iv) folgt (vgl.Satz B.1 und B.3 im Anhang). ¤

16

4 Halbleitersimulation

Neben den Elementen des RCL-Netzwerkes sollen zusatzlich Halbleiterele-mente betrachtet werden. Diese Bauelemente werden in das vorhandene Netz-werk eingefugt ohne die Anzahl n der Netzwerkknoten zu verandern. Wir be-trachten ausschließlich Elemente wie z.B. p-n-Dioden mit jeweils zwei Kontak-ten. Die Lage der Halbleiterelemente im Netzwerk wird durch die Inzidenzma-trix AH ∈ MIR(n−1, kH) beschrieben, wobei kH die Anzahl der Halbleiterbau-elemente ist.

Als mathematisches Modell fur die Halbleitersimulation betrachten wir im Ab-schnitt 4.1 zunachst das Drift-Diffusionsmodell [26]. Fur die Untersuchung derKopplungstruktur beschranken wir uns jedoch auf eine einfache, eindimensio-nale partielle Differentialgleichung, die im Abschnitt 4.2 vorgestellt wird.

4.1 Das Drift-Diffusionsmodell

Um das Verhalten von Halbleiterbauelementen numerisch untersuchen zu kon-nen, ist ein mathematisches Modell dieser Bauelemente notig. Die Poissonglei-chung

ε∆ψ = q · (n− p− C) (10)

beschreibt fur t ∈ I ⊂ IR das elektrostatische Potential ψ(y, t) in Abhangigkeitvon der Elektronendichte n(y, t), der Locherdichte p(y, t) und dem DotierprofilC(y) im Halbleitergebiet Ω ⊂ IRr, r ∈ 1, 2, 3. ε ist die Dielektrizitatskon-stante des Halbleitermaterials und mit q wird die Elementarladung bezeichnet.Uber die Kontinuitatsgleichungen fur die Elektronenstromdichte Jn(y, t) unddie Locherstromdichte Jp(y, t),

−q · ∂n

∂t+ div Jn = q ·R,

q · ∂p

∂t+ div Jp = −q ·R,

(11)

hangen n und p ihrerseits von ψ ab, denn nach dem Drift-Diffusions-Modellfur den Ladungstragertransport ist

Jn = q µn (UT∇n− n∇ψ),

Jp = −q µp (UT∇p + p∇ψ).(12)

Dabei beschreiben µn und µp die Beweglichkeit der Ladungstrager und UT =k Tq

ist die thermische Spannung bei Raumtemperatur mit der Boltzmann-Konstanten k. Die Temperatur T wird als konstant vorausgesetzt, so dass

17

auch UT eine Konstante ist. In (11) beschreibt R die Generation und Rekom-bination von Ladungstragern. Unterschiedliche Modelle wie zum Beispiel dieShockley-Read-Hall-Rekombinationen

R = R(n, p) =n p− n2

i

τp(n + ni) + τn(p + ni),

wurden untersucht [17, 26]. ni ist die intrinsische Ladungstragerdichte und τn

bzw. τp die Lebensdauer der Elektronen und Locher.

In Tabelle 1 sind die numerischen Werte der verwendeten Großen zusammen-gestellt [21].

Bezeichnung Bedeutung numerischer Wert

ε0 Dielektrizitatskonstante 8.8542·10−14[

AsV cm

]im Vakuum

ε Dielektrizitatskonstante 11.7·ε0

in Siliziumq Elementarladung 1.6022·10−19 [As]T Temperatur 300 [K]

k Boltzmannkonstante 1.38·10−23[

JK

]UT Thermische Spannung bei Raum-

temperatur UT = k Tq

0.0259 [V ]

ni intrinsische Ladungstragerdichte 1·1010[

1cm3

]

µn Beweglichkeit der Elektronen 500[

cm2

V s

]

µp Beweglichkeit der Locher 500[

cm2

V s

]τn Lebensdauer der Elektronen 1·10−7 [s]τp Lebensdauer der Locher 1·10−7 [s]

Tabelle 1: Numerische Werte der Modellkonstanten

Die Gleichungen (10)-(12) bilden ein nichtlineares System gekoppelter partiel-ler Differentialgleichungen. Dabei sind (11) parabolische PDEs, (10) ist jedochvon elliptischem Charakter. In Analogie zu den DAEs werden bereits solcheSysteme von einigen Autoren als Partielle Differential-Algebraische Gleichun-gen (PDAEs) bezeichnet [10]. Die elliptische Gleichung wird dabei als ent-artete parabolische Gleichung betrachtet. Eine Ortsdiskretisierung uberfuhrtdas Gleichungssystem in eine Algebro-Differentialgleichung [15].

Um die eindeutige Losbarkeit des Systems zu sichern, mussen geeignete Rand-bedingungen gefordert werden. Als Beispiel betrachten wir die zweidimensio-nale Diode aus Abbildung 2. An Kontakten zu Isolatoren (B-C, C-D, E-A)

18

gelten homogene Neumann-Randbedingungen

∂ψ

∂n= 0,

∂n

∂n= 0,

∂p

∂n= 0, (13)

wobei n die außere Einheitsnormale an das Simulationsgebiet Ω ist.

A B C

DE

Anode

Kathode

Abbildung 2: Skizze einer zweidimensionalen Diode

An den Ohmschen Kontakten (A-B, D-E) sind Dirichlet-Randbedingungen

ψ = ψbi + ψE, (14)

n =1

2

(C +

√C2 + 4n2

i ), p =1

2

(− C +√

C2 + 4n2i ) (15)

zu fordern. Die Gleichungen (15) ergeben sich aus der Uberlegung, dass anOhmschen Kontakten die Raumladung verschwindet und dass sich die La-dungstrager dort im thermodynamischen Gleichgewicht befinden [26]:

p− n + C = 0, n p = n2i .

ψE ist das von außen angelegte externe Potential und ψbi bezeichnet das ’built-in’-Potential des Halbleiters. Nach [26] ergibt sich ψbi aus dem Dotierprofilnaherungsweise zu

ψbi(x) = UT ln(C +

√C2 + 4n2

i

2 ni

). (16)

Das System (10)-(12) wird als Drift-Diffusionsmodell bezeichnet. Van Roos-broeck hat 1950 in [30] die grundlegende Form dieser Gleichungen angegeben.In [26] zeigt Selberherr, dass die Gleichungen (10), (11) aus den MaxwellschenGesetzen abgeleitet werden konnen. Zur Ableitung von (12) sind dagegen eineReihe zusatzlicher physikalischer Annahmen notig (vgl. [26], Seite 17).

19

Die angelegte externe Spannung ψE definiert die Randbedingungen fur ψ. DieAuswertung der Stromdichte J = Jn + Jp + JD liefert den Strom durch dasBauelement. Dabei ist

JD = −ε grad ψt

die Driftkomponente des Stromes. Diese Driftkomponente simuliert das kapa-zitive Verhalten des Halbleiters.

Im gesamten Halbleitergebiet Ω gilt zu jedem Zeitpunkt t

div J = 0, (17)

denn wenn (10) nach der Zeit t differenziert wird, folgt aus (11)

div J = div(Jn + Jp + JD) = q∂n

∂t+ q R− q

∂p

∂t− q R − q

∂n

∂t+ q

∂p

∂t= 0.

Hierbei wird vorausgesetzt, dass das Dotierprofil C zeitlich konstant ist. DerStrom iH durch die Diode aus Abbildung 2 ergibt sich also zu

iH(t) =

∫

O1

J(y, t) · n(y) dy = −∫

O2

J(y, t) · n(y) dy,

wenn n der Vektor der außeren Einheitsnormalen an Ω ist. O1 und O2 bezeich-nen die Ohmschen Kontakte A-B und D-E.

Falls das Simulationsgebiet Ω = (0, l) ⊂ IR1, l > 0, eindimensional ist, erhaltenwir fur den Kopplungsstrom

iH(t) = J(0, t) = −J(l, t)

= σ1

( −∇n(0,t)∇p(0,t)

)UT + σ1

( n(0,t)p(0,t)

)∂ψ

∂y(0, t)− ε

∂ψt

∂y(0, t).

(18)

mit der linearen Abbildung σ1(ab ) = −q (µn a + µpb).

4.2 Das PDE-Modell

Im Kapitel 5 soll das RCL-Netzwerk aus Kapitel 3 um Halbleiterbauelementeerweitert werden. Das so entstehende gekoppelte System aus DAEs und PDEssoll als Abstraktes Differential-Algebraisches System untersucht werden. DieADA-Systeme wurden vorerst lediglich fur den linearen Fall formuliert. Da-her sind wir an einem linearen Modell partieller Differentialgleichungen furden Halbleiter interessiert, um das ADAS-Konzept anwenden zu konnen (Ab-schnitt 2).

20

Im Rahmen der vorliegenden Arbeit kann folglich nicht das vollstandige, nicht-lineare Drift-Diffusionsmodell aus Abschnitt 4.1 untersucht werden. Es musszu einem linearen Modell ubergegangen werden, das dann notwendigerweisenur einige, fur die Kopplung wesentliche Elemente des Drift-Diffusionsmodellsaufgreifen kann.

Anstelle von (10)-(15) wird ausschließlich die partielle Differentialgleichung

ε ∆ψ(y, t) = f(y, t), y ∈ Ω, t ∈ I, (19)

mit einer rechten Seite f ∈ C1(I, L2(Ω)

)betrachtet. Dabei ist Ω = (0, l) ⊂ IR.

An den zwei Kontakten y = 0 und y = l werden Dirichlet-Randbedingungen

ψ(y, t) = ψbi(y) + ψE(y, t) y ∈ ∂Ω, t ∈ I, (20)

gefordert.

Im Wesentlichen wird also auf die Simulation der Ladungstrager verzichtet, sodass im Weiteren auch nicht von einem

”Kopplungsstrom“ gesprochen werden

kann. iH wird daher als”Kopplungsgroße“ bezeichnet.

Werden n und p als konstant angesehen, entfallt in (18) der von ∇n und ∇p

abhangige Term. Außerdem ist dann σ1 = σ1

( n(0,t)p(0,t)

) ∈ IR konstant.

Es wird versucht, bei der Definition von iH die verbleibenden vom Potential ψund von ψt abhangigen Komponenten mitzusimulieren. Sei

r : H2(Ω) → IR, ϕ 7→ rϕ =∂ϕ

∂y(0) (21)

der Restriktionsoperator, der den Gradienten von ϕ ∈ H2(Ω) am linken Inter-vallende auswertet. Als Kopplungsgroße betrachten wir

iH(t) = σ1 r ψ(·, t)− σ2 rψt(·, t) (22)

mit Konstanten σ1, σ2 ∈ IR. Es sei σ2 ≥ 0.

Damit rϕ definiert ist, muss die (schwache) Ableitung ∂ϕ∂y

Randwerte besitzen.

Es ist daher ∂ϕ∂y∈ H1(Ω) zu fordern. In (21) wird folglich H2(Ω) als Defini-

tionsbereich von r gewahlt.

Fur spatere Berechnungen wird im folgenden Lemma das Bild im r bestimmt.

Lemma 16 Sei l > 0 und Ω = (0, l) ⊂ IR. Dann ist r ein linearer, be-schrankter Operator mit im r = IR.

21

Beweis:

(i) Die Linearitat von r ist klar. Die Beschranktheit folgt aus der stetigenEinbettung H2(Ω) → C1(Ω) bei Ω ⊂ IR1 (Satz A.5 im Anhang), dennfur alle ϕ ∈ H2(Ω) ist

|rϕ| = |ϕ′(0)| ≤ maxy∈Ω

|ϕ(y)|+ maxy∈Ω

|ϕ′(y)|

= ‖ϕ‖C1(Ω) ≤ c ‖ϕ‖H2(Ω)

mit einer von ϕ unabhangigen Konstante c.

(ii) Es ist im r ⊂ IR. Sei a ∈ IR beliebig. Fur y ∈ Ω setzen wir ϕa(y) =a y

l(l− y). Dann ist ϕa ∈ H2(Ω) ∩H1

0 (Ω) mit r ϕa = ∂ϕa

∂y(0) = a, was die

Surjektivitat von r zeigt. ¤

Aus dem Beweis folgt, dass auch die Einschrankung r : H2(Ω) ∩H10 (Ω) → IR

surjektiv ist. Im Weiteren wird nicht zwischen r und r unterschieden. Ausden verwendeten Funktionenraumen wird stets deutlich werden, ob r oder dieEinschrankung r zu betrachten ist.

Wegen f ∈ C1(I, L2(Ω)

)ist die Abbildung

F : I → L2(Ω), t 7→ F (·, t) mit F (y, t) =1

ε

y∫

0

ξ∫

0

f(τ, t) dτ dξ

ein Element von C1(I, H2(Ω)

)(vgl. die Beispiele A.4 und A.8 im Anhang A).

Die zweimalige Integration von (19) nach dem Ort liefert die Darstellung

ψ(y, t) = F (y, t) + y r ψ(·, t) + ψbi(0) + ψE(0, t),

wobei wegen (20) ψ(l, t) = ψbi(l) + ψE(l, t) gelten muss. Folglich ist

rψ(·, t) = −1

l

(F (l, t) + ψbi(0) + ψE(0, t)− ψbi(l)− ψE(l, t)

).

Falls das externe Potential ψE(0, ·), ψE(l, ·) nach t differenzierbar ist, liegt ψin C1

(I, H2(Ω)

), so dass die Ableitung ψt definiert ist.

Die Abbildung ψbi stellt im Drift-Diffusionsmodell das vom Dotierprofil Cabhangige ’built-in’-Potential (16) dar. Hier ist ψbi zunachst eine beliebigeFunktion aus H2(Ω).

Die bisher getroffenen Vereinbarungen werden als Voraussetzung (V4’) zusam-mengefasst.

(V4’) Es sei Ω = (0, l), l > 0, f ∈ C1(I, L2(Ω)

)und ψbi ∈ H2(Ω).

Fur die Kopplungsgroße iH gelte σ2 ≥ 0.

22

4.3 Bemerkungen zum PDE-Modell

Das System (19), (20) ist das typisches Beispiel einer linearen, elliptischenpartiellen Differentialgleichung zweiter Ordnung. In [32] wird das System zumBeispiel als stationares Erhaltungsgesetz studiert. Dort ist ψ(y) die Tempera-tur eines Korpers Ω im Punkt y ∈ Ω.

Fur die reale Halbleitersimulation besitzt dieses Modell jedoch keine physi-kalische Relevanz. Der Stromfluss durch das Halbleiterelement kann nur beivorhandenen Ladungstragern simuliert werden. Ohne die Modellierung vonLadungstragern kann auch nicht von einem

”Kopplungsstrom“ gesprochen wer-

den.

Obwohl die Kontinuitatsgleichungen (11) vernachlassigt werden, ist das Ge-samtsystem aus Netzwerk-DAE und Halbleiter-PDE dennoch nicht entkop-pelt, da die angelegte externe Spannung ψE wie beim vollstandigen Drift-Diffusionsmodell die Randbedingungen fur ψ definiert. Das externe PotentialψE variiert mit der Zeit t, so dass fur unterschiedliche Zeitpunkte im Allge-meinen auch unterschiedliche Probleme (19), (20) zu losen sind.

Die vom Zeitpunkt t abhangige Losung ψ(·, t) bestimmt nun ihrerseits einenWert iH , der als Kopplungsgroße in die Netzwerksimulation zuruckfließt. Beider Definition (22) von iH wird die Abhangigkeit des Kopplungsstromes von∇ψ und ∇ψt beachtet.

Im vorliegenden PDE-Modell gibt es keine Entsprechung fur die Eigenschaft(17) des Drift-Diffusionsmodells. ∇ψ und ∇ψt werden am linken Rand des Si-mulationsgebietes ausgewertet. Es ist nicht gesichert, dass die Auswertung amrechten Rand zu den gleichen (absoluten) Werten fuhrt. Dennoch wird spateriH im Rahmen der Kirchhoffschen Knotengleichungen als Kopplungsgroße furbeide Kontakte des Halbleiters verwendet.

Trotz dieser starken Einschrankungen wird das PDE-Modell aus den folgendenGrunden dennoch untersucht:

(i) Das lineare Modell ermoglicht die Untersuchung des gekoppelten Systemsals Abstraktes Differential-Algebraisches System.

(ii) Die mit diesem Modell gewonnenen Aussagen konnen als Ausgangs-punkt fur weitere Untersuchungen mit umfangreicheren Modellen wiedem vollstandigen Drift-Diffusionsmodell dienen und

(iii) die gewonnenen Erkenntnisse konnen die Erweiterung des ADAS-Kon-zeptes auf den nichtlinearen Fall weiter motivieren.

Es ist zu erwarten, dass sich trotz der Unzulanglichkeiten des Modells qualita-tive Aussagen zur Kopplungsstruktur und zum ADAS-Index eventuell auch auf

23

den nichtlinearen Fall ubertragen lassen. Damit ist klar, dass in der vorliegen-den Arbeit die qualitative Untersuchung der Kopplung im Vordergrund stehenmuss. Eine analytische oder auch nur numerische Losung der gekoppelten Sys-teme ist nicht das Ziel der Betrachtungen. Selbst wenn uber Integration o.a.eine explizite Darstellung fur ψ moglich ist, wird hiervon bei der Formulierungdes Gesamtsystems kein Gebrauch gemacht, da dann ein entkoppeltes Systemfolgen wurde. Gerade die Kopplungsstruktur und die Auswirkungen auf denADAS-Index sollen aber im Vordergrund der weiteren Betrachtungen stehen(vgl. Beispiel 5 im Kapitel 2).

Das Modell (19), (20) ist zwar durch die Halbleiterphysik motiviert, es simu-liert jedoch kein reales Bauelement. Die hypothetischen, durch (19), (20) be-schriebenen Elemente werden daher als

”Halbleiter“ (in Anfuhrungsstrichen)

bezeichnet.

24

5 PDAEs als Abstrakte

Differential-Algebraische Systeme

Nachdem in den Kapiteln 3 und 4 die Netzwerk-DAE und ein PDE-Modellvorgestellt wurden, werden nun beide Modelle als gekoppeltes System betrach-tet. Solche gekoppelten Systeme partieller und differential-algebraischer Glei-chungen werden auch als Partial Differential-Algebraic Equations (PDAEs)bezeichnet.

Zunachst wird das PDAE-System in den Abschnitten 5.1 und 5.2 als Abstrak-tes Differential-Algebraisches System formuliert, dessen Index im Abschnitt5.3 untersucht wird. Dabei betrachten wir zunachst die Kopplungsgroße

iH(t) = σ1 r ψ(·, t),

bei der die”Driftkomponente“ rψt(·, t) vernachlassigt wird. Die Aussagen zum

ADAS-Index bei der Kopplung mittels

iH(t) = σ1 r ψ(·, t)− σ2 rψt(·, t)

nehmen jedoch den großeren Raum ein. Die verschiedenen Kopplungsgroßenentsprechen unterschiedlichen Detaillierungsgraden bei der Modellbildung. Imersten Fall bleiben kapazitive Effekte der

”Halbleiterelemente“ unberucksich-

tigt.

Die wesentlichen Ergebnisse werden in den Satzen 20 und 27 zusammengefasst.Dem Beweis von Satz 27 ist der Abschnitt 5.4 gewidmet. Im Anschluss daranwerden einige Beispiele betrachtet, die die wesentlichen Aussagen aber auchdie Beweise illustrieren.

5.1 Das gekoppelte System

Wir betrachten ein verbundenes RCL-Netzwerk mit n Knoten und der In-zidenzmatrix (ARACALAV AI). QC sei ein Projektor auf ker AT

C . Mit PC =In−1 −QC wird der zu QC komplementare Projektor bezeichnet.

Die Modifizierte Knotenanalyse (MNA) fuhrt auf die Algebro-Differentialglei-chung

A(Dx(t)

)′+ Bx(t) = q(t), t ∈ I = [t0, T ], (8)

aus Kapitel 3 mit proper formuliertem Hauptterm und x(t) =

(e(t)iV (t)iL(t)

)∈ IRm,

m = n−1+kL +kV . Wieder ist n die Anzahl der Knoten im Netzwerk und kL,

25

kV gibt die Anzahl der Induktivitaten bzw. Spannungsquellen an. Wir setzen(V1)-(V3) voraus (vgl. Abschnitt 3.4). Insbesondere hat (8) den Index 1.

Dieses Netzwerk wird durch”Halbleiterelemente“ mit jeweils zwei Kontakten

erganzt, die durch das Modell aus Abschnitt 4.2 beschrieben werden. Wir be-trachten kH solcher Elemente, die jeweils der Bedingung (V4’) genugen.

Fur k = 1, . . . , kH ist das System (19), (20)

ε ∆ψk(y, t) = fk(y, t), y ∈ Ω = (0, l), t ∈ I,

ψk(y, t) = ψbi,k(y) + ψE,k(y, t) y ∈ ∂Ω, t ∈ I,(23)

zu losen. Hier wurde angenommen, dass fur die Simulationsgebiete Ωk =(0, lk) = (0, l) gilt. Da (19), (20) eine lineare PDE ist, stellt dies keine Ein-schrankung der Allgemeinheit dar. Ωk kann durch eine lineare Transformationauf das Intervall (0, l) abgebildet werden. Mit den Bezeichnungen

ψ =

(ψ1

...ψkH

), ψE =

(ψE,1

...ψE,kH

), ψbi =

(ψbi

...ψbi

)und f =

(f1

...fkH

)

kann das obige System vektoriell als

ε ∆ψ(y, t) = f(y, t), y ∈ ΩkH = (0, l)kH , t ∈ I,

ψ(y, t) = ψbi(y) + ψE(y, t) y ∈ (∂Ω)kH , t ∈ I,(24)

geschrieben werden10. Der Laplace-Operator ∆ wird hier komponentenweiseangewandt.

Im Folgenden soll der Restriktionsoperator r ebenfalls komponentenweise aufVektoren operieren:

r : [H2(Ω)]kH → IRkH , ϕ 7→ rϕ =(∂ϕk

∂y(0)

)k=1,...,kH

. (25)

Mit r 1 : H2(Ω) → IR, r 1ϕ = ∂ϕ∂y

(0) wird die skalare Version bezeichnet.

Die Kopplungsgroße ist dann durch den Vektor

iH(t) = σ1 r ψ(·, t)− σ2 rψt(·, t), (26)

definiert, σ1, σ2 ∈ IR und σ2 ≥ 0.

10(24) ist lediglich eine verkurzte Schreibweise fur (23). Die Komponenten yi von y variie-ren fur i = 1, . . . , kH in Ω = (0, l). Wegen (23) sind Randbedingungen fur yi = 0 und yi = lvorzugeben. Randbedingungen werden folglich nur fur y ∈ (∂Ω)kH und nicht fur y ∈ ∂(ΩkH )formuliert.

26

Fur die weiteren Uberlegungen sind die tatsachlichen Zahlenwerte von σ1 undσ2 irrelevant. Bei der Formulierung von iH wird angenommen, dass fur jedes derkH ”

Halbleiterelemente“ die gleichen Werte σ1, σ2 ∈ IR gewahlt wurden. Dannkonnen σ1, σ2 weiterhin als Zahlen behandelt werden, was die Schreibweisevieler Formeln vereinfacht11.

In iH gehen ausschließlich die Ortsableitungen von ψ und ψt =(

∂ψk

∂t

)k=1,...,kH

ein, so dass iH nur von der Differenz ψE(0, t)−ψE(l, t) abhangt, die ihrerseitsaufgrund des Kirchhoffschen Spannungsgesetzes (KVL) durch die Knotenpo-tentiale e(t) gegeben ist:

ψE(0, t)− ψE(l, t) = AH e(t).

Wir konnen ohne Einschrankung annehmen, dass ψE(l, t) = 0 gilt. Dann istψE(0, t) = AT

H e(t) und das PDE-System nimmt die Form (27)-(29) an12.

ε ∆ψ(y, t) = f(y, t), y ∈ ΩkH = (0, l)kH , t ∈ I, (27)

ψ(0, t) = ψbi(0) + ATH e(t), t ∈ I, (28)

ψ(l, t) = ψbi(l), t ∈ I. (29)

Fur die Randbedingung (28) konnen wir mit E =(

AH00

)auch

ψ(0, t) = ψbi(0) + ATH e(t)) = ψbi(0) + ET x(t)

schreiben.

Die Kopplungsgroßen iH,k des Vektors iH mussen bei der Formulierung derKnotengleichungen beachtet werden, da sie als

”Kopplungsstrome“ einen Bei-

trag zum Kirchhoffschen Stromgesetz (KCL) liefern. Statt (8) muss die DAE

A(Dx(t)

)′+ Bx(t) + EiH = q(t), t ∈ I = [t0, T ], (30)

betrachtet werden.

Wie im Kapitel 2 angedeutet, werden die Randbedingungen (28), (29) in ADA-Systemen durch eine Homogenisierungsprozedur und die Wahl passender Funk-tionenraume beachtet. Sei dazu

(V5’) h ∈ H2(Ω) mit h(0) = 1 und h(l) = 0.

Es kann beispielsweise h(y) = 2l3

y3 − 3l2

y2 + 1 gewahlt werden. Fur dieseFunktion h gilt r 1h = ∂h

∂y(0) = 0. Wir setzen

11Defacto lassen wir nur Vektoren (σi,1, . . . , σi,kH)T = σi1 zu. Dies ist keine Beschrankung

der Allgemeinheit, da unter Verwendung der komponentenweisen Multiplikation ⊗ (vgl.Seite 48) alle Ergebnisse unverandert auch fur Vektoren (σi,1, . . . , σi,kH

)T bewiesen werdenkonnen.

12Mit ψ(0, t) ist hier ψ(0 · 1, t) gemeint. Analog ist ψ(l, t) und ψbi(0) als ψ(l · 1, t) undψbi(0 · 1) zu interpretieren.

27

ψ0(y, t) = ψbi(y) + h(y) ET x(t).

Fur

ψ(y, t) = ψ(y, t)− ψ0(y, t)

ist dann

ε∆ψ(y, t)=ε∆ψ(y, t)−ε∆ψ0(y, t)=f(y, t)−ε∆ψbi(y)−ε(∆h(y)

)ET x(t),

y ∈ Ω, t ∈ I, und

ψ(0, t) = ψ(l, t) = 0, t ∈ I.

Fur die Kopplungsgroße berechnen wir

iH(t) = σ1 rψ(·, t)−σ2 rψt(·, t)= σ1 r ψ(·, t)−σ2 r ψt(·, t)+σ1 rψ0(·, t)− σ2 rψ0,t(·, t)= σ1 r ψ(·, t)−σ2 r ψt(·, t)+σ1 ·

(r ψbi+(r 1h) ET x(t)

)−σ2 ·(r 1h) ET x′(t).

Wenn

rψbi = 0 und r 1h = 0

gefordert wird, nimmt iH auch fur ψ die bekannte Struktur

iH(t) = σ1 r ψ − σ2 r ψt

an. Die Voraussetzungen (V4’) und (V5’) werden daher wie folgt erganzt13:

(V4) Es sei Ω = (0, l), l > 0, f ∈ C1(I, L2(Ω)

)und ψbi ∈ H2(Ω) mit

r ψbi = 0. Fur die Kopplungsgroße iH gelte σ2 ≥ 0.

(V5) h ∈ H2(Ω) mit h(0) = 1, h(l) = 0 und r 1h =∂h

∂y(0) = 0.

Wenn wir ψ wieder mit ψ bezeichnen, nimmt das System (26)-(29) die Gestalt

ε ∆ψ(y, t) = f(y, t)−ε∆ψbi(y)−ε(∆h(y)

)ET x(t), y ∈ ΩkH , t ∈ I, (31)

ψ(0, t) = ψ(l, t) = 0, t ∈ I, (32)

iH(t) = σ1 rψ(·, t)− σ2 r ψt(·, t), (33)

an.

13(i) Im Drift-Diffusionsmodell ist ψbi das ’built-in’-Potential (16). Fur reale Halbleiterele-mente ist haufig das Dotierprofil C in der Nahe der Kontakte in Normalenrichtung konstant,so dass dann r ψbi = 0 folgt.(ii) Fur die angegebene Funktion h(y) = 2

l3 y3 − 3l2 y2 + 1 ist r 1h = 0 erfullt.

28

Bemerkung 17 Die Voraussetzung (V5) sichert, dass ∆h 6= 0 ∈ L2(Ω) ist.Diese Eigenschaft wird bei der Indexuntersuchung eine Rolle spielen.

Beweis: Ware ∆h = 0, so wurde fur alle y ∈ Ω

h(y) = h(0) + y∂h

∂y(0) = 1

gelten – ein Widerspruch zu (V5), denn es ist h(l) = 0. ¤

5.2 Das Abstrakte Differential-Algebraische System

Mit den Gleichungen (30), (31)-(33) haben wir das gekoppelte System ausNetzwerk-DAE und

”Halbleiter“-PDE vollstandig formuliert. Dieses PDAE-

System wird als Abstraktes Differential-Algebraisches System untersucht.

Wir betrachten dazu die folgenden linearen Raume

X = IRm × [H2(Ω) ∩H10 (Ω)]kH × IRkH ,

Y = IRm × [L2(Ω)]kH × IRkH ,

Z = IRn+kL × IRkH .

Die Normen

|y|2 =

(k∑

i=1

y2i

) 12

, y ∈ IRk,

‖g‖2 =

∫

Ω

|g(y)|2 dy

12

, g ∈ L2(Ω),

‖g‖2,2 =

(2∑

i=0

∥∥ ∂i

∂yig(y)

∥∥2

2

) 12

, g ∈ H2(Ω),

induzieren Normen ‖ · ‖X , ‖ · ‖Y , ‖ · ‖Z auf den Produktraumen X, Y und Z.Im Anhang A wird gezeigt, dass X, Y und Z Hilbert-Raume sind.

Mit den Operatoren

A : Z → Y, ( vw ) 7→ A ( v

w ) =

A 00 00

√σ2 ·IkH

(vw

)

D : X → Z,(

vϕw

)7→ D

(vϕw

)=

(D 0 00 −√σ2 ·r 0

)

vϕw

B : X → Y,(

vϕw

)7→ B

(vϕw

)=

B 0 E

ε(∆h(·))ET ε∆ 0

0 σ1 r −IkH

vϕw

29

ist (30), (31)-(33) aquivalent zu dem Abstrakten Differential-AlgebraischenSystem

A(Du(t))′

+ Bu(t) = q(t), t ∈ I = [t0, T ], (34)

wobei

u(t) =

(x(t)

ψ(·,t)iH(t)

)und q(t) =

(q(t)

f(·,t)−ε∆ψbi(·)0

)

ist. In (34) werden sowohl der Laplace-Operator als auch der Restriktionsope-rator r komponentenweise angewandt.

Bemerkung 18

(i) Bei der Definition des Restriktionsoperators r in (21) musste H2(Ω) alsDefinitionsbereich gewahlt werden. Die Operatoren D und B wenden r

auf die zweite Komponente von u(t) an. Fur diese Komponente mussdaher ebenfalls H2(Ω) gewahlt werden.

(ii) Die zweite Komponente ψ(·, t) von u(t) muss zu allen Zeitpunkten denhomogenen Randbedingungen (32), ψ(·, t) ∈ [H1

0 (Ω)]kH , genugen. DieRandbedingungen werden also schon durch die Wahl des Funktionen-raumes X = IRm × [H2(Ω) ∩H1

0 (Ω)]kH × IRkH beachtet14.

(iii) Bei der Formulierung von (34) wurde ausgenutzt, dass die Reihenfolge derAbleitungen ∂

∂trψ(·, t)=r ψt(·, t) vertauscht werden kann (Anhang A).

Offenbar ist A : Z → Y ein linearer, beschrankter und normal losbarer Opera-tor. Wegen der Linearitat von r ist auch D : X → Z ein linearer Operator. Dasich die Aussagen des Lemmas 16 sinngemaß auf die vektorielle Version von r

aus (25) ubertragen lassen, folgt aus der Beschranktheit und der Surjektivitatvon r die Beschranktheit und normale Losbarkeit von D. Daruber hinaus fuhrtdie Wahl des Funktionenraumes X auch zu einem beschrankten Operator B,da der Laplace-Operator ∆ : H2(Ω) → L2(Ω) beschrankt ist. Es ist

A(t) ∈ Lb(Z, Y ), D(t) ∈ Lb(X, Z), und B(t) ∈ Lb(X,Y ). (35)

Wir bestimmen

kerA = ker A× ker(√

σ2 ·IkH

)und imD = im D × im

(√σ2 ·r

).

14Statt [H2(Ω)∩H10 (Ω)]kH darf nicht [H2

0 (Ω)]kH gewahlt werden, da in diesem Fall r ϕ = 0fur alle ϕ ∈ [H2

0 (Ω)]kH folgen wurde (vgl. Anhang A).

30

Da der Hauptterm der Netzwerk-DAE (30) aufgrund von (V3) proper formu-liert ist, gilt

Z = kerA ⊕ imD,

denn

• σ2 = 0 ⇒ ker(√

σ2 ·IkH

)= IRkH , im

(√σ2 ·r

)= 0

• σ2 > 0 ⇒ ker(√

σ2 ·IkH

)= 0, im

(√σ2 ·r

)= im r = IRkH .

Es ist

AD =

A 00 00

√σ2 ·IkH

(D 0 00 −√σ2 ·r 0

)=

AD 0 00 0 00 −σ2 r 0

.

Wegen ker AD = ker D und im AD = im A folgt außerdem

kerAD= ker AD × ker(σ2r )× IRkH = ker D × ker(σ2r )× IRkH = kerD,imAD = im AD × 0 × im(σ2r ) = im A× 0 × im(σ2r ) = imA,

so dass der Hauptterm von (34) entsprechend Definition 2 ebenfalls properformuliert ist.

5.3 Aussagen zum ADAS-Index

Zur Bestimmung des ADAS-Index des Systems (34) wird die Operatorkette Gi,i ≥ 0, aus Kapitel 2 betrachtet. Nach Definition 4 ist zu prufen, ab welchemIndex µ der Operator Gµ injektiv und dicht losbar ist.

Da A, D und B nicht explizit von der Zeit abhangen, kann stets ein konstanterProjektor Wi mit kerWi = imGi gewahlt werden. Damit ist die in Definition 4geforderte stetige t-Abhangigkeit von Wi(t) = Wi stets gegeben.

Zunachst wollen wir festhalten, dass ein ADAS-Index µ < 1 fur das System(34) nicht zu erwarten ist.

Satz 19 Fur den ADAS-Index µ von (34) ist stets µ ≥ 1.

Beweis: Es ist G0 = AD =(

AD 0 00 0 00 −σ2 r 0

), so dass G0 wegen kerG0 =

ker(AD)× ker(σ2r )× IRkH nicht injektiv ist. Folglich muss µ ≥ 1 gelten. ¤

Seien W0 und Rσ2 Projektoren langs im AD bzw. langs im(σ2r ). Dann wird

durch W0 =(W0 0 0

0 I 00 0 Rσ2

)ein Projektor langs imG0 definiert.

31

Da der Restriktionsoperator r nach Lemma 16 surjektiv ist, kann

Rσ2 =

IkH

, falls σ2 = 0

0 , falls σ2 > 0

gewahlt werden. In jedem Fall ist dim imW0 > 0.

5.3.1 ADAS ohne Driftkomponente

Fur einen Moment wollen wir die Kopplungsgroße

iH(t) = σ1 r ψ(·, t)betrachten, die sich fur σ2 = 0 ergibt. In diesem Fall wird also die

”Driftkom-

ponente“ rψt(·, t) vernachlassigt. Wegen

G0 = AD =(

AD 0 00 0 00 0 0

)

und kerG0 = ker(AD)× [H2(Ω) ∩H10 (Ω)]kH × IRkH = kerG0 ist

Q0 =

(Q0 0 00 I 00 0 IkH

)

ein Projektor auf N0 = kerG0, wobei Q0 der Projektor auf ker(AD) aus Ab-schnitt 3.4 ist.

Der folgende Satz zeigt, dass im Fall σ2 = 0 der ADAS-Index stets µ = 1betragt.

Satz 20 Es gelte (V1)-(V5). Die Matrix l ·IkH− σ1E

T Q0G−11 E sei regular.

Falls σ2 = 0 ist, gilt fur den ADAS-Index µ von (34) stets µ = 1.

Beweis: Es ist G1 = AD + BQ0 =

G1 0 Eε(∆h

)ET Q0 ε ∆ 0

0 σ1 r −IkH

, denn

G1 = AD + BQ0. Die Operatoren ∆ : H2(Ω) → L2(Ω) und r : H2(Ω) → IRsind beschrankt. Daher ist G1 ∈ Lb(X,Y ), was kerG0 = kerG0 bedeutet. Es ist

(vϕw

)∈ N1 = kerG1 ⇔

G1v + Ew = 0ε(∆h

)ET Q0v + ε ∆ϕ = 0

σ1 r ϕ− w = 0

⇔

v = −G−11 Eσ1 r ϕ

∆ϕ = σ1 ·(∆h

)ET Q0G

−11 E rϕ

w = σ1 rϕ

, (36)

32

da G1 wegen (V2) regular ist. Es wurde (V5) vorausgesetzt. Wegen Bemer-kung 17 folgt ∆h 6= 0 ∈ L2(Ω), so dass die rechte Seite der PDE ∆ϕ(y) =σ1 ·

(∆h(y)

)ET Q0G

−11 Er ϕ im Allgemeinen nicht verschwindet. Wird diese

Gleichung bezuglich des Ortes y integriert, ergibt sich die Darstellung

∂ϕ

∂y(y) = σ1 ·

(∂h

∂y(y)− ∂h

∂y(0)

)ET Q0G

−11 E r ϕ +

∂ϕ

∂y(0),

wobei wegen (V5) ∂h∂y

(0) = 0 ist. Eine weitere Integration liefert

ϕ(y) = σ1 ·(h(y)− h(0)

)ET Q0G

−11 E rϕ + y

∂ϕ

∂y(0) + ϕ(0). (37)

Wieder wegen (V5) ist h(0) = 1. Da(

vϕw

)∈ X, also ϕ ∈ [H2(Ω) ∩ H1

0 (Ω)]kH

ist, gilt ϕ(0) = ϕ(l) = 0. Die Eigenschaft ϕ(l) = 0 bringt insbesondere

0 = −σ1ET Q0G

−11 E rϕ + l

∂ϕ

∂y(0) =

(l·IkH

− σ1ET Q0G

−11 E

)rϕ,

denn h(l) = 0. Da l·IkH−σ1E

T Q0G−11 E nach Voraussetzung regular ist, muss

r ϕ = 0 gelten. Mit (36) und (37) folgt somit auch ϕ = 0 sowie v = 0 undw = 0, was die Injektivitat von G1 zeigt.

Der Raum Y0 = IRm× [C∞(Ω)]kH×IRkH liegt dicht in Y (vgl. Anhang A, [14]).Um die dichte Losbarkeit von G1 nachzuweisen, genugt es daher Y0 ⊂ imG1 zuzeigen.

Sei also(

z1ηz2

)∈ Y0 beliebig.

(vT , ϕT , wT

)T ∈ X ist ein Urbildelement, wenn

G1

(vϕw

)=

(z1ηz2

)⇔

v = G−11

(z1 − E

(σ1 r ϕ− z2

))

∆ϕ = 1εη − (

∆h)ET Q0G

−11

(z1 − E

(σ1 rϕ− z2

))

w = σ1 rϕ− z2

.

Zweimalige Integration liefert

ϕ(y) =1

ε

y∫

0

ξ∫

0

η(τ) dτ dξ+(1−h(y)

)ET Q0G

−11

(z1−E

(σ1 rϕ−z2

))+y r ϕ,

so dass wegen ϕ(l) = 0

(l·IkH

− σ1ET Q0G

−11 E

)rϕ = −1

ε

l∫

0

ξ∫

0

η(τ) dτ dξ −ET Q0G−11

(z1 + Ez2

)

33

gelten muss. Wieder wegen der Regularitat von l · IkH− σ1E

T Q0G−11 E ist

r ϕ eindeutig bestimmt. Die obigen Darstellungen fur ϕ, v und w liefern dasgesuchte Urbildelement. ¤

Eine genauere Untersuchung der Struktur der Matrix ET Q0G−11 E im Lem-

ma 34 des Abschnitts 5.4 wird zeigen, dass sich ET Q0G−11 E als

ET Q0G−11 E = AT

HQCQV−CH−1MR

QTV−CQT

CAH

darstellen lasst. QC und QV−C sind Projektoren auf ker ATC bzw. ker AT

V QC ,und HMR

ist eine positiv definite Matrix.

Wir konnen daher den Satz 20 auch wie folgt formulieren:

Folgerung 21 Es gelte (V1)-(V5). Die Matrix

l·IkH− σ1A

THQCQV−CH−1

MRQT

V−CQTCAH

sei regular. Falls σ2 = 0 ist, gilt fur den ADAS-Index von (34) stets µ = 1. ¤

Bemerkung 22 Es sei σ1 6= 0. Die Matrix

l·IkH− σ1E

T Q0G−11 E = l·IkH

− σ1ATHQCQV−CH−1

MRQT

V−CQTCAH

ist genau dann regular, wenn lσ1

kein Eigenwert der Matrix ET Q0G−11 E =

ATHQCQV−CH−1

MRQT

V−CQTCAH ist.

Zur Uberprufung der geforderten Regularitat sind also die Eigenwerte der Ma-trix AT

HQCQV−CH−1MR

QTV−CQT

CAH zu bestimmen.

Es soll besonders darauf hingewiesen werden, dass diese Matrix die Strukturder Matrizen aus Lemma 10 aufweist. Es ist

ATHQCQV−CH−1

MRQT

V−CQTCAH = MT

HH−1MR

MH

mit MH = QTV−CQT

CAH und der positiv definiten Matrix H−1MR

.

Die im Satz 20 geforderte Regularitat der Matrix

l·IkH− σ1E

T Q0G−11 E

kann stets erreicht werden. Im Falle der Singularitat genugt es, die Lange ldes Simulationsgebietes geringfugig zu andern, so dass l

σ1kein Eigenwert der

Matrix ET Q0G−11 E = MT

HH−1MR

MH ist.

34

Abgesehen von diesen mathematischen Uberlegungen zeigt der Satz aber, dassauch das Simulationsgebiet in Form der Lange l Einfluss auf den ADAS-Indexhaben konnte. Die Eigenwerte der Matrix ET Q0G

−11 E hangen jedoch lediglich

von den Daten des ursprunglichen Netzwerkes (Q0, G−11 ) und von der Lage der

”Halbleiter“ im Netzwerk (E, AH) ab.

5.3.2 ADAS mit Driftkomponente

Von nun an gelte σ2 > 0, so dass die Kopplungsgroße

iH(t) = σ1 r ψ(·, t)− σ2 rψt(·, t)

die”Driftkomponente“ mitsimuliert.

Es zeigt sich, dass in diesem Fall der ADAS-Index 1 oder 2 betragen kann. EinADAS-Index µ ≥ 3 kann fur das System (34) ausgeschlossen werden.

Wird die”Driftkomponente“ rψt(·, t) mitsimuliert, hangt der Index des ADA-

Systems ausschließlich von der Topologie des Netzwerkes ab. Sowohl die Lageder Elemente im RCL-Netzwerk (insbesondere Kapazitaten und Spannungs-quellen) als auch die Lage der

”Halbleiterelemente“ haben Einfluss auf den

Index.

In Definition 7 wurde bereits der Begriff der Schleife in einem Netzwerk ein-gefuhrt. Um im Satz 25 einen Zusammenhang zwischen der Netzwerktopologieund dem ADAS-Index hergestellen zu konnen, mussen spezielle Schleifen ein-zelner Bauelemente untersucht werden.

Ist E ein Element im Netzwerk, so wollen wir Zweige, die dem Element Eentsprechen, als E-Zweige bezeichnen.

Definition 23 Eine Schleife heißt

(i) H-Schleife, wenn sie ausschließlich aus H-Zweigen besteht.

(ii) CH-Schleife, wenn sie mindestens einen H-Zweig enthalt und ansonstenausschließlich aus H- und C-Zweigen gebildet wird.

(iii) CVH-Schleife, wenn sie mindestens einen H-Zweig enthalt und ansonstenausschließlich aus H-, V - und C-Zweigen gebildet wird.

Jede H-Schleife ist also eine spezielle CH-Schleife, und jede CH-Schleife kannals CVH-Schleife interpretiert werden. Die Umkehrung ist im Allgemeinen je-doch nicht richtig.

35

Es ist G0 =(

AD 0 00 0 00 −σ2 r 0

)mit σ2 > 0, so dass W0 =

(W0 0 00 I 00 0 0

)ein Projektor

langs imG0 ist, wenn wieder W0 ein Projektor langs im(AD) ist. Offenbar istdim imW0 > 0.

Im Beweis von Satz 19 wurde kerG0 bereits bestimmt. Wegen σ2 6= 0 ist nun

kerG0 = ker(AD)× ker r × IRkH = kerG0 = N0.

Sei Qr 1 ein beliebiger Projektor auf ker r 1 ⊂ H2(Ω) ∩H10 (Ω), so ist

Qr =

0BB@ Qr 1 0 · · · 00 Qr 1 · · · 0

. . . . . . . . . . . . . . . . . . . . . . .0 0 · · · Qr 1

1CCA (38)

ein Projektor auf ker r ⊂ [H2(Ω) ∩H10 (Ω)]kH . Ist außerdem Q0 der Projektor

auf N0 = ker(AD) aus Abschnitt 3.4, so ist

Q0 =

(Q0 0 00 Qr 00 0 IkH

)

ein Projektor auf N0. Wir berechnen

G1 = AD + BQ0

=

G1 0 E

ε(∆h

)ET Q0 ε ∆Qr 0

0 σ1 r Qr−σ2 r −IkH

=

G1 0 E

ε(∆h

)ET Q0 ε ∆Qr 0

0 −σ2 r −IkH

,

denn rQr = 0. Wie im vorigen Abschnitt ist G1 : X → Y ein linearer be-schrankter Operator. Der Index des ADA-Systems ist genau dann 1, wenndieser Operator injektiv und dicht losbar ist.

Das folgende Lemma stellt mogliche Projektoren Qr 1 bereit.

Lemma 24 Es sei Ω = (0, l) und r 1 : H2(Ω) ∩ H10 (Ω) → IR sei der in (25)

definierte Restriktionsoperator (genauer: dessen Einschrankung).Wir fixieren η ∈ H2(Ω) ∩H1

0 (Ω) mit r 1η = 1 (z.B. ist η = yl(l − y) moglich).

Dann ist

Qr 1,η : H2(Ω) ∩H10 (Ω) → H2(Ω) ∩H1

0 (Ω),(Qr 1,ηϕ

)(y)=ϕ(y)−η(y)·r 1ϕ

ein Projektor auf ker r 1.

Beweis: Sei ϕ ∈ H2(Ω) ∩H10 (Ω) beliebig.

Wegen∂(Qr 1,η ϕ)

∂y(y) = ∂ϕ

∂y(y)− ∂η

∂y(y)·r 1ϕ ist

r 1

(Qr 1,ηϕ

)=

∂ϕ

∂y(0)− r 1ϕ = 0, (39)

36

denn ∂η∂y

(0) = r 1η = 1, so dass im(Qr 1,ηϕ

) ⊂ ker r 1 folgt.Ist andererseits ϕ ∈ ker r 1, so ist r 1ϕ = 0, was Qr 1,ηϕ = ϕ impliziert. Wirhaben also im Qr 1,η = ker r 1 gezeigt.

Ist wieder ϕ ∈ H2(Ω) ∩H10 (Ω) beliebig, so folgt wegen (39)

(Q2

r 1,ηϕ)(y) =

(Qr 1,ηϕ

)(y)− η(y)·r 1

(Qr 1,ηϕ

)=

(Qr 1,ηϕ

)(y)

fur jedes y ∈ Ω, so dass Qr 1,η ein Projektor ist. ¤

In Satz 19 wurde gezeigt, dass fur den ADAS-Index stets µ ≥ 1 gilt. DieseAussage kann nun verscharft werden.

Satz 25 Es gelte (V1)-(V5) und es sei σ2 > 0. Der ADAS-Index µ von (34)ist genau dann µ = 1, wenn das Netzwerk keine CVH-Schleife enthalt.

Das ursprungliche RCL-Netzwerk hat den Index 1. Wegen (V2) und Satz 14enthalt es daher keine CV -Schleifen und keine LI-Schnittmengen. Im Satz 25wird nun behauptet, dass sich der Index 1 genau dann auf das ADA-Systemubertragt, wenn durch die

”Halbleiterelemente“ keine CVH-Schleifen geschlos-

sen werden. In diesem Fall existieren auch keine LH- und LIH-Schnittmengen.

Falls spater auch Index-2-Netzwerke zugelassen werden sollen, ist die Nicht-existenz dieser Schnittmengen sicherlich als zusatzliche Bedingung zu fordern,um ein Index-1-Kriterium fur den ADAS-Index angeben zu konnen.

Der Beweis des Satzes 25 wird im Abschnitt 5.4 gefuhrt.

Sobald das Netzwerk eine CVH-Schleife enthalt, ist der ADAS-Index µ ≥ 2.Im Abschnitt 5.4 wird gezeigt, dass in diesem Fall dim imW1 > 0 fur jedenProjektor W1 langs imG1 gilt.

Der nachste Satz sichert, dass (34) niemals einen ADAS-Index µ ≥ 3 besitzenkann.

Satz 26 Es gelte (V1)-(V5) und es sei σ2 > 0.Fur den ADAS-Index von (34) gilt dann stets µ ≤ 2.

Dieser Satz wird ebenfalls im Abschnitt 5.4 bewiesen.

Wir haben damit Abstrakte Differential-Algebraische Systeme (34) aus derSchaltungsimulation vollstandig bzgl. ihres ADAS-Index klassifiziert.

Die Aussagen der Satze 19, 25 und 26 sollen daher nochmals im Satz 27 zu-sammengefasst werden.

37

Satz 27 Es gelte (V1)-(V5) und es sei σ2 > 0.

(i) Der ADAS-Index von (34) betragt stets µ ≥ 1.

(ii) Der ADAS-Index ist genau dann µ = 1, wenn das Netzwerk keine CVH-Schleifen enthalt.

(iii) In allen anderen Fallen betragt der ADAS-Index µ = 2. ¤

5.3.3 Zusammenfassung

Wir betrachten in Abbildung 3 erneut die Beispielschaltung aus der Einleitung.In der Praxis werden Halbleiterelemente wie hier die Diode durch Ersatzschalt-bilder angenahert. Der Detaillierungsgrad des Ersatzschaltbildes kann dabeiden Index des resultierenden Netzwerkes mit beeinflussen [6].

Abbildung 3: Approximation mit Ersatzschaltbildern

Wird die Halbleiterdiode lediglich als nichtlineare Stromquelle interpretiert,so enthalt das resultierende Netzwerk keine CV -Schleife. Folglich fuhrt dasErsatzschaltbild auf ein Index-1-Netzwerk (Satz 14). Sobald aber kapazitiveEffekte mitsimuliert werden, ergibt sich hier ein Index-2-Netzwerk, da der Kon-densator des Kompaktmodells eine CV -Schleife schließt.

Die Satze 20 und 27 zeigen ein ganz analoges Verhalten der ADA-Systemein der Schaltungssimulation. Wird auf die Simulation kapazitiver Effekte ver-zichtet (σ2 = 0), so betragt der ADAS-Index µ = 1, wenn nur die Matrix

38

l·IkH− σ1E

T Q0G−11 E regular ist (Satz 20). Bei Beachtung des

”Driftstromes“

(σ2 > 0) kann es jedoch zu einer Indexerhohung kommen (Satz 27).

Im Satz 20 musste die Regularitat der Matrix l·IkH−σ1E

T Q0G−11 E gefordert

werden. Im Fall σ2 = 0 konnten daher Eigenschaften des Bauelementes, hierin Form der Lange l, Einfluss auf den Index haben (Bemerkung 22).

Die Driftstromkomponente sollte jedoch schon aus physikalischen Uberlegun-gen stets mit beachtet werden (Kapitel 4), so dass auf weitere Untersuchungenhinsichtlich des Index bei Aufgaben mit σ2 = 0 verzichtet wurde.

Schließlich zeigt der Satz 27, dass im Fall σ2 > 0 der ADAS-Index ausschließlichvon der Topologie des Netzwerkes abhangt. Parameter des einzelnen Halblei-terelementes spielen keine Rolle fur den Index. Einzig die Lage der Elementeim Netzwerk ist von Bedeutung. Damit konnen vorhandene, effiziente topo-logische Index-Tests, wie sie z.B. im Schaltungssimulator Titan der InfineonTechnologies AG implementiert sind, problemlos auch fur ADA-Systeme in derSchaltungssimulation erweitert werden.

Der Satz 27 stellt eine naturliche Fortsetzung des Satzes 14 dar, so dass sichdas Konzept des ADAS-Index harmonisch in die Welt der Schaltungssimulationeingliedert. Ein ADAS-Index µ < 1 ist schon wegen dem Vorhandensein derquasi-stationaren Poisson-Gleichung ausgeschlossen. Die Falle µ = 1 und µ = 2stimmen mit den Erfahrungen aus der Modellierung mit Ersatzschaltbildernuberein.

Hinzu kommt lediglich die Tatsache, dass schon Schaltungen mit reinen CH-oder H-Schleifen zu einem ADAS-Index µ = 2 fuhren (vgl. Abbildung 7 imAbschnitt 5.5). Das oben vorgestellte, aus nichtlinearer Stromquelle und Kon-densator bestehende Ersatzschaltbild fur Halbleiterdioden wurde in diesem Fallnoch auf den Index 1 fuhren, da die Strome durch die zugehorigen Kondensato-ren bei der Modifizierten Knotenanalyse nicht zum Vektor x der Unbekanntenzahlen wurden.

Allerdings durfen Analogien zwischen dem hier vorgestellten PDAE-Systemund der realen Schaltungssimulation in Anbetracht der Unzulanglichkeiten desverwendeten PDE-Modells nicht uberstrapaziert werden. Untersuchungen mitdem vollstandigen Drift-Diffusionsmodell mussen zeigen, inwieweit sich derBegriff des ADAS-Index, dann fur nichtlineare Systeme, tatsachlich bewahrt.

39

5.4 Beweise der Hauptresultate

In diesem Abschnitt werden Beweise der Satze 25 und 26 geliefert. Im gesamtenAbschnitt wird daher (V1)-(V5) vorausgesetzt. Zunachst werden die Aussagenzur Existenz der fur den ADAS-Index relevanten CVH-Schleifen auf die Un-tersuchung von Kernen linearer Abbildungen zuruckgefuhrt (Lemma 29).

Damit hangt die Existenz von CVH-Schleifen von der Singularitat bestimmterMatrizen ab (Lemma 29 und 34), die ihrerseits wiederum die Injektivitat unddie dichte Losbarkeit der Operatoren Gi bestimmen.

Wir betrachten wieder das RCL-Netzwerk aus Abschnitt 3.2 mit n Knoten undder elementbezogenen Inzidenzmatrix (ARACALAV AI). Außerdem sei AH dieInzidenzmatrix der

”Halbleiterelemente“. Wie zuvor sei QC ∈ Lb(IR

n−1, IRn−1)ein Projektor auf ker AT

C . Weiter sei QV−C ∈ Lb(IRn−1, IRn−1) ein Projektor

auf ker ATV QC .

In endlichdimensionalen Raumen konnen lineare Abbildungen mit den sie dar-stellenden Matrizen identifiziert werden. Das heißt, wir betrachten QC undQV−C als Elemente von MIR(n− 1, n− 1).

Lemma 28 Fur die obigen Projektoren gilt

(i) im QC = ker ATC , ker QT

C = im AC ,

(ii) im QV−C = ker ATV QC , ker QT

V−C = im QTCAV .

Beweis:

(i) Nach Definition ist im QC = ker ATC . Wegen QT

CAC = (ATCQC)T = 0 ist

im AC ⊂ ker QTC . Die Satze B.1 bis B.4 aus dem Anhang zeigen

dim im AC = n− 1− dim(im AC)⊥ = n− 1− dim ker ATC

= n− 1− dim im QC = n− 1− rang(QC)

= n− 1− rang(QTC) = n− 1− dim im QT

C = dim ker QTC ,

so dass im AC = ker QTC folgt.

(ii) Wiederum gilt im QV−C = ker ATV QC nach Definition, was im QT

CAV ⊂ker QT

V−C impliziert. Wie bei (i) wird dim im QTCAV = dim ker QT

V−C be-rechnet, so dass ker QT

V−C = im QTCAV ist. ¤

Wir geben nun ein Kriterium fur die Existenz von CVH-Schleifen an, das mitTechniken aus [28] bewiesen wird:

40

Lemma 29 Das obige Netzwerk enthalt genau dann eine CVH-Schleife, wennker QT

V−CQTCAH 6= 0 ist.

Beweis:

(⇒) Das Netzwerk enthalte eine CVH-Schleife. Nach Lemma 9 sind die Spal-

ten von (ACAV AH) linear abhangig, d.h. es existiert(

abc

)∈ IRkC+kV +kH

mit ACa + AV b + AHc = 0. Da jede CVH-Schleife mindestens einenH-Zweig enthalt, kann c 6= 0 angenommen werden. Multiplikation mitQT

V−CQTC von links bringt QT

V−CQTCAHc = 0, also c ∈ ker QT

V−CQTCAH .

(⇐) Ist umgekehrt c ∈ ker QTV−CQT

CAH ein nicht-triviales Kernelement, so istQT

CAHc ∈ ker QTV−C = im QT

CAV (Lemma 28). Sei also b ∈ IRkV so, dassQT

CAHc = QTCAV b. Wieder wegen Lemma 28 ist AHc− AV b ∈ ker QT

C =im AC und es existiert a ∈ IRkC mit ACa + AV b − AHc = 0, so dass dieSpalten von (ACAV AH) linear abhangig sind. Die zugehorigen Zweigebilden nach Lemma 9 eine Schleife, die wegen c 6= 0 mindestens einenH-Zweig enthalt, so dass eine CVH-Schleife existiert. ¤

Die Matrix QTV−CQT

CAH wollen wir mit MH bezeichnen. Analog sei MR =QT

V−CQTCAR.

Lemma 30 Fur die Matrix MR = QTV−CQT

CAR gilt:

(i) QMR= In−1 −QCQV−C ist ein Projektor auf ker MT

R und

(ii) HMR= MRGMT

R + QTMR

QMRist positiv definit.

Dabei ist G = diag(G1, . . . , GkG) die Leitwertmatrix des RCL-Netzwerkes.

Beweis:

(i) Wir betrachten PC = In−1 − QC und PV−C = In−1 − QV−C . Wegenim PC ⊂ ker AT

V QC = im QV−C = ker PV−C ist PV−CPC = 0, so dassQV−CPC = PC gilt. Damit ist PMR

= In−1 −QMR= QCQV−C wegen

P 2MR

= (QCQV−C)2 = QCQV−C(In−1 − PC)QV−C

= QCQV−C −QCQV−CPCQV−C = QCQV−C = PMR

ein Projektor, so dass auch QMRein Projektor ist.

Wegen im QMR= ker PMR

ist noch ker MTR = ker QCQV−C zu zeigen.

41

Die Inklusion ker QCQV−C ⊂ ker MTR ist offensichtlich. Sei a ∈ ker MT

R

beliebig. Dann ist b = QCQV−Ca ∈ ker ATR. Wegen im QC = ker AT

C

und im QV−C = ker ATV QC ist b ∈ ker AT

C ∩ ker ATV . Aus (V2) folgt nun

QCQV−Ca = b = 0, denn b ∈ ker(ACARAV )T = 0. Somit gilt auchker MT

R ⊂ ker QCQV−C .

(ii) Fur beliebiges z ∈ IRn−1 ist

zT HMRz = (MT

Rz)T G(MTRz) + |QMR

z|22 ≥ 0,

da die Leitwertmatrix G nach (V1) positiv definit ist. Damit ist auch

zT HMRz = 0 ⇔ MT

Rz = 0 und |QMRz|2 = 0

⇔ z ∈ ker MTR = im QMR

und z ∈ ker QMR

⇔ z = 0,

da fur QMRals Projektor im QMR

∩ ker QMR= 0 gilt. ¤

Aufgrund der Voraussetzung (V2) ist die zur Netzwerk-DAE gehorende Matrix

G1 =

(ACCAT

C+ARGATRQC 0 AV

−ATLQC L 0

ATV QC 0 0

)

regular, so dass die Inverse G−11 existiert. Passend zur Blockstruktur von G1

schreiben wir

G−11 =

(X11 X12 X13X21 X22 X23X31 X32 X33

). (40)

Die Matrix HMRliefert eine Darstellung fur die Komponente X11 von G−1

1 .

Lemma 31 Mit den obigen Bezeichnungen ist

G−11 =

X11 0 X13

L−1ATLQCX11 L−1 L−1AT

LQCX13

XT13Q

TC 0 X33

(41)

mit

QCX11 = H−1MR

QTV−CQT

C = QCQV−CH−1MR

QTV−CQT

C = (QCX11)T . (42)

Beweis: Aus G1G−11 = Im = G−1

1 G1 folgt fur j = 1, . . . , 3(ACCAT

C + ARGATRQC

)X1j+AV X3j = δ1jIn−1 (43)

−ATLQCX1j+ LX2j = δ2jIkL

(44)

ATV QCX1j = δ3jIkV

(45)

42

und fur i = 1, . . . , 3 ergeben sich die Beziehungen

Xi1

(ACCAT

C + ARGATRQC

)−Xi2ATLQC+Xi3A

TV QC = δi1In−1 (46)

Xi2L = δi2IkL(47)

Xi1AV = δi3IkV, (48)

wobei δij das Kronecker-Symbol ist. Die Induktivitatsmatrix L ist nach (V1)regular. Daher liefert (47)

X12 = 0, X22 = L−1, X32 = 0,

so dass aus (44)

X21 = L−1ATLQCX11, X23 = L−1AT

LQCX13

folgt. Weiterhin bringen (45) und (48) fur i = j = 3

ATV QCX13 =IkV

=ATV XT

31 ⇒ AV

(QCX13 −XT

31

)=0 ⇒ X31 =XT

13QTC ,

denn nach Voraussetzung (V2) ist ker AV = 0.Nun sei j = 1. Aufgrund von (45) ist

X11 = QV−CX11.

Nach Multiplikation mit QTV−CQT

C von links lautet (43) somit

MRGMTRX11 = (QCQV−C)T = P T

MR,

wobei PMR= In−1 − QMR

= QCQV−C der komplementare Projektor zu QMR

ist. Das Lemma 30 liefert die Darstellung

MRGMTR = HMR

−QTMR

QMR= HMR

(In−1 −QMR) = HMR

PMR

mit der regularen, positiv definiten Matrix HMR. Es folgt

PMRX11 = H−1

MRP T

MR= H−1

MRQT

V−CQTC ,

so dass

QCX11 = QCQV−CX11 = PMRX11 = H−1

MRP T

MR= H−1

MRQT

V−CQTC

ist. Wegen P 2MR

= PMRist auch QCX11 = QCQV−CH−1

MRQT

V−CQTC .

Aus der Symmetrie der Leitwertmatrix G = diag(G1, . . . , GkG) folgt die Sym-

metrie von HMR. Damit H−1

MRsymmetrisch und es folgt QCX11 = (QCX11)

T .¤

Die Darstellung (41) von G−11 ist zusammen mit (42) ausreichend fur die wei-

teren Uberlegungen. Die Beziehungen (43)-(48) gestatten es jedoch, samtlicheKomponenten von G−1

1 zu bestimmen.

43

Lemma 32 Wir definieren die Matrizen

MV C = ATV QCQT

CAV , S =(In−1 − AV M−1

V CATV QCQT

C

),

MV CR = ATV QCQT

CAR, T =(In−1 − ARGAT

RH−1MR

QTV−CQT

C

),

HAC= ACCAT

C + QTCQC ,

wobei G wiederum die Leitwertmatrix ist. Dann gilt

X11 = H−1AC

S T + H−1MR

QTV−CQT

C , (49)

X31 = M−1V CAT

V QCQTC

(In−1 − ARGAT

RH−1MR

QTV−CQT

C

), (50)

X13 =(In−1 −

(H−1

ACS T + H−1

MRQT

V−CQTC

)ARGAT

R

)QCQT

CAV (M−1V C)T , (51)

X33 = M−1V CMV CRG

(MT

RH−1MR

MR −G−1)GT MT

V CR(M−1V C)T . (52)

Beweis: Wegen (V2) ist ker QTCAV = 0. Aus Lemma 10 folgt die Re-

gularitat von MV C = ATV QCQT

CAV , denn

ker(AT

V QCQTCAV

)= ker(QT

CAV )T In−1(QTCAV ) = ker QT

CAV = 0.

Damit ist insbesondere S korrekt definiert. Wird (43) fur j = 1 von links mitAT

V QCQTC multipliziert, folgt

MV CX31 = ATV QCQT

C

(In−1 − ARGAT

RQCX11

).

(50) ergibt sich, wenn zusatzlich (42) beachtet wird.

Sei PC = In−1 − QC . Nach einiger Rechnung folgt wegen ACCATC = HAC

PC

aus (43) fur j = 1 die Beziehung

HACPCX11 = ACCAT

CPCX11 = ACCATCX11

= In−1 − ARGATRQCX11 − AV X31 = S · T.

Da wie im Lemma 30 gezeigt werden kann, dass auch HACpositiv definit,

also regular ist, haben wir wieder wegen (42) fur X11 die Darstellung (49)gewonnen.

Wird die Gleichung (46) bei i = 3 von rechts mit QCQTCAV multipliziert, ergibt

sich mit (50)

X33MV C = −X31ARGATRQCQT

CAV

= M−1V CMV CRGMT

RH−1MR

MRGMTV CR −M−1

V CMV CRGMTV CR

= M−1V CMV CRG

(MT

RH−1MR

MR −G−1)GT MT

V CR

44

denn X32 = 0 (Lemma 31) und die Leitwertmatrix G ==diag(G1, . . . , GkG) ist

regular und symmetrisch. Da auch MV C regular und symmetrisch ist, wurdesomit auch (52) gezeigt.

Fur i = 1 liefert (46) bei Multiplikation mit QCQTCAV von rechts

X13MV C =(In−1 −X11ARGAT

R

)QCQT

CAV

=(In−1 −

(H−1

AC· S · T + H−1

MRQT

V−CQTC

)ARGAT

R

)QCQT

CAV ,

also (51). ¤

Bemerkung 33 Die im Lemma 32 definierten Matrizen

S =(In−1−AV M−1

V CATV QCQT

C

)und T =

(In−1−ARGAT

RH−1MR

QTV−CQT

C

)

sind Projektoren.

Beweis:

(i) (AV M−1V CAT

V QCQTC)2 = AV M−1

V CMV CM−1V CAT

V QCQTC = AV M−1

V CATV QCQT

C

zeigt, dass AV M−1V CAT

V QCQTC und damit auch S ein Projektor ist.

(ii) Wegen H−1MR

QTV−CQT

C = QCX11 und (46), (48) ist ARGATRH−1

MRQT

V−CQTC

aufgrund von

(ARGATRH−1

MRQT

V−CQTC)2 = ARGAT

RH−1MR

QTV−CQT

CARGATRQCX11

= ARGATRH−1

MRQT

V−CQTC(In−1 − AV X31 − ACCAT

CX11)

= ARGATRH−1

MRQT

V−CQTC − ARGAT

RQCX11AV X31

= ARGATRH−1

MRQT

V−CQTC

ein Projektor. Gleiches gilt somit auch fur T . ¤

Mit der Darstellung (41), (42) von G−11 aus Lemma 31 konnen wir nun die

Struktur der Matrix ET Q0G−11 E untersuchen, die schon im Satz 20 ein Rolle

gespielt hat. Wir zeigen, dass der Kern dieser Matrix mit ker QTV−CQT

CAH

ubereinstimmt.

Wie im Abschnitt 5.1 sei E =(

AH00

)und Q0 =

(QC 0 00 0 00 0 IkV

)sei der Projektor

auf N0 = ker(AD) aus Abschnitt 3.4.

Lemma 34 Es sei MH = QTV−CQT

CAH . Dann ist

(i) ET Q0G−11 E = AT

HQCQV−CH−1MR

QTV−CQT

CAH = MTHH−1

MRMH und

(ii) ET Q0G−11 E ist genau dann regular, wenn ker MH = 0 ist.

45

Beweis:

(i) Mit (41) zeigt eine einfache Matrixmultiplikation

ET Q0G−11 E = ( AT

H 0 0 )

(QC 0 00 0 00 0 IkV

) (X11 X12 X13X21 X22 X23X31 X32 X33

)(AH00

)

= ATHQCX11AH

die Behauptung, denn wegen (42) ist

QCX11 = QCQV−CH−1MR

QTV−CQT

C .

(ii) Mit HMRist auch H−1

MRpositiv definit (Lemma B.6 im Anhang). Die

Darstellung (i) und Lemma 10 liefern (ii). ¤

Wir haben nun alle Bezeichungen und Hilfsaussagen bereitgestellt, die einenBeweis des Satzes 25 ermoglichen. Zur besseren Orientierung wird der Satzhier nochmals formuliert.

Satz 25 Es gelte (V1)-(V5) und es sei σ2 > 0. Der ADAS-Index µ von (34)ist genau dann µ = 1, wenn das Netzwerk keine CVH-Schleife enthalt.

Beweis: Wir betrachten wieder den linearen Operator G1 : X → Y mitder Matrixdarstellung

G1 = AD + BQ0 =

G1 0 E

ε(∆h

)ET Q0 ε∆Qr 0

0 −σ2 r −IkH

und dem Definitionsbereich

DG1 = X = IRm × [H2(Ω) ∩H10 (Ω)]kH × IRkH .

Qr wird dabei mit Hilfe eines beliebigen Projektors Qr 1 auf ker r 1 wie in (38)definiert. Nach Definition 4 ist die Injektivitat und dichte Losbarkeit von G1

zu untersuchen.

(vϕw

)∈ kerG1 ⇔

G1v + Ew = 0ε(∆h

)ET Q0v + ε ∆Qrϕ = 0

−σ2 rϕ− w = 0

⇔

v = σ2G−11 E r ϕ

∆Qrϕ = −σ2

(∆h

)ET Q0G

−11 E r ϕ

w = −σ2 rϕ

(53)

46

Es wurde beachtet, dass G1 wegen (V2) regular ist. Aus (V5) folgt auch hier∆h 6= 0. Zweimalige Ortsintegration der zweiten Gleichung liefert

(Qrϕ

)(y) = −σ2

(h(y)− 1

)ET Q0G

−11 E rϕ ∀ y ∈ Ω, (54)

denn wegen ϕ ∈ [H2(Ω) ∩H10 (Ω)]kH und im Qr = ker r ist

(Qrϕ

)(0) = 0 und

∂(Qrϕ

)

∂y(0) = rQr = 0.

Hier wurde wie im Beweis von Satz 20 auf Seite 32 gerechnet. Außerdem muss(Qrϕ

)(l) = 0 gelten, so dass mit (54)

0 = σ2ET Q0G

−11 E rϕ (55)

folgt (beachte (V5)). Es wurde σ2 6= 0 vorausgesetzt. Somit besitzt (55) genaudann nur die triviale Losung r ϕ = 0, wenn ET Q0G

−11 E regular ist.

Das Lemma 35 wird sogar zeigen, dass G1 genau dann injektiv ist, wennET Q0G

−11 E regular ist.

Insgesamt ergibt sich dann mit den Lemmata 29 und 34

G1 injektiv ⇔ ET Q0G−11 E regular ⇔ ker MH = 0

⇔ es existiert keine CVH-Schleife.

Wie im Beweis von Satz 20 wird gezeigt, dass G1 im Falle der Regularitat vonET Q0G

−11 E dicht losbar ist, was den Beweis des Satzes beendet. ¤

Lemma 35 Es gelte (V1)-(V5) und es sei σ2 > 0. Qr 1 sei ein beliebiger Pro-jektor auf ker r 1. Wenn Qr wie in (38) definiert wird, gilt: G1 ist genau danninjektiv, wenn ker ET Q0G

−11 E = 0 ist.

Beweis:

(i) Es ist ker Qr 1 6= 0.Aus kerQr 1=0 wurde H2(Ω)∩H1

0 (Ω)=imQr 1⊕ker Qr 1=im Qr 1=kerr 1

folgen – ein Widerspruch, denn ϕa(y) = ayl(l− y) ∈ H2(Ω)∩H1

0 (Ω) liegtfur a 6= 0 nicht in ker r 1.

(ii) Sei ϕ∗ ∈ ker Qr 1 mit ϕ∗ 6= 0. Dann ist a∗ = rϕ∗ 6= 0.Ware ϕ∗ ∈ ker Qr 1 \ 0 mit rϕ∗ = 0, so wurde ϕ∗ ∈ ker Qr 1 ∩ ker r 1 =ker Qr 1 ∩ im Qr 1 = 0, also ϕ∗ = 0, folgen – Widerspruch.

47

(⇒) Sei G1 injektiv. Annahme: ker ET Q0G−11 E 6= 0.

Dann existiert a ∈ ker ET Q0G−11 E, a 6= 0. ⇒ ∃ λ ∈ IRkH , λ 6= 0, so

dass a = λ⊗ (a∗ ·1) = λ⊗ (r 1ϕ∗ ·1) = λ⊗ r (ϕ∗ ·1) = r (ϕ∗ ·λ). Hier ist

⊗ : IRkH × IRkH → IRkH , (a, b) 7→ a⊗ b =

(a1b1...

akHbkH

).

die komponentenweise Multiplikation auf IRkH . Fur ϕ∗ ∈ H2(Ω)∩H10 (Ω)

und λ ∈ IRkH ist (ϕ∗·λ)∈ [H2(Ω)∩H10 (Ω)]kH . Wegen (38) ist Qr (ϕ∗·λ) = 0.

Mit va = σ2G−11 E r

(ϕ∗ ·λ

)und wa = −σ2 r

(ϕ∗ ·λ

)ist

(va

ϕ∗·λwa

)∈ kerG1.

Im Widerspruch zur Injektivitat von G1 ist aber(

vaϕ∗·λwa

)6= 0. Folglich

muss doch ker ET Q0G−11 E = 0 gelten.

(⇐) Sei ker ET Q0G−11 E = 0.

(53) und (54) zeigen, dass aus (vT, ϕT, wT )T∈kerG1 sofort v=σ2G−11 E r ϕ,

ϕ ∈ ker Qr , w = −σ2 r ϕ und r ϕ ∈ ker ET Q0G−11 E = 0 folgt. Damit

ist rϕ = 0, v = 0, w = 0 und ϕ ∈ ker Qr ∩ ker r = ker Qr ∩ im Qr = 0,was (vT , ϕT , wT )T = 0 und damit die Injektivitat von G1 zeigt. ¤

Mit Satz 25 ist der Index-1-Fall vollstandig charakterisiert. Wir wollen im Fol-genden daher voraussetzen, dass ET Q0G

−11 E singular ist, also dass das Netz-

werk eine CVH-Schleife enthalt.

In diesem Fall existiert ein z ∈ IRkH , das nicht im Bild von ET Q0G−11 E ent-

halten ist. Da Ω beschrankt ist, liegt die konstante Abbildung Ψ(y) = 2εσ2

l2z

in [L2(Ω)]kH . Es ist∫ l

0

∫ ξ

0Ψ dτdξ = ε σ2 z. Damit folgt

(0Ψ0

)6∈ imG1, denn

andernfalls musste eine Folge(ϕn

)n

in [H2(Ω) ∩H10 (Ω)]kH existieren mit

ET Q0G−11 Erϕn → 1

εσ2

l∫

0

ξ∫

0

Ψ dτdξ = z,

was wegen z 6∈ im ET Q0G−11 E ausgeschlossen ist. Folglich gilt dim imW1 > 0

fur jeden Projektor W1 langs imG1.

Der Beweis von Satz 25 liefert eine Darstellung des Kernes von G1. Wir stellendazu einen Projektor QK auf

K = ϕ ∈ [H2(Ω) ∩H10 (Ω)]kH | rϕ ∈ ker ET Q0G

−11 E = ker MH

bereit.

48

Lemma 36 Sei Qr ein Projektor auf ker r ⊂ [H2(Ω) ∩ H10 (Ω)]kH und R sei

der eindeutig bestimmte Projektor mit

im R = ker MH und ker R = (ker MH)⊥.

Dann ist

QK : [H2(Ω) ∩H10 (Ω)]kH → [H2(Ω) ∩H1

0 (Ω)]kH , ϕ 7→ R + Qr (I − R)

ein Projektor auf K.

Beweis: Nach Lemma 34 ist ET Q0G−11 E = MT

HH−1MR

MH , so dass ausLemma 10 und Satz B.3

im ET Q0G−11 E = im MT

H = (ker MH)⊥ = ker R,

ker ET Q0G−11 E = ker MH = im R

folgt, denn H−1MR

ist positiv definit. R = RT ist ein Orthoprojektor, denn auch

fur RT ist im RT = (ker R)⊥ = ker MH und ker RT = (im R)⊥ = (ker MH)⊥.Bei der Definition von QK fassen wir den Projektor R ∈ Lb(IR

kH ) als Elementvon Lb

([H2(Ω) ∩H1

0 (Ω)]kH)

auf. Auch dort ist R ein Projektor. Fur QK gilt

Q2K =

(R + Qr (I − R)

)(R + Qr (I − R)

)

= R + R Qr (I − R) + Qr (I − R)−Qr R Qr (I − R)

= QK + R Qr (I − R)−Qr R Qr (I − R) = QK,

denn fur jedes ϕ ∈ [H2(Ω)∩H10 (Ω)]kH ist r R Qrϕ = R rQrϕ = 0 (r ist linear),

also R Qrϕ ∈ ker r und somit R Qr = Qr R Qr . Damit ist QK ein Projektor.

Es ist im QK = K zu zeigen.

Ist zunachst ϕ ∈ K, so folgt rϕ ∈ ker MH = im R und ϕ−Rϕ ∈ ker r = im Qr ,d.h. ϕ−Rϕ=Qr

(ϕ−Rϕ

). Damit ist QKϕ=Rϕ+Qr (I−R)ϕ=Rϕ+ϕ−Rϕ=ϕ,

also ϕ ∈ im QK.

Sei andererseits ϕ∈ im QK, d.h. ϕ=QKϕ, so ist rϕ=rQKϕ=r(R+Qr (I−R)

)ϕ

=r Rϕ=R rϕ ∈ ker MH , was ϕ ∈ K bedeutet. ¤

Folgerung 37 Es gelte (V1)-(V5) und es sei σ2 > 0. Dann ist

N1 = kerG1 =(

vϕw

)∈ X

∣∣∣ϕ ∈ ker Qr , r ϕ ∈ ker MH ,

v= σ2G−11 Erϕ

w= −σ2rϕ

.

und

Q1 =

(0 σ2G−1

1 ErQK 00 Pr QK 00 −σ2rQK 0

)

ist ein Projektor auf N1. Dabei ist Qr = I − Pr mit Hilfe eines beliebigenProjektors Qr 1 auf ker r 1 wie in (38) definiert und QK ist der Projektor auf Kaus Lemma 36.

49

Beweis: rϕ kann in (55) jeden Wert rϕ = a ∈ ker MH annehmen, dennaufgrund von Lemma 34 ist

ker ET Q0G−11 E = ker MH .

Mit (53) und (54) ist also

(vϕw

)∈ kerG1 ⇔

v = σ2G−11 E r ϕ

Qrϕ = 0w = −σ2 rϕ

rϕ ∈ ker MH

. (56)

Aus der Beschranktheit von G1 ∈ Lb(X,Y ) folgt die Abgeschlossenheit desKernes kerG1 = kerG1 und N1 = kerG1 lasst sich in der behaupteten Artdarstellen.

Wegen

r (Prϕ) = r (ϕ)− r (Qrϕ) = r (ϕ) ∀ ϕ ∈ [H2(Ω) ∩H10 (Ω)]kH (57)

ist QKPrQKϕ=PrQKϕ fur alle ϕ, denn QKϕ∈K impliziert PrQKϕ∈K=im QK.Damit ist PrQK ein Projektor, denn

(PrQK)2 = PrQKPrQK = PrPrQK = PrQK.

Dies zeigt, dass wegen (57) auch der oben angegebene Operator

Q1 : X → X

ein Projektor ist. Wir zeigen imQ1 = N1 = kerG1.

Sei(

vϕw

)∈ X beliebig. Dann ist

G1Q1

(vϕw

)=

(σ2ErQKϕ−σ2ErQKϕ

εσ2(∆h) ET Q0G−11 ErQKϕ+ε ∆Qr Pr QKϕ

−σ2 rPr QKϕ+σ2rQKϕ

)= 0,

denn es gilt (57), QrPr =0 und wegen rQKϕ∈ker MH ist ET Q0G−11 ErQKϕ=0.

Damit ist imQ1 ⊂ N1. Die umgekehrte Inklusion ist wegen (56) offensichtlich,

denn(

vϕw

)∈N1 impliziert ϕ∈K∩ker Qr =im QK∩im Pr , so dass ϕ=QKϕ=Prϕ

ist. ¤

50

Nun sind wir in der Lage, den Operator G2 zu berechnen, dessen Eigenschaftenden Index-2-Fall charakterisieren.

G2 = G1 + BP0Q1

= G1 +

B 0 E

ε(∆h

)ET ε ∆ 0

0 σ1 r −IkH

0 σ2P0G−11 Er QK 0

0 Pr QK 00 0 0

=

G1 0 E

ε(∆h

)ET Q0 ε∆Qr 0

0 −σ2 r −IkH

+

0 σ2BP0G−11 ErQK 0

0 εσ2 ·(∆h

)ET P0G

−11 ErQK + ε∆Pr QK 0

0 σ1 rQK 0

=

G1 σ2BP0G−11 Er QK E

ε(∆h

)ET Q0 ε∆

(Qr + Pr QK

)+ εσ2 ·

(∆h

)ET P0G

−11 ErQK 0

0 −σ2 r + σ1 rQK −IkH

. (58)

Satz 26 Es gelte (V1)-(V5) und es sei σ2 > 0.Fur den ADAS-Index von (34) gilt dann stets µ ≤ 2.

Beweis:(

vϕw

)∈ kerG2 ⇐⇒

v = σ2G−11 Erϕ−

(σ1Im + σ2G

−11 BP0

)G−1

1 Er QKϕ

∆(Qr+PrQK

)ϕ = −(

∆h)ET Q0v − σ2 ·

(∆h

)ET P0G

−11 ErQKϕ

w = −σ2 r ϕ + σ1 r QKϕ

⇐⇒

v = σ2G−11 Erϕ−

(σ1Im + σ2G

−11 BP0

)G−1

1 Er QKϕ

∆(Qr+PrQK

)ϕ = −σ2 ·

(∆h

)ET Q0G

−11 Erϕ

−σ2 ·(∆h

)ET P0G

−11 Er QKϕ

+(∆h

)ET Q0

(σ1Im + σ2G

−11 BP0

)G−1

1 ErQKϕ

w = −σ2 r ϕ + σ1 r QKϕ

Wegen rQK ∈ ker MH = ker ET Q0G−11 E ist ET Q0G

−11 ErQK = 0. Wir mussen

daher die folgende Gleichung untersuchen:

∆(Qr+PrQK

)ϕ =− σ2 ·

(∆h

)ET Q0G

−11 Er ϕ− σ2 ·

(∆h

)ET G−1

1 ErQKϕ

+ σ2 ·(∆h

)ET Q0G

−11 BP0G

−11 Er QKϕ.

Die Matrix G1 ist als G1 = AD + BQ0 definiert, was G1Q0 = BQ0, alsoG−1

1 BQ0 = Q0 zeigt. Folglich ist G−11 BP0 = G−1

1 B −Q0 und

∆(Qr+PrQK

)ϕ =− σ2 ·

(∆h

)ET Q0G

−11 Erϕ− σ2 ·

(∆h

)ET G−1

1 Er QKϕ

+ σ2 ·(∆h

)ETQ0G

−11 BG−1

1 Er QKϕ−σ2 ·(∆h

)ETQ0G

−11 ErQKϕ

=− σ2 ·(∆h

)ET Q0G

−11 Erϕ

− σ2 ·(∆h

)ET

(Im −Q0G

−11 B

)G−1

1 Er QKϕ,

51

denn wie oben ist ET Q0G−11 Er QKϕ = 0. Auch hier ∆h 6= 0 ∈ L2(Ω).

Wird diese Gleichung zweimal bezuglich des Ortes integriert, erhalten wir we-gen rQr = 0 und rPr = r((Qr+PrQK)ϕ

)(y) =− σ2 ·

(h(y)−1

)ET Q0G

−11 Er ϕ (59)

− σ2 ·(h(y)−1

)ET

(Im−Q0G

−11 B

)G−1

1 Er QKϕ + y r QKϕ,

so dass nun

0 =((Qr+PrQK)ϕ

)(l)

=σ2ET Q0G

−11 Er ϕ+σ2E

T(Im−Q0G

−11 B

)G−1

1 ErQKϕ+l rQKϕ (60)

folgt, denn es ist auch(Qr+PrQK

)ϕ ∈ [H2(Ω) ∩H1

0 (Ω)]kH .

Das nachste Lemma (Lemma 38) wird zeigen, dass aus der Beziehung (60)sofort r ϕ = 0 folgt.

Insbesondere ist daher ϕ ∈ K = im QK, was ϕ = QKϕ bedeutet. Folglich istauch rQKϕ = rϕ = 0. Die obige Darstellung von kerG2 und (59) bringen dannv = 0, w = 0 und

(Qr +PrQK

)ϕ = ϕ − PrPKϕ = 0 mit PK = I − QK. Wegen

ϕ = PrPKϕ ∈ im Pr = ker Qr folgt ϕ ∈ im Qr ∩ ker Qr = 0, also ϕ = 0, sodass G2 stets injektiv ist.

Die dichte Losbarkeit von G2 folgt aus der Surjektivitat der Abbildung

F :[H2(Ω) ∩H10 (Ω)]kH → IRkH ,

ϕ 7→ σ2ET Q0G

−11 Erϕ +

(l IkH

+ σ2ET(Im −Q0G

−11 B

)G−1

1 E)r QKϕ.

Die technischen Details werden im Lemma 40 ausgefuhrt. ¤

Lemma 38 Es gelte (V1)-(V5) und es sei σ2 > 0. Dann gilt die folgendeImplikation:

0 = σ2ET Q0G

−11 Erϕ + l rQKϕ

+σ2ET(Im−Q0G

−11 B

)G−1

1 ErQKϕ

⇒ rϕ = 0 (61)

Beweis: Wir setzen die Gultigkeit der Gleichung

0 = σ2ET Q0G

−11 Er ϕ+l r QKϕ+σ2E

T(Im−Q0G

−11 B

)G−1

1 ErQKϕ (62)

voraus und mussen rϕ = 0 zeigen. Wird (62) von links mit dem Vektor(rQKϕ)T multipliziert, folgt

0 = σ2 · (rQKϕ)T ET Q0G−11 Erϕ + l |rQKϕ|22

+ σ2 · (rQKϕ)T ET(Im −Q0G

−11 B

)G−1

1 Er QKϕ(63)

Die Struktur der hier auftretenden Matrizen wird genauer untersucht:

52

(i) Die Matrix ET Q0G−11 E = MT

HH−1MR

MH wurde im Lemma 34 betrachtet.Wegen QKϕ ∈ K ist rQKϕ ∈ ker MH , so dass

(rQKϕ)T ET Q0G−11 E = (rQKϕ)T MT

HH−1MR

MH

= (MHrQKϕ)T H−1MR

MH = 0

folgt.

(ii) Mit der Formulierung (40) fur G−11 wird

ET G−11 E = ( AT

H 0 0 )(

X11 X12 X13X21 X22 X23X31 X32 X33

)(AH00

)= AT

HX11AH

berechnet.

(iii) In Lemma 31 haben wir die Struktur von G−11 genauer untersucht. Wir

erhalten damit

ET Q0G−11 BG−1

1 E = ( ATH 0 0 )

(QC 0 00 0 00 0 IkV

)G−1

1 B(

X11 0 X13X21 X22 X23X31 0 X33

) (AH00

)

= ( ATHQC 0 0 )

(X11 0 X13X21 X22 X23X31 0 X33

) (ARGAT

R AL AV

−ATL 0 0

ATV 0 0

)( X11AHX21AHX31AH

)

= ( ATHQCX11 0 AT

HQCX13 )

(ARGAT

RX11AH+ALX21AH+AV X31AH

−ATLX11AH

ATV X11AH

)

= ATHQCX11ARGAT

RX11AH + ATHQCX13A

TV X11AH

+ATHQCX11ALX21AH + AT

HQCX11AV X31AH

Wegen (48) ist X11AV = 0. Außerdem wurde im Lemma 31

X21 = L−1ATLQCX11, QCX13 = XT

31 und QCX11 = (QCX11)T

gezeigt, so dass die Darstellung

ET Q0G−11 BG−1

1 E = ATH(QCX11)

T ARGATRX11AH

+ ATHXT

31ATV X11AH

+ ATH(QCX11)

T ALL−1ATLQCX11AH

gilt. Wird (43) aus Lemma 31 beachtet15,

(QCX11)T ARGAT

R =(ARGATRQCX11)

T

=(In−1 − AV X31 − ACCATCX11)

T

=In−1 −XT31A

TV −XT

11ACCATC ,

15Wir betrachten RCL Netzwerke, so dass die Kapazitatsmatrix C und die LeitwertmatrixG als Diagonalmatrizen symmetrisch sind.

53

folgt sogar

ET Q0G−11 BG−1

1 E = ATHX11AH − AT

HXT11ACCAT

CX11AH

+ ATH(QCX11)

T ALL−1ATLQCX11AH .

Mit (i)-(iii) nimmt (63) die Gestalt

0 = l |r QKϕ|22 + σ2 · (ATCX11AHr QKϕ)T C(AT

CX11AHrQKϕ)

− σ2 · (ATLQCX11AHr QKϕ)T L−1(AT

LQCX11AHrQKϕ)

an. Diese Gleichung kann weiter vereinfacht werden. Wegen der Eigenschaftr QKϕ ∈ ker MH und der Darstellung QCX11 = H−1

MRQT

V−CQTC aus Lemma 31

ist

QCX11AHrQKϕ = H−1MR

QTV−CQT

CAHrQKϕ = H−1MR

MHrQKϕ = 0.

Da nach Voraussetzung l > 0 und σ2 > 0 ist, folgt somit aus

0 = l |rQKϕ|22 + σ2 · (ATCX11AHrQKϕ)T C(AT

CX11AHrQKϕ)

wegen der positiven Definitheit der Kapazitatsmatrix C

ATCX11AHr QKϕ = 0 und rQKϕ = 0.

Die Gleichung (62) reduziert sich somit auf

0 = σ2ET Q0G

−11 Er ϕ,

was rϕ ∈ ker ET Q0G−11 E = ker MH bedeutet. Damit ist ϕ ein Element des

Unterraumes K = im QK von [H2(Ω) ∩ H10 (Ω)]kH . Wir haben also ϕ = QKϕ

gezeigt, was wegen

rϕ = r QKϕ = 0

den Beweis des Lemmas abschließt. ¤

Bemerkung 39 In Lemma 38 ist

Q = Q0G−11 B

der kanonische Projektor auf N0 [8]. Fur den Projektor P = Im −Q0G−11 B ist

also ker P = N0. Lemma 38 zeigt, dass die Implikation

0 = σ2ET Q0G

−11 Er ϕ +

(l·IkH

+ σ2ET PG−1

1 E)rQKϕ ⇒ rϕ = 0

stets richtig ist.

54

Lemma 40 Es gelte (V1)-(V5) und es sei σ2 > 0. Dann ist die Abbildung

F : [H2(Ω) ∩H10 (Ω)]kH → IRkH ,

ϕ 7→ σ2ET Q0G

−11 Erϕ +

(l IkH

+ σ2ET PG−1

1 E)rQKϕ

surjektiv. Dabei ist P = Im −Q0G−11 B der kanonische Projektor langs N0.

Beweis:

Wie im Lemma 36 sei R = RT der eindeutig bestimmte Projektor mit

im R = ker MH = ker ET Q0G−11 E und

ker R = (ker MH)⊥ = im ET Q0G−11 E.

Sei z ∈ IRkH beliebig. Wir geben ein ϕ ∈ [H2(Ω) ∩H10 (Ω)]kH mit Fϕ = z an.

(i) Nach Voraussetzung ist l > 0 und σ2 > 0, so dass

M = l·IkH+ σ2(A

TCX11AH)T C(AT

CX11AH)

wegen der positiven Definitheit der Kapazitatsmatrix C selbst positivdefinit ist. Mit Lemma 10 folgt

im RTMR = im RT = im R und ker RTMR = ker R,

so dass RTMR : im R → im R bijektiv ist. Damit existiert genau einz1 = Rz1 ∈ im R mit RTMRz1 = Rz.

(ii) Nach Konstruktion ist σ2ET Q0G

−11 E : ker R → ker R bijektiv und es

existiert genau ein z2 ∈ ker R mit

σ2ET Q0G

−11 Ez2 = (IkH

−R)z − (IkH−R)MRz1.

(iii) Sei z = z1 + z2. Dann ist Rz = z1 und (IkH− R)z = z2. Da aus den

Berechnungen im Beweis von Lemma 38 folgt, dass(l·IkH

+ σ2ET PG−1

1 E)R

=(l·IkH

+ σ2(ATCX11AH)T C(AT

CX11AH))R

− σ2(ATLQCX11AH)T L−1(AT

LQCX11AH)R = MR

gilt, erhalten wir bei Verwendung von (i) und (ii)

σ2ET Q0G

−11 Ez +

(l·IkH

+ σ2ET PG−1

1 E)Rz

= σ2ET Q0G

−11 E(IkH

−R)z + MRz

= σ2ET Q0G

−11 Ez2 + (IkH

−R)MRz1 + RMRz1

= (IkH−R)z + RTMRz1 = (IkH

−R)z + Rz = z.

55

(iv) Seien (b1, . . . , br) und (b1, . . . , bs) Basen von im R bzw. ker R. Da derRestriktionsoperator r : [H2(Ω) ∩ H1

0 (Ω)]kH → IRkH nach Lemma 16surjektiv ist, existieren Elemente ϕi, ϕj ∈ [H2(Ω) ∩H1

0 (Ω)]kH mit

rϕi = bi, i = 1, . . . , r, und r ϕj = bj, j = 1, . . . , s.

Ohne Einschrankung kann ϕj ∈ ker QK gewahlt werden, denn auch furPKϕj ist

rPKϕj = r(ϕj − Rϕj −Qr (I − R)ϕj

)

= r ϕj − Rr ϕj = bj − Rbj = bj.

Hier wurde die Darstellung von QK aus Lemma 36 genutzt.

Es existieren αi, αj ∈ IR, so dass z1 =∑r

i=1 αibi und z2 =∑s

j=1 αj bj.

Wir setzen ϕ =∑r

i=1 αiϕi +∑s

j=1 αjϕj. Dann ist ϕ ∈ [H2(Ω)∩H10 (Ω)]kH

mit

rϕ =r∑

i=1

αirϕi +s∑

j=1

αjr ϕj =r∑

i=1

αibi +s∑

j=1

αj bj = z1 + z2 = z

r QKϕ =r∑

i=1

αirQKϕi +s∑

j=1

αjrQKϕj =r∑

i=1

αirϕi = z1 = R z

Mit (iii) folgt schließlich die Behauptung: Fϕ = z. ¤

5.5 Beispiele

In diesem Abschnitt werden die Aussagen des Abschnitts 5.3.2 anhand einerkleinen Beispielschaltung verdeutlicht. Insbesondere werden die fur die Beweiserelevanten Matrizen explizit berechnet.

Die graphische Darstellung der Bauelemente ist Tabelle 2 zu entnehmen. Umeinfache Zahlenbeispiele zu erhalten, wird stets 1 als Kenngroße der Bauele-mente gewahlt.

Im Mittelpunkt soll die Beispielschaltung aus Abbildung 4 stehen. Generellwollen wir σ1 = σ2 = l = 1 voraussetzen. Die Knoten der Schaltung wurden mitarabischen Ziffern nummeriert. Das Netzwerk wird durch die Inzidenzmatrix

(1 0

−1 00 10 00 0︸ ︷︷ ︸

AC

1 0 00 1 00 −1 00 0 10 0 −1︸ ︷︷ ︸

AR

010

−10︸︷︷︸

AV

00001︸︷︷︸

AL

)

56

linearer zeitinvarianter Widerstand (G = 1 Ω−1)

linearer zeitinvarianter Kondensator (C = 1 F)

lineare zeitinvariante Induktivitat (L = 1 H)

unabhangige Spannungsquelle (V )

”Halbleiterelement“ (H)

Tabelle 2: Graphische Darstellung der Bauelemente

G

C

L

V

1

2

3

4

1

C

G

G2 5

3

1

1

2

1

Abbildung 4: Die Beispielschaltung.(DAE-Index 1)

G

C

L

V

1

2

3

4

1

C

G

G2 5

3

1

1

2

1

H1 H2

Abbildung 5: H1 und H2 schließen kei-ne CVH-Schleifen. (ADAS-Index 1)

beschrieben. Die modifizierte Knotenanalyse (MNA) fuhrt auf die Algebro-Differentialgleichung0BBBBBBB@

1 −1 0 0 0 0 0−1 1 0 0 0 0 0

0 0 1 0 0 0 00 0 0 0 0 0 00 0 0 0 0 0 00 0 0 0 0 1 00 0 0 0 0 0 0

1CCCCCCCA0BBBBBBB@

e1

e2

e3

e4

e5

iLiV

1CCCCCCCA′

+

0BBBBBBB@1 0 0 0 0 0 00 1 −1 0 0 0 10 −1 1 0 0 0 00 0 0 1 −1 0 −10 0 0 −1 1 1 00 0 0 0 −1 0 00 1 0 −1 0 0 0

1CCCCCCCA0BBBBBBB@

e1

e2

e3

e4

e5

iLiV

1CCCCCCCA =

0BBBBBBB@000000vs

1CCCCCCCA.

Auf die Angabe des Argumentes t wurde verzichtet. vs beschreibt die un-abhangige Spannungsquelle. Es bezeichne (e1, . . . , e5) die kanonische Basis16

im IR5. Wegen ker ATC = span

((e1 + e2), e4, e5

)ist

QC =

0BBB@1/2 1/2 0 0 01/2 1/2 0 0 00 0 0 0 00 0 0 1 00 0 0 0 1

1CCCA16Um Verwechslungen mit den Knotenpotentialen e zu vermeiden, werden die Basisvek-

toren mit aufrechten Buchstaben im Fettdruck bezeichnet. Es ist ei = (eij) mit eij = δij .

57

ein Projektor auf ker ATC . Mit den Bezeichnungen aus Abschnitt 3.4 folgt

Q0 =

0BBBBBBB@1/2 1/2 0 0 0 0 01/2 1/2 0 0 0 0 00 0 0 0 0 0 00 0 0 1 0 0 00 0 0 0 1 0 00 0 0 0 0 0 00 0 0 0 0 0 1

1CCCCCCCA und G1 =

0BBBBBBB@3/2 −1/2 0 0 0 0 0−1/2 3/2 0 0 0 0 1−1/2 −1/2 1 0 0 0 0

0 0 0 1 −1 0 −10 0 0 −1 1 0 00 0 0 0 −1 1 0

1/2 1/2 0 −1 0 0 0

1CCCCCCCA.

Das Netzwerk hat den Index µ = 1, denn G1 ist wegen det(G1) = −2 regular.Wir berechnen

ATV QC =

(1/2 1/2 0 −1 0

),

so dass

QV−C =

0BBB@5/6 −1/6 0 1/3 0

−1/6 5/6 0 1/3 00 0 1 0 0

1/3 1/3 0 1/3 00 0 0 0 1

1CCCAals Projektor auf ker AT

V QC gewahlt werden kann.

Beim Hinzufugen von”Halbleiterelementen“ zwischen den Knoten 1 und 5

bzw. 3 und 5 werden keine CVH-Schleifen geschlossen (Abbildung 5). Furdieses System liefert Satz 25 daher den ADAS-Index 1. Tatsachlich ist wegen

AH =

0BBB@1 00 00 −10 0

−1 1

1CCCA, E =

AH

00

,

die Matrix

ET Q0G−11 E =

(1 −1−1 3/2

)

regular. Die Matrix

MH =

0BBB@1/3 01/3 00 0

1/3 0−1 1

1CCCAbesitzt den gleichen Kern (Lemma 34). Es ist ker MH = 0.In Abbildung 6 schließt das Element H1 eine CVH-Schleife. Nach Satz 25betragt der ADAS-Index dieser Schaltung µ ≥ 2, so dass mit Satz 26 µ = 2folgt. Es ist

AH =

0BBB@1 00 00 −1

−1 00 1

1CCCA, ET Q0G−11 E =

(0 00 3/2

)und MH =

0BBB@0 00 00 00 00 1

1CCCA,

58

so dass der Kern ker MH = ker ET Q0G−11 E = span

(( 1

0 ))

nicht trivial ist.

Auch wenn H1 lediglich eine CH-Schleife schließt (Abbildung 7), betragt derADAS-Index µ = 2, denn wegen

AH =

0BBB@1 0

−1 00 −10 00 1

1CCCA und MH =

0BBB@0 00 00 00 00 1

1CCCAist auch hier ker MH nicht trivial.

Abschließend soll explizit nachgerechnet werden, dass der ADAS-Index derSchaltungen aus Abbildung 6 und 7 jeweils 2 betragt. Dazu ist die Gultigkeitder Implikation (61) aus Lemma 38 zu prufen.

In beiden Fallen ist ker MH = span(( 1

0 )). Es folgt

K = ϕ ∈ [H2(Ω) ∩H10 (Ω)]kH

∣∣ rϕ ∈ ker MH = ϕ = ( ϕ1

ϕ2 ) ∈ [H2(Ω) ∩H10 (Ω)]kH

∣∣ r 1ϕ2 = 0 .

r operiert komponentenweise auf Vektoren ϕ ∈ [H2(Ω)]kH . Zur leichterenOrientierung soll r 1 die skalare Version von r bezeichnen. Es sei Qr 1 ein Pro-jektor auf ker r 1. Dann ist

QK =(I 0

0 Qr 1

)

ein Projektor auf K.

Entsprechend der Bemerkung 39 gelte

0 = ET Q0G−11 Erϕ +

(I2 + ET PG−1

1 E)r QKϕ.

G

C

L

V

1

2

3

41

C

G

G2 5

3

1

1

2

1

H1

H2

Abbildung 6: H1 schließt eineCVH-Schleife. (ADAS-Index 2)

G

C

L

V

1

2

3

41

C

G

G2 5

3

1

1

2

1

H1

H2

Abbildung 7: H1 schließt CH-Schleife.(ADAS-Index 2)

59

fur die Schaltung aus Abbildung 6. Dann ist

0 = ET Q0G−11 Er ϕ +

(I2 + ET (I5 −Q0G

−11 B)G−1

1 E)rQKϕ

=(

0 00 3/2

)(r 1ϕ1

r 1ϕ2

)+

(2 −1/2

−1/2 −3/4

)(r 1ϕ1

r 1Qr 1ϕ2

)=

(2 r 1ϕ1

32r 1ϕ2 − 1

2r 1ϕ1

),

so dass r ϕ = ( r 1ϕ1r 1ϕ2 ) = 0 folgt. Satz 26 liefert den ADAS-Index 2 fur diese

Beispielschaltung.

Analog wird im Fall der Schaltung aus Abbildung 7 gerechnet. Auch hier ist

ET Q0G−11 E =

(0 00 3/2

)und

I2 + ET (I5 −Q0G−11 B)G−1

1 E =(

2 −1/2−1/2 −3/4

),

so dass der ADAS-Index 2 betragt.

60

6 Thesen

1. Es wurden PDAE-Netzwerkmodelle der elektrischen Schaltungssimula-tion untersucht. Diese Partial Differential-Algebraic Equations sind ge-koppelte Systeme aus Netzwerk- und Bauelementsimulation. Der Netz-werkteil wird durch die Algebro-Differentialgleichung (DAE)

A(Dx(t)

)′+ Bx(t) + EiH = q(t), t ∈ I = [t0, T ], (64)

modelliert. (64) habe den Index 1. Halbleiterbauelemente werden mitpartiellen Differentialgleichungen (PDEs) wie dem Drift-Diffusionsmo-dell beschrieben. Um das Gesamtsystem als Abstraktes Differential-Alge-braisches System behandeln zu konnen, wird zu dem linearen Modell

ε ∆ψ(y, t) = f(y, t), y ∈ ΩkH = (0, l)kH , t ∈ I,

ψ(y, t) = ψbi(y) + ψE(y, t) y ∈ (∂Ω)kH , t ∈ I,(65)

ubergegangen.

2. Die Formulierung der Kopplung von (64) und (65) uber Dirichlet-Rand-bedingungen

ψ(y, t) = ψbi(y) + ψE(y, t), y ∈ (∂Ω)kH , t ∈ I,

und die Kopplungsgroße

iH(t) = σ1 rψ(·, t)− σ2 r ψt(·, t). (66)

ist wesentlich fur die Untersuchung des Gesamtsystems. Es wird

ψ(0, t) = ψbi(0) + ET x(t), t ∈ I,

ψ(l, t) = ψbi(l), t ∈ I,

vorgeschlagen, da iH(t) nur von der Differenz ψE(0, t)−ψE(l, t)=ET x(t)abhangt.

3. Die gekoppelten Gleichungen (64), (65) und (66) konnen als AbstraktesDifferential-Algebraisches System

A(Du(t))′

+ Bu(t) = q(t), t ∈ I = [t0, T ], (67)

formuliert werden. Dabei sind A, D und B Operatoren in reellen Hilbert-Raumen. Fur (67) wird (V1)-(V5) vorausgesetzt (Seite 16 und 28).

61

4. Der ADAS-Index von (67) hangt vom Detaillierungsgrad der Modellie-rung ab.

1) iH(t) = σ1 rψ(·, t)Bei Vernachlassigung kapazitiver Effekte ist der ADAS-Index stets 1,wenn die Matrix l·IkH

− σ1ET Q0G

−11 E regular ist.

2) iH(t) = σ1 rψ(·, t)− σ2 r ψt(·, t)Werden kapazitiver Effekte mitsimuliert, hangt der ADAS-Indexausschließlich von der Topologie des Netzwerkes ab (vgl. 5.).

5. Bei Verwendung der Kopplungsgroße (66) mit σ2 > 0, also bei Beachtungkapazitiver Effekte, konnte gezeigt werden:

1) Der ADAS-Index betragt stets µ ≥ 1.

2) Der ADAS-Index ist genau dann µ = 1, wenn das Netzwerk keineCVH-Schleifen enthalt.

3) In allen anderen Fallen betragt der ADAS-Index µ = 2.

6. Die Ergebnisse aus 5. stimmen mit den Erfahrungen aus der Modellierungmit Kompaktmodellen uberein. Damit fugt sich das Konzept des ADAS-Index harmonisch in die Welt der Schaltungssimulation ein.

7. Untersuchungen des vollstandigen Drift-Diffusionsmodells im Rahmenvon nichtlinearen Abstrakten Differential-Algebraischen Systemen mus-sen zeigen, wie sich der ADAS-Index bei realen Schaltungen bewahrt.

62

A Funktionalanalytische Grundlagen

Definition A.1 (Produktraum) Seien (X, ‖ · ‖X) und (Y, ‖ · ‖Y ) normierteRaume. Dann heißt der normierte Raum (X × Y, ‖ · ‖X×Y ) mit

‖ · ‖X×Y : X × Y → IR, (x, y) 7→ ‖(x, y)‖X×Y = ‖x‖X + ‖y‖Y .

Produktraum von X und Y .

Sind(X, (·, ·)X

)und

(Y, (·, ·)Y

)Hilbert-Raume mit entsprechenden Skalarpro-

dukten (·, ·)X und (·, ·)Y , so ist(X × Y, (·, ·)X×Y

)mit dem Skalarprodukt

(·, ·)X×Y :(X × Y )× (X × Y ) → IR,((x1, y1), (x2, y2)

)X×Y

= (x1, x2)X + (y1, y2)Y

ebenfalls ein Hilbert-Raum.

Die Normeigenschaften von ‖ · ‖X×Y ergeben sich durch einfaches Einsetzenaus den Eigenschaften der Normen ‖ · ‖X und ‖ · ‖Y . Gleiches gilt fur dieEigenschaften des Skalarproduktes (·, ·)X×Y .

Offensichtlich kann die Definition des Produktraumes induktiv auf n Kompo-nenten erweitert werden: Sind (Xi, ‖ · ‖Xi

), i = 1, . . . , n, normierte Raume, soist (Z, ‖ · ‖Z) mit

Z = X1 × · · · ×Xn, ‖(x1, . . . , xn)‖Z =n∑

i=1

‖xi‖Xi

ebenfalls ein normierter Raum. Ist(Xi, (·, ·)Xi

)Hilbert-Raum, i = 1, . . . , n, so

ist auch(Z, (·, ·)Z

)ein Hilbertraum mit dem Skalarprodukt

((x1, . . . , xn), (x1, . . . , xn)

)Z

=n∑

i=1

(xi, xi)Xi.

Definition A.2 Sei Ω ⊂ IRm offen, m ≥ 1 und k ≥ 0.

(i) Ck(Ω) bezeichnet die Menge aller Abbildungen f : Ω → IR, die stetigepartielle Ableitungen bis zur Ordnung k (einschließlich) besitzen.

(ii) Ck(Ω) umfasst alle Abbildungen f ∈ Ck(Ω), fur die sich alle partiellenAbleitungen bis zur Ordnung k (einschließlich) stetig auf ganz Ω fort-setzen lassen.

(iii) C∞(Ω) =⋂∞

l=0 C l(Ω),

C∞c (Ω) = f ∈ C∞(Ω) | ∃K =K(f) ⊂ Ω kompakt mit f =0 auf Ω \K

63

Mit der Ableitung 0-ter Ordnung f (0) ist stets die Abbildung f selbst gemeint.

Die Bezeichnung Ck(Ω) ist nicht ganz glucklich, da Ω durch Ω nicht eindeutigbestimmt ist17. Aus dem Kontext wird jedoch stets klar, welche Menge Ωgemeint ist.

Mit der Abbildung

‖ · ‖Ck(Ω) : Ck(Ω) → IR, f 7→ ‖f‖Ck(Ω) =k∑

i=0

maxy∈Ω

|f (i)(y)|

wird(Ck(Ω), ‖ · ‖Ck(Ω)

)zu einem normierten Raum.

Definition A.3 Sei Ω ⊂ IRm offen, m ≥ 1, und f, g ∈ L1,loc(Ω). f heißtschwache Ableitung von g bzgl. des Multiindex α, wenn

∫Ω

gDαv dy = (−1)|α|∫Ω

fv dy ∀ v ∈ C∞c (Ω)

gilt. Dabei ist Dαv = ∂|α|v∂y

α11 ···∂yαm

m.

Hierbei bedeutet f ∈ L1,loc(Ω), dass f : Ω → IR messbar ist und außerdemf ∈ L1(K) fur jede prakompakte Menge K mit K ⊂ Ω gilt.

Schwache Ableitungen bzgl. α = (α1, . . . , αm) ∈ Zm, αi ≥ 0, werden ebenfallsmit dem Symbol Dαf =g bezeichnet. Fur 1 ≤ p < ∞ setzen wir

‖f‖p =(∫

Ω|f(y)|p dy

) 1p , ‖f‖k,p =

(∑0≤|α|≤k ‖Dαf‖p

p

) 1p.

Definition 41

(i) Lp(Ω) = f : Ω→IR | f messbar, ‖f‖p <∞ heißt Lebesgue-Raum und

(ii) W kp (Ω) = f ∈Lp(Ω) | Dαf ∈Lp(Ω) ∀ 0≤|α|≤k heißt Sobolev-Raum.

(Lp(Ω), ‖ · ‖p

)und

(W k

p (Ω), ‖f‖k,p

)sind normierte Raume, wenn Funktionen,

die auf Mengen vom Lebesgue-Maß null ubereinstimmen, miteinander identi-fiziert werden. Es ist W 0

p (Ω) = Lp(Ω). Mit den Skalarprodukten

(f, g)2 =∫Ω

fg dy und (f, g)k,2 =∑

0≤|α|≤k(Dαf, Dαg)2

werden(L2(Ω), (·, ·)2

)und

(W k

2 (Ω), (·, ·)k,2

)zu Hilbert-Raumen.

Fur W k2 (Ω) wird auch Hk(Ω) geschrieben.

Ist Ω ⊂ IRm ein Gebiet, so liegt der Raum C∞(Ω) fur 1 ≤ p < ∞ dicht inLp(Ω). Einen ausfuhrlichen Beweis geben Lieb/Loss in [14].

17Beispiel: fur Ω1 = (0, 1) und Ω2 = (0, 12 )∪ ( 1

2 , 1) ist zwar Ω1 6= Ω2 aber Ω1 = Ω2 = [0, 1]

64

Beispiel A.4 Sei Ω = (0, l) ⊂ IR1, l > 0. Fur f ∈ L2(Ω) definieren wir

g : Ω → IR, y 7→ g(y) =∫ y

0f(τ) dτ.

Mit der Holder-Ungleichung folgt

‖g‖22 =

∫ l

0

∣∣ ∫ y

0f(τ) dτ

∣∣2 dy ≤ ∫ l

0‖g‖2

1 dy ≤ ∫ l

0‖1‖2

2‖g‖22 dy = l2‖f‖2

2,

so dass auch g ein Element aus L2(Ω) ist. Wir zeigen, dass f (schwache) Ab-leitung von g ist. Fur beliebiges v ∈ C∞

c (Ω) ist

∫ l

0f(τ) v(τ) dτ =

∫ l

0f(τ)

∫ τ

0∂v∂y

(ξ) dξ dτ =∫

Gf(τ)∂v

∂y(ξ) dξ dτ.

Hier wird uber die Menge

G = (τ, ξ) | 0 ≤ τ ≤ l, 0 ≤ ξ ≤ τ = (τ, ξ) | 0 ≤ ξ ≤ l, ξ ≤ τ ≤ l integriert. Durch Umparametrisieren folgt also

∫ l

0f(τ) v(τ) dτ =

∫ l

0

∫ l

ξf(τ)∂v

∂y(ξ) dτ dξ =

∫ l

0

(∫ l

0f(τ) dτ−∫ ξ

0f(τ) dτ

)∂v∂y

(ξ) dξ

= g(l)(v(l)− v(0)

)− ∫ l

0g(ξ) ∂v

∂y(ξ) dξ = − ∫ l

0g(ξ) ∂v

∂y(ξ) dξ,

denn wegen v ∈ C∞c (Ω) ist v(0) = v(l) = 0. Damit ist ∂g

∂y= f ∈ L2(Ω), also

g ∈ H1(Ω).

Die Sobolevschen Einbettungssatze machen Aussagen uber die stetige Einbet-tung der Raume W k

p (Ω) in Raume stetig-differenzierbarer Funktionen. Daruberhinaus existiert ein Spuroperator γ, der verallgemeinerte Randwerte bestimmt.Satz A.5 fasst diese Eigenschaften fur p = 2 zusammen [32].

Satz A.5 Sei Ω ⊂ IRm ein beschranktes Gebiet mit Lipschitzstetigem Rand.Dann gilt:

(i) Fur k − j > m2

ist W k2 (Ω) stetig in Cj(Ω) eingebettet, W k

2 (Ω) → Cj(Ω),

d.h. zu jedem Element f ∈ W k2 (Ω) existiert genau ein f ∈ Cj(Ω), so

dass f und f sich nur auf einer Menge vom Maß null unterscheiden.Daruberhinaus ist ‖f‖Cj(Ω) ≤ c ‖f‖k,2 mit einer von f unabhangigenKonstanten c.

(ii) Es existiert genau ein linearer, stetiger Operator

γ : W 12 (Ω) → L2(∂Ω),

so dass γ(f) : ∂Ω → IR fur jedes f ∈ C1(Ω) die Randwerte von fbeschreibt. γ(f) heißt Spur von f ∈ W 1

2 (Ω) und γ wird als Spuroperatorbezeichnet.

Ist f ∈ W 12 (Ω) mit γ(f) = 0, so wird auch

”f = 0 auf ∂Ω“ geschrieben.

65

Definition A.6 Der Abschluss von C∞c (Ω) in W k

p (Ω) wird mit

W kp (Ω) be-

zeichnet. Fur den Raum

W k2 (Ω) wird das Symbol Hk

0 (Ω) verwendet.

Ist f ∈ Hk0 (Ω), k ≥ 1, so existiert eine Folge (fn) in C∞

c (Ω) mit fn → f . Wegenγ(fn) = 0 ∀ n folgt auch γ(f) = 0, d.h. f = 0 auf ∂Ω. Fur alle |α| ≤ k− 1 istsogar Dαf = 0 auf ∂Ω, denn Dαf ∈ H1

0 (Ω).

Ausfuhrlichere Informationen zu Sobolevraumen sowie Beweise des Satzes A.5konnen z.B. in [32] nachgelesen werden.

Definition A.7 Sei (X, ‖ · ‖) ein normierter Raum und I = [t0, T ] ⊂ IR einbeschranktes, abgeschlossenes Intervall. Fur f : I → X schreiben wir

(i) f ∈ C(I, X

), wenn gilt:

∀ t ∈ I ∀ ε > 0 ∃ δ > 0 : |t− t| ≤ δ ⇒ ‖f(t)− f(t)‖ ≤ ε,

(ii) f ∈ C1(I, X

), wenn gilt:

∀ t ∈ I existiert der Grenzwert ∂f∂t

(t) = limh→0

1

h

[f(t+h)−f(t)

]in X und

∂f∂t

ist ein Element von C(I, X

).

Fur t = t0 und t = T werden lediglich einseitige Grenzwerte betrachtet.

Beispiel A.8 Sei Ω = (0, l) ⊂ IR1, l > 0 und f ∈ C1(I, L2(Ω)

), d.h. fur jedes

t ∈ I ist f(·, t) ein Element von L2(Ω). Wir setzen

g : I → L2(Ω), t 7→ g(·, t) mit g(y, t) =∫ y

0f(τ, t) dτ.

Beispiel A.4 zeigt, dass sogar g(·, t) ∈ H1(Ω) fur alle t ∈ I gilt. Wir zeigeng ∈ C1

(I, H1(Ω)

).

Sei (hn) eine beliebige Folge reeller Zahlen, hn 6= 0 ∀ n, mit hn → 0. Sei nunt ∈ I und ε > 0 beliebig. Wegen f ∈ C1

(I, L2(Ω)

)existiert ein n0 ∈ IN mit

∥∥∥ 1

hn

[f(·, t + hn)− f(·, t)]− ∂f

∂t(·, t)

∥∥∥2≤ ε√

l2 + 1∀ n ≥ n0.

Aus f ∈ C1(I, L2(Ω)

)und dem Beispiel A.4 folgt auch, dass fur

g : I → H1(Ω), t 7→ g(·, t) mit g(y, t) =∫ y

0∂f∂t

(τ, t) dτ

g ∈ C(I, H1(Ω)

)gilt, denn ∂f

∂t∈ C

(I, L2(Ω)

). Fur alle n ≥ n0 ist dann

∥∥∥ 1

hn

[g(·, t+hn)−g(·, t)]−g(·, t)∥∥∥

2

1,2

=∥∥∥ 1

hn

[g(·, t+hn)−g(·, t)]−g(·, t)∥∥∥

2

2+

∥∥∥ 1

hn

[f(·, t+hn)−f(·, t)]−∂f

∂t(·, t)

∥∥∥2

2

≤ ∫ l

0

∥∥ 1

hn

[f(·, t+hn)−f(·, t)]−∂f

∂t(·, t)

∥∥2

1dy+

ε2

l2+1≤ l2

ε2

l2+1+

ε2

l2+1= ε2,

so dass 1hn

[g(·, t+hn)−g(·, t)] → g(·, t) in H1(Ω), also g ∈ C1(I, H1(Ω)

).

66

Lemma A.9 Fur jedes ψ ∈ C1(I, H2(Ω)

)gilt ∂

∂t∂∂y

ψ = ∂∂y

∂∂t

ψ.

Beweis: Sei t ∈ I beliebig fixiert. Wir zeigen ∂∂t

∂∂y

ψ(·, t) = ∂∂y

∂∂t

ψ(·, t).Sei dazu (hn)n eine beliebige Nullfolge in IR, hn 6= 0 ∀n.

(i) Wegen ψ ∈ C1(I, H2(Ω)

)ist

ϕt,n =1

hn

[ ∂

∂yψ(·, t + hn)

∂

∂yψ(·, t)]

fur jedes n ∈ IN die schwache Ableitung von

ψt,n =1

hn

[ψ(·, t + hn)− ψ(·, t)].

(ii) Aus ψ ∈ C1(I, H2(Ω)

)folgt auch ∂ψ

∂y∈ C1

(I, H1(Ω)

), denn wenn

(ψt,n

)n

Cauchyfolge in H2(Ω) ist, so ist wegen der Definition der H2-Norm(ϕt,n

)n

Cauchyfolge in H1(Ω). Damit gilt ϕt,n → ∂∂t

∂ψ∂y

(·, t) in H1(Ω).

(iii) Schließlich bedeutet ψ ∈ C1(I, H2(Ω)

)auch ψt,n → ∂ψ

∂t(·, t) in H2(Ω).

Die Eigenschaften (i)-(iii) sichern, dass ∂∂t

∂ψ∂y

(·, t) die schwache Ableitung von∂ψ∂t

(·, t) ist, also ∂∂t

∂∂y

ψ = ∂∂y

∂∂t

ψ. ¤

67

B Grundlagen aus der Linearen Algebra

Wir betrachten den linearen Raum (IRm, | · |2) mit dem Skalarprodukt

(x, y) =m∑

k=1

xkyk und |x|2 =√

(x, x).

Satz B.1 (Orthogonales Komplement) Sei K ⊂ IRm ein linearer Teil-raum mit orthogonalem Komplement K⊥ = x ∈ IRm | (x, y) = 0 ∀ y ∈ K .Dann ist IRm = K ⊕K⊥, also m = dim IRm = dim K + dim K⊥.

Beweis: x ∈ K ∩ K⊥ ⇒ 0 = (x, x) = |x|2 ⇒ x = 0 Damit istK ∩K⊥ = 0. Es ist noch IRm = K + K⊥ zu zeigen. Sei dazu P : IRm → Kdie orthogonale Projektion auf K. Fur beliebiges x ∈ IRm sei w = x− Px. DaP orthogonale Projektion ist, folgt (w, y) = 0 ∀ y ∈ K, d.h. w ∈ K⊥. Damitwurde gezeigt: x = Px + w ∈ K + K⊥. ¤

Satz B.2 (Hauptsatz uber lineare Abbildungen) Sei A : IRm → IRn ei-ne lineare Abbildung. Dann ist m = dim IRm = dim ker A + dim im A.

Beweis: z.B. in [7] ¤

Satz B.3 Sei A : IRm → IRn eine lineare Abbildung. Dann gilt

(i) ker A = (im AT )⊥ (ii) im A = (ker AT )⊥.

Beweis: Beachte AT : IRn → IRm. Zunachst zeigen wir”⊂“.

x ∈ kerA, y=AT z ∈ imAT bel.⇒ (x, y)=xT AT z=(Ax)T z=0 ⇒ x ∈ (imAT )⊥.

y=Ax∈ imA, z∈kerAT bel. ⇒ (y, z)=(Ax)T z=xT (AT z)=0 ⇒ y ∈ (kerAT )⊥.

Die Inklusionen”⊃“ folgen nun aus Dimensionsgrunden:

dim(im AT )⊥ ≥ dim ker A =Satz B.2

m− dim im A ≥ m− dim(ker AT )⊥

=Satz B.1

m−(n−dim ker AT ) =Satz B.2

m−(n−(n−dim im AT ))

= m− dim im AT =Satz B.1

dim(im AT )⊥

Wegen ker A ⊂ (im AT )⊥ und im A ⊂ (ker AT )⊥ aber dim ker A = dim(im AT )⊥

und dim im A = dim(ker AT )⊥ folgt die Behauptung. ¤

Satz B.4 (Spaltenrang=Zeilenrang) Sei A : IRm → IRn linear. Dann ist

rang(A) = rang(AT ).

68

Beweis: rang(A) = dim im A =Satz B.3

dim(ker AT )⊥ =Satz B.1

n−dim ker AT

=Satz B.2

dim im AT = rang(AT ) ¤

Definition B.5 Die Matrix M ∈ Lb(IRm) heißt positiv definit, wenn

∀ x ∈ IRm, x 6= 0 ⇒ xT Mx > 0.

Im Allgemeinen ist M nicht symmetrisch.

Lemma B.6 M ∈ Lb(IRm) sei positiv definit. Dann gilt:

(i) M ist regular (ii) M−1 ist selbst positiv definit.

Beweis:(i) x ∈ ker M ⇒ Mx = 0 ⇒ xT Mx = 0 ⇒ x = 0, d.h. ker M = 0(ii) Sei x ∈ IRm, x 6= 0. Dann ist y = M−1x 6= 0, denn M−1 ist regular. ⇒xT M−1x = (xT M−1x)T = xT (M−1)T x = xT (M−1)T MM−1x = yT My > 0. ¤

Schließlich ist noch das Lemma 10 aus [28] zu beweisen, von dem schon so oftGebrauch gemacht wurde.

Lemma 10 M ∈ Lb(IRm) sei positiv definit. Fur N ∈ Lb(IR

m, IRn) gilt dann

(i) ker NMNT = ker NT (ii) im NMNT = im N.

Beweis: (i)”⊂“ x ∈ ker NMNT ⇒ 0 = xT NMNT x = (NT x)T M(NT x)

⇒ NT x = 0 ⇒ x ∈ ker NT ;”⊃“ ist trivial

(ii)”⊂“ ist trivial. Mit dieser Inklusion folgt aber dim im N ≥ dim im NMNT

=Satz B.2

n − dim ker NMNT =(i)

n − dim ker NT = rang(NT ) =Satz B.4

rang(N) =

dim im N , so dass dim im N = dim im NMNT und mit”⊂“ auch (ii) gilt. ¤

69

C Verzeichnis der verwendeten Symbole

IN Menge der naturlichen ZahlenIRk linearer Raum der k-dimensionalen Spaltenvektoren, k ≥ 1

1 1 =( 1...

1

)∈ IRk, k ∈ IN

(a, b) ⊂ IR offenes Intervall: (a, b) = y ∈ IR | a < y < b[a, b] ⊂ IR abgeschlossenes Intervall: [a, b] = y ∈ IR | a ≤ y ≤ b∃

”es existiert“

∀”fur alle“

∈ x ∈ A ist ein Element der Menge A⊂ Teilmenge: A ⊂ B ⇔ ( ∀ x ∈ A ⇒ x ∈ B )

Insbesondere ist A = B zugelassen.\ Differenzmenge: x ∈ A \B ⇔ ( x ∈ A und x 6∈ B )∪ Vereinigung: x ∈ A ∪B ⇔ (x ∈ A oder x ∈ B)MIR(n,m) Menge der reellen n×m MatrizenLb(IR

m, IRn) Raum der linearen Abbildungen IRm → IRn

Mk k-faches kartesisches Produkt M × · · · ×M der Menge MM Abschluss der Menge M (in der entsprechenden Topologie)

δij Kronecker-Symbol mit δij =

1, i = j

0, i 6= j

Ik Einheitsmatrix Ik = (δij) ∈ MIR(k, k) der Dimension kAT Transponierte Matrix: A = (aij) ⇒ AT = (aji)ker A Kern von A: x ∈ ker A ⇔ Ax = 0im A Bild von A: y ∈ im A ⇔ ∃ x : Ax = yrang(A) Rang der Matrix A: rang(A) = dim im Adim U Dimension des Unterraumes Uspan(M) lineare Hulle der Menge M⊥ U⊥ ist das orthogonale Komplement des Unterraumes U⊕ direkte Summe von UnterraumenC∞

c (Ω)Lp(Ω)Hk(Ω)Hk

0 (Ω)

Funktionenraume, vgl. Anhang A

I Identitat: I(x) = x ∀ xDA Definitionsbereich des Operators A| · | Absolutbetrag im IR1.

| · |2 Euklidische Norm im IRk: |y|2 =(∑k

i=1 y2i

) 12, y ∈ IRk

‖ · ‖p , ‖ · ‖k,p Normen auf Funktionenraumen, vgl. Anhang A⊗ komponentenweise Multiplikation

⊗ : IRkH × IRkH → IRkH , (a, b) 7→ a⊗ b =

(a1b1...

akHbkH

)

70

∇, grad Nabla-Operator, Gradient: grad =(

∂∂y1

, . . . , ∂∂yk

)

∆ Laplace-Operator: ∆ =∑k

i=0∂2

∂y2i

Ω Simulationsgebiet der”Halbleiter“

ψ elektrostatisches Potential, Losung der”Halbleiter“-PDE

ϕ Element von [H2(Ω) ∩H10 (Ω)]kH

A, A Koeffizient der Netzwerk-DAE bzw. des ADASAE Inzidenzmatrix des Teilgraphen der Elemente EB, B Koeffizient der Netzwerk-DAE bzw. des ADASC KapazitatsmatrixD, D Koeffizient der Netzwerk-DAE bzw. des ADAS

E =(

AH00

)erweiterte Inzidenzmatrix der

”Halbleiterelemente“

f, F rechte Seite der PDE und deren OrtsintegralG, Gi, Gi Leitwertmatrix, Matrixkette der DAE bzw. Operatorkette

des ADASH, Hk(Ω) allg. Hilbert-Raum, Sobolev-RaumI IdentitatI beschranktes, abgeschlossenes Intervall in IR1

K Unterraum von H2(Ω) ∩H10 (Ω) zur Bestimmung von kerG1

L, Lp(Ω) Induktivitatsmatrix, Lebesgue-RaumMH , MR Matrizen MH = QT

V−CQTCAH und MR = QT

V−CQTCAR

Ni, Ni Kern von Gi bzw. Gi

O1, O2 Ohmsche KontaktePi, Pi Projektoren langs Ni bzw. Ni, Pi = Im −Qi, Pi = I − Qi

Qi, Qi Projektoren auf Ni bzw. Ni

QC , QV−C Projektoren auf ker ATC bzw. auf ker AT

V QC

r , R, R Restriktionsoperator, Projektor auf im D langs ker Abzw. auf imD langs kerA

R Projektor auf ker MH langs (ker MH)⊥

S, T Projektoren aus Lemma 32t Zeitargument t ∈ I

u(t) Losung des ADA-Systemsv Vektor in IRn−1

w, Wi Vektor in IRkH und Projektor langs imGi

x(t), X Losung der Netzwerk-DAEund X = IRm × [H2(Ω) ∩H1

0 (Ω)]kH × IRkH

y, Y Ortsvariable y ∈ Ω und Y = IRm × [L2(Ω)]kH × IRkH

z, Z Vektor in IRk und Z = IRn+kL × IRkH

71

Literatur

[1] Alı, G., Bartel, A., Gunther, M., Tischendorf, C.: An existence analysisfor elliptic partial differential-algebraic multiphysics models in electricalnetwork design, Preprint, Universitat Karlsruhe (TH) (2002)

[2] Brenan, K.E., Campbell, S.L., Petzold, L.R.: Numerical solution of initial-value problems in differential-algebraic equations, North-Holland, NewYork (1989)

[3] Engl, W., Laur, R., Dirks, H.: MEDUSA - A Simulator for Modular Cir-cuits, IEEE Transactions on Computer-Aided Design of Integrated Cir-cuits and Systems. Vol. CAD-1, No. 2., April 1982.

[4] Estvez Schwarz, D., Tischendorf, C.: Structural analysis of electric circuitsand consequences for MNA, International Journal of Circuit Theory andApplications, 28, 131-162 (2000)

[5] Feldmann, U., Gunther, M.: CAD-based electric-circuit modeling in indu-stry – I. Mathematical structure and index of network equations, Surveyon Mathematics for Industry, Springer (1999)

[6] Feldmann, U., Gunther, M.: The DAE-index in electric circuit simulation,Technische Universitat Munchen, TUM-M9319, (1993)

[7] Fischer, G.: Lineare Algebra, Vieweg Studium, Grundkurs Mathematik,(1997)

[8] Griepentrog, E., Marz, R.: Differential-algebraic equations and their nu-merical treatment, Teubner, Leipzig (1986)

[9] Gunther, M.: Partielle differential-algebraische Systeme in der numeri-schen Zeitbereichsanalyse elektrischer Schaltungen, Fortschritt-BerichteVDI, Reihe 20, Nr. 343, (2001)

[10] Gunther, M., Rentrop, P.: PDAE-Netzwerkmodelle in der elektrischenSchaltungssimulation, In: Analog ’99: 5. ITG/GMM-Diskussionssitzung

”Entwicklung von Analogschaltungen mit CAE-Methoden mit dem

Schwerpunkt Entwurfsmethodik und parasitare Effekte“, Hrsg.: W. John,Frankfurt (2001), 31-38

[11] Hairer, E., Wanner, G.: Solving ordinary differential equations II: stiffand differential algebraic problems., Springer, Berlin Heidelberg New YorkTokyo (1991)

[12] Kuh, E.S., Desoer, C.A.: Basic Circuit Theory, McGraw-Hill Book Com-pany (1969)

72

[13] Lamour, R., Marz, R., Tischendorf, C.: PDAEs and Furter Mixed Sys-tems as Abstract Differential Algebraic Systems, Preprint 11, Humboldt-Universitat Berlin, (2001)

[14] Lieb, E.H., Loss, M.: Analysis, Graduate Studies in Mathematics, Volume14, American Mathematical Society (1996)

[15] Lucht, W., Strehmel, K., Eichler-Liebenow, C.: Indexes and Special Dis-cretization Methods for Linear Partial Differential Algebraic Equations,BIT, Vol. 39, No. 3, pp. 484-512 (1999)

[16] Markowich, P.A.: The Stationary Semiconductor Device Equations, Sprin-ger, Wien New York (1986)

[17] Markowich, P.A., Ringhofer, C.A., Schmeiser, C.: Semiconductor Equati-ons, Springer, Wien New York (1990)

[18] Marz, R.: Numerical methods for differential algebraic equations, ActaNumer. I: 141-198 (1992)

[19] Marz, R.: The index of linear differential algebraic equations with properlystated leading term, Preprint 7, Humboldt-Universitat (2001), erscheint in

”Results in Mathematics“

[20] Mayaram, K.: CODECS: A Mixed-Level Circuit and Device Simulator,Memorandum No. UCB/ERL M88/71, Electronics Research Laboratory,College of Engineering, University of California, Berkeley (1988)

[21] Montrone, Francesco: Ein robustes adaptives Verfahren zur numerischenLosung partieller Differentialgleichungen bei elektronischen Bauelementenin drei Raumdimensionen, Dissertation an der TU Munchen (1995)

[22] Rienen, U.v., Gunther, M., Hecht, D. (Eds.): Scientific Computing in Elec-trical Engineering, Lecture Notes in Computational Science and Enginee-ring 18, Springer, Berlin Heidelberg New York (2001)

[23] Rotella, F.M.: Mixed Circuit and Device Simulation for Analysis, Design,and Optimization of Opto-Electronic, Radio Frequency, and High SpeedSemiconductor Devices, Dissertation, Department of Electrical Enginee-ring, Stanford University, (2000)

[24] Rotella, F.M., Yu, Z., Dutton, R.W.: PISCES-2ET and Its ApplicationSubsystems, PART III: Mixed-Mode Device/Circuit Simulator Integra-ted Circuits Laboratory Stanford, University Stanford, California 94305(1994)

[25] Schumilina, I.: Index-3-Algebro-Differentialgleichungen mit proper formu-liertem Hauptterm, Preprint 20, Humboldt-Universitat Berlin, (2001)

73

[26] Selberherr, S.: Analysis and Simulation of Semiconductor Devices,Springer, Wien New York, (1984)

[27] Strehmel, K., Weiner, R.: Numerik gewohnlicher Diffrentialgleichungen,Teubner Verlag (1995)

[28] Tischendorf, C.: Topological index calculation of differential-algebraicequations in circuit simulation, Surveys on Mathematics for Industry 8,Springer (1999)

[29] Tischendorf, C.: Numerische Simulation elektrischer Netzwerke, Vorle-sungsskript, Humboldt-Universitat zu Berlin, Sommersemester 2000

[30] Van Roosbroeck, W.V.: Theory of Flow of Electrons and Holes in Ger-manium and Other Semiconductors, Bell System Technical Journal 29, S.560-607, (1950)

[31] Vlach, J., Singhal, K.: Computer methods for circuit analysis and design,Van Nostrand Reinhold, New York (1983)

[32] Zeidler, E.: Nonlinear Functional Analysis and its Applications II/A –Linear Monotone Operators, Springer, New York Berlin Heidelberg (1990)

[33] Zielke, G.: Auflosung beliebiger linearer algebraischer Gleichungen durchBlockzerlegung, Beitrage zur Numerischen Mathematik 8 (1980), 181-199

74

Erklarung

Ich erklare, dass ich die vorliegende Arbeit selbststandigund nur unter Verwendung der angegebenen Literaturund Hilfsmittel angefertigt habe.

Berlin, 30. September 2002