Camille Schreck

Chris Wojtan

Coupling 3D Liquid Simulation with 2D Wave Propagation for Large Scale Water Surface Animation Using the Equivalent Sources Method

Figure 1: A dragon springing out of water. Our method couples a FLIP fluid simulation and 2D wave propagation. The dotted ellipse shows the limit between the two.

Introduction

Progress in fluid simulation has made water scenes in animated movies and games increasingly realistic. However, the simulation of very large scenes with convincing visual details is still extremely expensive and challenging. This paper proposes an efficient method for simulating large bodies of water like a sea by splitting the domain into regions of complex and simple motions: We simulate complex fluid motions like splashes, breaking waves, and droplets with a fully 3D Navier-Stokes solver, while we simulate the calmer areas of water with a 2D heightfield solver based on fundamental solutions of the Helmholtz equation [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF].

Unlike previous methods for coupling 3D and 2D fluid animation techniques, our approach avoids the use of a fixed grid for wave simulation, which allows us to simulate arbitrarily detailed waves and infinitely large bodies of water. This 2D solver is based on Airy wave theory, which faithfully reproduces the relative speeds of waves in any water depth, not just shallow water. We then couple the simulations together using state-of-the-art techniques for limiting artificial wave reflections where the 2D and 3D domains meet.

Despite our measures to minimize numerical artifacts when coupling the simulations together, we show that differences in numerical errors between wave simulations can to produce distracting numerical artifacts. To solve this problem, we introduce a novel systematic approach to compensate for numerical errors between wave simulations, by matching the effective wave speeds of a 3D fluid solver. Our contributions are the following:

• A method to couple any 3D liquid simulation with a 2D wave propagation method based on fundamental solutions. • A method to estimate the effective dispersion law for any 3D liquid simulator, which we use to match wave properties between different simulations and eliminate spurious numerical refraction artifacts.

State of the Art

3D Liquid Simulation

The motion of a liquid surface is typically described by the incompressible Euler Equations:

∂v v v ∂t + v v v • ∇v v v = -∇P + g g g (1) ∇ • v v v = 0 (2)
where v v v is the 3D velocity of the fluid, P is the pressure, and g g g is the acceleration due to gravity. Researchers in computer graphics solve these equations numerically. Advection is typically handled using a variant of semi-Lagrangian advection [START_REF] Stam | Stable fluids[END_REF] or Lagrangian particles [START_REF] Brackbill | FLIP: A method for adaptively zoned, particle-in-cell calculations of fluid flows in two dimensions[END_REF]. The incompressibility constraint in the second equation is enforced by treating pressure as a Lagrange multiplier and solving a Poisson equation [START_REF] Bridson | Fluid simulation for computer graphics[END_REF]. To track the motion of the surface, one can advect an implicit surface using the level set method [START_REF] Osher | Level set methods and dynamic implicit surfaces[END_REF], track an explicit surface mesh [START_REF] Wojtan | Liquid simulation with mesh-based surface tracking[END_REF], or reconstruct a surface from FLIP particles [START_REF] Yu | Reconstructing surfaces of particle-based fluids using anisotropic kernels[END_REF].

Simulating liquids in this manner has immense advantages: it creates highly detailed motions, including non-linear wave propagation and effects like splashes and overturning waves that cannot be captured by a heightfield alone. The main disadvantage of this approach is its computational complexity; detailed visual effects require very fine grids, which require long computation times and large amounts of memory.

To combat this computational complexity, many researchers proposed adaptive simulation methods [MWN*17] which gain speedups and reduce memory by locally reducing simulation detail where it may not be needed. These methods typically use adaptive spatial data structures like octrees [START_REF] Losasso | Simulating water and smoke with an octree data structure[END_REF], warping grids [IWT*18], or meshes [ATW13; BWHT07; WT08; KFCO06] instead of regular grids. Zhu et al. [ZLC*13] adaptively stretch out the grid cells away from a central point of interest, effectively extending the size of the computational domain and reducing detail away from the grid center. By adaptively reducing computational degrees of freedom, these methods all make trade-offs between visual detail and computational speed. For surface wave simulation, reducing simulation resolution has the unfortunate side effect of removing detailed water surface ripples with short wavelengths.

2D Water Wave Simulation

To reduce computational complexity without removing details at the surface, several researchers investigated the simulation of water waves on top of a flat fluid surface. One approach assumes that the fluid velocity is irrotational, allowing it to be expressed as the gradient of a velocity potential. This assumption significantly simplifies the equation so it can be solved using lower-dimensional degrees of freedom limited to the two-dimensional surface [Tes04; Tes14; CMT*16; JW15; JW17; JSM*18] or even the one-dimensional boundary of the fluid domain [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF]. Another approach to the problem assumes that the water is so shallow that the fluid velocity is effectively constant throughout the depth. This assumption gives rise to tall cell methods [IGLF06; CM11] and the shallow water equations (SWE) [LvdP02; WMT07; CM10; SBC*11] which again lead to a computational speed-up due to the reduction in dimensionality from 3D to 2D. The equivalent sources (i.e., "fundamental solutions") approach of Schreck et al. [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] also enables infinitely large, open fluid domains. The water wavelength is also independent from simulation degrees of freedom, so the method creates high frequency visual details without increasing the simulation resolution.

Despite their computational efficiency and increased visual detail, these specialized wave simulation approaches have significant limitations. They enforce particular constraints on the velocity of the fluid below the surface (e.g. potential flow or constant with respect to depth), so they are not as expressive as fully threedimensional simulations. The potential-flow methods follow Airy wave theory (also known as "linear" wave theory) [START_REF] Airy | Tides and waves[END_REF], so they cannot model non-linear wave propagation. These methods also assume the water surface is a heightfield, so they cannot model breaking waves or flows that make the surface exhibit interesting topology changes.

submitted to EUROGRAPHICS 2022.

Coupling 2D and 3D liquid simulations

Inspired by the strengths of both techniques, several researchers have combined 2D wave simulation methods with 3D fluid solvers. Many methods (ours included) take advantage of the efficiency of 2D solvers to simulate extremely large simulation domains by combining a 3D simulation in the area of interest with a 2D simulation further away. Thuerey et al. [START_REF] Thürey | Animation of open water phenomena with coupled shallow water and free surface simulations[END_REF] On a related note, some researchers proposed methods for coupling a 3D fluid simulation into an existing fluid environment. Nielsen & Bridson [NB11] propose using an existing fluid simulation as a guide shape for a new 3D fluid simulation; English et al. [START_REF] English | Chimera grids for water simulation[END_REF] propose a "chimera grid" strategy for matching fluid simulation variables at domain boundaries; Bojsen-Hansen & Wojtan [START_REF] Bojsen-Hansen | Generalized Non-Reflecting Boundaries for Fluid Re-Simulation[END_REF] use generalized perfectly matched layers to blend the 3D simulation into the surrounding environment; and Stomakhin et al. [START_REF] Stomakhin | Fluxed Animated Boundary Method[END_REF] describe how to compute fluid fluxes through the animated boundary of a fluid simulation domain.

Our method couples the 2D wave solver of Schreck et al. [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] to the boundary of a 3D fluid simulation. This particular surface wave solver has a number of unique advantages over ones used in the past: First of all, unlike the constant-speed waves produced by SWE and Lattice Boltzmann method wave simulations, [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] is based on Airy wave theory and produces waves with physical behaviors closer to those of the 3D fluid simulation. Secondly, contrary to the methods that use 2D grids or meshes to represent the surrounding water surface waves, our combined method allows fluid domains that are literally infinitely large. Finally, because [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] gives us full control over the dispersion behavior of the waves, so we can actually modify the wave solver's behavior to compensate for numerical errors in the 3D fluid simulation, practically eliminating visual artifacts at the domain boundary.

Technical background

3D liquid simulation

Our work couples a 3D fluid simulation to a 2D wave solver. We consider 3D fluid simulators based on a regular grid, which represents velocity and pressure fields as functions discretized on the grid, and also tracks the location of the liquid surface. We have implemented our technique into two different methods: the fluid solver implemented in the Houdini software by SideFX, which appears to use narrow-band FLIP [FAW*16] for advection and extracts the surface from the FLIP particles [START_REF] Yu | Reconstructing surfaces of particle-based fluids using anisotropic kernels[END_REF]; and the solver described by Bojsen-Hansen & Wojtan [START_REF] Bojsen-Hansen | Generalized Non-Reflecting Boundaries for Fluid Re-Simulation[END_REF], which uses semi-Lagrangian advection and represents the liquid surface using the level set method [START_REF] Osher | Level set methods and dynamic implicit surfaces[END_REF].

2D Water wave simulation

We use the method of Schreck et al. [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] for 2D water wave simulation. By assuming that water depth is constant and that waves are not sloped too steeply, they represent the water surface motion with Airy theory [START_REF] Airy | Tides and waves[END_REF]. The water surface can then be decomposed as a set of sinusoidal frequency components of wave number k -related to the wavelength λ by k= 2π λ . Each component is a function of the form a cos(kxωt), where a is the amplitude and ω is the angular frequency (ω = 2π f where f is the time frequency). The heightfield is a discrete sum of these frequency components:

η(x x x,t) = Re(u(x x x,t)) and u(x x x,t) = ∑ k p k (x x x)e -iω k t , (3)
where p k is the complex spatial field corresponding to wavenumber k. Each p k is a solution of the Helmholtz equation:

∇ 2 p + k 2 p = 0.
Water waves are dispersive, which means that their angular frequency depends on their wavenumber. This dispersion relation makes their phase speed cp depend on wavelength as well, via the relationship cp = ω k . In particular, gravity waves follow the theoretical relation ω k = gk tanh(kh) with g being the gravity constant and h being water depth. We refer interested readers to [START_REF] Johnson | A Modern Introduction to the Mathematical Theory of Water Waves[END_REF] for a more complete introduction to the linear theory of water surface waves.

Schreck et al. propose to use the Equivalent Sources Method (ESM) to compute the spatial components p sc

k of the field u sc scattered by the obstacle. The idea of ESM is to approximate a field by a set of n fundamental solutions:

p sc k (x x x) = n ∑ i=1 a k,i φ k (x x x -x x x i ) (4)
where the a k,i are complex amplitudes and φ k (x x xx x x i ) is a fundamental solution of the Helmholtz equation, or source, located at x x x i . More precisely, this source is a point emitting circular waves and can be explicitly computed as

φ k (x x x -x x x i ) = -i 4 H (2) 0 (k||x x x -x x x i ||) where H (2)
0 is the 0 th order Hankel function of the second kind. The sources are uniformly set on an offset curve inside the obstacle.

To compute how an incoming wave u in is scattered by an obstacle, Schreck et al. define the boundary condition at the border of the simulation domain:

u sc (x x x,t) = -u in (x x x,t) (5)
for all times t and for all x x x on the boundary. The functions e -iω k t k (defined in Equation 3) are orthogonal, so the boundary conditions are enforced separately for each wave number k. Combining Equation 3 and Equation 5 gives the boundary condition for each k: The height of the surface can finally be computed using Equation 3as

p sc k (x x x) = -p in k (x x x) ∀x
η(x x x,t) = Re ∑ k ∑ i a k,i φ k (y y y j -x x x i )e -iω k t + u in (x x x,t) . ( 8
)
To animate aperiodic wave motions, Schreck et al. replace each amplitude a k,i by a time-dependent function a k,i (t). The latter is a linear interpolation of a set of discrete amplitudes a h k,i = a k,i (t T ) , where t T is the time after T timesteps. In this case, p k (x x x,t) = ∑ i a k,i (tx x x/cg)φ k (xx i ), where cg = ∂ω ∂k is the wave's group speed. This method animates time-dependent waves at the cost of solving a system similar to Equation 7 at each time step.

Note that, unlike 2D wave simulations that discretize the heightfield onto a grid, this method analytically defines the wave heights everywhere on R 2 . This makes it trivial to represent infinitely large bodies of water and to represent high-frequency ripples without using an adaptive grid representation. In a sense, it maximizes the benefits of 2D wave simulations (efficient computation of large heightfield domains) and is a good candidate for coupling to a 3D simulation for locally enhanced details.

Objectives and Overview

The goal of this project is to propose a tool to efficiently simulate a detailed scene on an arbitrarily large water surface. For this we want to use a computationally expensive full 3D fluid simulation only where really needed and a cheaper 2D surface wave propagation using sources everywhere else (see Figure 2).

As waves propagate outward from a disturbance, the law of energy conservation causes their amplitudes to decrease, and these small-amplitude waves are well-described by linear wave theory. We assume that the 3D fluid simulation is large enough not only to capture all 3D non-linear effects, but also such that the waves leaving the 3D simulation domain can be approximated using Airy theory -the water depth and wavelengths are significantly larger than the wave height-(see Figure 3).

The surface in these regions can then be approximated using the equivalent sources method.

The 3D fluid simulation represents only part of a larger surface, so we need to prevent any erroneous internal reflections on the border of the simulation. We use either a perfectly matched layer (PML) method or an increasing viscosity to damp to waves before they reach the end of the simulated domain.

