Zeruko gorputzen dinamika

Indar zentrala eta kontserbakorra

Adibideak

Saiakuntza

Programazio-kodea

Orri honetan aztertuko dugu, nola mugitzen den gorputz bat, erakarpen-zentro batetik r1 distantziara kokatu, eta v1 abiaduraz jaurtitzen denean. Norabide erradialarekiko perpendikularki jaurtitzen da, gero eta abiadura handiagoaz, alboko irudiak erakusten duen bezala. Gorputz horren masa m da.

Kepler-en legeek zehazki deskribatzen dute planeten higidurak nolakoak diren eguzkiaren inguruan, baina ez dute azaltzen zein den higidura horren zergatia.

1.-Planetek orbita eliptikoak jarraitzen dituzte eta Eguzkia elipsearen foku batean dago.

2.-Eguzkitik planetaraino doan posizio-bektoreak denbora-tarte konstanteetan azalera konstantea ekortzen (edo estaltzen) du.

3.-Orbiten periodoen karratuak eta elipseen ardatzerdien kuboak proportzionalak dira.

Newton-en grabitazioaren legeak Kepler-en legeei azalpena ematen die, baina are gehiago, planeten ibilbideak denboraren menpe kalkulatzeko balio du. Orokorrean gorputz batek, erakarpen grabitatorioa jasaten duenean, ibilbide laua deskribatzen du, "konika" bat: elipsea, parabola edo hiperbola.

Aurreko atal batean komentatu da, indar grabitatorioa zentrala eta kontserbakorra denez, eguzkiaren eta planeta baten arteko erakarpen-indarrak lehen ordenako ekuazio diferentzial bi ematen dituela. Ekuazio diferentzial bi horiek koordenatu polarretan adierazten badira, eta denbora eliminatzen bada, ibilbidearen ekuazioa ematen du, izan ere konika baten ekuazioa.

Orri honetako programa interaktiboak bestela jokatzen du: planetaren posizioa jakinda, azelerazioaren X eta Y osagaiak kalkulatzen ditu, eta hasierako abiadura ezaguna denez, abiadura eta posizio berriak kalkulatzen ditu, integralak prozedura numerikoez ebatzita.

 Higidura ekuazioen soluzio numerikoa

Demagun m masadun partikula batek (planeta batek) indar erakarle bat jasaten duela puntu finko baterantz (Eguzkirantz). Suposatuko dugu Eguzkia finkoa dela, bere masa, M, partikularena baino askoz handiagoa delako.

Planetak jasaten duen F indarrak norabide erradiala du, eta eguzkiaren zentrorantz apuntatzen du uneoro. Grabitazio unibertsalaren legeak indar horren modulua ematen du:

Adierazpen horretako r da, partikulatik indar-zentroraino dagoen distantzia. Defini ditzagun koordenatu cartesiarrak eguzkia koordenatuen jatorrian kokatuta: Indarraren osagaiak honela idatz daitezke:

Newton-en bigarren legea aplikatuz, eta azelerazioa idazten bada posizioaren bigarren deribatu gisa, ekuazio diferentzial bi lortzen dira bigarren ordenakoak:

Hasierako baldintzak ezagututa (posizioa eta hasierako abiadura), ekuazio diferentzial bikote hori prozedura numerikoez ebatz daiteke, esaterako Runge-Kutta metodoa.

Eskalak

Ekuazio diferentzialak ebazten hasi aurretik, komenigarria da ordenagailuaren kalkuluetan zenbaki handiegiak edota txikiegiak saihestea. Horretarako ekuazioak berridatz daitezke:

Aukera dezagun unitate-sistema egokia: lehenik, luzera neurtzeko, unitate astronomikoa, alegia, Lurraren eta Eguzkiaren arteko distantzia: L=UA bat =1.496·10¹¹ m. Eta bigarrenik, denbora neurtzeko, urtea,  P= 1 urte= 365.26 egun =3.156·10⁷ s.

Unitate sistema horretan aldagaiak honela berridazten dira: x=Lx’, eta t=P·t’, Eta lehen ekuazio diferentziala honela:

Lurrak Eguzkiaren inguruan duen orbitaren ardatzerdi nagusia L da, eta periodoa P, bira bat osatzeko tardatzen duen denbora. M eguzkiaren masa denez, Kepler-en hirugarren legearen arabera:

Unitate-sistema berrian, ekuazio bikotea honela idazten da (primak kendu dira erosoago jarduteko):

Energiaren kontserbazio-printzipioa

Partikularen energia totala konstantea da uneoro, eta hasieran, t=0 aldiunean, hau da:

Partikularen energia totala E₀<0 bada, partikula harrapatuta geratuko da erakarlearen inguruan, eta aldiz, E₀≥0 bada, bere orbitak azkenean infinitura ihes egingo du.

Partikularen energia edozein t aldiunetan:

Eta unitate-sistema egokian:

v=v’·L/P,   x=x’·L,   y=y’·L,   d=d’·L

Etor gaitezen berriro hasierako unitateetara, eta defini dezagun e energia masa unitateko:

Ondorengo programa interaktiboak uneoro honako zatidura kalkulatzen du:

zatidura horri errorearen portzentajea deituko diogu. Kalkulu numerikoan zehar partikularen energia ez litzateke aldatu behar, baina aldatzen joango balitz, eta zatidura hori "bat" baino handiago izatera iritsiko balitz, programak kalkulua geldiaraziko du. Kalkulatutako orbita, orbita errealetik nahikoa urruntzen ari dela esan nahi du. Hori gerta daiteke, partikula erakarpen zentrotik oso hurbil pasatzen denean.

Indar zentral eta kontserbakorra

Irudiak erakusten du partikula bat orbita eliptiko batean mugitzen, eta erakarpen-zentroa elipsearen foku batean dago. Hurbiltze maximoko distantzia r₁ da, eta urruntze maximokoa r₂. Partikulak duen abiadura, bi posizio horietan, v₁ eta v₂ dira hurrenez hurren.  Partikulak posizio horietan duen energia eta momentu angeluarra idatz daitezke, eta biak konstanteak direnez:

Ekuazio bi horiek ebazteko bi aukera izaten dira:

• r₁ eta r₂ ezagututa, kalkulatu behar dira v₁ eta v₂.

Lehen ekuazioan v₁ bakandu eta bigarrenean ordezkatu. Hona emaitza:

Zenbait ikaslek egiten duten ohiko akatsa izaten da, v₁ eta v₂ kalkulatzea baina higidura zirkular uniformearen dinamikan oinarrituta, erakarpen-zentroa eta elipsearen kurbadura-zentroa nahasi egiten dituztelako. Dei diezaiogun R elipsearen kurbadura-erradioari, esaterako, hurbiltze maximoko posizioan. Orduan, higidura zirkular uniformearen dinamika honela idatz daiteke:

Emaitza hori ikusita R da, elipsearen ekuazioa koordenatu polarretan adierazten denean agertzen den d parametroa.

• Beste aukera izaten da, r₁ eta v₁ ezagututa, kalkulatu behar dira r₂ eta v₂.

Jaurtiketaren uneko r₁ posizioa eta v₁ abiadura ezagunak badira, partikularen energia idatz daiteke, eta konstantea denez, v₂ kalkula daiteke. Zenbait kalkulu burutu ondoren:

       (1)

Ihes-abiadura

Partikula bat erakarpen-zentrotik r₁ distantziara dagoenean, eta jaurti egiten bada (esaterako, erakarpen zentrotik urrunduz), badago jaurtitze-abiadura minimo bat, partikula infinituraino iristeko: dei diezaiogun v_e jaurtitze-abiadura minimo horri (e azpindizea erabili da, gaztelerazko "escape"-gatik).

Esaterako, Lurraren gainazaletik ihes egiteko abiadura minimoa hau da, v_e=11190.7 m/s. (r₁=6.37·10⁶ m, eta M=5.98·10²⁴ kg).

Izan ere, v₂ abiadura, urruntze maximoarena (1), ihes abiaduraren menpe idatz daiteke:

Eta momentu angeluarraren kontserbazioa aplikatuz urruntze maximoaren distantzia ere kalkula daiteke:

Higidura zirkularra

Partikula jaurtitzen bada, justu abiadura jakin batez (vc), orduan orbita zirkularra deskribatuko du r erradioaz:

Erraz egiazta daiteke (1) formulan, partikulak orbita zirkularra jarraitzen badu, orduan eta  v₂=v₁.

Lurraren inguruko satelite artifizialen kasuan: Jaurtitze-abiadura, v1, orbita zirkularraren vc abiadura baino txikiagoa bada, orduan jaurtiketaren puntua apogeoa da. Jaurtitze-abiadura, v1, orbita zirkularraren vc abiadura baino handiagoa bada, orduan jaurtiketaren puntua perigeoa da. Irudiko orbita gorria zirkularra da eta beste biak eliptikoak, baina jaurtiketaren puntua, eskumakoa, apogeoa da berdean eta perigeoa urdinean.

Adibideak:

1 adibidea. Lurra

Lurraren orbita ia zirkularra da eta bere erradioa r=1.496·10¹¹ m=1 UA. Higidura zirkular uniformearen dinamika aplikatuz, Lurraren abiadura kalkula daiteke. Eguzkiaren masa, M=1.98·10³⁰ kg, denez, lurraren abiadura v=29711.8 m/s. Lehen atalean azaldutako unitate sistemaren arabera:

• Luzera neurtzeko unitatea, L=1UA=1.496·10¹¹ m
• Denbora neurtzeko unitatea, P=urtebete=365.26 egun=3.156·10⁷ s.

abiadura unitate horietan: v'

v'= 6.268 UA/urte

