Python: calculation of required plate thickness of exposed type penstock against buckling pressure by Brent method

A critical buckling pressure of exposed type penstocks can be obtained by following equation.

$\begin{equation} p_k=\cfrac{2\cdot E_s}{1-\nu_s{}^2}\cdot \left(\cfrac{t}{D_0'}\right)^3 \end{equation}$

$p_k$	critical buckling pressure of penstock
$t$	plate thickness excluding corrosion allowance
$D_0'$	design external diameter
$E_s$	elastic modulus of penstock (=206,000MPa)
$\nu_s$	Poisson's modulus of penstock (=0.3)

usually, a corrosion allowance is set to $\epsilon=1.5 mm$ , and the design plate thickness including corrosion allowance $t_0$ and the design external diameter $D_0'$ can be calculated as follows.

$\begin{equation} t_0=t+\epsilon \qquad D_0'=D_0+t_0 \end{equation}$

Where, $D_0$ means the design internal diameter.

In case of exposed type penstocks, it is required to withstand the buckling pressure of 0.02 MPa which includes a safety factor of 1.5, and the plate thickness with this capacity shall be defined as a lower limit of plate thickness.

The equation to be solved by Brent method is shown below. The equation has a shape of $f=0$ .

$\begin{equation} f=p_k-\cfrac{2\cdot E_s}{1-\nu_s{}^2}\cdot \left(\cfrac{t}{D_0+t+\epsilon}\right)^3=0 \end{equation}$

Required two initial values for Brent method are given as t1=1.0, t2=50.0 in the program.

Program code is shown below.

import numpy as np
from scipy import optimize


def func(t,d0,eps,pk):
    Es=206000 # elastic modulus of steel
    po=0.3    # Poisson's ratio of steel
    f=pk-2*Es/(1-po**2)*(t/(d0+t+eps))**3
    return f


def main():
    d0=4000.0 # internal diameter of penstock
    eps=1.5   # corrosion alloowance
    pk=0.02   # critical buckling pressure
    t1=1.0    # initial value for Brent method
    t2=50.0   # initial value for Brent method
    tt=optimize.brentq(func,t1,t2,args=(d0,eps,pk))
    t0=np.ceil(tt+eps) # required plate thickness
    print('t + eps=',tt+eps)
    print('t0=',t0)


#==============
# Execution
#==============
if __name__ == '__main__': main()