We consider a fluid surface represented by domain D that can be divided into two subdomains D 3d and D 2d (see Figure 2). D 2d corresponds to the regions where the surface of the fluid can be represented by linear (Airy) waves and displayed on a heightfield. D 3d is the domain where more complex behavior is needed (3D water droplets or breaking waves, for example).

In Section 5, we explain how to ensure a seamless transition between D 2d and D 3d for linear waves radiating from D 3d . We build upon [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF] and use the state given by the 3D liquid simulation as boundary conditions for the 2D solver. We also discuss how to modify the boundary conditions of the 3D fluid simulation to prevent artificial internal reflections, and to send in new waves from the 2D simulation.

When a physical wave propagates from one medium to another, it can exhibit reflections and refractions. We observe similar effects when coupling two different simulation wave representations; the subtle differences in discretization affect the way that waves propagate and can lead to similar artificial wave reflections and refractions when a wave exits one simulation and enters another. This is a common problem with spatially adaptive discretizations, where high-frequency waves cannot be represented by a coarser discretization and cause an internal reflection [START_REF] Vichnevetsky | Propagation through numerical mesh refinement for hyperbolic equations[END_REF]. Artificial damping and dispersion from the 3D simulation can significantly change the behavior of water surface waves, especially when compared to the theoretically exact behaviors described by linear wave theory. These differences cause significant artificial refraction artifacts at the boundary between the 3D and 2D simulations, as illustrated in Figure 6. In Section 6, we explain how to estimate the effective wave behaviors from any 3D fluid simulation; we then use these imperfect behaviors in our 2D wave solver to remove artificial refractions.

Coupling method

To obtain a seamless transition between the two regions, we must take care to set boundary conditions of each solver so that waves cleanly propagate from one domain to the other. We first describe our novel method for transferring waves from the 3D simulation into the 2D simulation. Afterward, we describe wave-blending techniques for transferring waves back from 2D to 3D.

Boundary conditions for the 2D simulation

In Section 3, we explained how Schreck et al. approximate the waves dispersed by an obstacle. We use a similar approach to approximate the wave field radiating from a 3D domain with a set of 2D sources.

As explained in Section 4, we assume that the 3D surface respects Airy theory near the boundary and can be represented as a heightfield. Let's call u 3d the field representing the surface of the submitted to EUROGRAPHICS 2022.

water from the fluid simulation in D 3d and u 2d its approximate representation in D 2d using a set of n fundamental sources centered in {x x x i } i=1...n :

u 2d (x x x,t) = ∑ k n ∑ i=1 a i,k (t)Φ(x x x -x x x i )e -iω k t (9)
We then use a method similar as the one described in Section 3 to compute the complex amplitudes a i,k such that the transition between D 3d to D 2d is as seamless as possible. For this we simply set the 2D heightfield equal to the 3D surface:

u 2d (x x x,t) = u 3d (x x x,t)∀x x x ∈ δD, ( 10
)
where δD is the boundary between the two domains. By resolving the boundary condition separately for each wave number k, we have:

p 2d k (x x x,t) = p 3d k (x x x,t)∀x x x ∈ δD (11) with p 2d k (x x x,t) = ∑ n i=1 a k,i (t)φ k (x x x -x x x i ).
In our examples, we set δD to be an ellipse bounding the 3D simulation, though this choice of domain boundary is up to the user. The sources (red crosses in Figure 3) are placed at an offset distance from δD. To satisfy this boundary condition, we discretize δD as a set of m boundary points y y y j . As done in [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF], we uniformly sample the sources and boundary points with respect to the Nyquist sampling criterion. This gives us a linear system to solve for the amplitudes a k,i (t): This system requires a frequency decomposition of the 3D heightfield at each boundary point y y y j . Setting H j (t) = u 3d (y y y j ,t) = ∑ k p 3d k (y y y j ,t)e -iω k t , we apply a Fast Fourier Transform (FFT) over a small time window of size T = w∆t centered at time τ:

FFT ([H j (τ - T 2 ), . . . , H j (τ + T 2 )]) = [p 3d f0 (τ), . . . , p 3d fw (τ)], (13)
where [ f 0 , . . . , fw] is a range of discrete frequencies separated by ∆ f = 1/T . Note that the set of wave numbers used in Equation 9 are the ones corresponding to these frequencies: [k 0 , . . . , kw] with f ξ = 2πω(k ξ ). The precision of the spectral decomposition depends on the window size, so a compromise must be found between precision in frequency and in time. Our implementation uses w = 256 for wavelengths ranging between 0.13 and 33 meters. In practice we only keep a few of them that recover the wavelengths created by the simulation.

In our examples, we usually use a dozen of wavelengths ranging from 0.6m to 10m.

Boundary conditions for the 3D simulation

Next we discuss how to inject waves from the 2D region into the 3D fluid simulation. Though similar techniques already exist in the computer animation literature [BW16; SS17], we describe our approach here for completeness. Our results use two different 3D fluid simulation solvers, so this section will discuss how we enforce boundary conditions for both types of simulations.

For the semi-Lagrangian level set simulation, we use the generalized perfectly-matched layers approach of Bojsen-Hansen & Wojtan, treating our 2D wave solver as the background flow. For the narrow band FLIP simulation, we add a boundary layer with gradually increasing viscosity near the boundary, and we blend the positions and velocities of the FLIP particles close to the boundary. Though this gradual blending is not as efficient as the PML-based approach, we found it sufficient to inject the 2D waves and dampen artificial internal reflections.

In any case, the 3D boundary conditions require us to compute the 3D velocity and displacement fields induced by the 2D waves. For a planar wave η(x x x,t) = a cos(kxωt) (with x x x = (x, y, z) and z being the vertical direction), Airy wave theory gives us an expression for the velocity and particle displacement below the surface [START_REF] Stomakhin | Fluxed Animated Boundary Method[END_REF]:

vx(x x x,t) = ẋ(t) = aωe kz cos(kx -ωt) (14) vz(x x x,t) = ż(t) = aωe kz sin(kx -ωt) (15)
Similarly, the displacement can then be computed as:

ξx(x x x,t) = -ae kz sin(kx -ωt) (16) ξz(x x x,t) = ae kz cos(kx -ωt) ( 17
)
We adapt these formulae for a wave field defined by sources. For a surface defined by a source η(x x x,t) = a φ k (x x xx 0 )e -iωt centered in x x x 0 , we obtain the following fields:

vr(x x x,t) = ṙ(t) = aωe kz φ k (r)e -iωt (18) vz(x x x,t) = ż(t) = aωe kz φ k (r)e -i(ωt+ π 2 ) (19) ξr(x x x,t) = -ae kz φ k (r)e -i(ωt-π 2 ) (20) ξz(x x x,t) = ae kz φ k (r)e -iωt (21)
submitted to EUROGRAPHICS 2022.

