Note
Click here to download the full example code
Vertical flowrate¶
In this example VerticalFlowRate boundary condition is illustrated. Geometry
is defined so that every side is a wall. At the bottom we set two vertical flowrate
boundary conditions, the first being an inlet and the second an outlet.
import os, sys
import matplotlib.pylab as plt
import numpy as np
import freshkiss3d as fk
import freshkiss3d.extra.plots as fk_plt
os.system('rm -r outputs')
#sphinx_gallery_thumbnail_number = 2
0
Time loop:¶
simutime = fk.SimuTime(final_time=1, time_iteration_max=10000,
second_order=True)
Mesh:¶
dir_path = os.path.abspath(os.path.dirname(sys.argv[0]))
triangular_mesh = fk.TriangularMesh.from_msh_file(dir_path + '/inputs/simple_canal_2.msh')
fk_plt.plot_mesh(triangular_mesh)
Layers:¶
NL = 4
layer = fk.Layer(NL, triangular_mesh, topography=0.)
Primitives:¶
primitives = fk.Primitive(triangular_mesh, layer, free_surface=.5)
Tracer:¶
tracers = [fk.Tracer(triangular_mesh, layer, primitives, Tinit=0.0, label='temperature')]
Boundary conditions:¶
VerticalFlowRate boundary conditions consist of imposing vertical flowrate at the
bottom on designated cells. In case of inlet you can also impose a Tracer value.
If more than one cell is selected, imposed flowrate is shared based on cell area.
As any other boundary condition, there can be as much VerticalFlowRate as needed
as long as they are defined in a list.
slides = [fk.Slide(ref=r) for r in [1, 2, 3, 4]]
vertical_flowrate_in = fk.VerticalFlowRate(
times=[0., 10, 20],
flowrates=[0.1, 0.1, 0.1],
cells=[97],
tracers_inlet={'temperature':10.})
vertical_flowrate_out = fk.VerticalFlowRate(
times=[0., 10, 20],
flowrates=[-0.1, -0.1, -0.1],
cells=[94],
tracers_inlet={'temperature':10.})
vertical_flowrates = [vertical_flowrate_in, vertical_flowrate_out]
NUM_PARAMS={'ipres':True,
'flux_type':1,
'implicit_exchanges':True}
Writter:¶
NB_VTK = 10
vtk_writer = fk.VTKWriter(triangular_mesh, scheduler=fk.schedules(count=NB_VTK),
scale_h=2.)
Problem definition:¶
problem = fk.Problem(simutime, triangular_mesh, layer, primitives,
tracers=tracers,
slides=slides,
vertical_flowrates=vertical_flowrates,
numerical_parameters=NUM_PARAMS,
vtk_writer=vtk_writer)
===================================================================
| INITIALIZATION |
===================================================================
Problem size: Nodes=144, Layers=4, Triangles=240,
Iter = 0 ; Dt = 0.0000s ; Time = 0.00s ; ETA = 0.00s
Problem solving:¶
problem.solve()
===================================================================
| TIME LOOP |
===================================================================
Iter = 7 ; Dt = 0.0033s ; Time = 0.02s ; ETA = 0.84s
Iter = 13 ; Dt = 0.0033s ; Time = 0.04s ; ETA = 0.80s
Iter = 19 ; Dt = 0.0033s ; Time = 0.06s ; ETA = 0.79s
Iter = 25 ; Dt = 0.0033s ; Time = 0.08s ; ETA = 0.78s
Iter = 32 ; Dt = 0.0033s ; Time = 0.10s ; ETA = 0.76s
Iter = 38 ; Dt = 0.0033s ; Time = 0.12s ; ETA = 0.74s
Iter = 44 ; Dt = 0.0033s ; Time = 0.14s ; ETA = 0.70s
Iter = 50 ; Dt = 0.0033s ; Time = 0.16s ; ETA = 0.69s
Iter = 56 ; Dt = 0.0033s ; Time = 0.18s ; ETA = 0.69s
Iter = 63 ; Dt = 0.0033s ; Time = 0.21s ; ETA = 0.66s
Iter = 69 ; Dt = 0.0033s ; Time = 0.23s ; ETA = 5.84s
Iter = 75 ; Dt = 0.0033s ; Time = 0.25s ; ETA = 0.63s
Iter = 81 ; Dt = 0.0033s ; Time = 0.27s ; ETA = 0.62s
Iter = 88 ; Dt = 0.0033s ; Time = 0.29s ; ETA = 0.60s
Iter = 94 ; Dt = 0.0033s ; Time = 0.31s ; ETA = 0.58s
Iter = 100 ; Dt = 0.0033s ; Time = 0.33s ; ETA = 0.57s
Iter = 106 ; Dt = 0.0033s ; Time = 0.35s ; ETA = 0.56s
Iter = 113 ; Dt = 0.0033s ; Time = 0.37s ; ETA = 0.52s
Iter = 119 ; Dt = 0.0033s ; Time = 0.39s ; ETA = 0.52s
Iter = 125 ; Dt = 0.0033s ; Time = 0.41s ; ETA = 0.50s
Iter = 131 ; Dt = 0.0033s ; Time = 0.43s ; ETA = 0.48s
Iter = 137 ; Dt = 0.0033s ; Time = 0.45s ; ETA = 4.18s
Iter = 144 ; Dt = 0.0033s ; Time = 0.47s ; ETA = 0.45s
Iter = 150 ; Dt = 0.0033s ; Time = 0.49s ; ETA = 0.44s
Iter = 156 ; Dt = 0.0033s ; Time = 0.51s ; ETA = 0.41s
Iter = 162 ; Dt = 0.0033s ; Time = 0.53s ; ETA = 0.39s
Iter = 169 ; Dt = 0.0033s ; Time = 0.55s ; ETA = 0.38s
Iter = 175 ; Dt = 0.0033s ; Time = 0.57s ; ETA = 0.36s
Iter = 181 ; Dt = 0.0033s ; Time = 0.59s ; ETA = 0.35s
Iter = 187 ; Dt = 0.0033s ; Time = 0.61s ; ETA = 0.33s
Iter = 193 ; Dt = 0.0033s ; Time = 0.63s ; ETA = 0.30s
Iter = 200 ; Dt = 0.0033s ; Time = 0.66s ; ETA = 0.29s
Iter = 206 ; Dt = 0.0033s ; Time = 0.68s ; ETA = 0.27s
Iter = 212 ; Dt = 0.0033s ; Time = 0.70s ; ETA = 0.26s
Iter = 218 ; Dt = 0.0033s ; Time = 0.71s ; ETA = 0.25s
Iter = 225 ; Dt = 0.0033s ; Time = 0.74s ; ETA = 0.22s
Iter = 231 ; Dt = 0.0033s ; Time = 0.76s ; ETA = 0.21s
Iter = 237 ; Dt = 0.0033s ; Time = 0.78s ; ETA = 0.19s
Iter = 243 ; Dt = 0.0033s ; Time = 0.80s ; ETA = 0.17s
Iter = 249 ; Dt = 0.0033s ; Time = 0.82s ; ETA = 0.16s
Iter = 256 ; Dt = 0.0033s ; Time = 0.84s ; ETA = 0.14s
Iter = 262 ; Dt = 0.0033s ; Time = 0.86s ; ETA = 0.12s
Iter = 268 ; Dt = 0.0033s ; Time = 0.88s ; ETA = 0.10s
Iter = 274 ; Dt = 0.0033s ; Time = 0.90s ; ETA = 0.09s
Iter = 281 ; Dt = 0.0033s ; Time = 0.92s ; ETA = 0.07s
Iter = 287 ; Dt = 0.0033s ; Time = 0.94s ; ETA = 0.05s
Iter = 293 ; Dt = 0.0033s ; Time = 0.96s ; ETA = 0.03s
Iter = 299 ; Dt = 0.0033s ; Time = 0.98s ; ETA = 0.02s
Iter = 306 ; Dt = 0.0001s ; Time = 1.00s ; ETA = 0.00s
===================================================================
| END |
===================================================================
Problem.solve() completed in 1.0554919242858887s (wall time)
Plots:¶
Velocity norm is plotted for each layer:
vertex_labels = triangular_mesh.vertex_labels
boundary_labels = triangular_mesh.boundary_labels
x = np.asarray(triangular_mesh.triangulation.x)
y = np.asarray(triangular_mesh.triangulation.y)
trivtx = np.asarray(triangular_mesh.triangulation.trivtx)
NC = len(x)
index = range(NL)
velocity = np.zeros((NC))
fig, axes = plt.subplots(nrows=len(index)) #, ncols=2)
v = np.linspace(-.1, .1, 10, endpoint=True)
for L, ax in zip(reversed(index), axes.flat):
lab = 'L={}'.format(L)
ax.triplot(x, y, trivtx, color='k', lw=0.5)
# VELOCITY:
for C in range(NC):
velocity[C] = problem.primitives.W[C,L]
# velocity[C] = np.sqrt( problem.primitives.U[C,L]**2 + \
# problem.primitives.V[C,L]**2 + \
# problem.primitives.W[C,L]**2 )
im = ax.tricontourf(x, y, trivtx, velocity[:], v, cmap=plt.cm.bwr)
cbar = fig.colorbar(im, ax=ax, aspect=8, ticks=v)
cbar.set_label(lab, rotation=90)
axes[0].set_title('Vertical veloctiy (W)'.format(L))
for L in range(len(index)-1): plt.setp(axes[L].get_xticklabels(), visible=False)
for a, b, label in zip(x, y, vertex_labels):
if a < 0.45 and a > 0.35 and b < 0.3 and b > 0.2:
axes[len(index)-1].plot(a, b, marker='o', color='k', markersize=5)
if a < 1.8 and a > 1.7 and b < 0.3 and b > 0.2:
axes[len(index)-1].plot(a, b, marker='o', color='k', markersize=5)
axes[len(index)-1].annotate('VerticalFlowRate in', xy=(0.37, 0.25), xytext=(0.2, -0.3),
arrowprops=dict(facecolor='black', width=0.2, headwidth=5, shrink=0.1))
axes[len(index)-1].annotate('VerticalFlowRate out', xy=(1.73, 0.25), xytext=(1.3, -0.3),
arrowprops=dict(facecolor='black', width=0.2, headwidth=5, shrink=0.1))
plt.show()
Total running time of the script: ( 0 minutes 1.826 seconds)