An early version of the gas particle simulation program.
program Oscill;
uses
Crt, Graph;
const
Memory = 100;
Windows = 4;
Count = 49;
dt = 5e-3;
MaxStat = 120;
type
ResolutionPreference = (Lower, Higher);
ColorList = array [1..Windows] of integer;
TVector = record
x,y : Double;
end;
TParticle = record
acc,
act,
vel,
cur : TVector;
radius : Double;
underF : Boolean;
end;
var
Xmax,
Ymax,
ViewXmax,
ViewYmax : integer;
Colors: ColorList;
Ch: char;
BackColor:integer;
GraphDriver, GraphMode : integer;
MaxColors: word;
Particles: array[1..Count] of TParticle;
BoxXmax, BoxYmax, BoxXmin, BoxYmin: Double;
CurPage : Integer;
StatData : array[1..MaxStat, 1..Count] of Double;
EnableStat, StatChanged : Boolean;
CurStat : Integer;
function IntToStr(I: Real): String;
var
S: string[16];
begin
Str(I:16:8, S);
IntToStr := S;
end;
procedure ChangePage;
begin
if CurPage = 0 then
begin
CurPage := 1;
SetActivePage(1);
SetVisualPage(0);
end
else
begin
CurPage := 0;
SetActivePage(0);
SetVisualPage(1);
end;
end;
procedure FullPort;
{ Set the view port to the entire screen }
begin
SetViewPort(0, 0, Xmax, Ymax, ClipOn);
end; { FullPort }
procedure Frame;
begin
SetViewPort(0, 0, Xmax, Ymax-(TextHeight('M')+4)-1,ClipOn);
SetColor(MaxColors);
Rectangle(0, 0, Xmax-1, (Ymax-(TextHeight('M')+4)-1)-1);
SetViewPort(1, 1, Xmax-2, (Ymax-(TextHeight('M')+4)-1)-2,ClipOn);
end { Frame };
procedure MessageFrame(Msg:string);
begin
ChangePage;
FullPort;
SetColor(MaxColors);
SetTextStyle(DefaultFont, HorizDir, 1);
SetTextJustify(CenterText, TopText);
SetLineStyle(SolidLn, 0, NormWidth);
SetFillStyle(EmptyFill, 0);
Bar(0, Ymax-(TextHeight('M')+4), Xmax, Ymax);
Rectangle(0, Ymax-(TextHeight('M')+4), Xmax, Ymax);
OutTextXY(Xmax div 2, Ymax-(TextHeight('M')+2), Msg);
Frame;
ChangePage;
FullPort;
SetColor(MaxColors);
SetTextStyle(DefaultFont, HorizDir, 1);
SetTextJustify(CenterText, TopText);
SetLineStyle(SolidLn, 0, NormWidth);
SetFillStyle(EmptyFill, 0);
Bar(0, Ymax-(TextHeight('M')+4), Xmax, Ymax);
Rectangle(0, Ymax-(TextHeight('M')+4), Xmax, Ymax);
OutTextXY(Xmax div 2, Ymax-(TextHeight('M')+2), Msg);
Frame;
end;
procedure ShowStats;
var i: Integer;
Nrg: Double;
begin
Inc(CurStat);
for i:=1 to Count do
StatData[CurStat, i] := Sqrt(Particles[i].vel.x * Particles[i].vel.x + Particles[i].vel.y * Particles[i].vel.y);
MessageFrame(Concat('Current measure, # ', IntToStr(CurStat)));
end;
procedure DrawStatPlot;
const
countGaps = 30;
var maxVel: Double;
i, k, m, maxGap: Integer;
Gap: array[1..countGaps] of Integer;
begin
maxVel := 0;
for i:=1 to MaxStat do
for k:=1 to Count do
if maxVel < StatData[i, k] then maxVel := StatData[i, k];
for i:=1 to countGaps do Gap[i] := 0;
for i:=1 to MaxStat do
for k:=1 to Count do
for m:=1 to countGaps do
if (StatData[i,k] > maxVel/countGaps*(m - 1)) and (StatData[i,k] < maxVel/countGaps*(m))
then Gap[m] := Gap[m] + 1;
maxGap := 0;
for i:=1 to countGaps do
if maxGap < Gap[i] then maxGap := Gap[i];
ChangePage;
Graph.ClearViewPort;
Graph.Rectangle(Round(BoxXmin), Round(BoxYmin), Round(BoxXmax), Round(BoxYmax));
SetColor(4);
for i:=1 to countGaps do
Graph.Rectangle(
Round(BoxXmin + (BoxXmax - BoxXmin)/countGaps*(i-1)),
Round(BoxYmax),
Round(BoxXmin + (BoxXmax - BoxXmin)/countGaps*i),
Round(BoxYmax - (BoxYmax - BoxYmin)/maxGap*Gap[i]));
ChangePage;
end;
procedure WaitToGo;
var
Ch : char;
begin
MessageFrame('Press any key to continue... Esc aborts');
repeat until KeyPressed;
Ch := ReadKey;
if Ch = #27 then begin
CloseGraph;
Writeln('All done.');
Halt(1);
end
else
ClearViewPort;
MessageFrame('Press a key to stop action, Esc quits.');
end;
procedure CheckForUserInput;
begin
if KeyPressed then begin
Ch := ReadKey;
if Ch = 's' then ShowStats
else if Ch <> #27 then WaitToGo;
end;
end;
procedure DrawPoints;
var
i : Integer;
begin
ChangePage;
Graph.ClearViewPort;
SetColor(7);
Graph.Rectangle(Round(BoxXmin), Round(BoxYmin), Round(BoxXmax), Round(BoxYmax));
for i:=1 to Count do
if Particles[i].underF then
begin
SetColor(4);
SetFillStyle(SolidFill, 4);
Graph.FillEllipse(Round(Particles[i].cur.x), Round(Particles[i].cur.y), 3, 3);
end
else
begin
SetColor(9);
SetFillStyle(SolidFill, 9);
Graph.FillEllipse(Round(Particles[i].cur.x), Round(Particles[i].cur.y), 3, 3);
end;
end;
function DistFromTo(This, That : TParticle): Double;
var res: Double;
begin
res := Sqrt((This.cur.x - That.cur.x)*(This.cur.x - That.cur.x) + (This.cur.y - That.cur.y)*(This.cur.y - That.cur.y));
if res = 0 then DistFromTo := 1
else DistFromTo := res;
end;
function GetForce(This, That : TParticle): Double;
var res : Double;
begin
res := This.radius - DistFromTo(This, That);
if res > 0 then GetForce := res * 1e4
else GetForce := 0;
end;
function GetForceX(This, That : TParticle): Double;
begin
GetForceX := GetForce(This, That) * (This.cur.x - That.cur.x) / DistFromTo(This, That);
end;
function GetForceY(This, That : TParticle): Double;
begin
GetForceY := GetForce(This, That) * (This.cur.y - That.cur.y) / DistFromTo(This, That);
end;
procedure MainCycle(dt: Double);
var i,k : Integer;
wall: TParticle;
begin
for i:=1 to Count do
begin
Particles[i].acc.x := 0;
Particles[i].acc.y := 0;
end;
for i:=1 to Count do
begin
wall.radius := 3;
wall.cur.x := BoxXmin;
wall.cur.y := Particles[i].cur.y;
Particles[i].acc.x := Particles[i].acc.x + GetForceX(Particles[i], wall);
Particles[i].acc.y := Particles[i].acc.y + GetForceY(Particles[i], wall);
wall.cur.x := BoxXmax;
wall.cur.y := Particles[i].cur.y;
Particles[i].acc.x := Particles[i].acc.x + GetForceX(Particles[i], wall);
Particles[i].acc.y := Particles[i].acc.y + GetForceY(Particles[i], wall);
wall.cur.x := Particles[i].cur.x;
wall.cur.y := BoxYmin;
Particles[i].acc.x := Particles[i].acc.x + GetForceX(Particles[i], wall);
Particles[i].acc.y := Particles[i].acc.y + GetForceY(Particles[i], wall);
wall.cur.x := Particles[i].cur.x;
wall.cur.y := BoxYmax;
Particles[i].acc.x := Particles[i].acc.x + GetForceX(Particles[i], wall);
Particles[i].acc.y := Particles[i].acc.y + GetForceY(Particles[i], wall);
end;
EnableStat := True;
for i:=1 to Count do
for k:=1 to Count do
if i<>k then
begin
Particles[i].acc.x := Particles[i].acc.x + GetForceX(Particles[i], Particles[k]);
Particles[i].acc.y := Particles[i].acc.y + GetForceY(Particles[i], Particles[k]);
if (Particles[i].acc.x <> 0) or (Particles[i].acc.y <> 0) then Particles[i].underF := True
else Particles[i].underF := False;
EnableStat := EnableStat and not Particles[i].underF;
end;
StatChanged := StatChanged or not EnableStat;
for i:=1 to Count do
begin
Particles[i].vel.x := Particles[i].vel.x + (Particles[i].acc.x + Particles[i].acc.x)/2 * dt;
Particles[i].vel.y := Particles[i].vel.y + (Particles[i].acc.y + Particles[i].acc.y)/2 * dt;
Particles[i].cur.x := Particles[i].cur.x + Particles[i].vel.x * dt;
Particles[i].cur.y := Particles[i].cur.y + Particles[i].vel.y * dt;
Particles[i].act.x := Particles[i].acc.x;
Particles[i].act.y := Particles[i].acc.y;
end;
if EnableStat and StatChanged then
begin
ShowStats;
StatChanged := False;
end;
DrawPoints;
end;
procedure DoMagic;
var i : Integer;
test: TParticle;
begin
BoxXmin := ViewXmax/6*1;
BoxXmax := ViewXmax/6*5;
BoxYmin := ViewYmax/6*1;
BoxYmax := ViewYmax/6*5;
Randomize;
for i:=1 to Count do
begin
Particles[i].cur.x := (BoxXmax - BoxXmin) * Random + BoxXmin;
Particles[i].cur.y := (BoxYmax - BoxYmin) * Random + BoxYmin;
{ Particles[i].vel.x := 500 * (0.5 - Random);
Particles[i].vel.y := 500 * (0.5 - Random);}
Particles[i].vel.y := 170;
Particles[i].radius := 7;
end;
repeat
MainCycle(dt);
{ test.cur.x := Particles[1].cur.x + 1;
test.cur.y := Particles[1].cur.y + 1;
MessageFrame(IntToStr(Particles[1].acc.x));}
CheckForUserInput;
until (Ch = #27) or (CurStat = MaxStat);
if Ch<>#27 then DrawStatPlot;
WaitToGo;
end;
procedure TestGraphError(GraphErr: integer);
begin
if GraphErr <> grOk then begin
Writeln('Graphics error: ', GraphErrorMsg(GraphErr));
repeat until keypressed;
ch := readkey;
Halt(1);
end;
end;
procedure Init;
var
Err, I: integer;
StartX, StartY: integer;
Resolution: ResolutionPreference;
s: string;
begin
Resolution := Lower;
if paramcount > 0 then begin
s := paramstr(1);
if s[1] = '/' then
if upcase(s[2]) = 'H' then
Resolution := Higher;
end;
Ch := ' ';
GraphDriver := Detect;
DetectGraph(GraphDriver, GraphMode);
TestGraphError(GraphResult);
case GraphDriver of
CGA : begin
GraphDriver := CGA;
GraphMode := CGAC1;
end;
MCGA : begin
case GraphMode of
MCGAMed, MCGAHi: GraphMode := MCGAC1;
end;
end;
EGA : begin
If Resolution = Lower then
GraphMode := EGALo
else
GraphMode := EGAHi;
end;
EGA64 : begin
If Resolution = Lower then
GraphMode := EGA64Lo
else
GraphMode := EGA64Hi;
end;
PC3270 : begin
GraphDriver := CGA;
GraphMode := CGAC1;
end;
ATT400 : case GraphMode of
ATT400C1,
ATT400C2,
ATT400Med,
ATT400Hi :
begin
GraphMode := ATT400C1;
end;
end;
end;
InitGraph(GraphDriver, GraphMode, '');
TestGraphError(GraphResult);
SetGraphMode(1);
SetTextStyle(DefaultFont, HorizDir, 1);
SetTextJustify(CenterText, TopText);
MaxColors := GetMaxColor;
BackColor := 0;
Xmax := GetMaxX;
Ymax := GetMaxY;
ViewXmax := Xmax - 2;
ViewYmax := (Ymax - (TextHeight('M') + 4) - 1) - 2;
end; {init}
begin
Init;
Frame;
MessageFrame('Press a key to stop action, Esc quits.');
DoMagic;
CloseGraph;
RestoreCrtMode;
Writeln('The End.');
end.