From 7c541bf7dc0902935b247b4cbd700f2116d98a25 Mon Sep 17 00:00:00 2001
From: Paul Dechamps <paul.dechamps@uliege.be>
Date: Sun, 10 Nov 2024 17:54:47 +0100
Subject: [PATCH] (feat) Update agard wing validation case

---
 blast/validation/agardValidation.py | 76 ++++++++++++++++-------------
 1 file changed, 43 insertions(+), 33 deletions(-)

diff --git a/blast/validation/agardValidation.py b/blast/validation/agardValidation.py
index a97d4d4..bcc84bd 100644
--- a/blast/validation/agardValidation.py
+++ b/blast/validation/agardValidation.py
@@ -43,7 +43,7 @@ def cfgInviscid(nthrds, verb):
     'Verb' : verb, # verbosity
     # Model (geometry or mesh)
     'File' : os.path.dirname(os.path.dirname(os.path.dirname(os.path.abspath(__file__)))) + '/modules/dartflo/dart/models/agard445.geo', # Input file containing the model
-    'Pars' : {'xL': 7., 'yL': 3., 'zL': 6., 'xO': -3, 'zO': -3, 'msLeRt': 0.0028, 'msTeRt': 0.0056, 'msLeTp': 0.0018, 'msTeTp': 0.0036, 'msF': 1.0}, # parameters for input file model
+    'Pars' : {'xL': 22., 'yL': 22., 'zL': 22., 'xO': -11., 'zO': -11., 'msLeRt': 0.0056, 'msTeRt': 0.0056, 'msLeTp': 0.0036, 'msTeTp': 0.0036, 'msF': 1.0}, # parameters for input file model
     'Dim' : 3, # problem dimension
     'Format' : 'vtk', # save format (vtk or gmsh)
     # Markers
@@ -57,7 +57,7 @@ def cfgInviscid(nthrds, verb):
     'Upstream' : 'upstream',
     # Freestream
     'M_inf' : 0.96,     # freestream Mach number
-    'AoA' : 1,         # freestream angle of attack
+    'AoA' : 0.,         # freestream angle of attack
     # Geometry
     'S_ref' : 0.35,    # reference surface length
     'c_ref' : 0.47,      # reference chord length
@@ -71,16 +71,16 @@ def cfgInviscid(nthrds, verb):
     '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' : 50        # solver maximum number of iterations
