1.7.2.1. PitzDaily¶
[3]:
import sys, time
import numpy as np
seed = 1435
np.random.seed(seed)
# import problem classes
import Exeter_CFD_Problems as TestProblems
verbose = False
sys.argv = sys.argv[:1]
print('Demonstration of the PitzDialy test problem.')
# set up directories.
settings = {
'source_case': 'Exeter_CFD_Problems/data/PitzDaily/case_fine/',
'case_path': 'Exeter_CFD_Problems/data/PitzDaily/case_single/',
'boundary_files': ['Exeter_CFD_Problems/data/PitzDaily/boundary.csv'],
'fixed_points_files': ['Exeter_CFD_Problems/data/PitzDaily/fixed.csv']
}
# instantiate the problem object
prob = TestProblems.PitzDaily(settings)
# get the lower and upper bounds
lb, ub = prob.get_decision_boundary()
# generate random solutions satisfying the lower and upper bounds.
x = np.random.random((1000, lb.shape[0])) * (ub - lb) + lb
rand_x = []
for i in range(x.shape[0]):
if prob.constraint(x[i]): # check to see if the random solution is valid
rand_x.append(x[i])
# evaluate for a solution
print('Number of control points: ', prob.n_control)
print('Decision vector: ', rand_x[0])
print('Running simulation ...')
start = time.time()
res = prob.evaluate(rand_x[0], verbose=verbose)
print('Objective function value:', res)
print('Time taken:', time.time()-start, ' seconds.')
Demonstration of the PitzDialy test problem.
Number of control points: [5]
Decision vector: [ 0.11063283 -0.03252882 0.0333295 -0.02732676 0.07319029 0.00579979
0.08095974 -0.04046204 0.06019682 -0.02585734]
Running simulation ...
Objective function value: 0.203349406523
Time taken: 106.35813975334167 seconds.
[4]:
print(max(lb),min(lb))
-0.01 -0.05
[5]:
print(max(ub),min(ub))
0.287397 0.014
[6]:
x
[6]:
array([[ 5.78530692e-02, -3.31068377e-02, 5.50604920e-02, ...,
-4.71955674e-02, 1.88146599e-01, -1.32952546e-02],
[-5.67533117e-03, -7.48384848e-03, 1.38897278e-01, ...,
-3.25892854e-02, 2.22551915e-02, -7.50087530e-03],
[ 2.71721553e-02, -4.34923342e-02, 2.63459670e-01, ...,
-1.55789562e-02, 7.12139945e-02, -1.24998786e-02],
...,
[ 2.62482277e-01, 8.32653794e-03, 1.36335014e-01, ...,
-6.39175269e-03, 1.56653572e-01, -2.76899372e-02],
[ 2.32482592e-01, -2.66450403e-02, 1.01493483e-01, ...,
-2.29677976e-02, -8.32203134e-03, 1.86693078e-04],
[ 2.35917506e-02, -1.43564186e-02, -9.60306583e-03, ...,
-3.61561063e-02, 2.43942977e-01, -2.07163078e-03]])
[7]:
! ls
Exeter_CFD_Problems plotHeatExchangerTemperature.py
figures plotKaplanDuct.py
HeatExchanger.ipynb plotPitzDailyGrid.py
KaplanDuct.ipynb plotPitzDailyStreamlines.py
LICENSE plotPitzDailyVelocity.py
log PlyParser_FoamFileParser_parsetab.py
log2 __pycache__
log3 README.md
log4 sample_script.py
PitzDaily.ipynb U_centerplane.vtk
plotHeatExchanger.py
[8]:
! postProcess -case Exeter_CFD_Problems/data/PitzDaily/case_single -func sample > log 2>&1
[9]:
# I have run this case from cases_fine to see what the results are for non_optimised case
! postProcess -case Exeter_CFD_Problems/data/PitzDaily/case_non_optimised -func sample > log2 2>&1
[24]:
import plotPitzDailyVelocity
import importlib
importlib.reload(plotPitzDailyVelocity) # have to reload module if changes are made
%matplotlib inline
filename = \
'Exeter_CFD_Problems/data/PitzDaily/case_non_optimised/postProcessing/sample/1634/U_nearWall.vtk'
plotPitzDailyVelocity.plot(filename)
[25]:
filename_2 = \
'Exeter_CFD_Problems/data/PitzDaily/case_single/postProcessing/sample/0/U_nearWall.vtk'
plotPitzDailyVelocity.plot(filename_2)
1.7.2.2. The Triangulated Grid¶
[17]:
import plotPitzDailyGrid
import importlib
importlib.reload(plotPitzDailyGrid) # have to reload module if changes are made
%matplotlib inline
[18]:
! foamToVTK -case Exeter_CFD_Problems/data/PitzDaily/case_non_optimised > log3 2>&1
[19]:
grid = \
'Exeter_CFD_Problems/data/PitzDaily/case_non_optimised/VTK/case_non_optimised_1634.vtk'
[20]:
plotPitzDailyGrid.plot(grid)
[21]:
! foamToVTK -case Exeter_CFD_Problems/data/PitzDaily/case_single > log4 2>&1
[22]:
grid_2 = \
'Exeter_CFD_Problems/data/PitzDaily/case_single/VTK/case_single_0.vtk'
[23]:
plotPitzDailyGrid.plot(grid_2)
1.7.2.3. Turbulucid¶
[1]:
import warnings; warnings.simplefilter('ignore')
[2]:
%matplotlib inline
[3]:
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
from turbulucid import *
matplotlib.rcParams['figure.figsize'] = (15, 10)
[4]:
matplotlib.rcParams['axes.labelsize'] = 15
matplotlib.rcParams['xtick.labelsize'] = 15
matplotlib.rcParams['ytick.labelsize'] = 15
matplotlib.rcParams['text.usetex'] = True
[5]:
from os.path import join
import turbulucid
import importlib
importlib.reload(turbulucid) # have to reload module if changes are made
case = Case(join("./Exeter_CFD_Problems/data/PitzDaily/", \
"case_single","postProcessing","sample", "0", "U_nearWall.vtk"),pointData=True)
[6]:
print(case.fields)
['U']
[7]:
case.xlim[0]
[7]:
-0.025776666930589998
[8]:
case.ylim[0]
[8]:
-0.026246666349433335
[9]:
case.xlim[1]
[9]:
0.29517665818919003
[10]:
max(case.cellCentres[0])
[10]:
0.00020892665
[11]:
case.cellCentres
[11]:
array([[-1.9027714e-02, 2.0892665e-04],
[-1.9027714e-02, 1.0446333e-04],
[-2.0072391e-02, 3.1876573e-04],
...,
[ 2.6091382e-01, -1.9279622e-02],
[ 2.5608581e-01, -1.9622434e-02],
[ 2.5503725e-01, -1.9715633e-02]], dtype=float32)
[12]:
case['U']
[12]:
array([[ 9.8294477e+00, -8.0274665e-01, 5.5096077e-25],
[ 6.5214849e+00, -7.0288688e-01, 4.2253022e-25],
[ 9.9577675e+00, -2.6880601e-01, -3.0189468e-24],
...,
[ 1.3043121e+00, 1.0591840e-01, 0.0000000e+00],
[ 1.2289770e+00, 7.9705156e-02, -2.8791849e-19],
[ 1.2289770e+00, 7.9705156e-02, -2.8791849e-19]], dtype=float32)
[13]:
case["magU"] = np.linalg.norm(case["U"], axis=1)
[14]:
case["magU"]
[14]:
array([9.86217213, 6.55925417, 9.96139526, ..., 1.30860567, 1.2315588 ,
1.2315588 ])
[212]:
plot=plot_field(case, "magU", cmap='rainbow',colorbar=False)
plt.xlim([-0.025,0.2])
plt.ylim([-0.03,0.03])
plt.xlabel(r'$x \ (m)$')
plt.ylabel(r'$y \ (m)$')
bar=add_colorbar(plot, aspect=10, padFraction=1)
bar.set_label(r'$Velocity \ Magnitude \ (m/s)$')
[216]:
plot2=plot_vectors(case,'U', colorField='magU',cmap='rainbow', normalize=True, \
sampleByPlane=True,planeResolution=[20, 80], width=0.0015)
plt.xlim([-0.025,0.2])
plt.ylim([-0.03,0.03])
plt.xlabel(r'$x \ (m)$')
plt.ylabel(r'$y \ (m)$')
bar2=add_colorbar(plot2, aspect=10, padFraction=1)
bar2.set_label(r'$Velocity \ Magnitude \ (m/s)$')
[17]:
plot3=plot_streamlines(case, 'U',colorField='magU',cmap='rainbow',planeResolution=[100,100],density=2)
plt.xlim([-0.025,0.2])
plt.ylim([-0.03,0.03])
plt.xlabel(r'$x \ (m)$')
plt.ylabel(r'$y \ (m)$')
bar3=add_colorbar(plot3.lines, aspect=10, padFraction=1)
bar3.set_label(r'$Velocity \ Magnitude \ (m/s)$')
[19]:
case2 = Case(join("./Exeter_CFD_Problems/data/PitzDaily/", \
"case_non_optimised","postProcessing","sample", "1634", "U_nearWall.vtk"),pointData=True)
[20]:
case2["magU"] = np.linalg.norm(case2["U"], axis=1)
[22]:
plot4=plot_streamlines(case2, 'U',colorField='magU',cmap='rainbow',planeResolution=[90,90],density=2)
plt.xlim([-0.025,0.2])
plt.ylim([-0.03,0.03])
plt.xlabel(r'$x \ (m)$')
plt.ylabel(r'$y \ (m)$')
bar4=add_colorbar(plot4.lines, aspect=10, padFraction=1)
bar4.set_label(r'$Velocity \ Magnitude \ (m/s)$')
[ ]: