Program sluzy do rysowania krzywej modelowej predkosci radialnej gwiazdy, wywolanej istnieniem planet wokol tej gwiazdy. Wszystkie parametry trzeba wpisac w kodzie ( program o nic nas nie pyta po uruchomieniu ). Nalezy ustawic liczbe planet N ( w przykladzie ponizej sa dwie planety ), poczatkowy czas t, liczbe krokow czasowych LK, dlugosc kroku h. Nalezy takze podac parametry orbitalne V0,K,P,e,w,tau(w kodzie podalem przykladowe parametry i wyjanilem co oznaczaja). Nalezy takze podac nazwe piku wyjsciowego, do ktorego zostana zapisane punkty krzywej predkosci radialnej. W programie mozna zmieniac czas poczatkowy krzywej, liczbe i dlugosc kroku czasowego. ---------- <verbatim> / *modelowa_krzywa.c* / #include<stdio.h> #include<math.h> #define DZIEN 86400 / *liczba sekund w 1 dniu* / #define N 2 / *liczba planet* / #define B 0.0001 / *dokladnosc rozwiazania r-nia Keplera* / #define LK 25000 / *liczba krokow* / #define PI 3.1415926536 / *liczba pi* / float NR(float P, float tau, float t, float e, float b, float pi); main() { int i,k; float h,t,Vr,niu,E; float V0,K[N],P[N],e[N],w[N],tau[N]; FILE *wynik; char *nazwa; t=400*DZIEN; /* poczatkowy czas */ h=0.1*DZIEN; /* dlugosc kroku */ /*parametry orbitalne - tu dokonujemy zmian */ V0=-5.5; /* predkosc srodka masy */ K[0]=15; /* ampituda predkosci radialnej wywolana 1. planeta */ P[0]=1080*DZIEN; /* okres 1. planety */ e[0]=0.19; /* ekscentycznosc orbity 1. planety */ w[0]=81*PI/180; /* dlugosc perycentrum orbity 1. planety */ tau[0]=539*DZIEN; /* czas przejscia przez pericentrum 1. plnety */ K[1]=6; /* podobnie dla kolejnych planet */ P[1]=10.5262*DZIEN; /* -- */ e[1]=0.016; /* -- */ w[1]=53*PI/180; /* -- */ tau[1]=101.4*DZIEN; /* -- */ nazwa="HD114729_modelowa.dat"; /* nazwa pliku wyjsciowego */ wynik=fopen(nazwa,"w"); for (k=1; k<=LK; ++k) { Vr=V0; for (i=0; i<=N-1; ++i) { E=NR(P[i],tau[i],t,e[i],B,PI); niu=2*atan(sqrt((1+e[i])/(1-e[i]))*tan(E/2)); Vr=Vr+K[i]*(cos(w[i]+niu)+e[i]*cos(w[i])); } fprintf(wynik,"%10.4f%10.2f%10.2f\n",t/DZIEN,Vr,1); t=t+h; } fclose(wynik); return 0; } /*funkcja obliczajaca E w czasie t metoda Newtona-Raphsona IV rzedu z dokladnoscia b*/ float NR(float P, float tau, float t, float e, float b, float pi) { float MM,EE,f0,f1,f2,f3,d1,d2,d3; MM=2*pi*(t-tau)/P; EE=MM; f0=EE-e*sin(EE)-MM; f1=1-e*cos(EE); f2=e*sin(EE); f3=e*cos(EE); d1=-f0/f1; d2=-f0/(f1+d1*f2/2); d3=-f0/(f1+d2*f2/2+d2*d2*f3/6); while ((d3>b) || (d3<-b)) { EE=EE+d3; f0=EE-e*sin(EE)-MM; f1=1-e*cos(EE); f2=e*sin(EE); f3=e*cos(EE); d1=-f0/f1; d2=-f0/(f1+d1*f2/2); d3=-f0/(f1+d2*f2/2+d2*d2*f3/6); } return EE; } </verbatim> -- Main.CezaryMigaszewski - 26 Feb 2005
This topic: Main
>
TWikiUsers
>
CezaryMigaszewski
>
ModelowaKrzywa
Topic revision: revision 3 (raw view)
Copyright © CC-BY-SA by the contributing authors. All material on this collaboration platform is copyrighted under CC-BY-SA by the contributing authors unless otherwise noted.
Ideas, requests, problems regarding Foswiki?
Send feedback