From ea893bcc352ea627b119e28d2fd48e54ae41cae5 Mon Sep 17 00:00:00 2001
From: Paul Dechamps <paul.dechamps@uliege.be>
Date: Sun, 10 Nov 2024 17:59:14 +0100
Subject: [PATCH] (feat) Update API
API only works for 2D cases
---
blast/api/{core.py => blaster_api.py} | 103 +++++++-------
blast/tests/apiTest.py | 197 ++++++++++++++++++++++++++
2 files changed, 250 insertions(+), 50 deletions(-)
rename blast/api/{core.py => blaster_api.py} (53%)
create mode 100644 blast/tests/apiTest.py
diff --git a/blast/api/core.py b/blast/api/blaster_api.py
similarity index 53%
rename from blast/api/core.py
rename to blast/api/blaster_api.py
index 5811adf..3d34d6d 100644
--- a/blast/api/core.py
+++ b/blast/api/blaster_api.py
@@ -19,7 +19,9 @@
## Initialize blast computation
# Paul Dechamps
-def initBlast(cfg, icfg, iSolverName='DART'):
+AVAILABLE_SOLVERS = ['DART']
+
+def init_blaster(cfg, icfg, iSolverName='DART', task='analysis'):
"""
Inputs
------
@@ -55,60 +57,61 @@ def initBlast(cfg, icfg, iSolverName='DART'):
from blast.coupler import Coupler
import blast
- if iSolverName == 'DART':
+ if iSolverName not in AVAILABLE_SOLVERS:
+ raise RuntimeError('BLASTERAPI: Invalid inviscid solver name', iSolverName)
+ elif iSolverName == 'DART':
from blast.interfaces.dart.DartInterface import DartInterface as interface
- # Check viscous solver parameters.
- if 'Re' in cfg and cfg['Re'] > 0:
- _Re = cfg['Re']
- else:
- raise RuntimeError('Missing or invalid Reynolds number')
- if 'Minf' in cfg and cfg['Minf'] > 0:
- _Minf = cfg['Minf']
- else:
- _Minf = 0.1
- if 'CFL0' in cfg and cfg['CFL0'] > 0:
- _CFL0 = cfg['CFL0']
- else:
- _CFL0 = 1
- if 'Verb' in cfg and 0<= cfg['Verb'] <= 3:
- _verbose = cfg['Verb']
- else:
- _verbose = 1
- _xtrF = [-1, -1]
- if 'xtrF' in cfg:
- for i, xtr in enumerate(cfg['xtrF']):
- if not(0 <= xtr <= 1) and xtr != -1:
- raise RuntimeError("Incorrect forced transition location.")
- _xtrF[i] = xtr
- _span = 0
- _nSections = 1
-
- # Check coupler parameters.
- if 'couplIter' in cfg and cfg['couplIter'] > 0:
- __couplIter = cfg['couplIter']
- else:
- __couplIter = 150
- if 'couplTol' in cfg:
- __couplTol = cfg['couplTol']
- else:
- __couplTol = 1e-4
- if 'iterPrint' in cfg:
- __iterPrint = cfg['iterPrint']
- else:
- __iterPrint = 1
- if 'resetInv' in cfg:
- __resetInv = cfg['resetInv']
- else:
- __resetInv = False
+ # Viscous solver
+ if 'Re' not in cfg or cfg['Re'] <= 0.:
+ raise RuntimeError('BLASTERAPI: Missing or invalid Reynolds number')
+ _Re = cfg['Re']
+
+ _Minf = cfg.get('Minf', 0.1)
+ _cfl0 = cfg.get('CFL0', 1.)
+ _verbose = cfg.get('Verb', 0)
+ _xtrF = cfg.get('xtrF', [-1, -1])
+ _span = cfg.get('span', 1.)
+ _nSections = cfg.get('nSections', 1)
+
+ if _Minf <= 0:
+ raise RuntimeError('BLASTERAPI: Invalid Mach number', _Minf)
+ if _cfl0 <= 0:
+ raise RuntimeError('BLASTERAPI: Invalid CFL number', _cfl0)
+ for xtr in _xtrF:
+ if not(0 <= xtr <= 1) and xtr != -1:
+ raise RuntimeError("BLASTERAPI: Incorrect forced transition location.", xtr)
+ if _span < 0:
+ raise RuntimeError('BLASTERAPI: Invalid span', _span)
+ if _nSections <= 0:
+ raise RuntimeError('BLASTERAPI: Invalid number of sections', _nSections)
+
+ # Coupler
+ _couplIter = cfg.get('couplIter', 150)
+ _couplTol = cfg.get('couplTol', 1e-4)
+ _iterPrint = cfg.get('iterPrint', 1)
+ _resetInv = cfg.get('resetInv', False)
+
+ if _couplIter < 0:
+ raise RuntimeError('BLASTERAPI: Invalid number of coupling iterations', _couplIter)
+ if _couplTol <= 0:
+ raise RuntimeError('BLASTERAPI: Invalid coupling tolerance', _couplTol)
+ if _iterPrint < 0:
+ raise RuntimeError('BLASTERAPI: Invalid iteration print frequency', _iterPrint)
# Viscous solver object.
- vSolver = blast.Driver(_Re, _Minf, _CFL0, _nSections, _xtrF[0], _xtrF[1], _span, _verbose =_verbose)
+ vSolver = blast.Driver(_Re, _Minf, _cfl0, _nSections, _xtrF[0], _xtrF[1], _span, _verbose =_verbose)
# Solvers interface.
- solversAPI = interface(icfg, vSolver, cfg)
+ solverAPI = interface(icfg, vSolver, cfg)
# Coupler
- coupler = Coupler(solversAPI, vSolver, _maxCouplIter = __couplIter, _couplTol = __couplTol, _iterPrint = __iterPrint, _resetInv = __resetInv)
+ coupler = Coupler(solverAPI, vSolver, _maxCouplIter = _couplIter, _couplTol = _couplTol, _iterPrint = _iterPrint, _resetInv = _resetInv)
+
+ # Adjoint objects
+ cAdj = None
+ if task == 'optimization':
+ cAdj = blast.CoupledAdjoint(solverAPI.adjointSolver, vSolver)
return {'coupler': coupler,
- 'solversAPI': solversAPI,
- 'vSolver': vSolver}
+ 'isol': solverAPI,
+ 'vsol': vSolver,
+ 'adj': cAdj}
diff --git a/blast/tests/apiTest.py b/blast/tests/apiTest.py
new file mode 100644
index 0000000..3b5bb9b
--- /dev/null
+++ b/blast/tests/apiTest.py
@@ -0,0 +1,197 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+
+# Copyright 2022 University of Liège
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+# http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+
+# @author Paul Dechamps
+# @date 2024
+# Test the blaster adjoint implementation.
+
+# Imports.
+import numpy as np
+
+from fwk.wutils import parseargs
+import fwk
+from fwk.testing import *
+from fwk.coloring import ccolors
+
+import fwk
+from fwk.testing import *
+from fwk.coloring import ccolors
+import numpy as np
+import os
+
+from matplotlib import pyplot as plt
+from blast.api.blaster_api import init_blaster
+import blast.utils as viscUtils
+
+def cfgInviscid(nthrds, verb):
+ import os.path
+ # Parameters
+ return {
+ # Options
+ 'scenario' : 'aerodynamic',
+ 'task' : 'optimization',
+ 'Threads' : nthrds, # number of threads
+ 'Verb' : verb, # verbosity
+ # Model (geometry or mesh)
+ 'File' : os.path.dirname(os.path.abspath(__file__)) + '/../models/dart/n0012.geo', # Input file containing the model
+ 'Pars' : {'xLgt' : 50, 'yLgt' : 50, 'msF' : 10, 'msTe' : 0.01, 'msLe' : 0.01}, # parameters for input file model
+ 'Dim' : 2, # problem dimension
+ 'Format' : 'gmsh', # save format (vtk or gmsh)
+ # Markers
+ 'Fluid' : 'field', # name of physical group containing the fluid
+ 'Farfield' : ['upstream', 'farfield', 'downstream'], # LIST of names of physical groups containing the farfield boundaries (upstream/downstream should be first/last element)
+ 'Wing' : 'airfoil', # LIST of names of physical groups containing the lifting surface boundary
+ 'Wake' : 'wake', # LIST of names of physical group containing the wake
+ 'WakeTip' : 'wakeTip', # LIST of names of physical group containing the edge of the wake
+ 'Te' : 'te', # LIST of names of physical group containing the trailing edge
+ 'dbc' : False,
+ 'Upstream' : 'upstream',
+ # Freestream
+ 'M_inf' : 0.2, # freestream Mach number
+ 'AoA' : 2., # freestream angle of attack
+ # Geometry
+ 'S_ref' : 1., # reference surface length
+ 'c_ref' : 1., # reference chord length
+ 'x_ref' : .25, # reference point for moment computation (x)
+ 'y_ref' : 0.0, # reference point for moment computation (y)
+ 'z_ref' : 0.0, # reference point for moment computation (z)
+ # Numerical
+ 'LSolver' : 'SparseLu', # inner solver (Pardiso, MUMPS or GMRES)
+ 'G_fill' : 2, # fill-in factor for GMRES preconditioner
+ 'G_tol' : 1e-5, # tolerance for GMRES
+ 'G_restart' : 50, # restart for GMRES
+ 'Rel_tol' : 1e-6, # relative tolerance on solver residual
+ 'Abs_tol' : 1e-8, # absolute tolerance on solver residual
+ 'Max_it' : 20, # solver maximum number of iterations
+ }
+
+def cfgBlast(verb):
+ return {
+ 'Re' : 1e6, # Freestream Reynolds number
+ 'Minf' : 0.2, # Freestream Mach number (used for the computation of the time step only)
+ 'CFL0' : 1, # Inital CFL number of the calculation
+ 'Verb': verb, # Verbosity level of the solver
+ 'couplIter': 100, # Maximum number of iterations
+ 'couplTol' : 1e-2, # Tolerance of the VII methodology
+ 'iterPrint': 20, # int, number of iterations between outputs
+ 'resetInv' : True, # bool, flag to reset the inviscid calculation at every iteration.
+ 'sections' : [0], # List of sections for boundary layer calculation
+ 'xtrF' : [-1, -1], # Forced transition location
+ 'interpolator' : 'Matching', # Interpolator for the coupling
+ }
+
+def main():
+
+ tms = fwk.Timers()
+
+ tms['total'].start()
+
+ # Parse agrs
+ args = parseargs()
+ verb = args.verb
+ nthrds = args.k
+
+ icfg = cfgInviscid(nthrds, verb)
+ cfg = cfgBlast(verb)
+
+ obj = init_blaster(cfg, icfg, iSolverName='DART', task='optimization')
+ coupler, isol, vsol, adj = obj['coupler'], obj['isol'], obj['vsol'], obj['adj']
+
+ # Forward problem
+ coupler.run()
+
+ # Adjoint problem
+ isol.adjointSolver.run()
+ adj.run()
+
+ vSolution = viscUtils.getSolution(isol.sec, write=False)[0]
+ #plotSolution(vSolution)
+ tms['total'].stop()
+
+ # Make the api crash
+ params = {
+ 'Re': -1,
+ 'Minf': -1,
+ 'CFL0': -1,
+ 'couplIter': -1,
+ 'couplTol': -1,
+ 'iterPrint': -1
+ }
+
+ # Loop over the parameters
+ for param, value in params.items():
+ original = cfg[param]
+ try:
+ cfg[param] = value
+ obj = init_blaster(cfg, icfg, iSolverName='DART', task='optimization')
+ raise AssertionError(f'API initialized with {param} = ', cfg[param])
+ except AssertionError as e:
+ raise RuntimeError(e)
+ except:
+ print(ccolors.ANSI_GREEN + f'API crashed with {param} = ' + str(cfg[param]) + ccolors.ANSI_RESET)
+ cfg[param] = original
+
+ print('Statistics')
+ print(tms)
+ print(ccolors.ANSI_BLUE + 'PyTesting' + ccolors.ANSI_RESET)
+ tests = CTests()
+ tests.add(CTest('Cl', isol.getCl(), 0.226, 1e-2))
+ tests.add(CTest('Cd', isol.getCd()+vsol.Cdf, 0.00587, 1e-3))
+ tests.add(CTest('dCl_dAoA', adj.tdCl_AoA, 5.47, 1e-4))
+ tests.add(CTest('dCd_dAoA', adj.tdCd_AoA, 0.09531, 1e-3))
+ tests.run()
+
+ # eof
+ print('')
+
+def plotSolution(vSolution):
+ plt.figure()
+ plt.plot(vSolution['x'], vSolution['cf'], label='$c_f$')
+ plt.xlim([0, 1])
+ plt.xlabel('x/c')
+ plt.ylabel('$c_f$')
+ plt.legend()
+ plt.grid()
+ plt.savefig('cf.png')
+ plt.draw()
+
+ plt.figure()
+ plt.plot(vSolution['x'], vSolution['ctEq'], label='$ctEq$')
+ plt.xlim([0, 1])
+ plt.xlabel('x/c')
+ plt.ylabel('$ctEq$')
+ plt.legend()
+ plt.grid()
+ plt.savefig('cf.png')
+ plt.draw()
+
+ plt.figure()
+ plt.plot(vSolution['x'], vSolution['theta'], label='$\theta$')
+ plt.plot(vSolution['x'], vSolution['deltaStar'], label='$\delta^*$')
+ plt.plot(vSolution['x'], vSolution['theta']*vSolution['H'], label='$\delta^*from$')
+ plt.xlim([0, 1])
+ plt.xlabel('x/c')
+ plt.ylabel('$\delta^*$')
+ plt.legend()
+ plt.grid()
+ plt.savefig('cf.png')
+ plt.draw()
+
+
+if __name__ == '__main__':
+ main()
\ No newline at end of file
--
GitLab