diff options
| -rw-r--r-- | Physikunit.pas | 108 | ||||
| -rw-r--r-- | Plasmapropagation.lps | 103 |
2 files changed, 172 insertions, 39 deletions
diff --git a/Physikunit.pas b/Physikunit.pas index 5d8881a..27cca80 100644 --- a/Physikunit.pas +++ b/Physikunit.pas @@ -91,12 +91,14 @@ type tImpulsPunkt = class(tObject) // Besetzungsdichten verschiedener Teilchenspezies für ein Paar (r;p) private - raumpunkt: tRaumPunkt; - nachbarn: array[0..1,boolean] of tImpulsPunkt; // x- und px-Nachbarn (jeweil in Richtung - und +) + raumpunkt: tRaumPunkt; + nachbarn: array[0..1,boolean] of tImpulsPunkt; // x- und px-Nachbarn (jeweil in Richtung - und +) + istNullPkt: boolean; public werte: array of array[boolean] of extended; // Besetzungsdichten und deren Ableitungen für verschiedene Teilchenspezies iMGamma: array of extended; // 1/m/gamma (= v/p) - constructor create(rp: tRaumPunkt; linksRaum,linksImpuls: tImpulsPunkt; anzTeilchen: longint); + constructor create(anzTeilchen: longint); overload; // NullPunkt erstellen + constructor create(rp: tRaumPunkt; linksRaum,linksImpuls: tImpulsPunkt; anzTeilchen: longint); overload; destructor destroy; override; {$DEFINE LiKotImpulsPunktHeader} {$INCLUDE linearkombinationen.inc} @@ -106,6 +108,9 @@ type function gibMirWas(var defizit: extended; entfernung, teilchen,richtung: longint; positiv: boolean): boolean; inline; // ich habe jemanden erreicht procedure setzeNull; inline; procedure berechneAbleitungen(pX,iDX,iDPX: extended); inline; + function gradient(rtg,tlc: longint; iD,fak: extended): extended; inline; + property istNullPunkt: boolean read istNullPkt; + function gibNullPunkt: tImpulsPunkt; end; { tRaumPunkt } @@ -580,12 +585,37 @@ end; // tImpulsPunkt **************************************************************** -constructor tImpulsPunkt.create(rp: tRaumPunkt; linksRaum,linksImpuls: tImpulsPunkt; anzTeilchen: longint); +constructor tImpulsPunkt.create(anzTeilchen: longint); var i: longint; b: boolean; begin inherited create; + istNullPkt:=true; // als NullPunkt erstellen + raumPunkt:=nil; + fillchar(werte,sizeof(werte),#0); + setlength(werte,anzTeilchen); + for i:=0 to length(werte)-1 do + for b:=false to true do + werte[i,b]:=0; + fillchar(iMGamma,sizeof(iMGamma),#0); + setlength(iMGamma,length(werte)); + for i:=0 to length(iMGamma)-1 do + iMGamma[i]:=1; + nachbarn[0,false]:=self; + nachbarn[0,true]:=self; + nachbarn[1,false]:=self; + nachbarn[1,true]:=self; +end; + +constructor tImpulsPunkt.create(rp: tRaumPunkt; linksRaum,linksImpuls: tImpulsPunkt; anzTeilchen: longint); +var + i: longint; + b: boolean; + np: tImpulsPunkt; +begin + inherited create; + istNullPkt:=false; raumPunkt:=rp; fillchar(werte,sizeof(werte),#0); setlength(werte,anzTeilchen); @@ -601,8 +631,21 @@ begin nachbarn[1,false]:=linksImpuls; nachbarn[1,true]:=nil; for i:=0 to 1 do - if assigned(nachbarn[i,false]) then + if assigned(nachbarn[i,false]) and + not nachbarn[i,false].istNullPunkt then nachbarn[i,false].nachbarn[i,true]:=self; + np:=nil; + for i:=0 to 1 do + if assigned(nachbarn[i,false]) and + not assigned(np) then + np:=nachbarn[i,false].gibNullPunkt; + if not assigned(np) then + np:=tImpulsPunkt.create(anzTeilchen); + for i:=0 to 1 do + for b:=false to true do + if not assigned(nachbarn[i,b]) then begin + nachbarn[i,b]:=np; + end; end; destructor tImpulsPunkt.destroy; @@ -612,7 +655,8 @@ var begin for i:=0 to 1 do for b:=false to true do - if assigned(nachbarn[i,b]) then + if assigned(nachbarn[i,b]) and + not nachbarn[i,b].istNullPunkt then nachbarn[i,b].nachbarn[i,not b]:=nil; setlength(werte,0); setlength(iMGamma,0); @@ -726,7 +770,9 @@ begin ( raumpunkt.emFelder[efEX,false] + raumpunkt.matFelder[i,mfPY,false] * iMGamma[i] * raumpunkt.emFelder[efBZ,false]); werte[i,true]:= - tX * ( + gradient(0,i,iDX,tX) + + gradient(1,i,iDPX,tP); +(* tX * ( (nachbarn[0,true].werte[i,false] - nachbarn[0,false].werte[i,false]) + sign(tX)*raumpunkt.felder.teilchen[i].diffusion[0] * ( nachbarn[0,true].werte[i,false] + @@ -741,7 +787,7 @@ begin nachbarn[1,false].werte[i,false] - 2*werte[i,false] ) - )* iDPX/2; + )* iDPX/2; *) end; // Der jeweils zweite Summand entspricht (in etwa) einer thermischen @@ -749,6 +795,50 @@ begin // d(ln n)/d[xp] * d[XP] = "diffusion" gewährleisten. end; +function tImpulsPunkt.gradient(rtg,tlc: longint; iD,fak: extended): extended; +begin + + { + // Variante mit Diffusionsterm + result:= + fak * ((nachbarn[rtg,true].werte[tlc,false] - nachbarn[rtg,false].werte[tlc,false]) + + sign(fak)*raumpunkt.felder.teilchen[tlc].diffusion[rtg] * ( + nachbarn[rtg,true].werte[tlc,false] + + nachbarn[rtg,false].werte[tlc,false] - + 2*werte[tlc,false] + ) + ) * iD/2; + } + + result:= + fak * ( - 1/60 * nachbarn[rtg,false].nachbarn[rtg,false].nachbarn[rtg,false].werte[tlc,false] + + 3/20 * nachbarn[rtg,false].nachbarn[rtg,false].werte[tlc,false] + - 3/4 * nachbarn[rtg,false].werte[tlc,false] + + 3/4 * nachbarn[rtg,true].werte[tlc,false] + - 3/20 * nachbarn[rtg,true].nachbarn[rtg,true].werte[tlc,false] + + 1/60 * nachbarn[rtg,true].nachbarn[rtg,true].nachbarn[rtg,true].werte[tlc,false] + ) * iD; + + +// probier mal den Gradienten anders zu berechnen und zwar so: +// https://en.wikipedia.org/wiki/Finite_difference_coefficient +end; + +function tImpulsPunkt.gibNullPunkt: tImpulsPunkt; +var + i: longint; +begin + if istNullPunkt then begin + result:=self; + exit; + end; + result:=nil; + for i:=0 to 1 do + if assigned(nachbarn[i,false]) and + not assigned(result) then + result:=nachbarn[i,false].gibNullPunkt; +end; + // tRaumPunkt ****************************************************************** constructor tRaumPunkt.create(linkerNachbar: tRaumPunkt; derBesitzer: tFelder; teilchen,_aP: longint; _dP: extended); @@ -871,7 +961,7 @@ var i: longint; begin // df/dt = - v Nabla f - q(E + v/c x B) Nabla_p f - for i:=1 to aP-2 do // Rand wird ausgelassen! + for i:=0 to aP-1 do phasenraum[i].berechneAbleitungen((i-(aP div 2))*dP,iDX,1/dP); end; diff --git a/Plasmapropagation.lps b/Plasmapropagation.lps index 9a7b886..9330883 100644 --- a/Plasmapropagation.lps +++ b/Plasmapropagation.lps @@ -15,13 +15,14 @@ <Unit1> <Filename Value="Physikunit.pas"/> <IsPartOfProject Value="True"/> + <IsVisibleTab Value="True"/> <EditorIndex Value="1"/> - <TopLine Value="1577"/> - <CursorPos Y="714"/> - <FoldState Value=" T3mB0{F5 pjdmK0G[945jZ0N4]90jX0A2 piYj90C7]9ajL0v5]9fkW0E[j4Ejg082]93jH03013 T0%LPc071F311#"/> - <UsageCount Value="176"/> + <TopLine Value="611"/> + <CursorPos X="77" Y="819"/> + <FoldState Value=" T3mG0{F5 pjdmK0G[945jZ0L3]93jW0C31117 piXnD0F3]9bjH0J4[I4jR0Q6 pkMlL03013 T0)PPc071F]sQ"/> + <UsageCount Value="179"/> <Bookmarks Count="1"> - <Item0 Y="1247"/> + <Item0 Y="1337"/> </Bookmarks> <Loaded Value="True"/> </Unit1> @@ -30,16 +31,15 @@ <IsPartOfProject Value="True"/> <EditorIndex Value="-1"/> <CursorPos X="20" Y="8"/> - <UsageCount Value="138"/> + <UsageCount Value="141"/> </Unit2> <Unit3> <Filename Value="input.plap"/> <IsPartOfProject Value="True"/> - <IsVisibleTab Value="True"/> <EditorIndex Value="2"/> <TopLine Value="30"/> <CursorPos X="23" Y="62"/> - <UsageCount Value="137"/> + <UsageCount Value="140"/> <Loaded Value="True"/> <DefaultSyntaxHighlighter Value="None"/> </Unit3> @@ -49,14 +49,14 @@ <EditorIndex Value="-1"/> <TopLine Value="285"/> <CursorPos X="20" Y="323"/> - <UsageCount Value="85"/> + <UsageCount Value="88"/> </Unit4> <Unit5> <Filename Value="linearkombinationen.inc"/> <IsPartOfProject Value="True"/> <EditorIndex Value="-1"/> <TopLine Value="64"/> - <UsageCount Value="48"/> + <UsageCount Value="51"/> </Unit5> <Unit6> <Filename Value="input.epost"/> @@ -155,83 +155,126 @@ <UsageCount Value="10"/> </Unit19> </Units> - <JumpHistory Count="19" HistoryIndex="18"> + <JumpHistory Count="30" HistoryIndex="29"> <Position1> <Filename Value="Physikunit.pas"/> - <Caret Line="1683" Column="12" TopLine="1135"/> + <Caret Line="1229" Column="51" TopLine="1217"/> </Position1> <Position2> <Filename Value="Physikunit.pas"/> - <Caret Line="1685" TopLine="1198"/> + <Caret Line="637" Column="74" TopLine="608"/> </Position2> <Position3> <Filename Value="Physikunit.pas"/> - <Caret Line="1686" TopLine="1198"/> + <Caret Line="625" Column="14" TopLine="609"/> </Position3> <Position4> <Filename Value="Physikunit.pas"/> - <Caret Line="1687" TopLine="1198"/> + <Caret Line="972" Column="8" TopLine="204"/> </Position4> <Position5> <Filename Value="Physikunit.pas"/> - <Caret Line="1223" TopLine="1162"/> + <Caret Line="139" Column="3" TopLine="121"/> </Position5> <Position6> <Filename Value="Physikunit.pas"/> - <Caret Line="1227" TopLine="1162"/> + <Caret Line="656" Column="16" TopLine="632"/> </Position6> <Position7> <Filename Value="Physikunit.pas"/> - <Caret Line="1217" TopLine="1069"/> + <Caret Line="1666" TopLine="973"/> </Position7> <Position8> <Filename Value="Physikunit.pas"/> - <Caret Line="1229" Column="51" TopLine="1217"/> + <Caret Line="959" Column="74" TopLine="926"/> </Position8> <Position9> <Filename Value="Physikunit.pas"/> - <Caret Line="637" Column="74" TopLine="608"/> + <Caret Line="956" Column="5" TopLine="936"/> </Position9> <Position10> <Filename Value="Physikunit.pas"/> - <Caret Line="625" Column="14" TopLine="609"/> + <Caret Line="734" Column="22" TopLine="718"/> </Position10> <Position11> <Filename Value="Physikunit.pas"/> - <Caret Line="972" Column="8" TopLine="204"/> + <Caret Line="768" Column="27" TopLine="729"/> </Position11> <Position12> <Filename Value="Physikunit.pas"/> - <Caret Line="139" Column="3" TopLine="121"/> + <Caret Line="15"/> </Position12> <Position13> <Filename Value="Physikunit.pas"/> - <Caret Line="656" Column="16" TopLine="632"/> + <Caret Line="731" Column="24" TopLine="703"/> </Position13> <Position14> <Filename Value="Physikunit.pas"/> - <Caret Line="1666" TopLine="973"/> + <Caret Line="730" Column="28" TopLine="703"/> </Position14> <Position15> <Filename Value="Physikunit.pas"/> - <Caret Line="622" TopLine="713"/> + <Caret Line="765" Column="66" TopLine="729"/> </Position15> <Position16> <Filename Value="Physikunit.pas"/> - <Caret Line="30" Column="61" TopLine="46"/> + <Caret Line="632" Column="42" TopLine="603"/> </Position16> <Position17> <Filename Value="Physikunit.pas"/> - <Caret Line="959" Column="74" TopLine="926"/> + <Caret Line="644" Column="26" TopLine="612"/> </Position17> <Position18> <Filename Value="Physikunit.pas"/> - <Caret Line="956" Column="5" TopLine="936"/> + <Caret Line="112" Column="60" TopLine="82"/> </Position18> <Position19> <Filename Value="Physikunit.pas"/> - <Caret Line="30" Column="64" TopLine="75"/> + <Caret Line="630" Column="30" TopLine="609"/> </Position19> + <Position20> + <Filename Value="Physikunit.pas"/> + <Caret Line="610" Column="26" TopLine="587"/> + </Position20> + <Position21> + <Filename Value="Physikunit.pas"/> + <Caret Line="101" Column="71" TopLine="81"/> + </Position21> + <Position22> + <Filename Value="Physikunit.pas"/> + <Caret Line="839" Column="60" TopLine="802"/> + </Position22> + <Position23> + <Filename Value="Physikunit.pas"/> + </Position23> + <Position24> + <Filename Value="Physikunit.pas"/> + <Caret Line="610" Column="49" TopLine="587"/> + </Position24> + <Position25> + <Filename Value="Physikunit.pas"/> + <Caret Line="101" Column="61" TopLine="75"/> + </Position25> + <Position26> + <Filename Value="Physikunit.pas"/> + <Caret Line="628" Column="21" TopLine="523"/> + </Position26> + <Position27> + <Filename Value="Physikunit.pas"/> + <Caret Line="114" TopLine="88"/> + </Position27> + <Position28> + <Filename Value="Physikunit.pas"/> + <Caret Line="643" TopLine="588"/> + </Position28> + <Position29> + <Filename Value="Physikunit.pas"/> + <Caret Line="638" Column="33" TopLine="619"/> + </Position29> + <Position30> + <Filename Value="Physikunit.pas"/> + <Caret Line="640" Column="27" TopLine="620"/> + </Position30> </JumpHistory> </ProjectSession> </CONFIG> |