pyADflow API
- class adflow.pyADflow.ADFLOW(*args, **kwargs)[source]
Create the ADflow object.
- Parameters:
- commMPI intra comm
The communicator on which to create ADflow. If not given, defaults to MPI.COMM_WORLD.
- optionsdict
The list of options to use with ADflow. This keyword argument is NOT OPTIONAL. It must always be provided. It must contain, at least the ‘gridFile’ entry for the filename of the grid to load
- debugbool
Set this flag to true when debugging with a symbolic debugger. The MExt module deletes the copied .so file when not required which causes issues debugging.
- dtypestr
String type for float: ‘d’ or ‘D’. Not needed to be used by user.
- addActuatorRegion(fileName, axis1, axis2, familyName, thrust=0.0, torque=0.0, heat=0.0, relaxStart=None, relaxEnd=None, coordXfer=None)[source]
Add an actuator disk zone defined by the supplied closed surface in the plot3d file “fileName”. This surface defines the physical extent of the region over which to apply the source terms. The plot3d file may be multi-block but all the surface normals must point outside and no additional surfaces can be inside. Internally, we find all of the CFD volume cells that have their center inside this closed surface and apply the source terms over these cells. This surface is only used with the original mesh coordinates to mark the internal CFD cells, and we keep using these cells even if the geometric design and the mesh coordinates change. For example, the marked cells can lie outside of the original closed surface after a large design change, but we will still use the cells that were inside the surface with the baseline design.
axis1 and axis2 define the vector that we use to determine the direction of the thrust addition. Internally, we compute a vector by $axisVec = axis2-axis1$ and then we normalize this vector. When adding the thrust terms, the direction of the thrust is obtained by multiplying the thrust magnitude by this vector.
Optionally, the source terms in the actuator zone can be gradually ramped up as the solution converges. This continuation approach can be more robust but the users should be careful with the side effects of partially converged solutions. This behavior can be controlled by relaxStart and relaxEnd parameters. By default, the full magnitude of the source terms are added. relaxStart controls when the continuation starts ramping up. The value represents the relative convergence in a log10 basis. So relaxStart = 2 means the source terms will be inactive until the residual is decreased by a factor of 100 compared to free stream conditions. The source terms are ramped to the full magnitude at relaxEnd. E.g., a relaxEnd value of 4 would result in the full source terms after a relative reduction of 1e4 in the total residuals. If relaxStart is not provided, but relaxEnd is provided, then the relaxStart is assumed to be 0. If both are not provided, we do not do ramping and just activate the full source terms from the beginning. When this continuation is used, we internally ramp up the magnitude of the source terms monotonically to prevent oscillations; i.e., decrease in the total residuals increase the source term magnitudes, but an increase in the residuals do not reduce the source terms back down.
- Parameters:
- fileNamestr
Surface Plot 3D file (multiblock ascii) defining the closed region over which the integration is to be applied.
- axis1numpy array, length 3
The physical location of the start of the axis
- axis2numpy array, length 4
The physical location of the end of the axis
- familyNamestr
The name to be associated with the functions defined on this region.
- thrustscalar, float
The total amount of axial force to apply to this region, in the direction of axis1 -> axis2
- torquescalar, float
The total amount of torque to apply to the region, about the specified axis.
- heatscalar, float
The total amount of head added in the actuator zone with source terms
- relaxStartscalar, float
The start of the relaxation in terms of orders of magnitude of relative convergence
- relaxEndscalar, float
The end of the relaxation in terms of orders of magnitude of relative convergence
- coordXferfunction
A callback function that performs a coordinate transformation between the original reference frame and any other reference frame relevant to the current CFD case. This allows user to apply arbitrary modifications to the loaded plot3d surface. The call signature is documented in DVGeometry’s
addPointset
method.
- addArbitrarySlices(normals, points, sliceType='relative', groupName=None, sliceDir=None)[source]
Add slices that vary arbitrarily in space. This is a generalization of the
addSlices
routine, where we have the user specify a list of “normals” and “points” that define slice planes. Rest of the code is the same. This way users can add slices that follow the dihedral of a wing for example.- Parameters:
- normals3-d array or ndarray (nSlice, 3)
The normals of the slice directions. If an array of size 3 is passed, we use the same normal for all slices. if an array of multiple normals are passed, we use the individual normals for each point. in this case, the numbers of points and normals must match.
- pointsndarray (nSlice, 3)
Point coordinates that define a slicing plane along with the normals
- sliceDirndarray (nSlice, 3)
If this is provided, only the intersections that is in this direction starting from the normal is kept. This is useful if you are slicing a closed geometry and only want to keep one side of it. It is similar to the functionality provided in
addCylindricalSlices
.- sliceTypestr {‘relative’, ‘absolute’}
Relative slices are ‘sliced’ at the beginning and then parametrically move as the geometry deforms. As a result, the slice through the geometry may not remain planar. An absolute slice is re-sliced for every output so is always exactly planar and always at the initial position the user indicated.
- groupNamestr
The family to use for the slices. Default is None corresponding to all wall groups.
- addCylindricalSlices(pt1, pt2, nSlice=25, sliceBeg=0.0, sliceEnd=360.0, sliceType='relative', groupName=None)[source]
Add cylindrical projection slices. The cylindrical projection axis is defined from pt1 to pt2. Then we start from the lift direction and rotate around this axis until we are back at the beginning. The rotation direction follows the right hand rule; we rotate the plane in the direction of vectors from pt1 to pt2 crossed with the lift direction. This routine will not work as is if the cylinder axis is exactly aligned with the lift direction. If that use case is needed, consider using
addArbitrarySlices
.- Parameters:
- pt1numpy array, length 3
beginning point of the vector that defines the rotation axis
- pt2numpy array, length 3
end point of the vector that defines the rotation axis
- nSliceint
number of slices around a 360 degree full rotation.
- sliceBegfloat
beginning of the slices in the rotation direction in degrees
- sliceEndfloat
end of the slices in the rotation direction in degrees
- sliceTypestr {‘relative’, ‘absolute’}
Relative slices are ‘sliced’ at the beginning and then parametrically move as the geometry deforms. As a result, the slice through the geometry may not remain planar. An absolute slice is re-sliced for every output so is always exactly planar and always at the initial position the user indicated.
- groupNamestr
The family to use for the slices. Default is None corresponding to all wall groups.
- addFunction(funcName, groupName, name=None)[source]
Add a “new” function to ADflow by restricting the integration of an existing ADflow function by a section of the mesh defined by ‘groupName’. The function will be named ‘funcName_groupName’ provided that the ‘name’ keyword argument is not given. It is is, the function will be use that name. If necessary, this routine may be used to “change the name” of a function. For example,
>>> addFunction('cd', None, 'super_cd')
will add a function that is the same as ‘cd’, but called ‘super_cd’.
- Parameters:
- funcNamestr or list
The name of the built-in adflow function
- groupNamestr or list
The family or group of families to use for the function
- Namestr or list
An overwrite name.
- Returns:
- nameslist
The names if the functions that were added
- addFunctions(funcNames, groupNames, names=None)[source]
Add a series of new functions to ADflow. This is a vector version of the addFunction() routine. See that routine for more documentation.
- addIntegrationSurface(fileName, familyName, isInflow=True, coordXfer=None)[source]
Add a specific integration surface for performing massflow-like computations.
- Parameters:
- fileNamestr
Surface Plot 3D file (may have multiple zones) defining integration surface.
- familyNamestr
User supplied name to use for this family. It should not already be a name that is defined by the surfaces of the CGNS file.
- isInflowbool
Flag to treat momentum forces as if it is an inflow or outflow face. Default is True
- coordXferfunction
A callback function that performs a coordinate transformation between the original reference frame and any other reference frame relevant to the current CFD case. This allows user to apply arbitrary modifications to the loaded plot3d surface. The call signature is documented in DVGeometry’s
addPointset
method.
- addLiftDistribution(nSegments, direction, groupName=None)[source]
Add a lift distribution to the surface output.
- Parameters:
- nSegmentsint
Number of slices to use for the distribution. Typically 150-250 is sufficient
- directionstr {‘x’, ‘y’, ‘z’}
The normal of the slice direction. If you have z-axis ‘out the wing’ you would use ‘z’
- groupName: str
The family group to use for the lift distribution. Default of None corresponds to all wall surfaces.
- addSlices(direction, positions, sliceType='relative', groupName=None)[source]
Add parametric slice positions. Slices are taken of the wing at time addParaSlices() is called and the parametric positions of the intersections on the surface mesh are stored. On subsequent output, the position of the slice moves as the mesh moves/deforms. This effectively tracks the same location on the wing.
- Parameters:
- directionstr {‘x’, ‘y’, ‘z’}
The normal of the slice direction. If you have z-axis ‘out the wing’ you would use ‘z’
- positionsfloat or array
The list of slice positions along the axis given by direction.
- sliceTypestr {‘relative’, ‘absolute’}
Relative slices are ‘sliced’ at the beginning and then parametrically move as the geometry deforms. As a result, the slice through the geometry may not remain planar. An absolute slice is re-sliced for every output so is always exactly planar and always at the initial position the user indicated.
- groupNamestr
The family to use for the slices. Default is None corresponding to all wall groups.
- addUserFunction(funcName, functions, callBack)[source]
Add a new function to ADflow by combining existing functions in a user-supplied way. The allows the user to define a function such as L/D while only requiring a single adjoint solution.
- Parameters:
- funcNamestr
The name of user-supplied function
- functionslist of strings
The exisitng function strings which are the input to the callBack function
- callBackpython function handle
The user supplied function that will produce the user function. This routine must be complex-step safe for the purpose of computing sensitivities.
Examples
>>> ap = AeroProblem(....,evalFuncs=['L/D']) >>> CFDSolver = ADFLOW(options=....) >>> def LoverD(funcs): funcs['L/D'] = funcs['cl']/funcs['cd'] return funcs >>> CFDSolver.addUserFunction('L/D', ['cl','cd'], LoverD)
- checkPartitioning(nprocs)[source]
This function determines the potential load balancing for nprocs. The intent is this function can be run in serial to determine the best number of procs for load balancing. The grid is never actually loaded so this function can be run with VERY large grids without issue.
- Parameters:
- nProcsint
The number of procs check partitioning on
- Returns:
- loadInbalancefloat
The fraction of inbalance for cells. 0 is good. 1.0 is really bad
- faceInbalancefloat
The fraction of inbalance for faces. This gives an idea of communication code. 0 is god. 1.0 is really bad.
- computeJacobianVectorProductBwd(resBar=None, funcsBar=None, fBar=None, wDeriv=None, xVDeriv=None, xSDeriv=None, xDvDeriv=None, xDvDerivAero=None)[source]
This the main python gateway for producing reverse mode jacobian vector products. It is not generally called by the user by rather internally or from another solver. A mesh object must be present for the
xSDeriv=True
flag and a mesh and DVGeo object must be present forxDvDeriv=True
flag. Note that more than one of the specified return flags may be spcified. If more than one return is specified, the order of return is :(wDeriv, xVDeriv, XsDeriv, xDvDeriv, dXdvDerivAero)
.- Parameters:
- resBarnumpy array
Seed for the residuals (dwb in adflow)
- funcsBardict
Dictionary of functions with reverse seeds. Only nonzero seeds need to be provided. All other seeds will be taken as zero.
- fBarnumpy array
Seed for the forces (or tractions depending on the option value) to use in reverse mode.
- wDerivbool
Flag specifiying if the state (w) derivative (wb) should be returned
- xVDerivbool
Flag specifiying if the volume node (xV) derivative should be returned
- xSDerivbool
Flag specifiying if the surface node (xS) derivative should be returned.
- xDvDerivbool
Flag specifiying if the design variable (xDv) derviatives should be returned. This will include both geometric and aerodynamic derivatives
- xDvDerivAerobool
Flag to return just the aerodynamic derivatives. If this is True and xDvDeriv is False,*just* the aerodynamic derivatives are returned.
- Returns:
- wbar, xvbar, xsbar, xdvbar, xdvaerobararray, array, array, dict, dict
One or more of these are returned depending on the
*Deriv
flags provided.
- computeJacobianVectorProductBwdFast(resBar=None)[source]
This fast routine computes only the derivatives of the residuals with respect to the states. This is the operator used for the matrix-free solution of the adjoint system.
- Parameters:
- resBarnumpy array
Seed for the residuals (dwb in adflow)
- Returns:
- wbar: array
state derivative seeds
- computeJacobianVectorProductFwd(xDvDot=None, xSDot=None, xVDot=None, wDot=None, residualDeriv=False, funcDeriv=False, fDeriv=False, groupName=None, mode='AD', h=None, evalFuncs=None)[source]
This the main python gateway for producing forward mode jacobian vector products. It is not generally called by the user by rather internally or from another solver. A DVGeo object and a mesh object must both be set for this routine.
- Parameters:
- xDvDotdict
Perturbation on the geometric design variables defined in DVGeo.
- xSDotnumpy array
Perturbation on the surface
- xVDotnumpy array
Perturbation on the volume
- wDotnumpy array
Perturbation the state variables
- residualDerivbool
Flag specifiying if the residualDerivative (dwDot) should be returned
- funcDerivbool
Flag specifiying if the derviative of the cost functions (as defined in the current aeroproblem) should be returned.
- fderivbool
Flag specifiying if the derviative of the surface forces (tractions) should be returned
- groupNamestr
Optional group name to use for evaluating functions. Defaults to all surfaces.
- modestr [“AD”, “FD”, or “CS”]
Specifies how the jacobian vector products will be computed.
- hfloat
Step sized used when the mode is “FD” or “CS
- Returns:
- dwdot, funcsdot, fDotarray, dict, array
One or more of the these are return depending on the
*Deriv
flags
- computeStabilityParameters()[source]
run the stability derivative driver to compute the stability parameters from the time spectral solution
- evalFunctions(aeroProblem, funcs, evalFuncs=None, ignoreMissing=False)[source]
Evaluate the desired functions given in iterable object, ‘evalFuncs’ and add them to the dictionary ‘funcs’. The keys in the funcs dictionary will be have an
<ap.name>_
prepended to them.- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to to get the solution for
- funcsdict
Dictionary into which the functions are saved.
- evalFuncsiterable object containing strings
If not None, use these functions to evaluate.
- ignoreMissingbool
Flag to suppress checking for a valid function. Please use this option with caution.
- aeroProblem
Examples
>>> funcs = {} >>> CFDsolver(ap) >>> CFDsolver.evalFunctions(ap1, funcs, ['cl', 'cd']) >>> funcs >>> # Result will look like (if aeroProblem, ap1, has name of 'wing'): >>> # {'wing_cl':0.501, 'wing_cd':0.02750}
- evalFunctionsSens(aeroProblem, funcsSens, evalFuncs=None)[source]
Evaluate the sensitivity of the desired functions given in iterable object, ‘evalFuncs’ and add them to the dictionary ‘funcSens’. The keys in the ‘funcsSens’ dictionary will be have an
<ap.name>_
prepended to them.- Parameters:
- funcSensdict
Dictionary into which the function derivatives are saved.
- evalFuncsiterable object containing strings
The additional functions the user wants returned that are not already defined in the aeroProblem
Examples
>>> funcSens = {} >>> CFDsolver.evalFunctionsSens(ap1, funcSens, ['cl', 'cd'])
- getAdjoint(objective)[source]
Return the adjoint values for objective if they exist. Otherwise just return zeros
- getAdjointResNorms()[source]
Return the following adjoint residual norms: initRes Norm: Norm the adjoint RHS startRes Norm: Norm at the start of adjoint call (with possible non-zero restart) finalCFD Norm: Norm at the end of adjoint solve
- getAdjointStateSize()[source]
Return the number of ADJOINT degrees of freedom (states) that are on this processor. The reason this is different from getStateSize() is that if frozenTurbulence is used for RANS, the nonlinear system has 5+neq turb states per cell, while the adjoint still has 5.
- getConvergenceHistory(workUnitTime=None)[source]
Retrieve the convergence history from the fortran level.
This information is printed to the terminal during a run. It is iterTot, IterType, CFL, Step, Linear Res, and CPU time (if added as a monitor variable) and the data for each monitor variable.
- Parameters:
- workUnitTimefloat
The scaling factor specific to the processor. If provided and CPU time is a monitor variable (showcpu is true), the work units (a processor independent time unit) ~will be added to the returned dict too.
- Returns:
- convergeDictdict
A dictionary of arrays and lists. The keys are the data types. The indices of the arrays are the major iteration numbers.
Examples
>>> CFDsolver(ap) >>> CFDsolver.evalFunctions(ap1, funcs, ['cl', 'cd']) >>> hist = CFDSolver.getConvergenceHistory() >>> if MPI.COMM_WORLD.rank == 0: >>> with open(os.path.join(output_dir, "convergence_history.pkl"), "wb") as f: >>> pickle.dump(hist, f)
- getForces(groupName=None, TS=0)[source]
Return the forces on this processor on the families defined by groupName.
- Parameters:
- groupNamestr
Group identifier to get only forces cooresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- TSint
Spectral instance for which to get the forces
- Returns:
- forcesarray (N,3)
Forces (or tractions depending on that forceAsTractions options) on this processor. Note that N may be 0, and an empty array of shape (0, 3) can be returned.
- getHeatFluxes(groupName=None, TS=0)[source]
Return the heat fluxes for isothermal walls on the families defined by group name on this processor.
- Parameters:
- groupNamestr
Group identifier to get only heat fluxes cooresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- TSint
Spectral instance for which to get the fluxes.
- Returns:
- heatFluxesarray (N)
HeatFluxes on this processor. Note that N may be 0, and an empty array of shape (0) can be returned.
- getResNorms()[source]
Return the initial, starting and final Res Norms. Typically used by an external solver.
- getResidual(aeroProblem, res=None, releaseAdjointMemory=True)[source]
Return the residual on this processor. Used in aerostructural analysis
- getSolution(groupName=None)[source]
Retrieve the basic solution variables from the solver. This will return all variables defined in basicCostFunctions for the specified group. This is a lower level function than evalFunctions() which should be used for optimization.
- Parameters:
- groupNamestr
The family group on which to evaluate the functions.
- getSolverMeshIndices()[source]
Get the list of indices to pass to the mesh object for the volume mesh mapping
- getSpatialPerturbation(seed=314)[source]
This is is a debugging routine only. It is used only in regression tests when it is necessary to compute a consistent random spatial vector seed that is independent of per-processor block distribution.
- Parameters:
- seedinteger
Seed to use for random number. Only significant on root processor
- getSpatialSize()[source]
Return the number of degrees of spatial degrees of freedom on this processor. This is (number of nodes)*(number of spectral instances)*3
- getStatePerturbation(seed=314)[source]
This is is a debugging routine only. It is used only in regression tests when it is necessary to compute a consistent random state vector seed that is independent of per-processor block distribution. This routine is not memory scalable as a complete state vector is generated on each process.
- Parameters:
- seedinteger
Seed to use for random number. Only significant on root processor
- getStateSize()[source]
Return the number of degrees of freedom (states) that are on this processor. This is (number of states)*(number of cells)*(number of spectral instances)
- getSurfaceConnectivity(groupName=None, includeZipper=True, includeCGNS=False)[source]
Return the connectivity dinates at which the forces (or tractions) are defined. This is the complement of getForces() which returns the forces at the locations returned in this routine.
- Parameters:
- groupNamestr
Group identifier to get only forces cooresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- includeCGNSbool
Whether or not this function should return the indices of the CGNS blocks that each face patch belongs to. Zipper mesh patches will have cgnsBlockID = -1.
- getSurfacePerturbation(seed=314)[source]
This is is a debugging routine only. It is used only in regression tests when it is necessary to compute a consistent random surface perturbation seed that is independent of per-processor block distribution.
- Parameters:
- seedinteger
Seed to use for random number. Only significant on root processor
- getSurfacePoints(groupName=None, includeZipper=True, TS=0)[source]
Return the coordinates for the surfaces defined by groupName.
- Parameters:
- groupNamestr
Group identifier to get only coordinates cooresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- TSint
The time spectral instance to use for the forces.
- getUniqueSpatialPerturbationNorm(dXv)[source]
This is is a debugging routine only. It is used only in regression tests when it is necessary to compute the norm of a spatial perturbuation on meshes that are split. This will unique-ify the nodes and accumulate onto the unique nodes thus giving the same norm independent of the block splits. Again, this routine is not memory scalable and should only be used for debugging purposes.
- Parameters:
- dXvnumpy vector
Spatial perturbation of size getSpatialSize()
- globalAdjointPreCon(inVec, outVec)[source]
This function is ONLY used as a preconditioner to the global Aero-Structural ADJOINT system. This computes outVec = M^(-1)*inVec where M^(-1) is the approximate inverse application of the preconditing matrix.
- Parameters:
- inVecarray
inVec must be size self.getAdjointStateSize()
- Returns:
- outVecarray
Preconditioned vector
- globalNKPreCon(inVec, outVec)[source]
This function is ONLY used as a preconditioner to the global Aero-Structural system. This computes outVec = M^(-1)*inVec where M^(-1) is the approximate inverse application of the preconditing matrix.
- Parameters:
- inVecarray
inVec must be size self.getStateSize()
- Returns:
- outVecarray
Preconditioned vector
- mapVector(vec1, groupName1, groupName2, vec2=None, includeZipper=True)[source]
This is the main workhorse routine of everything that deals with families in ADflow. The purpose of this routine is to convert a vector ‘vec1’ (of size Nx3) that was evaluated with ‘groupName1’ and expand or contract it (and adjust the ordering) to produce ‘vec2’ evaluated on groupName2.
A little ascii art might help. Consider the following “mesh” . Family ‘fam1’ has 9 points, ‘fam2’ has 10 pts and ‘fam3’ has 5 points. Consider that we have also also added two additional groups: ‘f12’ containing ‘fam1’ and ‘fam2’ and a group ‘f23’ that contains families ‘fam2’ and ‘fam3’. The vector we want to map is ‘vec1’. It is length 9+10. All the ‘x’s are significant values.
The call: mapVector(vec1, ‘f12’, ‘f23’)
will produce the “returned vec” array, containing the significant values from ‘fam2’, where the two groups overlap, and the new values from ‘fam3’ set to zero. The values from fam1 are lost. The returned vec has size 15.
fam1 fam2 fam3 |---------+----------+------| |xxxxxxxxx xxxxxxxxxx| <- vec1 |xxxxxxxxxx 000000| <- returned vec (vec2)
It is also possible to pass in vec2 into this routine. For that case, the existing values in the array will not be kept. In the previous examples, the values corresponding to fam3 will retain their original values.
- Parameters:
- vec1Numpy array
Array of size Nx3 that will be mapped to a different family set.
- groupName1str
The family group where the vector vec1 is currently defined
- groupName2str
The family group where we want to the vector to mapped into
- vec2Numpy array or None
Array containing existing values in the output vector we want to keep. If this vector is not given, the values will be filled with zeros.
- Returns:
- vec2Numpy array
The input vector maped to the families defined in groupName2.
- propagateUncertainty(aeroProblem, evalFuncs=None, UQDict=None)[source]
Use the first order second moment method to predict output uncertainties for the current solution.
- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to solve
- evalFuncsiterable object containing strings
The functions the user wants the uncertainty of
- UQDictdict
Dictionary containing the mean and std dev. of the input parameters that are providing the uncertain input.
- aeroProblem
- resetANKCFL(aeroProblem)[source]
Resets the ANK CFL of the aeroProblem to the value set by “ANKCFL0” option. During the first ANK iteration that follows, this will be overwritten by the ANK solver itself with the CFL Min logic based on relative convergence. This is useful if keeping a really high ANK CFL is not desired.
- Parameters:
- aeroProblem
AeroProblem
The aeroproblem whose ANK CFL will be reset.
- aeroProblem
- resetAdjoint(obj)[source]
Reset an adjoint ‘obj’ that is stored in the current aeroProblem. If the adjoint does not yet exist, nothing is done
- Parameters:
- objstr
String identifing the objective.
- resetFlow(aeroProblem, releaseAdjointMemory=True)[source]
Reset the flow after a failure or for a complex step derivative operation.
- Parameters:
- aeroProblem
AeroProblem
The aeroproblem with the flow information we would like to reset the flow to.
- aeroProblem
- rootSetOption(name, value, reset=False)[source]
Set ADflow options from the root proc. The option will be broadcast to all procs when __call__ is invoked.
- Parameters:
- namestr
The name of the option
- valueAny
The value of the option
- resetBool
If True, we reset all previously-set rootChangedOptions.
See also
- saveAdjointMatrix(baseFileName)[source]
Save the adjoint matrix to a binary petsc file for possible future external testing
- Parameters:
- basefileNamestr
Filename to use. The Adjoint matrix, PC matrix(if it exists) and RHS will be written
- setAdjoint(adjoint, objective=None)[source]
Sets the adjoint vector externally. Used in coupled solver
- setAeroProblem(aeroProblem, releaseAdjointMemory=True)[source]
Set the supplied aeroProblem to be used in ADflow.
- Parameters:
- aeroProblem
AeroProblem
The supplied aeroproblem to be set.
- releaseAdjointMemorybool, optional
Flag to release the adjoint memory when setting a new aeroproblem, by default True
- aeroProblem
- setDisplacements(aeroProblem, dispFile)[source]
This function allows the user to perform aerodynamic analysis/optimization while using a fixed set of displacements computed from a previous structural analysis. Essentially this allows the jig shape to designed, but performing analysis on the flying shape.
- Parameters:
- aeroProblem
AeroProblem
The AP object that the displacements should be applied to.
- dispFilestr
The file contaning the displacements. This file should have been obtained from TACS
- aeroProblem
Notes
The fixed set of displacements do not affect the sensitivities, since they are fixed and not affected by any DVs.
Also, in the case where the current surface mesh was not used to generate the displacements file, a nearest neighbor search is used to apply the displacements.
- setMesh(mesh)[source]
Set the mesh object to the aero_solver to do geometric deformations
- Parameters:
- meshMBMesh or USMesh object
The mesh object for doing the warping
- setResNorms(initNorm=None, startNorm=None, finalNorm=None)[source]
Set one of these norms if not None. Typlically used by an external solver
- setRotationRate(rotCenter, rotRate, cgnsBlocks=None)[source]
Set the rotational rate for the grid:
rotCenter: 3-vectorThe center of rotation for steady motion rotRate: 3-vector or aeroProblem: If it is a 3-vector, take the rotations to about x-y-z, if it is an aeroProblem with p,q,r, defined, use that to compute rotations. cgnsBlocks: The list of blocks to set. NOTE: This must be in 1-based ordering!
- setStates(states)[source]
Set the states on this processor. Used in aerostructural analysis and for switching aeroproblems
- setSurfaceCoordinates(coordinates, groupName=None)[source]
Set the updated surface coordinates for a particular group.
- Parameters:
- coordinatesnumpy array
Numpy array of size Nx3, where N is the number of coordinates on this processor. This array must have the same shape as the array obtained with getSurfaceCoordinates()
- groupNamestr
Name of family or group of families for which to return coordinates for.
- setTargetCp(CpTargets, groupName=None, TS=0)[source]
Set the CpTarget distribution for am inverse design problem.
- Parameters:
- CpSurfnumpy array
Array of CP targets to set for the surface. This size must correpsond to the size of the surface obtained using the same groupName.
- groupNamestr
Group identifier to set only CPtargets corresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- TSint
Time spectral instance to set.
- setWallTemperature(temperature, groupName=None, TS=0)[source]
Set the temperature of the isothermal walls.
- Parameters:
- temperaturenumpy array
Dimensional temperature to set for wall. This size must correpsond to the size of the heat flux obtained using the same groupName.
- groupNamestr
Group identifier to set only temperatures corresponding to the desired group. The group must be a family or a user-supplied group of families. The default is None which corresponds to all wall-type surfaces.
- TSint
Time spectral instance to set.
- solveAdjointForRHS(inVec, relTol=None)[source]
Solve the adjoint system with an arbitary RHS vector.
- Parameters:
- inVecnumpy array
Array of size w
- Returns:
- outVecnumpy array
Solution vector of size w
- solveCL(aeroProblem, CLStar, alpha0=None, CLalphaGuess=None, delta=0.5, tol=0.001, autoReset=False, maxIter=20, relaxCLa=1.0, relaxAlpha=1.0, L2ConvRel=None, stopOnStall=False, writeSolution=False, workUnitTime=None, updateCutoff=1e-16)[source]
This is a simple secant method search for solving for a fixed CL. This really should only be used to determine the starting alpha for a lift constraint in an optimization.
- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to solve
- CLStarfloat
The desired target CL
- alpha0angle (deg)
Initial guess for secant search (deg). If None, use the value in the aeroProblem
- CLalphaGuessfloat or None
The user can provide an estimate for the lift curve slope in order to accelerate convergence. If the user supply a value to this option, it will not use the delta value anymore to select the angle of attack of the second run. The value should be in 1/deg.
- deltaangle (deg)
Initial step direction for secant search
- tolfloat
Desired tolerance for CL
- autoResetbool
Flag to reset flow between solves. The Euler NK method has issues when only the alpha is changed (Martois effect we think). This will reset the flow after each solve which solves this problem. Not necessary (or desired) when using the RK solver. The useCorrection option is the preferred way of addressing this issue with the NK/ANK methods. If solver still struggles after restarts with the correction, this option can be used.
- maxIterint
Maximum number of secant solver iterations.
- relaxCLafloat
Relaxation factor for the CLalpha update. Value of 1.0 will give the standard secant solver. A lower value will damp the update to the CL alpha. Useful when the CL output is expected to be noisy.
- relaxAlphafloat
Relaxation factor for the alpha update. Value of 1.0 will give the standard secant solver. A lower value will damp the alpha update. Useful to prevent large changes to alpha.
- L2ConvRelfloat or list or None
Temporary relative L2 convergence for each iteration of the CL solver. If the option is set to None, we don’t modify the L2ConvergenceRel option. If a float is provided, we use this value as the target relative L2 convergence for each iteration. If a list of floats is provided, we iterate over the values and set the relative L2 convergence target separately for each iteration. If the list length is less than the maximum number of iterations we can perform, we keep using the last value of the list until we run out of iterations.
- stopOnStallbool
Flag to determine if we want to stop early if the solver fails to achieve the relative or total L2 convergence. This is useful when the solver is expected to converge to the target L2 tolerances at each call.
- writeSolutionbool
Flag to enable a call to self.writeSolution as the CL solver is finished.
- workUnitTimefloat or None
Optional parameter that will be passed to getConvergenceHistory after each adflow call.
- updateCutofffloat
Required change in alpha to trigger the clalpha update with the Secant method. If the change in alpha is smaller than this option we don’t even bother changing clalpha and keep using the clalpha value from the last iteration. This prevents clalpha from getting bad updates due to noisy cl outputs that result from small alpha changes.
- aeroProblem
- Returns:
- resultsDictdictionary
Dictionary that contains various results from the cl solve. The dictionary contains the following data:
- convergedbool
Flag indicating cl solver convergence. True means both the CL solver target tolerance AND the overall target L2 convergence is achieved. False means either one or both of the tolerances are not reached.
- iterationsint
Number of iterations ran.
- l2convergencefloat
Final relative L2 convergence of the solver w.r.t. free stream residual. Useful to check overall CFD convergence.
- alphafloat
Final angle of attack value
- clfloat
Final CL value.
- clstarfloat
The original target CL. returned for convenience.
- errorfloat
Error in cl value.
- clalphafloat
Estimate to clalpha used in the last solver iteration
- timefloat
Total time the solver needed, in seconds
- historylist
List of solver convergence histories. Each entry in the list contains the convergence history of the solver from each CL Solve iteration. The history is returned using the “getConvergenceHistory” method.
- solveDirectForRHS(inVec, relTol=None)[source]
Solve the direct system with an arbitary RHS vector.
- Parameters:
- inVecnumpy array
Array of size w
- Returns:
- outVecnumpy array
Solution vector of size w
- solveErrorEstimate(aeroProblem, funcError, evalFuncs=None)[source]
Evaluate the desired function errors given in iterable object ‘evalFuncs’, and add them to the dictionary ‘funcError’. The keys in the funcError dictionary will be have an
<ap.name>_
prepended to them.- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to get the error for
- funcErrordict
Dictionary into which the function errors are saved. We define error to be \(\epsilon = f^\ast - f\), where \(f^\ast\) is the converged solution and \(f\) is the unconverged solution.
- evalFuncsiterable object containing strings
If not None, use these functions to evaluate.
- aeroProblem
Examples
>>> CFDsolver(ap) >>> funcsSens = {} >>> CFDSolver.evalFunctionsSens(ap, funcsSens) >>> funcError = {} >>> CFDsolver.solveErrorEstimate(ap, funcError) >>> # Result will look like (if aeroProblem, ap, has name of 'wing'): >>> print(funcError) >>> # {'wing_cl':0.00085, 'wing_cd':0.000021}
- solveSep(aeroProblem, sepStar, nIter=10, alpha0=None, delta=0.1, tol=0.001, expansionRatio=1.2, sepName=None)[source]
This is a safe-guarded secant search method to determine the alpha that yields a specified value of the separation sensor. Since this function is highly nonlinear we use a linear search to get the bounding range first.
- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to solve
- sepStarfloat
The desired target separation sensor value
- nIterint
Maximum number of iterations
- alpha0angle (deg) or None
Initial guess. If none, use what is in aeroProblem.
- deltaangle (deg)
Initial step. Only the magnitude is significant
- tolfloat
Desired tolerance for sepSensor
- sepNamestr or None
User supplied function to use for sep sensor. May be a user-added group function.
- aeroProblem
- Returns:
- None, but the correct alpha is stored in the aeroProblem
- solveTargetFuncs(aeroProblem, funcDict, tol=0.0001, nIter=10, Jac0=None)[source]
Solve the an arbitrary set of function-dv sets using a Broyden method.
- Parameters:
- aeroProblem
AeroProblem
The aerodynamic problem to be solved
- funcDictdict
Dictionary of function DV pairs to solve: {‘func’:{‘dv’:str, ‘dvIdx’:idx,’target’:val,’initVal’:val,’initStep’:val}} func : Name of function that is being solved for dv : design variable that has dominant control of this function value dvIdx : index into dv array if dv is not scalar target : target function value initVal : initial design variable value initStep : initial step for this dv in when generating finite difference starting jacobian
- tolfloat
Tolerance for the L2 norm of function error from target values
- nIterint
Maximum number of iterations.
- Jac0nxn numpy array
Initial guess for the func-dv Jacobian. Usually obtained from a previous analysis and saves n function evaluations to produce the initial Jacobian.
- aeroProblem
- solveTrimCL(aeroProblem, trimFunc, trimDV, dvIndex, CLStar, trimStar=0.0, alpha0=None, trim0=None, da=0.001, deta=0.01, tol=0.0001, nIter=10, Jac0=None, liftFunc='cl')[source]
Solve the trim-Cl problem using a Broyden method.
- Parameters:
- ASProblemASProblem instance
The aerostructural problem to be solved
- trimFuncstr
Solution variable to use for trim. Usually ‘cmy’ or ‘cmz’
- trimDVstr
Dame of DVGeo design variable to control trim
- dvIndexint
Index of the trimDV function to use for trim
- CLStarfloat
Desired CL value
- trimStarfloat
Desired trimFunc value
- alpha0float or None
Starting alpha. If None, use what is in the aeroProblem
- trim0float or None
Starting trim value. If None, use what is in the DVGeo object
- dafloat
Initial alpha step for Jacobian
- detafloat
Initial stet in the ‘eta’ or trim dv function
- tolfloat
Tolerance for trimCL solve solution
- nIterint
Maximum number of iterations.
- Jac02x2 numpy array
Initial guess for the trim-cl Jacobian. Usually obtained from a previous analysis and saves two function evaluations to produce the initial Jacobian.
- liftFuncstr
Solution variable to use for lift. Usually ‘cl’ or a custom function created from cl.
- writeActuatorRegions(fileName, outputDir=None)[source]
Debug method that writes the cells included in actuator regions to a tecplot file. This routine should be called on all of the procs in self.comm. The output can then be used to verify that the actuator zones are defined correctly.
- Parameters:
- fileNamestr
Name of the output file
- outputDirstr
output directory. If not provided, defaults to the output directory defined with the aero_options
- writeForceFile(fileName, TS=0, groupName=None, cfdForcePts=None)[source]
This function collects all the forces and locations and writes them to a file with each line having: X Y Z Fx Fy Fz. This can then be used to set a set of structural loads in TACS for structural only optimization
Like the getForces() routine, an external set of forces may be passed in on which to evaluate the forces. This is only typically used in an aerostructural case.
- writeLiftDistributionFile(fileName)[source]
Evaluate and write the lift distibution to a tecplot file.
- Parameters:
- fileNamestr
File of lift distribution. Should have .dat extension.
- writeMeshFile(fileName)[source]
Write the current mesh to a CGNS file. This call isn’t used normally since the volume solution usually contains the grid
- Parameters:
- fileNamestr
Name of the mesh file
- writeSlicesFile(fileName)[source]
Evaluate and write the defined slice information to a tecplot file.
- Parameters:
- fileNamestr
Slice file. Should have .dat extension.
- writeSolution(outputDir=None, baseName=None, number=None, writeSlices=True, writeLift=True)[source]
This is a generic shell function that potentially writes the various output files. The intent is that the user or calling program can call this file and ADflow write all the files that the user has defined. It is recommended that this function is used along with the associated logical flags in the options to determine the desired writing procedure
- Parameters:
- outputDirstr
Use the supplied output directory
- baseNamestr
Use this supplied string for the base filename. Typically only used from an external solver.
- numberint
Use the user supplied number to index solution. Again, only typically used from an external solver.
- writeSlicesbool
Flag to determine if the slice files are written, if we have any slices added.
- writeLiftbool
Flag to determine if the lift files are written, if we have any lift distributions added.
- writeSurfaceSensitivity(fileName, func, groupName=None)[source]
Write a tecplot file of the surface sensitivity. It is up to the use to make sure the adjoint already computed before calling this function.
- Parameters:
- fileNamestr
String for output filename. Should end in .dat for tecplot.
- funcstr
ADFlow objective string for the objective to write.
- groupNamestr
Family group to use. Default to all walls if not given (None)
- writeSurfaceSolutionFile(fileName)[source]
Write the current state of the surface flow solution to a CGNS file. Keyword arguments:
- Parameters:
- fileNamestr
Name of the file. Should have .cgns extension.
- writeSurfaceSolutionFileTecplot(fileName)[source]
Write the current state of the surface flow solution to a teclot file.
- Parameters:
- fileNamestr
Name of the file. Should have .plt extension.
- writeVolumeSolutionFile(fileName, writeGrid=True)[source]
Write the current state of the volume flow solution to a CGNS file. This is a lower level routine; Normally one should call
writeSolution()
.- Parameters:
- fileNamestr
Name of the file. Should have .cgns extension.
- writeGridbool
Flag specifying whether the grid should be included or if links should be used. Always writing the grid is recommended even in cases when it is not strictly necessary. Note that if
writeGrid = False
the volume files do not contain any grid coordinates rendering the file useless if a separate grid file was written out and is linked to it.