diff --git a/baseclasses/__init__.py b/baseclasses/__init__.py index 9ecf32e..9dbaee9 100644 --- a/baseclasses/__init__.py +++ b/baseclasses/__init__.py @@ -1,24 +1,21 @@ -__version__ = "1.8.0" +__version__ = "1.9.0" from .problems import ( AeroProblem, - TransiProblem, - StructProblem, AeroStructProblem, + EngineProblem, + FieldPerformanceProblem, + FluidProperties, + FuelCase, + ICAOAtmosphere, + LGProblem, MissionProblem, MissionProfile, MissionSegment, + StructProblem, + TransiProblem, WeightProblem, - FuelCase, - FluidProperties, - ICAOAtmosphere, - EngineProblem, - FieldPerformanceProblem, - LGProblem, ) - -from .solvers import BaseSolver, AeroSolver - -from .utils import getPy3SafeString - +from .solvers import AeroSolver, BaseSolver from .testing import BaseRegTest, getTol +from .utils import getPy3SafeString, tecplotIO diff --git a/baseclasses/problems/pyWeight_problem.py b/baseclasses/problems/pyWeight_problem.py index 4c84392..5b931c7 100644 --- a/baseclasses/problems/pyWeight_problem.py +++ b/baseclasses/problems/pyWeight_problem.py @@ -4,8 +4,13 @@ Holds the weightProblem class for weightandbalance solvers. """ -import numpy as np import copy +from pathlib import Path + +import numpy as np + +from ..utils import TecplotFEZone, TecplotOrderedZone, writeTecplot +from ..utils.tecplotIO import ZoneType try: from pygeo import geo_utils @@ -68,9 +73,9 @@ def addComponents(self, components): # *components? """ # Check if components is of type Component or list, otherwise raise Error - if type(components) == list: + if isinstance(components, list): pass - elif type(components) == object: + elif isinstance(components, object): components = [components] else: raise Error("addComponents() takes in either a list of or a single component") @@ -132,7 +137,7 @@ def setSurface(self, surf): """ - if type(surf) == list: + if isinstance(surf, list): self.p0 = np.array(surf[0]) self.v1 = np.array(surf[1]) self.v2 = np.array(surf[2]) @@ -194,24 +199,28 @@ def writeSurfaceTecplot(self, fileName): File name for tecplot file. Should have a .dat extension. """ - f = open(fileName, "w") - f.write('TITLE = "weight_problem Surface Mesh"\n') - f.write('VARIABLES = "CoordinateX" "CoordinateY" "CoordinateZ"\n') - f.write("Zone T=%s\n" % ("surf")) - f.write("Nodes = %d, Elements = %d ZONETYPE=FETRIANGLE\n" % (len(self.p0) * 3, len(self.p0))) - f.write("DATAPACKING=POINT\n") - for i in range(len(self.p0)): - points = [] - points.append(self.p0[i]) - points.append(self.p0[i] + self.v1[i]) - points.append(self.p0[i] + self.v2[i]) - for i in range(len(points)): - f.write(f"{points[i][0]:f} {points[i][1]:f} {points[i][2]:f}\n") + # Build the FETriangle data array + dataArrays = np.zeros((len(self.p0) * 3, 3), dtype=float) + dataArrays[::3] = self.p0 + dataArrays[1::3] = self.p0 + self.v1 + dataArrays[2::3] = self.p0 + self.v2 + data = {"CoordinateX": dataArrays[:, 0], "CoordinateY": dataArrays[:, 1], "CoordinateZ": dataArrays[:, 2]} + # Create the connectivity + conn = np.zeros((len(self.p0), 3), dtype=int) for i in range(len(self.p0)): - f.write("%d %d %d\n" % (3 * i + 1, 3 * i + 2, 3 * i + 3)) + conn[i, :] = [3 * i + 1, 3 * i + 2, 3 * i + 3] - f.close() + # Create the single zone + zones = [TecplotFEZone("surf", data, conn, zoneType=ZoneType.FETRIANGLE)] + + writeTecplot( + fileName, + title="weight_problem Surface Mesh", + zones=zones, + datapacking="POINT", + precision="SINGLE", + ) def writeTecplot(self, fileName): """ @@ -362,9 +371,9 @@ def addFuelCases(self, cases): """ # Check if case is a single entry or a list, otherwise raise Error - if type(cases) == list: + if isinstance(cases, list): pass - elif type(cases) == object: + elif isinstance(cases, object): cases = [cases] else: raise Error("addFuelCases() takes in either a list of or a single fuelcase") @@ -447,7 +456,7 @@ def _getComponentKeys(self, include=None, exclude=None, includeType=None, exclud if includeType is not None: # Specified a list of component types to include - if type(includeType) == str: + if isinstance(includeType, str): includeType = [includeType] weightKeysTmp = set() for key in weightKeys: @@ -457,21 +466,21 @@ def _getComponentKeys(self, include=None, exclude=None, includeType=None, exclud if include is not None: # Specified a list of compoents to include - if type(include) == str: + if isinstance(include, str): include = [include] include = set(include) weightKeys.intersection_update(include) if exclude is not None: # Specified a list of components to exclude - if type(exclude) == str: + if isinstance(exclude, str): exclude = [exclude] exclude = set(exclude) weightKeys.difference_update(exclude) if excludeType is not None: # Specified a list of compoent types to exclude - if type(excludeType) == str: + if isinstance(excludeType, str): excludeType = [excludeType] weightKeysTmp = copy.copy(weightKeys) for key in weightKeys: @@ -490,48 +499,42 @@ def writeMassesTecplot(self, filename): filename: str filename for writing the masses. This string will have the - .dat suffix appended to it. + # .dat suffix appended to it if it does not already have it. """ - fileHandle = filename + ".dat" - f = open(fileHandle, "w") nMasses = len(self.nameList) - f.write('TITLE = "%s: Mass Data"\n' % self.name) - f.write('VARIABLES = "X", "Y", "Z", "Mass"\n') locList = ["current", "fwd", "aft"] + zones = [] for loc in locList: - f.write('ZONE T="%s", I=%d, J=1, K=1, DATAPACKING=POINT\n' % (loc, nMasses)) - - for key in self.components.keys(): + dataArray = np.zeros((nMasses, 4), dtype=float) + for i, key in enumerate(self.components.keys()): CG = self.components[key].getCG(loc) mass = self.components[key].getMass() - x = np.real(CG[0]) - y = np.real(CG[1]) - z = np.real(CG[2]) - m = np.real(mass) - - f.write(f"{x:f} {y:f} {z:f} {m:f}\n") - - # end - f.write("\n") - # end - - # textOffset = 0.5 - # for loc in locList: - # for name in self.nameList: - # x= np.real(self.componentDict[name].CG[loc][0]) - # y= np.real(self.componentDict[name].CG[loc][1]) - # z= np.real(self.componentDict[name].CG[loc][2])+textOffset - # m= np.real(self.componentDict[name].W) - - # f.write('TEXT CS=GRID3D, HU=POINT, X=%f, Y=%f, Z=%f, H=12, T="%s"\n'%(x,y,z,name+' '+loc)) - # # end - - # # end - - f.close() - return + dataArray[i, 0] = CG[0] + dataArray[i, 1] = CG[1] + dataArray[i, 2] = CG[2] + dataArray[i, 3] = mass + + data = { + "CoordinateX": dataArray[:, 0], + "CoordinateY": dataArray[:, 1], + "CoordinateZ": dataArray[:, 2], + "Mass": dataArray[:, 3], + } + + zones.append(TecplotOrderedZone(loc, data)) + + # Create the path with the .dat extension + filePath = Path(filename).with_suffix(".dat") + + writeTecplot( + filePath, + title=f"{self.name}: Mass Data", + zones=zones, + datapacking="POINT", + precision="SINGLE", + ) def writeProblemData(self, fileName): """ diff --git a/baseclasses/solvers/pyAero_solver.py b/baseclasses/solvers/pyAero_solver.py index 27c3ae6..32cdfee 100644 --- a/baseclasses/solvers/pyAero_solver.py +++ b/baseclasses/solvers/pyAero_solver.py @@ -13,7 +13,8 @@ # Extension modules # ============================================================================= from .BaseSolver import BaseSolver -from ..utils import CaseInsensitiveDict, Error +from ..utils import CaseInsensitiveDict, Error, TecplotFEZone, writeTecplot +from ..utils.tecplotIO import ZoneType # ============================================================================= # AeroSolver Class @@ -21,7 +22,6 @@ class AeroSolver(BaseSolver): - """ Abstract Class for Aerodynamic Solver Object """ @@ -207,24 +207,25 @@ def writeTriangulatedSurfaceTecplot(self, fileName, groupName=None, **kwargs): """ [p0, v1, v2] = self.getTriangulatedMeshSurface(groupName, **kwargs) if self.comm.rank == 0: - f = open(fileName, "w") - f.write('TITLE = "%s Surface Mesh"\n' % self.name) - f.write('VARIABLES = "CoordinateX" "CoordinateY" "CoordinateZ"\n') - f.write("Zone T=%s\n" % ("surf")) - f.write("Nodes = %d, Elements = %d ZONETYPE=FETRIANGLE\n" % (len(p0) * 3, len(p0))) - f.write("DATAPACKING=POINT\n") - for i in range(len(p0)): - points = [] - points.append(p0[i]) - points.append(p0[i] + v1[i]) - points.append(p0[i] + v2[i]) - for i in range(len(points)): - f.write(f"{points[i][0]:f} {points[i][1]:f} {points[i][2]:f}\n") + dataArray = np.zeros((len(p0) * 3, 3), dtype=float) + dataArray[::3] = p0 + dataArray[1::3] = p0 + v1 + dataArray[2::3] = p0 + v2 + data = {"CoordinateX": dataArray[:, 0], "CoordinateY": dataArray[:, 1], "CoordinateZ": dataArray[:, 2]} + conn = np.zeros((len(p0), 3), dtype=int) for i in range(len(p0)): - f.write("%d %d %d\n" % (3 * i + 1, 3 * i + 2, 3 * i + 3)) + conn[i, :] = [3 * i + 1, 3 * i + 2, 3 * i + 3] - f.close() + zones = [TecplotFEZone("surf", data, conn, zoneType=ZoneType.FETRIANGLE)] + + writeTecplot( + fileName, + title=f"{self.name} Surface Mesh", + zones=zones, + datapacking="POINT", + precision="SINGLE", + ) def checkSolutionFailure(self, aeroProblem, funcs): """Take in a an aeroProblem and check for failure. Then append the diff --git a/baseclasses/utils/__init__.py b/baseclasses/utils/__init__.py index a79b2be..5e0737f 100644 --- a/baseclasses/utils/__init__.py +++ b/baseclasses/utils/__init__.py @@ -1,8 +1,9 @@ -from .containers import CaseInsensitiveSet, CaseInsensitiveDict +from .containers import CaseInsensitiveDict, CaseInsensitiveSet from .error import Error -from .utils import getPy3SafeString, pp, ParseStringFormat -from .fileIO import writeJSON, readJSON, writePickle, readPickle, redirectIO, redirectingIO +from .fileIO import readJSON, readPickle, redirectingIO, redirectIO, writeJSON, writePickle from .solverHistory import SolverHistory +from .tecplotIO import TecplotFEZone, TecplotOrderedZone, TecplotZone, readTecplot, writeTecplot +from .utils import ParseStringFormat, getPy3SafeString, pp __all__ = [ "CaseInsensitiveSet", @@ -17,5 +18,10 @@ "redirectIO", "redirectingIO", "SolverHistory", + "TecplotZone", + "TecplotFEZone", + "TecplotOrderedZone", + "writeTecplot", + "readTecplot", "ParseStringFormat", ] diff --git a/baseclasses/utils/tecplotIO.py b/baseclasses/utils/tecplotIO.py new file mode 100644 index 0000000..3f009d7 --- /dev/null +++ b/baseclasses/utils/tecplotIO.py @@ -0,0 +1,1840 @@ +import re +import struct +from abc import ABC, abstractmethod +from enum import Enum +from pathlib import Path +from typing import Any, Dict, Generic, List, Literal, TextIO, Tuple, TypeVar, Union + +import numpy as np +import numpy.typing as npt + +T = TypeVar("T", bound="TecplotZone") + + +# ============================================================================== +# ENUMS +# ============================================================================== +class ZoneType(Enum): + """Tecplot finite element zone types""" + + UNSET = -1 + ORDERED = 0 + FELINESEG = 1 + FETRIANGLE = 2 + FEQUADRILATERAL = 3 + FETETRAHEDRON = 4 + FEBRICK = 5 + + +class DataPacking(Enum): + """Tecplot data packing formats""" + + BLOCK = 0 + POINT = 1 + + +class VariableLocation(Enum): + """Grid location of the variable data""" + + NODE = 0 + CELL_CENTER = 1 + NODE_AND_CELL_CENTER = 2 + + +class DataPrecision(Enum): + """Tecplot data precision""" + + SINGLE = 6 + DOUBLE = 12 + + +class BinaryDataPrecisionCodes(Enum): + """Binary data precision codes""" + + SINGLE = 1 + DOUBLE = 2 + + +class DTypePrecision(Enum): + """Numpy data types for single and double precision""" + + SINGLE = np.float32 + DOUBLE = np.float64 + + +class FileType(Enum): + """Tecplot file types""" + + FULL = 0 + GRID = 1 + SOLUTION = 2 + + +class SectionMarkers(Enum): + """Tecplot section markers""" + + ZONE = 299.0 # V11.2 marker + DATA = 357.0 + + +class BinaryFlags(Enum): + """Binary boolean flags""" + + NONE = -1 + FALSE = 0 + TRUE = 1 + + +class StrandID(Enum): + """Strand ID default codes""" + + PENDING = -2 + STATIC = -1 + + +class SolutionTime(Enum): + """Solution time default codes""" + + UNSET = -1 + + +class Separator(Enum): + """Separator characters""" + + SPACE = " " + COMMA = "," + TAB = "\t" + NEWLINE = "\n" + CARRIAGE_RETURN = "\r" + + +# ============================================================================== +# DATA STRUCTURES +# ============================================================================== +class TecplotZone: + """Base class for Tecplot zones.""" + + def __init__( + self, + name: str, + data: Dict[str, npt.NDArray], + solutionTime: float = 0.0, + strandID: int = -1, + ): + """Create a tecplot zone object. + + Parameters + ---------- + name : str + The name of the zone. + data : Dict[str, npt.NDArray] + A dictionary of variable names and their corresponding data. + solutionTime : float, optional + The solution time of the zone, by default 0.0 + strandID : int, optional + The strand id of the zone, by default -1 + """ + self.name = name + self.data = data + self.solutionTime = solutionTime + self.strandID = strandID + self.zoneType: Union[str, ZoneType] = ZoneType.UNSET + self._validateName() + self._validateData() + self._validateSolutionTime() + self._validateStrandID() + + @property + def variables(self) -> List[str]: + return list(self.data.keys()) + + @property + def shape(self) -> Tuple[int, ...]: + return self.data[self.variables[0]].shape + + @property + def nNodes(self) -> int: + return np.multiply.reduce(self.shape) + + def _validateName(self) -> None: + """Check that the zone name is a valid string. + + Raises + ------ + TypeError + If the zone name is not a valid string. + """ + if not isinstance(self.name, str): + raise TypeError("Zone name must be a string.") + + def _validateData(self) -> None: + """Check that the data is a valid dictionary and the values + are numpy arrays that have the same shape. + + Raises + ------ + TypeError + If the data is not a dictionary or the values are not numpy + arrays. + ValueError + If the variables do not have the same shape. + """ + if not isinstance(self.data, dict): + raise TypeError("Data must be a dictionary.") + + for val in self.data.values(): + if not isinstance(val, np.ndarray): + raise TypeError("Data values must be numpy arrays.") + + if not all(self.data[var].shape == self.shape for var in self.variables): + raise ValueError("All variables must have the same shape.") + + def _validateSolutionTime(self) -> None: + """Check that the solution time is a valid float. + + Raises + ------ + TypeError + If the solution time is not a float. + ValueError + If the solution time is less than zero. + """ + if not isinstance(self.solutionTime, float): + raise TypeError("Solution time must be a float.") + + if self.solutionTime < 0.0: + raise ValueError("Solution time must be greater than or equal to zero.") + + def _validateStrandID(self) -> None: + """Check that the strand ID is a valid integer. + + Raises + ------ + TypeError + If the strand ID is not an integer. + """ + if not isinstance(self.strandID, int): + raise TypeError("Strand ID must be an integer.") + + +class TecplotOrderedZone(TecplotZone): + """Tecplot ordered zone. These zones do not contain connectivity + information because the data is ordered in an (i, j, k) grid. + """ + + def __init__( + self, + name: str, + data: Dict[str, npt.NDArray], + solutionTime: float = 0.0, + strandID: int = -1, + ): + """To create a tecplot ordered zone object: + + .. code-block:: python + + # --- Example usage --- + # Create the data + nx, ny = 10, 10 + X = np.random.rand(nx, ny) + Y = np.random.rand(nx, ny) + + # Create the zone + zone = TecplotOrderedZone( + "OrderedZone", + {"X": X, "Y": Y}, + solutionTime=0.0, + strandID=-1, + ) + + Parameters + ---------- + name : str + The name of the zone. + data : Dict[str, npt.NDArray] + A dictionary of variable names and their corresponding data. + zoneType : Union[str, ZoneType], optional + The type of the zone, by default ZoneType.ORDERED + solutionTime : float, optional + The solution time of the zone, by default 0.0 + strandID : int, optional + The strand id of the zone, by default -1 + """ + super().__init__(name, data, solutionTime=solutionTime, strandID=strandID) + self.zoneType = ZoneType.ORDERED + + @property + def iMax(self) -> int: + return self.shape[0] + + @property + def jMax(self) -> int: + return self.shape[1] if len(self.shape) > 1 else 1 + + @property + def kMax(self) -> int: + return self.shape[2] if len(self.shape) > 2 else 1 + + +class TecplotFEZone(TecplotZone): + """Tecplot finite element zone. These zones contain connectivity + information to describe the elements in the zone. The type of + element is determined by the shape of the connectivity array and + the ``tetrahedral`` flag. The connectivity array is 0-based. + + The following shapes correspond to the following element types, + where n is the number of elements: + + - ``(n, 2)``: FELINESEG + - ``(n, 3)``: FETRIANGLE + - ``(n, 4)``, FEQUADRILATERAL + - ``(n, 4)``, FETETRAHEDRON + - ``(n, 8)``: FEBRICK + """ + + def __init__( + self, + zoneName: str, + data: Dict[str, npt.NDArray], + connectivity: npt.NDArray, + zoneType: Union[str, ZoneType], + solutionTime: float = 0.0, + strandID: int = -1, + ): + """To create a tecplot finite element zone object: + + .. code-block:: python + + # --- Example usage --- + # Create node coordinates + x = np.linspace(0, 1, nx + 1) + nodes = np.column_stack((x, x**2)) + + # Create element connectivity + connectivity = np.column_stack((np.arange(nx), np.arange(1, nx + 1))) + + # Create the zone + zone = TecplotFEZone( + "FEZone", + {"x": nodes[:, 0], "y": nodes[:, 1]}, + connectivity, + zoneType="FELINESEG", + solutionTime=0.0, + strandID=-1, + ) + + Parameters + ---------- + zoneName : str + The name of the zone. + data : Dict[str, npt.NDArray] + A dictionary of variable names and their corresponding data. + connectivity : npt.NDArray + The connectivity array that describes the elements in the + zone. + zoneType : Union[str, ZoneType] + The type of the zone. Can be a string that matches an entry + in the ZoneType enum or the ZoneType enum itself. + solutionTime : float, optional + The solution time of the zone, by default 0.0 + strandID : int, optional + The strand id of the zone, by default -1 + """ + super().__init__(zoneName, data, solutionTime=solutionTime, strandID=strandID) + self._connectivity = connectivity + self.zoneType = zoneType + self._validateZoneType() + self._validateConnectivity() + self._uniqueIndices = np.unique(self._connectivity.flatten()) + self._uniqueConnectivity = self._remapConnectivity() + + @property + def connectivity(self) -> npt.NDArray: + return self._connectivity + + @connectivity.setter + def connectivity(self, value: npt.NDArray) -> None: + self._connectivity = value + self._validateConnectivity() + self._uniqueIndices = np.unique(self._connectivity.flatten()) + self._uniqueConnectivity = self._remapConnectivity() + + @property + def nElements(self) -> int: + return self.connectivity.shape[0] + + @property + def triConnectivity(self) -> npt.NDArray: + if self.zoneType == ZoneType.FETRIANGLE: + return self.connectivity + elif self.zoneType == ZoneType.FEQUADRILATERAL: + return np.row_stack((self.connectivity[:, [0, 1, 2]], self.connectivity[:, [0, 2, 3]])) + else: + raise TypeError(f"'triConnectivity' not supported for {self.zoneType.name} zone type.") + + @property + def uniqueData(self) -> Dict[str, npt.NDArray]: + return {var: self.data[var][self._uniqueIndices] for var in self.variables} + + @property + def uniqueConnectivity(self) -> npt.NDArray: + return self._uniqueConnectivity + + def _validateZoneType(self) -> None: + supportedZones = [zone.name for zone in ZoneType if zone.name != "ORDERED"] + if isinstance(self.zoneType, str): + if self.zoneType.upper() not in supportedZones: + raise ValueError("Invalid zone type.") + self.zoneType = ZoneType[self.zoneType.upper()] + elif isinstance(self.zoneType, ZoneType): + if self.zoneType.name not in supportedZones: + raise ValueError("Invalid zone type.") + else: + raise ValueError("Invalid zone type.") + + def _validateConnectivity(self) -> None: + if self.zoneType == ZoneType.FELINESEG: + assert self.connectivity.shape[1] == 2, "Connectivity shape does not match zone type." + elif self.zoneType == ZoneType.FETRIANGLE: + assert self.connectivity.shape[1] == 3, "Connectivity shape does not match zone type." + elif self.zoneType == ZoneType.FEQUADRILATERAL: + assert self.connectivity.shape[1] == 4, "Connectivity shape does not match zone type." + elif self.zoneType == ZoneType.FETETRAHEDRON: + assert self.connectivity.shape[1] == 4, "Connectivity shape does not match zone type." + elif self.zoneType == ZoneType.FEBRICK: + assert self.connectivity.shape[1] == 8, "Connectivity shape does not match zone type." + else: + # Prior validation step should ensure we don't reach this point + # but raise an error just in case. + raise TypeError("Invalid zone type.") + + def _remapConnectivity(self) -> npt.NDArray: + # Create a remapping array that's as large as the maximum index + remap = np.full(self._uniqueIndices.max() + 1, -1, dtype=np.int64) + remap[self._uniqueIndices] = np.arange(len(self._uniqueIndices)) + return remap[self._connectivity] + + +# ============================================================================== +# ASCII WRITERS +# ============================================================================== +def writeArrayToFile( + arr: npt.NDArray[np.float64], + handle: TextIO, + maxLineWidth: int = 4000, + precision: int = 6, + separator: Separator = Separator.SPACE, +) -> None: + """ + Write a 2D numpy array to a file using numpy.array_str with custom + formatting. + + Parameters + ---------- + arr : npt.NDArray[np.float64] + 2D numpy array to write. + file : TextIO + A file-like object (e.g., opened with 'open()') to write to. + maxLineWidth : int, optional + Maximum width of each line in characters, by default 4000. + precision : int, optional + Number of decimal places for floating-point numbers, by default + 6. + separator : Literal[" ", ",", "\\t", "\\n", "\\r"], optional + Separator to use between elements, by default Separator.SPACE + + Raises + ------ + ValueError + If the input array is not 2-dimensional. + """ + if arr.ndim != 2: + raise ValueError("Input must be a 2D numpy array") + + for row in arr: + rowStr = np.array2string( + row, + max_line_width=maxLineWidth, + separator=separator.value, + formatter={"float_kind": lambda x: f"{x:.{precision}E}"}, + threshold=np.inf, + ) + cleanedRowStr = rowStr.strip("[]").strip().replace("\n ", "\n") # Remove brackets and extra spaces + handle.write(cleanedRowStr + "\n") + + +class TecplotZoneWriterASCII(Generic[T], ABC): + def __init__( + self, + zone: T, + datapacking: Literal["BLOCK", "POINT"], + precision: Literal["SINGLE", "DOUBLE"], + separator: Separator = Separator.SPACE, + ) -> None: + """Abstract base class for writing Tecplot zones to ASCII files. + + Parameters + ---------- + zone : T + The Tecplot zone to write. + datapacking : Literal["BLOCK", "POINT"] + The data packing format. BLOCK is row-major, POINT is + column-major. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + """ + self.zone = zone + self.datapacking = DataPacking[datapacking].name + self.fmtPrecision = DataPrecision[precision].value + self.separator = separator + + @abstractmethod + def writeHeader(self, handle: TextIO): + pass + + @abstractmethod + def writeFooter(self, handle: TextIO): + pass + + def writeData(self, handle: TextIO): + data = np.stack([self.zone.data[var] for var in self.zone.variables], axis=-1) + + if self.datapacking == "POINT": + data = data.reshape(-1, len(self.zone.variables)) + else: + data = data.reshape(-1, len(self.zone.variables)).T + + writeArrayToFile(data, handle, maxLineWidth=4000, precision=self.fmtPrecision, separator=self.separator) + + +class TecplotOrderedZoneWriterASCII(TecplotZoneWriterASCII[TecplotOrderedZone]): + def __init__( + self, + zone: TecplotOrderedZone, + datapacking: Literal["BLOCK", "POINT"], + precision: Literal["SINGLE", "DOUBLE"], + separator: Separator = Separator.SPACE, + ) -> None: + """Writer for Tecplot ordered zones in ASCII format. + + Parameters + ---------- + zone : TecplotOrderedZone + The ordered zone to write. + datapacking : Literal["BLOCK", "POINT"] + The data packing format. BLOCK is row-major, POINT is + column-major. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + separator : Separator, optional + Separator to use between elements. The Separator + is an enum defined in :meth:`Separator `, + by default Separator.SPACE + """ + super().__init__(zone, datapacking, precision, separator) + + def writeHeader(self, handle: TextIO): + """Write the zone header to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + # Write the zone header + zoneString = f'ZONE T="{self.zone.name}"' + zoneString += f", I={self.zone.iMax}" + + if self.zone.jMax > 1: + zoneString += f", J={self.zone.jMax}" + + if self.zone.kMax > 1: + zoneString += f", K={self.zone.kMax}" + + # Write the strand ID and solution time + if self.zone.strandID != StrandID.STATIC.value and self.zone.strandID != StrandID.PENDING.value: + # ASCII format does not support the -1 or -2 strand IDs + # So we only write the strand ID if it is not -1 + zoneString += f", STRANDID={self.zone.strandID}" + + # Only write the solution time if it is set + if self.zone.solutionTime != SolutionTime.UNSET.value: + zoneString += f", SOLUTIONTIME={self.zone.solutionTime}" + + zoneString += f", DATAPACKING={self.datapacking}\n" + + handle.write(zoneString) + + def writeFooter(self, handle: TextIO): + """Write the zone footer to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + handle.write("\n") + + +class TecplotFEZoneWriterASCII(TecplotZoneWriterASCII[TecplotFEZone]): + def __init__( + self, + zone: TecplotFEZone, + datapacking: Literal["BLOCK", "POINT"], + precision: Literal["SINGLE", "DOUBLE"], + separator: Separator = Separator.SPACE, + ) -> None: + """Writer for Tecplot finite element zones in ASCII format. + + Parameters + ---------- + zone : TecplotFEZone + The finite element zone to write. + datapacking : Literal["BLOCK", "POINT"] + The data packing format. BLOCK is row-major, POINT is + column-major. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + separator : Separator, optional + Separator to use between elements. The Separator is an enum + defined in :meth:`Separator `, + by default Separator.SPACE + """ + super().__init__(zone, datapacking, precision, separator) + + def writeHeader(self, handle: TextIO): + """Write the zone header to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + # Write the zone header + zoneString = f'ZONE T="{self.zone.name}"' + zoneString += f", DATAPACKING={self.datapacking}" + + # Write the node and element information + zoneString += f", NODES={np.multiply.reduce(self.zone.shape):d}" + zoneString += f", ELEMENTS={self.zone.nElements:d}" + zoneString += f", ZONETYPE={self.zone.zoneType.name}" + + # Write the strand ID and solution time + if self.zone.strandID != StrandID.STATIC.value and self.zone.strandID != StrandID.PENDING.value: + # ASCII format does not support the -1 or -2 strand IDs + # So we only write the strand ID if it is not -1 + zoneString += f", STRANDID={self.zone.strandID}" + + # Only write the solution time if it is set + if self.zone.solutionTime != SolutionTime.UNSET.value: + zoneString += f", SOLUTIONTIME={self.zone.solutionTime}\n" + + handle.write(zoneString) + + def writeFooter(self, handle: TextIO): + """Write the zone footer to the file. This includes the + connectivity information. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + connectivity = self.zone.connectivity + 1 + # Get the max characters in the connectivity + maxChars = len(str(connectivity.max())) + + np.savetxt(handle, connectivity, fmt=f"%{maxChars}d") + + handle.write("\n") + + +class TecplotWriterASCII: + def __init__( + self, + title: str, + zones: List[TecplotZone], + datapacking: Literal["BLOCK", "POINT"], + precision: Literal["SINGLE", "DOUBLE"], + separator: Separator = Separator.SPACE, + ) -> None: + """Writer for Tecplot files in ASCII format. + + Parameters + ---------- + title : str + The title of the Tecplot file. + zones : List[TecplotZone] + A list of Tecplot zones to write. + datapacking : Literal["BLOCK", "POINT"] + The data packing format. BLOCK is row-major, POINT is + column-major. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + separator : Separator, optional + Separator to use between elements. The Separator + is an enum defined in :meth:`Separator `, + by default Separator.SPACE + """ + self.title = title + self.zones = zones + self.datapacking = DataPacking[datapacking].name + self.precision = DataPrecision[precision].name + self.separator = separator + self._validateVariables() + + def _validateVariables(self) -> None: + """Check that all zones have the same variables.""" + if not all(set(self.zones[0].variables) == set(zone.variables) for zone in self.zones): + raise ValueError("All zones must have the same variables.") + + def _writeVariables(self, handle: TextIO) -> None: + """Write the variable names to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + variables = [f'"{var}"' for var in self.zones[0].variables] + variableString = ", ".join(variables) + handle.write(f"VARIABLES = {variableString}\n") + + def _writeZone(self, handle: TextIO, zone: TecplotZone) -> None: + """Write a Tecplot zone to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + zone : TecplotZone + The zone to write. + + Raises + ------ + ValueError + If the zone type is invalid. + """ + if isinstance(zone, TecplotOrderedZone): + writer = TecplotOrderedZoneWriterASCII(zone, self.datapacking, self.precision, self.separator) + elif isinstance(zone, TecplotFEZone): + writer = TecplotFEZoneWriterASCII(zone, self.datapacking, self.precision, self.separator) + else: + raise ValueError("Invalid zone type.") + + writer.writeHeader(handle) + writer.writeData(handle) + writer.writeFooter(handle) + + def write(self, filename: Union[str, Path]) -> None: + """Write the Tecplot file to disk. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object. + """ + with open(filename, "w") as handle: + handle.write(f'TITLE = "{self.title}"\n') + self._writeVariables(handle) + for zone in self.zones: + self._writeZone(handle, zone) + + +# ============================================================================== +# BINARY WRITERS +# ============================================================================== +def _writeInteger(handle: TextIO, value: int) -> None: + """Write an integer to a binary file as int32. + + Parameters + ---------- + handle : TextIO + The file handle to write to. + value : int + Integer value to write. + """ + handle.write(struct.pack("i", value)) + + +def _writeFloat32(handle: TextIO, value: float) -> None: + """Write a float to a binary file as float32. + + Parameters + ---------- + handle : TextIO + The file handle to write to. + value : float + Float value to write. + """ + handle.write(struct.pack("f", value)) + + +def _writeFloat64(handle: TextIO, value: float) -> None: + """Write a float to a binary file as float64. + + Parameters + ---------- + handle : TextIO + The file handle to write to. + value : float + Float value to write. + """ + handle.write(struct.pack("d", value)) + + +def _writeString(handle: TextIO, value: str) -> None: + """Write a string to a binary file as a string. + + Parameters + ---------- + handle : TextIO + The file handle to write to. + value : str + String value to write. + """ + for char in value: + asciiValue = ord(char) + handle.write(struct.pack("i", asciiValue)) + + handle.write(struct.pack("i", 0)) + + +class TecplotZoneWriterBinary(Generic[T], ABC): + def __init__( + self, + title: str, + zone: T, + precision: Literal["SINGLE", "DOUBLE"], + ) -> None: + """Abstract base class for writing Tecplot zones to binary + files. + + Parameters + ---------- + title : str + The title of the Tecplot file. + zone : T + The Tecplot zone to write. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + """ + self.title = title + self.zone = zone + self.datapacking = "BLOCK" + self.precision = DataPrecision[precision].name + + def _writeCommonHeader(self, handle: TextIO) -> None: + """Write the common header information for all zones. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + # Write the zone marker + _writeFloat32(handle, SectionMarkers.ZONE.value) # Write the zone marker + _writeString(handle, self.zone.name) # Write the zone name + _writeInteger(handle, BinaryFlags.NONE.value) # Write the parent zone + _writeInteger(handle, self.zone.strandID) # Write the strand ID + _writeFloat64(handle, self.zone.solutionTime) # Write the solution time + _writeInteger(handle, BinaryFlags.NONE.value) # Write the default color + _writeInteger(handle, self.zone.zoneType.value) # Write the zone type + _writeInteger(handle, DataPacking.BLOCK.value) # Data Packing (Always block for binary) + _writeInteger(handle, VariableLocation.NODE.value) # Specify the variable location + _writeInteger(handle, BinaryFlags.FALSE.value) # Are raw 1-1 face neighbors supplied + + @abstractmethod + def writeHeader(self, handle: TextIO): + pass + + def writeData(self, handle: TextIO): + """Write the zone data to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + # Get the data into a single array + data = np.stack([self.zone.data[var] for var in self.zone.variables], axis=-1) + + # Flatten the data such that each variable is a row + data = data.reshape(-1, len(self.zone.variables)).T + + _writeFloat32(handle, SectionMarkers.ZONE.value) # Write the zone marker + + # Write the variable data format for each variable + for _ in range(len(self.zone.variables)): + _writeInteger(handle, BinaryDataPrecisionCodes[self.precision].value) + + _writeInteger(handle, BinaryFlags.FALSE.value) # No passive variables + _writeInteger(handle, BinaryFlags.FALSE.value) # No variable sharing + _writeInteger(handle, BinaryFlags.NONE.value) # No connectivity sharing + + # Write the min/max values for the variables + for i in range(len(self.zone.variables)): + _writeFloat64(handle, data[i, ...].min()) + _writeFloat64(handle, data[i, ...].max()) + + # Write the data using the specified data format (single or double) + data.astype(DTypePrecision[self.precision].value).tofile(handle) + + @abstractmethod + def writeFooter(self, handle: TextIO): + pass + + +class TecplotOrderedZoneWriterBinary(TecplotZoneWriterBinary[TecplotOrderedZone]): + def __init__( + self, + title: str, + zone: TecplotOrderedZone, + precision: Literal["SINGLE", "DOUBLE"], + ) -> None: + """Writer for Tecplot ordered zones in binary format. + + Parameters + ---------- + title : str + The title of the Tecplot file. + zone : TecplotOrderedZone + The ordered zone to write. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + """ + super().__init__(title, zone, precision) + + def writeHeader(self, handle: TextIO): + """Write the zone header to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + self._writeCommonHeader(handle) # Write the common header information + + # --- Specific to Ordered Zones --- + _writeInteger(handle, self.zone.iMax) # Write the I dimension + _writeInteger(handle, self.zone.jMax) # Write the J dimension + _writeInteger(handle, self.zone.kMax) # Write the K dimension + _writeInteger(handle, BinaryFlags.FALSE.value) # No aux data + + def writeFooter(self, handle: TextIO): + """Write the zone footer to the file. This is not used for + ordered zones. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + pass + + +class TecplotFEZoneWriterBinary(TecplotZoneWriterBinary[TecplotFEZone]): + def __init__( + self, + title: str, + zone: TecplotFEZone, + precision: Literal["SINGLE", "DOUBLE"], + ) -> None: + """Writer for Tecplot finite element zones in binary format. + + Parameters + ---------- + title : str + The title of the Tecplot file. + zone : TecplotFEZone + The finite element zone to write. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + """ + super().__init__(title, zone, precision) + + def writeHeader(self, handle: TextIO): + """Write the zone header to the file. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + self._writeCommonHeader(handle) # Write the common header information + + # --- Specific to FE Zones --- + _writeInteger(handle, self.zone.nNodes) # Write the number of nodes + _writeInteger(handle, self.zone.nElements) # Write the number of elements + _writeInteger(handle, 0) # iCellDim (future use, set to 0) + _writeInteger(handle, 0) # jCellDim (future use, set to 0) + _writeInteger(handle, 0) # kCellDim (future use, set to 0) + _writeInteger(handle, BinaryFlags.FALSE.value) # No aux data + + def writeFooter(self, handle: TextIO): + """Write the zone footer to the file. This includes the + connectivity information. + + Parameters + ---------- + handle : TextIO + The file handle. + """ + self.zone.connectivity.astype("int32").tofile(handle) + + +class TecplotWriterBinary: + def __init__( + self, + title: str, + zones: List[TecplotZone], + precision: Literal["SINGLE", "DOUBLE"], + ) -> None: + """Writer for Tecplot files in binary format. + + This writer only supports files formatted using the format + designated by magic number ``#!TDV112``. + + See the Tecplot 360 User's Manual for more information on the + binary file format and the magic number. + + Parameters + ---------- + title : str + The title of the Tecplot file. + zones : List[TecplotZone] + A list of Tecplot zones to write. + precision : Literal["SINGLE", "DOUBLE"] + The floating point precision to write the data. + """ + self._magicNumber = b"#!TDV112" + self.title = title + self.zones = zones + self.precision = precision + self._checkVariables() + + def _checkVariables(self) -> None: + """Check that all zones have the same variables.""" + if not all(set(self.zones[0].variables) == set(zone.variables) for zone in self.zones): + raise ValueError("All zones must have the same variables.") + + def _getZoneWriter(self, zone: TecplotZone) -> TecplotZoneWriterBinary: + """Get the appropriate zone writer based on the zone type. + + Parameters + ---------- + zone : TecplotZone + The Tecplot zone to write. + + Returns + ------- + TecplotZoneWriterBinary + The appropriate zone writer object. + + Raises + ------ + ValueError + If the zone type is invalid. + """ + if isinstance(zone, TecplotOrderedZone): + return TecplotOrderedZoneWriterBinary(self.title, zone, self.precision) + elif isinstance(zone, TecplotFEZone): + return TecplotFEZoneWriterBinary(self.title, zone, self.precision) + else: + raise ValueError("Invalid zone type.") + + def write(self, filename: Union[str, Path]) -> None: + """Write the Tecplot file to disk. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object. + """ + with open(filename, "wb") as handle: + handle.write(self._magicNumber) # Magic number + _writeInteger(handle, 1) # Byte order + _writeInteger(handle, FileType.FULL.value) # Full filetype + _writeString(handle, self.title) # Write the title + _writeInteger(handle, len(self.zones[0].variables)) # Write the number of variables + + for var in self.zones[0].variables: + _writeString(handle, var) # Write the variable names + + # Write the zone headers + for zone in self.zones: + writer = self._getZoneWriter(zone) + writer.writeHeader(handle) + + # Write the data marker + _writeFloat32(handle, SectionMarkers.DATA.value) + + # Write the data and footer for each zone + for zone in self.zones: + writer = self._getZoneWriter(zone) + writer.writeData(handle) + writer.writeFooter(handle) + + +# ============================================================================== +# ASCII READERS +# ============================================================================== +def readArrayData(filename: str, iCurrent: int, nVals: int, nVars: int) -> Tuple[npt.NDArray, int]: + """Read array data from a Tecplot ASCII file. + + Parameters + ---------- + filename : str + The filename of the Tecplot file. + iCurrent : int + The current line number in the file. + nVals : int + The number of values to read. + nVars : int + The number of variables in the data. + + Returns + ------- + Tuple[npt.NDArray, int] + The data and the number of lines read. + """ + with open(filename, "r") as handle: + lines = handle.readlines() + + pattern = r"[\s,\t\n\r]+" # Separator pattern + + data = [] + nLines = 0 + while len(data) < nVals * nVars: + line = lines[iCurrent].strip() # Remove leading/trailing whitespace + vals = [float(val) for val in re.split(pattern, line) if val] # Split the line into values + data.extend(vals) # Add the values to the data list + iCurrent += 1 # Increment the line number + nLines += 1 + + data = np.array(data) + + return data, nLines + + +class TecplotASCIIReader: + def __init__(self, filename: Union[str, Path]) -> None: + """Reader for Tecplot files in ASCII format. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object. + """ + self.filename = filename + + def _readZoneHeader(self, lines: List[str], iCurrent: int) -> Tuple[Dict[str, Any], int]: + """Read the zone header information from a line in a Tecplot file. + + Parameters + ---------- + line : str + The line containing the zone header information. + + Returns + ------- + Dict[str, Any] + A dictionary containing the parsed zone header information. + """ + # Get all the header lines into a single string + header = [] + + # Loop until the line starts with a number which denotes the start of the data section + exitPattern = re.compile(r"^\s*\d") + while not exitPattern.match(lines[iCurrent]): + header.append(lines[iCurrent].strip("\n")) + iCurrent += 1 + + # Join the header lines into a single string + headerString = ", ".join(header) + + # Use regex to parse the header information + zoneNameMatch = re.search(r'(zone t)\s*=\s*[\'""]?([^\'""\n,]+)[\'""]?(?=[,\n]|$)', headerString, re.IGNORECASE) + zoneName = zoneNameMatch.group(2) if zoneNameMatch else None + + zoneTypeMatch = re.search(r"zonetype\s*=\s*(\w+)", headerString, re.IGNORECASE) + zoneType = zoneTypeMatch.group(1) if zoneTypeMatch else "ORDERED" + + datapackingMatch = re.search(r"datapacking\s*=\s*(\w+)", headerString, re.IGNORECASE) + datapacking = datapackingMatch.group(1) if datapackingMatch else None + + nNodesMatch = re.search(r"nodes\s*=\s*(\d+)", headerString, re.IGNORECASE) + nNodes = int(nNodesMatch.group(1)) if nNodesMatch else None + + nElementsMatch = re.search(r"elements\s*=\s*(\d+)", headerString, re.IGNORECASE) + nElements = int(nElementsMatch.group(1)) if nElementsMatch else None + + iMaxMatch = re.search(r"i\s*=\s*(\d+)", headerString, re.IGNORECASE) + iMax = int(iMaxMatch.group(1)) if iMaxMatch else 1 + + jMaxMatch = re.search(r"j\s*=\s*(\d+)", headerString, re.IGNORECASE) + jMax = int(jMaxMatch.group(1)) if jMaxMatch else 1 + + kMaxMatch = re.search(r"k\s*=\s*(\d+)", headerString, re.IGNORECASE) + kMax = int(kMaxMatch.group(1)) if kMaxMatch else 1 + + solutionTimeMatch = re.search(r"solutiontime\s*=\s*(\d+\.\d+)", headerString, re.IGNORECASE) + solutionTime = float(solutionTimeMatch.group(1)) if solutionTimeMatch else 0.0 + + strandIDMatch = re.search(r"strandid\s*=\s*(\d+)", headerString, re.IGNORECASE) + strandID = int(strandIDMatch.group(1)) if strandIDMatch else -1 + + headerDict = { + "zoneName": zoneName, + "zoneType": zoneType, + "datapacking": datapacking, + "nNodes": nNodes, + "nElements": nElements, + "iMax": iMax, + "jMax": jMax, + "kMax": kMax, + "solutionTime": solutionTime, + "strandID": strandID, + } + + return headerDict, iCurrent + + def _readOrderedZoneData( + self, iCurrent: int, variables: List[str], zoneHeaderDict: Dict[str, Any] + ) -> Tuple[TecplotOrderedZone, int]: + """Read the data for an ordered Tecplot zone. + + Parameters + ---------- + iCurrent : int + The current line number in the file. + variables : List[str] + The list of variable names. + zoneHeaderDict : Dict[str, Any] + The zone header information. + + Returns + ------- + Tuple[TecplotOrderedZone, int] + The ordered zone object and the number of lines read. + """ + iMax = zoneHeaderDict["iMax"] + jMax = zoneHeaderDict["jMax"] + kMax = zoneHeaderDict["kMax"] + nNodes = iMax * jMax * kMax + shape = (iMax, jMax, kMax, len(variables)) + + if zoneHeaderDict["datapacking"] == "POINT": + # Point data is column-major + nodalData, nodeOffset = readArrayData(self.filename, iCurrent, nNodes, len(variables)) + nodalData = nodalData.reshape(shape, order="C").squeeze() + else: + # Block data is row-major + nodalData, nodeOffset = readArrayData(self.filename, iCurrent, nNodes, len(variables)) + nodalData = nodalData.reshape(shape, order="F").squeeze() + + data = {var: nodalData[..., i] for i, var in enumerate(variables)} + zone = TecplotOrderedZone( + zoneHeaderDict["zoneName"], + data, + solutionTime=zoneHeaderDict["solutionTime"], + strandID=zoneHeaderDict["strandID"], + ) + + return zone, nodeOffset + + def _readFEZoneData( + self, iCurrent: int, variables: List[str], zoneHeaderDict: Dict[str, Any] + ) -> Tuple[TecplotFEZone, int]: + """Read the data for a finite element Tecplot zone. + + Parameters + ---------- + iCurrent : int + The current line number in the file. + variables : List[str] + The list of variable names. + zoneHeaderDict : Dict[str, Any] + The zone header information. + + Returns + ------- + Tuple[TecplotFEZone, int] + The finite element zone object and the number of lines read. + """ + nNodes = zoneHeaderDict["nNodes"] + nElements = zoneHeaderDict["nElements"] + + if zoneHeaderDict["datapacking"] == "POINT": + # Point data is column-major + nodalData, nodeOffset = readArrayData(self.filename, iCurrent, nNodes, len(variables)) + nodalData = nodalData.reshape(nNodes, len(variables), order="C") + else: + # Block data is row-major + nodalData, nodeOffset = readArrayData(self.filename, iCurrent, nNodes, len(variables)) + nodalData = nodalData.reshape(nNodes, len(variables), order="F") + + connectivity = np.loadtxt(self.filename, skiprows=iCurrent + nodeOffset, max_rows=nElements, dtype=int) + + # Check if the nodal data is 1D + if nodalData.ndim == 1: + nodalData = nodalData.reshape(-1, len(variables)) + + if connectivity.ndim == 1: + connectivity = connectivity.reshape(nElements, -1) + + data = {var: nodalData[..., i] for i, var in enumerate(variables)} + zone = TecplotFEZone( + zoneHeaderDict["zoneName"], + data, + connectivity - 1, + zoneType=zoneHeaderDict["zoneType"], + solutionTime=zoneHeaderDict["solutionTime"], + strandID=zoneHeaderDict["strandID"], + ) + + return zone, nodeOffset + nElements + + def _readZoneData(self, lines: List[str], iLine: int, variables: List[str]) -> Tuple[TecplotZone, int]: + """Read the data for a Tecplot zone. + + Parameters + ---------- + lines : List[str] + The list of lines in the Tecplot file. + iLine : int + The current line number in the file. + variables : List[str] + The list of variable names. + + Returns + ------- + Tuple[TecplotZone, int] + The Tecplot zone object and the number of lines read. + """ + zoneHeaderDict, iLine = self._readZoneHeader(lines, iLine) + + if zoneHeaderDict["zoneType"] == "ORDERED": + zone, iOffset = self._readOrderedZoneData(iLine, variables, zoneHeaderDict) + else: + zone, iOffset = self._readFEZoneData(iLine, variables, zoneHeaderDict) + + return zone, iLine + iOffset + + def read(self) -> Tuple[str, List[TecplotZone]]: + """Read the Tecplot file and return the title and zones. + + Returns + ------- + Tuple[str, List[TecplotZone]] + The title of the Tecplot file and a list of Tecplot zones. + + Raises + ------ + ValueError + If the file is not a valid Tecplot file. + ValueError + If the title is missing. + """ + with open(self.filename, "r") as handle: + lines = handle.readlines() + + zones = [] + + # Get the title + title = re.search(r'title\s*=\s*"(.*)"', lines[0], re.IGNORECASE) + if title is None: + raise ValueError("Tecplot file must have a title on the first line.") + title = title.group(1) + + # Get the variable names + variables = re.findall(r'"([^"]*)"', lines[1]) + + iLine = 2 + while iLine < len(lines): + zone, iLine = self._readZoneData(lines, iLine, variables) + zones.append(zone) + + # Skip any empty lines + while iLine < len(lines) and not lines[iLine].strip(): + iLine += 1 + + return title, zones + + +# ============================================================================== +# BINARY READERS +# ============================================================================== +class TecplotBinaryReader: + def __init__(self, filename: Union[str, Path]) -> None: + """Reader for Tecplot files in binary format. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object + """ + self.filename = filename + self._nVariables = 0 + self._variables = [] + + def _readString(self, handle: TextIO) -> str: + """Read a string from a binary file that is null-terminated. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + + Returns + ------- + str + The string read from the file. + """ + result = [] + while True: + data = handle.read(4) + integer = struct.unpack_from("i", data, 0)[0] + + if integer == 0: + break + + result.append(chr(integer)) + + return "".join(result) + + def _readInteger(self, handle: TextIO, offset: int = 0) -> int: + """Read an integer from a binary file as int32. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + int + The integer read from the file. + """ + return int(np.fromfile(handle, dtype=np.int32, count=1, offset=offset)[0]) + + def _readIntegerArray(self, handle: TextIO, nValues: int, offset: int = 0) -> npt.NDArray[np.int32]: + """Read an array of integers from a binary file as int32. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + nValues : int + The number of values to read. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + npt.NDArray[np.int32] + The integer array read from the file. + """ + return np.fromfile(handle, dtype=np.int32, count=nValues, offset=offset) + + def _readFloat32(self, handle: TextIO, offset: int = 0) -> float: + """Read a float from a binary file as float32. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + float + The float read from the file. + """ + return float(np.fromfile(handle, dtype=np.float32, count=1, offset=offset)[0]) + + def _readFloat32Array(self, handle: TextIO, nValues: int, offset: int = 0) -> npt.NDArray[np.float32]: + """Read an array of floats from a binary file as float32. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + nValues : int + The number of values to read. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + npt.NDArray[np.float32] + The float array read from the file. + """ + return np.fromfile(handle, dtype=np.float32, count=nValues, offset=offset) + + def _readFloat64(self, handle: TextIO, offset: int = 0) -> float: + """Read a float from a binary file as float64. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + float + The float read from the file. + """ + return float(np.fromfile(handle, dtype=np.float64, count=1, offset=offset)[0]) + + def _readFloat64Array(self, handle: TextIO, nValues: int, offset: int = 0) -> npt.NDArray[np.float64]: + """Read an array of floats from a binary file as float64. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + nValues : int + The number of values to read. + offset : int, optional + The offset (in bytes) from the file's current position, + by default 0 + + Returns + ------- + npt.NDArray[np.float64] + The float array read from the file. + """ + return np.fromfile(handle, dtype=np.float64, count=nValues, offset=offset) + + def _readOrderedZone(self, handle: TextIO, zoneName: str, strandID: int, solutionTime: float) -> TecplotOrderedZone: + """Read an ordered Tecplot zone header from a binary file. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + zoneName : str + The name of the zone. + strandID : int + The strand ID. + solutionTime : float + The solution time. + + Returns + ------- + TecplotOrderedZone + The ordered Tecplot zone object. + """ + iMax = self._readInteger(handle) + jMax = self._readInteger(handle) + kMax = self._readInteger(handle) + + return TecplotOrderedZone( + zoneName, + {var: np.zeros((iMax, jMax, kMax)).squeeze() for _, var in enumerate(self._variables)}, + solutionTime=solutionTime, + strandID=strandID, + ) + + def _readFEZone( + self, handle: TextIO, zoneName: str, zoneType: int, strandID: int, solutionTime: float + ) -> TecplotFEZone: + """Read a finite element Tecplot zone header from a binary file. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + zoneName : str + The name of the zone. + zoneType : int + The zone type. + strandID : int + The strand ID. + solutionTime : float + The solution time. + + Returns + ------- + TecplotFEZone + The finite element Tecplot zone object. + + Raises + ------ + ValueError + If the zone type is invalid. + """ + nNodes = self._readInteger(handle) + nElements = self._readInteger(handle) + iCellDim = self._readInteger(handle) # NOQA: F841 + jCellDim = self._readInteger(handle) # NOQA: F841 + kCellDim = self._readInteger(handle) # NOQA: F841 + + if zoneType == ZoneType.FELINESEG.value: + connectivity = np.zeros((nElements, 2), dtype=int) + elif zoneType == ZoneType.FETRIANGLE.value: + connectivity = np.zeros((nElements, 3), dtype=int) + elif zoneType == ZoneType.FEQUADRILATERAL.value: + connectivity = np.zeros((nElements, 4), dtype=int) + elif zoneType == ZoneType.FETETRAHEDRON.value: + connectivity = np.zeros((nElements, 4), dtype=int) + elif zoneType == ZoneType.FEBRICK.value: + connectivity = np.zeros((nElements, 8), dtype=int) + else: + raise ValueError("Invalid zone type.") + + return TecplotFEZone( + zoneName, + {var: np.zeros(nNodes) for _, var in enumerate(self._variables)}, + connectivity, + zoneType=ZoneType(zoneType), + solutionTime=solutionTime, + strandID=strandID, + ) + + def _readZoneHeader(self, handle: TextIO) -> Union[TecplotOrderedZone, TecplotFEZone]: + """Read a Tecplot zone header from a binary file. + + Parameters + ---------- + handle : TextIO + The file handle to read from. + + Returns + ------- + Union[TecplotOrderedZone, TecplotFEZone] + The Tecplot zone object, either ordered or finite element. + """ + zoneName = self._readString(handle) + parentZone = self._readInteger(handle) # NOQA: F841 + strandID = self._readInteger(handle) + solutionTime = self._readFloat64(handle) + defaultColor = self._readInteger(handle) # NOQA: F841 + zoneType = self._readInteger(handle) + datapacking = self._readInteger(handle) # NOQA: F841 + variableLocation = self._readInteger(handle) # NOQA: F841 + rawFaceNeighbors = self._readInteger(handle) # NOQA: F841 + + if zoneType == ZoneType.ORDERED.value: + zone = self._readOrderedZone(handle, zoneName, strandID, solutionTime) + else: + zone = self._readFEZone(handle, zoneName, zoneType, strandID, solutionTime) + + return zone + + def read(self) -> Tuple[str, List[TecplotZone]]: + """Read the Tecplot file and return the title and zones. + + Returns + ------- + Tuple[str, List[TecplotZone]] + The title of the Tecplot file and a list of Tecplot zones. + + Raises + ------ + ValueError + If the file is not a valid Tecplot file. + """ + file = open(self.filename, "rb") + file.seek(0, 2) + fileSize = file.tell() + file.seek(0) + + magic = file.read(8).decode("utf-8") + if magic != "#!TDV112": + raise ValueError("Invalid Tecplot binary file version.") + + byteOrder = self._readInteger(file) # NOQA: F841 + filetype = self._readInteger(file) # NOQA: F841 + title = self._readString(file) + self._nVariables = self._readInteger(file) + self._variables = [self._readString(file) for _ in range(self._nVariables)] + zones: List[Union[TecplotOrderedZone, TecplotFEZone]] = [] + + # Read all the zone headers + while True: + marker = self._readFloat32(file) + + if marker == SectionMarkers.ZONE.value: + # Initialize zone from the header + zone = self._readZoneHeader(file) + zones.append(zone) + if marker == SectionMarkers.DATA.value: + break + + # Read the data for each zone + izone = 0 + while file.tell() < fileSize: + zoneMarker = self._readFloat32(file) # NOQA: F841 + dataFormats = [self._readInteger(file) for _ in range(self._nVariables)] + passiveVariables = self._readInteger(file) # NOQA: F841 + variableSharing = self._readInteger(file) # NOQA: F841 + connSharing = self._readInteger(file) # NOQA: F841 + minMaxArray = self._readFloat64Array(file, 2 * self._nVariables).reshape(self._nVariables, 2) # NOQA: F841 + + if isinstance(zones[izone], TecplotOrderedZone): + iMax = zones[izone].iMax + jMax = zones[izone].jMax + kMax = zones[izone].kMax + + for i in range(self._nVariables): + if dataFormats[i] == BinaryDataPrecisionCodes.SINGLE.value: + readData = self._readFloat32Array(file, iMax * jMax * kMax).reshape(iMax, jMax, kMax).squeeze() + else: + readData = self._readFloat64Array(file, iMax * jMax * kMax).reshape(iMax, jMax, kMax).squeeze() + + zones[izone].data[self._variables[i]][...] = readData + + if isinstance(zones[izone], TecplotFEZone): + nNodes = zones[izone].nNodes + nElements = zones[izone].nElements + + for i in range(self._nVariables): + if dataFormats[i] == BinaryDataPrecisionCodes.SINGLE.value: + readData = self._readFloat32Array(file, nNodes) + else: + readData = self._readFloat64Array(file, nNodes) + + zones[izone].data[self._variables[i]][...] = readData + + connectivitySize = zones[izone].connectivity.size + connectivity = self._readIntegerArray(file, connectivitySize).reshape(nElements, -1) + zones[izone].connectivity = connectivity + + izone += 1 + + file.close() + + return title, zones + + +# ============================================================================== +# PUBLIC FUNCTIONS +# ============================================================================== +def writeTecplot( + filename: Union[str, Path], + title: str, + zones: List[TecplotZone], + datapacking: Literal["BLOCK", "POINT"] = "POINT", + precision: Literal["SINGLE", "DOUBLE"] = "SINGLE", + separator: Separator = Separator.SPACE, +) -> None: + """Write a Tecplot file to disk. The file format is determined by the + file extension. If the extension is .plt, the file will be written in + binary format. If the extension is .dat, the file will be written in + ASCII format. + + .. note:: + + - ASCII files can be written with either BLOCK or POINT data packing. + - Binary files are always written with BLOCK data packing. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object. + title : str + The title of the Tecplot file. + zones : List[TecplotZone] + A list of Tecplot zones to write + datapacking : Literal["BLOCK", "POINT"], optional + The data packing format. BLOCK is row-major, POINT is + column-major, by default "POINT" + precision : Literal["SINGLE", "DOUBLE"], optional + The floating point precision to write the data, by default + "SINGLE" + separator : Separator, optional + The separator to use when writing ASCII files. The Separator + is an enum defined in :meth:`Separator `, + by default Separator.SPACE + + Raises + ------ + ValueError + If the file extension is invalid. + + Examples + -------- + .. code-block:: python + + from baseclasses import tecplotIO as tpio + import numpy as np + + nx, ny, nz = 10, 10, 10 + X, Y, Z = np.meshgrid(np.random.rand(nx), np.random.rand(ny), np.random.rand(nz), indexing="ij") + data = {"X": X, "Y": Y, "Z": Z} + zone = tpio.TecplotOrderedZone("OrderedZone", data) + + # Write the Tecplot file in ASCII format + tpio.writeTecplot("ordered.dat", "Ordered Zone", [zone], datapacking="BLOCK", precision="SINGLE", separator=tpio.Separator.SPACE) + + # Write the Tecplot file in binary format + tpio.writeTecplot("ordered.plt", "Ordered Zone", [zone], precision="SINGLE") + + """ + filepath = Path(filename) + if filepath.suffix == ".plt": + writer = TecplotWriterBinary(title, zones, precision) + writer.write(filepath) + elif filepath.suffix == ".dat": + writer = TecplotWriterASCII(title, zones, datapacking, precision, separator) + writer.write(filename) + else: + raise ValueError("Invalid file extension. Must be .plt (binary) or .dat (ASCII).") + + +def readTecplot(filename: Union[str, Path]) -> Tuple[str, List[Union[TecplotOrderedZone, TecplotFEZone]]]: + """Read a Tecplot file from disk. The file format is determined by + the file extension. If the extension is .plt, the file will be read + in binary format. If the extension is .dat, the file will be read in + ASCII format. + + Parameters + ---------- + filename : Union[str, Path] + The filename as a string or pathlib.Path object. + + Returns + ------- + Tuple[str, List[Union[TecplotOrderedZone, TecplotFEZone]]] + The title of the Tecplot file and a list of Tecplot zones. + + Raises + ------ + ValueError + If the file extension is invalid. + + Examples + -------- + .. code-block:: python + + from baseclasses import tecplotIO as tpio + + # Read a Tecplot file in ASCII format + title, zones = tpio.readTecplot("ordered.dat") + + # Read a Tecplot file in binary format + title, zones = tpio.readTecplot("ordered.plt") + """ + filepath = Path(filename) + if filepath.suffix == ".plt": + reader = TecplotBinaryReader(filename) + title, zones = reader.read() + elif filepath.suffix == ".dat": + reader = TecplotASCIIReader(filename) + title, zones = reader.read() + else: + raise ValueError("Invalid file extension. Must be .plt (binary) or .dat (ASCII).") + + return title, zones diff --git a/tests/input/airfoil_000_slices.dat b/tests/input/airfoil_000_slices.dat new file mode 100644 index 0000000..a9f7eef --- /dev/null +++ b/tests/input/airfoil_000_slices.dat @@ -0,0 +1,285 @@ + Title = "ADflow Slice Data" +Variables = "CoordinateX" "CoordinateY" "CoordinateZ" "XoC" "YoC" "ZoC" "Density" "VelocityX" "VelocityY" "CoefPressure" "Mach" "SkinFrictionMagnitude" "SkinFrictionX" +Zone T= "Slice_0001 None Absolute Regular z = 0.50000000" + Nodes = 142 Elements= 138 ZONETYPE=FELINESEG + DATAPACKING=POINT + 7.823774E-01 1.691034E-02 5.000000E-01 7.823466E-01 1.690968E-02 0.000000E+00 8.769829E-01 1.730982E-02 -1.106485E-03 -8.302255E-02 0.000000E+00 2.180030E-03 2.175599E-03 + 7.622641E-01 1.816350E-02 5.000000E-01 7.622340E-01 1.816278E-02 0.000000E+00 8.795820E-01 1.661022E-02 -9.950567E-04 -7.566586E-02 0.000000E+00 2.147108E-03 2.143266E-03 + 7.400762E-01 1.943594E-02 5.000000E-01 7.400470E-01 1.943518E-02 0.000000E+00 8.840950E-01 1.568833E-02 -9.324715E-04 -6.295219E-02 0.000000E+00 2.082592E-03 2.078989E-03 + 7.157196E-01 2.090906E-02 5.000000E-01 7.156914E-01 2.090823E-02 0.000000E+00 8.812874E-01 1.467928E-02 -1.014514E-03 -7.089337E-02 0.000000E+00 2.002484E-03 1.997728E-03 + 6.891330E-01 2.295514E-02 5.000000E-01 6.891058E-01 2.295424E-02 0.000000E+00 8.696365E-01 1.400634E-02 -1.156615E-03 -1.038830E-01 0.000000E+00 1.963482E-03 1.956759E-03 + 6.602989E-01 2.550903E-02 5.000000E-01 6.602729E-01 2.550802E-02 0.000000E+00 8.576853E-01 1.397531E-02 -1.201179E-03 -1.380027E-01 0.000000E+00 2.009371E-03 2.001965E-03 + 6.292542E-01 2.810639E-02 5.000000E-01 6.292294E-01 2.810529E-02 0.000000E+00 8.473740E-01 1.434909E-02 -1.171540E-03 -1.677053E-01 0.000000E+00 2.110641E-03 2.103706E-03 + 5.961002E-01 3.072295E-02 5.000000E-01 5.960767E-01 3.072174E-02 0.000000E+00 8.339816E-01 1.472225E-02 -1.155850E-03 -2.060563E-01 0.000000E+00 2.211227E-03 2.204499E-03 + 5.610108E-01 3.344136E-02 5.000000E-01 5.609887E-01 3.344005E-02 0.000000E+00 8.182221E-01 1.495282E-02 -1.125149E-03 -2.509287E-01 0.000000E+00 2.288083E-03 2.281663E-03 + 5.242375E-01 3.611316E-02 5.000000E-01 5.242168E-01 3.611174E-02 0.000000E+00 8.046896E-01 1.497465E-02 -1.054999E-03 -2.893953E-01 0.000000E+00 2.327903E-03 2.322113E-03 + 4.861092E-01 3.873006E-02 5.000000E-01 4.860900E-01 3.872853E-02 0.000000E+00 7.989295E-01 1.474552E-02 -9.504222E-04 -3.059910E-01 0.000000E+00 2.321257E-03 2.316348E-03 + 4.470274E-01 4.112937E-02 5.000000E-01 4.470098E-01 4.112775E-02 0.000000E+00 7.977514E-01 1.429271E-02 -8.107427E-04 -3.098850E-01 0.000000E+00 2.269757E-03 2.266117E-03 + 4.074550E-01 4.316362E-02 5.000000E-01 4.074389E-01 4.316192E-02 0.000000E+00 7.788086E-01 1.374141E-02 -6.616893E-04 -3.640119E-01 0.000000E+00 2.191907E-03 2.189512E-03 + 3.678986E-01 4.480381E-02 5.000000E-01 3.678841E-01 4.480205E-02 0.000000E+00 7.199728E-01 1.368958E-02 -5.040894E-04 -5.307141E-01 0.000000E+00 2.182102E-03 2.180701E-03 + 3.288871E-01 4.595042E-02 5.000000E-01 3.288741E-01 4.594861E-02 0.000000E+00 6.318774E-01 1.547437E-02 -3.337712E-04 -7.804645E-01 0.000000E+00 2.449687E-03 2.449031E-03 + 2.909459E-01 4.656835E-02 5.000000E-01 2.909344E-01 4.656652E-02 0.000000E+00 5.566506E-01 2.033488E-02 -1.234061E-04 -9.951881E-01 0.000000E+00 3.177438E-03 3.177259E-03 + 2.545710E-01 4.660779E-02 5.000000E-01 2.545610E-01 4.660596E-02 0.000000E+00 5.244056E-01 2.735012E-02 2.374738E-04 -1.089817E+00 0.000000E+00 4.198999E-03 4.198697E-03 + 2.202045E-01 4.605987E-02 5.000000E-01 2.201958E-01 4.605806E-02 0.000000E+00 5.277605E-01 3.342428E-02 8.403877E-04 -1.083680E+00 0.000000E+00 5.024009E-03 5.022085E-03 + 1.882140E-01 4.494129E-02 5.000000E-01 1.882066E-01 4.493952E-02 0.000000E+00 5.442207E-01 3.669597E-02 1.644182E-03 -1.038860E+00 0.000000E+00 5.379277E-03 5.373338E-03 + 1.588780E-01 4.328696E-02 5.000000E-01 1.588718E-01 4.328525E-02 0.000000E+00 5.608540E-01 3.811115E-02 2.573374E-03 -9.919544E-01 0.000000E+00 5.423373E-03 5.410392E-03 + 1.323788E-01 4.116000E-02 5.000000E-01 1.323735E-01 4.115838E-02 0.000000E+00 5.760369E-01 3.947691E-02 3.681957E-03 -9.485354E-01 0.000000E+00 5.426116E-03 5.401896E-03 + 1.088015E-01 3.862899E-02 5.000000E-01 1.087972E-01 3.862747E-02 0.000000E+00 5.923225E-01 4.177502E-02 5.127649E-03 -9.021378E-01 0.000000E+00 5.520156E-03 5.478123E-03 + 8.814071E-02 3.576593E-02 5.000000E-01 8.813724E-02 3.576452E-02 0.000000E+00 6.122635E-01 4.500841E-02 7.060907E-03 -8.457289E-01 0.000000E+00 5.696456E-03 5.626584E-03 + 7.031156E-02 3.264592E-02 5.000000E-01 7.030879E-02 3.264464E-02 0.000000E+00 6.376871E-01 4.874236E-02 9.594187E-03 -7.739780E-01 0.000000E+00 5.896192E-03 5.784006E-03 + 5.516427E-02 2.934346E-02 5.000000E-01 5.516210E-02 2.934231E-02 0.000000E+00 6.694468E-01 5.250913E-02 1.278610E-02 -6.841842E-01 0.000000E+00 6.073830E-03 5.899995E-03 + 4.249863E-02 2.594613E-02 5.000000E-01 4.249695E-02 2.594511E-02 0.000000E+00 7.073704E-01 5.578594E-02 1.658521E-02 -5.764983E-01 0.000000E+00 6.200159E-03 5.941386E-03 + 3.203145E-02 2.253495E-02 5.000000E-01 3.203019E-02 2.253406E-02 0.000000E+00 7.516866E-01 5.776977E-02 2.081023E-02 -4.499086E-01 0.000000E+00 6.239092E-03 5.867480E-03 + 2.331672E-02 1.909245E-02 5.000000E-01 2.331580E-02 1.909169E-02 0.000000E+00 8.025574E-01 5.712975E-02 2.494108E-02 -3.034481E-01 0.000000E+00 6.103029E-03 5.589531E-03 + 1.591438E-02 1.552938E-02 5.000000E-01 1.591375E-02 1.552877E-02 0.000000E+00 8.601468E-01 5.237826E-02 2.801593E-02 -1.361041E-01 0.000000E+00 5.656655E-03 4.982834E-03 + 9.658338E-03 1.178916E-02 5.000000E-01 9.657958E-03 1.178869E-02 0.000000E+00 9.328105E-01 4.272251E-02 2.900778E-02 7.562705E-02 0.000000E+00 4.796930E-03 3.963653E-03 + 4.674145E-03 7.870585E-03 5.000000E-01 4.673961E-03 7.870275E-03 0.000000E+00 1.044568E+00 2.814464E-02 2.629317E-02 3.982327E-01 0.000000E+00 3.457744E-03 2.529153E-03 + 1.281341E-03 3.839618E-03 5.000000E-01 1.281290E-03 3.839466E-03 0.000000E+00 1.206823E+00 1.132702E-02 1.599617E-02 8.615307E-01 0.000000E+00 1.678967E-03 9.668870E-04 + 2.352898E-18 -5.130713E-19 5.000000E-01 0.000000E+00 0.000000E+00 0.000000E+00 1.290470E+00 6.340877E-03 -6.986974E-03 1.099159E+00 0.000000E+00 1.229277E-03 4.491090E-04 + 1.281240E-03 -3.074255E-03 5.000000E-01 1.281190E-03 -3.074134E-03 0.000000E+00 1.160849E+00 2.173498E-02 -2.338417E-02 7.311508E-01 0.000000E+00 2.482573E-03 1.645839E-03 + 4.673777E-03 -5.479147E-03 5.000000E-01 4.673593E-03 -5.478931E-03 0.000000E+00 9.857179E-01 3.948160E-02 -2.127343E-02 2.324595E-01 0.000000E+00 3.598765E-03 3.162125E-03 + 9.657579E-03 -7.462679E-03 5.000000E-01 9.657198E-03 -7.462385E-03 0.000000E+00 9.183535E-01 4.666479E-02 -1.538191E-02 3.750595E-02 0.000000E+00 4.131688E-03 3.918618E-03 + 1.000039E+00 4.335176E-04 5.000000E-01 9.999999E-01 4.335005E-04 0.000000E+00 9.178565E-01 3.218372E-02 -5.635738E-03 3.712361E-02 0.000000E+00 1.164162E-03 9.562320E-04 + 9.998942E-01 4.463465E-04 5.000000E-01 9.998548E-01 4.463289E-04 0.000000E+00 9.532969E-01 3.939878E-02 -6.167425E-03 1.362680E-01 0.000000E+00 1.434919E-03 1.429344E-03 + 9.996793E-01 4.653449E-04 5.000000E-01 9.996399E-01 4.653266E-04 0.000000E+00 9.334895E-01 2.687128E-02 -1.948203E-03 8.086840E-02 0.000000E+00 1.185505E-03 1.180902E-03 + 9.993613E-01 4.934446E-04 5.000000E-01 9.993219E-01 4.934251E-04 0.000000E+00 9.375248E-01 3.241384E-02 -2.952429E-03 9.193923E-02 0.000000E+00 1.512094E-03 1.506228E-03 + 9.988916E-01 5.349283E-04 5.000000E-01 9.988522E-01 5.349072E-04 0.000000E+00 9.376982E-01 3.108600E-02 -2.733864E-03 9.228407E-02 0.000000E+00 1.663165E-03 1.656720E-03 + 9.981996E-01 5.960037E-04 5.000000E-01 9.981603E-01 5.959802E-04 0.000000E+00 9.374421E-01 2.603351E-02 -2.293258E-03 9.144503E-02 0.000000E+00 1.650891E-03 1.644504E-03 + 9.971840E-01 6.855608E-04 5.000000E-01 9.971447E-01 6.855338E-04 0.000000E+00 9.371065E-01 2.062899E-02 -1.817270E-03 9.033775E-02 0.000000E+00 1.582297E-03 1.576190E-03 + 9.957013E-01 8.161055E-04 5.000000E-01 9.956620E-01 8.160733E-04 0.000000E+00 9.366509E-01 1.628945E-02 -1.435026E-03 8.885614E-02 0.000000E+00 1.522780E-03 1.516925E-03 + 9.935539E-01 1.004757E-03 5.000000E-01 9.935147E-01 1.004718E-03 0.000000E+00 9.356150E-01 1.325044E-02 -1.167588E-03 8.572860E-02 0.000000E+00 1.495616E-03 1.489896E-03 + 9.904791E-01 1.273998E-03 5.000000E-01 9.904401E-01 1.273948E-03 0.000000E+00 9.332168E-01 1.141662E-02 -9.996309E-04 7.876507E-02 0.000000E+00 1.497252E-03 1.491572E-03 + 9.861469E-01 1.651517E-03 5.000000E-01 9.861080E-01 1.651452E-03 0.000000E+00 9.305094E-01 1.060698E-02 -9.212272E-04 7.092268E-02 0.000000E+00 1.509788E-03 1.504120E-03 + 9.821027E-01 2.001929E-03 5.000000E-01 9.820640E-01 2.001851E-03 0.000000E+00 9.285245E-01 1.051754E-02 -9.105222E-04 6.511007E-02 0.000000E+00 1.526705E-03 1.521040E-03 + 9.778355E-01 2.369477E-03 5.000000E-01 9.777970E-01 2.369384E-03 0.000000E+00 9.265370E-01 1.082341E-02 -9.310436E-04 5.930848E-02 0.000000E+00 1.554483E-03 1.548785E-03 + 9.733502E-01 2.753336E-03 5.000000E-01 9.733118E-01 2.753227E-03 0.000000E+00 9.243767E-01 1.130566E-02 -9.654937E-04 5.304624E-02 0.000000E+00 1.590739E-03 1.584986E-03 + 9.686497E-01 3.152840E-03 5.000000E-01 9.686116E-01 3.152715E-03 0.000000E+00 9.222332E-01 1.184954E-02 -1.005224E-03 4.685132E-02 0.000000E+00 1.628497E-03 1.622689E-03 + 9.637331E-01 3.567703E-03 5.000000E-01 9.636952E-01 3.567562E-03 0.000000E+00 9.200808E-01 1.240052E-02 -1.044457E-03 4.064388E-02 0.000000E+00 1.665345E-03 1.659492E-03 + 9.585926E-01 3.998262E-03 5.000000E-01 9.585548E-01 3.998105E-03 0.000000E+00 9.179270E-01 1.293799E-02 -1.081959E-03 3.443700E-02 0.000000E+00 1.701687E-03 1.695793E-03 + 9.532109E-01 4.445688E-03 5.000000E-01 9.531734E-01 4.445513E-03 0.000000E+00 9.157741E-01 1.344964E-02 -1.116836E-03 2.823230E-02 0.000000E+00 1.737908E-03 1.731976E-03 + 9.475603E-01 4.912065E-03 5.000000E-01 9.475229E-01 4.911872E-03 0.000000E+00 9.136108E-01 1.392595E-02 -1.148354E-03 2.199604E-02 0.000000E+00 1.773831E-03 1.767864E-03 + 9.416026E-01 5.400229E-03 5.000000E-01 9.415655E-01 5.400017E-03 0.000000E+00 9.114136E-01 1.436534E-02 -1.176172E-03 1.566132E-02 0.000000E+00 1.809263E-03 1.803266E-03 + 9.352941E-01 5.913251E-03 5.000000E-01 9.352572E-01 5.913018E-03 0.000000E+00 9.091619E-01 1.477763E-02 -1.200471E-03 9.170598E-03 0.000000E+00 1.844320E-03 1.838304E-03 + 9.285925E-01 6.453568E-03 5.000000E-01 9.285559E-01 6.453313E-03 0.000000E+00 9.068578E-01 1.518004E-02 -1.221647E-03 2.530993E-03 0.000000E+00 1.879184E-03 1.873170E-03 + 9.214652E-01 7.022170E-03 5.000000E-01 9.214289E-01 7.021894E-03 0.000000E+00 9.045332E-01 1.558881E-02 -1.240143E-03 -4.165837E-03 0.000000E+00 1.913765E-03 1.907782E-03 + 9.138872E-01 7.618947E-03 5.000000E-01 9.138512E-01 7.618646E-03 0.000000E+00 9.022295E-01 1.601212E-02 -1.256568E-03 -1.080105E-02 0.000000E+00 1.947607E-03 1.941684E-03 + 9.058310E-01 8.244078E-03 5.000000E-01 9.057953E-01 8.243753E-03 0.000000E+00 8.999812E-01 1.644836E-02 -1.271797E-03 -1.727622E-02 0.000000E+00 1.980050E-03 1.974206E-03 + 8.972635E-01 8.898770E-03 5.000000E-01 8.972281E-01 8.898419E-03 0.000000E+00 8.978009E-01 1.688845E-02 -1.286974E-03 -2.355606E-02 0.000000E+00 2.010552E-03 2.004790E-03 + 8.881456E-01 9.585664E-03 5.000000E-01 8.881106E-01 9.585287E-03 0.000000E+00 8.956373E-01 1.731879E-02 -1.303141E-03 -2.978710E-02 0.000000E+00 2.039106E-03 2.033411E-03 + 8.784302E-01 1.030898E-02 5.000000E-01 8.783955E-01 1.030858E-02 0.000000E+00 8.933518E-01 1.772203E-02 -1.319517E-03 -3.636127E-02 0.000000E+00 2.066522E-03 2.060875E-03 + 8.680511E-01 1.107397E-02 5.000000E-01 8.680169E-01 1.107353E-02 0.000000E+00 8.908439E-01 1.807475E-02 -1.330781E-03 -4.356164E-02 0.000000E+00 2.094151E-03 2.088557E-03 + 8.569008E-01 1.188539E-02 5.000000E-01 8.568670E-01 1.188492E-02 0.000000E+00 8.883060E-01 1.834260E-02 -1.327444E-03 -5.084494E-02 0.000000E+00 2.122633E-03 2.117151E-03 + 8.448105E-01 1.274679E-02 5.000000E-01 8.447772E-01 1.274629E-02 0.000000E+00 8.862150E-01 1.848015E-02 -1.302759E-03 -5.686511E-02 0.000000E+00 2.150373E-03 2.145099E-03 + 8.315513E-01 1.366267E-02 5.000000E-01 8.315185E-01 1.366213E-02 0.000000E+00 8.846171E-01 1.844526E-02 -1.262203E-03 -6.148111E-02 0.000000E+00 2.173458E-03 2.168446E-03 + 8.168677E-01 1.464627E-02 5.000000E-01 8.168355E-01 1.464569E-02 0.000000E+00 8.824909E-01 1.822369E-02 -1.220767E-03 -6.754388E-02 0.000000E+00 2.187738E-03 2.182920E-03 + 8.005336E-01 1.572371E-02 5.000000E-01 8.005021E-01 1.572309E-02 0.000000E+00 8.791559E-01 1.783739E-02 -1.178774E-03 -7.694478E-02 0.000000E+00 2.190819E-03 2.186112E-03 + 7.823774E-01 1.691034E-02 5.000000E-01 7.823466E-01 1.690968E-02 0.000000E+00 8.769829E-01 1.730982E-02 -1.106485E-03 -8.302255E-02 0.000000E+00 2.180030E-03 2.175599E-03 + 8.568334E-01 -1.615150E-03 5.000000E-01 8.567996E-01 -1.615086E-03 0.000000E+00 9.418312E-01 2.095896E-02 3.658859E-04 1.005031E-01 0.000000E+00 2.420404E-03 2.420054E-03 + 8.679829E-01 -1.431087E-03 5.000000E-01 8.679487E-01 -1.431030E-03 0.000000E+00 9.425060E-01 2.084891E-02 3.425828E-04 1.024651E-01 0.000000E+00 2.410123E-03 2.409815E-03 + 8.783611E-01 -1.270333E-03 5.000000E-01 8.783265E-01 -1.270283E-03 0.000000E+00 9.431044E-01 2.064718E-02 3.181916E-04 1.041987E-01 0.000000E+00 2.401830E-03 2.401560E-03 + 8.880758E-01 -1.129756E-03 5.000000E-01 8.880408E-01 -1.129712E-03 0.000000E+00 9.436463E-01 2.039254E-02 2.932791E-04 1.057626E-01 0.000000E+00 2.394905E-03 2.394672E-03 + 8.971929E-01 -1.007254E-03 5.000000E-01 8.971576E-01 -1.007214E-03 0.000000E+00 9.441251E-01 2.012229E-02 2.683172E-04 1.071401E-01 0.000000E+00 2.389094E-03 2.388895E-03 + 9.057598E-01 -9.011546E-04 5.000000E-01 9.057241E-01 -9.011191E-04 0.000000E+00 9.445325E-01 1.986380E-02 2.438155E-04 1.083097E-01 0.000000E+00 2.384381E-03 2.384213E-03 + 9.138154E-01 -8.099012E-04 5.000000E-01 9.137794E-01 -8.098693E-04 0.000000E+00 9.448696E-01 1.963215E-02 2.202656E-04 1.092763E-01 0.000000E+00 2.380714E-03 2.380574E-03 + 9.213928E-01 -7.319502E-04 5.000000E-01 9.213565E-01 -7.319213E-04 0.000000E+00 9.451449E-01 1.943061E-02 1.980562E-04 1.100646E-01 0.000000E+00 2.377825E-03 2.377711E-03 + 9.285195E-01 -6.657734E-04 5.000000E-01 9.284829E-01 -6.657472E-04 0.000000E+00 9.453682E-01 1.925154E-02 1.773417E-04 1.107029E-01 0.000000E+00 2.375301E-03 2.375208E-03 + 9.352206E-01 -6.099201E-04 5.000000E-01 9.351837E-01 -6.098960E-04 0.000000E+00 9.455488E-01 1.907820E-02 1.577954E-04 1.112177E-01 0.000000E+00 2.372875E-03 2.372802E-03 + 9.415286E-01 -5.631533E-04 5.000000E-01 9.414915E-01 -5.631311E-04 0.000000E+00 9.456946E-01 1.888898E-02 1.387161E-04 1.116305E-01 0.000000E+00 2.370603E-03 2.370545E-03 + 9.474858E-01 -5.244859E-04 5.000000E-01 9.474484E-01 -5.244653E-04 0.000000E+00 9.458081E-01 1.866310E-02 1.195926E-04 1.119468E-01 0.000000E+00 2.368645E-03 2.368602E-03 + 9.531360E-01 -4.931011E-04 5.000000E-01 9.530984E-01 -4.930817E-04 0.000000E+00 9.458866E-01 1.838528E-02 1.004076E-04 1.121561E-01 0.000000E+00 2.366900E-03 2.366869E-03 + 9.585172E-01 -4.682654E-04 5.000000E-01 9.584794E-01 -4.682470E-04 0.000000E+00 9.459277E-01 1.804749E-02 8.156609E-05 1.122485E-01 0.000000E+00 2.364747E-03 2.364727E-03 + 9.636574E-01 -4.492731E-04 5.000000E-01 9.636194E-01 -4.492554E-04 0.000000E+00 9.459282E-01 1.764847E-02 6.362841E-05 1.122142E-01 0.000000E+00 2.360922E-03 2.360909E-03 + 9.685736E-01 -4.354232E-04 5.000000E-01 9.685354E-01 -4.354061E-04 0.000000E+00 9.458838E-01 1.719666E-02 4.712843E-05 1.120444E-01 0.000000E+00 2.353896E-03 2.353889E-03 + 9.732737E-01 -4.260205E-04 5.000000E-01 9.732353E-01 -4.260037E-04 0.000000E+00 9.457942E-01 1.672631E-02 3.215995E-05 1.117490E-01 0.000000E+00 2.343695E-03 2.343692E-03 + 9.777586E-01 -4.203871E-04 5.000000E-01 9.777201E-01 -4.203705E-04 0.000000E+00 9.456877E-01 1.633673E-02 2.027334E-05 1.114240E-01 0.000000E+00 2.336193E-03 2.336192E-03 + 9.820255E-01 -4.178789E-04 5.000000E-01 9.819868E-01 -4.178624E-04 0.000000E+00 9.453731E-01 1.623640E-02 9.331972E-06 1.105504E-01 0.000000E+00 2.346100E-03 2.346099E-03 + 9.860694E-01 -4.178996E-04 5.000000E-01 9.860305E-01 -4.178832E-04 0.000000E+00 9.451589E-01 1.681199E-02 -1.387892E-06 1.100018E-01 0.000000E+00 2.380906E-03 2.380906E-03 + 9.904012E-01 -4.203243E-04 5.000000E-01 9.903622E-01 -4.203077E-04 0.000000E+00 9.455695E-01 1.858514E-02 -7.756270E-06 1.112473E-01 0.000000E+00 2.423344E-03 2.423344E-03 + 9.934757E-01 -4.234474E-04 5.000000E-01 9.934366E-01 -4.234308E-04 0.000000E+00 9.457769E-01 2.197832E-02 -1.742003E-05 1.119648E-01 0.000000E+00 2.466043E-03 2.466041E-03 + 9.956230E-01 -4.262769E-04 5.000000E-01 9.955837E-01 -4.262601E-04 0.000000E+00 9.452957E-01 2.719731E-02 -3.534340E-05 1.107957E-01 0.000000E+00 2.527701E-03 2.527698E-03 + 9.971056E-01 -4.285279E-04 5.000000E-01 9.970663E-01 -4.285110E-04 0.000000E+00 9.447832E-01 3.426382E-02 -5.186382E-05 1.095708E-01 0.000000E+00 2.612713E-03 2.612709E-03 + 9.981211E-01 -4.302056E-04 5.000000E-01 9.980818E-01 -4.301887E-04 0.000000E+00 9.445112E-01 4.283458E-02 -6.408236E-05 1.090059E-01 0.000000E+00 2.697456E-03 2.697452E-03 + 9.988130E-01 -4.314106E-04 5.000000E-01 9.987737E-01 -4.313936E-04 0.000000E+00 9.445384E-01 5.106542E-02 -1.014609E-04 1.092347E-01 0.000000E+00 2.707156E-03 2.707152E-03 + 9.992827E-01 -4.322566E-04 5.000000E-01 9.992433E-01 -4.322396E-04 0.000000E+00 9.449224E-01 5.442339E-02 2.618544E-05 1.105295E-01 0.000000E+00 2.505741E-03 2.505737E-03 + 9.996007E-01 -4.328422E-04 5.000000E-01 9.995613E-01 -4.328251E-04 0.000000E+00 9.389288E-01 4.678001E-02 -5.834615E-04 9.422103E-02 0.000000E+00 2.013439E-03 2.013435E-03 + 9.998156E-01 -4.332438E-04 5.000000E-01 9.997762E-01 -4.332267E-04 0.000000E+00 9.652188E-01 6.378814E-02 3.115090E-03 1.674000E-01 0.000000E+00 2.308594E-03 2.308590E-03 + 9.999607E-01 -4.335176E-04 5.000000E-01 9.999213E-01 -4.335005E-04 0.000000E+00 9.194356E-01 4.924938E-02 7.209547E-03 3.962806E-02 0.000000E+00 1.592428E-03 1.449650E-03 + 9.999738E-01 -2.890117E-04 5.000000E-01 9.999344E-01 -2.890003E-04 0.000000E+00 8.834192E-01 1.656528E-04 -3.672523E-02 -6.044585E-02 0.000000E+00 1.647927E-03 -1.488262E-04 + 9.999869E-01 -1.445059E-04 5.000000E-01 9.999475E-01 -1.445002E-04 0.000000E+00 9.200216E-01 -6.768082E-03 -5.688034E-02 4.276053E-02 0.000000E+00 2.152794E-03 -1.944214E-04 + 1.000000E+00 5.460349E-14 5.000000E-01 9.999606E-01 5.460185E-14 0.000000E+00 9.146963E-01 -1.246480E-03 -8.292641E-03 2.920251E-02 0.000000E+00 9.568096E-04 -2.843208E-05 + 1.000013E+00 1.445059E-04 5.000000E-01 9.999737E-01 1.445002E-04 0.000000E+00 9.191401E-01 2.078684E-03 3.861437E-02 4.127059E-02 0.000000E+00 1.463157E-03 1.321394E-04 + 1.000026E+00 2.890117E-04 5.000000E-01 9.999868E-01 2.890003E-04 0.000000E+00 8.910295E-01 5.427128E-03 3.024639E-02 -3.787454E-02 0.000000E+00 1.366723E-03 1.234304E-04 + 1.000039E+00 4.335176E-04 5.000000E-01 9.999999E-01 4.335005E-04 0.000000E+00 9.178565E-01 3.218372E-02 -5.635738E-03 3.712361E-02 0.000000E+00 1.164162E-03 9.562320E-04 + 9.657579E-03 -7.462679E-03 5.000000E-01 9.657198E-03 -7.462385E-03 0.000000E+00 9.183535E-01 4.666479E-02 -1.538191E-02 3.750595E-02 0.000000E+00 4.131688E-03 3.918618E-03 + 1.591313E-02 -9.124254E-03 5.000000E-01 1.591250E-02 -9.123894E-03 0.000000E+00 9.025181E-01 4.806517E-02 -1.097433E-02 -1.093925E-02 0.000000E+00 4.324674E-03 4.213519E-03 + 2.331488E-02 -1.055139E-02 5.000000E-01 2.331396E-02 -1.055098E-02 0.000000E+00 8.965480E-01 4.685386E-02 -7.923268E-03 -3.005983E-02 0.000000E+00 4.370616E-03 4.308254E-03 + 3.202893E-02 -1.182037E-02 5.000000E-01 3.202767E-02 -1.181991E-02 0.000000E+00 8.896745E-01 4.398983E-02 -5.632399E-03 -5.053972E-02 0.000000E+00 4.312881E-03 4.277367E-03 + 4.249529E-02 -1.297016E-02 5.000000E-01 4.249361E-02 -1.296965E-02 0.000000E+00 8.832212E-01 4.035761E-02 -3.862412E-03 -6.914083E-02 0.000000E+00 4.199755E-03 4.180288E-03 + 5.515994E-02 -1.399362E-02 5.000000E-01 5.515776E-02 -1.399307E-02 0.000000E+00 8.797259E-01 3.665546E-02 -2.528761E-03 -7.906228E-02 0.000000E+00 4.066228E-03 4.056307E-03 + 7.030603E-02 -1.485197E-02 5.000000E-01 7.030326E-02 -1.485139E-02 0.000000E+00 8.795595E-01 3.316821E-02 -1.550091E-03 -7.943075E-02 0.000000E+00 3.917887E-03 3.913398E-03 + 8.813378E-02 -1.550534E-02 5.000000E-01 8.813031E-02 -1.550473E-02 0.000000E+00 8.819338E-01 3.000358E-02 -8.520976E-04 -7.254173E-02 0.000000E+00 3.757152E-03 3.755484E-03 + 1.087929E-01 -1.592600E-02 5.000000E-01 1.087886E-01 -1.592538E-02 0.000000E+00 8.859163E-01 2.726018E-02 -3.688623E-04 -6.105432E-02 0.000000E+00 3.597488E-03 3.597062E-03 + 1.323683E-01 -1.609599E-02 5.000000E-01 1.323631E-01 -1.609536E-02 0.000000E+00 8.908008E-01 2.498222E-02 -4.091464E-05 -4.697405E-02 0.000000E+00 3.452021E-03 3.451967E-03 + 1.588655E-01 -1.600922E-02 5.000000E-01 1.588593E-01 -1.600859E-02 0.000000E+00 8.962510E-01 2.312065E-02 1.747640E-04 -3.126929E-02 0.000000E+00 3.323192E-03 3.323082E-03 + 1.881992E-01 -1.568312E-02 5.000000E-01 1.881918E-01 -1.568250E-02 0.000000E+00 9.019354E-01 2.158914E-02 3.096616E-04 -1.490025E-02 0.000000E+00 3.206942E-03 3.206618E-03 + 2.201872E-01 -1.514454E-02 5.000000E-01 2.201785E-01 -1.514395E-02 0.000000E+00 9.074908E-01 2.033350E-02 3.860694E-04 1.083300E-03 0.000000E+00 3.101968E-03 3.101426E-03 + 2.545510E-01 -1.444234E-02 5.000000E-01 2.545410E-01 -1.444177E-02 0.000000E+00 9.126535E-01 1.933639E-02 4.179903E-04 1.591620E-02 0.000000E+00 3.011081E-03 3.010400E-03 + 2.909230E-01 -1.363838E-02 5.000000E-01 2.909115E-01 -1.363785E-02 0.000000E+00 9.170070E-01 1.857688E-02 4.230656E-04 2.841757E-02 0.000000E+00 2.935557E-03 2.934819E-03 + 3.288612E-01 -1.277407E-02 5.000000E-01 3.288483E-01 -1.277356E-02 0.000000E+00 9.203483E-01 1.802618E-02 4.164445E-04 3.801695E-02 0.000000E+00 2.874478E-03 2.873731E-03 + 3.678697E-01 -1.188383E-02 5.000000E-01 3.678552E-01 -1.188336E-02 0.000000E+00 9.228653E-01 1.765932E-02 4.040552E-04 4.524603E-02 0.000000E+00 2.826600E-03 2.825878E-03 + 4.074229E-01 -1.099858E-02 5.000000E-01 4.074069E-01 -1.099815E-02 0.000000E+00 9.247671E-01 1.744995E-02 3.883316E-04 5.070356E-02 0.000000E+00 2.789766E-03 2.789092E-03 + 4.469923E-01 -1.014522E-02 5.000000E-01 4.469746E-01 -1.014482E-02 0.000000E+00 9.260519E-01 1.736175E-02 3.758366E-04 5.440127E-02 0.000000E+00 2.759775E-03 2.759145E-03 + 4.860709E-01 -9.318470E-03 5.000000E-01 4.860518E-01 -9.318103E-03 0.000000E+00 9.268147E-01 1.735837E-02 3.730613E-04 5.662453E-02 0.000000E+00 2.732195E-03 2.731580E-03 + 5.241963E-01 -8.507160E-03 5.000000E-01 5.241756E-01 -8.506825E-03 0.000000E+00 9.275234E-01 1.742258E-02 3.766485E-04 5.869957E-02 0.000000E+00 2.705529E-03 2.704913E-03 + 5.609667E-01 -7.720093E-03 5.000000E-01 5.609446E-01 -7.719789E-03 0.000000E+00 9.284599E-01 1.754804E-02 3.811533E-04 6.142411E-02 0.000000E+00 2.679941E-03 2.679326E-03 + 5.960533E-01 -6.967906E-03 5.000000E-01 5.960298E-01 -6.967632E-03 0.000000E+00 9.295134E-01 1.772769E-02 3.856693E-04 6.448824E-02 0.000000E+00 2.655437E-03 2.654827E-03 + 6.292047E-01 -6.257197E-03 5.000000E-01 6.291799E-01 -6.256950E-03 0.000000E+00 9.305617E-01 1.795580E-02 3.901036E-04 6.754712E-02 0.000000E+00 2.632162E-03 2.631559E-03 + 6.602470E-01 -5.594318E-03 5.000000E-01 6.602210E-01 -5.594098E-03 0.000000E+00 9.315928E-01 1.822770E-02 3.939439E-04 7.056301E-02 0.000000E+00 2.610339E-03 2.609749E-03 + 6.890788E-01 -4.983729E-03 5.000000E-01 6.890517E-01 -4.983533E-03 0.000000E+00 9.325991E-01 1.853630E-02 3.977322E-04 7.351174E-02 0.000000E+00 2.589769E-03 2.589193E-03 + 7.156633E-01 -4.425241E-03 5.000000E-01 7.156351E-01 -4.425067E-03 0.000000E+00 9.335855E-01 1.887325E-02 4.018549E-04 7.640634E-02 0.000000E+00 2.570030E-03 2.569468E-03 + 7.400180E-01 -3.918171E-03 5.000000E-01 7.399888E-01 -3.918017E-03 0.000000E+00 9.345871E-01 1.923058E-02 4.056590E-04 7.934534E-02 0.000000E+00 2.550773E-03 2.550227E-03 + 7.622041E-01 -3.461570E-03 5.000000E-01 7.621741E-01 -3.461434E-03 0.000000E+00 9.356022E-01 1.959702E-02 4.094466E-04 8.232463E-02 0.000000E+00 2.531354E-03 2.530823E-03 + 7.823159E-01 -3.051539E-03 5.000000E-01 7.822851E-01 -3.051419E-03 0.000000E+00 9.366271E-01 1.995737E-02 4.138259E-04 8.533382E-02 0.000000E+00 2.511033E-03 2.510515E-03 + 8.004707E-01 -2.683917E-03 5.000000E-01 8.004391E-01 -2.683811E-03 0.000000E+00 9.377133E-01 2.029576E-02 4.167475E-04 8.851532E-02 0.000000E+00 2.489863E-03 2.489360E-03 + 8.168035E-01 -2.357672E-03 5.000000E-01 8.167713E-01 -2.357579E-03 0.000000E+00 9.388784E-01 2.059302E-02 4.142201E-04 9.191304E-02 0.000000E+00 2.468922E-03 2.468444E-03 + 8.314859E-01 -2.072892E-03 5.000000E-01 8.314531E-01 -2.072810E-03 0.000000E+00 9.400241E-01 2.082186E-02 4.039887E-04 9.524550E-02 0.000000E+00 2.449731E-03 2.449290E-03 + 8.447441E-01 -1.827142E-03 5.000000E-01 8.447108E-01 -1.827070E-03 0.000000E+00 9.410212E-01 2.095083E-02 3.870335E-04 9.814534E-02 0.000000E+00 2.433453E-03 2.433057E-03 + 8.568334E-01 -1.615150E-03 5.000000E-01 8.567996E-01 -1.615086E-03 0.000000E+00 9.418312E-01 2.095896E-02 3.658859E-04 1.005031E-01 0.000000E+00 2.420404E-03 2.420054E-03 + 1 2 + 2 3 + 3 4 + 4 5 + 5 6 + 6 7 + 7 8 + 8 9 + 9 10 + 10 11 + 11 12 + 12 13 + 13 14 + 14 15 + 15 16 + 16 17 + 17 18 + 18 19 + 19 20 + 20 21 + 21 22 + 22 23 + 23 24 + 24 25 + 25 26 + 26 27 + 27 28 + 28 29 + 29 30 + 30 31 + 31 32 + 32 33 + 33 34 + 34 35 + 35 36 + 37 38 + 38 39 + 39 40 + 40 41 + 41 42 + 42 43 + 43 44 + 44 45 + 45 46 + 46 47 + 47 48 + 48 49 + 49 50 + 50 51 + 51 52 + 52 53 + 53 54 + 54 55 + 55 56 + 56 57 + 57 58 + 58 59 + 59 60 + 60 61 + 61 62 + 62 63 + 63 64 + 64 65 + 65 66 + 66 67 + 67 68 + 68 69 + 69 70 + 70 71 + 72 73 + 73 74 + 74 75 + 75 76 + 76 77 + 77 78 + 78 79 + 79 80 + 80 81 + 81 82 + 82 83 + 83 84 + 84 85 + 85 86 + 86 87 + 87 88 + 88 89 + 89 90 + 90 91 + 91 92 + 92 93 + 93 94 + 94 95 + 95 96 + 96 97 + 97 98 + 98 99 + 99 100 + 100 101 + 101 102 + 102 103 + 103 104 + 104 105 + 105 106 + 106 107 + 108 109 + 109 110 + 110 111 + 111 112 + 112 113 + 113 114 + 114 115 + 115 116 + 116 117 + 117 118 + 118 119 + 119 120 + 120 121 + 121 122 + 122 123 + 123 124 + 124 125 + 125 126 + 126 127 + 127 128 + 128 129 + 129 130 + 130 131 + 131 132 + 132 133 + 133 134 + 134 135 + 135 136 + 136 137 + 137 138 + 138 139 + 139 140 + 140 141 + 141 142 diff --git a/tests/input/airfoil_000_surf.plt b/tests/input/airfoil_000_surf.plt new file mode 100644 index 0000000..cea1695 Binary files /dev/null and b/tests/input/airfoil_000_surf.plt differ diff --git a/tests/test_tecplotIO.py b/tests/test_tecplotIO.py new file mode 100644 index 0000000..210a0db --- /dev/null +++ b/tests/test_tecplotIO.py @@ -0,0 +1,646 @@ +import itertools +import tempfile +import unittest +from itertools import product +from pathlib import Path +from typing import List, Tuple + +import numpy as np +import numpy.typing as npt +from parameterized import parameterized + +from baseclasses.utils import TecplotFEZone, TecplotOrderedZone, readTecplot, writeTecplot +from baseclasses.utils.tecplotIO import Separator, ZoneType + +# --- Save tempfile locally or in a temp directory --- +SAVE_TEMPFILES = False +SAVE_DIR = Path(__file__).parent / "tmp_tecplot" if SAVE_TEMPFILES else None +if SAVE_DIR is not None: + SAVE_DIR.mkdir(exist_ok=True) + +# --- Define test matrices for Ordered and FE Zones --- +TEST_CASES_ORDERED = [ + tc + for tc in itertools.product( + [(10,), (100,), (10, 10), (100, 100), (10, 10, 10)], + ["SINGLE", "DOUBLE"], + ["POINT", "BLOCK"], + [".dat", ".plt"], + [Separator.SPACE, Separator.COMMA, Separator.TAB, Separator.NEWLINE, Separator.CARRIAGE_RETURN], + ) +] +# Add a single stress test case for each datapacking +TEST_CASES_ORDERED.append(((100, 100, 100), "DOUBLE", "BLOCK", ".dat", Separator.SPACE)) +TEST_CASES_ORDERED.append(((100, 100, 100), "DOUBLE", "POINT", ".dat", Separator.COMMA)) + +# filter out cases that use POINT datapacking with binary '.plt' extensions +TEST_CASES_ORDERED = [tc for tc in TEST_CASES_ORDERED if not (tc[2] == "POINT" and tc[3] == ".plt")] + +TEST_CASES_FE = [ + tc + for tc in itertools.product( + [ZoneType.FELINESEG, ZoneType.FETRIANGLE, ZoneType.FEQUADRILATERAL, ZoneType.FETETRAHEDRON, ZoneType.FEBRICK], + ["SINGLE", "DOUBLE"], + ["POINT", "BLOCK"], + [".dat", ".plt"], + [Separator.SPACE, Separator.COMMA, Separator.TAB, Separator.NEWLINE, Separator.CARRIAGE_RETURN], + ) +] +# filter out cases that use POINT datapacking with binary '.plt' extensions +TEST_CASES_FE = [tc for tc in TEST_CASES_FE if not (tc[2] == "POINT" and tc[3] == ".plt")] + + +def createBrickGrid(ni: int, nj: int, nk: int) -> Tuple[npt.NDArray, npt.NDArray]: + """Create a 3D grid of hexahedral 'brick' elements. + + Parameters + ---------- + ni : int + The number of elements in the i-direction. + nj : int + The number of elements in the j-direction. + nk : int + The number of elements in the k-direction. + + Returns + ------- + Tuple[npt.NDArray, npt.NDArray] + A tuple containing the node coordinates and element connectivity. + """ + x = np.linspace(0, 1, ni) + y = np.linspace(0, 1, nj) + z = np.linspace(0, 1, nk) + x, y, z = np.meshgrid(x, y, z) + nodal_data = np.column_stack((x.flatten(), y.flatten(), z.flatten())) + + connectivity = [] + for i in range(nk - 1): + for j in range(nj - 1): + for k in range(ni - 1): + index = i * (ni * nj) + j * ni + k + connectivity.append( + [ + index, + index + 1, + index + ni + 1, + index + ni, + index + (ni * nj), + index + (ni * nj) + 1, + index + (ni * nj) + ni + 1, + index + (ni * nj) + ni, + ] + ) + + return nodal_data, np.array(connectivity) + + +def createTetGrid(ni: int, nj: int, nk: int) -> Tuple[npt.NDArray, npt.NDArray]: + """Create a 3D grid of tetrahedral elements. + + Parameters + ---------- + ni : int + The number of elements in the i-direction. + nj : int + The number of elements in the j-direction. + nk : int + The number of elements in the k-direction. + + Returns + ------- + Tuple[npt.NDArray, npt.NDArray] + A tuple containing the node coordinates and element connectivity. + """ + # Create node coordinates + x = np.linspace(0, 1, ni + 1) + y = np.linspace(0, 1, nj + 1) + z = np.linspace(0, 1, nk + 1) + xx, yy, zz = np.meshgrid(x, y, z) + nodes = np.column_stack((xx.flatten(), yy.flatten(), zz.flatten())) + + # Create element connectivity + connectivity = [] + for k in range(nk): + for j in range(nj): + for i in range(ni): + # Get the eight corners of the hexahedron + n0 = i + j * (ni + 1) + k * (ni + 1) * (nj + 1) + n1 = n0 + 1 + n2 = n1 + (ni + 1) + n3 = n2 - 1 + n4 = n0 + (ni + 1) * (nj + 1) + n5 = n4 + 1 + n6 = n5 + (ni + 1) + n7 = n6 - 1 + + # Basic subdivision (always add these) + connectivity.extend( + [ + [n0, n1, n3, n7], + [n0, n1, n4, n7], + [n1, n2, n3, n7], + [n1, n2, n6, n7], + [n1, n4, n5, n7], + [n1, n5, n6, n7], + ] + ) + + # Additional tetrahedra for boundary faces + if i == 0: + connectivity.append([n0, n3, n4, n7]) + if i == ni - 1: + connectivity.append([n1, n2, n5, n6]) + if j == 0: + connectivity.append([n0, n1, n4, n5]) + if j == nj - 1: + connectivity.append([n2, n3, n6, n7]) + if k == 0: + connectivity.append([n0, n1, n2, n3]) + if k == nk - 1: + connectivity.append([n4, n5, n6, n7]) + + connectivity = np.array(connectivity) + + return nodes, connectivity + + +# Create a finite element mesh with connectivity +def createQuadGrid(ni: int, nj: int) -> Tuple[npt.NDArray, npt.NDArray]: + """Create a 2D grid of quadrilateral elements. + + Parameters + ---------- + ni : int + The number of elements in the i-direction. + nj : int + The number of elements in the j-direction. + + Returns + ------- + Tuple[npt.NDArray, npt.NDArray] + A tuple containing the node coordinates and element connectivity. + """ + # Create node coordinates + x = np.linspace(0, 1, ni + 1) + y = np.linspace(0, 1, nj + 1) + xx, yy = np.meshgrid(x, y) + nodes = np.column_stack((xx.flatten(), yy.flatten(), np.zeros_like(xx.flatten()))) + + # Create element connectivity + connectivity = [] + for j in range(nj): + for i in range(ni): + n1 = i + j * (ni + 1) + n2 = n1 + 1 + n3 = n2 + ni + 1 + n4 = n3 - 1 + connectivity.append([n1, n2, n3, n4]) + + connectivity = np.array(connectivity) + + return nodes, connectivity + + +def createTriGrid(ni: int, nj: int) -> Tuple[npt.NDArray, npt.NDArray]: + """Create a 2D grid of triangular elements. + + Parameters + ---------- + ni : int + The number of elements in the i-direction. + nj : int + The number of elements in the j-direction. + + Returns + ------- + Tuple[npt.NDArray, npt.NDArray] + A tuple containing the node coordinates and element connectivity. + """ + # Create node coordinates + x = np.linspace(0, 1, ni + 1) + y = np.linspace(0, 1, nj + 1) + xx, yy = np.meshgrid(x, y) + nodes = np.column_stack((xx.flatten(), yy.flatten(), np.zeros_like(xx.flatten()))) + + # Create element connectivity + connectivity = [] + for j in range(nj): + for i in range(ni): + n1 = i + j * (ni + 1) + n2 = n1 + 1 + n3 = n2 + ni + 1 + connectivity.append([n1, n2, n3]) + + n1 = i + j * (ni + 1) + n2 = n1 + ni + 1 + n3 = n2 + 1 + connectivity.append([n1, n2, n3]) + + connectivity = np.array(connectivity) + + return nodes, connectivity + + +def createLineSegGrid(ni: int) -> Tuple[npt.NDArray, npt.NDArray]: + """Create a 1D grid of line segments. + + Parameters + ---------- + ni : int + The number of elements in the i-direction. + + Returns + ------- + Tuple[npt.NDArray, npt.NDArray] + A tuple containing the node coordinates and element connectivity. + """ + # Create node coordinates + x = np.linspace(0, 1, ni + 1) + nodes = np.column_stack((x, x**2, np.zeros_like(x))) + + # Create element connectivity + connectivity = np.column_stack((np.arange(ni), np.arange(1, ni + 1))) + + return nodes, connectivity + + +class TestTecplotIO(unittest.TestCase): + N_PROCS = 1 + + def setUp(self): + thisDir = Path(__file__).parent + self.externalFileAscii = thisDir / "input" / "airfoil_000_slices.dat" + self.externalFileBinary = thisDir / "input" / "airfoil_000_surf.plt" + + def assertOrderedZonesEqual(self, zones1: List[TecplotOrderedZone], zones2: List[TecplotOrderedZone]): + self.assertEqual(len(zones1), len(zones2)) + for zone1, zone2 in zip(zones1, zones2): + self.assertEqual(zone1.name, zone2.name) + self.assertEqual(zone1.shape, zone2.shape) + self.assertListEqual(zone1.variables, zone2.variables) + self.assertEqual(zone1.iMax, zone2.iMax) + self.assertEqual(zone1.jMax, zone2.jMax) + self.assertEqual(zone1.kMax, zone2.kMax) + self.assertEqual(zone1.solutionTime, zone2.solutionTime) + self.assertEqual(zone1.strandID, zone2.strandID) + + def assertFEZonesEqual(self, zones1: List[TecplotFEZone], zones2: List[TecplotFEZone]): + self.assertEqual(len(zones1), len(zones2)) + for zone1, zone2 in zip(zones1, zones2): + self.assertEqual(zone1.name, zone2.name) + self.assertEqual(zone1.shape, zone2.shape) + self.assertListEqual(zone1.variables, zone2.variables) + self.assertEqual(zone1.nElements, zone2.nElements) + self.assertEqual(zone1.nNodes, zone2.nNodes) + self.assertEqual(zone1.solutionTime, zone2.solutionTime) + self.assertEqual(zone1.strandID, zone2.strandID) + + def test_orderedZone(self): + # Create a 3D grid of shape (nx, ny, nz, 3) + nx, ny, nz = 10, 10, 10 + X = np.random.rand(nx, ny, nz) + Y = np.random.rand(nx, ny, nz) + Z = np.random.rand(nx, ny, nz) + + # Create a Tecplot zone + zone = TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z}) + + self.assertEqual(zone.name, "Grid") + self.assertEqual(zone.shape, (nx, ny, nz)) + self.assertEqual(len(zone.variables), 3) + self.assertListEqual(zone.variables, ["X", "Y", "Z"]) + self.assertEqual(zone.iMax, nx) + self.assertEqual(zone.jMax, ny) + self.assertEqual(zone.kMax, nz) + self.assertEqual(zone.solutionTime, 0.0) + self.assertEqual(zone.strandID, -1) + + msg = "Zone name must be a string" + with self.assertRaises(TypeError, msg=msg): + TecplotOrderedZone(123, {"X": X, "Y": Y, "Z": Z}, solutionTime=0.0, strandID=-1) + + msg = "Solution time must be a float" + with self.assertRaises(TypeError): + TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z}, solutionTime="1.0", strandID=-1) + + msg = "Solution time must be greater than or equal to zero" + with self.assertRaises(ValueError): + TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z}, solutionTime=-1.0, strandID=1) + + msg = "Strand ID must be an integer" + with self.assertRaises(TypeError): + TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z}, solutionTime=0.0, strandID="1") + + msg = "Data values must be numpy arrays." + with self.assertRaises(TypeError): + TecplotOrderedZone( + "Grid", {"X": X.tolist(), "Y": Y.tolist(), "Z": Z.tolist()}, solutionTime=0.0, strandID=-1 + ) + + msg = "Data must be a dictionary." + with self.assertRaises(TypeError): + TecplotOrderedZone("Grid", X, solutionTime=0.0, strandID=-1) + + msg = "All variables must have the same shape." + with self.assertRaises(ValueError): + TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z[:-1]}, solutionTime=0.0, strandID=-1) + + def test_FEZone(self): + # Create a quad grid + nx, ny = 10, 10 + nodes, connectivity = createQuadGrid(nx, ny) + + solutionTime = 10.0 + strandID = 4 + # Create a Tecplot zone with Quad elements + zone = TecplotFEZone( + "QuadGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1]}, + connectivity, + zoneType=ZoneType.FEQUADRILATERAL, + solutionTime=solutionTime, + strandID=strandID, + ) + + self.assertEqual(zone.name, "QuadGrid") + self.assertEqual(zone.shape, ((nx + 1) * (ny + 1),)) + self.assertEqual(len(zone.variables), 2) + self.assertListEqual(zone.variables, ["X", "Y"]) + self.assertEqual(zone.nElements, connectivity.shape[0]) + self.assertEqual(zone.nNodes, (nx + 1) * (ny + 1)) + self.assertEqual(zone.solutionTime, solutionTime) + self.assertEqual(zone.strandID, strandID) + self.assertEqual(zone.zoneType, ZoneType.FEQUADRILATERAL) + + nx, ny, nz = 10, 10, 10 + nodes, connectivity = createTetGrid(nx, ny, nz) + + # Create a Tecplot zone with Tet elements + zone = TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + solutionTime=solutionTime, + strandID=strandID, + ) + + self.assertEqual(zone.name, "TetGrid") + self.assertEqual(zone.shape, ((nx + 1) * (ny + 1) * (nz + 1),)) + self.assertEqual(len(zone.variables), 3) + self.assertListEqual(zone.variables, ["X", "Y", "Z"]) + self.assertEqual(zone.nElements, connectivity.shape[0]) + self.assertEqual(zone.nNodes, (nx + 1) * (ny + 1) * (nz + 1)) + self.assertEqual(zone.solutionTime, solutionTime) + self.assertEqual(zone.strandID, strandID) + self.assertEqual(zone.zoneType, ZoneType.FETETRAHEDRON) + + msg = "Zone name must be a string" + with self.assertRaises(TypeError, msg=msg): + TecplotFEZone( + 123, + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + ) + + msg = "Solution time must be a float" + with self.assertRaises(TypeError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + solutionTime="1.0", + ) + + msg = "Solution time must be greater than or equal to zero" + with self.assertRaises(ValueError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + solutionTime=-1.0, + ) + + msg = "Strand ID must be an integer" + with self.assertRaises(TypeError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + strandID="1", + ) + + msg = "Data values must be numpy arrays." + with self.assertRaises(TypeError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0].tolist(), "Y": nodes[:, 1].tolist(), "Z": nodes[:, 2].tolist()}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + ) + + msg = "Data must be a dictionary." + with self.assertRaises(TypeError): + TecplotFEZone("TetGrid", nodes, connectivity, zoneType=ZoneType.FETETRAHEDRON) + + msg = "All variables must have the same shape." + with self.assertRaises(ValueError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:-1, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + ) + + msg = "Connectivity shape does not match zone type." + with self.assertRaises(AssertionError): + TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETRIANGLE, + ) + + @parameterized.expand( + TEST_CASES_ORDERED, + name_func=lambda f, n, p: parameterized.to_safe_name(f"{f.__name__}_{p[0]}"), + ) + def test_ReadWriteOrderedZones( + self, shape: Tuple[int, ...], precision: str, datapacking: str, ext: str, separator: Separator + ): + zones: List[TecplotOrderedZone] = [] + + if len(shape) == 1: + X = np.linspace(0, 1, shape[0]) + Y = np.zeros_like(X) + Z = np.zeros_like(X) + elif len(shape) == 2: + x = np.linspace(0, 1, shape[0]) + y = np.linspace(0, 1, shape[1]) + + X, Y = np.meshgrid(x, y) + Z = np.zeros_like(X) + elif len(shape) == 3: + x = np.linspace(0, 1, shape[0]) + y = np.linspace(0, 1, shape[1]) + z = np.linspace(0, 1, shape[2]) + + X, Y, Z = np.meshgrid(x, y, z, indexing="ij") + + zone = TecplotOrderedZone("Grid", {"X": X, "Y": Y, "Z": Z}) + zones.append(zone) + + title = "ASCII ORDERED ZONES TEST" + prefix = f"ORDERED_{shape}_{precision}_{datapacking}_" + with tempfile.NamedTemporaryFile(suffix=ext, delete=not SAVE_TEMPFILES, dir=SAVE_DIR, prefix=prefix) as tmpfile: + writeTecplot(tmpfile.name, title, zones, datapacking=datapacking, precision=precision, separator=separator) + titleRead, zonesRead = readTecplot(tmpfile.name) + + if ext == ".dat": + with open(tmpfile.name, "r") as f: + lines = f.readlines() + + for line in lines: + self.assertTrue(len(line) <= 4000) + + self.assertEqual(titleRead, title) + self.assertOrderedZonesEqual(zones, zonesRead) + + @parameterized.expand( + TEST_CASES_FE, + name_func=lambda f, n, p: parameterized.to_safe_name(f"{f.__name__}_{p[0]}"), + ) + def test_ReadWriteFEZones( + self, zoneType: ZoneType, precision: str, datapacking: str, ext: str, separator: Separator + ): + zones: List[TecplotFEZone] = [] + + if zoneType == ZoneType.FELINESEG: + rawDataList = [createLineSegGrid(nx) for nx in range(2, 10)] + elif zoneType == ZoneType.FETRIANGLE: + rawDataList = [createTriGrid(nx, ny) for nx, ny in product(range(2, 10), range(2, 10))] + elif zoneType == ZoneType.FEQUADRILATERAL: + rawDataList = [createQuadGrid(nx, ny) for nx, ny in product(range(2, 10), range(2, 10))] + elif zoneType == ZoneType.FETETRAHEDRON: + rawDataList = [createTetGrid(nx, ny, nz) for nx, ny, nz in zip(range(2, 10), range(2, 10), range(2, 10))] + elif zoneType == ZoneType.FEBRICK: + rawDataList = [createBrickGrid(nx, ny, nz) for nx, ny, nz in zip(range(2, 10), range(2, 10), range(2, 10))] + + for i, (nodes, connectivity) in enumerate(rawDataList): + zone = TecplotFEZone( + f"Grid_{i}_{nodes.shape}", + {"X": nodes[..., 0], "Y": nodes[..., 1], "Z": nodes[..., 2]}, + connectivity, + zoneType=zoneType, + ) + zones.append(zone) + + title = f"ASCII {zoneType.name} ZONES TEST" + prefix = f"{zoneType.name}_{precision}_{datapacking}_" + with tempfile.NamedTemporaryFile(suffix=ext, delete=not SAVE_TEMPFILES, dir=SAVE_DIR, prefix=prefix) as tmpfile: + writeTecplot(tmpfile.name, title, zones, datapacking=datapacking, precision=precision, separator=separator) + titleRead, zonesRead = readTecplot(tmpfile.name) + + if ext == ".dat": + with open(tmpfile.name, "r") as f: + lines = f.readlines() + + for line in lines: + self.assertTrue(len(line) <= 4000) + + self.assertEqual(titleRead, title) + self.assertFEZonesEqual(zones, zonesRead) + + def test_ASCIIReadWriteExternal(self): + try: + title, zones = readTecplot(self.externalFileAscii) + + for zone in zones: + conn = zone.connectivity + uniqueConn = zone.uniqueConnectivity + data = zone.data + uniqueData = zone.uniqueData + + for variable in zone.variables: + np.testing.assert_allclose(data[variable][conn], uniqueData[variable][uniqueConn]) + except Exception as e: + self.fail(f"Reading external ASCII file {self.externalFileAscii} failed with error: {e}") + + def test_BinaryReadWriteExternal(self): + try: + title, zones = readTecplot(self.externalFileBinary) + + for zone in zones: + conn = zone.connectivity + uniqueConn = zone.uniqueConnectivity + data = zone.data + uniqueData = zone.uniqueData + + for variable in zone.variables: + np.testing.assert_allclose(data[variable][conn], uniqueData[variable][uniqueConn]) + + except Exception as e: + self.fail(f"Reading external binary file {self.externalFileBinary} failed with error: {e}") + + def test_TriToTriConn(self): + # Create a fe ordered tri zone + ni, nj = 10, 10 + nodes, connectivity = createTriGrid(ni, nj) + zone = TecplotFEZone( + "TriGrid", {"X": nodes[:, 0], "Y": nodes[:, 1]}, connectivity, zoneType=ZoneType.FETRIANGLE + ) + + triConn = zone.triConnectivity + + np.testing.assert_array_equal(triConn, zone.connectivity) + + def test_QuadToTriConn(self): + # Create a fe ordered quad zone + ni, nj = 10, 10 + nodes, connectivity = createQuadGrid(ni, nj) + zone = TecplotFEZone( + "QuadGrid", {"X": nodes[:, 0], "Y": nodes[:, 1]}, connectivity, zoneType=ZoneType.FEQUADRILATERAL + ) + + triConn = zone.triConnectivity + + self.assertEqual(triConn.shape, (zone.nElements * 2, 3)) + + def test_TriConnBadZoneType(self): + # Create an line seg zone + ni = 10 + nodes, connectivity = createLineSegGrid(ni) + zone = TecplotFEZone("LineSeg", {"X": nodes[:, 0], "Y": nodes[:, 1]}, connectivity, zoneType=ZoneType.FELINESEG) + + with self.assertRaises(TypeError): + zone.triConnectivity + + # Create a tet zone + ni, nj, nk = 10, 10, 10 + nodes, connectivity = createTetGrid(ni, nj, nk) + zone = TecplotFEZone( + "TetGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FETETRAHEDRON, + ) + + with self.assertRaises(TypeError): + zone.triConnectivity + + # Create a brick zone + ni, nj, nk = 10, 10, 10 + nodes, connectivity = createBrickGrid(ni, nj, nk) + zone = TecplotFEZone( + "BrickGrid", + {"X": nodes[:, 0], "Y": nodes[:, 1], "Z": nodes[:, 2]}, + connectivity, + zoneType=ZoneType.FEBRICK, + ) + + with self.assertRaises(TypeError): + zone.triConnectivity