Astronomija

Energijos taupymas Barneso-Huto algoritme

Energijos taupymas Barneso-Huto algoritme


We are searching data for your request:

Forums and discussions:
Manuals and reference books:
Data from registers:
Wait the end of the search in all databases.
Upon completion, a link will appear to access the found materials.

Parašiau kodą, kuris atitinka Barneso-Huto gravitacinės dinamikos algoritmą. Viskas atrodo gana gerai, išskyrus tai, kad kai apskaičiuoju bendrą energiją $$ E_ {tot} equiv E_k + E_u = const $$, gaunu tokį rezultatą:

ir jūs galite pamatyti, kad nors tai atrodo pastovus didžiąją laiko dalį, pradžioje jis vis labiau praeina. Kode neradau klaidos, kuri galėtų būti atsakinga už tai. Ar gali būti kitas paaiškinimas, galintis tai suprasti?


SUPERBOX - efektyvus kodas be susidūrimų galaktikos dinamikai

Pateikiame „Superbox“ - dalelių tinklelio kodą su didelės skiriamosios galios tinkleliais ir NGP (artimiausio tinklo taško) jėgos skaičiavimo schemą, pagrįstą antrosiomis potencialo išvestinėmis. „Superbox“ įdiegia greitą mažai saugomą FFT algoritmą, suteikiantį galimybę dirbti su milijonais dalelių staliniuose kompiuteriuose. Bandymo skaičiavimai rodo energijos ir kampinio impulso išsaugojimą vienai daliai per 10 5 per perėjimo laiką. Tinklelio ir skaitinio atsipalaidavimo poveikis išlieka nereikšmingas, net jei šie skaičiavimai apima Hablo evoliucijos laiką. Kadangi antriniai tinkleliai seka atskirų galaktikų trajektorijas, kodas leidžia labai tiksliai išspręsti sąveiką galaktikų grupėse, pvz., Didelio greičio susidūrimus tarp elipsinių galaktikų ir nykštukinių galaktikų potvynio ardymą. Puikus susitarimas pasiekiamas lyginant su tiesiogine suma N- kūno kodas, veikiantis specialios paskirties „Grape 3“ aparatinėje įrangoje. Palydovinių galaktikų orbitinis skilimas dėl dinaminės trinties, gautos naudojant „Superbox“, sutinka su Chandrasekhar & # x27s gydymu, kai Coulombo logaritmas lnΛ≈1.5.


Pavadinimas: periodiškai greitas daugiapolis metodas: optimalios žaliosios funkcijos aproksimavimas

Greitas daugiapolis metodas (FMM) laikosi periodinių ribinių sąlygų „savaime“, jei jis naudoja periodinę „Green“ funkciją, kad apskaičiuotų daugiapolio išplėtimą kiekvieno FMM aštuonmedžio mazgo sąveikos zonoje. Galima apibrėžti „optimalią“ Žaliąją funkciją tokiam metodui, kurio rezultatas yra skaitinis sprendimas, kuris konverguoja į lygiavertį „Dalelių-tinklelio“ tirpalą pakankamai aukštų daugybinių polių riboje. Galima išvesti atskirą funkcinę optimalios „Green“ funkcijos lygtį, tačiau ji praktiškai nėra naudinga, nes jos sprendimo būdai nėra žinomi. Vietoj to, šiame dokumente pateikiamas optimalios žaliosios funkcijos apytikslis dydis, kuris yra tikslus, lygus geriau nei 1e-3 LMAX normoje ir 1e-4 L2 normoje, kad būtų praktiškai naudingas kelių polių skaičius. Tokia apytiksliai optimali „Green“ funkcija suteikia praktišką būdą FMM įgyvendinti su periodinėmis ribinėmis sąlygomis „savaime“, nereikalaujant apskaičiuoti grotelių sumų ar pasikliauti hibridiniais FMM-PM metodais.


Medienos modifikavimas: mechaninės savybės ir pritaikymas

4.2. Kompozicinės struktūros įtaka higromechaniniam elgesiui

Sutankėjusių dribsnių pasiskirstymas ir savybės lemia plokštės keliamąją galią. Orientuotos medžio drožlių plokštės gali padidinti medinių kompozitų tvirtumą ir standumą. Nepaisant to, tokių plokščių deformacija yra labai didelė (McNatt 1973). Dėl to medienos kompozicinių plokščių kokybė priklauso nuo dribsnių savybių ir išsidėstymo (Geimer 1982).

Paprastai mediniai kompozitai neturi simetrijos vidurio plokštumos atžvilgiu, jų savybės gali būti asimetrinės arba antisimetrinės (Brauns ir Rocens 1994a). Kompozitas analizuojamas kaip lygiagrečių elementarių sluoksnių sistema. Sluoksnio kryptis k yra ϕ ( k ) (8 pav.). Vietinė koordinačių sistema siejama su pagrindinėmis elementinių sluoksnių kryptimis. Siekiant apytiksliai įvertinti įtempio būseną, modelyje naudojami nuo įtempio gradiento priklausomi jėgos ir poros įtempiai.

8 paveikslas. Daugiasluoksnis medžio kompozito modelis.

Padermių ir įtempių priklausomybė nuo drėgmės kiekio pokyčių, ΔW, pagrindinio sluoksnio ašyse galima parašyti taip:

Eqns specifinės drėgmės padermės β ii ′ ir suvaržantys įtempiai r jj ′. (9) ir (10) yra susijusios su kompozicinės medžiagos standumu (Tsai ir Hahn 1980, Brauns ir Rocens 2004).

Daugiasluoksnis modelis ir modelis, pagrįstas laminato analogija (Halpin ir kt. 1971) gali būti naudojami sluoksniuotų ir pluoštinių kompozitų elgesiui nustatyti. Nesimetriškas drėgmės pasiskirstymas sukelia tiesinį išsiplėtimą ir higromechaninį deformaciją. Kreiviai, atsirandantys dėl medienos higrodeformacijos, buvo apskaičiuoti naudojant Braunso ir Rocenso (1994a, 2004 ir žr. Medienos ir betono elastingumas: higromechaniniai efektai) sukurtą metodą.

Norėdami parodyti dalelių išlyginimo poveikį higromechaninėms ir mechaninėms savybėms, laminatai su elementariais sluoksniais, skirtingai orientuoti, pavyzdžiui, į mašinos kryptį (x1) (Brauns ir Rocens 1997). Bet kurio laminato vidutinis kampas (neatsižvelgiant į ženklą) apibrėžiamas kaip φ - ir santykinis išlyginimas ξ gali būti parašyta taip:

apibrėžiant atsitiktinę lentą su nuline lygiuote.

Šis modelis buvo pritaikytas mediniams kompozitams, susidedantiems iš fanerų, riboto ilgio dribsnių ir klijų (Brauns ir Rocens 1991). Atsižvelgiama į tankinimo poveikį gaminant medžiagą (Brauns ir Rocens 1994b). Analizuojant naudojamos tankinto beržo lukšto techninės charakteristikos: elastingumo moduliai E1= 16 500 ir E2= 700 (MPa) šlyties moduliai G 12 = 900, G 13 = 1500 ir G23= 300 (MPa) Puasono koeficientai ν 21 = 0,45, ν 31 = 0,34 ir ν 32 = 0,30. Elastingumo ir šlyties modulių drėgmės pokyčio korekcijos yra: α 1 = 250, α 2 = 25, α 12 = 25, α 13 = 30 ir α23= 20 (MPa). Specifinės drėgmės padermės esant tam tikroms pradinėms sąlygoms yra: β 1 = 0,00003, β 2 = 0,002, β3=0.003 (%) − 1 .

Medinių kompozitų vidurinės plokštumos deformacijos ir kreivės priklauso nuo drėgmės kiekio ir pasiskirstymo. Higroskopinė faneros deformacija pateikta 9 pav. Modelyje yra septyni vienodo storio 1,4 mm lukštai. Dalinis faneros (beržo) tankis yra 660 kg m - 3. Rezultatai, lyginantys tankintos ir įprastos medienos savybes, iliustruoja asimetrinio netiesinio drėgmės pasiskirstymo ir antisimetrinės struktūros poveikį kreivumui esant vienodam drėgmės kiekiui.

9 paveikslas. Storio patinimas (ɛ3) iš faneros: 1, įprasta mediena 2, tankinto medžio išlinkimai: 3, κ1 4, κ2 (nesimetriškas drėgmės pasiskirstymas) iškrypimas: 5, κ6× 10 (antisimetrinė struktūra). Po Brauno J, Rocens K 1997 Higromechaninis medinių kompozitų elgesys. Wood Sci. Technol. 31, 193–204.

