Ajalooliselt on arvutusteadus enamasti piirdunud teadlaste ja doktorantide valdkonnaga. Kuid aastate jooksul - võib-olla suurema tarkvarakogukonna teadmata - meile teaduslikud arvutiprotsessid on vaikselt koostanud avatud lähtekoodiga koostööraamatukogusid, mis tegelevad valdava enamusega raskete raskustega . Tulemuseks on see, et matemaatiliste mudelite rakendamine on nüüd lihtsam kui kunagi varem ja kuigi valdkond ei pruugi olla veel massitarbimiseks valmis, on eduka rakendamise latt drastiliselt langetatud. Uue arvutusliku koodibaasi nullist väljatöötamine on tohutu ettevõtmine, mida mõõdetakse tavaliselt aastates, kuid need avatud lähtekoodiga teaduslikud arvutusprojektid võimaldavad nende arvutusvõimete suhteliselt kiiresti võimendamiseks kasutada hõlpsasti kättesaadavaid näiteid.
Kuna teadusliku arvutamise eesmärk on anda teaduslik ülevaade looduses eksisteerivatest reaalsetest süsteemidest, esindab teadusharu arvutitarkvara lähenemisviisi reaalsuseks muutmise tipptaset. Tegelikku maailma nii täpsuse kui ka eraldusvõimega väga kõrgel tasemel jäljendava tarkvara valmistamiseks tuleb kasutada keerukat diferentsiaalmatemaatikat, mis nõuab teadmisi, mida ülikoolide, rahvuslaborite või ettevõtete teadus- ja arendusosakondade seintest harva leidub. Lisaks sellele tekitavad olulised arvulised väljakutsed millal püüdes kirjeldada reaalse maailma pidevat, lõpmatult väikest struktuuri nullide ja üksikute diskreetse keele kasutamine. Hoolika arvulise teisendamise ammendav pingutus on vajalik algoritmide muutmiseks, mis on nii arvutuslikult jälgitavad kui ka sisukad tulemused. Teisisõnu, teaduslik arvutus on raske töö.
Mulle isiklikult meeldib see FEniCSi projekt , kasutades seda oma lõputöö jaoks, ja näitan oma kallutatust selle õpetuse jaoks meie koodi näite valimisel. (On ka teisi väga kvaliteetseid projekte nagu DUNE mida võiks samuti kasutada.)
FEniCS-i kirjeldatakse ise kui „koostööprojekti automatiseeritud teadusliku arvutamise uuenduslike kontseptsioonide ja tööriistade väljatöötamiseks, pöörates erilist tähelepanu diferentsiaalvõrrandite automatiseeritud lahendamisele piiratud elementide meetoditega.” See on võimas raamatukogu tohutu hulga probleemide ja teaduslike arvutirakenduste lahendamiseks. Selle kaasautorite hulka kuuluvad Simula uurimislabor, Cambridge'i ülikool, Chicago ülikool, Baylori ülikool ja KTH Kuninglik Tehnoloogiainstituut, kes on selle viimase kümnendi jooksul ühiselt hindamatuks ressursiks ehitanud ( vaadake FEniCS-koodisoojust ).
Hämmastav on see, kui palju vaeva on FEniCSi teek meid kaitsnud. Et mõista aimu nende teemade hämmastavast sügavusest ja ulatusest, mida projekt katab vaadata nende avatud lähtekoodiga õpikut , kus 21. peatükis võrreldakse isegi erinevaid lõplike elementide skeeme kokkusurumatute voogude lahendamiseks.
Lava taga on projekt integreerinud meie jaoks suure hulga avatud lähtekoodiga teadusarvutuste raamatukogusid, mis võivad huvi pakkuda või otseselt kasutada. Need hõlmavad mingis kindlas järjekorras projekte, mida FEniCSi projekt kutsub:
See projekti integreeritud välispakettide loend annab meile aimu selle päritud võimalustest. Näiteks võimaldab MPI integreeritud tugi skaleerida kaugtöötajaid arvutusklastri keskkonnas (st see kood töötab super- või sülearvutis).
Samuti on huvitav märkida, et lisaks teaduslikule arvutusele on palju rakendusi, mille jaoks neid projekte saaks kasutada, sealhulgas rahaline modelleerimine, pilditöötlus, optimeerimisprobleemid ja võib-olla isegi videomängud. Näiteks oleks võimalik luua videomäng, mis kasutab mõnda neist avatud lähtekoodiga algoritmidest ja tehnikatest, et lahendada voolav kahemõõtmeline vedelik, näiteks ookeani / jõe hoovused, millega mängija suhtleb (võib-olla proovige ja proovige purjetada vahelduva tuule- ja veevooluga paadiga).
Püüan siin anda maitset sellele, mida arvulise mudeli väljatöötamine hõlmab, näidates, kuidas ühes neist avatud lähtekoodiga teekidest on välja töötatud ja rakendatud põhiline arvutusliku vedeliku dünaamika skeem - antud juhul FEniCSi projekt . FEnICS pakub API-sid nii Pythonis kui ka C ++ -s. Selle näite puhul kasutame nende Pythoni API-d.
Arutleme mõne üsna tehnilise sisu üle, kuid eesmärk on lihtsalt anda maitse sellest, mida selline teadusliku arvutuskoodi väljatöötamine tähendab, ja sellest, kui palju tänapäevased avatud lähtekoodiga tööriistad meie heaks teevad. Loodetavasti aitame selle käigus demüstifitseerida teadusliku arvutamise keerukat maailma. (Pange tähele, et Lisa on ette nähtud üksikasjalik teave kõigi matemaatiliste ja teaduslike aluste kohta neile, keda see üksikasjalik tase huvitab.)
VASTUVÕTU: Nende lugejate jaoks, kellel on teadusliku arvutustarkvara ja rakenduste taust vähe või puudub, võivad selle näite osad tekitada teile järgmise tunde:
klientide hankimine konsultandiks
Kui jah, siis ärge heitke meelt. Peamine väljavõte on siin see, kuivõrd olemasolevad avatud lähtekoodiga projektid võivad paljusid neist ülesannetest oluliselt lihtsustada.
Seda silmas pidades alustame FEnICSi vaatamisest kokkusurumatute Navier-Stokesi demo . See demo modelleerib L-kujulise painde kaudu voolava kokkupressimatu vedeliku, näiteks torustiku, rõhku ja kiirust.
Lingitud demolehel olev kirjeldus annab koodi käivitamiseks vajalike toimingute suurepärase ja kokkuvõtliku seadistuse ning soovitan teil sellega kiiresti tutvuda. Kokkuvõtteks võib öelda, et demo lahendab kiiruse ja rõhu läbi painde kokkusurumatud vooluvõrrandid . Demo käivitab aja jooksul voolava vedeliku lühikese simulatsiooni, animeerides tulemusi vastavalt sellele. See saavutatakse toru ruumi esindava võrgu seadistamise ja toru abil Lõplike elementide meetod kiiruse ja rõhu arvuline lahendamine võrgusilma igas punktis. Seejärel kordame ajas, ajakohastades kiirus- ja rõhuvälju, kasutades uuesti meie käsutuses olevaid võrrandeid.
Demo töötab hästi nagu praegu, kuid muudame seda veidi. Demo kasutab Chorin lõhenes , kuid me kasutame inspireeritud pisut erinevat meetodit Kim ja Moin selle asemel, mis loodetavasti on stabiilsem. See lihtsalt nõuab, et muudaksime konvektiivsete ja viskoossete terminite lähendamiseks kasutatud võrrandit, kuid selleks peame salvestama eelmise ajaetapi kiirusvälja ja lisama värskendusvõrrandile kaks täiendavat mõistet, mis kasutavad seda eelmist teave täpsema arvulise lähenduse saamiseks.
Tehkem see muudatus. Esiteks lisame uue Function
seadistamise vastu. See on objekt, mis tähistab abstraktset matemaatilist funktsiooni, näiteks vektori või skalaarvälja. Nimetame seda un1
see salvestab eelmise kiirusvälja meie funktsiooniruumis
V
.
... # Create functions (three distinct vector fields and a scalar field) un1 = Function(V) # the previous time step's velocity field we are adding u0 = Function(V) # the current velocity field u1 = Function(V) # the next velocity field (what's being solved for) p1 = Function(Q) # the next pressure field (what's being solved for) ...
Järgmisena peame simulatsiooni igas etapis muutma viisi, kuidas „esialgset kiirust” uuendatakse. See väli tähistab ligikaudset kiirust järgmisel etapil, kui rõhku ignoreeritakse (selles punktis pole rõhk veel teada). Siinkohal asendame Chorini split-meetodi uuema Kimi ja Moini astmelise meetodiga. Teisisõnu muudame välja avaldist F1
:
Asenda:
# Tentative velocity field (a first prediction of what the next velocity field is) # for the Chorin style split # F1 = change in the velocity field + # convective term + # diffusive term - # body force term F1 = (1/k)*inner(u - u0, v)*dx + inner(grad(u0)*u0, v)*dx + nu*inner(grad(u), grad(v))*dx - inner(f, v)*dx
Koos:
# Tentative velocity field (a first prediction of what the next velocity field is) # for the Kim and Moin style split # F1 = change in the velocity field + # convective term + # diffusive term - # body force term F1 = (1/k)*inner(u - u0, v)*dx + (3.0/2.0) * inner(grad(u0)*u0, v)*dx - (1.0/2.0) * inner(grad(un1)*un1, v)*dx + (nu/2.0)*inner(grad(u+u0), grad(v))*dx - inner(f, v)*dx
nii, et demo kasutab nüüd kiirendatud vahevälja lahendamiseks meie uuendatud meetodit, kui ta kasutab F1
Lõpuks veenduge, et värskendaksime iga kordusetapi lõpus eelmist kiirusvälja un1
... # Move to next time step un1.assign(u0) # copy the current velocity field into the previous velocity field u0.assign(u1) # copy the next velocity field into the current velocity field ...
Nii et järgmine on meie FEniCS CFD demo täielik kood koos muudatustega:
'''This demo program solves the incompressible Navier-Stokes equations on an L-shaped domain using Kim and Moin's fractional step method.''' # Begin demo from dolfin import * # Print log messages only from the root process in parallel parameters['std_out_all_processes'] = False; # Load mesh from file mesh = Mesh('lshape.xml.gz') # Define function spaces (P2-P1) V = VectorFunctionSpace(mesh, 'Lagrange', 2) Q = FunctionSpace(mesh, 'Lagrange', 1) # Define trial and test functions u = TrialFunction(V) p = TrialFunction(Q) v = TestFunction(V) q = TestFunction(Q) # Set parameter values dt = 0.01 T = 3 nu = 0.01 # Define time-dependent pressure boundary condition p_in = Expression('sin(3.0*t)', t=0.0) # Define boundary conditions noslip = DirichletBC(V, (0, 0), 'on_boundary && (x[0] 0.5 - DOLFIN_EPS))') inflow = DirichletBC(Q, p_in, 'x[1] > 1.0 - DOLFIN_EPS') outflow = DirichletBC(Q, 0, 'x[0] > 1.0 - DOLFIN_EPS') bcu = [noslip] bcp = [inflow, outflow] # Create functions un1 = Function(V) u0 = Function(V) u1 = Function(V) p1 = Function(Q) # Define coefficients k = Constant(dt) f = Constant((0, 0)) # Tentative velocity field (a first prediction of what the next velocity field is) # for the Kim and Moin style split # F1 = change in the velocity field + # convective term + # diffusive term - # body force term F1 = (1/k)*inner(u - u0, v)*dx + (3.0/2.0) * inner(grad(u0)*u0, v)*dx - (1.0/2.0) * inner(grad(un1)*un1, v)*dx + (nu/2.0)*inner(grad(u+u0), grad(v))*dx - inner(f, v)*dx a1 = lhs(F1) L1 = rhs(F1) # Pressure update a2 = inner(grad(p), grad(q))*dx L2 = -(1/k)*div(u1)*q*dx # Velocity update a3 = inner(u, v)*dx L3 = inner(u1, v)*dx - k*inner(grad(p1), v)*dx # Assemble matrices A1 = assemble(a1) A2 = assemble(a2) A3 = assemble(a3) # Use amg preconditioner if available prec = 'amg' if has_krylov_solver_preconditioner('amg') else 'default' # Create files for storing solution ufile = File('results/velocity.pvd') pfile = File('results/pressure.pvd') # Time-stepping t = dt while t Programmi käivitamine näitab voogu küünarnuki ümber. Käivitage teaduslik arvutuskood ise, et seda näha! Lõpliku kaadri ekraanid on esitatud allpool.
Suhtelised rõhud simulatsiooni lõpus olevas paindes, suurusjärgu järgi suurendatud ja värvitud (mittemõõtmelised väärtused):

Suhtelised kiirused simulatsiooni lõpus olevas paindes, kui vektorglüüfid on suurusjärgu järgi muudetud ja värvitud (mittemõõtmelised väärtused).

Nii et oleme teinud olemasoleva demo, mis juhtub rakendama meie omaga väga sarnast skeemi üsna hõlpsalt, ja muutnud seda paremaks lähenduseks, kasutades eelmise ajaetapi teavet.
Siinkohal võite mõelda, et see on tühine muudatus. See oli ja suuresti see ongi point. See avatud lähtekoodiga teaduslik arvutusprojekt on võimaldanud meil kiiresti rakendada muudetud numbrilist mudelit, muutes nelja koodirida. Sellised muudatused võivad suurtes uurimiskoodeksites võtta kuid.
Projektil on palju muid demosid, mida saaks kasutada lähtepunktina. Neid on isegi mitmeid avatud lähtekoodiga rakendused projektile, mis rakendab erinevaid mudeleid.
Järeldus
Teaduslik arvutus ja selle rakendused on tõepoolest keerulised. Sellest ei saa mööda. Kuid nagu nii paljudes valdkondades üha enam paika peab, võib olemasolevate avatud lähtekoodiga tööriistade ja projektide üha kasvav maastik oluliselt lihtsustada seda, mis muidu oleks äärmiselt keeruline ja tüütu programmeerimisülesanne. Ja võib-olla on aeg isegi lähedal, kus teaduslik arvutus muutub piisavalt kättesaadavaks, et leida ennast hõlpsasti kasutamiseks väljaspool teadlaskonda.
LISA: Teaduslikud ja matemaatilised alused
Huviliste jaoks on siin ülal meie arvutusliku vedeliku dünaamika juhendi tehnilised alused. Järgnev on väga kasulik ja ülevaatlik kokkuvõte teemadest, mida tavaliselt käsitletakse kümmekonna kraadiõppe kursuse jooksul. Teema sügavast mõistmisest huvitatud kraadiõppurid ja matemaatilised tüübid võivad leida, et see materjal on üsna köitev.
Vedeliku mehaanika
'Modelleerimine' on üldiselt mingi reaalsüsteemi jada lähendustega lahendamine. Mudel hõlmab sageli pidevaid võrrandeid, mis ei sobi arvuti rakendamiseks, ja seetõttu tuleb seda arvuliste meetoditega veelgi lähendada.
Vedeliku mehaanika puhul alustame seda juhendit põhivõrranditest, Navier-Stokesi võrranditest, ja kasutame seda CFD-skeemi väljatöötamiseks.
The Navier-Stokesi võrrandid on seeria osalised diferentsiaalvõrrandid (PDE) mis kirjeldavad vedeliku voogu väga hästi ja on seega meie lähtepunktiks. Neid saab tuletada a-st läbi visatud massi-, impulss- ja energiasäästuseadustest Reynoldsi transporditeoreem ja kandideerimine Gaussi teoreem ja tuginedes Stoke'i hüpoteesile. Võrrandid nõuavad pidevat eeldust, kus eeldatakse, et meil on piisavalt vedeliku osakesi, et anda statistilisi omadusi nagu temperatuur, tihedus ja kiirus. Lisaks on vaja eeldada lineaarset suhet pinna pingetensori ja deformatsioonitensori vahel, sümmeetriat pingetensoris ja isotroopse vedeliku vahel. Oluline on teada eeldusi, mida selle arengu käigus teeme ja pärime, et saaksime hinnata rakendatavust saadud koodis. Navier-Stokesi võrrandid Einsteini tähistus ilma pikema jututa:
Massi säilitamine:
Hoogu säilitamine:

Energia säästmine:

kus hälbiv stress on:
Ehkki füüsilises maailmas reguleeritakse enamikku vedeliku voogudest väga üldiselt, pole neist otsest erilist kasu. Võrrandite teadaolevaid täpseid lahendusi on suhteliselt vähe ja a 1 000 000 dollarit Millenniumi auhind on olemas kõigile, kes suudavad lahendada olemasolu ja sujuvuse probleem . Oluline osa on see, et meil on lähtepunkt oma mudeli väljatöötamiseks, tehes komplekssuse vähendamiseks hulga eeldusi (need on klassikalise füüsika kõige raskemad võrrandid).
mida tähendab personali suurendamine
Asjade “lihtsuse” hoidmiseks kasutame oma valdkonnapõhiseid teadmisi, et teha vedelikule kokkusurumatu oletus ja eeldame püsivaid temperatuure, nii et energiavõrrandi säilitamine, millest saab soojusvõrrand, ei ole vajalik (lahti ühendatud). Nüüd on meil kaks võrrandit, endiselt PDE-d, kuid oluliselt lihtsamad, lahendades siiski suure hulga reaalseid probleeme.
Järjepidevuse võrrand
Hoogvõrrandid

Siinkohal on meil nüüd kena matemaatiline mudel kokkusurumatute vedelike voogude jaoks (näiteks väikese kiirusega gaasid ja vedelikud, näiteks vesi). Nende võrrandite otsene käsitsi lahendamine pole lihtne, kuid see on tore, kuna saame lihtsate probleemide jaoks 'täpseid' lahendusi. Nende võrrandite kasutamine huvipakkuvate probleemide lahendamiseks, näiteks tiiva kohal voolav õhk või mõnda süsteemi läbiv vesi, eeldab, et lahendame need võrrandid arvuliselt.
erinevus s corp c corp
Numbrilise skeemi ülesehitamine
Keerukamate probleemide lahendamiseks arvuti abil on vaja meetodit meie kokkusurumatute võrrandite numbriliseks lahendamiseks. Osaliste diferentsiaalvõrrandite või isegi diferentsiaalvõrrandite arvuline lahendamine pole triviaalne. Kuid meie selles juhendis toodud võrranditel on eriline väljakutse (üllatus!). See tähendab, et peame lahendama impulsivõrrandid, hoides lahuse lahknemise vabana, nagu järjepidevus nõuab. Lihtne ajaintegratsioon millegi sarnase kaudu Runge-Kutta meetod on raskendatud, kuna järjepidevuse võrrandil pole ajalist tuletist.
Võrrandite lahendamiseks pole õiget või isegi parimat meetodit, kuid toimivaid võimalusi on palju. Aastakümnete jooksul on selle probleemi lahendamiseks leitud mitmeid lähenemisviise, näiteks ümbersõnastamine keerise ja voogfunktsiooni osas, kunstliku kokkusurutatavuse sisseviimine ja operaatori lõhestamine. Chorin (1969) ja seejärel Kim ja Moin (1984, 1990) sõnastasid väga eduka ja populaarse murdeastmelise meetodi, mis võimaldab meil võrrandeid integreerida, lahendades survevälja otse, mitte kaudselt. Murdosa astme meetod on üldine meetod võrrandite lähendamiseks nende operaatorite jagamise teel, antud juhul piki rõhku jagades. Lähenemine on suhteliselt lihtne ja samas jõuline, motiveerides selle valikut siin.
Esiteks peame võrrandid numbriliselt diskrimineerima ajas, et saaksime astuda ühest ajahetkest teise. Kim ja Moin (1984) järgimisel kasutame otsesõnu teist järku Adams-Bashforthi meetod konvektiivterminite puhul kaudne implitsiit Vända-Nicholsoni meetod viskoossete terminite korral ajaline tuletise jaoks lihtne lõplik erinevus, jättes samas rõhugradiendi tähelepanuta. Need valikud ei ole kaugeltki ainsad lähendused, mida saab teha: nende valik on osa skeemi ülesehitamisest, kontrollides skeemi numbrilist käitumist.

