Merge pull request #25 from pkck28/main

Updated documentation and few fixes

blackbox/base/airfoilbaseclass.py

@@ -6,7 +6,6 @@
 from pyDOE2 import lhs
 from baseclasses import AeroProblem
 from scipy.io import savemat
-from scipy import integrate
 from mpi4py import MPI
 
 comm = MPI.COMM_WORLD
@@ -272,11 +271,14 @@ def getObjectives(self, x: np.ndarray) -> tuple:
         pkgdir = sys.modules["blackbox"].__path__[0]
 
         # Setting filepath based on the how alpha is treated alpha
-        if self.options["alpha"] == "explicit":
-            filepath = os.path.join(pkgdir, "runscripts/airfoil/runscript_airfoil.py")
+        if type(self).__name__ == "AirfoilCSTMultipoint":
+            filepath = os.path.join(pkgdir, "runscripts/airfoil/runscript_airfoil_cst_mp.py")
         else:
-            # filepath = os.path.join(pkgdir, "runscripts/airfoil/runscipt_airfoil_cst_opt.py")
-            filepath = os.path.join(pkgdir, "runscripts/airfoil/runscript_airfoil_rf.py")
+            if self.options["alpha"] == "explicit":
+                filepath = os.path.join(pkgdir, "runscripts/airfoil/runscript_airfoil.py")
+            else:
+                # filepath = os.path.join(pkgdir, "runscripts/airfoil/runscipt_airfoil_cst_opt.py")
+                filepath = os.path.join(pkgdir, "runscripts/airfoil/runscript_airfoil_rf.py")
 
         # Copy the runscript to analysis directory
         shutil.copy(filepath, "{}/{}/runscript.py".format(directory, self.genSamples+1))
@@ -288,7 +290,7 @@ def getObjectives(self, x: np.ndarray) -> tuple:
             self._writeCoords(coords=points, filename="deformedAirfoil.dat")
 
         if self.options["plotAirfoil"]:
-            self._plotAirfoil(self.plt, self.coords, points)
+            self._plotAirfoil(self.plt, self.origCoords, points)
 
         if self.parametrization == "FFD" and self.options["writeDeformedFFD"]:
             self.DVGeo.writePlot3d("deformedFFD.xyz")
@@ -332,7 +334,16 @@ def getObjectives(self, x: np.ndarray) -> tuple:
             filehandler.close()
 
             # Calculate the area
-            output["area"] = integrate.simpson(points[:,0], points[:,1], even="avg")
+            x = points[:,0]
+            y = points[:,1]
+
+            area = 0.0
+            N = len(x)
+            j = N - 1
+            for i in range(0,N):
+                area += (x[j] + x[i]) * (y[j] - y[i])
+                j = i
+            output["area"] = abs(area)/2.0
 
             if self.options["getFlowFieldData"]:
                 # Reading the cgns file
@@ -405,7 +416,7 @@ def getAirfoil(self, x: np.ndarray) -> np.ndarray:
 
         # If no geometric design variable is present, then return the original airfoil
         if self.DVGeo.getNDV() == 0:
-            return self.coords[:,0:2]
+            return self.origCoords[:,0:2]
 
         # Creating dictionary from x
         newDV = {}
@@ -414,6 +425,10 @@ def getAirfoil(self, x: np.ndarray) -> np.ndarray:
             loc = loc.reshape(-1,)
             newDV[dv] = x[loc]
 
+        if type(self).__name__ == "AirfoilCSTMultipoint":
+            for ap in self.options["aeroProblem"]:
+                ap.setDesignVars(newDV)
+
         if self.parametrization == "FFD" and self.options["fixLETE"]:
 
             # Adjusting LE FFD points
@@ -459,9 +474,13 @@ def calculateArea(self, x: np.ndarray) -> float:
         x = points[:,0]
         y = points[:,1]
 