with r = |x x xx x x 0 |.

The waves coming toward the 3D domain can be for example coming from reflection on an obstacle (represented by another set of sources as described in [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF]), or background planar wind waves. When injecting back these waves -represented by additional background heightfield η bg -into the fluid simulation, we need to account for them by enforcing the boundary condition η 2d (x x x,t) + η bg (x x x,t) = η 3d (x x x,t). Figure 7 shows how this 2D-to-3D coupling looks for circular waves entering a level-set simulation.

Matching wave behaviors

Every numerical method for simulating fluids adds some amount of numerical error to the solution. When coupling two simulations which use two different numerical methods, or even coupling two simulations using the same method but with different resolutions, visible artifacts will appear. These errors usually result in different wave speeds or damping rates, known as numerical dispersion and numerical viscosity. If the two simulations have different apparent values for the dispersion relation ω(k), waves will visibly refract at the boundary between the simulations.

In our case, the 2D wave solver has a hard-coded dispersion relation, so it actually exhibits no numerical dispersion errors. On the contrary, in the 3D simulation, the grid discretization significantly affects the behavior of the fluid -the effective dispersion law depends on the parameters of the simulation. Instead of trying to make the 3D simulation more accurate, we propose to intentionally add errors to the hard-coded 2D dispersion law, so that it better matches the 3D simulation. We do this by estimating the effective dispersion relation of the 3D fluid, and then using this estimated law in place of the analytic ω k in ESM. Rather than analytically trying to derive the dominant error for each specific method, we chose an easy to run method which is based on numerical experiments and agnostic to the kind of 3D solver.

We use the following experiment to estimate a 3D fluid simulation's frequency response f = 2πω to a given wavelength λ = 2 π k : we initialize a standing wave of wavelength λ centered within a a virtual wave tank of size l = 2λ. We then begin the simulation and measure the height H(t) of the point in the center of the tank over time. Next, we apply an FFT on H(t) and identify the frequency of the spectrum with the highest amplitude to compute f and ω = f /2π. Figure 4 illustrates this experimental setup.

By repeating the experiment for different wavelengths we can plot the dispersion law ω = ωs(k) and compare it to the theoretical one ω = gk tanh(kh). Figure 5 plots our findings for both the FLIP simulation and levelset simulators for a fixed simulation resolution and time step ∆t = 0.033. For the levelset simulation we use a cell size of 0.0625m, and for FLIP we use a particle separation of 0.075 with the Houdini parameters particle radius scale equal to 1.2 and a grid scale equal to 2. Figure 5 shows the dispersion laws we obtained for these two simulations compared to the theoretical gravity waves dispersion.

The FLIP simulation substantially increases the wave speed and adds noise. This discrepancy is illustrated in Figure 6 (a) where the difference between the dispersion laws between the 3D FLIP simulation (left) and 2D wave simulation (right) leads to sharp wave refraction at the interface. To compensate for this increased frequency, we found that simply re-scaling the theoretical law by a constant coefficient d ≈ 1.4 is sufficient to practically eliminate artificial wave refractions, as illustrated in Figure 6 (b).

Although the level-set simulation is close enough to theory for the refraction to be far less noticeable, it still exhibits slightly slower waves than Airy wave theory predicts. This difference in wave speed can create inconsistencies along the border when injecting 2D waves into the 3D simulation, resulting in spurious waves as shown in Figure 7 (left).To correct this, we directly use data from the experiments to associate a more accurate frequency (ω = 3.13) to the injected wavelength (λ = 6m).

Although our examples apply simple approaches to match the sources dispersion law with the data from the 3D simulation, any data-fitting approach could work; our only requirement is that ωs(k) be differentiable, so we can compute the group speed ∂ω

Post-processing details

Surfaces reconstructed from particles tend to exhibit highfrequency noise, which makes the 3D FLIP simulations look substantially different from the perfectly smooth analytical solutions given by the 2D wave solver. To make the surfaces look similar, we add Perlin noise [START_REF] Perlin | An image synthesizer[END_REF] to the 2D wave surfaces:

η noise (x x x,t) = Re(a(x x x)e -iωt )
submitted to EUROGRAPHICS 2022. where a(x x x) is a complex amplitude whose real and imaginary parts are two noise functions with the same parameters. We also linearly blend the surface functions between the 2D and 3D domains, as illustrated by the purple interpolation curve in Figure 3.

Figure 8 shows an example with and without these postprocessing steps.

Results

Please see our supplementary video for animations of our results. We first verify our coupling method in the absence of complicating factors by replacing the 3D fluid simulation with a 2D heightfield with a known analytical solution. Figure 9 shows how we generated a few circular ripples within the dashed ellipse, and solved for the fundamental solutions which best matched the behavior at the simulation boundary.

Figure 10 illustrates the consequences of different dispersion laws, comparing the deep water dispersion law ω = gk tanh(kh) to the shallow water one ω = k √ gh. The two laws exhibit very different behaviors; the velocity of the waves is similar for all frequencies for shallow water, but they are dispersive in deeper water. Deep water waves also have a phase velocity that differs from their group velocity, causing waves to fade away at the edge of the group (see the video). This is not the case for shallow water for which the group velocity and phase velocity are equal.

Our method is practically agnostic to the type of 3D fluid simulation used. To illustrate this benefit, we show examples using two different kinds of simulation. Figure 11 shows coupling with a semi-Lagrangian levelset method to simulate a few drops falling into water with ripples propagating far away from the original 3D domain. Figure 1 shows a FLIP simulation with a dragon springing out of water, with 2D waves propagating outward into an unbounded fluid domain.

Our method's choice of 2D wave propagation gives an analytical solution defined everywhere in the 2D domain, so the surface can be adaptively computed at render-time to represent very large scenes as well as more detailed close-ups. Figure 12 shows a rotating water wheel with the camera zooming out to view waves propagating far away. Figure 13 shows a zoom on the boundary between the 2D and 3D domain of the example shown in Figure 1, using a finer grid to obtain a smooth and detailed close-up view.

Finally, we compare our coupling method with a full 3d fluid simulation in Figure 14. The visual results are very similar, but our method takes only a fraction of the time to compute: our coupling method takes around 30 seconds per frame (20s for generating the 3D simulation at the center and 11s to compute the 2D wave propagation) and the full simulation takes more than 2 minutes.

Table 1 shows the times for our examples as well as the number of sources, boundary points, and resolution of the display grid. Note that, as described in [START_REF] Schreck | Fundamental solutions for water wave animation[END_REF], we are using a different number of sources for each wavelength. For example in Figure 14 (left), the biggest wavelength is 45 meters and the corresponding number of sources is 16. The smallest wavelength is 1.5m with 566 sources.

