Real-time landscape-size convective clouds simulation and rendering

: This paper presents an e(cid:30)cient, physics-based procedural model for the real-time animation and visualization of cumulus clouds at landscape size. We couple a coarse Lagrangian model of air parcels with a procedural ampli(cid:28)cation using volumetric noise. Our Lagrangian model draws an aerology i.e., the atmospheric physics of hydrostatic atmosphere with thermodynamics transforms, augmented by a model of mixing between parcels and environment. In addition to the particle-particle interactions, we introduce particle-implicit environment interactions. In contrast to the usual (cid:29)uid simulation, we thus do not need to sample the transparent environment, a key property for real-time e(cid:30)ciency and scalability to large domains. Inheriting from the high-level physics of aerology, we also validate our simulation by comparing it to predictive diagrams, and we show how the user can easily control key aspects of the result such as the cloud base and top altitude. Our model is thus fast, physical and controllable.


Introduction
Clouds are ubiquitous in outdoor naturalistic scenes, especially cumulus clouds which make interesting and conspicuous patterns.Thus they form an important component of many outdoor graphics applications, such as video games and ight simulators.Since these applications are becoming increasingly realistic, covering large and distributed landscapes, a similar trend is expected from their clouds which is currently not the case.Moreover, clouds could be one of the secondary features in these applications.Therefore, their animation and rendering should take only a small fraction of the total frame time.
The requirements are contradictory.Purely procedural models like [1,2,3] cannot convey the animated look of cumulus clouds which is a complex and evolving volumetric phenomena.On the other hand, traditional uid simulation methods like computational uid dynamics (CFD) simulations are not aordable in real-time, especially at landscape scale or with thin 3D details.
The evolving cloud simulation through huge Eulerian solvers is challenging even for physicists, due to the fact that small-scale and large ranges are required at the same time causing the numerical blending to aect the thermodynamics.In the scope of real-time graphics, cloud simulations do exist, for example based on Eulerian physics [4].However, they address only very limited volumes or show very coarse resolutions.One of the big disadvantages both with CFD and Eulerian methods is that most of the sky volume is lled with invisible air that needs to be stored and computed.Cloud simulation has also been proposed using Lagrangian framework [5] but such methods are dependent on a high particle count to achieve realistic resolution.
We propose a hybrid physics-inspired procedural method for the real-time simulation and rendering of detailed clouds.The main contributions of our approach are: • We introduce smoothed-particle hydrodynamics (SPH) based Lagrangian framework for the cloud physics, operating at the macro-level on large particles called thermals making it computationally inexpensive.• To avoid sampling the ubiquitous clear sky around clouds, we extend the Lagrangian formalism to account for the hybrid particle-implicit environment which is the key for real-time and scalability.For closure, we also introduce a new formulation of cloud-environment mixing (entrainment/detrainment). • The use of Lagrangian formalism, aerology physics and separation of environment from the dynamics brings more robustness, oers meaningful high-level handles for the user control on environment prole and evolution, cloud base and top.• For visualization, we treat these thermals as procedural units containing hypertexture.The volumetric texture inside a thermal operates at the micro-level and is produced based on its macro physical state.The model lightness makes possible the rst real-time, landscape-size convective cloud interactive simulation, with rendering capturing small-scale details.It is at the same time physically validated and controllable, which is important for applications aiming at the contradictory requirements of natural look, predictable aspect and performance.

