Despite the fact that computational fluid dynamics (CFD) software is now
(relatively) fast and freely available, it is still amazingly difficult to use.
Inaccessible software imposes a significant entry barrier on students and
junior engineers, and even senior researchers spend less time developing
insights and more time on software issues. Lily Pad was developed as an initial
attempt to address some of these problems. The goal of Lily Pad is to lower the
barrier to CFD by adopting simple high-speed methods, utilising modern
programming features and environments, and giving immediate visual feed-back to
the user. The resulting software focuses on the fluid dynamics instead of the
computation, making it useful for both education and research. LilyPad is open
source and available online at https://github.com/weymouth/lily-pad for all use
under the MIT license

Weymouth, G.D.

arXiv

Southampton (e-Prints Soton)

Lily Pad: Towards Real-time InteractiveComputational Fluid DynamicsGabriel D. Weymouth, University of Southampton, UKG.D.Weymouth@soton.ac.uk1 IntroductionAmazing advances have been made in both computing resources and computational methodologies in thelast 15 years - enabling first-principle predictions of ship science flows which are fast enough to be usefulin engineering. In addition, open source projects such as OpenFOAM [9] provide a free and adaptablealternative to standard engineering software. However, despite the fact that computational fluid dynamics(CFD) software is now (relatively) fast and freely available, it is still amazingly difficult to use. CFDsoftware is complex, fussy, opaque, and poorly designed [14], even compared to other scientific softwaresuch as the free and open source R-project for statistical analysis [13].1 Inaccessible software imposesa significant entry barrier on students and junior engineers, and even senior researchers spend less timedeveloping insights and more time on software issues [5].Lily Pad [1] was developed as an initial attempt to address some of these problems. The goal of LilyPad is to lower the barrier to CFD by adopting simple high-speed methods, utilizing modern programmingfeatures and environments, and giving immediate visual feed-back to the user. The resulting softwarefocuses on the fluid dynamics instead of the computation, making it useful for both education andresearch.It would be impossible to genuinely achieve this goal unless Lily Pad was genuine CFD solver. There-fore the full two-dimensional Navier-Stokes equations are solved and the exact body boundary conditionsare applied. However, most of the complications plaguing CFD are avoided by using the Boundary DataImmersion Method (BDIM [17, 7]) to immerse solid bodies into the fluid domain. Despite this simplicity(or perhaps because of it), BDIM is very accurate, it has some nice analytic properties, and has beenextensively validated.Lily Pad is written in Processing [12], a programming language initially based on Java which promotes“software literacy within the visual arts and visual literacy within technology”. Lily Pad utilizes theProcessing Development Environment, which integrates the writing, testing, and usage modes into asingle platform. Specifically, a short script is written on the front-page (Fig 1a), the code is executed bypressing a button, and the flow is simultaneously simulated and visualized in a pop-up window (Fig 1b).Fig 1 is a complete working example, the default simulation of the flow around a circular cylinder. Andit is interactive - the circle can be dragged around the window and the flow is updated in real-time.While the default flow in Fig 1 is at low resolution and features fairly simple physics, the method iscompletely extendable, and Lily Pad simulations have already been published in the Journal of Fluid Me-chanics [15, 11] and have been used to initiate and test many additional research projects, as summarizedin the examples section. First, the numerical method and implementation in Lily Pad is summarized inmore detail.2 Flow Solver MethodologySolving the incompressible Navier-Stokes equation numerically presents a few irreducible difficulties,namely time integration of the finite-scale nonlinear momentum equations and inversion of a Poissonequation for the pressure. However, the majority of the complexity in CFD software comes not fromthe fluid equations, but from coupling those equations to the irregular data at domain boundaries.The adoption of boundary-fitted meshes to achieve this coupling increases the memory requirements,computational effort and algorithmic complexity of the numerical method, as well as presenting conceptualand modeling difficulties to the user.1To say nothing of software that people actually enjoy using.(a) (b)Figure 1: Lily Pad example in the Processing Development Environment: (a) complete example script, and (b)real-time interactive visualization.Lily Pad utilizes fast and robust methods for dealing with all of these issues, starting with the adoptionof the Boundary Data Immersion Method for coupling the fluid and solid equations. The basic idea ofBDIM (and all Immersed Boundary methods since [10]) is to adjust the equations of motion to accountfor the interaction of the fluid and solid and solve those equations on a simple numerical grid. In LilyPad this is a uniform Cartesian grid. This avoids adopting questionable physical simplifications (suchas treating fluids as particles, with issues documented in [6]) or hiding the numerical complexity in ablack box (such as automatic mesh refinement). Instead, the analytic governing equations themselves areadjusted, allowing the numerical method to remain as simple and efficient as possible.The full details are in [7], but briefly, the analytic BDIM equations are obtained by the convolutionof two governing equations, each valid over their own region of the total domain. In Lily Pad this is abody with prescribed velocity immersed in a 2D viscous fluid. Using a symmetric kernel with width  forthe convolution and keeping up to second order terms gives~u = µ0 ~f + (1− µ0)~b+ µ1 ∂∂n(~f −~b)(1)~∇ · ~u = (1− µ0)~∇ ·~b− µ1 ∂∂n(~∇ ·~b)(2)where ~f and ~b are the update equations for the velocity in the fluid and body domains, and where µ0and µ1 are the zeroth and first moments of the kernel function with respect to the fluid domain. Theconvoluted velocity ~u update equation (1), takes the form of a mixing equation, with the fluid equationused in the fluid domain (where µ0 = 1, µ1 = 0) the body equation used in the body domain (whereµ0 = µ1 = 0) and a combination of the two within  of the boundary. The velocity divergence equation (2)is used to solve for the pressure required to satisfy the divergence condition in both domains, which enablesvolume-changing bodies found in biological flow problems to be modeled. Three important points: (a)the normal derivative term in the velocity equation (active in the boundary region) enables accurateboundary layer simulation including separation prediction [7]; (b) the presence of the kernel momentsresults in a variable coefficient Poisson equation for the pressure which implicitly enforces the correctpressure boundary conditions and significantly improves force prediction [17]; and (c) the body velocityequation ~b can be easily coupled to the state of the fluid, enabling prescribed or predicted fluid-structureinteraction simulations [15]. BDIM is therefore capable of fluid-structure interaction predictions using auniform Cartesian grid, and this is leveraged in Lily Pad to greatly simplify the software.Lily Pad uses implicit Large Eddy Simulation (ILES, [8]) to model the scales below the finite gridspacing in order to avoid the standard turbulence model selection and parameterization issues. Depsiteits simplicity (or, again, perhaps because of it) ILES has been shown to have greater effective subgridresolution than explicit LES models [2], but like all LES models it requires explicit time-stepping andis difficult to scale to high Reynolds number flows. Either a flux-limited QUICK or a Semi-Lagrangianconvection scheme can be used. The Semi-Lagrangian scheme is universally stable, enabling larger timesteps at the cost of reduced the accuracy near the body where the flow gradients are large. Additionally,because Lily Pad is two-dimensional, the results at intermediate Reynolds numbers will tend to havefeatures of two-dimensional turbulence, which can be quite distinct from their three-dimensional coun-terparts [4]. Within these limitations, the choice of ILES enables a wide variety of interesting flows to besimulated in Lily Pad quickly and easily.The variable coefficient pressure equation is solved using a Multi-Grid method with Jacobi smoothing.The unit Cartesian grid, sparse symmetric matrix, and relatively small time steps make this extremelyefficient; the residual is typically driven below the tolerance level in only a single V-cycle, making deter-mination of the pressure less expensive than evaluating the momentum equations.3 Implementation and UsageLily Pad is written in Processing in an object-oriented programming style, using classes of objects to setup and solve the CFD system. The core classes are as follows:• Field contains the data structures and method for handling scalar fields such as pressure. Thisincludes methods for differentiation, interpolation2 and advection of the field, as well as settingboundary conditions and displaying the scalar values on-screen.• VectorField contains the methods for vector fields such as the fluid velocity, made up of twoField components. A projection operator (using Multi-Grid) is defined here in addition to vectordifferentiation and the convection/diffusion operators.• Body is the parent class of all the solid geometries which are used in Lily Pad. Methods to transform(translate, rotate), query (force?, moment?, distance?) and display the body are defined.• BDIM is the high-level class used on the front-page to set-up and solve the general solid-fluid inter-action problem described by equations 1 and 2.The default script in Fig 1 shows how these classes and methods are used in Lily Pad. The Processingfront-page is broken into four parts: variable declaration, the setup method, the draw method, andadditional problem specific methods. In the default script three objects are declared:BDIM flow;CircleBody body;FloodPlot flood;where CircleBody extends the parent Body; defining the geometry to be a circle and supplying somesimplified methods such as the distance function. FloodPlot is a visualization class which makes niceflood plots to extend Field.display. Within setup we have:int n = (int)pow(2,6); // number of grid pointssize(400,400); // display window sizeWindow view = new Window(n,n);body = new CircleBody(n/3,n/2,n/8,view); // define geom2The unit background grid makes differentiation and interpolation operations nearly trivial.flow = new BDIM(n,n,1.5,body); // define flow using BDIMflood = new FloodPlot(view);flood.setLegend("vorticity",-.75,.75);The first line sets the number of points in the domain. This needs to be a large power of two (Javauses the pow function for exponents) to allow sub-division in the multi-grid solver.3 Next the size of thepop-up display is defined in number of pixels. The third defines a Window object to map the numericalgrid to a portion of the display. Any number of windows may be defined and used to plot different ‘views’of the simulation.In the next block, the body is constructed. The arguments define the size of the circle, and its locationin the simulation. This process is indefinitely extendable - any number of additional bodies of any sizeor shape can be defined to form a composite Body object. Note that the arguments are all given in termsof the number of grid cells, ie the diameter is n/8 = 8 meaning the circle is 8 grid cells wide. In Lily Padall distances are given in terms of the unit grid spacing since it is the only predetermined length scale.This has the pedagogical value of making all dimensions into statements of “resolution”. Since D is only8, we can safely assume this simulation is under-resolved. If we increase the number of points, this willincrease the resolution by the same factor, etc.Next, the flow is constructed using the number of points, time step, and body object as arguments.The default BDIM object uses a unit background flow ~U = (1, 0), but any velocity field can be defined andused as an initial condition. Since the distance scale is the grid spacing, the time step ∆t is equal to theCFL number. If ∆t is set to zero, adaptive time stepping is used with ∆t = (max~x |u|+ 3ν)−1 to ensureexplicit stability of the convective and diffusive terms. The kinematic viscosity ν is an optional argumentto the BDIM constructor and its value needs to be set using the grid length scale. For example to matchReD = UD/ν = 1000 the grid-based Reynolds number is Reg = U/ν = ReD/D = 125 and ν = 0.008.Again, there is pedagogical value in this approach - ν is so small that the physical viscosity is likely tobe drowned out by the ILES subgrid damping. Finally, the FloodPlot object is constructed by passingit a view, and then setting the color legend text and limits.It is worth emphasizing that the flow solver has been set-up in four lines of code. In standard CFDsoftware, setting-up the fluid domain requires tremendous time and effort and entails the explicit orimplicit specification of thousands of numerical parameters (the grid!). In Lily Pad, the fluid domainset-up is a non-issue and the only effective free numerical parameter is the ‘resolution’. 4The Processing draw method is looped indefinitely and is used to iterate the update and display calls:body.update();flow.update(body); flow.update2();flood.display(flow.u.vorticity());body.display();In this code, the body is updated, taking into account any changes input by the user. Alternatively,the body can be updated using prescribed motion or fluid-structure interaction equations. The updatedbody is used to update the conditions for the flow, and the flow is integrated in separate predictor andcorrector steps. Next the flood plot is used to display the vorticity of the flow velocity5 and the bodyis displayed on top. The display commands have been optimized for speed, and the educational value ofwatching the flow develop in real-time can hardly be overstated.4 Advanced Examples and Next StepsThis default example gives a basic introduction to the grammar used in a Lily Pad script, but one ofthe software’s strengths is how easily the methods are extended to different flow types. The most basicextension is to change the body’s shape and dynamics, such as in [15, 11] or a flapping foil energyextractor (Fig 2a). Increasing the number of bodies is also very simple, such as an array of 20 offshorerisers (Fig 2b). Combining multiple bodies with independent dynamics extends the potential application:from blind fish swimming past obstacles (Fig 2c) and drag reducing spinning cylinders (Fig 2e) to tandem3Multi-grid methods are especially well suited to object oriented programming, enabling a complete implementation inless than 90 lines of code.4The kernel width  is not a free parameter. Previous studies have shown that  = 2 is optimal.5Simple chaining of commands is another advantage of the object oriented coding style.(a) (b)(c) (d)(e) (f)Figure 2: Lily Pad example applications courtesy of (a) Amanda Persichetti, (b) Chen Yongxin, (c) AudreyMaertens, (d) Luke Muscutt, (e) James Schulmeister, and (f) Jacob Izraelevitz.swimming flippers of a turtle (Fig 2f) or even a Plesiosaur (Fig 2d). Extensions can also be made to thevisualization and user interface, as seen on the Lily Pad website [1].However, there are still limitations in this approach. The fundamental drawback is being restricted totwo-dimensional flow. A three-dimensional version, called Lotus, is in development but 3D simulationsincrease the computational cost by many orders of magnitude. Lotus attempts to compensate for thisby implementing multiprocessor versions (MPI and CUDA) and fast programming and visualizationlanguages (Fortran/Paraview) but is still much slower than real-time. Advances in 2D/3D models mayhelp bridge the gap, at least in the initial design and simulation set-up stages. High Reynolds numbersimulations will also require modeling advances, particularly in the development of near-wall modelsappropriate for Immersed Boundary methods.Nevertheless, two-dimensional and moderate Reynolds number simulations still have great value,especially in education. All the results in Fig 2 are student projects, and most were developed bystudents without formal training in CFD. Lily Pad’s ability to quickly and easily go from concept tosimulation has helped all of these students add depth to their studies. In addition, Lily Pad has beenused to develop stand-alone demo applications which are excellent for education and outreach [3, 16]. Assuch, Lily Pad has at least partially achieved its goal of lowering the barrier to CFD, and is a first steptowards more general real-time interactive solvers.References[1] Lily pad: Real-time two-dimensional fluid dynamics simulations in processing. http://github.com/weymouth/lily-pad. Accessed: 2015-09-29.[2] Stuart Dalziel Andrew Aspden, Nikos Nikiforakis and John B. Bell. Analysis of implicit les methods.Communications in Applied Mathematics and Computational Science, 3(1):103–126, 2008.[3] CNN MainSail. How has foiling made boats much faster? http://edition.cnn.com/videos/tv/2015/03/11/spc-mainsail-design-special-b.cnn/video/playlists/intl-mainsail. Ac-cessed: 2015-09-29.[4] Y Couder and C Basdevant. Experimental and numerical study of vortex couples in two-dimensionalflows. Journal of Fluid Mechanics, 173:225–251, 1986.[5] Derek Groen, Xiaohu Guo, James A Grogan, Ulf D Schiller, and James M Osborne. Softwaredevelopment practices in academia: a case study comparison. arXiv preprint arXiv:1506.05272,2015.[6] Areti Kiara, Kelli Hendrickson, and Dick KP Yue. Sph for incompressible free-surface flows. part i:Error analysis of the basic assumptions. Computers & Fluids, 86:611–624, 2013.[7] Audrey P Maertens and Gabriel D Weymouth. Accurate cartesian-grid simulations of near-bodyflows at intermediate reynolds numbers. Computer Methods in Applied Mechanics and Engineering,283:106–129, 2015.[8] LG Margolin. Finite scale theory: The role of the observer in classical fluid flow. Mechanics ResearchCommunications, 2013.[9] OpenFOAM Foundation. Openfoam. http://www.openfoam.org. Accessed: 2015-09-29.[10] Charles S. Peskin. Flow patterns around heart valves: A numerical method. Journal of Computa-tional Physics, 10(2):252 – 271, 1972.[11] Delyle T. Polet, David E. Rival, and Gabriel D. Weymouth. Unteady dynamics of rapid perchingmanoeuvres. Journal of Fluid Mechanics, 767(323–341), 2015.[12] Processing Foundation. Processing. http://www.processing.org. Accessed: 2015-09-29.[13] R Core Team. R: A language and environment for statistical computing. http://www.R-project.org. Accessed: 2015-09-29.[14] Damian Rouson, Jim Xia, and Xiaofeng Xu. Scientific software design: the object-oriented way.Cambridge University Press, 2011.[15] G. D. Weymouth. Chaotic rotation of a towed elliptical cylinder. Journal of Fluid Mechanics,743:385–398, 3 2014.[16] Gabriel D Weymouth. Interactive research examples. http://web.mit.edu/weymouth/www. Ac-cessed: 2015-09-29.[17] Gabriel D. Weymouth and Dick K.-P. Yue. Boundary data immersion method for cartesian-gridsimulations of fluid-body interaction problems. Journal of Computational Physics, 230(16):6233–6247, July 2011.

Lily Pad: towards real-time interactive computational fluid dynamics

Abstract

Similar works

Full text

Available Versions

Southampton (e-Prints Soton)