# Cremonesi et al benchmark

### From Gerris

Revision as of 01:16, 29 March 2011Delauxs (Talk | contribs) ← Previous diff |
Revision as of 01:28, 29 March 2011Delauxs (Talk | contribs) Next diff → |
||

Line 25: |
Line 25: | ||

Here the slide is modeled as a [http://en.wikipedia.org/wiki/Bingham_plastic Bingham fluid]. The setup is well described in (Cremonesi et al., 2011) except for the initial position of the slide that I could not find in the papers. Fortunately, I managed to digitize the initial position of the slide from Figure 10-1 and 11-1 of (Cremonesi et al., 2011). | Here the slide is modeled as a [http://en.wikipedia.org/wiki/Bingham_plastic Bingham fluid]. The setup is well described in (Cremonesi et al., 2011) except for the initial position of the slide that I could not find in the papers. Fortunately, I managed to digitize the initial position of the slide from Figure 10-1 and 11-1 of (Cremonesi et al., 2011). | ||

- | Two configurations are tested here. First the case of a slide with no initial velocity (initial x coordinate of the leading edge of the slide = -0.5 m with respect to the shoreline), and then that of a slide with an initial velocity of 3.153 m/s (initial x coordinate of the leading edge of the slide = 0.84 m with respect to the shoreline). | + | Two configurations are tested here. First the case of a slide with no initial velocity (initial x coordinate of the leading edge of the slide = -0.5 m with respect to the shoreline), and then that of a slide with an initial velocity of 3.153 m/s (initial x coordinate of the leading edge of the slide = 0.84 m with respect to the shoreline). |

+ | |||

+ | The experiments are quasi two-dimensional and therefore as (Cremonesi et al., 2011) we start by doing the two-dimensional simulations. | ||

+ | |||

+ | == Results == | ||

+ | |||

+ | Those are preliminary results obtained using the viscosity and yield stress that (Cremonesi et al., 2011) had to use so that their Bingham fluid would collide into the water at the same velocity as that recorded in the experiments. | ||

+ | |||

+ | The results are not impressive with both amplitude and timing significantly wrong. | ||

+ | |||

+ | |||

+ | |||

+ | == Simulation file == | ||

+ | |||

+ | <source lang=gfs> | ||

+ | |||

+ | # Maximum resolution | ||

+ | Define LEVEL 6 | ||

+ | |||

+ | # Minimum background resolution | ||

+ | Define LMIN 2 | ||

+ | |||

+ | # Initial velocity of the slide | ||

+ | Define U0 3.153 | ||

+ | |||

+ | # Distance between the leading edge of the slide and the shoreline | ||

+ | Define X0 0.84 | ||

+ | |||

+ | # Bingham fluid | ||

+ | Define MODEL 3 | ||

+ | |||

+ | 13 17 GfsSimulationMoving GfsBox GfsGEdge { x = -0.2 y = 0.2 } { | ||

+ | |||

+ | # Sets the duration of the simulation to 3 time units | ||

+ | Time { end = 4 dtmax = 0.01} | ||

+ | |||

+ | # We force the time-step to be as low as 0.005 for the first 5 time-steps | ||

+ | Event { step = 0.00001 iend = 5 } | ||

+ | |||

+ | # The CFL tolerance is set to 0.4 | ||

+ | AdvectionParams { cfl = 0.4 } | ||

+ | |||

+ | Global { | ||

+ | #define r_w 1000.0 /* Density of water */ | ||

+ | #define r_a 10.0 /* Density of air */ | ||

+ | #define r_s 1500.0 /* Density of solid */ | ||

+ | |||

+ | #define mu_w 0.0001 /* viscosity of water */ | ||

+ | #define mu_a 0.0000178 /* viscosity of air */ | ||

+ | #define mu_l 75. /* viscosity of liquid slide */ | ||

+ | #define mu_s 1000. /* viscosity of solid slide */ | ||

+ | |||

+ | double slide (double x, double y, double z) { | ||

+ | |||

+ | return MIN(MIN(-(y - (x + 2.*X0+sqrt(2.)*0.6 )), (y - (x + 2.*X0))),-(y + (x - sqrt(2)*0.236))); | ||

+ | } | ||

+ | |||

+ | double RHO (double tw, double ts) { | ||

+ | return 1./(ts*1./r_s+(1.-ts)*((1-tw)*1./r_w + tw*1./r_a)); | ||

+ | } | ||

+ | |||

+ | double MU (double tw, double ts, double d2) { | ||

+ | double ty = 35.; | ||

+ | double m, mu, N; | ||

+ | |||

+ | switch (MODEL) { | ||

+ | case 0: /* Newtonian */ | ||

+ | ty = 0.; N = 1.; break; | ||

+ | case 1: /* Power-law (shear-thinning) */ | ||

+ | ty = 0.; N = 0.5; break; | ||

+ | case 2: /* Herschel-Bulkley */ | ||

+ | N = 0.5; break; | ||

+ | case 3: /* Bingham */ | ||

+ | N = 1.; break; | ||

+ | } | ||

+ | if (d2 > 0.) | ||

+ | m = ty/(2.*d2) + mu_l*exp ((N - 1.)*log (d2)); | ||

+ | else { | ||

+ | if (ty > 0. || N < 1.) m = mu_s; | ||

+ | else m = N == 1. ? mu_l : 0.; | ||

+ | } | ||

+ | mu = MIN (m, mu_s); | ||

+ | |||

+ | return 1./(ts*1./mu+(1.-ts)*((1-tw)*1./mu_w + tw*1./mu_a)); | ||

+ | } | ||

+ | } | ||

+ | |||

+ | |||

+ | |||

+ | # Physical length of a box | ||

+ | PhysicalParams { L = 1. } | ||

+ | |||

+ | # Slope | ||

+ | Solid (y + x) | ||

+ | |||

+ | # Air/Water VOF tracer | ||

+ | VariableTracerVOF T | ||

+ | |||

+ | # Slide VOF tracer | ||

+ | VariableTracerVOF TS | ||

+ | |||

+ | # Initial air/water interface | ||

+ | InitFraction T (y) | ||

+ | |||

+ | # Initial slide position | ||

+ | InitFraction TS (slide(x,y,z)) | ||

+ | |||

+ | # Sets initial slope, air/water interface and slope refinement | ||

+ | RefineSolid LEVEL | ||

+ | RefineSurface { return x > -1. ? LEVEL : LMIN;} (y) | ||

+ | RefineSurface LEVEL (slide(x,y,z)) | ||

+ | |||

+ | # We lower the tolerance on the Poisson solver | ||

+ | ApproxProjectionParams { tolerance = 1e-5 } | ||

+ | ProjectionParams { tolerance = 1e-5 } | ||

+ | |||

+ | # Sets the initial velocity of the slide | ||

+ | Init { } { U = (TS*U0/sqrt(2.)) V = -(TS*U0/sqrt(2.)) } | ||

+ | |||

+ | # Density | ||

+ | PhysicalParams { alpha = 1./RHO(T,TS) } | ||

+ | |||

+ | # The osition of the VOF interface is stored in X and Y | ||

+ | VariablePosition Y T y | ||

+ | VariablePosition X T x | ||

+ | |||

+ | # Gravity | ||

+ | Source {} V -9.81 | ||

+ | |||

+ | Init { istep = 1 } { SEB = MU(T, TS, D2) } | ||

+ | |||

+ | # Refinement if the air water interface | ||

+ | AdaptFunction { istart = 1 istep = 1 } { | ||

+ | cmax = 0 | ||

+ | maxlevel = (x < 4 ? LEVEL : LEVEL - 2) | ||

+ | minlevel = (T < 1. ? LMIN:0) | ||

+ | } (T > 0 && T < 1) | ||

+ | |||

+ | # High refinement inside the slide | ||

+ | AdaptFunction { istart = 1 istep = 1 } { | ||

+ | cmax = 0 | ||

+ | maxlevel = (LEVEL - 1) | ||

+ | } ( TS ) | ||

+ | |||

+ | # Refinement of the slide's interface | ||

+ | AdaptFunction { istart = 1 istep = 1 } { | ||

+ | cmax = 0 | ||

+ | maxlevel = (x < 4 ? LEVEL : LEVEL - 2) | ||

+ | minlevel = (T < 1. ? LMIN:0) | ||

+ | } (TS > 0 && TS < 1) | ||

+ | |||

+ | # Refinement using the vorticity criterion | ||

+ | AdaptVorticity { istep = 1 } { | ||

+ | cmax = 1.0 | ||

+ | maxlevel = LEVEL | ||

+ | minlevel = LMIN | ||

+ | } | ||

+ | |||

+ | # Wave gauges | ||

+ | OutputLocation { istep = 1 } { awk -f distance.awk > probe1 } probe1x { interpolate = 0 } | ||

+ | OutputLocation { istep = 1 } { awk -f distance.awk > probe2 } probe2x { interpolate = 0 } | ||

+ | OutputLocation { istep = 1 } { awk -f distance.awk > probe3 } probe3x { interpolate = 0 } | ||

+ | OutputLocation { istep = 1 } { awk -f distance.awk > probe4 } probe4x { interpolate = 0 } | ||

+ | |||

+ | # Output control statistics on the simulation | ||

+ | OutputTime { istep = 1 } stderr | ||

+ | OutputBalance { istep = 1 } stderr | ||

+ | OutputProjectionStats { istep = 1 } stderr | ||

+ | |||

+ | # Output simulation files every 0.1 s | ||

+ | OutputSimulation { step = 0.1 end = 3} sim-%g.gfs | ||

+ | |||

+ | OutputSimulation { step = 0.5 start = 1.5} sim-%g.gfs | ||

+ | |||

+ | # Viscosity of the fluids | ||

+ | SourceViscosity {} (MU(T,TS,D2)) { | ||

+ | beta = 1 | ||

+ | } | ||

+ | |||

+ | # Generates a movie | ||

+ | GModule gfsview | ||

+ | OutputView { step = 2e-2 } { ppm2mpeg > movie.mpg } { width = 800 height = 600 } view2D.gfv | ||

+ | |||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | bottom = Boundary { BcDirichlet U 0 } | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | GfsBox { | ||

+ | top = BoundaryOutflow | ||

+ | } | ||

+ | 1 2 right | ||

+ | 2 3 right | ||

+ | 3 4 right | ||

+ | 4 5 right | ||

+ | 5 6 right | ||

+ | 7 8 right | ||

+ | 8 9 right | ||

+ | 9 10 right | ||

+ | 10 11 right | ||

+ | 11 12 right | ||

+ | 1 7 top | ||

+ | 2 8 top | ||

+ | 3 9 top | ||

+ | 4 10 top | ||

+ | 5 11 top | ||

+ | 6 12 top | ||

+ | 7 13 left | ||

+ | |||

+ | </source> |

## Revision as of 01:28, 29 March 2011

This benchmark for simulations of landslide tsunamis generated by a deformable slide using Gerris is based on these articles:

*Zweifel, A., Zuccala, D., Gatti, D.* - **Comparison between computed and experimentally generated impulse waves**

- Journal of Hydraulic Engineering 133(2):208--216, 2007
- Bibtex

*Cremonesi, M., Frangi, A., Perego, U.* - **A Lagrangian finite element approach for the simulation of water-waves induced by landslides**

- Computers and Structures In press, 2011
- Bibtex

Here the slide is modeled as a Bingham fluid. The setup is well described in (Cremonesi et al., 2011) except for the initial position of the slide that I could not find in the papers. Fortunately, I managed to digitize the initial position of the slide from Figure 10-1 and 11-1 of (Cremonesi et al., 2011).

Two configurations are tested here. First the case of a slide with no initial velocity (initial x coordinate of the leading edge of the slide = -0.5 m with respect to the shoreline), and then that of a slide with an initial velocity of 3.153 m/s (initial x coordinate of the leading edge of the slide = 0.84 m with respect to the shoreline).

The experiments are quasi two-dimensional and therefore as (Cremonesi et al., 2011) we start by doing the two-dimensional simulations.

## Results

Those are preliminary results obtained using the viscosity and yield stress that (Cremonesi et al., 2011) had to use so that their Bingham fluid would collide into the water at the same velocity as that recorded in the experiments.

The results are not impressive with both amplitude and timing significantly wrong.

## Simulation file

# Maximum resolution

Define LEVEL 6

# Minimum background resolution

Define LMIN 2

# Initial velocity of the slide

Define U0 3.153

# Distance between the leading edge of the slide and the shoreline

Define X0 0.84

# Bingham fluid

Define MODEL 3

13 17 GfsSimulationMoving GfsBox GfsGEdge { x = -0.2 y = 0.2 } {

# Sets the duration of the simulation to 3 time units

Time { end = 4 dtmax = 0.01}

# We force the time-step to be as low as 0.005 for the first 5 time-steps

Event { step = 0.00001 iend = 5 }

# The CFL tolerance is set to 0.4

AdvectionParams { cfl = 0.4 }

Global {

#define r_w 1000.0 /* Density of water */

#define r_a 10.0 /* Density of air */

#define r_s 1500.0 /* Density of solid */

#define mu_w 0.0001 /* viscosity of water */

#define mu_a 0.0000178 /* viscosity of air */

#define mu_l 75. /* viscosity of liquid slide */

#define mu_s 1000. /* viscosity of solid slide */

double slide (double x, double y, double z) {

return MIN(MIN(-(y - (x + 2.*X0+sqrt(2.)*0.6 )), (y - (x + 2.*X0))),-(y + (x - sqrt(2)*0.236)));

}

double RHO (double tw, double ts) {

return 1./(ts*1./r_s+(1.-ts)*((1-tw)*1./r_w + tw*1./r_a));

}

double MU (double tw, double ts, double d2) {

double ty = 35.;

double m, mu, N;

switch (MODEL) {

case 0: /* Newtonian */

ty = 0.; N = 1.; break;

case 1: /* Power-law (shear-thinning) */

ty = 0.; N = 0.5; break;

case 2: /* Herschel-Bulkley */

N = 0.5; break;

case 3: /* Bingham */

N = 1.; break;

}

if (d2 > 0.)

m = ty/(2.*d2) + mu_l*exp ((N - 1.)*log (d2));

else {

if (ty > 0. || N < 1.) m = mu_s;

else m = N == 1. ? mu_l : 0.;

}

mu = MIN (m, mu_s);

return 1./(ts*1./mu+(1.-ts)*((1-tw)*1./mu_w + tw*1./mu_a));

}

}

# Physical length of a box

PhysicalParams { L = 1. }

# Slope

Solid (y + x)

# Air/Water VOF tracer

VariableTracerVOF T

# Slide VOF tracer

VariableTracerVOF TS

# Initial air/water interface

InitFraction T (y)

# Initial slide position

InitFraction TS (slide(x,y,z))

# Sets initial slope, air/water interface and slope refinement

RefineSolid LEVEL

RefineSurface { return x > -1. ? LEVEL : LMIN;} (y)

RefineSurface LEVEL (slide(x,y,z))

# We lower the tolerance on the Poisson solver

ApproxProjectionParams { tolerance = 1e-5 }

ProjectionParams { tolerance = 1e-5 }

# Sets the initial velocity of the slide

Init { } { U = (TS*U0/sqrt(2.)) V = -(TS*U0/sqrt(2.)) }

# Density

PhysicalParams { alpha = 1./RHO(T,TS) }

# The osition of the VOF interface is stored in X and Y

VariablePosition Y T y

VariablePosition X T x

# Gravity

Source {} V -9.81

Init { istep = 1 } { SEB = MU(T, TS, D2) }

# Refinement if the air water interface

AdaptFunction { istart = 1 istep = 1 } {

cmax = 0

maxlevel = (x < 4 ? LEVEL : LEVEL - 2)

minlevel = (T < 1. ? LMIN:0)

} (T > 0 && T < 1)

# High refinement inside the slide

AdaptFunction { istart = 1 istep = 1 } {

cmax = 0

maxlevel = (LEVEL - 1)

} ( TS )

# Refinement of the slide's interface

AdaptFunction { istart = 1 istep = 1 } {

cmax = 0

maxlevel = (x < 4 ? LEVEL : LEVEL - 2)

minlevel = (T < 1. ? LMIN:0)

} (TS > 0 && TS < 1)

# Refinement using the vorticity criterion

AdaptVorticity { istep = 1 } {

cmax = 1.0

maxlevel = LEVEL

minlevel = LMIN

}

# Wave gauges

OutputLocation { istep = 1 } { awk -f distance.awk > probe1 } probe1x { interpolate = 0 }

OutputLocation { istep = 1 } { awk -f distance.awk > probe2 } probe2x { interpolate = 0 }

OutputLocation { istep = 1 } { awk -f distance.awk > probe3 } probe3x { interpolate = 0 }

OutputLocation { istep = 1 } { awk -f distance.awk > probe4 } probe4x { interpolate = 0 }

# Output control statistics on the simulation

OutputTime { istep = 1 } stderr

OutputBalance { istep = 1 } stderr

OutputProjectionStats { istep = 1 } stderr

# Output simulation files every 0.1 s

OutputSimulation { step = 0.1 end = 3} sim-%g.gfs

OutputSimulation { step = 0.5 start = 1.5} sim-%g.gfs

# Viscosity of the fluids

SourceViscosity {} (MU(T,TS,D2)) {

beta = 1

}

# Generates a movie

GModule gfsview

OutputView { step = 2e-2 } { ppm2mpeg > movie.mpg } { width = 800 height = 600 } view2D.gfv

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

bottom = Boundary { BcDirichlet U 0 }

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

GfsBox {

top = BoundaryOutflow

}

1 2 right

2 3 right

3 4 right

4 5 right

5 6 right

7 8 right

8 9 right

9 10 right

10 11 right

11 12 right

1 7 top

2 8 top

3 9 top

4 10 top

5 11 top

6 12 top

7 13 left