Informatika és tudomány

Tényleg nem létezik?
  • sublimiter
    #1344
    1.215947e-04
    2.432175e-04
    2.000232e+00

    Lol, a kapott szog ketszerese a newtoninak, pedig csak a fenysebesseg kiszamolasahoz hasznaltam az egyik Schwarzschild egyenletet, mit nemreg linkeltem.

    A progi a fenti egyenlettel szamolja a feny iranyvaltozasat, tehat fenytorest szammol.
    dpx2 = dpx *c2*c2/(c1*c1);
    Ennyi.

    Vagyis nem egeszen. Tomeggel rendelkezo testre fenyszeru belso mozgasokat kell szamolni. Ez hazi feladat. Az eredmeny tobb, mint meglepo. Nem en fogom leirni a megoldast.


    #include <stdio.h>
    #include <stdlib.h>
    #include <math.h>




    typedef long double skalar;



    struct vektor
    {
    skalar x,y,z;

    vektor() {x=0;y=0;z=0;};
    vektor(int x_,int y_,int z_) {x=x_;y=y_;z=z_;};
    vektor(skalar x_,skalar y_,skalar z_) {x=x_;y=y_;z=z_;};

    vektor operator + (vektor v) {vektor v_;v_.x=x+v.x;v_.y=y+v.y;v_.z=z+v.z;return v_;};
    vektor operator - (vektor v) {vektor v_;v_.x=x-v.x;v_.y=y-v.y;v_.z=z-v.z;return v_;};
    vektor operator * (skalar s) {vektor v_;v_.x=x*s;v_.y=y*s;v_.z=z*s;return v_;};
    vektor operator / (skalar s) {vektor v_;v_.x=x/s;v_.y=y/s;v_.z=z/s;return v_;};
    };

    skalar skalar_szorzat(vektor v1,vektor v2) { return (v1.x*v2.x + v1.y*v2.y + v1.z*v2.z);}
    skalar hossz(vektor v1) {return sqrtl(skalar_szorzat(v1,v1));}
    vektor normalt(vektor v1) { return (v1/hossz(v1));}
    vektor vektor_szorzat(vektor va,vektor vb)
    {
    vektor v_;

    v_.x=(va.y*vb.z)-(vb.y*va.z);
    v_.y=(va.z*vb.x)-(vb.z*va.x);
    v_.z=(va.x*vb.y)-(vb.x*va.y);

    return v_;
    };










    skalar g,c,m,rs,dr,rx,ry,dt;
    vektor p,dp,ddp;
    int lepes;


    void newton_gravity()
    {
    skalar r;

    while(p.x<rx)
    {
    p = p + dp*dt ;

    ddp = normalt(p);
    r = hossz(p);

    ddp = normalt(p)*(-g*m/(r*r));
    dp = dp + ddp*dt;
    p = p + ddp*dt*dt/2;
    }
    }
    void refract_gravity()
    {
    skalar r,dpx,dpy,dpx2,dpy2,c1,c2;
    vektor nx,ny;

    while(p.x<rx)
    {
    r = hossz(p);
    c1 = c*(1-g*m/(2*c*c*r)) / pow(1+g*m/(2*c*c*r),3.0);

    p = p + dp*dt ;

    r = hossz(p);
    c2 = c*(1-g*m/(2*c*c*r)) / pow(1+g*m/(2*c*c*r),3.0);//uj fenysebesseg

    ny = normalt(p);//beesesi meroleges
    nx.x = ny.y;
    nx.y = -ny.x;
    nx.z = ny.z;

    dpx = skalar_szorzat(nx,dp);
    dpy = skalar_szorzat(ny,dp);


    dpx2 = dpx *c2*c2/(c1*c1);
    dpy2 = sqrtl(c2*c2 - dpx2*dpx2);
    if(dpy<0) dpy2 = -dpy2;

    dp = (nx*dpx2 + ny*dpy2);
    }
    }



    void alapallapot()
    {
    g = (skalar)6.67428e-11;
    c = (skalar)2.99792458e8;
    m = (skalar)1.9891e30;

    rs = (skalar)2*m*g/(c*c);


    ry = 1392000e3;//sun R
    rx = 150e9;
    lepes = 2000000;
    dt = (rx*2/c)/lepes;


    p = vektor(-rx,ry,(skalar)0.0);
    dp = vektor((skalar)c,(skalar)0.0,(skalar)0.0);
    }

    int main()
    {
    alapallapot();
    newton_gravity();
    skalar fi = -atanl((p.y-ry)/rx)*180/M_PI;
    printf("%Le \n",fi);


    alapallapot();
    refract_gravity();
    skalar fi2 = -atanl((p.y-ry)/rx)*180/M_PI;
    printf("%Le \n",fi2);
    printf("%Le \n",fi2/fi);


    return 0;
    }