Ryšys tarp linijinio išsiplėtimo, storio išsipūtimo ir dribsnių išlyginimo, kai drėgmės kiekis padidėja 12%, parodytas 10 pav. Norint pagerinti lentos našumą, pageidautina formuoti lentas su santykiniu dribsnių išlyginimu ξ = 0 - 40%, kuriai dribsnių išlyginimo poveikis išilginėms išsiplėtimo deformacijoms ɛ1 ir ɛ2 yra nereikšmingas, o storis-patinimas ɛ3 pasiekia maksimumą. 10 paveiksle taip pat parodomi ryšiai tarp dribsnių išlyginimo ir specifinės koncentruotos reakcijos F* lenkiant (esant 1 cm įlinkiui) ir kritiniam įtempiui σ 1 cr suspaudžiant šarnyrinę kvadratinę plokštę.

10 paveikslas. Dribsnių išlyginimo poveikis plokščių patinimuisi ir mechaninėms savybėms, esant 12% drėgmės kiekiui: 1, ɛ1 2, ɛ2 3, ɛ3 (5) 4, F * (kN) 5, σ 1 kr (MPa). Po Brauno J, Rocens K 1997 Higromechaninis medinių kompozitų elgesys. Wood Sci. Technol. 31, 193–204.


SPACE / MATH - konstantos ir skaičiavimų lygtys

Šį sąrašą iš pradžių sudarė Dale'as Greeris. Būtų dėkingi papildymai.

Skaičiai skliaustuose yra apytiksliai, kurie bus naudojami daugeliui mėlynojo dangaus tikslų.

„Unix“ sistemos teikia „vienetų“ programą, naudingą konvertuojant tarp skirtingų sistemų (metrinė / anglų, CGS / MKS ir kt.)

7726 m / s (8000) - Žemės orbitos greitis 300 km aukštyje 3075 m / s (3000) - Žemės orbitos greitis 35786 km (geosinchroninis) 6371 km (6400) - vidutinis Žemės spindulys 6378 km (6400) - Pusiaujo žemės spindulys 1738 km (1700) - Vidutinis Mėnulio spindulys 5,974e24 kg (6e24) - Žemės masė 7,348e22 kg (7e22) - Mėnulio masė 1,989e30 kg (2e30) - Saulės masė 3.986e14 m ^ 3 / s ^ 2 (4e14) - Žemės masės gravitacijos pastoviųjų laikų masė 4.903e12 m ^ 3 / s ^ 2 (5e12) - Mėnulio gravitacinių pastoviųjų laikų masė 1.327e20 m ^ 3 / s ^ 2 ( Saulės gravitacinių pastoviųjų laikų masė 384401 km (4e5) - Vidutinis atstumas nuo Žemės iki Mėnulio 1,496e11 m (15e10) - Vidutinis Žemės ir Saulės atstumas (astronominis vienetas)

1 megatonas (MT) TNT = apie 4,2 e15 J arba energijos ekvivalentas apie 0,05 kg (50 g) medžiagos. Nuoroda: JR Williamsas, „Daiktų energijos lygis“, Karinių oro pajėgų specialiųjų ginklų centras (ARDC), Kirtlando oro pajėgų bazė, Naujoji Meksika, 1963 m. Taip pat žr. „Branduolinių ginklų padariniai“, sudaryti S. Glasstone'o ir PJ Dolano. , paskelbtą JAV gynybos departamento (gauti iš GPO).

Kur d yra atstumas, v yra greitis, a yra pagreitis, t yra laikas. Papildomų specializuotų lygčių galima rasti:

Pastoviam pagreičiui d = d0 + vt + .5at ^ 2 v = v0 + esant v ^ 2 = 2ad

Pagreitis cilindre (kosminėje kolonijoje ir kt.), Kurio spindulys r ir sukimosi periodas t:

Apskritinėms Keplerio orbitoms, kur: Vc = apskritos orbitos greitis Vesc = pabėgimo greitis M = bendra orbituojančių ir orbituojamų kūnų masė G = gravitacijos konstanta (apibrėžta toliau) u = G * M (galima išmatuoti daug tiksliau nei G arba M ) K = -G * M / 2 / ar = orbitos spindulys (matuojamas nuo sistemos masės centro) V = orbitos greitis P = orbitos periodas a = pusiau didelė orbitos ašis

Vc = sqrt (M * G / r) Vesc = sqrt (2 * M * G / r) = sqrt (2) * Vc V ^ 2 = u / a P = 2 pi / (Sqrt (u / a ^ 3) ) K = 1/2 V ** 2 - G * M / r (energijos taupymas)

Ekscentriškos orbitos periodas yra toks pat kaip apskritos orbitos, turinčios tą pačią pusiau pagrindinę ašį, periodas.

Greičio pokytis, reikalingas kampo phi plokštumai pakeisti žiedinėje orbitoje:

delta V = 2 kvrt (GM / r) nuodėmė (phi / 2)

Energija masei m įkišti į apskritą orbitą (neatsižvelgiama į sukimosi greitį, kuris šiek tiek sumažina energiją).

GMm (1 / Re - 1 / 2Rcirc) Re = žemės spindulys Rcirc = apskritos orbitos spindulys.

Klasikinė raketinė lygtis, kur dv = greičio pokytis Isp = specifinis variklio impulsas Ve = išmetimo greitis x = reakcijos masė m1 = raketos masė, išskyrus reakcijos masę g = 9,80665 m / s ^ 2

Ve = Isp * g dv = Ve * ln ((m1 + x) / m1) = Ve * ln ((galutinė masė) / (pradinė masė))

Reliatyvistinė raketos lygtis (pastovus pagreitis)

t (be pagreičio) = c / a * sinh (a * t / c) d = c ** 2 / a * (cosh (a * t / c) - 1) v = c * tanh (a * t / c)

Reliatyvistinė raketa su išmetimo greičiu Ve ir masės santykiu MR:

t (be pagreičio) = c / a * sinh (Ve / c * ln (MR)) d = c ** 2 / a * (cosh (Ve / C * ln (MR)) - 1) v = c * tanh ( Ve / C * ln (MR))

Konvertuojant iš paralakso į atstumą:

d (parsekais) = 1 / p (lanko sekundėmis) d (astronominiais vienetais) = 206265 / p

Įvairūs f = ma - jėga yra masė ir pagreitis w = fd - darbas (energija) yra jėga ir atstumas

Atmosferos tankis kinta kaip exp (-mgz / kT), kur z yra aukštis, m - molekulinė masė oro kg, g - vietinis sunkio pagreitis, T - temperatūra, k - Bolztmanno konstanta. Žemėje iki 100 km,

kur d yra tankis, d0 yra tankis 0 km atstumu, yra maždaug teisinga, taigi

d @ 12 km (40000 pėdų) = d0 * .18 d @ 9 km (30000 pėdų) = d0 * .27 d @ 6 km (20000 pėdų) = d0 * .43 d @ 3 km (10000 pėdų) = d0 *. 65

Atmosferos skalės aukštis Sauso praleidimo greitis (km esant emisijos lygiui) (K / km) ------------------------- ------- ------- Žemė 7,5 9,8 Marsas 11 4,4 Venera 4,9 10,5 Titanas 18 1,3 Jupiteris 19 2,0 Saturnas 37 0,7 Uranas 24 0,7 Neptūnas 21 0,8 Tritonas 8 1

Titiuso-Bode dėsnis, skirtas apytiksliam planetų atstumui:

R (n) = 0,4 + 0,3 * 2 ^ N astronominiai vienetai

Tai gana gerai tinka Merkurijui (N = -begalybė), Venerai (N = 0), Žemei (N = 1), Marsui (N = 2), Jupiteriui (N = 4), Saturnui (N = 5), Uranui ( N = 6) ir Plutonas (N = 7).

6.62618e-34 Js (7e-34) - Plancko konstantos „h“ 1.054589e-34 Js (1e-34) - Plancko konstantos / (2 * PI), „h bar“ 1.3807e-23 J / K ( 1.4e-23) - Boltzmanno konstanta „k“ 5.6697e-8 W / m ^ 2 / K (6e-8) - Stephan-Boltzmann konstanta „sigma“ 6.673e-11 N m ^ 2 / kg ^ 2 (7e -11) - Niutono gravitacinė konstanta „G“ 0,0029 m K (3e-3) - Vienos konstanta „sigma (W)“ 3,827e26 W (4e26) - Saulės šviesumas 1370 W / m ^ 2 (1400) - - Saulės pastovumas (intensyvumas esant 1 AG) 6,96e8 m (7e8) - Saulės spindulys 1738 km (2e3) - Mėnulio spindulys 299792458 m / s (3e8) - šviesos greitis vakuume "c" 9,46053e15 m (1e16) - šviesos metai 206264.806 AU (2e5) - 3.2616 šviesmečiai (3) - -> parsek 3.0856e16 m (3e16) - /

Juodosios skylės spindulys (dar vadinamas Schwarzschildo spinduliu):

2GM / c ^ 2, kur G yra Niutono Grav konstanta, M yra BH masė, c yra šviesos greitis

Dalykai, kuriuos reikia pridėti (kažkas jų ieško!) Pagrindiniai raketos skaičiai ir lygtys Aerodinaminiai dalykai Energija, skirta įterpti svarą į orbitą arba pagreitėti iki žvaigždžių greičio. Ne žiediniai atvejai?

APSKAIČIAVIMŲ ATLIKIMAS IR DUOMENŲ FORMATŲ AIKTINIMAS

KOSMOSO ERDVĖS ORBITŲ IR TRAJEKTORIŲ SKAIČIAVIMAS

Tinklalapyje dažnai rekomenduojamos nuorodos yra šios:

„Astrodinamikos pagrindai“ Rogeris Bate'as, Donaldas Muelleris, Jerry'as White'as 1971 m., „Dover Press“, 455 p. 8,95 USD (JAV) (minkštas viršelis). ISBN 0-486-60061-0

NASA kosminių skrydžių vadovai (datuojami 1960 m.) SP-33 orbitinio skrydžio vadovas (3 dalys) SP-34 mėnulio skrydžio vadovas (3 dalys) SP-35 planetinių skrydžių vadovas (9 dalys)

Jų galima rasti universitetų aeronautikos bibliotekose arba užsisakyti per JAV vyriausybę. „Printing Office“ (GPO), nors jiems užsisakyti tikriausiai reikėtų daugiau informacijos.

M. A. Minovitch, _ Balistinių tarpplanetinių trajektorijų nustatymas ir charakteristikos veikiant daugybei planetų atrakcijų_, techninė ataskaita 32-464, reaktyvinių variklių laboratorija, Pasadena, Kalifornija, 1963 m.

Pavadinimas sako viską. Pradedama nuo pagrindų ir tęsiama. Labai gerai. Jame yra papildomas straipsnis:

M. Minovitch, _ Didelių planetų perubacijų panaudojimas projektuojant giluminę kosminę saulės zondą ir ne ekliptikos trajektorijas_, techninė ataskaita 32-849, JPL, Pasadena, Kalifornija, 1965 m.

Pirmiausia turite perskaityti pirmąjį, kad iš tikrųjų suprastumėte šį. Jame yra _trumpa_ santrauka, jei galite rasti tik antrąją.

Norėdami gauti šių ataskaitų, susisiekite su JPL.

„Erdvėlaivio požiūrio dinamika“, Peteris C. Hughesas 1986, Johnas Wiley ir sūnūs.

„Dangaus mechanika: kompiuterio vadovas praktikui“, Lawrence G. Taffas (Wiley-Interscience, Niujorkas, 1985).

Pradeda nuo pagrindų (2 kūno problema, koordinatės) ir veikia nustatant orbitą, trikdžius ir diferencines korekcijas. Taffas taip pat trumpai aptaria žvaigždžių dinamiką, įskaitant trumpą n-kūno problemų aptarimą.

PLANETARINĖS POZICIJOS SKAIČIAVIMAS

"Aiškinamasis astronominio almanacho priedas" (pataisytas leidimas), Kennethas Seidelmannas, University Science Books, 1992. ISBN 0-935702-68-7. Kietais viršeliais 65 USD.

Gili matematika visiems AA algoritmams ir lentelėms.

Suteikia serijas, kad būtų galima apskaičiuoti tikslumą iki 1 lanko minutės + arba - 300 metų laikotarpiui. Plutonas yra įtrauktas, tačiau teigiama, kad jo tikslumas yra tik apie 15 lanko minučių.

_ Daugiametis interaktyvus kompiuterių almanachas (MICA), pagamintas JAV karinio jūrų laivyno observatorijos. Galioja 1990–1999 m. 55 USD (80 USD už JAV ribų). Galima įsigyti „IBM“ (užsakymo Nr. PB93-500163HDV) arba „Macintosh“ (užsakymo nr. PB93-500155HDV). Iš NTIS prekybos centro (703) -487-4650. Manau, kad tai skirta pakeisti USNO interaktyvų kompiuterių efemerį.

_Interaktyvus kompiuterinis efemeris (iš JAV jūrų observatorijos), platinamas IBM ir PC diskeliuose, 35 USD („Willmann-Bell“). Apima 1800-2049 m.

„Planetos programos ir lentelės nuo -4000 iki +2800“, Bretagnon & Simon 1986, Willmann-Bell.

Diskelius galima įsigyti atskirai.

„Dangaus mechanikos pagrindai“ (2-asis leidimas), J.M.A. Danby 1988 m., Willmannas-Bellas.

Geras pagrindinis tekstas. Apima „BASIC“ programas, atskirą diskelių rinkinį galima įsigyti atskirai.

"Astronominės skaičiuoklių formulės" (4 leidimas), J. Meeus 1988, Willmann-Bell.

„Astronominiai algoritmai“, J. Meeus 1991, Willmann-Bell.

Jei aktyviai naudosite vieną iš „Astronominių formulių skaičiuoklėms“ leidimų, norėsite jį pakeisti „Astronomijos algoritmais“. Ši nauja knyga yra labiau orientuota į kompiuterius, o ne į skaičiuotuvus, ir joje pateikiamos planetų judėjimo formulės, pagrįstos šiuolaikiniu „Jet Propulsion Laboratory“, JAV karinio jūrų laivyno observatorijos ir „Bureau des Longitudes“ darbu. Ankstesnės knygos buvo sukurtos remiantis formulėmis, kurios dažniausiai buvo sukurtos praėjusiame amžiuje.

Algoritmai prieinami atskirai diskelyje.

„Praktinė astronomija su jūsų skaičiuokle“ (3-asis leidimas), P. Duffett-Smith, 1988, Kembridžo universiteto leidykla.

„Orbitos mėgėjams su mikrokompiuteriu“, D. Tattersfield 1984, Stanley Thornes, Ltd.

Įtraukia programų pavyzdžius į BASIC.

„Orbitos II mėgėjams“, D. Tattersfieldas 1987 m., Johnas Wiley ir Sonsas.

„Astronomija / mokslinė programinė įranga“ - dalinamųjų programų, viešosios nuosavybės ir komercinės programinės įrangos katalogas, skirtas IBM ir kitiems kompiuteriams. Astronomijos programinė įranga apima planetariumo modeliavimą, efemerių generatorius, astronomines duomenų bazes, saulės sistemos modeliavimą, palydovinio stebėjimo programas, dangaus mechanikos simuliatorius ir kt.

„Andromeda Software, Inc.“ P.O. Box 605 Amherst, NY 14226-0605

KRATERINIŲ DIMETRŲ APSKAIČIAVIMAS IŠ ĮŽEMĮ ĮTAKANČIŲ ASTEROIDŲ

Astrogeologas Gene Shoemaker siūlo šią formulę, pagrįstą branduolinių bandymų sukelto krateriavimo tyrimais.

(1 / 3.4) D = S S c K W: kraterio skersmuo km g p f n

(1/6) S = (g / g): kūnų, išskyrus g e t Žemę, gravitacijos korekcijos koeficientas, kur g = 9,8 m / s ^ 2, o g yra tikslinio kūno paviršiaus e t sunkis. Šis mastelis nurodomas mėnulio krateriams ir gali atitikti kitus kūnus.

(1 / 3.4) S = (p / p): tikslinio tankio p korekcijos koeficientas, pattp = 1,8 g / cm ^ 3 aliumui Jangle U kraterio vietoje, p = 2,6 g / cm ^ 3 vidutiniam uolienai žemyniniai skydai.

C: kraterio žlugimo koeficientas, 1 - kraterių skersmens, 1,3 - didesnių kraterių (Žemėje).

(1 / 3,4) K: 0,74 km / (kT TNT ekvivalentas) n empiriškai nustatytas iš „Jangle U“ branduolinio bandymo kraterio.

3 2 22 W = pi * d * delta * V / (12 * 4,185 * 10): sviedinio kinetinė energija MT TNT ekvivalentu, nurodytu skersmeniu d, greičiu v ir sviedinio tankio delta CGS vienetais. maždaug 3 g / cm ^ 3 delta yra gana gera asteroidui.

Žemę kertantiems asteroidams gali būti naudojamas RMS greitis V = 20 km / sek.

Remiantis šiomis prielaidomis, kūnas, Arizonoje sukūręs „Barringer Meteor Crater“ (1,13 km skersmens), būtų buvęs apie 40 metrų skersmens.

Apskritai galima naudoti (po Gehrels, 1985):

Asteroidas Objektų skaičius Poveikio tikimybė Poveikio energija, kaip Hirosimos bombos skersmuo (km) (smūgiai per metus)

10 10 10^-8 10^9 1 1 000 10^-6 10^6 0.1 100 000 10^-4 10^3

priimant paprastus mastelio dėsnius. Manoma, kad Hirosimos sprogimas yra 0,013 MT TNT ekvivalentas arba maždaug 5 * 10 ^ 13 džaulių.

Gehrels, T. 1985 Asteroidai ir kometos. _Fizika šiandien_ 38, 32–41. [puiki bendra pasauliečių temos apžvalga]

Shoemaker, E. M., 1983 m. Žemės asteroidas ir kometa. _ Met. Šv. Žemės planeta. Sci._ 11, 461-494. [labai ilgas ir gana techniškas, tačiau išsamus dalyko nagrinėjimas]

Batsiuvys, E. M., J. G. Williamsas, E. F. Helinas ir R. F. Wolfe 1979 m. Žemę kertantys asteroidai: orbitos klasės, susidūrimo su Žeme dažnis ir kilmė. Publikacijoje _Asteroids_, T. Gehrels, red., P. 253-282, Arizonos universiteto leidykla, Tuksonas.

Cunningham, C. J., 1988 m. - „Įvadas į asteroidus: kita siena“ (Richmond: Willman-Bell, Inc.) [apima visus asteroidų tyrimų aspektus ir yra puikus įvadas į temą įvairaus lygio žmonėms. Jame taip pat yra labai platus informacinis sąrašas, apimantis iš esmės visą šios srities informacinę medžiagą.]

ŽEMĖLAPIO PROJEKCIJOS IR Sferinė trigonometrija

Du lengvai randami žemėlapių projekcijų šaltiniai yra „Encyclopaedia Britannica“ (ypač senesni leidimai) ir pamoka, rodoma _Graphics Gems_ (Academic Press, 1990). Pastarasis buvo parašytas turint omenyje ekspozicijos paprastumą ir tinkamumą skaitmeniniam skaičiavimui (taip pat rodomos sferinės trig formulės, kaip ir skaitmeniniu būdu suplanuoti pavyzdžiai).

John Snyderio USGS leidinyje „Žemėlapio projekcijos - darbo vadovas“, USGS profesionaliame dokumente 1395. Daugiau nei jums kada nors rūpėjo žinoti apie žemėlapio projekcijas. Jame pateikiami išsamūs 32 projekcijų aprašymai su istorija, ypatybėmis, projekcijų formulėmis (abiem sferinėms žemė ir elipsoidinė žemė) ir skaitiniai bandymo atvejai. Tai tvarkinga knyga, verta visų 382 puslapių. Tai yra 20 USD.

Taip pat galbūt norėsite Snyderio ir Philipo Voxlando bendraautorio tomo „Žemėlapių projekcijų albumas“, USGS Professional Paper 1453. Čia pateikiama mažiau informacijos apie 130 projekcijų ir variantų. Formulės yra gale, pavyzdiniai siužetai priekyje. 14 USD, 250 puslapių.

Galite užsisakyti šiuos 2 būdus. Pigus, lėtas būdas yra tiesiai iš USGS: Žemės mokslo informacijos centro, JAV geologijos tarnybos, 507 nacionalinio centro, Restono, VA 22092. (800) -USA-MAPS. Jie gali jums pasiūlyti kainą ir pasakyti, kur siųsti pinigus. Tikimasi 6–8 savaičių trukmės.

Daug greitesnis būdas (apie 1 savaitę) yra „Timely Discount Topos“, (303) -469-5022, 9769 W. 119th Drive, Suite 9, Broomfield, CO 80021. Paskambinkite jiems ir pasakykite, ko norite. Jie nurodys kainą, jūs atsiųsite čekį, tada jie eis į „USGS“ klientų aptarnavimo skyrių ir pasiims jį jums. Pridėkite maždaug 3–4 USD aptarnavimo mokestį ir pristatymą.

(Gal labiau prieinamas) žemėlapių straipsnis yra:

R. Milleris ir F. Reddy, „Pasaulio žemėlapis Paskalyje“, Baitas V12 # 14, 1987 m. Gruodžio mėn

Turi „Turbo Pascal“ procedūras, skirtas penkioms įprastoms žemėlapio projekcijoms. Demo programą „CARTOG.PAS“ ir nedidelius (6000 taškų) pakrantės duomenis galima rasti „CompuServe“, „GEnie“ ir daugelyje „BBS“.

Kai kurios sferinės trignometrijos nuorodos yra šios:

_Sferinė astronomija_, W.M. Smart, Kembridžo U. Spauda, ​​1931 m.

_ Sferinės astronomijos sąvadas, S. Newcombas, Doveris, 1960 m.

_Sferinė astronomija_, R.M. Green, Cambridge U. Press., 1985 („Smart“ atnaujinimas).

_Sferinė astronomija_, E Woolard ir G.Clemence, „Academic Press“, 1966 m.

VEIKSMINGAI ATLIKTI N-KŪNO MODULIAVIMUS

„Kompiuterinis modeliavimas naudojant daleles“ R. W. Hockney ir J. W. Eastwoodas (Adam Hilger Bristol ir Philadelphia 1988)

„Greitas potencialių laukų įvertinimas dalelių sistemose“, L. Greengard MIT Press, 1988.

Proveržio O (N) modeliavimo metodas. Buvo lygiagreti.

L. Greengard ir V. Rokhlin, „Greitas dalelių modeliavimo algoritmas“, Journal of Computational Physics, 73: 325-348, 1987.

„O (N) algoritmas trimačiams N kūno modeliavimams“, MSEE darbas, Feng Zhao, MIT AILab techninė ataskaita 995, 1987

Apima O (N ^ 2) FORTRAN kodą, kurį parašė šios srities pradininkas Aarsethas.

Hierarchiniai (N log N) medžių metodai aprašyti šiuose dokumentuose:

A. W. Appel, „Efektyvi daugelio kūno modeliavimo programa“, SIAM mokslinių ir statistinių skaičiavimų leidinys, t. 6, p. 85, 1985 m.

Barnes & Hut, „Hierarchinis O (N log N) jėgos skaičiavimo algoritmas“, Gamta, V324 # 6096, 1986 m. Gruodžio 4–10 d.

L. Hernquist, „Hierarchiniai N kūno metodai“, Kompiuterinės fizikos komunikacijos, t. 48, p. 107, 1988 m.

ATSIŽVELGIANT Į „FITS“ Vaizdo formatą

Jei jums tiesiog reikia ištirti FITS vaizdus, ​​naudokite paketą ppm (žr. DUK apie comp.graphics), jei norite konvertuoti juos į pageidaujamą formatą. Norėdami gauti daugiau informacijos apie formatą ir kitą programinę įrangą, kad galėtumėte ją skaityti ir rašyti, skaitykite DUK sci.astro.fits.

TRYS Dimensijų žvaigždės / galaktikos koordinatės

Norėdami sugeneruoti astronominių objektų 3D koordinates, pirmiausia gaukite astronominę duomenų bazę, kurioje nurodomas teisingas objektų pakilimas, deklinacija ir paralaksas. Konvertuokite paralaksą į atstumą naudodami DUK 6 dalyje pateiktą formulę, konvertuokite RA ir deklinaciją į vienetinės sferos koordinates (išsamesnės informacijos apie tai rasite šiame skyriuje anksčiau pateiktose nuorodose apie planetų padėtį ir sferinę trignometriją) ir skalę Atstumas.

Dvi šiam tikslui naudingos duomenų bazės yra „Yale Bright Star“ katalogas (šaltiniai išvardyti DUK 3 skyriuje) arba „Žvaigždžių katalogas per 25 parsekus nuo saulės“.


Šis elementas yra PBL (probleminis mokymasis) užsiėmimas, skirtas įžanginei fizikai, susijusiai su impulso, trinties jėgos ir kinetinės energijos išsaugojimu. Šis scenarijus susijęs su automobilio avarija, kai mažą automobilį plačiai smūgiuoja transporto priemonė, kurios masė viršija dvigubai. Studentai turi užduotį nustatyti, ar kuris nors vairuotojas užsiima neapdairiu vairavimu. Norėdami išspręsti problemą, studentai bendradarbiauja nustatydami trinties koeficientą važiuojamojoje dalyje, kiekvienos transporto priemonės greitį avarijos metu ir transporto priemonės greitį prieš stabdant. Jie naudos strategijų derinį: darbo ir energijos teoremą ir kinematines lygtis.

Patikimai laikydamiesi PBL metodo, studentai persiųs informaciją, kad atskirtų naudingą nuo nesusijusių duomenų, patys suras trūkstamą informaciją ir paskui taikys fiziką ieškodami sprendimų.

Šiame šaltinyje yra atsispausdinamas mokinio vadovas ir slaptažodžiu apsaugotas mokytojo vadovas su sprendimais ir patarimais instruktoriams. (Norint naudotis mokytojo vadovu, reikia užsiregistruoti pas autorius.) ŽR. ŠIAME PUSLAPYJE SUSIJUSIUS PUNKTUS, pateikiančius nuorodą į visą tų pačių autorių PBL pratybų rinkinį.


5 Daugelio polių priėmimo kriterijaus optimizavimas

Turint omenyje patobulintus klaidų įvertinimus, pagaliau galima apsvarstyti praktišką FMM įgyvendinimą siekiant didelio tikslumo. Pagrindiniai šiame kontekste kylantys klausimai yra šie:

Ką pasirinkti plėtros centrams z ir s?

Kada laikyti dvi ląsteles gerai atskirtomis?

Kokia išplėtimo tvarka p naudoti?

Galimi atsakymai į šiuos klausimus turi įtakos skaičiavimo kainai ir apytikslio tikslumui. Taigi tam tikram tikslumo tikslui yra optimalus pasirinkimas visiems šiems parametrams, atsižvelgiant į minimalų procesoriaus laiko (ir atminties) suvartojimą. Šiame skyriuje nagrinėjami šios problemos algoritminiai aspektai, t. Y. Pasirinkimas z ir s ir daugiapolio priėmimo kriterijaus funkcinė forma. Parametrų derinimas (daugiapolio priėmimo kriterijaus ir išplėtimo eiliškumo) siekiant minimalių skaičiavimo pastangų tam tikram tikslumui yra kito skyriaus tema.

Stebėtina, kad šis optimalaus pasirinkimo klausimas z ir s o kelių polių priėmimo kriterijus nebuvo daug ištirtas. Vietoj to, įgyvendinant daugiapolius metodus, dažnai naudojama viena iš dviejų paprastų strategijų. Medžio kodas paprastai naudoja fiksuotą tvarką p ir išsiplėtimas, sutelktas į ląstelių masės centrus, tuo tarpu dvi ląstelės laikomos gerai atskirtomis, jei tenkinamas paprastas geometrinis daugiapolio priėmimo kriterijus (5), kad θ crit kontroliuoja tikslumą.

Kita vertus, naudojant tradicinius FMM, plėtros centrai z ir s abu yra laikomi geometriniais ląstelių centrais, o dvi ląstelės laikomos gerai atskirtomis, kai tik išsiplėtimas suartėja, ir tai atitinka θ crit = 1. Naudojant hierarchinius kubinius tinklelius (vietoj adaptyvaus medžio), tai įgyvendinama sąveikaujant tik su kaimyninėmis ląstelėmis tame pačiame tinklelio lygyje, kurių pagrindinės ląstelės yra kaimynės (pvz., Cheng et al. [1999]). Tikslumą kontroliuoja tik išplėtimo tvarka p.

5.1 Plėtimosi centrų pasirinkimas z ir s

Kiek man žinoma, visi esami FMM diegimai naudoja tą pačią padėtį daugiapoliams ir potencialiems plėtimosi centrams, t. Y. Z = s, kiekvienai ląstelei. Tradiciniam FMM jie yra lygūs geometriniams ląstelių centrams. Tai turi daug galimų sąveikos krypčių r ˆ naudą, ypač kai θ crit = 1, kuriai būtų galima iš anksto apskaičiuoti koeficientus Θ n m (r ˆ). Tačiau šie koeficientai yra apskaičiuojami skrendant dažnai greičiau nei ieškant lentelės. Be to, atsižvelgiant į 4 paveikslą, θ crit = 1 atrodo netinkamas dideliam tikslumui.

Tiesą sakant, apribojimas z = s sumažina laisvę, taigi ir potencialą optimizuoti metodą. Nepaisant to, kai siekiama mažo tikslumo, z = s = z com, ląstelių masės centrai, pasirinkimas turi tam tikrų pranašumų. Pirma, dipoliai išnyksta, o žemos eilės daugiapoliai būna beveik minimalūs. Antra, jei naudojant abipusę algoritmo versiją (kai sąveika A → B ir B → A atliekama vienu metu), skaičiavimo sąnaudos sumažėja, o apytikslės jėgos tiksliai tenkina trečiąjį Niutono dėsnį, ty F ab + F ba = 0 (Dehnen [2002 ]).

Tačiau praktiškai iš tokio tikslaus Newtono dėsnio paklusnumo nėra jokios naudos, nes visas impulsas nėra tiksliai išsaugotas dėl integracijos klaidų, kylančių dėl to, kad dalelės turi atskirus laiko žingsnius. Be to, tokiu atveju nukrypimo nuo tikslaus impulso išsaugojimo laipsnis neatspindi tikrųjų sukauptų jėgos paklaidų. Taikant bendresnį metodą, apytikslės jėgos nukryps nuo idealios F ab + F ba = 0 dydžiu, palyginamu su jų faktinėmis jėgos paklaidomis, o bendro impulso neišsaugojimas šiek tiek rodo sukauptą jėgos paklaidų poveikį ( taip pat žr. B.3 priedą).

5.1.1 Potencialaus plėtros centro pasirinkimas s

4 skirsnio rezultatai, ypač funkcinė E A → B forma (13) lygtyje, siūlo pasirinkti potencialius išsiplėtimo centrus s tokie, kad susidarę kriauklės spinduliai ρ s, taigi ir įvertintos sąveikos paklaidos, būtų minimalios. Taigi, s = z ses, centras mažiausia gaubianti sfera. Rasti mažiausią uždarą sferą rinkiniui n taškų sudėtingumas yra O (n). Tai darant kiekvienai kriaukle esančiai celei, tektų sumokėti visas O (N ln N) sąnaudas ir būtų nepaprastai brangu.

Vietoj to, aš naudoju tikslų apytikslį rezultatą, surasdamas kiekvienai ląstelei mažiausią sferą, apimančią jos dukterų ląstelių sferas. Tai patiria visas O (N) sąnaudas ir yra įgyvendinama per Kompiuterinės geometrijos algoritmų biblioteką (http://www.cgal.org, Fischer ir kt. [2013]), naudojant Matoušek ir kt. Algoritmą. ([1996]).

5.1.2 Daugiapolio išplėtimo centro pasirinkimas z

Kaip jau minėta pirmiau, z = z com nustatymas turi tam tikrų pranašumų mažiems išplėtimo užsakymams p. Tačiau esant dideliems išplėtimo užsakymams, didesnės eilės daugiapoliai tampa vis svarbesni, o tai rodo, kad z = z ses gali būti geresnis pasirinkimas. Norėdamas įvertinti šių metodų privalumus, pakartojau 4 skirsnio eksperimentus abiem metodais ir palyginau gautas didžiausias absoliučios ir santykinės jėgos paklaidas, patirtas tas pats ląstelių → ląstelių sąveika (kuriai du metodai suteikia skirtingą θ).

Radau, kad abiejų metodų klaidos yra labai panašios, kai vidutinis kvadratinis nuokrypis yra ∼ 0,15 dex, tačiau labai mažas vidutinis nuokrypis. Esant p ≲ 8, pastebima tikslesnių jėgų tendencija z = z com, tuo tarpu esant p smaller 8 gaunamos mažesnės klaidos, kai z = z ses. Ši tendencija yra tiesiog pasekmė, kai P k yra mažesnė, kai z = z com nei z = z ses esant mažam k ir didesnis aukštyje k. Tai kartu su patobulintais klaidų įvertinimais (14) taip pat paaiškina, kad (sąveikai A → B) z = z com teikia tikslesnes jėgas, jei ρ z, A & lt ρ s, B, o z = z ses yra linkę būti tiksliau, jei ρ z, A & gt ρ s, B.

5.2 Paprastas FMM įgyvendinimas

Let us first experiment with an implementation that uses the simple multipole-acceptance criterion (5) and a fixed expansion order p. This is the standard choice for the tree code and as such implemented in many gravity solvers used in astrophysics. The computational costs of such an implementation roughly scale as p α / θ crit 3 with α ∼ 2.3 , because the number of interactions increases as θ crit − 3 for large N, while the cost for one is ∝ p 2.3 . Together with the simple error estimate (11), this means that jei one aims each FMM interaction to satisfy δ a / a < ϵ , then the minimum cost for fixed ϵ occurs for

Thus, the optimal opening angle is independent of p. The accuracy is then controlled by the expansion order, requiring p ≳ 16 for ϵ < 10 − 8 (according to Figure 4). The computational costs rise roughly like | ln ϵ | α with decreasing ϵ.

I applied the FMM method with z = z ses , expansion order p = 8 , and θ crit = 0.4 to N = 10 7 equal-mass particles drawn from a Plummer sphere. Figure 7 plots the resulting distributions of absolute (top), relative (middle), and scaled (bottom) acceleration errors. All three distributions are mono-modal, but very wide, much wider than those obtain from GPU-based direct summation (Figure 3). In particular, there are extended tails towards very large relative or scaled errors, containing only ≲ 1% of the particles but reaching up to 1000 times the median error. These tails are due to particles at large radii and, for the relative errors only, also at small radii (see the discussion in Section 3.1).

Acceleration errors for naive FMM. Similar to Figure 3, but for N = 10 7 particles and accelerations obtained by FMM using expansion order p = 8 and multipole-acceptance criterion ρ z , A + ρ s , B < r θ crit (equation (5)) with θ crit = 0.4 . Bins are 0.01 dex wide.

There are two main effects responsible for these properties of the error distributions. First, errors from a single FMM interaction follow a distribution with variance of 1-2 dex. The maximum errors reported in Section 4 only occur for particles near the edges and corners of the sink cell, while most have smaller errors. Moreover, the force errors due to FMM interactions of the same sink cell with source cells in opposing directions tend to partially cancel rather than add up. Both explain why the median errors reported in Figure 7 are much smaller than the maximum relative error incurred by a single cell → cell interaction, which according to Figure 4 is ∼ 10 −4 .

More important is a second effect: the final force errors are not the sum of the relative errors of individual FMM interactions, which are controlled by the simple multipole-acceptance criterion, but of their absolute errors δa. Since, according to equation (11), δ a ∼ θ p M A / r 2 ∼ θ p + 2 M A / 4 ρ z , A 2 , the FMM interactions with cells of large surface density M / ρ z 2 dominate the error budget. In fact, the particles at very large radii have δ a / a ≈ δ a / f ∼ 10 − 4 , exactly as expected from a few FMM interactions with near maximal errors.

5.3 Towards better multipole-acceptance criteria

This discussion suggests that multipole-acceptance criteria which balance the absolute force errors of individual FMM interactions are preferable. When working with the simple estimators (11) or the error bound (10), this leads to critical opening angles which depend on the properties of the interacting cells, such as their mass or surface density.

Such an approach can indeed be made to work (Dehnen [2002]), but the aim here is to go beyond that and use the improved error estimates (14). This results in the multipole-acceptance criteria

with the aim to obtain δ a / a ≲ ϵ and δ a / f ≲ ϵ , respectively.

The black histograms in Figure 8 show the error distributions resulting from these criteria, when the values for a b and f b used in equation (16) have been taken from the direct-summation comparison run. The distributions for δ a / a in the left and δ a / f in both panels are remarkably narrow with a median error ∼ ϵ as targeted, a steep truncation towards large errors, and a maximum error ∼ 10 ϵ . The tail of large δ a / a in the right panel is due to particles at small radii, for which a ≪ f such that criterion (16b) allows large δ a / a .

Acceleration errors for improved FMM. Same as Figure 7 but for the multipole-acceptance criterion (16a) with ϵ = 2 × 10 − 7 (left) or (16b) with ϵ = 10 − 7 (right). The values for a ir f are either taken from the direct-summation run (black), or obtained by low-oder FMM (red, see Section 5.4.2). In all cases the computational effort is similar to that of the FMM run shown in Figure 7 (since a ≤ f criterion (16a) gives a tighter constraint than (16b) and hence requires larger ϵ for the same computational effort).

The difference between these error distributions and those shown in Figure 7 and resulting from the simple geometric multipole-acceptance criterion (5) is remarkable. While the median errors are comparable, the criteria (16) do not produce extended tails of large errors of the quantity controlled ( δ a / a in left and δ a / f in the right panels of Figure 8), and the maximum errors are more than 2 orders of magnitude smaller. What is more, the tails towards small errors have also been somewhat reduced, indicating that the improved criterion avoids overly accurate individual FMM interactions.

This improvement has been achieved without increasing the overall computational effort, but by carefully considering the error contribution from each approximated interaction.

5.4 Practical multipole-acceptance criteria

In a real application one has, of course, no a priori knowledge of a b or f b for any particle and must instead use something else in the multipole-acceptance criteria (16). In some situations, a suitable scale can be gleaned from the properties of the system modelled. For example, if simulating a star cluster of known mass profile M ( r ) and centre x 0 , one may simply use a b ∼ G M ( r b ) r b − 2 with r b = | x b − x 0 | . I now consider other options.

5.4.1 Using accelerations from the previous time step

Employing the accelerations a b from the previous time step in equation (16a) requires no extra computations. However, it means that the gravity solver is not self-contained, but requires some starter to get the initial accelerations.

Also, using information from the previous time step subtly introduces an artificial arrow of time into the simulation, because δ a new < ϵ a old implies δ a new / a new < ϵ a old / a new . Hence, a particle moving in a direction of increasing acceleration has, on average, smaller δ a / a than when moving in the opposite direction, or in reversed time. However, the time integration methods currently employed almost exclusively in N-body simulations of collisional stellar dynamics are irreversible and introduce their own arrow of time. This suggests, that the additional breach of time symmetry by the magnitude (not the direction) of the force error may not be a serious problem in practice. Footnote 4

5.4.2 Estimating a b or f b using low-order FMM

As Section 4 has shown, the error estimate E ˜ A → B used in the multipole-acceptance criteria (16) still has significant uncertainty, and using highly accurate values for a b or f b in equation (16) is unnecessary. Instead, rough estimates should suffice. Such estimates can be obtained via a low-order FMM. This amounts to running the FMM twice: once with a simple multipole-acceptance criterion to obtain rough estimates for a b or f b , and then again using the sophisticated criteria (16) employing the results of the first run.

The acceleration scale f (defined in equation (4)) is similar to the gravitational potential (1), except that its Greens function is | r | − 2 . This implies that it too can be estimated using FMM, albeit not using an explicitly harmonic formulation.

I implemented both options, estimating a arba f via FMM, using the lowest possible order ( p = 0 for f and p = 1 for gravity - recall that a = ∇ Ψ is approximated at one order lower than the potential Ψ) and multipole-acceptance criterion θ < 1 . To this end, I use s = z = z com and a mutual version of the dual tree walk. The resulting estimates for f or a = | a | have rms relative errors of ∼ 15%. The additional computational effort is still much smaller than that of the high-accuracy approximation of gravity itself, though estimating f is faster because it is a scalar rather than a vector and because no square-root needs to be calculated.

The distributions of acceleration errors resulting from using these estimates in equation (16) are shown in red in Figure 8. They are only very slightly worse than those in black, which have been obtained using the exact values of a b and f b in equation (16).


Turinys

Advantages Edit

  • By construction, SPH is a meshfree method, which makes it ideally suited to simulate problems dominated by complex boundary dynamics, like free surface flows, or large boundary displacement.
  • The lack of a mesh significantly simplifies the model implementation and its parallelization, even for many-core architectures. [4][5]
  • SPH can be easily extended to a wide variety of fields, and hybridized with some other models, as discussed in Modelling Physics.
  • As discussed in section on weakly compressible SPH, the method has great conservation features.
  • The computational cost of SPH simulations per number of particles is significantly less than the cost of grid-based simulations per number of cells when the metric of interest is related to fluid density (e.g., the probability density function of density fluctuations). [6] This is the case because in SPH the resolution is put where the matter is.

Limitations Edit

  • Setting boundary conditions in SPH such as inlets and outlets [7] and walls [8] is more difficult than with grid-based methods. In fact, it has been stated that "the treatment of boundary conditions is certainly one of the most difficult technical points of the SPH method". [9] This challenge is partly because in SPH the particles near the boundary change with time. [10] Nonetheless, wall boundary conditions for SPH are available [8][10][11]
  • The computational cost of SPH simulations per number of particles is significantly larger than the cost of grid-based simulations per number of cells when the metric of interest is not (directly) related to density (e.g., the kinetic-energy spectrum). [6] Therefore, overlooking issues of parallel speedup, the simulation of constant-density flows (e.g., external aerodynamics) is more efficient with grid-based methods than with SPH.

Fluid dynamics Edit

Smoothed-particle hydrodynamics is being increasingly used to model fluid motion as well. This is due to several benefits over traditional grid-based techniques. First, SPH guarantees conservation of mass without extra computation since the particles themselves represent mass. Second, SPH computes pressure from weighted contributions of neighboring particles rather than by solving linear systems of equations. Finally, unlike grid-based techniques, which must track fluid boundaries, SPH creates a free surface for two-phase interacting fluids directly since the particles represent the denser fluid (usually water) and empty space represents the lighter fluid (usually air). For these reasons, it is possible to simulate fluid motion using SPH in real time. However, both grid-based and SPH techniques still require the generation of renderable free surface geometry using a polygonization technique such as metaballs and marching cubes, point splatting, or 'carpet' visualization. For gas dynamics it is more appropriate to use the kernel function itself to produce a rendering of gas column density (e.g., as done in the SPLASH visualisation package).

One drawback over grid-based techniques is the need for large numbers of particles to produce simulations of equivalent resolution. In the typical implementation of both uniform grids and SPH particle techniques, many voxels or particles will be used to fill water volumes that are never rendered. However, accuracy can be significantly higher with sophisticated grid-based techniques, especially those coupled with particle methods (such as particle level sets), since it is easier to enforce the incompressibility condition in these systems. SPH for fluid simulation is being used increasingly in real-time animation and games where accuracy is not as critical as interactivity.

Recent work in SPH for fluid simulation has increased performance, accuracy, and areas of application:

  • B. Solenthaler, 2009, develops Predictive-Corrective SPH (PCISPH) to allow for better incompressibility constraints [12]
  • M. Ihmsen et al., 2010, introduce boundary handling and adaptive time-stepping for PCISPH for accurate rigid body interactions [13]
  • K. Bodin et al., 2011, replace the standard equation of state pressure with a density constraint and apply a variational time integrator [14]
  • R. Hoetzlein, 2012, develops efficient GPU-based SPH for large scenes in Fluids v.3 [15]
  • N. Akinci et al., 2012, introduce a versatile boundary handling and two-way SPH-rigid coupling technique that is completely based on hydrodynamic forces the approach is applicable to different types of SPH solvers [16]
  • M. Macklin et al., 2013 simulates incompressible flows inside the Position Based Dynamics framework, for bigger timesteps [17]
  • N. Akinci et al., 2013, introduce a versatile surface tension and two-way fluid-solid adhesion technique that allows simulating a variety of interesting physical effects that are observed in reality [18]
  • J. Kyle and E. Terrell, 2013, apply SPH to Full-Film Lubrication [19]
  • A. Mahdavi and N. Talebbeydokhti, 2015, propose a hybrid algorithm for implementation of solid boundary condition and simulate flow over a sharp crested weir [20]
  • S. Tavakkol et al., 2016, develop curvSPH, which makes the horizontal and vertical size of particles independent and generates uniform mass distribution along curved boundaries [21]
  • W. Kostorz and A. Esmail-Yakas, 2020, propose a general, efficient and simple method for evaluating normalization factors near piecewise-planar boundaries [11]
  • Colagrossi et al., 2019, study flow around a cylinder close to a free-surface and compare with other techniques [1]

Astrophysics Edit

Smoothed-particle hydrodynamics's adaptive resolution, numerical conservation of physically conserved quantities, and ability to simulate phenomena covering many orders of magnitude make it ideal for computations in theoretical astrophysics. [22]

Simulations of galaxy formation, star formation, stellar collisions, [23] supernovae [24] and meteor impacts are some of the wide variety of astrophysical and cosmological uses of this method.

SPH is used to model hydrodynamic flows, including possible effects of gravity. Incorporating other astrophysical processes which may be important, such as radiative transfer and magnetic fields is an active area of research in the astronomical community, and has had some limited success. [25] [26]

Solid mechanics Edit

Libersky and Petschek [27] [28] extended SPH to Solid Mechanics. The main advantage of SPH in this application is the possibility of dealing with larger local distortion than grid-based methods. This feature has been exploited in many applications in Solid Mechanics: metal forming, impact, crack growth, fracture, fragmentation, etc.

Another important advantage of meshfree methods in general, and of SPH in particular, is that mesh dependence problems are naturally avoided given the meshfree nature of the method. In particular, mesh alignment is related to problems involving cracks and it is avoided in SPH due to the isotropic support of the kernel functions. However, classical SPH formulations suffer from tensile instabilities [29] and lack of consistency. [30] Over the past years, different corrections have been introduced to improve the accuracy of the SPH solution, leading to the RKPM by Liu et al. [31] Randles and Libersky [32] and Johnson and Beissel [33] tried to solve the consistency problem in their study of impact phenomena.

Dyka et al. [34] [35] and Randles and Libersky [36] introduced the stress-point integration into SPH and Ted Belytschko et al. [37] showed that the stress-point technique removes the instability due to spurious singular modes, while tensile instabilities can be avoided by using a Lagrangian kernel. Many other recent studies can be found in the literature devoted to improve the convergence of the SPH method.

Recent improvements in understanding the convergence and stability of SPH have allowed for more widespread applications in Solid Mechanics. Other examples of applications and developments of the method include:

  • Metal forming simulations. [38]
  • SPH-based method SPAM (Smoothed Particle Applied Mechanics) for impact fracture in solids by William G. Hoover. [39]
  • Modified SPH (SPH/MLSPH) for fracture and fragmentation. [40]
  • Taylor-SPH (TSPH) for shock wave propagation in solids. [41]
  • Generalized coordinate SPH (GSPH) allocates particles inhomogeneously in the Cartesian coordinate system and arranges them via mapping in a generalized coordinate system in which the particles are aligned at a uniform spacing. [42]

Interpolations Edit

Kernel functions commonly used include the Gaussian function, the quintic spline and the Wendland C 2 > kernel. [44] The latter two kernels are compactly supported (unlike the Gaussian, where there is a small contribution at any finite distance away), with support proportional to h . This has the advantage of saving computational effort by not including the relatively minor contributions from distant particles.

Although the size of the smoothing length can be fixed in both space and time, this does not take advantage of the full power of SPH. By assigning each particle its own smoothing length and allowing it to vary with time, the resolution of a simulation can be made to automatically adapt itself depending on local conditions. For example, in a very dense region where many particles are close together, the smoothing length can be made relatively short, yielding high spatial resolution. Conversely, in low-density regions where individual particles are far apart and the resolution is low, the smoothing length can be increased, optimising the computation for the regions of interest.

Discretization of governing equations Edit

For particles of constant mass, differentiating the interpolated density ρ i > with respect to time yields

The other important equation for a compressible inviscid fluid is the Euler equation for momentum balance:

Similarly to continuity, the task is to define a discrete gradient operator in order to write

which has the property of being skew-adjoint with the divergence operator above, in the sense that

this being a discrete version of the continuum identity

This property leads to nice conservation properties. [45]

Notice also that this choice leads to a symmetric divergence operator and antisymmetric gradient. Although there are several ways of discretizing the pressure gradient in the Euler equations, the above antisymmetric form is the most acknowledged one. It supports strict conservation of linear and angular momentum. This means that a force that is exerted on particle i by particle j equals the one that is exerted on particle j by particle i including the sign change of the effective direction, thanks to the antisymmetry property ∇ W i j = − ∇ W j i =- abla W_> .

Variational principle Edit

The above SPH governing equations can be derived from a least action principle, starting from the Lagrangian of a particle system:

When applied to the above Lagrangian, it gives the following momentum equation:

Pluging the SPH density interpolation and differentiating explicitly ∂ ρ j ∂ r i >>_>>> leads to

Time integration Edit

From the work done in the 80's and 90's on numerical integration of point-like particles in large accelerators, appropriate time integrators have been developed with accurate conservation properties on the long term they are called symplectic integrators. The most popular in the SPH literature is the leapfrog scheme, which reads for each particle i :

Other symplectic integrators exist (see the reference textbook [48] ). It is recommended to use a symplectic (even low-order) scheme instead of a high order non-symplectic scheme, to avoid error accumulation after many iterations.

Integration of density has not been studied extensively (see below for more details).

Symplectic schemes are conservative but explicit, thus their numerical stability requires stability conditions, analogous to the Courant-Friedrichs-Lewy condition (see below).

Boundary techniques Edit

In case the SPH convolution shall be practiced close to a boundary, i.e. closer than s · h , then the integral support is truncated. Indeed, when the convolution is affected by a boundary, the convolution shall be split in 2 integrals,

where B(r) is the compact support ball centered at r , with radius s · h , and Ω(r) denotes the part of the compact support inside the computational domain, Ω ∩ B(r) . Hence, imposing boundary conditions in SPH is completely based on approximating the second integral on the right hand side. The same can be of course applied to the differential operators computation,

Several techniques has been introduced in the past to model boundaries in SPH.

Integral neglect Edit

The most straightforward boundary model is neglecting the integral,

such that just the bulk interactions are taken into account,

This is a popular approach when free-surface is considered in monophase simulations. [49]

The main benefit of this boundary condition is its obvious simplicity. However, several consistency issues shall be considered when this boundary technique is applied. [49] That's in fact a heavy limitation on its potential applications.

Fluid Extension Edit

Probably the most popular methodology, or at least the most traditional one, to impose boundary conditions in SPH, is Fluid Extension technique. Such technique is based on populating the compact support across the boundary with so-called ghost particles, conveniently imposing their field values. [50]

Along this line, the integral neglect methodology can be considered as a particular case of fluid extensions, where the field, A , vanish outside the computational domain.

The main benefit of this methodology is the simplicity, provided that the boundary contribution is computed as part of the bulk interactions. Also, this methodology has been deeply analyzed in the literature. [51] [50] [52]

On the other hand, deploying ghost particles in the truncated domain is not a trivial task, such that modelling complex boundary shapes becomes cumbersome. The 2 most popular approaches to populate the empty domain with ghost particles are Mirrored-Particles [53] and Fixed-Particles. [50]

Boundary Integral Edit

The newest Boundary technique is the Boundary Integral methodology. [54] In this methodology, the empty volume integral is replaced by a surface integral, and a renormalization:

with nj the normal of the generic j-th boundary element. The surface term can be also solved considering a semi-analytic expression. [54]

Hydrodynamics Edit

Weakly compressible approach Edit

Another way to determine the density is based on the SPH smoothing operator itself. Therefore, the density is estimated from the particle distribution utilizing the SPH interpolation. To overcome undesired errors at the free surface through kernel truncation, the density formulation can again be integrated in time. [54]

The weakly compressible SPH in fluid dynamics is based on the discretization of the Navier–Stokes equations or Euler equations for compressible fluids. To close the system, an appropriate equation of state is utilized to link pressure p and density ρ . Generally, the so-called Cole equation [55] (sometimes mistakenly referred to as the "Tait equation") is used in SPH. It reads

Usually the weakly-compressible schemes are affected by a high-frequency spurious noise on the pressure and density fields. [57] This phenomenon is caused by the nonlinear interaction of acoustic waves and by fact that the scheme is explicit in time and centered in space . [58]

Through the years, several techniques have been proposed to get rid of this problem. They can be classified in three different groups:

  1. the schemes that adopt density filters,
  2. the models that add a diffusive term in the continuity equation,
  3. the schemes that employ Riemann solvers to model the particle interaction.
Density filter technique Edit

The schemes of the first group apply a filter directly on the density field to remove the spurious numerical noise. The most used filters are the MLS (Moving Least Squares) and the Shepard filter [57] which can be applied at each time step or every n time steps. The more frequent is the use of the filtering procedure, the more regular density and pressure fields are obtained. On the other hand, this leads to an increase of the computational costs. In long time simulations, the use of the filtering procedure may lead to the disruption of the hydrostatic pressure component and to an inconsistency between the global volume of fluid and the density field. Further, it does not ensure the enforcement of the dynamic free-surface boundary condition.

Diffusive term technique Edit

A different way to smooth out the density and pressure field is to add a diffusive term inside the continuity equation (group 2) :

The first schemes that adopted such an approach were described in Ferrari [59] and in Molteni [56] where the diffusive term was modeled as a Laplacian of the density field. A similar approach was also used in Fatehi and Manzari . [60]

In Antuono et al. [61] a correction to the diffusive term of Molteni [56] was proposed to remove some inconsistencies close to the free-surface. In this case the adopted diffusive term is equivalent to a high-order differential operator on the density field. [62] The scheme is called δ-SPH and preserves all the conservation properties of the SPH without diffusion (e.g., linear and angular momenta, total energy, see [63] ) along with a smooth and regular representation of the density and pressure fields.

In the third group there are those SPH schemes which employ numerical fluxes obtained through Riemann solvers to model the particle interactions. [64] [65] [66]

Riemann solver technique Edit

Since the above discretization is very dissipative a straightforward modification is to apply a limiter to decrease the implicit numerical dissipations introduced by limiting the intermediate pressure by [67]

where the limiter is defined as

Incompressible approach Edit

Viscosity modelling Edit

In general, the description of hydrodynamic flows require a convenient treatment of diffusive processes to model the viscosity in the Navier–Stokes equations. It needs special consideration because it involves the laplacian differential operator. Since the direct computation does not provide satisfactory results, several approaches to model the diffusion have been proposed.

Introduced by Monaghan and Gingold [68] the artificial viscosity was used to deal with high Mach number fluid flows. It reads

The artificial viscosity also has shown to improve the overall stability of general flow simulations. Therefore, it is applied to inviscid problems in the following form

It is possible to not only stabilize inviscid simulations but also to model the physical viscosity by this approach. To do so

For low Reynolds numbers the viscosity model by Morris [69] was proposed.

Additional physics Edit

Multiphase extensions Edit

Astrophysics Edit

Often in astrophysics, one wishes to model self-gravity in addition to pure hydrodynamics. The particle-based nature of SPH makes it ideal to combine with a particle-based gravity solver, for instance tree gravity code, [70] particle mesh, or particle-particle particle-mesh.

Solid mechanics and fluid-structure interaction (FSI) Edit

Total Lagrangian formulation for solid mechanics Edit

To discretize the governing equations of solid dynamics, a correction matrix B 0 ^<0>> [71] [72] is first introduced to reproducing rigid-body rotation as


Conservation of Linear Momentum and Energy

A pendulum is allowed to "crash" into a bar, dramatically altering its motion, but energy is conserved as is evidenced by the return swing.

Requiring little or no skill and the aid of a carpenter's square, one shot sinks two balls into the pockets.

Currently we provide two different versions of this experiment. One version uses balls rolling down a long, inclined track (roughly 2 meters), and the other version uses balls that swing on a pendulum.

Trasa

.

Same as previous except that mass ratio of balls is 1:3 (softball:basketball) leaving basketball dead and softball four times the height.

A tennis ball/basketball combination is dropped to the floor . the tennis ball theoretically bounces to nine times the original drop height.

What it shows:

Demonstration of elastic collisions between metal balls to show conservation of momentum and energy.

How it works:

Newton's Cradle (less affectionately known as Newton's Balls) consists of six rigid balls hanging in a row with bifilar suspension. The balls hang so that they just barely touch their neighbor.

Various initial conditions can be employed. A single ball displaced will collide with the remaining four, sending the ball at the far end off. Same idea for two or three balls. Four balls, and only the first two will stop the center two.

What it shows:

Using conservation of energy, calculate the height from which Barney must jump so that his head just barely kisses the floor at the bottom of his bungee cord jump. Then verify by experiment. Oops . hate when that happens! It turns out that it's not so simple and there are important details that must be taken into account.

How it works:

Barney (the friendly pink dinosaur) is "sandbagged" (with a 5 kg weight, duct-taped around his waist) and suspended from the sky-hook by a 3.1 meter-long (unstretched) spring. The spring constant has been measured.

Use conservation of energy to predict the height the arrow will reach.

What it shows:

When the string of a bow and arrow is pulled from equilibrium, the elastic potential energy in the bow is converted to kinetic energy of the arrow when the string is released. When the arrow.

A toy car rolling down a loop-the-loop track demonstrates the minimum height it must start at to successfully negotiate the loop.

What it shows:

For an object to move in a vertical circle, its velocity must exceed a critical value vc=(Rg) 1/2 , where R is the radius of the circle and g the acceleration due to gravity. This ensures that, at the top of the loop, the centripetal force balances the body's weight. This can be shown using a toy car on a looped track.

How it works:

The car is released from the top of a ramp and runs down a slope towards.

Faith in the conservation of energy is tested by taking the demonstrator's nose to task.

What it shows:

The principle of conservation of energy ensures that a pendulum released at a particular amplitude will not exceed that amplitude on the return swing. A lecturer's faith in their subject is put to the test using a 50lb (22.7kg) iron ball.

How it works:

Technique is very important here. The best method to employ is to stand with your back against the blackboard with your head also touching the board. This ensures that you don't lean forward after release.

What it shows:

Two cars have the same mass and same spring bumper. When given a push and allowed to collide with a wall, one car bounces off with only a small reduction in speed ("elastic" collison) whereas the other car comes nearly to a complere stop ("inelastic" collision).

How it works:

There are two impulse cars made of identical materials and have the same mass. The car that models an elastic collision has all its lead sinkers securely attached to the frame so that they can't move. In contrast, the car that models an inelastic collision has the lead sinkers.

Radio controlled car moves one way while road moves the other.

What it shows:

We tell our students that, when a car drives down the road, the road and the Earth move in the opposite direction, albeit imperceptibly. This demonstration is a realization of that concept, made possible (and perceptible) by the fact that the road is ne attached to the Earth.


Žiūrėti video įrašą: Implementation of TLD Algorithm (Gruodis 2022).