+    'Max_it' : 75        # solver maximum number of iterations
     }
 
 def cfgBlast(verb):
     return {
-        'Re' : 6.7e6,       # Freestream Reynolds number
+        'Re' : 5.96e5,       # Freestream Reynolds number
         'Minf' : 0.96,     # Freestream Mach number (used for the computation of the time step only)
         'CFL0' : 1,         # Inital CFL number of the calculation
 
-        'sections' : np.linspace(0.026, 0.7, 15),
+        'sections' : np.linspace(0.026, 0.73, 15),
         'writeSections': np.linspace(0.01, 0.76, 15),
         'Sym':[0.],
         'span': 0.762,
@@ -88,7 +88,7 @@ def cfgBlast(verb):
         'rbftype': 'linear',
         'smoothing': 1e-8,
         'degree': 0,
-        'neighbors': 10,
+        'neighbors': 6,
         'saveTag': 5,
 
         'Verb': verb,       # Verbosity level of the solver
@@ -102,6 +102,7 @@ def cfgBlast(verb):
 
 def main():
     # Timer.
+    import os
     tms = fwk.Timers()
     tms['total'].start()
 
@@ -109,19 +110,21 @@ def main():
     icfg = cfgInviscid(args.k, args.verb)
     vcfg = cfgBlast(args.verb)
 
-    AoAVec = [-1, 0., 1]
-    sfxVec = ['_AoAminus1deg', '_AoA0deg', '_AoA1deg']
-    for i in range(3):
+    #AoAVec = [-.01, 0.0, .01]
+    #sfxVec = ['_a1-', '_a0', '_a1']
+    AoAVec = [0.0]
+    sfxVec = ['_a0']
+    for i in range(len(AoAVec)):
         vcfg['sfx'] = sfxVec[i]
         icfg['AoA'] = AoAVec[i]
 
-        parsViscous = {'nLe': 15, 'nMid': 30, 'nTe': 7, 'nSpan': 60, 'nWake': 30,
-                    'progLe': 1.07, 'progMid': 1.0,  'progTe': 1.0, 'progSpan': 1.0, 'progWake': 1.15}
+        parsViscous = {'nLe': 30, 'nMid': 45, 'nTe': 10, 'nSpan': 60, 'nWake': 25,
+                    'progLe': 1.1, 'progMid': 1.0,  'progTe': 1.0, 'progSpan': 1.0, 'progWake': 1.2}
         vMsh = viscUtils.mesh(os.path.dirname(os.path.dirname(os.path.abspath(__file__))) + '/models/dart/agard445_visc.geo', parsViscous)
         vcfg['vMsh'] = vMsh
 
         tms['pre'].start()
-        coupler, iSolverAPI, vSolver = viscUtils.initBlast(icfg, vcfg)
+        coupler, isol, vsol = viscUtils.initBlast(icfg, vcfg)
         tms['pre'].stop()
 
         print(ccolors.ANSI_BLUE + 'PySolving...' + ccolors.ANSI_RESET)
@@ -132,17 +135,11 @@ def main():
         # Display results.
         print(ccolors.ANSI_BLUE + 'PyRes...' + ccolors.ANSI_RESET)
         print('      Re        M    alpha       Cl       Cd      Cdp      Cdf       Cm')
-        print('{0:6.1f}e6 {1:8.2f} {2:8.1f} {3:8.4f} {4:8.4f} {5:8.4f} {6:8.4f} {7:8.4f}'.format(vcfg['Re']/1e6, iSolverAPI.getMinf(), iSolverAPI.getAoA()*180/math.pi, iSolverAPI.getCl(), vSolver.Cdt, vSolver.Cdp, vSolver.Cdf, iSolverAPI.getCm()))
+        print('{0:6.1f}e6 {1:8.2f} {2:8.1f} {3:8.4f} {4:8.4f} {5:8.4f} {6:8.4f} {7:8.4f}'.format(vcfg['Re']/1e6, isol.getMinf(), isol.getAoA()*180/math.pi, isol.getCl(), vsol.Cdt, vsol.Cdp, vsol.Cdf, isol.getCm()))
 
         # Write results to file.
-        vSolution = viscUtils.getSolution(vSolver)
-
-        # Write results to file.
-        for iSec in range(len(iSolverAPI.cfg['EffSections'])):
-            vSolution = viscUtils.getSolution(vSolver, iSec)
-            viscUtils.write(vSolution, vSolver.getRe(), sfx=sfxVec[i]+'slice'+str(iSec))
-        vSolution['Cdt_int'] = vSolver.Cdf + iSolverAPI.getCd()
-        tms['total'].stop()
+        isol.save(sfx='_viscous'+sfxVec[i])
+        vSolution = viscUtils.getSolution(isol.sec, write=True, toW='all')
 
         print(ccolors.ANSI_BLUE + 'PyTiming...' + ccolors.ANSI_RESET)
         print('CPU statistics')
@@ -150,30 +147,43 @@ def main():
         print('SOLVERS statistics')
         print(coupler.tms)
 
-    cps = []
+    cps_v = []
+    cps_i = []
+    import os
+    # Get file dir
+    dir = os.path.dirname(os.path.dirname(os.path.dirname(os.path.abspath(__file__)))) + '/workspace/blast_validation_agardValidation/'
+
     for s in sfxVec:
-        cps_s = []
+        cps_s_v = []
+        cps_s_i = []
         for i in range(len(vcfg['writeSections'])):
-            cps_s.append(np.loadtxt('/Users/pauldechamps/lab/softwares/blaster/workspace/blast_validation_agardValidation/agard445_viscous'+s+'_slice_'+str(i)+'.dat', delimiter=',', skiprows=1))
-        cps.append(cps_s)
+            cps_s_v.append(np.loadtxt(dir+'agard445_viscous'+s+'_slice_'+str(i)+'.dat', delimiter=',', skiprows=1))
+            cps_s_i.append(np.loadtxt(dir+'agard445_inviscid'+s+'_slice_'+str(i)+'.dat', delimiter=',', skiprows=1))
+        cps_v.append(cps_s_v)
+        cps_i.append(cps_s_i)
     # Plotting
     from matplotlib import pyplot as plt
     for i in range(len(vcfg['writeSections'])):
         plt.figure(i)
         for j in range(len(sfxVec)):
-            plt.plot(cps[j][i][:,3], cps[j][i][:,4], label=sfxVec[j].replace('_', ' '))
+            plt.plot(cps_i[j][i][:,3], cps_i[j][i][:,4], label=sfxVec[j].replace('_', ' ')+'_invicid')
+        for j in range(len(sfxVec)):
+            plt.plot(cps_v[j][i][:,3], cps_v[j][i][:,4], label=sfxVec[j].replace('_', ' ')+'_viscous')
         plt.legend()
+        plt.gca().invert_yaxis()
         plt.title('Section ' + str(i))
         plt.xlabel('x')
         plt.ylabel('cp')
     plt.show()
-    # Test solution
-    print(ccolors.ANSI_BLUE + 'PyTesting...' + ccolors.ANSI_RESET)
-    tests = CTests()
-    tests.add(CTest('Cl', iSolverAPI.getCl(), 0.069, 5e-2))
-    tests.add(CTest('Cd', vSolution['Cdt_int'], 0.00498, 1e-3, forceabs=True))
-    tests.add(CTest('Iterations', len(aeroCoeffs), 8, 0, forceabs=True))
-    tests.run()
+
+    # Test solution only for AoA = 0.
+    if len(AoAVec) == 1 and AoAVec[0] == 0.0:
+        print(ccolors.ANSI_BLUE + 'PyTesting...' + ccolors.ANSI_RESET)
+        tests = CTests()
+        tests.add(CTest('Cl', isol.getCl(), 0.0, 1e-3))
+        tests.add(CTest('Cd', isol.getCd()+vsol.Cdf, 0.00564, 1e-3, forceabs=True))
+        tests.add(CTest('Iterations', len(aeroCoeffs['Cl']), 5, 0, forceabs=True))
+        tests.run()
 
     # eof
     print('')
-- 
GitLab