Related Work
Much work exists on uid simulation in the graphics literature, most of which is out of reach of real-time constraints for cloud animation or rendering.Moreover, these methods often account for a very limited volume since transparent regions generally cost as much as the visible ones.There do exist a few multiscale simulation schemes which make use of multigrid and nite elements but the performance is still far from real-time.Moreover, cloud physics includes thermodynamics which is a crucial ingredient (and constraint) since advected scalars are reactive; clouds are not just smoke.
The existing approaches for cloud animation in computer graphics follow dierent strategies.In [6] the cloud evolution is simulated using cellular automaton where the dynamics are specied using simple transition rules.[7] proposes a method for animating clouds surrounding the earth at planetary level based on high and low pressure regions.Mostly [8,9,4] account for thermodynamics and are thus directly relevant to us.However, they come with all the problems of Eulerian simulations like limited domain, coarse resolution, harsh initial conditions, energy damping, numerical dissipation, lack of control, etc. Eulerian approaches also suer a lot from low resolution since articial blending modies the physics, damping energy, etc. misrepresenting important characteristics such as vorticity and free boundaries.
Lagrangian approaches like SPH [10] are considered more reliable on moving features, even at low resolutions.[5] recently proposed a particle-based method for cloud simulation based on position-based uids [11].Their approach is based on physics and thermodynamics and relies on sampling air in a limited domain which involves merging and splitting SPH particles in transparent and cloudy regions respectively.Consequently, a high particle count is necessitated wherein physics computation takes up the bulk of computational time.and visualization is achieved as a secondary, oine step.
A key issue with sampling the transparent air is that equations are only a part of the problem; initial conditions, state of the whole atmosphere and evolving boundary conditions like humidity and temperature on the oor are equally crucial.It is therefore, advisable to separate inputs comprising the tropospheric environment (for which the hydrostatic hypothesis is very reasonable), from simulated outputs the clouds.
In terms of the control of uids, the existing works aim at artistic control of uid animation or pure proceduralism [12,13].This is totally dierent to controlling a few high-level, visible characteristics emerging from the simulation such as nebulosity, size and altitude range of clouds.Contrary to these methods, we rely on the high-level physics to simply govern these emerging characteristics.Some work on real-time visualization of clouds is also based on 2D textural cloud layer or sprite-based clouds.The rendering part in our technique could benet by using one or more of these methods.The reader is referred to [14] for a more detailed insight into the rendering methods for clouds.Rising air parcels are subjected to the laws of thermodynamics, buoyancy and mixing with the surrounding environment, which we briey describe here.
Thermodynamics: Humidity is characterized by the mixing ratio ω which is the weight of vapor in grams per kilogram of dry air.The saturated mixing ratio ω s is the limit beyond which the extra vapor condenses into liquid water, releasing latent heat L(= 2260kJ/kg) per gram of condensed mass.The temperature change in a parcel with mass m due to condensation (or vaporization) of ∆m of vapor is considering adiabatic conditions.Here C p = 1.006 kJ/kg.K is the specic heat of air and α is a constant (similar to [5]) which is empirically set to 0.3.
The saturated mixing ratio is obtained from the pressure P (which changes with altitude h) and the saturated vapor pressure e s (which depends only on the temperature T in Kelvin): e s = 6.11 × 10 where P 0 ≈ 1013.25 hPascal.
Temperature, buoyancy and stability: Two key ingredients here are the atmospheric prole (temperature and humidity with altitude) and the initial conditions on the ground where the thermal bubbles start from.With the help of these parameters, one can obtain the altitude of cloud base known as the condensation level where the saturation starts (ω = ω s ) and the top where the buoyancy vanishes and the parcel reaches a stable equilibrium.To ease the comparisons and computations, atmospheric physicists introduced several alternate notions of temperature.
While potential temperature Θ = T P0 P 0.286 and virtual temperature T v = T (1 + 0.61ω) account for the change of pressure and moisture respectively on parcel, virtual potential temperature Θ v accounts for the combined eect of both.
With these denitions, Archimedes' buoyancy force per unit mass can be simplied to (by using ideal gas law P ∝ ρT ).
Parcel-environment interactions: Mixing has a profound eect on the parcels since it modies water and air content as well as the temperature.In this context, entrainment is engulng relatively dry environment air in the cloudy parcel while detrainment is leaking of cloud mass into the environment.A commonly accepted formulation often assumes the form δm mδz = ε − δ where δm is the change in the mass of parcel, z its altitude and ε and δ are entrainment and detrainment constants respectively [15].Parcels are also subjected to friction or drag force f airDrag = − 1 2 C D ρv 2 A where v is the relative velocity of parcel, A its reference surface and C D is the drag coecient.
From Laplace's law P V γ = constant (γ ≈ 1.4 for air) and accounting for possible mass ux and pressure change, the parcel's volume is estimated using SPH: Since our method is based on the discretization of uid physics with particles, we need to introduce the SPH formalism.The particles inuence each others' movements through pressure and viscosity forces.In SPH a weighting kernel W within a given support radius sr is associated to the particles and pressure is solved through a state law around equilibrium.Density is computed using the neighborhood set over support radius (summation over j) for each particle where m i refers to the mass of i th particle and x i its position.Pressure follows the spring mass formulation, where uid rest density is given by ρ 0 and K is the stiness constant.
The pressure force is derived from the pressure using Navier-Stokes equation and the viscosity force (or the drag between particles) using relative velocity between neighboring particles (µ is viscosity coecient) We refer the reader to [16] for the weighting kernels employed in our method and for more in-depth survey of SPH. 4 Model Overview We adopt a two-scale model, with large particles called parcels acting as the physics units and simultaneously containing volumetric texture.Through this decoupling between the macro and micro levels, a physics-inspired procedural method is made possible.We introduce these levels in the following and cover them in detail in section 5 and section 6 respectively.
1. Macro level : On the physics front, our particle solver inherits most of the SPH uid simulation formalism with the dierence that our parcels are large and correspond to the thermals in [12], with very little overlapping.Unlike existing methods, the environment prole is modelled implicitly and consists of piecewise linear potential temperature and humidity.This enables us to have particles only for thermals and clouds and not for the whole sky, thereby avoiding expensive physics and value storage for transparent regions.
For closure, we introduce a model of particle-implicit environment interaction, hybridised with SPH.
2. Micro level : On the visualization front, our amplication relies upon on-the-y procedural hypertexture generation within the parcels and a GPU volume ray-tracer.We evaluate a per-parcel lling ratio factor used to clamp a portion of the volume to account for the internal heterogeneity occurring based on the amount of condensed water content.This is inspired by the fact that condensation rst occurs on the top of the parcel and at cloud decay lower average droplet density corresponds to a fractal rather than a continuous haze.We rely on the classical Lagrangian model for the particle-particle interaction.The inter-parcel interaction is captured by applying repulsive pressure force f pressure (Equation 9) that prevents parcels from getting too close to each other and the viscosity force f viscous (Equation 10) that accounts for the drag between the parcels.Besides the SPH forces (in our model), cloud parcels are also subjected to the thermodynamic rules.Thus in the following, we describe the interaction of parcels with the implicit environment, inter-parcel mixing and the buoyant forces acting on parcels, thereby completing the total parcel physics leading to cloud formation.

Interactions with the implicit environment
In our model, particles cover only a part of the atmosphere which still must represent a continuum of all physical interactions.We also therefore, need to estimate the contribution of the implicit (i.e., unsampled) environment on the area of parcels exposed to the atmosphere.Thus to these particle-particle interactions, we have to add forces corresponding to the equivalent terms in particle-environment interactions, f atmP res and f airDrag respectively.Atmospheric force f atmP res is applied to the exposed surface area of each parcel which is in direct contact with the atmosphere.The exposed area (denoted by S) is determined by coarsely tessellating each parcel into quadrilaterals and determining the set of quadrilaterals that do not lie inside any other parcel.Instead of querying all the parcels, a quadrilateral needs to perform this test only on the SPH neighbors of the parcel it belongs to.Atmospheric pressure is obtained using Equation 2 and is applied to each tessellated quadrilateral along its normal weighted by the area.f atmP res needs an additional weighting factor, which we obtain through calibration based on our experiments.

5.2
Mixing with the environment The parcels interact with the surrounding environment through the process of entrainment and detrainment.The total amount of entrained air in a parcel during time interval ∆t is a function of both its exposed area S and the velocity, which in essence gives the swept volume of the parcel through the external relatively dry air.Even if a parcel is globally motionless, smallscale turbulence pursues the mixing; otherwise clouds would not decay.To account for this, we thus add to the macroscopic drag a small turbulence-related term δv t .Furthermore, fractal nature of clouds makes the real cloud-air contact surface larger than one for the macroscopic sampling primitives.For this we assume proportionality coecients, K e for entrainment and K d for detrainment, which is the surface multiplier to account for subgrid fractal surface.Entrainment entails the mixing of external air ∆m a e (Equation 11), which also carries with it the water vapor ∆m w e (Equation 12) already present in it (computed using mixing ratio), with the parcel.The total mass gain ∆m e is obtained by summing these two quantities (Equation 13).To incorporate it in our formulation, we modify the standard entrainment and detrainment constants to ε = ε r and δ = δ r, where r is the radius of the parcel.Here ε is the entrainment rate taken from [15] and A is the exposed surface area obtained from the computed S. The entrained mass comprising both dry air and vapor is at a dierent temperature than the parcel.The equilibrium temperature Θ is computed by mixing the inux with the existing parcel contents m, see Equation 14.
Subsequently, the quantity of detrained content ∆m d (Equation 15) is determined using a similar expression as for entrainment albeit with a few changes.ρ atm in Equation 11 is substituted with the average density ρ in the parcel to compute ∆m a d .Assuming a uniform mixing of vapor in the parcel, an equal proportion of water mass ∆m w d is ejected together with ∆m a d .Furthermore, ε is also replaced by δ, the detrainment constant.The total mass detrained from a parcel is given by Equation 15.

Inter-parcel mixing
Much as the parcels exchange mass with the external environment through their exposed surface, a similar process takes between the neighboring parcels through the contact surface.The interparcel mass exchange ∆m p is obtained using a simple diusion process from high concentration parcel to lower ones surrounding it.With this, the combined mixing occurring to a parcel could be summed up as:

Buoyant force
In addition to aforementioned forces, buoyant force f buoy is added to each parcel which is a key component in parcel rising and cloud formation.With the computation of updated potential temperature Θ (Equation 14) inside the parcel, buoyant force follows from Equation 5based on the temperature dierence with the surrounding environment.With this, the total force acting on a parcel is: The micro level deals with the procedural generation of cloud details on top of the existing macro structure of parcels.The hypertexture is computed on-the-y, and is controlled by a few macro-level handles such as particle radius, condensation level plane and lling ratio which are themselves derived from the ongoing parcel physics.This is detailed in subsection 6.1.We employ GPU-based ray-tracing to eciently render our clouds as explained in subsection 6.2.

Volumetric cloud texture
From the local list of parcels intersecting a given point in space, we employ the hypertexture framework of [17] to compute the union of clamped hypertextured parcels (see Figure 2): • The pseudo-distance is computed as ∪ i (sphere i ∩half Space i ) where the sphere corresponds to the parcel and the half-space corresponds to clamping by the condensation level point plane.
Since our SPH kernels have low super-imposition, we use spheres larger than these kernels for nice looking unions.Note that we can stop the ∪ i loop as soon as it saturates to 1. • The sphere and half-space shell thickness is a user-tunable constant, but is modulated by the parcel lling ratio f to account for the cloud decay so that the noise texture gets sponge-like and erodes with the average density proportional to f .The density is then computed based on the noise and distance function (where k is a constant): The factor f estimates the percentage of a parcel lled with condensed (opaque) material.The idea is that the lling is not homogeneous.Condensation starts at condensation level (denoted by half-space) and grows making the parcel more opaque while decaying carves a sponge-like parcel radius ratio relative to its initial value.
texture.Therefore, in practice an average droplet density corresponds to a mix of clear and opaque fractions.This relative volume of opaque fraction is the lling ratio.We estimate it as f = ω cloud −ωs V air 1kg .ql0 where V air 1kg is the volume of 1kg of air at parcel altitude and q l0 is the droplet density in the opaque fraction.

