Transcription
Open-source CFD software: FeatFlowThe FEA(S)T groupsInstitute of Applied Mathematics, LS IIIDortmund University of Technology, Germanymatthias.moeller@math.tu-dortmund.deLake Tahoe, January 7, 20091
Overview12345FeatFlowFeatFlow 1.xFeatFlow 2.0PreprocessingDeViSoR Grid3DExample applicationPoisson equationPostprocessingGeneral Mesh ViewerUnConventional HPCFEAST2
FeatFlow 1.xFinite Element Analysis Toolbox Flow solverHigh performance unstructured finite element package for thenumerical solution of the incompressible Navier-Stokes equationsbased on the finite element packages Feat2D and Feat3Dwritten in Fortran 77 (and some C routines) by Stefan Turekdesigned for education, scientific research and industrial applicationsfull source-code and user manuals are available onlinemany extensions are not included in the official releaseVisit the FeatFlow homepagehttp://www.featflow.de3
FeatFlow 1.xThe Virtual Album of Fluid //www.featflow.de/albumTheoretical background4
Discretization techniquesIncompressible Navier-Stokes equationsut ν u u · u p f , · u 0,in Ω (0, T ]Spatial discretization techniquesnonconforming rotated multilinear finite elements for upiecewise constant pressure approximation for pSamarskij upwind or streamline diffusion stabilizationTemporal discretization techniquesimplicit one-step-θ-scheme (Backward Euler, Crank-Nicolson)implicit fractional-step-θ scheme (second-order accurate)adaptive time-stepping based on local discretization error5
Solution techniquesDiscretized incompressible Navier-Stokes equationsGiven un , g and k, solve for u un 1 and p pn 1[M θkN (u)]u kBp g,BT u 0where g [M θ1 kN (un )]un θ2 kf n 1 θ3 kf nNonlinear/linear solution strategiescoupled fixed point defect correction method CC2D/CC3Dnonlinear discrete projection scheme PP2D/PP3Dlinear multigrid techniques with adaptive step-length controlILU/SOR or Vanka-like block Gauß-Seidel smoother/solver6
FeatFlow 2.0Finite Element Analysis Toolbox Flow solverThe modern successor of FeatFlow 1.x for the numericalsolution of flow problems by the finite element methodmodular object-oriented design by Michael Köster et al.written in Fortran 95 (kernel applications)external libraries in F77/C (BLAS, UMFPACK, LAPACK)designed for education of students and scientific researchdetailed in-place documentation of the source-codeOfficial release not yet available; get the ALPHA snapshothttp://www.featflow.de/download/Featflow2 2.0ALPHA.tar.gz7
PrerequisitesUnix/Linux and Mac OS Xcompatible C and F95 compilerGCC and G95 version 0.91Intel R C /FortranCompilers for LinuxSun Studio C, C andFortran CompilersGNU make utilityMakefiles are providedfor all applicationsWindows XP, VistaMicrosoft R VisualStudio2003, 2005 or 2008Project files are providedfor all applicationsIntel R C /FortranCompilers for WindowsCygwinTMenvironmentGeneral Mesh Viewer (GMV)8
PrerequisitesUnix/Linux and Mac OS Xcompatible C and F95 compilerGCC and G95 version 0.91Intel R C /FortranCompilers for LinuxSun Studio C, C andFortran CompilersGNU make utilityMakefiles are providedfor all applicationsWindows XP, VistaMicrosoft R VisualStudio2003, 2005 or 2008Project files are providedfor all applicationsIntel R C /FortranCompilers for WindowsCygwinTMenvironmentGeneral Mesh Viewer (GMV)Linux is used for this workshop8
Getting FeatFlow 2.0Unpack the downloaded archive file tar xvzf Featflow2 2.0ALPHA.tar.gzChange into the base directory cd Featflow2 ; robjectGlobals.sparcreadme.txtGlobals.x86Rules apps f90.mkGlobals.x86 64 Rules apps.mkkernelRules libs.mklibrariesVERSIONSMakefilematlab9
Getting FeatFlow 2.0Unpack the downloaded archive file tar xvzf Featflow2 2.0ALPHA.tar.gzChange into the base directory cd Featflow2 ; robjectGlobals.sparcreadme.txtGlobals.x86Rules apps f90.mkGlobals.x86 64 Rules apps.mkkernelRules libs.mklibrariesVERSIONSMakefilematlab9
Getting FeatFlow 2.0Unpack the downloaded archive file tar xvzf Featflow2 2.0ALPHA.tar.gzChange into the base directory cd Featflow2 ; robjectGlobals.sparcreadme.txtGlobals.x86Rules apps f90.mkGlobals.x86 64 Rules apps.mkkernelRules libs.mklibrariesVERSIONSMakefilematlab9
Building FeatFlow 2.0Top-level build options make [ALT xxx] [ID yyy] target Some values for targethelpidallappslibs-print additional help and further optionprint out settings for current IDcompile all libraries and application modulescompile all application modulescompile all libraries10
Building FeatFlow 2.0Top-level build options make [ALT xxx] [ID yyy] target Some values for targethelpidallappslibs-print additional help and further optionprint out settings for current IDcompile all libraries and application modulescompile all application modulescompile all librariesSome available make modifiersALT xxx - specify alternative ID-xxx to useID yyy- override the autodetected architecture ID by yyy10
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID conceptrun on x86 64 GNU/Linux make idMachine-ID (Barracuda) : pc64-core2-linuxCompilers to be used:C compiler:C compiler:Fortran compiler:F-Library archiver:C-Library archiver:Flags to be used:OPTFLAGS OPTFLAGSC OPTFLAGSCPP OPTFLAGSF OPTFLAGSDEBUG OPTFLAGSCDEBUG OPTFLAGSCPPDEBUG OPTFLAGSFDEBUG FCFLAGS CCFLAGS CPPFLAGS BUILDLIB BLASLIB LAPACKLIB LDLIBS LDFLAGS /usr/bin/gcc/usr/bin/g /usr/local/g95/32bit integers/0.91/bin/g95/usr/bin/ar/usr/bin/ar-O3 -m64 -ffast-math -fexpensive-optimizations -fprefetch-loop-arrays -mmmx -msse -msse2 -msse3-g-O0 -g -Wall -fbounds-check -ftrace full-pipe -fmod -march nocona-pipe -march nocona-pipe -march noconafeat3d feat2d sysutils umfpack2 amd umfpack4 minisplib lapack blas(standard BLAS, included in installation package)(standard LAPACK, included in installation package)11
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID concept make ID pc-core2-linux idrun on x86 64 GNU/LinuxMachine-ID (Barracuda) : pc-core2-linuxCompilers to be used:C compiler:C compiler:Fortran compiler:F-Library archiver:C-Library archiver:Flags to be used:OPTFLAGS OPTFLAGSC OPTFLAGSCPP OPTFLAGSF OPTFLAGSDEBUG OPTFLAGSCDEBUG OPTFLAGSCPPDEBUG OPTFLAGSFDEBUG FCFLAGS CCFLAGS CPPFLAGS BUILDLIB BLASLIB LAPACKLIB LDLIBS LDFLAGS /usr/bin/gcc/usr/bin/g /usr/local/g95/32bit integers/0.91/bin/g95/usr/bin/ar/usr/bin/ar-O3 -m32 -ffast-math -fexpensive-optimizations -fprefetch-loop-arrays -mmmx -msse -msse2 -msse3-g-O0 -g -Wall -fbounds-check -ftrace full-pipe -fmod -march nocona-pipe -march nocona-pipe -march noconafeat3d feat2d sysutils umfpack2 amd umfpack4 minisplib lapack blas(standard BLAS, included in installation package)(standard LAPACK, included in installation package)11
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID concept make ALT ifc idrun on x86 64 GNU/LinuxMachine-ID (Barracuda) : pc64-core2-linux-ifcCompilers to be used:C compiler:C compiler:Fortran compiler:F-Library archiver:C-Library archiver:Flags to be used:OPTFLAGS OPTFLAGSC OPTFLAGSCPP OPTFLAGSF OPTFLAGSDEBUG OPTFLAGSCDEBUG OPTFLAGSCPPDEBUG OPTFLAGSFDEBUG FCFLAGS CCFLAGS CPPFLAGS BUILDLIB BLASLIB LAPACKLIB LDLIBS LDFLAGS ar/usr/local/intel/fce/10.1.021/bin/xiar-O3 -ipo -xT-g-traceback-traceback-warn all -check all -traceback-cm -fpe0 -vec-report0 -module-vec-report0-vec-report0feat3d feat2d sysutils umfpack2 amd umfpack4 minisplib lapack blas(standard BLAS, included in installation package)(standard LAPACK, included in installation package)-lsvml11
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID concept make [ALT xxx] [ID yyy] target compiler settings are defined in the global configuration filesGlobal.[alpha,ia64,mac,power,sparc,x86,x86 64]new compilers and/or architectures can be easily includedspecial purpose settings, e.g. pc64-opteron-linux-ifclarge11
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID concept make [ALT xxx] [ID yyy] target compiler settings are defined in the global configuration filesGlobal.[alpha,ia64,mac,power,sparc,x86,x86 64]new compilers and/or architectures can be easily includedspecial purpose settings, e.g. pc64-opteron-linux-ifclargeBuilding the Poisson example application cd applications/poisson make. after some time .Done, poisson-pc64-core2-linux is ready.11
Building FeatFlow 2.0, cont’d.Understanding the ALT/ID concept make [ALT xxx] [ID yyy] target compiler settings are defined in the global configuration filesGlobal.[alpha,ia64,mac,power,sparc,x86,x86 64]new compilers and/or architectures can be easily includedspecial purpose settings, e.g. pc64-opteron-linux-ifclargeBuilding the Poisson example application cd applications/poisson make debugturn on debugging facilities of the compiler. after some time .Done, poisson-pc64-core2-linux is ready.11
Preprocessing: Grid generation12
Geometric multigrid approach1Construct initial coarse grid in the external preprocessing step2Generate hierarchy of regularly refined meshes in the application13
Geometric multigrid approach1Construct initial coarse grid in the external preprocessing step2Generate hierarchy of regularly refined meshes in the 213110519616142Two-level ordering strategyadopt all coordinates from coarser grid levels unchangedintroduce new coordinates at edge/face/cell midpoints13
Coarse grid descriptionFile format is described in the Feat2D/Feat3D manualsSupported element types: triangles, quads (2D) and hexahedra (3D)14
Coarse grid descriptionFile format is described in the Feat2D/Feat3D manualsSupported element types: triangles, quads (2D) and hexahedra (3D)Domain triangulation is specified in TRI file (2D/3D)coordinate values of vertices in the interiorparameter values of vertices at the boundaryelements in terms of their corner nodesfirst two lines are treated as header/comments!14
Coarse grid descriptionFile format is described in the Feat2D/Feat3D manualsSupported element types: triangles, quads (2D) and hexahedra (3D)Domain triangulation is specified in TRI file (2D/3D)coordinate values of vertices in the interiorparameter values of vertices at the boundaryelements in terms of their corner nodesfirst two lines are treated as header/comments!Boundary parametrization is specified in PRM file (only 2D)each boundary component is described by p [0, pmax ]the ‘interior’ is located left to the boundary do not mix up (counter-)clockwise orientationsupported boundary types: lines, (arcs of) circles14
DeViSoR Grid3DCoarse grid generator for FeatFlow and Feastwritten in Java OpenGL and published under the GPLavailable at http://www.feast.uni-dortmund.desend requests, bug reports to devisor@featflow.deon-line help system and self-contained tutorial includedUnpack the downloaded archive file unzip grid-3.0.21.zipStart the application cd grid-3.0.21 java -jar grid3d.jar15
Alternative grid generatorsGiD – the personal pre- and post-processorevaluation version is available at http://gid.cimne.upc.esfully automatic structured and unstructured coarse grid generatorsupports triangular, quadrilateral, tetrahedral, hexahedral elementsprovides an effective easy-to-use and geometric user interfaceGiD2Feat – set of tools to convert GiD meshes to PRM/TRI files16
Example application: Poisson equation17
Possion equationChange into the application source directory cd applications/poisson/src; lspoissonXd callback.f90 poisson.f90poissonXd methodYYY.f90Open the application main source file emacs poisson.f9018
Possion equationChange into the application source directory cd applications/poisson/src; lspoissonXd callback.f90 poisson.f90poissonXd methodYYY.f90Open the application main source file emacs poisson.f90Initialization and finalizationsystem init()storage init(999, 100)storage done()initialize system-wide settingsinitialize storage managementfinalize storage management18
Possion equation, cont’d.Sample problem: u fin Ω (0, 1)2 ,u 0 on Ωright hand sidef (x, y) 32(x(1 x) y(1 y))analytical solutionu(x, y) 16x(1 x)y(1 y)Open the demonstration module file emacs poisson2d method0 simple.f90contains the corresponding subroutineincludes all required kernel modulesprovides detailed step-by-step tutorialOpen the callback function module file emacs poisson2d callback.f9019
Possion equation, cont’d.Grid generationboundary read prmtria readTriFile2Dtria quickRefine2LevelOrderingtria initStandardMeshFromRawread boundary parametrizationread domain triangulationperform regular refinementgenerate data structuresSpatial discretizationspdiscr initBlockDiscr2Dspdiscr initDiscr simpleprepare block discretizationinitialize spatial discretizationbilf createMatrixStructurecreate scalar matrix structurebilf buildMatrixScalardiscretize the bilinear formlinf buildVectorScalardiscretize the linear form/r.h.s.20
Possion equation, cont’d.Dirichlet boundary conditionsboundary createRegionbcasm newDirichletBConRealBDmatfil discreteBCvecfil discreteBCrhsvecfil discreteBCsolspecify boundary segmentcalculate discrete b.c.’sset b.c.’s in system matrixset b.c.’s in right hand sideset b.c.’s in solution vectorLinear BiCGStab solverlinsol initBiCGStablinsol setMatricesinitialize linear solverattach system matrixlinsol initStructure, linsol initDatalinsol solveAdaptivelysolve linear system21
Possion equation, cont’d.Solution outputucd startGMVucd addVariableVertexBaseducd addVariableElementBaseducd writestart export on GMV formatadd vertex-based solution dataadd cell-based solution datawrite solution data to fileClean-up and finalizationXXX release, XXX releaseYYYfree allocated memoryXXX donestop sub-systemGeneral naming convention of subroutinesabbreviatedModulefile NameOfSubroutine22
Postprocessing: Visualization of the solution23
General Mesh Viewer3D scientific visualization tooldeveloped at Los Alamos National Lab by Frank Ortegaavailable at http://www-xdiv.lanl.gov/XCM/gmv/supported OS: UNIX/Linux, Mac OS X, Windows (Cygwin)unstructured meshes in 2D/3Dcutlines, cutplanes, cutspheresvertex-based, cell-based data setscontour, vector plotsAlternative: importer for ParaViewbased on development CVS-version24
Advanced topics: Multigrid, Mesh Adaptation, . . .25
Possion equation, revisitedSample problem: u fin Ω (0, 1)2 ,u 0 on Ωright hand sidef (x, y) 32(x(1 x) y(1 y))analytical solutionu(x, y) 16x(1 x)y(1 y)Open the demonstration module file emacs poisson2d method1 mg.f90linear geometric multigrid solverJacobi or ILU(0) smootherdirect coarse grid solver (UMFPACK)26
Anisotropic diffusionExample: applications/anisotropicdiffusion.f90 cos θ sin θk10cos θ sin θD sin θcos θ0 k2 sin θ cos θ · (D u) 0u 0 u 2k1 100,k2 1,in Ωon Γ0on Γ1θ π6linear finite elements, h 1/36Galerkin fails: umin 0.0553h-adaptation: umin 0.006827
Anisotropic diffusionExample: applications/anisotropicdiffusion.f90 cos θ sin θk10cos θ sin θD sin θcos θ0 k2 sin θ cos θ27
Mesh adaptationconformal mesh refinement based on red-green strategyvertex-locking algorithm for mesh re-coarsening procedurenodal generation function stores birth certificatesprovides complete characterization of elementsyoungest node corresponds to refinement levelis required for the vertex-locking algorithmstate function for element characterizationlocal two-level ordering strategy28
Mesh adaptationconformal mesh refinement based on red-green strategyvertex-locking algorithm for mesh re-coarsening procedurenodal generation function stores birth certificatesprovides complete characterization of elementsyoungest node corresponds to refinement levelis required for the vertex-locking algorithmstate function for element characterizationlocal two-level ordering strategyThis has been addressed in the talk on Monday28
Mesh adaptationconformal mesh refinement based on red-green strategyvertex-locking algorithm for mesh re-coarsening procedurenodal generation function stores birth certificatesprovides complete characterization of elementsyoungest node corresponds to refinement levelis required for the vertex-locking algorithmstate function for element characterizationlocal two-level ordering strategy28
Mesh genealogy, revisitedTriangulation Tm (Em , Vm ), m 0, 1, 2, . . . consists ofEm {Ωk : k 1, . . . , NE }and Vm {vi : i 1, . . . , NV }nodal generation function g : Vm N0 is defined recursivelyg(vi ) : 0max g(vj ) 1vj Γkl max g(vj ) 1vj Ωkif vi V0if vi Γkl : Ω̄k Ω̄lif vi Ωk \ Ωk29
Characterization of elementsIdea I: State function s : Em Z (in MSB representation)Set Bit[0] to 1 for quadrilateral, otherwise set it to zeroSet Bit[k 1.4] to 1 if both endpoints of edge k have same ageIf no two endpoints have same age, then find local position kof the youngest vertex, set Bit[k] to 1 and negate the stateIdea II: Define local ordering strategy within each element a priorielement characterizationtriangle/quadrilateral from T0green quadrilateralred quadrilateralinner red triangle’other’ trianglegreen triangleelement state0/13, 5, 9, 11,17, 217, 13, 19, 25142, 4, 8-8, -4, -230
Characterization of elementsIdea I: State function s : Em Z (in MSB representation)Set Bit[0] to 1 for quadrilateral, otherwise set it to zeroSet Bit[k 1.4] to 1 if both endpoints of edge k have same ageIf no two endpoints have same age, then find local position kof the youngest vertex, set Bit[k] to 1 and negate the stateIdea II: Define local ordering strategy within each element a priorielement characterizationtriangle/quadrilateral from T0green quadrilateralred quadrilateralinner red triangle’other’ trianglegreen triangleelement state0/13, 5, 9, 11,17, 217, 13, 19, 25142, 4, 8-8, -4, -230
Characterization of elementsIdea I: State function s : Em Z (in MSB representation)Set Bit[0] to 1 for quadrilateral, otherwise set it to zeroSet Bit[k 1.4] to 1 if both endpoints of edge k have same ageIf no two endpoints have same age, then find local position kof the youngest vertex, set Bit[k] to 1 and negate the stateIdea II: Define local ordering strategy within each element a priorielement characterizationtriangle/quadrilateral from T0green quadrilateralred quadrilateralinner red triangle’other’ trianglegreen triangleelement state0/13, 5, 9, 11,17, 217, 13, 19, 25142, 4, 8-8, -4, -230
Characterization of elementsIdea I: State function s : Em Z (in MSB representation)Set Bit[0] to 1 for quadrilateral, otherwise set it to zeroSet Bit[k 1.4] to 1 if both endpoints of edge k have same ageIf no two endpoints have same age, then find local position kof the youngest vertex, set Bit[k] to 1 and negate the stateIdea II: Define local ordering strategy within each element a priorielement characterizationtriangle/quadrilateral from T0green quadrilateralred quadrilateralinner red triangle’other’ trianglegreen triangleelement state0/13, 5, 9, 11,17, 217, 13, 19, 25142, 4, 8-8, -4, -230
Characterization of elementsIdea I: State function s : Em Z (in MSB representation)Set Bit[0] to 1 for quadrilateral, otherwise set it to zeroSet Bit[k 1.4] to 1 if both endpoints of edge k have same ageIf no two endpoints have same age, then find local position kof the youngest vertex, set Bit[k] to 1 and negate the stateIdea II: Define local ordering strategy within each element a priorielement characterizationtriangle/quadrilateral from T0green quadrilateralred quadrilateralinner red triangle’other’ trianglegreen triangleelement state0/13, 5, 9, 11,17, 217, 13, 19, 25142, 4, 8-8, -4, -2dca4312b#el 3d1a1#el 5#el 411c#el 2b#el 130
Storage managementConcept of dynamically allocatable arrays in Fortran 9xinteger, dimension(:), allocatable :: e array of size n dynamicallydeallocate dynamically allocated array31
Storage managementConcept of dynamically allocatable arrays in Fortran 9xinteger, dimension(:), allocatable :: Iarrayallocate array of size n allocate dynamically allocated arrayConcept of dynamically allocatable arrays in FeatFlowstorage init, storage donestorage newinitialize/finalize storageallocate new memory storagestorage reallocre-allocate memory storagestorage freedeallocate memory storagestorage getbase XXXaccess memory storage31
Storage management, cont’d.Supported data: 8/16/32/64 integer, SP/DP real, logical, stringsMemory storage is accessible via integer handle ihandlePointer to the memory is assigned via storage getbase XXXImplementation detailsmemory storage handling is internally mapped to de/allocatedifferent storage management can be implemented without changeshandles can be easily passed to subroutines and functionsderived types typically consist of a set of handles scalar data32
Storage management, cont’d.Derived type for triangulation structureSupported data: 8/16/32/64 integer, SP/DP real, logical, stringstype t triangulationMemorystorage is accessible via integer handle ihandleinteger :: ndim 0Pointerto the ::memoryvia storage getbase XXXintegerNVT is assigned0integer :: NMT 0integer:: NAT 0Implementationdetailsinteger :: NEL 0memoryhandling is internally mapped to de/allocate. . storage.integer::managementh DvertexCoords ST NOHANDLEdifferentstoragecan be implementedwithout changesinteger :: h IverticesAtElement ST NOHANDLEhandlescan be::easilypassed to subroutines and functionsintegerh IneighboursAtElementST NOHANDLEendtypederived types typically consist of a set of handles scalar data32
Example: Using handlesCreate new storagestorage new (scall, sname, Isize, ctype, ihandle, &cinitNewBlock, rheap )ctype {ST DOUBLE,ST SINGLE,ST INTx,ST CHAR,ST LOGICAL}cinitNewBlock {ST NEWBLOCK ZERO,ST NEWBLOCK NOINIT}Accessing storage, e.g. 2-dimensional double arrayreal(DP), dimension(:,:), pointer :: p Darraystorage getbase double2D (ihandle, p Darray, rheap )Releasing storagestorage free (ihandle)33
UnConventional High Performance Computing34
FeastFinite Element Analysis & Solution ToolsHigh performance finite element package for the efficientsimulation of large scale problems on heterogeneous hardwarewritten in Fortran 90 and C (MPI communication)CFD (Stokes, Navier-Stokes), CSM (Elasticity)macro-wise domain decomposition approachfast and robust parallel multigrid methodsunconventional hardware as FEM co-processorsVisit the Feast homepagehttp://www.feast.uni-dortmund.de35
Hardware-oriented NumericsScalable Recursive Clustering (ScaRC) solverCombine domain decomposition and cascaded multigrid methodsGlobally unstructuredLocally structuredhide anisotropies to improve robustness36
UnConventional HPCCell multi-core processor (PS3)7 synergistic processing units@ 3.2 GHz, 218 GFLOPS/sMemory @ 3.2 GHzGraphic processing unit (GPU)128 parallel scalar processors@ 1.35 GHz, 350 GFLOP/sGDDR3 memory @ 900 MHz37
Future visionUnified Feat Feast finite element packageDecompose the domain into globally unstructured macro-cellsUse generalized tensor product grid or unstructured mesh locallyReuse FEM co-processors in the FeatFlow packageEnable special features, (e.g. h-adaptation) per macro-cellOn-line resources and additional materialFeatFlow project:http://www.featflow.deFeast project:http://www.feast.uni-dortmund.de38
1 Open-sourceCFDsoftware:FeatFlow TheFEA(S)Tgroups Institute of Applied Mathematics, LS III Dortmund University of Technology, Germany matthias.moeller@math.tu-dortmund.de
