pyBlock

class pygeo.pyBlock.pyBlock(initType, fileName=None, FFD=False, symmPlane=None, kmax=4, **kwargs)[source]

pyBlock represents a collection of pySpline Volume objects.

It performs several functions including fitting and point projections. The actual b-spline volumes are of the pySpline Volume type.

Parameters
initTypestr

Initialization type. Only ‘plot3d’ is currently available.

fileNamestr

Filename of the plot3d file to be loaded. Should have a .fmt or .xyz extension. Also must be in ASCII format.

FFDbool

Flag to indicate that this object is to be created as an FFD. When this is true, no fitting is performed; The coordinates in the plot 3d file explicitly become the control points and uniform (and symmetric) knot vectors are assumed everywhere. This ensures a seamless FFD.

symmPlane{“x”, “y”, or “z”}

if a coordinate direciton is provided, the code will duplicate the FFD in the mirroring direction.

kmaxint

maximum order of the splines used for the underlying formulation. Default is a 4th order spline in each direction if the dimensions allow.

attachPoints(coordinates, ptSetName, interiorOnly=False, embTol=1e-10, nIter=100, eps=1e-12)[source]

Embed a set of coordinates into the volumes. This is the main high level function that is used by DVGeometry when pyBlock is used as an FFD.

Parameters
coordinatesarray, size (N,3)

The coordinates to embed in the object

ptSetNamestr

The name given to this set of coordinates.

interiorOnlybool

Project only points that lie fully inside the volume

embTolfloat

Tolerance on the distance between projected and closest point. Determines if a point is embedded or not in the FFD volume if interiorOnly is True.

epsfloat

Physical tolerance to which to converge Newton search

nIterint

Maximum number of Newton iterations to perform. The default of 100 should be sufficient for points that actually lie inside the volume, except for pathological or degenerate FFD volumes.

calcdPtdCoef(ptSetName)[source]

Calculate the (fixed) derivative of a set of embedded points with respect to the b-spline coefficients. This derivative consists of the b-spline basis functions

Parameters
ptSetNamestr

The name of the point set to use.

doConnectivity(fileName=None, nodeTol=0.0001, edgeTol=0.0001, greedyReorder=False)[source]

This function is used if a separate fitting topology is required for non-FFD creations. The sequence of calls is given in the examples section.

Parameters
fileNamestr

If the filename exists, read in the topology. Otherwise, write a an initial default topology

nodeTolfloat

Tolerance for co-incidient nodes

edgeTolfloat

Tolerance for co-incidient mid points of edges

greedyReorderbool

Flag to reorder numbering in a greedy form.

fitGlobal(greedyReorder=False)[source]

Determine the set of b-spline coefficients that best fits the set of volumes in the global sense. This is required for non-FFD creation.

Parameters
greedyReorderbool

Flag to compute ordering of initial mesh in a greedy ordering sense.

getAttachedPoints(ptSetName)[source]

Return all the volume points for an embedded volume with name ptSetName.

Parameters
ptSetNamestr

Name of a point set added with attachPoints()

Returns
coordinatesnumpy array (Nx3)

The coordinates of the embedded points. If a mask was used, only the points corresponding to the indices in mask will be non-zero in the array.

getBounds()[source]

Determine the extents of the set of volumes

Returns
xMinarray of length 3

Lower corner of the bounding box

xMaxarray of length 3

Upper corner of the bounding box

printConnectivity()[source]

Print the connectivity information to the screen

projectPoints(x0, checkErrors, embTol, eps, nIter)[source]

Project a set of points x0, into any one of the volumes. It returns the the volume ID, u, v, w, D of the point in volID or closest to it.

This is still technically a inefficient brute force search, but it uses some heuristics to give a much more efficient algorithm. Basically, we use the volume the last point was projected in as a ‘good guess’ as to what volume the current point falls in. This works since subsequent points are usually close together. This will not help for randomly distributed points.

Parameters
x0array of points (Nx3 array)

The list or array of points to use

checkErrorsbool

Flag to print out the error is points have not been projected to the tolerance defined by embTol.

See also

attachPoints

description of the other parameters

writePlot3d(fileName)[source]

Write the grid to a plot3d file. This isn’t efficient as it used ASCII format. Only useful for quick visualizations

Parameters
fileNameplot3d file name.

Should end in .xyz

writePlot3dCoef(fileName)[source]

Write the coefficients of the volumes to a plot3d file.

Parameters
fileNameplot3d file name.

Should end in .fmt

writeTecplot(fileName, vols=True, coef=True, orig=False, volLabels=False, edgeLabels=False, nodeLabels=False)[source]

Write a tecplot visualization of the pyBlock object.

Parameters
fileNamestr

Filename of tecplot file. Should have a .dat extension

volsbool. Default is True

Flag to write interpolated volumes

coefbool. Default is True

Flag to write spline control points

origbool. Default is True

Flag to write original data (if it exists)

volLabels: bool. Default is True

Flag to write volume labels in a separate tecplot file; filename is derived from the supplied fileName.

edgeLabels: bool. Default is False

Flag to write edge labels in a separate tecplot file; filename is derived from the supplied fileName.

nodeLabels: bool. Default is False

Flag to write node labels in a separate tecplot file; filename is derived from the supplied fileName.