Orbita bat osorik betetzeko denbora P=1 urte.

Ondorengo programa interaktiboan honako datuak sartuz gero: x=1.0 eta v_y=6.27, orduan eguzkiaren inguruan orbita zirkularra ateratzen da.

2 adibidea. Martitz

Hona hemen Martitzen datuak: Ardatzerdi nagusia: a=1.524 UA, Eszentrikotasuna: ε=0.093.

• Eguzkiraino hurbiltze maximoa (perihelioa) r₁=a-c=a-εa=1.382 UA=2.068·10¹¹ m

• Eguzkiraino urruntze maximoa (afelioa) r₂=a+c=a+εa=1.666 UA=2.492·10¹¹ m

Distantzia maximo eta minimoak ezagututa, alegia, r₁ eta r₂ perihelioko abiadura kalkula daiteke, v₁:

v₁=26420.7 m/s=5.573 UA/urte

Ondorengo programa interaktiboan honako datuak sartzen badira: x=1.382 eta v_y=5.573, orduan, bereziki, Martitzek Eguzkiaren inguruan burutzen duen elipsea lortzen da, eta orbita bat osorik betetzeko denbora P=1.86 urte.

 Saiakuntza

Aurretik marraztutako ibilbideak, Ezabatu botoiaz garbitzen dira, nahi denean.

Ariketa:

6.-Idatz bedi taulan, orbita osoa betetzeko tardatzen duen denbora, P, urteetan, eta egiazta bedi Kepler-en hirugarren legea taularen azken zutabean (P²/a³ zatidura gutxi gora behera konstantea izan behar da).

public abstract class RungeKutta {
	double h;
RungeKutta(double h){
	this.h=h;
}
void setPaso(double dt){
	this.h=dt;
}
public void resolver(Estado e){
	//variables auxiliares
	double k1, k2, k3, k4;
	double l1, l2, l3, l4;
	double q1, q2, q3, q4;
	double m1, m2, m3, m4;
//estado inicial
	double x=e.x;
	double y=e.y;
	double vx=e.vx;
	double vy=e.vy;
	double t=e.t;

	k1=h*vx;
	l1=h*f(x, y, vx, vy, t);
	q1=h*vy;
	m1=h*g(x, y, vx, vy, t);

	k2=h*(vx+l1/2);
	l2=h*f(x+k1/2, y+q1/2, vx+l1/2, vy+m1/2, t+h/2);
	q2=h*(vy+m1/2);
	m2=h*g(x+k1/2, y+q1/2, vx+l1/2, vy+m1/2, t+h/2);

	k3=h*(vx+l2/2);
	l3=h*f(x+k2/2, y+q2/2, vx+l2/2, vy+m2/2, t+h/2);
	q3=h*(vy+m2/2);
	m3=h*g(x+k2/2, y+q2/2, vx+l2/2, vy+m2/2, t+h/2);

	k4=h*(vx+l3);
	l4=h*f(x+k3, y+q3, vx+l3, vy+m3, t+h);
	q4=h*(vy+m3);
	m4=h*g(x+k3, y+q3, vx+l3, vy+m3, t+h);

	x+=(k1+2*k2+2*k3+k4)/6;
	vx+=(l1+2*l2+2*l3+l4)/6;
	y+=(q1+2*q2+2*q3+q4)/6;
	vy+=(m1+2*m2+2*m3+m4)/6;
	t+=h;

//estado final
	e.x=x;
	e.y=y;
	e.vx=vx;
	e.vy=vy;
	e.t=t;
}
abstract public double f(double x, double y, double vx, double vy, double t);
abstract public double g(double x, double y, double vx, double vy, double t);

}

public class Planeta extends RungeKutta{
Planeta(double h){
	super(h);
}
public double f(double x, double y, double vx, double vy, double t){
	double r=Math.sqrt(x*x+y*y);
	double z=-4*Math.PI*Math.PI*x/(r*r*r);
	return z;
}
public double g(double x, double y, double vx, double vy, double t){
	double r=Math.sqrt(x*x+y*y);
	double z=-4*Math.PI*Math.PI*y/(r*r*r);
	return z;
}
public double energia(double x, double y, double vx, double vy){
	double r=Math.sqrt(x*x+y*y);
	double z=(vx*vx+vy*vy)/2-4*Math.PI*Math.PI/r;
	return z;
}
}

//Objetos de la clase Planeta

	Estado estado=new Estado(0.0, x, 0.0, 0.0, vy); //estado inicial
	Planeta planeta=new Planeta(0.005);
 	eInicial=planeta.energia(estado.x, estado.y, estado.vx, estado.vy);

//rutina que calcula la trayectoria paso a paso
	planeta.resolver(valor);
	double energia=planeta.energia(estado.x, estado.y, estado.vx, estado.vy);
	double error=Math.abs((energia-eInicial)*100/eInicial);
	if(error>1.0){
		//detiene el movimiento
	}