Note
Click here to download the full example code
Height source¶
This case is showing the ability of freshkiss3d to apply source terms on water height to simulate rain or infiltration effects.
import os, sys
import argparse
import numpy as np
import random as rd
import matplotlib.pylab as plt
import freshkiss3d as fk
import freshkiss3d.extra.plots as fk_plt
os.system('rm -r outputs')
parser = argparse.ArgumentParser()
parser.add_argument('--nographics', action='store_true')
args = parser.parse_args()
#sphinx_gallery_thumbnail_number = 3
Time loop:¶
FINAL_TIME = 10.
simutime = fk.SimuTime(final_time=FINAL_TIME, time_iteration_max=20000, second_order=True)
create_figure_scheduler = fk.schedules(times=[0., 4., 8.])
Mesh:¶
dir_path = os.path.abspath(os.path.dirname(sys.argv[0]))
triangular_mesh = fk.TriangularMesh.from_msh_file(dir_path + '/inputs/square2.mesh')
if not args.nographics:
fk_plt.plot_mesh(triangular_mesh)
Layers:¶
NL = 2
layer = fk.Layer(NL, triangular_mesh, topography=0.)
Primitives:¶
Initial height is set to \(10 cm\).
primitives = fk.Primitive(triangular_mesh, layer, height=0.1)
Height source term:¶
ExternalEffects objects are called in problem.solve() at each time step. In
this example external effects are HeightSource.
There can be as much external effects as needed as long as they are stored in a
dictionary type.
Let us consider two cases:
In the first case we define continuous height source over time which means that the source given is either the rate of infiltration (if negative) or the rain intensity (if positive) both in m/s.
In the second case we define sparse height source over time to simulate water droplets falling on free surface at specified times. Water height source unit is set to m.
Note
Note that in both cases (continuous or sparse over time), the source can be either
homogeneous or heterogeneous in space i.e. user can set either a constant set of sources
with the sources parameter or a function (or set of function) of x, y with
the sources_funct parameter.
CASE = 2
Case 1: In this case rain intensity is set to \(2 cm/s\) and infiltration rate to \(1 cm/s\). With an initial height of \(10 cm\) and a final time of \(10 s\), height should double by the end of the simulation. In this case the source is continuous in time and homogeneous in space.
if CASE == 1:
rain = fk.HeightSource(
times=[1., FINAL_TIME],
sources=[0.02, 0.02],
source_unit='m/s')
infiltration = fk.HeightSource(
times=[1., FINAL_TIME],
sources=[-0.01, -0.01],
source_unit='m/s')
external_effects = {"rain": rain, "infiltration": infiltration}
Case 2: water droplets are defined by a gaussian function with random position and random time of fall. In this case, the source is sparse in time and heterogeneous in space.
if CASE == 2:
TG = triangular_mesh.triangulation
x_min = min(TG.x); x_max = max(TG.x)
y_min = min(TG.y); y_max = max(TG.y)
def function_generator():
def water_drop_source(x, y, x_0=rd.uniform(x_min, x_max),
y_0=rd.uniform(y_min, y_max)):
return .003*np.exp(-10.*((x-x_0)**2+(y-y_0)**2))
return water_drop_source
#randomized water drops:
drop_number = 200
water_drop_times = np.sort(np.array([rd.uniform(0.,FINAL_TIME) for t in range(drop_number)]))
water_drop_functions = [function_generator() for t in range(drop_number)]
rain_sparse = fk.HeightSource(
times=water_drop_times,
sources_funct=water_drop_functions,
source_unit='m')
external_effects = {"sparse rain": rain_sparse}
Boundary conditions:¶
slides = [fk.Slide(ref=r) for r in [1,2,3,4]]
Writter:¶
vtk_writer = fk.VTKWriter(triangular_mesh, scheduler=fk.schedules(count=10),
scale_h=20.)
Problem definition:¶
problem = fk.Problem(simutime, triangular_mesh,
layer, primitives,
slides=slides,
external_effects=external_effects,
vtk_writer=vtk_writer,
custom_funct={'plot':fk_plt.plot_freesurface_3d},
custom_funct_scheduler=create_figure_scheduler)
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=3303, Layers=2, Triangles=6408,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
Problem solving:¶
problem.solve()
if not args.nographics:
plt.show()
===================================================================
| TIME LOOP |
===================================================================
Iter = 6 ; Dt = 0.0346s ; Time = 0.21s ; ETA = 2.91s
Iter = 12 ; Dt = 0.0346s ; Time = 0.42s ; ETA = 2.83s
Iter = 18 ; Dt = 0.0345s ; Time = 0.62s ; ETA = 4.28s
Iter = 24 ; Dt = 0.0345s ; Time = 0.83s ; ETA = 4.13s
Iter = 30 ; Dt = 0.0345s ; Time = 1.04s ; ETA = 4.02s
Iter = 36 ; Dt = 0.0345s ; Time = 1.24s ; ETA = 2.40s
Iter = 42 ; Dt = 0.0345s ; Time = 1.45s ; ETA = 3.55s
Iter = 48 ; Dt = 0.0343s ; Time = 1.66s ; ETA = 2.32s
Iter = 54 ; Dt = 0.0344s ; Time = 1.86s ; ETA = 2.47s
Iter = 60 ; Dt = 0.0345s ; Time = 2.07s ; ETA = 3.55s
Iter = 66 ; Dt = 0.0345s ; Time = 2.28s ; ETA = 61.65s
Iter = 72 ; Dt = 0.0346s ; Time = 2.48s ; ETA = 2.27s
Iter = 77 ; Dt = 0.0346s ; Time = 2.66s ; ETA = 3.32s
Iter = 83 ; Dt = 0.0346s ; Time = 2.86s ; ETA = 3.20s
Iter = 89 ; Dt = 0.0346s ; Time = 3.07s ; ETA = 2.84s
Iter = 95 ; Dt = 0.0346s ; Time = 3.28s ; ETA = 3.03s
Iter = 101 ; Dt = 0.0346s ; Time = 3.49s ; ETA = 2.72s
Iter = 107 ; Dt = 0.0346s ; Time = 3.69s ; ETA = 1.74s
Iter = 113 ; Dt = 0.0344s ; Time = 3.90s ; ETA = 2.55s
Iter = 119 ; Dt = 0.0343s ; Time = 4.11s ; ETA = 1.60s
Iter = 125 ; Dt = 0.0344s ; Time = 4.31s ; ETA = 1.70s
Iter = 131 ; Dt = 0.0344s ; Time = 4.52s ; ETA = 2.53s
Iter = 137 ; Dt = 0.0345s ; Time = 4.73s ; ETA = 2.17s
Iter = 142 ; Dt = 0.0345s ; Time = 4.90s ; ETA = 2.06s
Iter = 148 ; Dt = 0.0345s ; Time = 5.11s ; ETA = 2.03s
Iter = 154 ; Dt = 0.0345s ; Time = 5.31s ; ETA = 1.93s
Iter = 160 ; Dt = 0.0345s ; Time = 5.52s ; ETA = 1.23s
Iter = 166 ; Dt = 0.0345s ; Time = 5.73s ; ETA = 1.28s
Iter = 172 ; Dt = 0.0345s ; Time = 5.93s ; ETA = 1.11s
Iter = 178 ; Dt = 0.0345s ; Time = 6.14s ; ETA = 1.57s
Iter = 184 ; Dt = 0.0345s ; Time = 6.35s ; ETA = 1.47s
Iter = 190 ; Dt = 0.0345s ; Time = 6.56s ; ETA = 1.44s
Iter = 196 ; Dt = 0.0342s ; Time = 6.76s ; ETA = 0.92s
Iter = 202 ; Dt = 0.0343s ; Time = 6.97s ; ETA = 0.86s
Iter = 208 ; Dt = 0.0344s ; Time = 7.17s ; ETA = 1.19s
Iter = 214 ; Dt = 0.0345s ; Time = 7.38s ; ETA = 1.10s
Iter = 219 ; Dt = 0.0345s ; Time = 7.55s ; ETA = 1.04s
Iter = 225 ; Dt = 0.0345s ; Time = 7.76s ; ETA = 0.94s
Iter = 231 ; Dt = 0.0345s ; Time = 7.97s ; ETA = 0.57s
Iter = 237 ; Dt = 0.0345s ; Time = 8.17s ; ETA = 0.77s
Iter = 243 ; Dt = 0.0344s ; Time = 8.38s ; ETA = 0.44s
Iter = 249 ; Dt = 0.0344s ; Time = 8.59s ; ETA = 0.58s
Iter = 255 ; Dt = 0.0345s ; Time = 8.79s ; ETA = 0.49s
Iter = 261 ; Dt = 0.0345s ; Time = 9.00s ; ETA = 0.45s
Iter = 267 ; Dt = 0.0345s ; Time = 9.21s ; ETA = 0.35s
Iter = 273 ; Dt = 0.0345s ; Time = 9.41s ; ETA = 0.24s
Iter = 279 ; Dt = 0.0345s ; Time = 9.62s ; ETA = 0.15s
Iter = 285 ; Dt = 0.0345s ; Time = 9.83s ; ETA = 0.07s
Iter = 290 ; Dt = 0.0331s ; Time = 10.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 6.636334180831909s (wall time)
Total running time of the script: ( 0 minutes 8.768 seconds)