Discussion about previous work: The recent article from Huang et al. [HQT*21] achieves impressive results by coupling a FLIP simulation with a Boundary Element Method (BEM). Our approach uses fundamental solutions instead of BEM to couple the 2D and 3D domain, leading to significantly different practical consequences. Although both methods assume linear waves outside of the 3D domain, our method is strictly limited to height field waves, while the BEM used by Huang et al. is not. On the other hand, our method is distinct from meshing-based approaches in that it can handle infinite domains and can simulate arbitrarily high-frequency waves without any spatial adaptivity in the simulation domain. Our work also offers preliminary insights into how to address the coupling errors at simulation boundaries, which may be applicable to methods like Huang et al. Figure 9: A controlled experiment to verify our approach: Our coupling method produces 2D waves that radiate away from another 2D domain within the dashed ellipse.

Discussion and Conclusion

We first note that our method only approximates the outgoing waves from a 3D simulation, so, despite our efforts to reduce visual artifacts, there will still be errors when the two simulations do not exactly match up on the boundary. One potential error source is the finite time source on our FFT used to estimate the 2D wave amplitudes; the choice of time window trades off between accuracy in spatial and temporal frequencies (also known as the Fourier uncertainty principle). The 2D wave solver also requires a choice of how densely to discretize the boundary and the space of wavenumbers; insufficient discretization here can also cause mismatches in the outgoing waves.

We chose to match the wave speeds across simulation boundaries by approximating the effective dispersion relation of the 3D solver. However, damping can also differ across simulations. Our initial investigations show that this is not nearly as strong of an effect as the mismatched wave speeds, but it is possible to develop a method for matching the effective damping rate in the future.

Finally, we note that our implementation does not take full advantage of the equivalent sources method. We used a parallel imsubmitted to EUROGRAPHICS 2022. plementation on CPU but due to constraints in implementing the method within Houdini, we did not parallelize the 2D wave solver on the GPU. We believe the performance of our method would be much faster with a more optimal implementation.

To conclude, we presented a method that adds radiating 2D waves to a localized 3D fluid simulation by coupling it with the equivalent sources method. The method allows efficient simulation of large open domains with arbitrarily high spatial frequencies, and it produces more physically plausible wave behaviors than other 2D solvers which assume constant wave speeds. Finally, we introduced a novel method for discovering and matching the effective dispersion behaviors of an arbitrary 3D fluid simulation, leading to a practical elimination of refraction artifacts at the boundary between simulations.

Figure 2 :

 2 Figure 2: Side view of our configuration: waves radiating from the 3d simulation are propagated as 2d wave on the surface.

Figure 3 :

 3 Figure 3: Top: side view of the 3d fluid simulation and 2d heightfield. Bottom: Top view of the configuration. Boundary points (purple dots) are placed on the border region of the 3d fluid simulation where the waves can be approximated as linear. The sources (red crosses) are placed at an offset distance from them. The surface of the fluid and the 2d heightfield are linearly interpolated (interpolation function in purple) to smooth the transition (between the dotted lines (top) and in the gray region (bottom).

  ∂k for Schreck et al. [SHW19]'s wavelet propagation.

Figure 4 :

 4 Figure 4: Set up of our simulation of a standing wave. The tank measures exactly two times the wavelength λ. We record the height in the middle of the tank.

Figure 5 :

 5 Figure 5: The recorded dispersion laws for a FLIP simulation and a levelset simulation compared to the theoretical relation of gravity waves.

Figure 6 :

 6 Figure 6: A FLIP simulation on the left side of the dotted line coupled with fundamental sources using the theoretical gravity wave dispersion law (a) and a corrected law to match better the recorded law shown in Figure 5 (b). The left picture shows a spurious refraction effect due to the difference of speed of the waves between the 3d fluid simulation and the 2d wave propagation. This effect is corrected on the right by using the new law.

Figure 7 :

 7 Figure 7: Injecting the waves from a source inside a levelset simulation using Bojsen-Hansen & Wojtan's [BW16] method. Left: The difference of velocity between the theoretical dispersion law used by the source and that of the simulation causes spurious waves. Right: Correcting the law of the source removes these artifacts.

  submitted to EUROGRAPHICS 2022.

Figure 8 :Table 1 :

 81 Figure 8: Our method with both linear blending and noise (left), without the linear blending (center-left), without noise (center-right) and without any post-processing (right). # wavelengths # sources # bound points Surface grid resolution time / frame (2D) time / frame (3D) Figure 1 17 7274 2000 450x600 141s 222s Figure 6 7 142 600 300x150 0.6s 3.7s Figure 9 7 140 300 200x200 0.4s -Figure 11 6 368 900 300x400 2.4s 1.2s Figure 12 18 3839 1200 150x200 6.8s 26s Figure 14 (left) 10 2352 2000 300x300 11s 20s Figure 14 (right) -----138s Table 1: Number of wavelengths, number of sources and boundary points, resolution of the 2D grid for our examples, as well as the time per frame for computing the 2D waves and the time per frame for the original 3D fluid simulation.

Figure 10 :

 10 Figure 10: Left: a full levelset simulation of one drop in the water in a deep pool. Right: a levelset simulation (in the center ring) coupled with 2D wave propagation using shallow water dispersion (top) and deep water dispersion (bottom). The dispersion law used by the 2D simulator has a significant affect on the resulting wave pattern and coupling, and should be chosen to match the 3D simulation.

Figure 11 :

 11 Figure 11: Drops falling into water. Our coupling method used with a levelset simulation (inside the dashed ellipse).

Figure 12 :

 12 Figure 12: Zooming out from a water wheel while still displaying waves radiating from the FLIP simulation in the center.

Figure 13 :

 13 Figure13: A close-up view of the example in Figure1, using the same wave sources but a finer resolution grid for the displayed surface.

Figure 14 :

 14 Figure 14: An armadillo being cruelly dipped into water. (left) Our coupling method with a FLIP simulation at the center (dashed ellipse) and 2d waves everywhere else. (right) Large FLIP simulation.

  it with Equation 4, a set of m linear equations with the amplitudes a k,i i=1..n as unknown is obtained for each k:

on the boundary. (6) By discretizing Equation 6 at a set of m points y y y j j=1..m and submitted to EUROGRAPHICS 2022. combining n ∑ i=1 a k,i φ k (y y y jx x x i ) = p in k (y y y j )∀y y y j . (7) Schreck et al. then solve this using a least squares optimization.

Acknowledgments: We wish to thank the anonymous reviewers and the members of the Visual Computing Group at IST Austria and MFX Team at INRIA for their valuable feedback. This research was supported by the Scientific Service Units (SSU) of IST Austria through resources provided by Scientific Computing. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme under grant agreement No. 638176.


Paul Dorbec
email: paul.dorbec@unicaen.fr

Fatma Kaci
email: fatma.kaci@univ-biskra.dz

Disjoint Maximal Independent Sets in Graphs and Hypergraphs

In this paper, we consider the question of the existence of disjoint maximal independent sets (mis) in graphs and hypergraphs. The question was raised in the 1970's independently by C. Berge and C. Payan. They considered the question of characterizing the graphs that admit disjoint mis, and in particular whether regular graphs do. In this paper, we are interested in the existence of disjoint mis in a graph or in its complement, motivated by the fact that most constructions of graphs that do not admit disjoint mis are such that their complement does. We prove that there are disjoint mis in a graph or its complement whenever the graph has diameter at least three or has chromatic number at most four. We also define a graph of chromatic number 5 and diameter 2 which does not admit disjoint mis nor its complement.

As our work was first motivated by a more recent work on disjoint mis in hypergraphs by Acharya (2010), we also consider the question of the existence of disjoint mis in hypergraphs. We answer a question by C. Jose and Z. Tuza (2009), proving that there exists balanced k-connected hypergraphs admitting no disjoint mis.

Introduction

We first introduce common definitions to graphs and hypergraphs, classical on graphs, that will be used throughout the rest of the paper. Specificities of the definitions for hypergraphs are discussed in Section 3

In a graph or a hypergraph H = (V, E), the neighborhood N (x) of a vertex x of V is the set of all the vertices that belong to a common edge with x. More generally, the neighborhood of a subset S ⊆ V of vertices is defined as N (S) = x∈S N (x). A subset S of V is an independent set of H if no two vertices of S belong to a same edge. An independent set S is maximal (inclusionwise) if there exist no superset S S that is also independent. Observe that a maximal independent set of H is also a dominating set of H, i.e. every vertex in H has a neighbor in the set S. Maximal independent sets are also known as independent dominating sets, which have been extensively studied in the literature, we refer to the survey by W. Goddard and M. Henning [START_REF] Goddard | Independent domination in graphs: A survey and recent results[END_REF].

In this paper, we are interested in the existence of disjoint maximal independent sets (disjoint mis) in graphs and hypergraphs. Observe first that all graphs do not contain disjoint mis. An example of a graph with no disjoint mis can be made by adding an isolated vertex to any connected graph. A connected example is depicted in Figure 1: any mis of G 2 contains at least three among the four degree 1 vertices, so any two mis intersect in at least two degree 1 vertices. Another known infinite family of graphs admitting no disjoint mis is formed by the coronas of odd cycles. Those graphs are obtained by adding a degree one vertex adjacent to every vertex of an odd cycle (see Figure 2). In these graphs, any maximal independent set contains at most half the vertices of the odd cycle, and thus more than half the degree one vertices. Therefore any two mis intersect in at least one of the degree one vertices. Based on this observation, O. Schaudt [START_REF] Schaudt | On disjoint maximal independent sets in graphs[END_REF] gave a polynomial time algorithm that either computes disjoint mis in a graph or returns an induced corona of an odd cycle. Note however that having no induced corona of a cycle is not a necessary condition to have disjoint mis. For example, the graph G 1 in Figure 1 does contain an induced corona of a 3-cycle, as well as disjoint mis. Actually, deciding whether a graph admits disjoint mis is a NP-complete problem, as shown by M. Henning et al. in [START_REF] Michael | Remarks about disjoint dominating sets[END_REF].

G 1 G 2
C. Berge proposed in conferences (as mentionned in [START_REF] Bondy | Graph theory with applications[END_REF][START_REF] Cockayne | Disjoint independent dominating sets in graphs[END_REF]) the following conjecture, which was independently proposed by C. Payan in [START_REF] Payan | Sur une classe de problèmes de couverture[END_REF]: Though, C. Payan [START_REF] Payan | A counter-example to the conjecture "every nonempty regular simple graph contains two disjoint maximal independent sets[END_REF] has found a family of counter-examples to the above conjecture, the smallest graph of which contains 630 vertices, all of degree 280.

However, the conjecture was proven to be true for (n -k)-regular graphs where 1 ≤ k ≤ 7 (see [START_REF] Cockayne | Disjoint independent dominating sets in graphs[END_REF]), k ≤ 10 (see [START_REF] Payan | Sur une classe de problèmes de couverture[END_REF]) and for k < 2 √ 2n -2 (see [START_REF] Erdős | Disjoint cliques and disjoint maximal independent sets of vertices in graphs[END_REF]). Morever C. Payan proved also the conjecture in particular for claw-free regular graphs [START_REF] Payan | Coverings by minimal transversals[END_REF].

2 About disjoint MIS or disjoint maximal cliques.

We consider here the question of the existence of disjoint mis in a graph or in its complement, or equivalently the question of a graph admitting disjoint mis or disjoint maximal cliques (i.e. disjoint maximal complete subgraphs). Observe that the simple examples of graphs not admitting disjoint mis shown above all admit disjoint maximal cliques. First, if a non trivial graph contains an isolate vertex v, then v itself is a maximal complete subgraph disjoint from any other clique in the graph. In the family of coronas of odd cycles, there are many maximal cliques formed by a vertex of degree one with its single neighbor, and all are disjoint. Similarly, we do not detail C. Payan's construction here, but it does contain disjoint maximal cliques.

The study of existence of disjoint mis in a graph and in its complement was already initiated by E. Cockayne and S.Hedetniemi in [START_REF] Cockayne | Disjoint independent dominating sets in graphs[END_REF], though with a different perspective, since Conjecture 1 was not disproved yet. They conjectured that if a graph admits disjoint mis, then its complement does too, which is disproved by the examples mentionned above.

This section is thus motivated by the following question, which we explore further :

Question 2 Let G be a non trivial graph. Is it true that either G or G admits disjoint mis? Equivalently, is it true that either G admits disjoint mis, or G admits disjoint maximal cliques?

About the diameter

One first direction to consider the question is in regard of the diameter of the graph. Note that the only graphs of diameter one are complete graphs, where each vertex is a mis in the graph. Thus non trivial graphs with diameter one admit disjoint mis. For larger diameters, we get the following:

Theorem 3 If G is a graph of diameter at least 3, then G contains disjoint maximal cliques.
Proof: Let G be a graph of diameter at least 3. Consider two vertices u and v at distance at least three in G, and two maximal cliques, one containing u and the other containing v. Since u and v are at distance more than two, they have no common neighbors, and the two cliques may not share a vertex. They are thus disjoint.

2

Observe in particular that the above theorem covers the (easy) case when the graph is not connected.

Corollary 4 If a graph admits neither disjoint mis, nor disjoint maximal cliques, it must be of diameter two.

About the chromatic number

Our next result uses the chromatic number of the graph as a parameter. Recall that a proper coloring of a graph G is an assignment of colors to the vertices so that no adjacent vertices receive the same color. The chromatic number χ(G) of the graph G is the minimum number of colors in a proper coloring of G. We consider the following greedy algorithm for producing a proper coloring of a graph: take the vertices in some order, and assign to each vertex the least color not used among its already colored neighbors.

Observation 5 For every graph, there exists an ordering of the vertices so that the greedy algorithm produces an optimal coloring.

Proof: Let G be a graph with chromatic number k, and c be an optimal coloring of G, with colors 1 to k. Consider any order of the vertices which begins with all vertices receiving color 1 in c, then all vertices receiving color 2, and so on until color k.

Let c be the coloring of G obtained by applying the above greedy algorithm with this ordering. Note that for a vertex v, the color c (v) may differ from c(v) (for example if c(v) = 2 but v has no neighbor colored 1). However, we can ensure inductively that for all v, c (v) ≤ c(v): by induction, all the already colored neighbors of a vertex v received colors less than c(v), so the least available color is at most c(v). Therefore, the coloring obtained by the greedy algorithm from that ordering uses at most k colors.

2

We proceed with a series of lemmas that will allow us to show that any graph with chromatic number at most 4 admits disjoint mis or disjoint cliques.

Lemma 6 If G is a connected graph with chromatic number χ(G) = 2, then G contains disjoint mis.

Proof: Consider a two-coloring of G, and denote by V 1 and V 2 the set of vertices colored 1 and 2, respectively (in other words, V 1 and V 2 are the partite sets of the bipartite graph G). Obviously, V 1 and V 2 are independent sets. Since the graph is connected, every vertex colored 2 has a neighbor colored 1, and reciprocally. So V 1 and V 2 are also maximal, they thus are disjoint mis.

2

Lemma 7 If G is a connected graph with chromatic number χ(G) = 3, then G contains disjoint mis or disjoint maximal cliques.

Proof: Consider an ordering of the vertices so that the greedy coloring of G uses three colors, and denote V 1 , V 2 and V 3 the sets of vertices attributed the corresponding colors by the greedy algorithm. Note that by definition, every vertex in

V 2 or V 3 has a neighbor in V 1 , So V 1 is necessarily a mis. Suppose first that every vertex in V 1 is adjacent to at least one vertex in V 2 .
Then V 2 is a dominating set in G, and thus V 1 and V 2 are disjoint mis, as required. Suppose next that every vertex in V 1 is adjacent to at least one vertex in V 3 . Then V 3 dominates V 1 , and it may be completed into a mis S using only vertices from V 2 . Then V 1 and S are disjoint mis, as required.

We now assume that both the above assumptions are false, i.e. there exists a vertex u in V 1 that has no neighbor in V 2 , and a vertex w in V 1 that has no neighbor in V 3 .

Since neither u nor w are isolated, there exists a vertex u ∈ V 3 adjacent to u and a vertex w ∈ V 2 adjacent to w. Since u has no neighbor in V 1 ∪ V 2 and u has no neighbor in V 3 , the complete subgraph {u, u } is maximal.

Similarly, since w has no neighbor in V 1 ∪ V 3 and w has no neighbor in V 2 , the set {w, w } is a maximal clique. Therefore, the sets {u, u } and {w, w } form disjoint maximal cliques.

2

Lemma 8 If G is a connected graph with chromatic number χ(G) = 4, then G contains disjoint mis or disjoint cliques.

Proof: Consider an ordering of the vertices so that the greedy coloring of G uses three colors, and denote V 1 , V 2 , V 3 and V 4 the sets of vertices attributed the corresponding colors by the greedy algorithm.

Similarly as in the proof of the previous lemma, if there exists t ∈ {2, 3, 4}, such that every vertex of V 1 has a neighbor in V t , then V 1 and a mis containing V t form disjoint mis. We thus assume that for every t ∈ {2, 3, 4}, there exists a vertex of V 1 that has no neighbor in V t .

Suppose some vertex u ∈ V 1 has no neighbor in two other partite sets (e.g. V 2 and V 3 ), and denote by V t the only partite set where u has neighbors (in the example, V 4 ). Since u is not isolated in G, u necessarily has a neighbor u in V t . Since u has neighbors only in V t , and V t is an independent set, {u, u } is a maximal clique in G. Let w ∈ V 1 denote a vertex with no neighbor in V t . Then any maximal clique C w containing w contains neither u nor u , and so C w and {u, u } are disjoint maximal cliques. Now suppose that neither of the above conditions hold, and thus that there are three distinct vertices u, v and w in V 1 that do not have neighbors in V 2 , V 3 and V 4 respectively. Let C u , C v and C w denote three maximal cliques containing u, v and w respectively. Note that by our assumption, each of C u , C v and C w contain at most three vertices. If any two such cliques are disjoint, the graph G contains disjoint maximal cliques, so let us assume they pairwise intersect.

Since u has no neighbor in V 2 and v has no neighbor in V 3 , C u and C v may only intersect in some vertex from V 4 , say w . Similarly, C u and C w must intersect in some vertex v in V 3 , and C v and C w intersect in some vertex u in V 2 . A picture of the current identified subgraph is provided in Figure 3. Note that since u has no neighbor in V 2 , v has no neighbor in V 3 and w has no neighbor in V 4 , we necessarily have C u = {u, v , w }, C v = {u , v, w } and C w = {u , v , w}. Moreover, {u, u }, {v, v } and {w, w } are independent sets. Let I u be a maximal independent set containing both u and u , and I v a maximal independent set containing both v and v . If I u and I v are disjoint, the lemma holds, thus assume by way of contradiction that there exist a vertex z in the intersection I u ∩ I v , and let C z be a maximal clique containing z.

V 1 V 2 V 3 V 4 u v w u v w
If C z does not contain w , then since u and v are not incident to z, C z and C u are disjoint, and the lemma holds.

Thus assume C z contains w . Then necessarily, w / ∈ C z . But then since the vertices u and v are not adjacent to z either, C w and C z are disjoint, and the lemma holds.

2

As a corollary, we get the following result.

Theorem 9 If G is a connected graph with chromatic number χ(G) ≤ 4, then G contains disjoint mis or disjoint cliques.

So we now can infer that the only graphs that may have neither disjoint mis nor disjoint maximal cliques have diameter two and chromatic number at least 5. Actually, we found a family of split graphs that have this property, which we describe now:

Let K p be a complete graph with vertex set {x 1 , x 2 , ..., x p }, where p > 4. G p is a graph obtained from K p by connecting any vertex z ij of Z = {z ij | i, j ∈ {1, 2, . . . , p} i < j} to any vertex x r of K p , where r / ∈ {i, j}. The graph G p has p + p 2 vertices, diameter 2 and chromatic number p (see Figure 4 for p = 5). Proof: In the graph G p , a maximal independent set may contain at most one vertex of the complete subgraph K p . It is thus either the set Z, or a set

S i = {x i } ∪ {z kl | i ∈ {k, l}} for some 1 ≤ i ≤ p.
No two of these independent sets are disjoint for p ≥ 4. In particular, S i and S j intersect in z i,j . Thus the graph G p admits no disjoint mis. On the other hand, a maximal clique in G p contains at most one vertex in Z, and thus is either K p or a clique C ij of the form {z ij } ∪ {x k | k / ∈ {i, j}} for 1 ≤ i < j ≤ p. Any two such cliques intersect in at least p -4 vertices, and thus no two are disjoint for p ≥ 4. This concludes the proof. 2

This implies that both above results are best possible in the way they are defined. However, the question of characterization of graphs with no disjoint mis nor disjoint maximal cliques remains open.

Hypergraphs

The study of disjoint mis in hypergraphs was started by B. Acharya in [START_REF] Acharya | Domination in hypergraphs II. New directions[END_REF] where he asks for a characterization of hypergraphs having disjoint mis [1, Problem 2]. Recall that we here use for definition of an independent set in a hypergraph a set where no two vertices belong to the same edge. With this definition, a set in a hypergraph is independent if and only if it corresponds to an independent set in the graph made from the hypergraph by replacing every edge by a complete subgraph. Though, we are interested in a conjecture from B. Jose and Z. Tuza [START_REF] Bibin | Hypergraph domination and strong independence[END_REF] about balanced hypergraphs, for which the equivalence may not be used.

A cycle in a hypergraph is an alternating sequence of vertices and edges (x 1 , E 1 , x 2 , E 2 , . . . , x k , E k ), with all edges and vertices distinct, such that x i , x i+1 ∈ E i for all 1 ≤ i < k, with the special case of E k containing both x k and x 1 . A hypergraph is balanced if every odd cycle uses an edge containing three vertices in the cycle (see Figure 5). As for graphs, the existence of disjoint mis in hypergraphs is not verified in general. B. Acharya conjectured that every balanced hypergraph admits disjoint mis. B. Jose and Z. Tuza disproved that conjecture [START_REF] Bibin | Hypergraph domination and strong independence[END_REF], and suggested that the question might still be relevant for k-connected balanced hypergraphs [START_REF] Bibin | Hypergraph domination and strong independence[END_REF]Problem 4], that is hypergraphs where there are no subset of k -1 vertices whose removal would disconnect the hypergraph. Here, we answer to this question by the negative proposing for any k an infinite family of balanced k-connected hypergraphs H k,t that does not have disjoint mis.

