Demonstration Functions

To aid in exploration we include a variety of test functions that can be used in conduction with the parameter space dimension reduction techniques. We divide these into two classes: those that are simple formulas and those that are the results of the numerical approximation of differential equations. The former are inexpensive to evaluate and equipped with analytic gradients; the latter can be expensive, may not have gradients, and may contain computational noise.

We provide two different ways to access these demo functions: a psdr.Function() interface working on the normalized domain and a low-level access to the underlying function and domain methods.

Formula-based Test Functions

Borehole Function

class psdr.demos.Borehole(domain='deterministic', dask_client=None)

The borehole test function

A function implementing the borehole test function [VLSE_borehole]. This function has the form

\[f(r_w, r, T_u, H_u, T_l, H_l, L, K_w) = \frac{ 2\pi T_u (H_u - H_l)}{ \ln(r/r_w) \left( 1 + \frac{2 L T_u}{\ln(r/r_w) r_w^2 K_w} + \frac{T_u}{T_l} \right) }\]

where the input variables have the domain

Variable	Interpretation
\(r_w\in [0.05, 0.15]\)	radius of borehole (m)
\(r \in [100, 50 \times 10^3]\)	radius of influence (m)
\(T_u \in [63070,115600]\)	trasmissivity of upper aquifer (m^2/yr)
\(H_u \in [ 990, 1110]\)	potentiometric head of upper aquifer (m)
\(T_l \in [63.1, 116]\)	transmissivity of lower aquifer (m^2/yr)
\(H_l \in [700, 820]\)	potentiometric head of lower aquifer (m)
\(L \in [1120, 1680]\)	length of borehole (m)
\(K_w \in [9855, 12045]\)	hydraulic conductivity of borehole (m/yr)

An alternative to this deterministic domain is an uncertain domain where \(r_w \sim \mathcal{N}(0.10, 0.0161812)\) and \(\log r \sim \mathcal{N}(7.71, 1.0056)\) and the remainder come from a uniform distribution on the domain previously specified.

Parameters:

• domain (['deterministic', 'uncertain']) – Which domain to use when constructing the function
• dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[VLSE_borehole]

Virtual Library of Simulation Experiments, Borehole Function https://www.sfu.ca/~ssurjano/borehole.html

Golinski Gearbox

class psdr.demos.GolinskiGearbox(dask_client=None)

The Golinski Gearbox Optimization Test Problem

This test problem originally descibed by Golinski [Gol70] and subsequently appearing in, for example, [Ray03] and [MDO], seeks to design a gearbox (speedreducer) to minimize volume subject to a number of constraints.

Here we take our formulation following [Ray03], which reduces the original 25 constraints to 11 by removing redundancy. We also shift these constraints such that they are satisfied if the values returned are negative and the constraints are violated if their value is positive.

Parameters:dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[MDO]

Langley Research Center: Multidisciplinary Optimization Test Suite, Golinski’s Speed Reducer

[Gol70]

“Optimal Synthesis Problems Solved by Means of Nonlinear Programming and Random Methods”, Jan Golinski, J. Mech. 5, 1970, pp.287–309. https://doi.org/10.1016/0022-2569(70)90064-9

[Ray03]

(1, 2) “Golinski’s Speed Reducer Problem Revisited”, Tapabrata Ray, AIAA Journal, 41(3), 2003, pp 556–558 https://doi.org/10.2514/2.1984

Nowacki Beam

class psdr.demos.NowackiBeam

The Nowacki Beam Problem

In this example, we seek to design a beam with fixed cross-sectional breath and height with respect to two objectives subject to five constraints described in [FSK08].

The two parameters for this model and their ranges are given below.

Variable	Interpretation
\(b \in [0.005, 0.050]\)	breadth (m)
\(h \in [0.020, 0.250]\)	height (m)

The following are the seven quantities of interest associated with this problem.

Function name	Formula	Role
cross-sectional area	\(a = bh\)	objective: minimize
bending stress	\(\sigma = 6FL/(bh^2)\)	objective: minimize
tip deflection	\(\delta = FL^3/(3EI_Y) \quad I_Y = bh^3/12\)	constraint: \(\delta \le 0.005\)
bending stress	\(\sigma_B = 6FL/(bh^2)\)	constraint: \(\sigma_B \le \sigma_Y\)
shear stress	\(\tau = 3F/(2bh)\)	constraint: \(\tau \le \sigma_Y/2\)
height to breadth ratio	\(\rho = h/b\)	constraint: \(\rho \le 10\)
tip force for buckling	\(F_T = (4/L^2)\sqrt{G I_T E I_Z/(1-\nu)}\)	constraint: \(F_T \ge 2F\)

where \(I_T = (b^3h+bh^3)/12\) and \(I_Z = b^3h/12\). The constants that appear above are taken from those values for mild steel:

Constant	Value
tip force	\(F=5\times 10^3\) N
beam length	\(L=1.5\) m
yield stress	\(\sigma_Y = 240\times 10^6\) Pa
Young’s modulus	\(E=216.62\times 10^9\) Pa
Poisson’s ratio	\(\nu=0.27\)
shear modulus	\(G=86.65\times 10^9\) Pa

This example originates in [Now80].

[FSK08]

Engineering Design via Surrogate Modeling: A Practical Guide Aleander I. J. Forrester, Andras Sobester, and Andy J. Keane Wiley, 2008

[Now80]

Modling of Design Decisions for CAD Horst Nowacki In: Computer Aided Design Modeling, J. Encarncao (editor) DOI:10.1007/BFb0040161

OTL Circuit Function

class psdr.demos.OTLCircuit(dask_client=None)

The OTL circuit test function

The OTL Circuit function “models an output transformerless push pull circuit” [VLSE]

\[\begin{split}f(R_{b1}, R_{b2}, R_f, R_{c1}, R_{c2}, \beta) :=& \frac{ \beta (V_{b1} + 0.74)(R_{c2} + 9)}{\beta(R_{c2} + 9) + R_f} + \frac{11.35 R_f}{\beta(R_{c2} + 9) + R_f} \\ &+\frac{0.74 R_f \beta (R_{c2} + 9)}{R_{c1}(\beta(R_{c2} + 9) + R_f)}, \\ \text{where } V_{b1} :=& \frac{12 R_{b2}}{R_{b1} + R_{b2}}\end{split}\]

Variable	Interpretation
\(R_{b1} \in [50, 150]\)	resistance b1 (K-Ohms)
\(R_{b2} \in [25, 75]\)	resistance b2 (K-Ohms)
\(R_{f} \in [0.5, 3]\)	resistance f (K-Ohms)
\(R_{c1} \in [1.2, 2.5]\)	resistance c1 (K-Ohms)
\(R_{c1} \in [0.25, 1.2]\)	resistance c2 (K-Ohms)
\(\beta \in [50, 300]\)	current gain (Amperes)

Parameters:dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[VLSE]

Virtual Library of Simulation Experiments, OTL Circuit https://www.sfu.ca/~ssurjano/otlcircuit.html

Piston Function

class psdr.demos.Piston(dask_client=None)

Piston test function

The Piston function “models the circular motion of a piston within a cylinder” [VLSE_piston].

\[\begin{split}f(M, S, V_0, k, P_0, T_a, T_0) :=& 2\pi \sqrt{\frac{M}{k+ S^2 \frac{P_0V_0}{T_0} \frac{T_a}{V^2}}}, \text{ where } \\ V &= \frac{S}{2k}\left( \sqrt{A^2 + 4k \frac{P_0V_0}{T_0}T_a} - A \right) \\ A &= P_0 S + 19.62M - \frac{k V_0}{S}\end{split}\]

Variable	Interpretation
\(M \in [30, 60]\)	piston weight (kg)
\(S \in [0.005, 0.020]\)	piston surface area (m^2)
\(V_0 \in [0.002, 0.010]\)	initial gas volume (m^3)
\(k \in [1000, 5000]\)	spring coefficient (N/m)
\(P_0 \in [90\times 10^3, 110 \times 10^3]\)	atmospheric pressure (N/m^2)
\(T_a \in [290, 296]\)	ambient temperature (K)
\(T_0 \in [340, 360]\)	filling gas temperature (K)

Parameters:dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[VLSE_piston]

Virtual Library of Simulation Experiments, Piston Function https://www.sfu.ca/~ssurjano/piston.html

Robot Arm Function

class psdr.demos.RobotArm(dask_client=None)

Robot Arm test function

A test function that measures the distance of a four segment arm from the origin [VLSE_robot].

\[\begin{split}f(\theta_1, \theta_2, \theta_3, \theta_4, L_1, L_2, L_3, L_4) &:= \sqrt{u^2 + v^2}, \text{ where }\\ u &:= \sum_{i=1}^4 L_i \cos \left( \sum_{j=1}^i \theta_j\right) \\ v &:= \sum_{i=1}^4 L_i \sin \left( \sum_{j=1}^i \theta_j\right)\end{split}\]

Variable	Interpretation
\(\theta_i \in [0, 2\pi]\)	angle of the ith arm segment
\(L_i \in [0,1]\)	length of the ith arm segment

Parameters:dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[VLSE_robot]

Virtual Library of Simulation Experiments, Robot Arm Function https://www.sfu.ca/~ssurjano/robot.html

Wing Weight Function

class psdr.demos.WingWeight

The wing weight test function

This function models the weight of a wing based on several design parameters [VLSE_wing]:

\[\begin{split}&f(S_w, W_{fw}, A, \Lambda, q, \lambda, t_c, N_z, W_{dg}, W_p) := \\ &\quad 0.036 S_w^{0.758} W_{fw}^{0.0035} \left( \frac{A}{\cos^2 \Lambda} \right)^{0.6} q^{0.006} \lambda^{0.04} \left( \frac{100 t_c}{\cos \Lambda} \right)^{-0.3} (N_z W_{dg})^{0.49} + S_w W_p\end{split}\]

Variable	Interpretation
\(S_w \in [150, 200]\)	wing area (ft^2)
\(W_{fw} \in [220, 300]\)	weight of fuel in the wing (lb)
\(A \in [6,10]\)	aspect ratio
\(\Lambda \in [-10,10]\)	quarter-chord sweep (degrees)
\(q \in [16,45]\)	dynamic pressure at cruise (lb/ft^2)
\(\lambda \in [0.5,1]\)	taper ratio
\(t_c\in [0.08,0.18]\)	aerfoil thickness to chord ratio
\(N_z \in [2.5, 6]\)	ultimate load factor
\(W_{dg} \in [1700, 2500]\)	flight design gross weight (lb)
\(W_p \in [0.025, 0.08]\)	paint weight (lb/ft^2)

References

[VLSE_wing]

Virtual Library of Simulation Experiments, Wing Weight Function https://www.sfu.ca/~ssurjano/wingweight.html

Model-based Test Functions

OpenAeroStruct

class psdr.demos.OpenAeroStruct

A test problem from OpenAeroStruct

A test problem using OpenAeroStruct similar to that described in [JHM18].

References

[JHM18]

Open-source coupled aerostructural optimization using Python. John P. Jasa, John T. Hwang, and Joaquim R. R. A. Martins. Structural and Multidisciplinary Optimization (2018) 57:1815-1827 DOI: 10.1007/s00158-018-1912-8

MULTI-F

class psdr.demos.MULTIF(truncate=1e-07, level=0, su2_maxiter=5000, workdir=None, dask_client=None, **kwargs)

An interface to the MULTI-F multiphysics jet nozzle model test problem.

The function describes a multiphysics model of a jet nozzle [FMAA18]. The domain consists of 96 design parameters and 40 uncertain variables and the output returns the mass, thrust, and many different failure constraints. This uses a Docker image to encapsulate the dependencies for MULTI-F; to use this function, install Docker and pull the multif image

>> docker pull jeffreyhokanson/multif:v25

This function encodes multiple different fidelity levels: in order of increasing fidelity and approximate one-core cost

Level	Physics	Mesh	Mechanics	Runtime (sec)
0	1d Non-ideal	N/A	linear	20
1	1d Non-ideal	N/A	nonlinear	205
2	2d Euler	Coarse	linear	82
3	2d Euler	Medium	linear	256
4	2d Euler	Fine	linear	713
5	3d Euler	Coarse	linear	1083
6	3d Euler	Medium	linear	3123
7	3d Euler	Fine	linear	13350
8	2d RANS	Coarse	linear
9	2d RANS	Medium	linear
10	2d RANS	Fine	linear
11	3d RANS	Coarse	linear	80690
12	3d RANS	Medium	linear	175888
13	3d RANS	Fine	linear	408512

Note at the present time the 2d RANS levels are buggy and not recommended for use.

Parameters:

• truncate (float in [0,1]; default 1e-7) – If non-zero, truncate the random domains by the specified probability. Truncation is necessary to provide a bounded domain for use with most sampling strategies.
• level (int in [0, 13]; default 0) – What level of fidelity to run (see description above)
• su2_maxiter (int; default: 5000) – Number of iterations used to solve for the fluid flow in SU2
• workdir (str; default: None) – If defined, this specifies the location where temporary files are written while solving the problem
• keep_data (bool; default: False) – If True, do not delete the work files; i.e., use this to preserve the flow solution for later use
• verbose (bool; default: False) – If True, print the output of running MULTI-F to stdout
• dask_client (dask.distributed.Client or None) – If specified, allows distributed computation with this function.

References

[FMAA18]

Reliability-Based Design Optimization of a Supersonic Nozzle Richard W. Fenrich, Victorien Menier, Philip Avery, and Juan J. Alonso, 6th European Conference on Computational Mechanics http://www.eccm-ecfd2018.org/admin/files/filePaper/p437.pdf

NACA0012

class psdr.demos.NACA0012(n_lower=10, n_upper=10, fraction=0.01, dask_client=None, **kwargs)

The lift and drag of a NACA0012 airfoil perturbed by Hicks-Henne bump functions

This function computes the lift and drag of a NACA0012 airfoil perturbed by Hicks-Henne bump functions. The lift and drag are computed using [SU2] using Euler flow at Mach number 0.8 and angle of attack 1.25.

This example slightly modifies the configuration file that appears in the SU2 quickstart folder: (inv_NACA0012.cfg)[https://github.com/su2code/SU2/blob/master/QuickStart/inv_NACA0012.cfg].

Parameters:

• n_lower (int, optional (default: 10)) – Number of bump functions for the lower surface
• n_upper (int, optional (default: 10)) – Number of bump functions for the upper surface
• fraction (float, optional (default:0.01)) –

References

[SU2]

https://su2code.github.io/