Vahekiirust saab nüüd integreerida, kuid see ignoreerib rõhu panust ja on nüüd lahknev (kokkusurumatus nõuab, et see oleks lahus vaba). Ülejäänud operaator on vajalik, et viia meid järgmise aja sammuni.

kus
on mingi skalaar, mille peame leidma, mille tulemuseks on erinev vaba liikumiskiirus. Me võime leida
korrigeerimisetapi lahknemise abil,

kus esimene termin on järjepidevuse järgi null, saades a Poissoni võrrand skalaarvälja jaoks, mis tagab solenoidse (divergentse vaba) kiiruse järgmisel ajahetkel.

Nagu näitavad Kim ja Moin (1984),
ei ole operaatori lõhenemise tagajärjel täpselt rõhk, kuid selle saab leida
.
Selles õpetuse punktis läheb meil üsna hästi, oleme diskrimineerinud valitsevad võrrandid ajaliselt, et saaksime neid integreerida. Nüüd peame operaatorid ruumiliselt diskrimineerima. Selle saavutamiseks on olemas mitmeid meetodeid Lõplike elementide meetod , Lõpliku helitugevuse meetod , ja Lõpliku erinevuse meetod , näiteks. Kimi ja Moini (1984) originaalloomingus jätkavad nad lõpliku erinevuse meetodit. Meetod on soodne oma suhtelise lihtsuse ja arvutustõhususe poolest, kuid kannatab keeruka geomeetria korral, kuna see nõuab a struktureeritud võrk .
Lõplike elementide meetod (FEM) on mugav valik selle üldisuse jaoks ja selle kasutamisel on mõned väga toredad avatud lähtekoodiga projektid. Eelkõige tegeleb see reaalsete geomeetriatega ühes, kahes ja kolmes mõõtmes, skaalad masinaklastrite väga suurte probleemide jaoks ja seda on suhteliselt lihtne kasutada kõrgekvaliteediliste elementide jaoks. Tavaliselt on meetod neist kolmest aeglasem, kuid see annab meile probleemide osas kõige suurema läbisõidu, seega kasutame seda siin.
Isegi FEM-i rakendamisel on palju valikuid. Siin me kasutame Galerkin FEM . Seejuures valeme võrrandid kaalutud jääkkujundisse, korrutades need igaüks testi funktsiooniga
vektorite jaoks ja
skalaarvälja jaoks ja domeeni kaudu integreerimine
. Seejärel teostame osalise integreerimise kõigi kõrgekvaliteediliste tuletisinstrumentide abil Stoke'i teoreem või Divergentsiteoreem . Seejärel püstitame variatsiooniprobleemi, saades soovitud CFD-skeemi.



Nüüd on meil rakendamiseks kena matemaatiline skeem 'mugavas' vormis, loodetavasti on see mõnevõrra mõistetav, mida sinna jõudmiseks vaja oli (palju matemaatikat ja suurepäraste uurijate meetodeid, mida me üsna palju kopeerime ja näpistame).