緩和係数の最適化2021年6月20日 | |
はじめにOptuna と Bayesian Optimization を用いて OpenFOAM のケースの緩和係数の最適化を試みる。 使用バージョン
計算ケース計算ケースには simpleFoam のチュートリアルケースである pitzDaily を用いる。 Optunaインストール$ pip install optuna テストtest-op.py
import optuna
def objective(trial):
x = trial.suggest_uniform("x", 0, 10)
return (x**2 - 2)**2
study = optuna.create_study()
#optuna.logging.set_verbosity(optuna.logging.INFO)
optuna.logging.set_verbosity(optuna.logging.WARNING)
study.optimize(objective, n_trials=100)
study.best_params
実行 $ python test-op.py 2 の平方根の値が出てくれば OK。 緩和係数の最適化run-op.py
import sys
import subprocess as sp
import re
import optuna
def run_command_with_output(command):
proc = sp.Popen(command,
shell=True, stdout=sp.PIPE, stderr=sp.PIPE)
while True:
line = proc.stdout.readline()
if not line:
break
line = line.decode().rstrip("\n")
yield line
if proc.wait() != 0:
while True:
line = proc.stderr.readline()
if not line:
break
line = line.decode().rstrip("\n")
yield line
raise OSError("Command execution faied: %s" % command)
def write_fvSolution(ap, au, at):
with open("system/fvSolution", "w") as f:
f.write("""FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSolution;
}
solvers
{
p
{
solver GAMG;
tolerance 1e-06;
relTol 0.1;
smoother GaussSeidel;
}
"(U|k|epsilon|omega|f|v2)"
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-05;
relTol 0.1;
}
}
SIMPLE
{
nNonOrthogonalCorrectors 0;
consistent yes;
residualControl
{
p 1e-2;
U 1e-3;
"(k|epsilon|omega|f|v2)" 1e-3;
}
}
relaxationFactors
{
fields
{
""")
f.write(" p %f;\n" % ap)
f.write(""" }
equations
{
""")
f.write(" U %f;\n" % au)
f.write(' ".*" %f;\n' %at)
f.write(""" }
}
""")
def run(ap, au, at):
print(ap, au, at)
write_fvSolution(ap, au, at)
time = 0
for line in run_command_with_output("simpleFoam"):
r = re.search("Foam::error", line)
if r is not None:
return 1000
r = re.search("ExecutionTime = (.+?) ", line)
if r is not None:
time = float(r.group(1))
return time
def objective(trial):
ap = trial.suggest_uniform("ap", 0.1, 1)
au = trial.suggest_uniform("au", 0.1, 1)
at = trial.suggest_uniform("at", 0.1, 1)
return run(ap, au, at)
study = optuna.create_study()
optuna.logging.set_verbosity(optuna.logging.WARNING)
study.optimize(objective, n_trials=100)
print(study.best_params)
ここで、圧力の緩和係数を ap、速度の緩和係数を au、それ以外 (乱流統計量) の緩和係数を at としている。それぞれ探索範囲は 0.1〜1 としている。 実行 $ python run-op.py 結果
{'ap': 0.6937429666961807, 'au': 0.9501927541779454, 'at': 0.7819762799229192}
Bayesian Optimizationインストール$ pip install bayesian-optimization テストtest-bo.py
from bayes_opt import BayesianOptimization
def objective(x):
return -(x**2 - 2)**2
optimizer = BayesianOptimization(
f=objective,
pbounds={'x': (0, 2)},
random_state=1
)
optimizer.maximize(init_points=3, n_iter=20)
print(optimizer.max)
緩和係数の最適化run-bo.py
import sys
import subprocess as sp
import re
from bayes_opt import BayesianOptimization
def run_command_with_output(command):
proc = sp.Popen(command,
shell=True, stdout=sp.PIPE, stderr=sp.PIPE)
while True:
line = proc.stdout.readline()
if not line:
break
line = line.decode().rstrip("\n")
yield line
if proc.wait() != 0:
while True:
line = proc.stderr.readline()
if not line:
break
line = line.decode().rstrip("\n")
yield line
raise OSError("Command execution faied: %s" % command)
def write_fvSolution(ap, au, at):
with open("system/fvSolution", "w") as f:
f.write("""FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object fvSolution;
}
solvers
{
p
{
solver GAMG;
tolerance 1e-06;
relTol 0.1;
smoother GaussSeidel;
}
"(U|k|epsilon|omega|f|v2)"
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-05;
relTol 0.1;
}
}
SIMPLE
{
nNonOrthogonalCorrectors 0;
consistent yes;
residualControl
{
p 1e-2;
U 1e-3;
"(k|epsilon|omega|f|v2)" 1e-3;
}
}
relaxationFactors
{
fields
{
""")
f.write(" p %f;\n" % ap)
f.write(""" }
equations
{
""")
f.write(" U %f;\n" % au)
f.write(' ".*" %f;\n' %at)
f.write(""" }
}
""")
def run(ap, au, at):
write_fvSolution(ap, au, at)
time = 0
for line in run_command_with_output("simpleFoam"):
r = re.search("Foam::error", line)
if r is not None:
return 1000
r = re.search("ExecutionTime = (.+?) ", line)
if r is not None:
time = float(r.group(1))
return time
def objective(ap, au, at):
return -run(ap, au, at)
optimizer = BayesianOptimization(
f=objective,
pbounds={"ap": (0.1, 1), "au": (0.1, 1), "at": (0.1, 1)},
random_state=1
)
optimizer.maximize(init_points=10, n_iter=100)
print(optimizer.max)
結果
{'target': -5.57, 'params': {'ap': 0.9821591542366661, 'at': 0.8978768367274595, 'au': 0.9 428528721941154}}
比較オリジナル (ap = 1、au = at = 0.9) ExecutionTime = 5.97 s ClockTime = 6 s SIMPLE solution converged in 280 iterations Optuna (ap = 0.694、au = 0.950、at = 0.782) ExecutionTime = 5.2 s ClockTime = 5 s SIMPLE solution converged in 215 iterations Bayesian Optimization (ap = 0.982、au = 0.943、at = 0.898) ExecutionTime = 5.4 s ClockTime = 5 s SIMPLE solution converged in 217 iterations | |
| PENGUINITIS | |