Ray-tracing clouds
The real-time volumetric details generated from subsection 6.1 are visualized using volume ray casting.Our simple physics engine can easily be coupled to a reasonably minimal volume raycasting fragment shader.To this end, we have used the uniform buer object (UBO) which allows us to pass arrays to the fragment shader.We transfer the entire grid together with the list of condensed parcels and the related attributes (parcel size, lling ratio, plane position) as UBO to the GPU. 7

Results and discussion
The proposed method was implemented using C++, OpenGL version 4.0 and tested on an Intel Core 3.4 GHz machine with 16 GB RAM and an NVidia GeForce GTX 970 graphics card on a display window of 730 × 545 pixels.
The number of parcels is user chosen and can vary from single to many.Our model allows realistic, landscape-scale cloud eects by using just few macro primitives.The physics part therefore, is a lightweight computation which can be carried out on the CPU itself.However, this is not a limitation of our method as in case where larger number of macro primitives are desired, physics computation could be easily and eciently transferred to GPU as laid out in [18].We also cast simplied shadow of clouds on the OpenGL terrain.This is achieved by shooting rays from the light source towards the ground to collect cloud opacities in a 2D texture in the rst pass and blending them with the scene at high resolution in the second pass.
The parcel tessellation for atmospheric pressure and mixing computation is kept coarse (18×9) in our case which provides a good balance between accuracy and computational eort.We set the environment dimensions to 10km × 10km × 10km where the simulation and visualization is carried out, although it is exible and user modiable (and could as well be irregular).In terms of Lagrangian parameters, density of air is taken to be around ρ 0 = 1.28kg/m 3 , stiness constant K = 100 (since air is compressible) and a low viscosity coecient µ ≈ 0.01 is chosen.Single-parcel kinematics is veried against aerology diagram for the following important parameters: • Lapse rates: As computed directly from Figure 4a, the dry adiabatic lapse rate comes out to be around −10 C/km (slope between points A and B), whereas the moist adiabatic lapse rate is about −5 C/km (average slope between points B and D).These values are close to the actual gures, as given in [19].Air parcels that get saturated as they rise cool at a smaller rate than the dry adiabatic lapse rate due to the heat released by the condensation of water vapor (B onwards).The wet adiabatic lapse rate varies with the altitude which is why we take into account the average value.
• Equilibrium altitude: In the case of no mixing with the atmosphere, the parcel reaches to a maximum altitude and comes to an equilibrium height after a series of damping oscillations, as is visible in Figure 4a, see also Figure 3a.The maximum altitude lies close to the point where parcel temperature intersects the atmospheric temperature prole.The parcel continues moving further up till its velocity is reversed by the negative buoyancy force.However, mixing makes a parcel lose its water content, thereby changing the internal dynamics.During this decaying process the parcel undergoes oscillatory motion too as it alternates between rising and sinking buoyant force [20], which is also validated by our simulation as shown in Figure 4b and Figure 3a.Note that while the parcel follows the same temperature trajectory during its oscillations in absence of mixing (Figure 4a), the temperature prole is altered with each oscillation if the mixing is enabled (Figure 4b).
• Mass and radius: The variation of parcel water mass and radius relative to its initial values is as demonstrated in Figure 3b and Figure 3c respectively.The initial mass gain with mixing enabled is explained by the fact that as per the model from [15], detrainment begins at condensation level while entrainment occurs throughout.It is for the same reason that the parcel begins to gain mass towards the end of its lifetime as it reaches the condensation level again.Eventually, the parcel exhausts its water content signicantly by leaking it to the atmosphere.The parcel radius grows as it ascends due to the decrease in atmospheric pressure.However, in the case of mixing with atmosphere the parcel radius shrinks due to the combined eect of water loss and descent.