The hypergraphs H k,t = (V, E) are constructed as follows: let

V = {u i,j | 1 ≤ i ≤ k + 1, 1 ≤ j ≤ t} and E = {X i | 1 ≤ i ≤ k + 1, i = 2} ∪ {Y j , Z j | 1 ≤ j ≤
t} where X i = {u i,1 , . . . , u i,t }, Y j = {u 2,j , . . . , u k+1,j } and Z j = {u 1,j , u 2,j }. Figure 6 illustrates that construction.

Proposition 11 Let k ∈ N, t ≥ 2k + 1. the hypergraph H k,t is balanced, kconnected, and does not have disjoint mis.

Proof: We first prove that H k,t is balanced. Note that any vertex of H k,t belongs to precisely two edges, one in We now prove that the hypergraph H k,t is k-connected, by giving at least k vertex disjoint walks between any two non-adjacent vertices u i,j and u i ,j (1 ≤ i ≤ i ≤ k, 1 ≤ j, j ≤ t). Note first that if i = i , then both are equal to 2 (or the corresponding vertices would be adjacent). For each case, we propose the (at least k) following vertex disjoint walks :

E 1 = {Y j | 1 ≤ j ≤ t} ∪ {X 1 } and one in E 2 = {Z j | 1 ≤ j ≤ t} ∪ {X i | 3 ≤ i ≤ k + 1}.
X 1 X 3 X 4 X 5 X k+1 Y 1 Y 2 Y 3
• if i = i = 2, {(u 2,j , u α,j , u α,j , u 2,j ), 1 ≤ α ≤ k + 1, α = 2}

