-
#1208
A legegyszerűbben úgy lehet fényelhajlást számolni, ha két nem túl távoli pontot ugyanazon irányba mozgatunk különböző sebességgel, és az végpontokból a tangens irányt újraszámoljuk.
Az ábrán egy zölddel rajzolt ellipszis pálya van newtoni gravitációval számolva és a hozzá, és a hozzá tartozó pályamenti sebesség görbe, és pirossal a fénytörés modellel számolt pálya.
Ja, hogy csak piros van? Persze, mivel pontosan fedi a zöldet, emiatt se a zöld pálya sem a zöld sebességgörbe nem látszik.
A lassabb fény ellenére ahogy közeledünk a tömeghez, a tér gradiense egyre nagyobb lesz, emiatt nagyobb a 'gyorsulás'. Habár egyenlőre csak sebességváltozást számol a program.
![]()
A fénytöréses modell kis körpályákat számol, mintha a test csavarvonalban mindig fénysebességgel haladna.
A modell függ a fény sebességétől, amíg a newtoni nem.Ki lehet próbálni, ha T*c*0.99-el számol a program, már jelentős az eltérés a newtoni pályától.
A newtoni szimuláció egyszerüsített, de teljesen helyes is ugyan ezt a pályát adja.Ilyen skálán nem számottevő a különbség.
typedef long double float1;
float1 g,c,m, Rs,r1,v1,scale, dr, dtime,geom_unit;
float3 position, speed;
int n;
void setup()
{
g=6.67428e-11;
c=2.99792e8;
m=1.9891e30;
dr=1e-4;
Rs=(float1)2.0*m*g/(c*c);
r1=Rs*40.0;
v1=sqrtl(m*g/r1)*0.8;
dtime=2e-7;
n=400000;
scale=400.0/r1;
position=float3((float1)r1,(float1)0.0L,(float1)0.0L);
speed=float3(v1*(float1)0.3L,v1*(float1)0.7L,(float1)0.0L);
}
void newton()
{
setup();
for(int j=0;j<n;j++) {
position=position + speed*dtime;
float1 r=sqrtl(dot(position,position));
float3 T=position/r;
float1 ac=dtime*m*g/(r*r);
speed=speed - T*ac;
float3 screen=position*scale;
pixel(500+(int)(screen.x),300+(int)(screen.y),0x00ff00);
screen.x=(float1)(j/100);
screen.y=sqrtl(dot(speed,speed))/400000.0;
pixel((int)(screen.x),900-(int)(screen.y),0x00aa00);
}
}
void refraction()
{
int u=12;
setup();
for(int j=0;j<n;j++) {
position=position + speed*dtime;
float dtime2=dtime/(float1)u;
for(int i=0;i<u;i++) {
float1 alpha=(float1)i*M_PI*2.0/(float1)u;
float3 T=float3((float1)cos(alpha),(float1)sin(alpha),(float1)0.0),N;
N.x= T.y; N.y=-T.x; N.z=0.0;
float3 p1=position - N*dr + N*dr*2.0;
float3 p2=position + N*dr + N*dr*2.0;
float1 c1=c*((float1)1.0 - Rs/sqrtl(dot(p1,p1)));
float1 c2=c*((float1)1.0 - Rs/sqrtl(dot(p2,p2)));
p1=p1 + (T*c1)*dtime2;
p2=p2 + (T*c2)*dtime2;
N=p2-p1;
N=N/sqrtl(dot(N,N));
T.x= -N.y; T.y=N.x; T.z=0.0;
speed=speed + T*c;
}
float3 screen=position*scale;
pixel(500+(int)(screen.x),300+(int)(screen.y),0xff0000);
screen.x=(float1)(j/100);
screen.y=sqrtl(dot(speed,speed))/400000.0;
pixel((int)(screen.x),900-(int)(screen.y),0xaa0000);
}
}