User control
The condensation level altitude and the nature of clouds obtained from parcels is dependent on the initial mixing ratio and the temperature prole.The cloud top is governed by the temperature inversion, the altitude at which the temperature begins to increase and hence the clouds parcel cease to rise further up.Though the default implicit environment is predened in our system, our tool (Figure 5) gives the user an easy and intuitive handle to: • alter the temperature prole graphically • select condensation level • choose mixing ratio of the atmosphere The system simulates the physics with the user requirements and the clouds are generated if the state is physically plausible.1: Performance split-up between the physics part on the CPU and the rendering part on the GPU together with the frame rates per second (fps).
Figure 1 illustrates the complete life cycle of clouds as parcels from their birth to death.Clouds often die out in order of the sequence of their creation and also depending on their exposure to the outside dry air.

Performance analysis
We have tried upto 85 dynamically produced parcels where the simulation together with the rendering continues to be interactive.The performance analysis of our system is summarised in Table 1, where the percentages are given relative to the total time spent in the system per frame.In there we also compare our performance with and without computation of normals (for voxel shading) using nite dierence method (FDM) for rendering.Our experiments suggest that rendering without nite dierences normal is signicantly faster than with normals.The only expensive part in parcel physics is the computation of area exposed to implicit environment on the tessellated sphere.On the other hand, hypertexture production and raycasting based visualization takes a larger share of overall computational time as compared to physics.
In Figure 6, a realistic and distributed cloud cover is obtained with just 14 parcels.Further, interesting cloud patterns are obtained by taking the union of density parcels in the overlapping regions as shown in Figure 7.In addition to the above mentioned experiments, being a Lagrangian-based technique our method can easily deal with the eects of wind, parcel-obstacle collision which could inuence the cloud movement and formation.Furthermore, the atmospheric temperature and humidity prole can be varied as a function of time of the day to inuence cloud formation.These aspects of our method are demonstrated in the accompanying video.We have presented an ecient, physics-based procedural hybrid method for real-time landscapesize simulation and visualization of cumulus clouds.The use of Lagrangian framework at the macro parcel level makes the physics computationally inexpensive.The introduction of the implicit atmosphere eliminates the need to sample transparent air and accounts for the parcelair interaction and mass exchange.The real-time visualization of clouds is achieved by raytracing the procedural amplication of volumetric noise on the GPU.Our model is controllable by the user and predictive, and can be validated against the areology diagram for basic parameters.
For future work, we would like to add more complex physics in order to permit billowing of cumulus clouds,and study hierarchical Lagrangian simulation for even larger scale scalability.

Figure 1 :
Figure1: Cloud life cycle demonstrated using 75 parcels resulting from the space-time distributed emission from the ground with probabilities depending on the local wetness and temperature.Clouds condense when reaching the condensation level altitude, rise up to equilibrium altitude, then decay and die due to the eects of entrainment and detrainment.

Figure 2 :
Figure 2: Left: hypertexture in the shell (between the blue and the red curves) of sphere ∪ half-space above condensation level.right: with lling ratio f =50%.

Figure 3 :
Figure 3: Attributes of a parcel varying with time in our simulation, with mixing enabled and disabled: (left) parcel position, (mid) parcel water mass ratio relative to its initial value, (right)

8
Conclusions and future work

Figure 5 :
Figure 5: Our interactive tool allows the user to experiment with intuitive parameters to generate plausible evolving elds of cumulus clouds.

Figure 6 :Figure 7 :
Figure 6: Parcel physics and hypertexture generated with just 14 parcels leading to fast and yet, realistic visualization.