-        # Calculate the area using simpson's rule
-        # Note: x and y are both flipped here
-        area = integrate.simpson(x, y, even='avg')
+        area = 0.0
+        N = len(x)
+        j = N - 1
+        for i in range(0,N):
+            area += (x[j] + x[i]) * (y[j] - y[i])
+            j = i
+        area = abs(area)/2.0
 
         return area
 
@@ -523,9 +542,18 @@ def _checkOptions(self, defaultOptions: list, requiredOptions: list, options: di
                     self._error("Second entry in \"numCST\" is less than 1.")
 
         ############ Validating aeroProblem
-        if "aeroProblem" in userProvidedOptions:
-            if not isinstance(options["aeroProblem"], AeroProblem):
-                self._error("\"aeroProblem\" attribute is not an aeroproblem.")
+        if type(self).__name__ == "AirfoilCSTMultipoint":
+            if "aeroProblem" in userProvidedOptions:
+                if not isinstance(options["aeroProblem"], list):
+                    self._error("\"aeroProblem\" attribute is not a list.")
+
+                for ap in options["aeroProblem"]:
+                    if not isinstance(ap, AeroProblem):
+                        self._error("\"aeroProblem\" attribute doesn't contain AeroProblem object.")
+        else:
+            if "aeroProblem" in userProvidedOptions:
+                if not isinstance(options["aeroProblem"], AeroProblem):
+                    self._error("\"aeroProblem\" attribute is not an aeroproblem.")
 
         ############ Validating solverOptions
         if "solverOptions" in userProvidedOptions:

blackbox/physics_models/airfoil/airfoil_cst.py

@@ -55,6 +55,12 @@ def __init__(self, options):
         # Updating/Appending the default option list with user provided options
         self._setOptions(options)
 
+        # Overiding/set some solver options
+        self.options["solverOptions"]["printAllOptions"] = False
+        self.options["solverOptions"]["printIntro"] = False
+        self.options["solverOptions"]["outputDirectory"] = "."
+        self.options["solverOptions"]["numberSolutions"] = False
+
         # Raise an error if pyvista is not installed
         if self.options["getFlowFieldData"]:
             if msg_pyvista != None:
@@ -83,26 +89,26 @@ def __init__(self, options):
             os.system("mkdir {}".format(directory))
 
         # Read the coordinate file
-        self.coords = readCoordFile(self.options["airfoilFile"])
+        self.origCoords = readCoordFile(self.options["airfoilFile"])
 
         # Some validation for coordinate file
-        if self.coords[0,0] != self.coords[-1,0]:
+        if self.origCoords[0,0] != self.origCoords[-1,0]:
             self._error("The X coordinate of airfoil doesn't start and end at same point.")
-        elif self.coords[0,1] != self.coords[-1,1]:
+        elif self.origCoords[0,1] != self.origCoords[-1,1]:
             self._error("The Y coordinate of airfoil doesn't start and end at same point.")
 
         # Initializing the parametrization object
         self.DVGeo = DVGeometryCST(self.options["airfoilFile"], numCST=self.options["numCST"], comm=comm)
         self.parametrization = "CST"
 
         # Adding pointset to the parametrization
-        self.coords = np.hstack(( self.coords, np.zeros((self.coords.shape[0], 1)) ))
-        self.DVGeo.addPointSet(self.coords, "airfoil")
+        self.origCoords = np.hstack(( self.origCoords, np.zeros((self.origCoords.shape[0], 1)) ))
+        self.DVGeo.addPointSet(self.origCoords, "airfoil")
 
         # Checking the number of points at trailing edge for blunt TE
         # Only two are allowed for CST. Otherwise, meshing will have problem.
         if not self.DVGeo.sharp:
-            if len(np.where(self.coords[1:-1,0] == self.coords[0,0])[0]) > 1:
+            if len(np.where(self.origCoords[1:-1,0] == self.origCoords[0,0])[0]) > 1:
                 self._error("There are more than two points in the trailing edge.")
 
         # Some initializations which will be used later