信頼性評価における限界状態関数法をサポートするためのユーティリティ機能をGitHub上からパッケージとして提供する．パッケージ名はUtility，モジュール名LimitStateである．
モジュール内のクラスLSFMを継承したクラスをユーザが定義することにより、限界状態関数法を容易に利用できる。LSFMのメソッドを以下に示す。
| LSFMのメソッド | 内容 |
|---|---|
| RFn() | Rackvitz Fiessler法による設計点探索 |
| GetBeta() | 信頼性指標βの取得 |
| GetAlpha() | 感度ベクトルの取得 |
| GetPOF() | 破損確率の取得 |
| GetDP() | 設計点の取得 |
| GetConv() | RF法の収束回数の取得 |
| GetPSF() | 部分安全係数の取得 |
RSモデルへの適用例
# RSモデル
# 限界状態関数クラスを新規に定義して利用する場合
class GRS(ls.LSFM): #LSFMから継承する
def __init__(self,n,muX,sigmmaX,dist):
super().__init__(n,muX,sigmmaX,dist)
def gcalc(self,X): #g値を計算する仮想関数
R=X[0]
S=X[1]
g=R-S #限界状態関数
return g #必ずg値を戻す
def dGdXcalc(self,X): #gの微分ベクトルを計算する仮想関数
R=X[0]
S=X[1]
dGdX=[0]*len(X) #Xと同じ大きさのリストを作る
dGdX[0] =1.0
dGdX[1] =-1.0
return dGdX #必ず微分ベクトルを戻す
dist=['normal','normal'] #分布形はweibull,normal,lognormal,gumbel,uniformのいずれか
muX=[200,100]
sigmmaX=[10,20]
n=2
rs=GRS(n,muX,sigmmaX,dist) #インスタンスの生成
rs.RFn()
rs.GetBeta()
#
この結果，以下の出力を確認する．このように、LSFMクラスを継承した上で、__init__関数、仮想関数のgcalc、dGdXcalc関数を定義するだけで、限界状態関数法の計算が可能となる
限界状態関数や変数を文字列で定義し，平均値，標準偏差を指定する場合．最も簡単な利用法．
注意:変数として使用可能な文字は小文字のみです．
g="r-s"
var =['r','s']
muX =[200,100]
sigmmaX=[10,20]
dist =["normal", "normal"]
aa =ls.GeneralG(g,var,len(var),muX,sigmmaX,dist)
aa.RFn()
aa.GetBeta()
この結果，以下の出力を確認する．
4.47213595499958
下記文献の検証
"Probability, Reliability and Statistical Methods
in Engineering Design"
Achintya Haldar & Sankaran Mahadevan
P.218 Table 7.5, 7.6
# 限界状態関数を文字列で与えて利用する場合
g="As*fy*d*(1.0-eta*As*fy/b/d/fcd)-M"
var =["As","fy","fcd","b","d","eta","M"]
muX =[1.56, 47.7, 3.5, 8.0, 13.2, 0.59, 326.25]
covX =[0.036, 0.15, 0.21, 0.045, 0.086, 0.05, 0.17]
sigmmaX = list(np.array(covX)*np.array(muX))
dist =["normal", "normal", "normal", "normal" ,"normal" ,"normal" ,"normal" ]
aa =ls.GeneralG(g,var,len(var),muX,sigmmaX,dist)
aa.RFn()
a1=aa.GetBeta()
dist[6]="lognormal"
aa =ls.GeneralG(g,var,len(var),muX,sigmmaX,dist)
aa.RFn()
a2=aa.GetBeta()
dist=["lognormal", "lognormal", "lognormal", "lognormal" ,"lognormal" ,"lognormal" ,"normal"]
aa =ls.GeneralG(g,var,len(var),muX,sigmmaX,dist)
aa.RFn()
a3=aa.GetBeta()
dist=["lognormal", "lognormal", "lognormal", "lognormal" ,"lognormal" ,"lognormal" ,"lognormal"]
aa =ls.GeneralG(g,var,len(var),muX,sigmmaX,dist)
aa.RFn()
a4=aa.GetBeta()
[a1,a2,a3,a4]
この結果，以下の信頼性指標の出力を確認する．
[3.8330281658128675, 3.7612536656241864, 4.3876843910772845, 4.090647398473586]
例題3と同じ問題について直接限界状態関数を定義して利用する方法。PYTHONの数式処理機能を利用すると便利。
例題3よりは、処理速度が早くなる。
from sympy import *
import numpy as np
class Gmulti(ls.LSFM):
def __init__(self,n,muX,sigmmaX,dist):
As,fy,fcd,b,d,eta,M = symbols('As fy fcd b d eta M')
g=As*fy*d*(1.0-eta*As*fy/b/d/fcd)-M
#あらかじめ微分の文字列を数式処理で作っておく
self.gg=str(g)
#親クラスの初期化は、self.gg定義後に行うこと
super().__init__(n,muX,sigmmaX,dist)
self.d0=str(diff(g,As))
self.d1=str(diff(g,fy))
self.d2=str(diff(g,fcd))
self.d3=str(diff(g,b))
self.d4=str(diff(g,d))
self.d5=str(diff(g,eta))
self.d6=str(diff(g,M))
def gcalc(self,X):
As=X[0]
fy=X[1]
fcd=X[2]
b=X[3]
d=X[4]
eta=X[5]
M=X[6]
g=eval(self.gg)
return g
def dGdXcalc(self,X):
As=X[0]
fy=X[1]
fcd=X[2]
b=X[3]
d=X[4]
eta=X[5]
M=X[6]
dGdX=[0]*len(X)
#微分値の評価
dGdX[0]=eval(self.d0)
dGdX[1]=eval(self.d1)
dGdX[2]=eval(self.d2)
dGdX[3]=eval(self.d3)
dGdX[4]=eval(self.d4)
dGdX[5]=eval(self.d5)
dGdX[6]=eval(self.d6)
return dGdX
dist =["normal", "normal", "normal", "normal" ,"normal" ,"normal" ,"normal" ]
muX =[1.56, 47.7, 3.5, 8.0, 13.2, 0.59, 326.25]
covX =[0.036, 0.15, 0.21, 0.045, 0.086, 0.05, 0.17]
sigmmaX = list(np.array(covX)*np.array(muX))
n=7
rs =Gmulti(n,muX,sigmmaX,dist)
rs.RFn()
a1=rs.GetBeta()
dist[6]="lognormal"
rs=Gmulti(n,muX,sigmmaX,dist)
rs.RFn()
a2=rs.GetBeta()
dist=["lognormal", "lognormal", "lognormal", "lognormal" ,"lognormal" ,"lognormal" ,"normal"]
rs=Gmulti(n,muX,sigmmaX,dist)
rs.RFn()
a3=rs.GetBeta()
dist=["lognormal", "lognormal", "lognormal", "lognormal" ,"lognormal" ,"lognormal" ,"lognormal"]
rs=Gmulti(n,muX,sigmmaX,dist)
rs.RFn()
a4=rs.GetBeta()
[a1,a2,a3,a4]
この結果，以下の信頼性指標の出力を確認する．
[3.8330281658128675, 3.7612536656241864, 4.3876843910772845, 4.090647398473586]
例題1と同じ問題について辞書型データでデータを与えて信頼性評価を行う。例題1やこの例題では、設計点探索のアルゴリズムとしてRackwitz Fiessler法を用いている。この場合、探索の開始点を指定する必要があるが、デフォールトとしては基準化空間の原点が指定されるようになっている。原点以外の任意点を指定する場合の例題は例題6に示してある。
# RSモデル
# 限界状態関数クラスを新規に定義して利用する場合
class GRS(ls.LSFM): #LSFMから継承する
def __init__(self,n,muX,sigmmaX,dist):
super().__init__(n,muX,sigmmaX,dist)
def gcalc(self,X): #g値を計算する仮想関数
R=X[0]
S=X[1]
g=R-S #限界状態関数
return g #必ずg値を戻す
def dGdXcalc(self,X): #gの微分ベクトルを計算する仮想関数
R=X[0]
S=X[1]
dGdX=[0]*len(X) #Xと同じ大きさのリストを作る
dGdX[0] =1.0
dGdX[1] =-1.0
return dGdX #必ず微分ベクトルを戻す
# main program #
#入力データは辞書形式で表現
data={"r":{"mean":200,"cov":0.05,"dist":"normal"},#辞書形式で与える入力データ
"s":{"mean":100,"cov":0.2,"dist":"normal"} }
n,muX,sigmmaX,dist=ls.dict2data(data)#辞書型データの入力データへの変換
rs=GRS(n,muX,sigmmaX,dist)
print('val=',rs.gcalc(muX))#初期点に対する限界状態関数値
rs.RFn()
print('beta=',rs.GetBeta())#信頼性指標の取得
print('alpha=',rs.GetAlpha())#感度の取得
print('POF=',rs.GetPOF())#破損確率の取得
例題5と同じであるが、スタート点を変更する例題。例題5やこの例題では、設計点探索のアルゴリズムとしてRackwitz Fiessler法を用いている。この場合、探索の開始点を指定する必要があるが、デフォールトとしては基準化空間の原点(平均点)が指定されるようになっている。しかし、場合によっては、原点以外の任意の座標点を開始点として指定したいこともある。以下にその例題を示す。### 1 ###と### 2 ###が変更点である
class GRS(ls.LSFM): #LSFMから継承する
def __init__(self,n,muX,sigmmaX,dist):
super().__init__(n,muX,sigmmaX,dist)
def gcalc(self,X): #g値を計算する仮想関数
R=X[0]
S=X[1]
g=R-S #限界状態関数
return g #必ずg値を戻す
def dGdXcalc(self,X): #gの微分ベクトルを計算する仮想関数
R=X[0]
S=X[1]
dGdX=[0]*len(X) #Xと同じ大きさのリストを作る
dGdX[0] =1.0
dGdX[1] =-1.0
return dGdX #必ず微分ベクトルを戻すam #
#入力データは辞書形式で表現
data={"r":{"mean":200,"cov":0.05,"dist":"normal"},#辞書形式で与える入力データ
"s":{"mean":100,"cov":0.2,"dist":"normal"} }
n,muX,sigmmaX,dist=ls.dict2data(data)#辞書型データの入力データへの変換
rs=GRS(n,muX,sigmmaX,dist)
print('val=',rs.gcalc(muX))#muXに対する限界状態関数値
### 1 ### まず探索開始点のリストを定義。例えば[100,100]の場合について例示。
nstart=[100,100]
### 2 ### RFnメソッドの引数startに'Coordinate'を指定し、Xstartにnstartを指定する
rs.RFn(start='Coordinate',Xstart=nstart)#dataに対する信頼性評価
#start='Origin'のとき、原点からスタート、start='Coordinate'のとき、Xstartで与える座標からスタートする。デフォールトは'Origin'
print('beta=',rs.GetBeta())#信頼性指標の取得
print('alpha=',rs.GetAlpha())#感度の取得
print('POF=',rs.GetPOF())#破損確率の取得
設計点探索のアルゴリズムとしてRackwitz Fiessler法を用いている場合、local minimumに収束してしまう場合がある。このような事態を防ぐためには、探索のスタート点を広域から、徐々に絞り込んでいけばよい。最適化プログラムOptunaを使った例題を以下に示す。限界状態関数は例題6と同じものを使う。Optuna適用のために、まずclass Objectを用意する。ただし、この探索方法は一般に相当の時間を要する。local minimumへの収束の懸念がある場合にのみ行う。
class GRS(ls.LSFM): #LSFMから継承する
def __init__(self,n,muX,sigmmaX,dist):
super().__init__(n,muX,sigmmaX,dist)
def gcalc(self,X): #g値を計算する仮想関数
R=X[0]
S=X[1]
g=R-S #限界状態関数
return g #必ずg値を戻す
def dGdXcalc(self,X): #gの微分ベクトルを計算する仮想関数
R=X[0]
S=X[1]
dGdX=[0]*len(X) #Xと同じ大きさのリストを作る
dGdX[0] =1.0
dGdX[1] =-1.0
return dGdX #必ず微分ベクトルを戻すam #
#入力データは辞書形式で表現
data={"r":{"mean":200,"cov":0.05,"dist":"normal"},#辞書形式で与える入力データ
"s":{"mean":100,"cov":0.2,"dist":"normal"} }
n,muX,sigmmaX,dist=ls.dict2data(data)#辞書型データの入力データへの変換
import numpy as np
import optuna
class Objective:
def __init__(self,data):
# 変数X,yの初期化
self.data=data
def __call__(self, trial):
n,muX,sigmmaX,dist=ls.dict2data(self.data)
rs=GRS(n,muX,sigmmaX,dist)
### スタート点の調査範囲の指定
r = trial.suggest_uniform('r', 50,300)
s = trial.suggest_uniform('s', 50,300)
start=[r,s]
rs.RFn(start='Coordinate',Xstart=start)
return rs.GetBeta() #信頼性指標を戻す
# main program #
#入力データは辞書形式で表現
data={"r":{"mean":200,"cov":0.05,"dist":"normal"},#辞書形式で与える入力データ
"s":{"mean":100,"cov":0.2,"dist":"normal"} }
optuna.logging.disable_default_handler() #optunaの出力を抑制する
objective = Objective(data)
study = optuna.create_study(direction='minimize')
study.optimize(objective, timeout=60)
# 最小の信頼性指標と、その値を与えたスタート点の出力
print('params:', study.best_params,',best_value:',study.best_value)
出力例
params: {'r': 107.86532411893185, 's': 126.05989917926681} ,best_value: 4.47213595499958