Note
Click here to download the full example code
Wind effect¶
In this example WindEffect
is illustrated. At the top of a cavity we set up
wind velocity and orientation. A circular flow is generated in the cavity.
import os, sys
import matplotlib.pylab as plt
import numpy as np
import freshkiss3d as fk
import freshkiss3d.extra.plots as fk_plt
#sphinx_gallery_thumbnail_number = 2
Wind effect:¶
The wind friction shear stress on top layer can be expressed as:
In freshkiss3d the friction coefficient Cd (empirical parameter) can be user-defined or chosen amongst a set of existing law (by default: Imberger & Patterson)
Imberger & Patterson (1984)
\[\begin{split}\begin{cases} C_d = 1,124.10^{-3} \quad & \text{if} \quad U_{wind} \leq 4 m/s \\ C_d = (0,94 + 0,041 U_{wind}) 10^{-3} \quad & \text{else} \end{cases}\end{split}\]Stefan & Ford (1975)
\[\begin{split}\begin{cases} C_d = 0,5.10^{-3} \sqrt{U_{wind}} \quad & \text{if} \quad U_{wind} \leq 15 m/s \\ C_d = 2,6.10^{-3} \quad & \text{else} \end{cases}\end{split}\]Coantic (1978)
\[C_d = (1 + 0.03 U_{wind}).10^{-3}\]
cases = ['Imberger_Patterson', 'Stefan_Ford', 'Coantic', 'Custom']
def custom_friction_law(v):
return (0.75 + 0.067*v)*1.e-3
In this example we compare all laws for a wind intensity of 20 m/s. Fluid dynamic viscosity is set to 1Pa.s and free surface to 2m.
Mesh:¶
dir_path = os.path.abspath(os.path.dirname(sys.argv[0]))
TG, vertex_labels, boundary_labels = fk.read_msh(dir_path + '/inputs/simple_canal_2.msh')
x = np.asarray(TG.x)
y = np.asarray(TG.y)
trivtx = np.asarray(TG.trivtx)
x *= 10
y *= 10
triangular_mesh = fk.TriangularMesh(TG, vertex_labels, boundary_labels)
fk_plt.plot_mesh(triangular_mesh)
Cases set-up:¶
WindEffect
class need a velocity (in m/s) and orientation (in Azimuth degree).
Both can be set to change over time. friction_law
and friction_law_funct
are optional arguments.
NL = 20
PHY_PARAMS={'friction_coeff':1.,
'horizontal_viscosity':0.01,
'vertical_viscosity':0.01}
simutime_l = []
layer_l = []
primitives_l = []
wind_effect_l = []
external_effects_l = []
problem_l = []
for C, law in enumerate(cases):
simutime_l.append(fk.SimuTime(final_time=5., time_iteration_max=10000,
second_order=False))
layer_l.append(fk.Layer(NL, triangular_mesh, topography=0.))
primitives_l.append(fk.Primitive(triangular_mesh, layer_l[C], free_surface=2.0))
if law == 'Custom':
wind_effect_l.append(fk.WindEffect(velocity=[20, 20, 20, 20],
orientation=[270., 270., 270., 270.],
times=[0., 3., 4., 5.],
friction_law_funct=custom_friction_law))
else:
wind_effect_l.append(fk.WindEffect(velocity=[20, 20, 20, 20],
orientation=[270., 270., 270., 270.],
times=[0., 3., 4., 5.],
friction_law=law))
external_effects_l.append({"wind": wind_effect_l[C]})
slides = [fk.Slide(ref=r) for r in [1, 2, 3, 4]]
Problem:¶
for C, law in enumerate(cases):
print(" ")
print(" ")
print(" SOLVING CASE: ")
print(" ")
print(" ")
problem_l.append(fk.Problem(simutime_l[C], triangular_mesh,
layer_l[C], primitives_l[C],
slides=slides,
physical_parameters=PHY_PARAMS,
external_effects=external_effects_l[C]))
problem_l[C].solve()
SOLVING CASE:
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=144, Layers=20, Triangles=240,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
===================================================================
| TIME LOOP |
===================================================================
Iter = 7 ; Dt = 0.0164s ; Time = 0.11s ; ETA = 0.59s
Iter = 13 ; Dt = 0.0164s ; Time = 0.21s ; ETA = 0.58s
Iter = 19 ; Dt = 0.0164s ; Time = 0.31s ; ETA = 0.57s
Iter = 25 ; Dt = 0.0164s ; Time = 0.41s ; ETA = 0.56s
Iter = 32 ; Dt = 0.0164s ; Time = 0.52s ; ETA = 0.54s
Iter = 38 ; Dt = 0.0164s ; Time = 0.62s ; ETA = 0.54s
Iter = 44 ; Dt = 0.0164s ; Time = 0.72s ; ETA = 0.53s
Iter = 50 ; Dt = 0.0164s ; Time = 0.82s ; ETA = 0.51s
Iter = 57 ; Dt = 0.0164s ; Time = 0.93s ; ETA = 0.50s
Iter = 63 ; Dt = 0.0164s ; Time = 1.03s ; ETA = 0.49s
Iter = 69 ; Dt = 0.0164s ; Time = 1.13s ; ETA = 0.47s
Iter = 75 ; Dt = 0.0164s ; Time = 1.23s ; ETA = 0.46s
Iter = 81 ; Dt = 0.0164s ; Time = 1.33s ; ETA = 0.46s
Iter = 88 ; Dt = 0.0164s ; Time = 1.44s ; ETA = 0.43s
Iter = 94 ; Dt = 0.0164s ; Time = 1.54s ; ETA = 0.43s
Iter = 100 ; Dt = 0.0164s ; Time = 1.64s ; ETA = 0.42s
Iter = 106 ; Dt = 0.0164s ; Time = 1.74s ; ETA = 0.40s
Iter = 113 ; Dt = 0.0164s ; Time = 1.85s ; ETA = 0.40s
Iter = 119 ; Dt = 0.0164s ; Time = 1.95s ; ETA = 0.38s
Iter = 125 ; Dt = 0.0164s ; Time = 2.05s ; ETA = 0.37s
Iter = 131 ; Dt = 0.0164s ; Time = 2.15s ; ETA = 0.35s
Iter = 137 ; Dt = 0.0164s ; Time = 2.25s ; ETA = 0.34s
Iter = 144 ; Dt = 0.0164s ; Time = 2.36s ; ETA = 0.32s
Iter = 150 ; Dt = 0.0164s ; Time = 2.46s ; ETA = 0.31s
Iter = 156 ; Dt = 0.0164s ; Time = 2.56s ; ETA = 0.29s
Iter = 162 ; Dt = 0.0164s ; Time = 2.65s ; ETA = 0.28s
Iter = 169 ; Dt = 0.0164s ; Time = 2.77s ; ETA = 0.27s
Iter = 175 ; Dt = 0.0164s ; Time = 2.87s ; ETA = 0.26s
Iter = 181 ; Dt = 0.0164s ; Time = 2.97s ; ETA = 0.25s
Iter = 187 ; Dt = 0.0164s ; Time = 3.06s ; ETA = 0.24s
Iter = 194 ; Dt = 0.0164s ; Time = 3.18s ; ETA = 0.22s
Iter = 200 ; Dt = 0.0164s ; Time = 3.28s ; ETA = 0.21s
Iter = 206 ; Dt = 0.0164s ; Time = 3.38s ; ETA = 0.20s
Iter = 212 ; Dt = 0.0164s ; Time = 3.47s ; ETA = 0.19s
Iter = 218 ; Dt = 0.0164s ; Time = 3.57s ; ETA = 0.18s
Iter = 225 ; Dt = 0.0164s ; Time = 3.69s ; ETA = 0.16s
Iter = 231 ; Dt = 0.0164s ; Time = 3.78s ; ETA = 0.15s
Iter = 237 ; Dt = 0.0164s ; Time = 3.88s ; ETA = 0.14s
Iter = 243 ; Dt = 0.0164s ; Time = 3.98s ; ETA = 0.13s
Iter = 250 ; Dt = 0.0164s ; Time = 4.10s ; ETA = 0.11s
Iter = 256 ; Dt = 0.0164s ; Time = 4.19s ; ETA = 0.11s
Iter = 262 ; Dt = 0.0164s ; Time = 4.29s ; ETA = 0.09s
Iter = 268 ; Dt = 0.0164s ; Time = 4.39s ; ETA = 0.07s
Iter = 275 ; Dt = 0.0164s ; Time = 4.50s ; ETA = 0.06s
Iter = 281 ; Dt = 0.0164s ; Time = 4.60s ; ETA = 0.05s
Iter = 287 ; Dt = 0.0164s ; Time = 4.70s ; ETA = 0.04s
Iter = 293 ; Dt = 0.0164s ; Time = 4.80s ; ETA = 0.02s
Iter = 300 ; Dt = 0.0164s ; Time = 4.91s ; ETA = 0.01s
Iter = 306 ; Dt = 0.0040s ; Time = 5.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 0.6217336654663086s (wall time)
SOLVING CASE:
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=144, Layers=20, Triangles=240,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
===================================================================
| TIME LOOP |
===================================================================
Iter = 7 ; Dt = 0.0164s ; Time = 0.11s ; ETA = 0.58s
Iter = 13 ; Dt = 0.0164s ; Time = 0.21s ; ETA = 0.57s
Iter = 19 ; Dt = 0.0164s ; Time = 0.31s ; ETA = 0.57s
Iter = 25 ; Dt = 0.0164s ; Time = 0.41s ; ETA = 0.55s
Iter = 32 ; Dt = 0.0164s ; Time = 0.52s ; ETA = 0.54s
Iter = 38 ; Dt = 0.0164s ; Time = 0.62s ; ETA = 0.53s
Iter = 44 ; Dt = 0.0164s ; Time = 0.72s ; ETA = 0.52s
Iter = 50 ; Dt = 0.0164s ; Time = 0.82s ; ETA = 0.51s
Iter = 57 ; Dt = 0.0164s ; Time = 0.93s ; ETA = 0.49s
Iter = 63 ; Dt = 0.0164s ; Time = 1.03s ; ETA = 0.48s
Iter = 69 ; Dt = 0.0164s ; Time = 1.13s ; ETA = 0.47s
Iter = 75 ; Dt = 0.0164s ; Time = 1.23s ; ETA = 0.46s
Iter = 81 ; Dt = 0.0164s ; Time = 1.33s ; ETA = 0.44s
Iter = 88 ; Dt = 0.0164s ; Time = 1.44s ; ETA = 0.44s
Iter = 94 ; Dt = 0.0164s ; Time = 1.54s ; ETA = 0.43s
Iter = 100 ; Dt = 0.0164s ; Time = 1.64s ; ETA = 0.42s
Iter = 106 ; Dt = 0.0164s ; Time = 1.74s ; ETA = 0.40s
Iter = 113 ; Dt = 0.0164s ; Time = 1.85s ; ETA = 0.39s
Iter = 119 ; Dt = 0.0164s ; Time = 1.95s ; ETA = 0.37s
Iter = 125 ; Dt = 0.0164s ; Time = 2.05s ; ETA = 0.36s
Iter = 131 ; Dt = 0.0164s ; Time = 2.15s ; ETA = 0.35s
Iter = 138 ; Dt = 0.0164s ; Time = 2.26s ; ETA = 0.34s
Iter = 144 ; Dt = 0.0164s ; Time = 2.36s ; ETA = 0.32s
Iter = 150 ; Dt = 0.0164s ; Time = 2.46s ; ETA = 0.31s
Iter = 156 ; Dt = 0.0164s ; Time = 2.56s ; ETA = 0.29s
Iter = 162 ; Dt = 0.0164s ; Time = 2.65s ; ETA = 0.28s
Iter = 169 ; Dt = 0.0164s ; Time = 2.77s ; ETA = 0.27s
Iter = 175 ; Dt = 0.0164s ; Time = 2.87s ; ETA = 0.26s
Iter = 181 ; Dt = 0.0164s ; Time = 2.96s ; ETA = 0.26s
Iter = 187 ; Dt = 0.0164s ; Time = 3.06s ; ETA = 0.24s
Iter = 194 ; Dt = 0.0164s ; Time = 3.18s ; ETA = 0.22s
Iter = 200 ; Dt = 0.0164s ; Time = 3.27s ; ETA = 0.21s
Iter = 206 ; Dt = 0.0164s ; Time = 3.37s ; ETA = 0.20s
Iter = 212 ; Dt = 0.0164s ; Time = 3.47s ; ETA = 0.18s
Iter = 219 ; Dt = 0.0164s ; Time = 3.59s ; ETA = 0.17s
Iter = 225 ; Dt = 0.0164s ; Time = 3.68s ; ETA = 0.16s
Iter = 231 ; Dt = 0.0164s ; Time = 3.78s ; ETA = 0.15s
Iter = 237 ; Dt = 0.0164s ; Time = 3.88s ; ETA = 0.14s
Iter = 244 ; Dt = 0.0164s ; Time = 3.99s ; ETA = 0.12s
Iter = 250 ; Dt = 0.0164s ; Time = 4.09s ; ETA = 0.11s
Iter = 256 ; Dt = 0.0164s ; Time = 4.19s ; ETA = 0.10s
Iter = 262 ; Dt = 0.0164s ; Time = 4.29s ; ETA = 0.09s
Iter = 269 ; Dt = 0.0164s ; Time = 4.40s ; ETA = 0.07s
Iter = 275 ; Dt = 0.0164s ; Time = 4.50s ; ETA = 0.06s
Iter = 281 ; Dt = 0.0164s ; Time = 4.60s ; ETA = 0.05s
Iter = 287 ; Dt = 0.0164s ; Time = 4.70s ; ETA = 0.04s
Iter = 293 ; Dt = 0.0164s ; Time = 4.80s ; ETA = 0.02s
Iter = 300 ; Dt = 0.0164s ; Time = 4.91s ; ETA = 0.01s
Iter = 306 ; Dt = 0.0076s ; Time = 5.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 0.6171374320983887s (wall time)
SOLVING CASE:
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=144, Layers=20, Triangles=240,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
===================================================================
| TIME LOOP |
===================================================================
Iter = 7 ; Dt = 0.0164s ; Time = 0.11s ; ETA = 0.60s
Iter = 13 ; Dt = 0.0164s ; Time = 0.21s ; ETA = 0.58s
Iter = 19 ; Dt = 0.0164s ; Time = 0.31s ; ETA = 0.58s
Iter = 25 ; Dt = 0.0164s ; Time = 0.41s ; ETA = 0.56s
Iter = 32 ; Dt = 0.0164s ; Time = 0.52s ; ETA = 0.55s
Iter = 38 ; Dt = 0.0164s ; Time = 0.62s ; ETA = 0.54s
Iter = 44 ; Dt = 0.0164s ; Time = 0.72s ; ETA = 0.53s
Iter = 50 ; Dt = 0.0164s ; Time = 0.82s ; ETA = 0.52s
Iter = 57 ; Dt = 0.0164s ; Time = 0.93s ; ETA = 0.50s
Iter = 63 ; Dt = 0.0164s ; Time = 1.03s ; ETA = 0.49s
Iter = 69 ; Dt = 0.0164s ; Time = 1.13s ; ETA = 0.47s
Iter = 75 ; Dt = 0.0164s ; Time = 1.23s ; ETA = 0.46s
Iter = 81 ; Dt = 0.0164s ; Time = 1.33s ; ETA = 0.46s
Iter = 88 ; Dt = 0.0164s ; Time = 1.44s ; ETA = 0.43s
Iter = 94 ; Dt = 0.0164s ; Time = 1.54s ; ETA = 0.42s
Iter = 100 ; Dt = 0.0164s ; Time = 1.64s ; ETA = 0.41s
Iter = 106 ; Dt = 0.0164s ; Time = 1.74s ; ETA = 0.40s
Iter = 113 ; Dt = 0.0164s ; Time = 1.85s ; ETA = 0.38s
Iter = 119 ; Dt = 0.0164s ; Time = 1.95s ; ETA = 0.37s
Iter = 125 ; Dt = 0.0164s ; Time = 2.05s ; ETA = 0.36s
Iter = 131 ; Dt = 0.0164s ; Time = 2.15s ; ETA = 0.35s
Iter = 137 ; Dt = 0.0164s ; Time = 2.25s ; ETA = 0.33s
Iter = 144 ; Dt = 0.0164s ; Time = 2.36s ; ETA = 0.32s
Iter = 150 ; Dt = 0.0164s ; Time = 2.46s ; ETA = 0.32s
Iter = 156 ; Dt = 0.0164s ; Time = 2.56s ; ETA = 0.31s
Iter = 162 ; Dt = 0.0164s ; Time = 2.65s ; ETA = 0.29s
Iter = 169 ; Dt = 0.0164s ; Time = 2.77s ; ETA = 0.27s
Iter = 175 ; Dt = 0.0164s ; Time = 2.87s ; ETA = 0.29s
Iter = 181 ; Dt = 0.0164s ; Time = 2.97s ; ETA = 0.27s
Iter = 187 ; Dt = 0.0164s ; Time = 3.06s ; ETA = 0.26s
Iter = 194 ; Dt = 0.0164s ; Time = 3.18s ; ETA = 0.25s
Iter = 200 ; Dt = 0.0164s ; Time = 3.28s ; ETA = 0.24s
Iter = 206 ; Dt = 0.0164s ; Time = 3.38s ; ETA = 0.23s
Iter = 212 ; Dt = 0.0164s ; Time = 3.47s ; ETA = 0.20s
Iter = 218 ; Dt = 0.0164s ; Time = 3.57s ; ETA = 0.20s
Iter = 225 ; Dt = 0.0164s ; Time = 3.69s ; ETA = 0.18s
Iter = 231 ; Dt = 0.0164s ; Time = 3.78s ; ETA = 0.17s
Iter = 237 ; Dt = 0.0164s ; Time = 3.88s ; ETA = 0.14s
Iter = 243 ; Dt = 0.0164s ; Time = 3.98s ; ETA = 0.13s
Iter = 250 ; Dt = 0.0164s ; Time = 4.10s ; ETA = 0.12s
Iter = 256 ; Dt = 0.0164s ; Time = 4.19s ; ETA = 0.11s
Iter = 262 ; Dt = 0.0164s ; Time = 4.29s ; ETA = 0.10s
Iter = 268 ; Dt = 0.0164s ; Time = 4.39s ; ETA = 0.09s
Iter = 275 ; Dt = 0.0164s ; Time = 4.51s ; ETA = 0.07s
Iter = 281 ; Dt = 0.0164s ; Time = 4.60s ; ETA = 0.06s
Iter = 287 ; Dt = 0.0164s ; Time = 4.70s ; ETA = 0.04s
Iter = 293 ; Dt = 0.0164s ; Time = 4.80s ; ETA = 0.03s
Iter = 299 ; Dt = 0.0164s ; Time = 4.90s ; ETA = 0.01s
Iter = 306 ; Dt = 0.0033s ; Time = 5.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 0.6513912677764893s (wall time)
SOLVING CASE:
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=144, Layers=20, Triangles=240,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
===================================================================
| TIME LOOP |
===================================================================
Iter = 7 ; Dt = 0.0164s ; Time = 0.11s ; ETA = 0.69s
Iter = 13 ; Dt = 0.0164s ; Time = 0.21s ; ETA = 0.60s
Iter = 19 ; Dt = 0.0164s ; Time = 0.31s ; ETA = 0.57s
Iter = 25 ; Dt = 0.0164s ; Time = 0.41s ; ETA = 0.57s
Iter = 32 ; Dt = 0.0164s ; Time = 0.52s ; ETA = 0.54s
Iter = 38 ; Dt = 0.0164s ; Time = 0.62s ; ETA = 0.52s
Iter = 44 ; Dt = 0.0164s ; Time = 0.72s ; ETA = 0.52s
Iter = 50 ; Dt = 0.0164s ; Time = 0.82s ; ETA = 0.51s
Iter = 57 ; Dt = 0.0164s ; Time = 0.93s ; ETA = 0.50s
Iter = 63 ; Dt = 0.0164s ; Time = 1.03s ; ETA = 0.49s
Iter = 69 ; Dt = 0.0164s ; Time = 1.13s ; ETA = 0.47s
Iter = 75 ; Dt = 0.0164s ; Time = 1.23s ; ETA = 0.46s
Iter = 81 ; Dt = 0.0164s ; Time = 1.33s ; ETA = 0.45s
Iter = 88 ; Dt = 0.0164s ; Time = 1.44s ; ETA = 0.43s
Iter = 94 ; Dt = 0.0164s ; Time = 1.54s ; ETA = 0.42s
Iter = 100 ; Dt = 0.0164s ; Time = 1.64s ; ETA = 0.41s
Iter = 106 ; Dt = 0.0164s ; Time = 1.74s ; ETA = 0.40s
Iter = 113 ; Dt = 0.0164s ; Time = 1.85s ; ETA = 0.38s
Iter = 119 ; Dt = 0.0164s ; Time = 1.95s ; ETA = 0.37s
Iter = 125 ; Dt = 0.0164s ; Time = 2.05s ; ETA = 0.36s
Iter = 131 ; Dt = 0.0164s ; Time = 2.15s ; ETA = 0.34s
Iter = 138 ; Dt = 0.0164s ; Time = 2.26s ; ETA = 0.33s
Iter = 144 ; Dt = 0.0164s ; Time = 2.36s ; ETA = 0.33s
Iter = 150 ; Dt = 0.0164s ; Time = 2.46s ; ETA = 0.31s
Iter = 156 ; Dt = 0.0164s ; Time = 2.56s ; ETA = 0.30s
Iter = 162 ; Dt = 0.0164s ; Time = 2.65s ; ETA = 0.29s
Iter = 169 ; Dt = 0.0164s ; Time = 2.77s ; ETA = 0.27s
Iter = 175 ; Dt = 0.0164s ; Time = 2.87s ; ETA = 0.26s
Iter = 181 ; Dt = 0.0164s ; Time = 2.97s ; ETA = 0.25s
Iter = 187 ; Dt = 0.0164s ; Time = 3.06s ; ETA = 0.24s
Iter = 194 ; Dt = 0.0164s ; Time = 3.18s ; ETA = 0.22s
Iter = 200 ; Dt = 0.0164s ; Time = 3.28s ; ETA = 0.21s
Iter = 206 ; Dt = 0.0164s ; Time = 3.37s ; ETA = 0.20s
Iter = 212 ; Dt = 0.0164s ; Time = 3.47s ; ETA = 0.18s
Iter = 219 ; Dt = 0.0164s ; Time = 3.59s ; ETA = 0.19s
Iter = 225 ; Dt = 0.0164s ; Time = 3.69s ; ETA = 0.20s
Iter = 231 ; Dt = 0.0164s ; Time = 3.78s ; ETA = 0.16s
Iter = 237 ; Dt = 0.0164s ; Time = 3.88s ; ETA = 0.15s
Iter = 243 ; Dt = 0.0164s ; Time = 3.98s ; ETA = 0.14s
Iter = 250 ; Dt = 0.0164s ; Time = 4.09s ; ETA = 0.12s
Iter = 256 ; Dt = 0.0164s ; Time = 4.19s ; ETA = 0.11s
Iter = 262 ; Dt = 0.0164s ; Time = 4.29s ; ETA = 0.09s
Iter = 268 ; Dt = 0.0164s ; Time = 4.39s ; ETA = 0.09s
Iter = 275 ; Dt = 0.0164s ; Time = 4.50s ; ETA = 0.07s
Iter = 281 ; Dt = 0.0164s ; Time = 4.60s ; ETA = 0.05s
Iter = 287 ; Dt = 0.0164s ; Time = 4.70s ; ETA = 0.04s
Iter = 293 ; Dt = 0.0164s ; Time = 4.80s ; ETA = 0.03s
Iter = 300 ; Dt = 0.0164s ; Time = 4.91s ; ETA = 0.01s
Iter = 306 ; Dt = 0.0054s ; Time = 5.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 0.6409673690795898s (wall time)
Plots:¶
Vertical velocity is plotted for max layer case as well as x velocity profiles for each discretization.
plt.rcParams["figure.figsize"] = [12, 8]
plt.style.use('seaborn-white')
plt.rcParams['font.family'] = 'serif'
plt.rcParams['font.serif'] = 'Ubuntu'
plt.rcParams['font.monospace'] = 'Ubuntu Mono'
plt.rcParams['font.size'] = 10
plt.rcParams['axes.labelsize'] = 8
plt.rcParams['axes.labelweight'] = 'bold'
plt.rcParams['axes.titlesize'] = 10
plt.rcParams['xtick.labelsize'] = 8
plt.rcParams['ytick.labelsize'] = 8
plt.rcParams['legend.fontsize'] = 8
plt.rcParams['figure.titlesize'] = 12
plt.rcParams["axes.grid"] = True
fig = plt.figure()
colors = ['blue', 'purple', 'green', 'red']
markers = ['s', '*', 'o', 'v']
def plot_U_profile(cell, ax, problem_l,char):
for C in range(len(cases)):
u = np.zeros((len(cases), NL))
z = np.zeros((len(cases), NL))
for L in range(NL):
u[C, L] = problem_l[C].primitives.U[cell,L]
z[C, L] = 0.5*(2.0/NL) + L*(2.0/NL)
ax.plot(u[C, :],z[C, :], color=colors[C], marker=markers[C],
label="{}".format(cases[C]))
ax.set_title("Velocity profile on node {}".format(char), fontstyle='italic')
ax.legend(loc=4, shadow=True, title="Friction law")
ax.set_xlabel('U (m/s)')
ax.set_ylabel('z (m)')
#Velocity profile
cell = 84
ax = plt.subplot2grid((3,4), (1,0), colspan=2, rowspan=2)
plot_U_profile(cell, ax, problem_l, 'A')
#Friction coeff comparisons
ax1 = plt.subplot2grid((3,4), (1,2), colspan=2, rowspan=2)
v_wind = np.linspace(0,30,100)
def Im(v): #'Imberger_Patterson'
if v <= 4.:
return 1.124*1.e-3
else:
return (0.94 + 0.041*v)*1.e-3
def St(v): #'Stefan_Ford'
if v <= 15.:
return 0.5*(np.sqrt(v))*1.e-3
else:
return 2.6*1.e-3
def Co(v): #'Coantic'
return (1.+0.03*v)*1.e-3
C1 = np.zeros(len(v_wind))
C2 = np.zeros(len(v_wind))
C3 = np.zeros(len(v_wind))
C4 = np.zeros(len(v_wind))
for C in range(len(C1)):
C1[C] = Im(v_wind[C])
C2[C] = St(v_wind[C])
C3[C] = Co(v_wind[C])
C4[C] = custom_friction_law(v_wind[C])
ax1.plot(v_wind,C1,color=colors[0], label='Imberger_Patterson')
ax1.plot(v_wind,C2,color=colors[1], label='Stefan_Ford')
ax1.plot(v_wind,C3,color=colors[2], label='Coantic')
ax1.plot(v_wind,C4,color=colors[3], label='User defined')
ax1.set_title("Comparison of friction coefficients", fontstyle='italic')
ax1.legend(loc=4, shadow=True, title="Friction law")
ax1.set_xlabel('Wind velocity (m/s)')
ax1.set_ylabel('Cd (m)')
#Mesh and surface vertical velocity
tri = triangular_mesh.triangulation
ax2 = plt.subplot2grid((3,4), (0,0), colspan=4)
ax2.triplot(tri.x, tri.y, tri.trivtx, color='k', lw=0.5)
im = ax2.tricontourf(tri.x, tri.y, tri.trivtx,
problem_l[len(cases)-1].primitives.W[:,int(NL/2)], 20,cmap=plt.cm.jet)
ax2.set_title('Vertical veloctiy on layer: {}'
.format(int(NL/2)))
cbar = fig.colorbar(im, ax=ax2, aspect=8)
cbar.set_label('W (m/s)', y=1.1, labelpad=-25, rotation=0)
ax2.set_xlabel('x (m)')
ax2.set_ylabel('y (m)')
ofs=0.05
#cell
ax2.plot(tri.x[cell], tri.y[cell], marker='o', color='k', markersize=5)
ax2.text(tri.x[cell]+ofs, tri.y[cell]+ofs, "A", color='k', fontsize=15)
plt.tight_layout()
plt.show()
/builds/builds/y49yFuNK/0/freshkiss3d/freshkiss3d/examples/simulations/example_wind.py:156: MatplotlibDeprecationWarning: The seaborn styles shipped by Matplotlib are deprecated since 3.6, as they no longer correspond to the styles shipped by seaborn. However, they will remain available as 'seaborn-v0_8-<style>'. Alternatively, directly use the seaborn API instead.
plt.style.use('seaborn-white')
Total running time of the script: ( 0 minutes 3.183 seconds)