• if i = 2 and i = 2, {(u i,j , u i,β , u i ,β , u i ,j ), 1 ≤ β ≤ t}

• if i = 2, i ≥ 3, {(u i,j , u α,j , u α,j , [u 2,j , ]u i ,j ), 1 ≤ α ≤ k + 1, α = 2}

• if i = 1, i = 2, {(u 1,j , u 1,α , u 2,α , u α,α , u α,j , u 2,j ) | 1 ≤ α ≤ k + 1, α = 2}

Let us finally show that H k,t has no disjoint mis. We denote by U 2 the set of vertices {u 2,j | 1 ≤ j ≤ t}. Let S be a maximal independent set of H k,t . Note first that S \ U 2 contains at most k vertices, or 2 would be in a same edge X i . For any j such that {u i,j | 1 ≤ i ≤ k + 1} ∩ (S \ U 2 ) = ∅, by maximality of S, u 2,j ∈ S. Therefore, S contains at least t -k vertices in U 2 . Now, since t ≥ 2k + 1, any two maximal independent sets intersect in some vertex u 2,j , and H k,t contains no disjoint M IS.

2

Figure 1 :

 1 Figure 1: Red and blue vertices form disjoint mis of G 1 . The graph G 2 does not admit disjoint mis.

Figure 2 :

 2 Figure 2: Coronas of odd cycle have no disjoint mis.

Conjecture 1 (

 1 C. Berge and C. Payan, 1970's) Every non empty regular simple graph contains disjoint mis.

Figure 3 :

 3 Figure 3: Situation in the proof of Lemma 8

Figure 4 :Proposition 10

 410 Figure 4: The graph G 5 which has neither disjoint mis nor disjoint maximal cliques.

6 Figure 5 :

 65 Figure 5: A balanced odd cycle in a hypergraph

Figure 6 :

 6 Figure 6: The hypergraph H k,t

  Therefore, any cycle C of H k,t uses alternating edges among the two sets E 1 and E 2 . Thus |C ∩ E 1 | is equal to |C ∩ E 2 | and C is an even cycle.