soncno obsevanje in klimatske spremembe po milankoviCevem modelu ŽIGA SMIT Fakulteta za matematiko in fiziko, Univerza v Ljubljani Institut Jozef Stefan PACS: 92.60.Ry, 92.30.Bc Osnova klimatskih modelov je toplota, ki jo Zemlja prejme od Sonca. Prispevek prikaze nazorno izpeljavo enacb, ki vodijo do dnevnega obseva, to je kolicine toplote, ki v poljubnem dnevu leta pade na enoto zemeljske povrsine v obliki krogle pri izbrani geografski sirini. Z analiticnimi priblizki v najnizjem redu je upostevano spreminjanje nagiba zemeljske vrtilne osi in spreminjanje ekscentricnosti Zemljinega tira. Z modelom je izracunana relativna kolicina toplote za obdobje würmske poledenitve (do 160 tisoc let nazaj). Primerjava z drugimi racuni in izmerjenimi klimatskimi spremembami kaze, da preprosti racun uspesno reproducira osnovno obliko klimatske krivulje. SOLAR IRRADIATION AND CLIMATIC CHANGES ACCORDING TO THE MILANKOVITCH MODEL Climatic models are essentially based on the energy, received by Earth from the Sun. The article shows methodic derivation of equations that lead to the daily insolation, i. e. the quantity of energy that hits a given area of spherically-shaped earth surface at given geographical amplitude on a given day. Analytical approximations are proposed for the time variation of the tilt of the geographical axis with respect to ecliptic and of the eccentricity of the Earth's orbit. The relative quantity of energy is calculated for the period of Wurm glace age (up to 160 thousand years in the past). A comparison with other calculations and measured climatic variations shows that the simple model reproduces well the main features of the climatic curve. Uvod Hitro spreminjanje klimatskih razmer v zadnjem času prav gotovo oživlja zanimanje za klimatske modele. Že sredi 18. stoletja so z opazovanjem geoloških pojavov zaznali pojav ledenih dob. Konec 19. stoletja so domnevali, da so nihanja hladnih in toplih obdobij posledica sprememb v sončnem obsevanju. Z idejo se je leta 1911 začel ukvarjati srbski matematik in naravoslovec Milutin Milankovič (1879-1958). Do podrobnosti jo je izdelal v avstro-ogrskem vojaskem ujetnistvu v blizini Budimpeste in jo leta 1920 objavil v knjigi, ki jo je izdala tedanja Jugoslovanska akademija znanosti in umetnosti v Zagrebu, izsla pa je pri zalozbi Gauthier-Villars v Parizu. Druga Milankovičeva knjiga je izsla v Beogradu leta 1941. V Milankovicevem modelu je klima odvisna od dnevne količine soncne toplote, ki posije na izbrani del zemeljske povrSine. Toplota je odvisna od kota, pod katerim padajo soncni zarki, in od dolzine dneva. Ti kolicini pa se dolgotrajno spreminjata zaradi precesije in opletanja zemeljske osi ter zaradi spreminjanja eliptičnosti zemeljskega tira. Milankovič je domneval, da je glavni vzrok za spreminjanje klime in nastop ledenih dob spreminjanje nagiba zemeljske osi, ki poteka s periodo 41 tisoč let. Toda ledene dobe nastopajo s periodo okoli 100 tisoč let, kar je nato kar za nekaj desetletij zmanjsalo zanimanje za Milankovičev model. Znova so ga ozivila globokomorska in ledeniska vrtanja, s katerimi je bilo mogoče rekonstruirati klimatske spremembe 500-800 tisoč let v preteklost. J. D. Hays, J. Imbrie in N. J. Shačkleton so v ključnem delu Variation in the Earth's Orbit: Pačema-ker of the Iče Ages leta 1976 z numeričnimi postopki določili osnovne periode pri spreminjanju klime in pokazali, da se ujemajo s spreminjanjem prečesije in nagiba zemeljske vrtilne osi. Vendar pa so te spremembe le priblizno periodične, ker gibanje Zemlje moti gravitačijsko polje tezjih planetov. Naloga je zato pogosto obrnjena: iz izmerjenih nihanj klime ugotavljajo variačije v Zemljini orbiti. Klimatska nihanja v preteklosti se kazejo v izotopski sestavi kisika v globokomorskih in ledenih vrtinah ter v sestavi in fizikalnih lastnostih se-dimentov. Merjenja teh pojavov so tudi uveljavljena metoda datiranja v geoloskih in starejsih arheoloskih obdobjih. Postopek, ko časovno skalo, ki jo dobimo iz klimatskih in geoloskih opazovanj, uskladimo s parametri po Milankovičevem modelu, imenujemo orbitalna uskladitev (orbital tuning). Ceprav vzroke za danasnje narasčanje temperatur isčemo v spremembah atmosfere, ki jih povzroča človek, pa je Milankovičev model se vedno osnova za računanje dolgoročnih klimatskih sprememb. S primerjavo podatkov za količino vpadle toplote pozimi in poleti ter pri različnih geografskih sirinah dobimo občutek, kako velike spremembe teh količin so pomembne. Izračun dnevnega obseva, to je količine toplote, ki pade na izbrani del zemeljske povrsine v enem dnevu, je zanimiv geometrijski problem, ki ga lahko z matematično nadarjenimi dijaki resimo ze v srednji soli. Namen prispevka pa je razviti model do te mere, da reprodučira osnovno obliko klimatske krivulje za obdobje würmske poledenitve, to je do 160 tisoč let nazaj. V ta namen predlaga analitične priblizke za časovno odvisnost eksčentričnosti zemeljskega tira in opletanje zemeljske vrtilne osi. Model Klima na Zemlji je odvisna od gostote energijskega toka s Sonča in od lastnosti atmosfere: odbojnosti za sončno svetlobo in absorptivnosti (emisivnosti) Slika 1. Vpadni kot sončnih žarkov na severni polobli poleti in pozimi opoldne. v infrardečem območju. Pri nasem modelu se bomo omejili na izračun dnevnega obseva, to je integrala gostote energijskega toka j (obsevanosti) v času od sončnega vzhoda do zahoda: Q f — ^ / j čos 6 did (1) Pri tem je 6 kot, ki ga sončni Zarki oklepajo z normalo na izbrani del zemeljske povrsine; 6 se spreminja zaradi vrtenja Zemlje in njenega gibanja okoli Sonča, saj zemeljska vrtilna os ni pravokotna na ekliptiko, ampak oklepa s pravokotničo nanjo kot ^q. Zaradi gibanja Zemlje po elipsi se med letom spreminja tudi j, ki pada s kvadratom razdalje Sonče-Zemlja. Za sevanje Sonča in s tem povezano sončno konstanto privzemimo, da se ne spreminja. Spremembe r od enega dneva do drugega pa so majhne, tako da lahko j v (1) postavimo pred integral. Za računanje kota 6 je ugodno, da začnemo meriti letni čas ti v trenutku, ko je Zemlja v pomladisču in je zemeljska os pravokotna na smer sončnih zarkov. Zemeljska os v splosnem oklepa s smerjo proti Sonču kot A: čos A = sin sin wltl, (2) pri čemer je wl kotna hitrost pri krozenju Zemlje okoli Sonča; wl = 2n/Tl, kjer je Tl perioda enega leta. Kot 6 zelimo podati kot funkčijo dnevnega časa td in letnega časa tl. Izberemo koordinatni sistem, pri katerem os z sovpada z zemeljsko vrtilno osjo. Osi x in y lezita v ekvatorialni ravnini; pri tem os y postavimo tako, da smer proti Sonču lezi v ravnini yz (slika 1). Sončni zarki oklepajo tedaj z osjo y kot Ö; vrednosti ö so v letni poloviči leta pozitivne, v zimski pa negativne. Kota ö in A se razlikujeta za n/2 (A = n/2 - 5), tako da velja sin ö = sin öo sin ti. (3) Smer proti Soncu je tedaj s = (0, - cos ö, sin ö). (4) Smer normale na izbrani toCki zemeljskega povrsja pri geografski sirini ^ je podana z vektorjem n = (cos wdtd cos — sin wdtd, cos sin ^). (5) Pri tem merimo dnevni cas td tako, da n ob casu td = 0 lezi v ravnini xz; Ud = 2n/Td, kjer je Td dolzina enega dneva. Iz enacb (4) in (5) dobimo cosQ: cos 6 = cos ö cos ^ sin udtd + sin ö sin (6) Iz zveze (6) lahko enostavno ocenimo, koliko je dan zaradi nagnjene zemeljske osi zjutraj in zvecer daljsi ali krajsi kot ob enakonocju. Za mejo svetlo-temno velja cos 6 = 0 in odtod — sin udtd = tan ö tan ^ = sin a, (7) pri cemer je a podaljsanje ali skrajsanje dneva v kotnih enotah. (V izvirni Milankovicevi izpeljavi namesto kota a nastopa kot n/2 + a, ki ustreza polovicni dolzini dneva v kotnih enotah.) Med arkticno nocjo in arkticnim dnevom je absolutna vrednost tan ö tan ^ vecja od ena. Arkticna noc ali dan nastopita pri geografski sirini ^ = n/2— | ö |, najmanjsa geografska sirina, pri kateri poznajo ta pojav, pa je ^ = n/2— | ö0 |. Zaradi loma svetlobe v ozracju (astronomske refrakcije) je dolzina dneva nekaj daljsa, kot jo napove enacba (6), pri arkticnem dnevu, ki se vlece nekaj mesecev, lahko tudi za nekaj dni. V enem dnevu na povrsino S pade toplota (1): TT+a Q = j f cos 6 d(udtd) = -j fcos u cos ^ cos a +(n + ^ sin ö sin ^ . Su^y \2 J —a (8) Pri tem smo upostevali, da se j v enem dnevu zaradi gibanja Zemlje po rahlo elipticnem tiru ne spremeni. V splosnem nas zanimajo relativne spremembe dnevnega obseva oziroma vpadle toplote Q. Za izhodi sce vzemimo vrednost Q0 pri ö = 0, torej koli-cino toplote, ki jo dobimo pri orientaciji zemeljske osi pravokotno na smer sončnih žarkov. Zaradi enostavnosti vzemimo referenčno vrednost toka j pri r = a, ko je oddaljenost Zemlje od Sonca enaka daljsi polosi elipse. Vrednost a je sicer vecja od povprecne razdalje Zemlja-Sonce, vendar pomeni boljso konstanto gibanja: med drugim nastopa v tretjem Keplerjevem zakonu in izrazu za energijo planeta. Relativno spreminjanje toplote Q, ki se na izbrano zemeljsko povrsino S izseva v enem dnevu, je tako Q / a V- Qo ^ / /n \ \ (^rJ (^cos Ö cos a ^ a + ^ sin S tan (9) Pri racunanju si pomagamo z izrazi za enacbo elipse v parametricni obliki. V tem zapisu so cas in koordinate planeta podani kot funkcija brezdimen-zijskega parametra ki se spreminja podobno kot polarni kot od 0 do 2n, vendar nima enostavne geometrijske ponazoritve. Za r velja (glej Landau, enacba 15.10): r = a(1 - ecos{), (10) pri cemer je e ekscentricnost elipse, parameter { pa je povezan s casom t' z zvezo 2n e - e sin e = = 2-ti. (11) Ti Pri majhnih vrednostih e je e priblizno enak polarnemu kotu, ki ga opise zveznica Sonce-Zemlja z zacetno lego v periheliju. Cas t', ki prav tako meri potovanje Zemlje okoli Sonca, ima torej v periheliju vrednost nic. Pri prakticnem racunanju resimo enacbo (11) numericno z iteracijo. Razmerje Q/Q0 zelimo izracunati za poljuben koledarski dan t v letu. V ta namen nam zadostuje, da pozabimo na prestopna leta in za dolzino leta vzamemo Tl = 365 dni. Po koledarju nastopi pomlad navadno 21. marca, tako da je od 1. januarja do pomladisca se 80 dni. Perihelij nastopi okoli 3. januarja. Pri casu tl v enacbi (3) in casu t' v enacbi (11) upostevamo casovna premika: tl = t - 80d; t' = t - 3d (12) V naslednjem koraku upostevamo dolgorocno spreminjanje Q/Q0, ki je opazno sele v obdobjih vec tisoc let. Zemeljska os v prostoru ni stalna, ampak zaradi precesije opleta po plascu stozca, ki je orientiran pravokotno na ekliptiko in ima odprtino 2S0. Ucinek tega gibanja je navidezno vrtenje elipse zemeljskega tira. Pri racunanju r iz enacb (10, 11) moramo torej pri casu t' upostevati dodatni fazni premik, ki je odvisen od geoloskega casa T; tega lahko merimo tako v preteklost kot v prihodnost: t' = t - 3d + TlT; (13) Tp Slika 2. Parametri Zemljinega tira, izračunana krivulja relativnega dnevnega obseva in dejansko spreminjanje temperature (prirejeno po Global Warming Art). Obmocje relativnega obsevanja do 160 tisoč let je označeno s kvadratkom. pri tem smo s Tp označili precesijski obhodni cas. Precesijski čas Tp glede na zvezde je 26 ka; s ka označujemo tisoče let. Cas, kot ga vidimo na Zemlji, je odvisen tudi od vrtenja glavne osi Zemljinega tira glede na zvezde in od drugih gibanj. Tudi tu je opaznih več harmoničnih komponent s periodami 19, 22 in 24 ka. V računu smo privzeli Tp =23 ka, kar se ujema z vrednostjo, ki so jo iz meritev na vrtini določili Hays in sodelavči. Drugi dve gibanji, ki vplivata na Q/Qo, sta spreminjanje nagiba zemeljske vrtilne osi glede na ekliptiko (kot 50) in spreminjanje eksčentričnosti elipse zemeljskega tira. Pri opisu teh dveh parametrov si pomagamo z analitičnimi priblizki. Spreminjanje nagiba zemeljske osi med skrajnima legama 22,1° in 24,5° poteka s periodo 41 ka. Trenutna vrednost je 23,44° in se zmanjsuje. V nasih računih podamo časovno odvisnost z analitičnim priblizkom T- 0,8 ka\ ^0 = 23,3° - 1,2° si^ 2n 41 ka / (14) Spreminjanje eksčentričnosti zemeljskega tira je kvaziperiodično zaradi vpliva masivnih planetov Jupitra in Saturna. Najmočnejse harmonične komponente imajo periode 95 ka, 125 ka in 413 ka. Analitično aproksimačijo, ki dovolj dobro opise spreminjanje e za zadnjih 400 ka, poisčemo z nastavkom e = A0 Ai sin ^o T + ^i 2n i=1 T,, (15) Slika 3. Spreminjanje nagiba zemeljske osi po enacbi (14) (levo). Spreminjanje ekscen-tricnosti Zemljinega tira po enacbi (16) (desno). Za prispevek s periodo 413 ka je amplituda 0,012. Minimum doseze pri času 400 ka v preteklosti, zato za izberemo -13 ka -1/4 x 413 ka. Preostale amplitude A in fazne čase ^ določimo z zvezami: 1. Največja vrednost e je 0,053, najmanjsa 0,005. 2. Danasnja vrednost e je 0,017. 3. Največja vrednost e je bila pred 200 ka. Iz pogoja (3) dobimo in z zahtevo, da morata biti ustrezni sinusni funkčiji v (15) enaki ena. Nato z upostevanjem pogojev (1) in (2) resimo sistem enačb za preostale tri amplitude A. Eksčentričnost je tako podana s semiempiričnim priblizkom / T + 33,75 ka\ e = 0,029+ 0,0062 sin 2n- ^^ , 95 ka / T - 13 ka \ - 0,012 čos 2n- + 0,0058 siM 2n / T - 18,75 ka 125 ka V 413 ka / (16) Rezultati in diskusija Model smo uporabili za prikaz nihanj v sončnem obsevanju za čas würmske poledenitve (do 160 tisoč let nazaj). Rezultate računa smo primerjali z vrednostmi, kot jih najdemo na spletu (slika 2). Kot kaze slika 3, z enačbo (14) pravilno opisemo periodo in začetni potek nihanja zemeljske osi, ne zajamemo pa rahlega spreminjanja amplitude. Spreminjanje eksčentričnosti o 0 50.4 0.3 45^ poleti pozimi O 20 40 60 SO 100 120 140 160 ap (ka) Slika 4. Izračunano relativno obsevanje za 45° severne Sirine. zemeljskega tira do 400 tisoč let nazaj pa zelo dobro ponazorimo s priblizkom (16). Enačba (16) napove enako strukturo tudi za obdobje med 400 in 800 tisoč leti v preteklost, vendar slika 2 kaze, da se oblika funkčije pri bolj oddaljenih obdobjih nekoliko popači. Krivulji relativnega sončnega obseva sta izračunani za dve geografski sirini: 45° (slika 4), kamor priblizno sodijo nasi kraji, in za 65° (slika 5). Sirine okoli 65° naj bi bile namreč odločilne za nastanek ledenih dob. Te nastopijo tedaj, kadar deli zemeljske povrsine dobijo poleti premalo toplote, da bi se stalil ves led, ki je ostal od zime. Sliko 5 lahko primerjamo s podrobnostmi na sliki 2: kot vidimo, nas račun dokaj dobro reprodučira posamezne strukture: vrh z visokimi temperaturami pred 10 tisoč leti, sedlast vrh z visokimi temperaturami med 30 in 45 tisoč leti, visoke temperature pred 75 in 100 tisoč leti, močno ohladitev pred 115 tisoč leti in zelo visoke temperature pred 125 tisoč leti (ta vrh navadno označujejo kot 5e ali 5,5 in ustreza eemianskemu interglačialu). Ohladitve nastopajo pred 20, 60 in 115 tisoč leti in stejejo kot poledenitve würm III, II in I. Od teh naj bi bila najhladnejsa würm I (pred 115 tisoč leti), vendar dejanske klimatske meritve kazejo, da so bile tedaj temperature visje kot pri drugih dveh. Sliki 4 in 5 pa pokazeta, da na izbrani del zemeljske povrsine vedno posije največ toplote poleti in da je ravno ta toplota odločilna za klimatske spremembe. Sklep V prispevku smo prikazali preprost model za spremljanje klimatskih nihanj po Milankovičevem predlogu. Model se omeji na bistvene pojave, ki zajemajo relativno velikost dnevnega obseva in spreminjanje osnovnih astronom- 2 9' 2 3 2 7 26. 2,5. O 2,3 2,2. 0 05 0,00. 65° poleti pOZitni SO BO 100 BP (ka) Slika 5. Izračunano relativno obsevanje za 65° severne širine. skih parametrov. Računanje po predstavljenih izrazih je mogoče opraviti s preprostim računalniskim programom ali čelo z navadnim kalkulatorjem. Seveda pa se mora braleč zavedati, da je spreminjanje klime odvisno se od mnogih drugih pojavov, ki niso zajeti v sedanjem modelu: povrsine Zemlje v obliki elipsoida, odbojnosti tal in vplivov atmosfere, kamor spadajo tako absorpčija sončnega sevanja kot učinki tople grede, in prenosov toplote z vetrovi in morskimi tokovi. Del teh pojavov je v model vključil ze Milankovič v prvi izdaji leta 1920. Z modelom lahko pogledamo tudi v prihodnost. Račun kaze, da bi v naslednjih nekaj tisoč letih morali pričakovati zmanjsanje temperature. Ce se to ne bo zgodilo, bodo za to krivi toplogredni pojavi. LITERATURA [1] M. J. Aitken, Science-based dating in archaeology, Longman, London 1990. [2] J. D. Hays, J. Imbrie in N. J. Shackleton, Variation in the Earth's Orbit: Pacemaker of the Ice Ages, Science 194 (1976), 1121-1132. [3] L. D. Landau in E. M. Lifshitz, Mechanics, Pergamon Press 1976. [4] M. Milankovitčh, Theorie Mathematique des Phenomenes Thermiques Produits par la Radiation Solaire, Gauthier-Villars, Paris 1920. [5] M. Milankovitčh, Kanon der Erdbestrahlung und seine Andwendung auf das Eiszeiten-problem, Srpska kraljevska akademija, Beograd 1941. [6] Wikipedia; gesla: Milankovitčh cycles, Milutin Milankovitčh http://www.dmfa-zaloznistvo.si/