Nozzle
In [1]:
Copied!
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline
import flowdyn.mesh as mesh
import flowdyn.xnum as xnum
import flowdyn.integration as tnum
import flowdyn.modelphy.euler as euler
import flowdyn.modeldisc as modeldisc
#
plt.rcParams['figure.dpi'] = 120
plt.rcParams['savefig.dpi'] = 120
plt.rcParams["animation.html"] = "jshtml" # for matplotlib 2.1 and above, uses JavaScript
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline
import flowdyn.mesh as mesh
import flowdyn.xnum as xnum
import flowdyn.integration as tnum
import flowdyn.modelphy.euler as euler
import flowdyn.modeldisc as modeldisc
#
plt.rcParams['figure.dpi'] = 120
plt.rcParams['savefig.dpi'] = 120
plt.rcParams["animation.html"] = "jshtml" # for matplotlib 2.1 and above, uses JavaScript
In [2]:
Copied!
nx = 200
meshsim = mesh.unimesh(ncell=nx, length=10.)
def S(x):
return 1.-.5*np.exp(-.5*(x-5.)**2)
plt.plot(meshsim.centers(), S(meshsim.centers())) ; plt.ylim(0,1)
model = euler.nozzle(sectionlaw=S)
bcL = { 'type': 'insub', 'ptot': 1.1, 'rttot': 1. }
bcR = { 'type': 'outsub', 'p': 1. }
monitors = {
'residual':{ 'name':'L2 residuals', 'frequency': 5 },
'data_average':{ 'data': 'mach', 'name':'Mach average', 'frequency': 1 }
}
rhs = modeldisc.fvm(model, meshsim, xnum.muscl(xnum.vanalbada), bcL=bcL, bcR=bcR)
solver = tnum.rk3ssp(meshsim, rhs, monitors=monitors)
# computation
#
nsol = 100 # number of intermediate resultats (for animation)
endtime = 100. # final physical time
cfl = 1.
finit = rhs.fdata_fromprim([ 1., 0.1, 1. ]) # rho, u, p
#finit = fsol[-1] ; finit.set_time(0)
fsol = solver.solve(finit, cfl, np.linspace(0., endtime, nsol+1))
solver.show_perf()
nx = 200
meshsim = mesh.unimesh(ncell=nx, length=10.)
def S(x):
return 1.-.5*np.exp(-.5*(x-5.)**2)
plt.plot(meshsim.centers(), S(meshsim.centers())) ; plt.ylim(0,1)
model = euler.nozzle(sectionlaw=S)
bcL = { 'type': 'insub', 'ptot': 1.1, 'rttot': 1. }
bcR = { 'type': 'outsub', 'p': 1. }
monitors = {
'residual':{ 'name':'L2 residuals', 'frequency': 5 },
'data_average':{ 'data': 'mach', 'name':'Mach average', 'frequency': 1 }
}
rhs = modeldisc.fvm(model, meshsim, xnum.muscl(xnum.vanalbada), bcL=bcL, bcR=bcR)
solver = tnum.rk3ssp(meshsim, rhs, monitors=monitors)
# computation
#
nsol = 100 # number of intermediate resultats (for animation)
endtime = 100. # final physical time
cfl = 1.
finit = rhs.fdata_fromprim([ 1., 0.1, 1. ]) # rho, u, p
#finit = fsol[-1] ; finit.set_time(0)
fsol = solver.solve(finit, cfl, np.linspace(0., endtime, nsol+1))
solver.show_perf()
cpu time computation (3944 it) : 5.704s 7.23 µs/cell/it
In [3]:
Copied!
# Figure / Plot of final results
monitor_names = ['data_average', 'residual']
fig, ax = plt.subplots(1, 2, figsize=(10,4))
ax[0].grid(linestyle='--', color='0.5')
monitors['data_average']['output'].plot_it(ax=ax[0])
ax[1].grid(linestyle='--', color='0.5')
monitors['residual']['output'].semilogplot_it(ax=ax[1])
# Figure / Plot of final results
monitor_names = ['data_average', 'residual']
fig, ax = plt.subplots(1, 2, figsize=(10,4))
ax[0].grid(linestyle='--', color='0.5')
monitors['data_average']['output'].plot_it(ax=ax[0])
ax[1].grid(linestyle='--', color='0.5')
monitors['residual']['output'].semilogplot_it(ax=ax[1])
In [4]:
Copied!
# Figure / Plot of final results
varname=['pressure', 'mach']
line = varname[:]
fig, ax = plt.subplots(1, 2, figsize=(10,4))
for i, var in enumerate(varname):
ax[i].set_ylabel(var) ; ax[i].set_ylim(0., 1.1*np.max(fsol[-1].phydata(var)))
ax[i].grid(linestyle='--', color='0.5')
line[i], = fsol[-1].plot(var, 'k-', axes=ax[i])
# Figure / Plot of final results
varname=['pressure', 'mach']
line = varname[:]
fig, ax = plt.subplots(1, 2, figsize=(10,4))
for i, var in enumerate(varname):
ax[i].set_ylabel(var) ; ax[i].set_ylim(0., 1.1*np.max(fsol[-1].phydata(var)))
ax[i].grid(linestyle='--', color='0.5')
line[i], = fsol[-1].plot(var, 'k-', axes=ax[i])
In [5]:
Copied!
import matplotlib.animation as anim
#
def animate(k):
#i = min(k, nsol)
for i, var in zip(range(len(varname)),varname):
fsol[k].set_plotdata(line[i], var)
return line
ani = anim.FuncAnimation(fig=fig, func=animate, frames=range(nsol+1), interval=100)#, blit=True)
ani
import matplotlib.animation as anim
#
def animate(k):
#i = min(k, nsol)
for i, var in zip(range(len(varname)),varname):
fsol[k].set_plotdata(line[i], var)
return line
ani = anim.FuncAnimation(fig=fig, func=animate, frames=range(nsol+1), interval=100)#, blit=True)
ani
Out[5]: