Logiciels libres pour le calcul scientifique : pourquoi et comment ? | LIEGE CREATIVE, 24.11.16

Christophe Geuzaine, Professeur ordinaire, ULg - Sciences Appliquées

Jeudi 24 novembre

Logiciels libres pour le calcul scientifique : pourquoi et comment ?

LIEGE CREATIVE, en partenariat avec :

FreeSoftwareforScientificComputing:WhyandHow?

ChristopheGeuzaineUniversityofLiège– MontefioreInstitute

LiègeCreative,November242016

Casestudies:Gmsh &ONELAB

http://gmsh.info |http://onelab.info

2

SomeBackground

• IamafullprofessorattheUniversityofLiège,whereIleadateamofabout20peopleattheintersectionofappliedmathematics,scientificcomputingandengineeringphysics

• Ourresearchinterests:modeling,analysis,algorithmdevelopment,andsimulationforproblemsarisinginvariousareasofengineeringandscience

• Currently:electromagnetics,biomedicalproblems,geophysicalwaves

3

SomeBackground

• Wewritealotofscientificcomputingsoftware,mostlysolversforpartialdifferentialequations

• Twocodesreleasedasfreesoftware:• Gmsh:meshgeneratorwithbuilt-inCADengineandpost-processor(http://gmsh.info)

• GetDP:generalfiniteelementsolverusingmixedfiniteelements(http://getdp.info)

• Thesearelongtermefforts… celebrating20yearsin2016-2017

• RecentlycombinedwithONELAB:alightweightinterfacetoFEsoftware (http://onelab.info)

FreeSoftwareforScientificComputing:WhyandHow?

5

WhatisFreeSoftware?

• “Freesoftware”meanssoftwarethatrespectsusers'freedomandcommunity

• Roughly,itmeansthattheusershavethefreedomtorun,copy,distribute,study,changeandimprovethesoftware

• Thus,“freesoftware”isamatterofliberty,notprice:thinkof“free”asin“freespeech,”notasin“freebeer”

6

WhatisFreeSoftware?

Aprogramisfreesoftwareif,asauser,youhave1. thefreedomtoruntheprogramasyouwish,for

anypurpose2. thefreedomtostudyhowtheprogramworks,and

changeitsoitdoesyourcomputingasyouwish3. thefreedomtoredistributecopiessoyoucanhelp

yourneighbor4. thefreedomtodistributecopiesofyourmodified

versionstoothers.Accesstothesourcecode(“opensource”)isapreconditionfor2and4.

7

WhatisScientificComputing?

Scientificcomputing(or“computationalscience”)=multidisciplinaryfieldthatusesadvancedcomputingcapabilitiestounderstandandsolvecomplexproblems• Appliedmathematics:gettingresultsoutof

applicationareas(science,engineering,humanities)• Numericalanalysis:resultsneedtobecorrectlyand

efficientlycomputable• Computing:thealgorithmsneedtobeimplemented

onmodernhardware

8

WhatisScientificComputing?

• Increasingdemandforscientificcomputing• Industryandacademia• Smart+X

• Coupledwithanexplosionofavailablecomputationalpowerandstorage(massivelydistributed/parallel/cloud)

• Numerousunderlyingmathematical/numericalconcepts/methodsfor• Analytics• Modelling&simulation• Optimization&control

9

WhyUse FreeSoftwareforScientificComputing?

• It’sfree!• Youcanuseiton10,000CPUsifyouwant

• It’sopen!• Youcanstudyallthedetails,modifywhatyouwant(ifyouhavetheskills),anduseitforever

• Qualityhasincreaseddramaticallyoverlast15years• thankstotheinternetandtotheavailabilityofhigh-quality,freeoperatingsystemsandtools(Linux,GCC,MPI,Python,Git,…)

• It’showscienceworks:• “standingontheshouldersofgiants”– reuse&extendexistingsoftwareforfasterandbetterresults

10

WhyUse FreeSoftwareforScientificComputing?

• Asaresult,freeprofessionalqualitytoolsforscientificcomputingarenowavailable• Examplesforengineering:OpenFOAM,Salome,R,Octave,deal.II,Paraview,Gmsh,…

• Similarforscienceorhumanitiesproblems• Downsides?• Lackofdocumentation/support• Heterogeneityofthetools• Lackof“expert”layer/userinterface• Usuallytakesmoretimeandefforttogetfirstresults• Thisiswheregoodcommercialtoolsshine

11

WhyDevelop FreeSoftwareforScientificComputing?

• Philosophicalreasons:• Openscience(opendata,opensource,openmethodology,openpeer-review,openaccess,openeducationalresources)forresponsibleandreproducibleresearch

• Publicresearchfunding:shouldgivebacktosociety• Practicalreasons:• Increasevisibility,numberofusers• Allowseasyco-operationwithfellowscientists/engineers

• Leadstobettercode

12

How?

• Asadeveloper,youchoosealicenseforyourcode:• Basedoncopyright(sometimeswithpatentclauses)• Variouschoices,frommorepermissive(“BSD”-type)tomorerestrictive(“GNUGPL”-type:“copyleft”)

• Differencesbetweenlicensescanbequitesubtle• Andyouimplementadevelopmentprocess,usually

tofosterco-operation(todaye.g.http://github.com)• Sourcecoderepositorywithversiontracking• Continuousintegrationtools• Bugtracking• Forums,mailinglists,wikis,…

13

IsitIncompatiblewithMakingMoney?

No

Butitdoesmakeitdifficulttomakealot ofmoney:• Sellservice/support;smallmargins,butstimulates

(high-quality)jobcreation• Addclosedlayersontop,withclearaddedvalue

(e.g.expertlayer)• Multiplelicensing(e.g.forintegrationintoclosed-

sourcesoftware)

Casestudy:Gmsh

15

Gmsh

Anopensourcefiniteelementmeshgenerator,withbuilt-inCADandpost-processingfacilities

Let’shavealook

16

Gmsh

Interactivedemoshowingtheuserinterfaceandtheconstructionofmeshesusing.geoand.stepfiles

17

Gmsh – ABitofHistory

• Startedin1996(“after-hours”,withJ.-F.Remacle)• Firstpublic(binary)releasein1998• OpensourcedunderGNUGPLin2003• Gmsh 2.0in2006• IJNMEpaperandswitchtoCMake in2009• Pythonwrappersin2011• Curvilinearmeshing:2012–• Quad-Hexmeshing:2012–• Homology:2013–• Parallelmeshing:2015–

18


• Funding:• Projectstartedasahobby(to“scratchanitch”),nosupportandnofundinguntil2006

• EDFfundedOpenCASCADE integrationin2006• Morerecently:• EUandRWprojects• SeveralPhDtheses&post-docs• Industrialusers• Commercialsoftwarevendors

• Noneoftherecentdevelopmentssignificantlyaffectthekernel

19


• Today,Gmsh represents• abouthalfamillionlinesofC++code• stillonly2coredevelopers(about100withrepositorywriteaccess)

• about1000peopleonmailinglists• about5000downloadsperweek(80%Windows)• about400citationsperyear– theGmsh paperiscitedabout2000times

• Gmsh oneofthemostpopularopensourcefiniteelementmeshgeneratorintheworld

2000 2005 2010 2015Year

0

100

200

300

400

500

600

700

References

"Gmsh" keyword on Google Scholar

20

Gmsh – StrategicChoices

• Designgoals:fast,lightanduser-friendly• Writtenin(simple)C++• GUIs:FLTK(desktop),UIKit (iOS),Android,web(soon!)• OpenGLgraphics• Highlyportable(OSes&compilers)• Easytodistribute:zerodependenciesoninstallation

• Handlingofnumerous(>40)thirdpartylibraries• BuildsystembasedonCMake – everythingisoptional(e.g.-DENABLE_OCC=0)

• Somelibsintegratedandredistributeddirectlyingmsh/contrib (BAMG,Metis4,Taucs,Concorde,…)

21


• Licensing• Gmsh isdistributedundertheGNUGeneralPublicLicenseV2orlater,withexceptionstoallowforeasierlinkingwithexternallibraries

• SpecificlicenseexemptionsgiventoGmsh bydevelopersofMetis,Concorde,BAMG,Tetgen/BR,…

• Somethirdpartylibrariesareunfortunatelynotcompatiblewiththeopensourceversion

• Wedouble-licenseforcommercial“embedders”

22


• Communityinfrastructure• Ourown(notonGitHub)toallowpublic/privateparts,forbothcommercialclientsorforspecificthirdpartyinterfaces(e.g.Parasolid)

• ContinuousintegrationandautomatednightlybuildsonWindows,MacandLinux(usingCMake/CTest/CDash onvirtualmachines)

• Timelineofchanges,wikiandbugtracking(usingTrac)• Twomailinglistsandawebsite(http://gmsh.info)withcompleteuserdocumentationandusertutorial

• Scientificaspectsofalgorithmsdetailedinjournalpapers

23

Gmsh – BasicConcepts

• Gmsh canbeusedat3levels• User:throughtheGUI(whichtranslatesmostactionsintoscriptingcommandsin“.geo”textfiles)

• Advanceduser:throughthededicated“.geo”language

• Developer:throughtheC++orPythonAPI• Maincharacteristic:allalgorithmsarewrittenin

termsofabstractCADentities,usinga“BoundaryREPresentation”(BREP)approach

24


• Gmsh interfacesdirectlywithmultipleCADkernels(Gmsh “.geo”,OpenCASCADE,ACIS,…)

• Gmsh nevertranslatesfromonekerneltoanother:• Read/writedirectlytoeachCADkernelusingownAPI• AllalgorithmswritteninCAD-agnosticmanner• Allowsmulti-kernelgeometries

25


Thegoalistodealwithverydifferentunderlyingdatarepresentationsinatransparentmanner

6

COMPUMAG–SYDNEY 2011, STATIC FIELDS AND QUASI-STATIC FIELDS (III), PA10.2, CMP–305 2

Fig. 1. Mesh used for the numerical experiments [20]: a) full head (300,000 nodes, 27 tissues), b) grey matter, c) white matter.

matter �G

(!) are modeled within a probabilistic framework,as functions of the random variable !. Therefore j

Avg

, jmax

,j99�perc

, eAvg

, emax

and e99�perc

are random as well. Inparticular, by using the maximum entropy principle [22] wemodel (arbitrarily) �

G

(!) and �W

(!) as independent randomvariables, uniformly distributed:

�G

(!) ⇠ U([0.0753 ; 0.5155]) (S/m) (1)�W

(!) ⇠ U([0.0533 ; 0.3020]) (S/m) (2)

C. The non intrusive approach

As the conductivities of the brain and the cerebellum aretwo independent random variables of finite variance, we canexpand them as a truncated series of order p

in

in the bi-dimensional Hermite polynomials of a random gaussian vector⇠(!) = (⇠1(!), ⇠2(!)), known as Hermite chaos polynomials[18]:

�G

(!) ⇡P

inX

i=0

�Gi

i

(⇠(!)) (3)

�W

(!) ⇡P

inX

i=0

�Wi

i

(⇠(!)) (4)

where �Gi

and �Wi

are scalar values that depend on theprobabilistic law of the conductivities, P

in

= Cp

in

2+p

in

is thenumber of bi-dimensional polynomials of order less than p

in

,and

i

is the ith bi-dimensional Hermite polynomial. Tosolve the stochastic problem, we use an approach based ona polynomial chaos decomposition of both the conductivityand the induced fields [18]. We assume the conductivities tobe of finite variance, with no assumption on the shape of theprobabilistic distribution.

The values of the induced fields—the average currentdensity in the brain j

Avg

(!) = jAvg

(⇠(!))—are computedby the finite element method from any couple of values(�

G

(⇠(!)),�W

(⇠(!))). The average density belongs to aspace that can be spanned by the polynomials (⇠(!)) andthus written as a truncated series to an order p

out

:

jAvg

(!) =P

outX

m=0

jAvg

m

m

(⇠(!)). (5)

To compute the value of the unknown real coefficients jAvg

m

,we use the orthogonality properties of the Hermite polynomi-als:

jAvg

m

=E[j

Avg

(!) m

(⇠(!))]

E[ m

(⇠(!))2], (6)

where E[·] is the mathematical expectation. The denominatorcan be computed analytically. The integral in the numeratoris computed by means of a Hermite Gauss integration schemewith d integration points [18]:

E[jAvg

(!) m

(⇠(!))] ⇡dX

i=1

...dX

j=1

wi,j

(jAvg

((t1, t2)i,j)) m

((t1, t2)i,j), (7)

with (t1, t2)i,j the i, j-th Gauss point and wi,j

the associatedweight in the bi-dimensional Cartesian rule. The deterministicproblem must thus be computed d2 times, with the conductiv-ity evaluated through (3) and (⇠1(!), ⇠2(!)) = (t1,= t2)i,j ,i, j = 1, . . . , d.

III. RESULTS AND DISCUSSION

The non intrusive method is governed by three parame-ters: p

in

, pout

and d; pin

is linked to the precision on the

26

Gmsh – RecentDevelopments

• Reparametrization andquad-meshing:2012–• Curvilinearmeshing:2012–• Homologysolver:2013–• Hexmeshing:2014–• Parallelmeshing:2015–

27

Gmsh – Reparametrization andQuad-MeshingApplied'&'Computa/onal'Electromagne/cs(ACE)

C.'Geuzaine';'SysMod'seminar';'Nov'18th'2011 59

Quad&meshing&of&FalconFalcon

Applied'&'Computa/onal'Electromagne/cs(ACE)

C.'Geuzaine';'SysMod'seminar';'Nov'18th'2011 58

Quad&meshing&of&FalconFalcon

Parametrizations of the surfaces of the Falcon aircraft in the {u,v}plane.

E.Marchandise,J.-F.Remacle andC.Geuzaine.Optimalparametrizationsforsurfaceremeshing.EngineeringwithComputers,December2012,pp.1-20J.-F.Remacle,J.Lambrechts,B.Seny,E.Marchandise,A.Johnen andC.Geuzaine.Blossom-Quad:anon-uniformquadrilateralmeshgeneratorusingaminimumcostperfectmatchingalgorithm.InternationalJournalforNumericalMethodsinEngineering89,pp.1102-1119,2012

28

Gmsh – CurvilinearMeshing

A.Johnen,J.-F.Remacle andC.Geuzaine.Geometricvalidityofcurvilinearfiniteelements.JournalofComputationalPhysics233,pp.359-372,2013A.Johnen,J.-F.Remacle andC.Geuzaine.Geometricvalidityofhigh-Ordertriangularfiniteelements.EngineeringwithComputers30(3),pp.375-382,2014T.Toulorge,C.Geuzaine,J.-F.Remacle,J.Lambrechts.Robustuntanglingofcurvilinearmeshes.JournalofComputationalPhysics254,pp.8-26,2013

Curvilinear meshing

References

A. Johnen, J.-F. Remacle and C. Geuzaine.Geometric validity of curvilinear finite elements. Journal of Computational Physics, 233, pp.359-372, 2013.

T. Toulorge, C. Geuzaine, J.-F. Remacle, J. Lambrechts.

Robust untangling of curvilinear meshes. Journal of Computational Physics 254, pp. 8-26, 2013.

29

Gmsh – (Co)homologySolver

M.Pellikka,S.Suuriniemi,L.Kettunen andC.Geuzaine.Homologyandcohomologycomputationinfiniteelementmodeling.SIAMJournalonScientificComputing35(5),pp.1195-1214,2013

14 M. Pellikka, S. Suuriniemi, L. Kettunen, and C. Geuzaine

Fig. 7.1: Induction heating geometry and surface mesh of Mc.

Fig. 7.2: A representation of the basis of the homology space H1(Ma).

Fig. 7.3: A representation of the basis of the cohomology space H1(Ma).









30

Gmsh – Hex-DominantMeshing

T.C.Baudouin,J.-F.Remacle,E.Marchandise,F.Henrotte,C.Geuzaine.Afrontalapproachtohex-dominantmeshgeneration.AdvancedModelingandSimulationinEngineering2014.

Hex-dominant meshing

Name #Vertices%Hex %Prism % Pyr % Tet cpu(s)CUBO 133,436 89.74 4.02 4,20 2.02 247BLADES 133,678 83.65 5.62 6.75 3.98 268CV745 102,946 78.55 7.63 8.81 5.91 225

Baudouin, T. C., Remacle, J. F., Marchandise, E., Henrotte, F., Geuzaine, C. (2014). A frontal

approach to hex-dominant mesh generation. Advanced Modeling and Simulation in Engineering

Sciences, 1(1), 1.

31

Gmsh – ParallelMeshing

J.-F.Remacle,V.Bertrand,C.Geuzaine.ATwo-LevelMultithreadedDelaunayKernel.24thInternationalMeshingRoundtable,ProcediaEngineeringVolume124,pp.6-17,2015

Multi-threaded meshing

104 105 106 107 108

2

4

6

8

10

12

14

16

18

20

Number of Points N

Spee

dup

M=2M=4M=8M=16M=32

Student Version of MATLAB

Less than 100 GB were necessary for 109 tetsIt took 142.8 sec. for generating a 109 mesh (M= 32,M

2

= 8)This corresponds to a rate of ' 7⇥106 tets/sec.

32

Gmsh – WastheMovetoFreeSoftwareUseful?

• Yes!• (Alot)moreusers,leadingtomorerobustcode• Dramaticallyimprovedinternationalvisibility• Newcollaborations,inunexpectedareas(oilandgas,neuroscience,puremath,…)

• Allowedusto“changedivision”(i.e.,reachworldclass)

• Downsides?• Veryfewsignificant externalcodecontributions(topic/codeisjusttoocomplicated)

• Lostsomecommercialopportunities

Casestudy:ONELAB

34

ONELAB

Free,lightweightinterfacetofiniteelementsoftware

Let’shavealook

35

ONELAB

Interactivedemo,openingaGetDP .profile,doingacalculation(e.g.benchmarks/superconductors),andshowingsomepost-processing

whilethesimulationisrunning.

36

ONELAB– Context

Economic• Growingimportanceofscientificcomputingineducationandindustry

• Prohibitivecostofcommercialpackagesforasignificantsubsetofpotentialusers(SMEs,education,occasionaluse)

Scientific• Highqualityoffree/open-sourcesoftwaredevelopedinuniversitiesandresearchcenters

• Sometimesaheadofcommercialequivalents

Practical• Nouser-friendlyinterfaceand/orpoordocumentationformostopensourceFiniteElementAnalysis(FEA)codes

37

ONELAB– Objective

DevelopaplatformforintegratingfreeFiniteElementAnalysis(FEA)software:• allowingtheintegration(byco-simultation)ofany

open-sourcecode,whatevertheircharacteristics• withanintuitiveGUIallowingnoviceuserstoget

startedandguidedintotheworldofFEmodeling—butwiththepossibilitytoconstructsophisticated,upgradable,multi-code,multi-platformscriptsforthespecializeduser

• andwiththepossibilitytoconstructbotheducation-andbusiness-specifictools

38

ONELAB– Design

• Don’treimplement,interfacetheexisting!• Makeitassmall,lightweightandaseasytomaintain

aspossible(nosolver-dependentcodeintheinterface)

• Makeiteasytoprovidetemplates,withinteractiveparametermodification

• ONELABrole=datacentralization,modificationandredispatching

Issuesofcompletenessandconsistencyoftheparametersetarecompletelydealtwithonthesolverside

39

ONELAB– Design

HeterogeneityfreeFEcodes?• Abstractinterfaceforparameterexchange(CAD,

mesh,physicalproperties,constraints,codedriver)

Missing“expert”layer,top-downvalidationanddocumentation?• Developmentanddocumentationoftemplates:• Model:blackbox,parameterizable viaabstractinterface

• Meta-model:setofmodels+selectionlogic

40

ONELAB– Design

Client-server:• Clients:CADkernels,meshers,solvers,post-

processors• Server:Gmsh (currently)+database

Abstractinterface:• Theserverhasnoaprioriknowledgeoftheclients

(nometa-languageorexchangefileformat)• Theserverdoesnotwriteinputfilesfor(native)

clients:theclientcommunicateswiththeservertodefinewhatinformationshouldbeexchanged

41

ONELAB– Implementation

Abstractinterfacetophysicalproperties,constraints&codedrivers:• ReferenceserverinC++,clientsinC++orPython• ExchangeofparametersthroughTCP/IP,Unixsocketsorin-memory

“Native”clientsuseC++orPythondirectly• Gmsh,GetDP,Octave,Fenics,BEM++,…

“Non-native”clientsusePython,byinstrumentingtheirinputfiles• Elmer,OpenFOAM,Code_Aster,Abaqus,Gnuplot,FreeFem,…

42

ONELAB– GiveitaTry

• AvailableonLinux,Mac,Windows,iOS(iPhone/iPad)andAndroid

• Downloadfromhttp://onelab.info,or searchforONELABontheAppStore ortheGooglePlayStore

43

Conclusion

• Freesoftwareforscientificcomputing• Growingneed/useinbothindustryandacademia• Complementarytoclosed/commercialtools

• Mypersonal(academic)experience:freesoftwareisthewaytogo• Philosophicallyandscientifically(openscience)• Pragmatically(bettercode,increasedvisibility)• MakingmycodesfreeclearlyledtomanyopportunitiesthatIwouldn’thavehadotherwise

• ...andcommercialendeavorsarestillpossible

44

[email protected]