ryujimiyaの日記

C#を使って数値解析したい

時間領域FEMによるLamb wave弾性波導波路伝達問題ーPMLを用いる場合ー

1. はじめに

これまで弾性波プレート導波路の伝達問題を周波数領域FEMで解いてきました。

今回時間領域FEM(Time Domain Finite Element Method)でLamb wave弾性波導波路の伝達問題を解いています。

解析端の境界条件にはPML(Perfectly Matched Layers)を用いています。

 

2. 定式化

pdfにまとめました。

周波数領域FEMのPMLとちがって少し複雑です。畳み込み積分(convolution)を用いたPML定式化になります。

Time Domain FEM Formulations for Lamb Wave Elastic Plate
Waveguide Discontinuity Problem Using Perfectly Matched
Layers(PML)

http://starlightparade.usamimi.info/ivyfem/doc/ElasticLambWaveguidePMLTD.pdf?p=0

 

3. 計算結果

ポアソン比ν=0.31とします。ヤング率Eや密度ρは反射、透過特性には無関係なので、E=1、ρ=1としておきます。

導波路に亀裂が入ったときの反射、透過特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

導波路幅をWとすると、図の幅d=W/2になります。

亀裂の深さhは、h/d = 0.8としています。

また亀裂の角度は2θ=10°としています。

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

入射波は素のガウシアンパルスを用いています。

不連続部の変位(displacement)ベクトル分布の時間変化 (ux(t), uy(t))

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

の計算結果は次のようになりました。

PML終端でスプリアス反射が見られたので物理領域に影響しないようにPMLは十分な厚さをとるひつようがありました。

反射、透過係数はux(t)、uy(t)をフーリエ変換(FFTです)して周波数特性を求めました。

f:id:ryujimiya:20200908065114j:plain

f:id:ryujimiya:20200908065134j:plain

f:id:ryujimiya:20200908065144j:plain

f:id:ryujimiya:20200908065157j:plain

f:id:ryujimiya:20200908065209j:plain

f:id:ryujimiya:20200908065219j:plain

 

4. まとめ

PMLを用いた時間領域FEMをLamb wave弾性波導波路伝達問題に適用しました。




SH wave(Shear Horizontal Wave)弾性波プレート導波路の伝達問題ーPMLを用いた場合とABC(吸収境界条件)を用いた場合ー

1. はじめに

前にSH wave弾性波導波路の伝達問題を計算しましたが、そのときポート境界条件は固有モード展開を用いました。

本記事ではポート境界条件にPMLを装荷した場合と、ABC(吸収境界条件, Absorbing Boundary Conditions)を用いた場合の2つの定式化をして計算しています。

 

2.定式化

pdfにまとめました。

〇PML

FEM Formulations for SH(Shear Horizontal) Wave Elastic Plate
Waveguide Discontinuity Problem Using Perfectly Matched
Layers(PML)

http://starlightparade.usamimi.info/ivyfem/doc/ElasticSHWaveguidePML.pdf?p=0

〇ABC

SH Wave(Shear Horizontal wave) Elastic Plate Waveguides
Discontinuity Problem Using Absorbing Boundary
Conditions(ABC)

http://starlightparade.usamimi.info/ivyfem/doc/ElasticSHWaveguideABC.pdf?p=0

 

3. 計算結果

材料

ポアソン比ν=0.31

ヤング率E=1

密度ρ=1

入力導波路の幅の50%の出力導波路をdown stepで接続したときの反射特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

入力導波路幅をW=W1とすると、図の幅d=W(1/2)になります。

図は導波路をmid-planeで半分にしたものです。

down step 

出力導波路幅W2=(1/2)W1

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

不連続部の変位(displacement)の分布(サーモグラフィー) Re(uz)

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

入力導波路の固有モードの分散特性、分布

の計算結果は次のようになりました。

 

3.1. PMLの場合の計算結果

f:id:ryujimiya:20200907183725j:plain

f:id:ryujimiya:20200907183743j:plain

f:id:ryujimiya:20200907183758j:plain

f:id:ryujimiya:20200907183812j:plain

f:id:ryujimiya:20200907183824j:plain

f:id:ryujimiya:20200907183836j:plain

f:id:ryujimiya:20200907183848j:plain

f:id:ryujimiya:20200907183904j:plain

 

3.2. ABCの場合の計算結果

f:id:ryujimiya:20200907183930j:plain

f:id:ryujimiya:20200907183942j:plain

f:id:ryujimiya:20200907183954j:plain

f:id:ryujimiya:20200907184006j:plain

f:id:ryujimiya:20200907184022j:plain

f:id:ryujimiya:20200907184033j:plain

f:id:ryujimiya:20200907184045j:plain

f:id:ryujimiya:20200907184056j:plain

 

4. まとめ

SH wave弾性波プレート導波路の伝達問題に対して、PMLを用いた定式化とABCを用いた定式化の2つの定式化を行い、down-stepの反射周波数特性を計算しました。

 

吸収境界条件(ABC,Absorbing Boundary Conditions)を用いたLamb wave弾性波プレート導波路の伝達問題の計算

1. はじめに

これまでにLamb waveの弾性波プレート導波路の伝達問題(transmission problem, 不連続問題 discontinuity problem, 散乱問題 scattering problem)について半無限長導波路を表現するために固有モード展開を用いる方法とPML(Perfectly Matched Layers)を用いる方法を試しました。

本記事ではABC(吸収境界条件、Absorbing Boundary Conditions)を用いる方法を試しています。

 

2. 定式化

不連続部に接続する2ポート境界にABCを課します。内部領域に入射面および反射、透過波を観測する参照面を設けます。

定式化をpdfにまとめました。

Lamb Wave Elastic Plate Waveguides Discontinuity Problem
Using Absorbing Boundary Conditions(ABC)

http://starlightparade.usamimi.info/ivyfem/doc/ElasticLambWaveguideABC.pdf?p=0

Lamb waveのABCはどうやらないようなので、地震波(seismic wave, P wave and S wave)のABCを用いました。

f:id:ryujimiya:20200814090839j:plain

3. 計算結果

ポアソン比ν=0.31とします。ヤング率Eや密度ρは反射、透過特性には無関係なので、E=1、ρ=1としておきます。

導波路に亀裂が入ったときの反射、透過特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

導波路幅をWとすると、図の幅d=W/2になります。

亀裂の深さhは、h/d = 0.8としています。

また亀裂の角度は2θ=10°としています。

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

不連続部の変位(displacement)ベクトル分布 Re(ux, -j uy)

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

入力導波路の固有モードの分散特性、分布

の計算結果は次のようになりました。

f:id:ryujimiya:20200814091126j:plain

f:id:ryujimiya:20200814091137j:plain

f:id:ryujimiya:20200814091148j:plain

f:id:ryujimiya:20200814091201j:plain

f:id:ryujimiya:20200814091211j:plain

f:id:ryujimiya:20200814091220j:plain

f:id:ryujimiya:20200814091230j:plain

f:id:ryujimiya:20200814091241j:plain

ABCの境界を、入射面あるいはport2から十分な距離をとると、Lamb waveのABCでないにもかかわらず良好な結果が得られました。

 

4. まとめ

ABCを用いたLamb wave弾性波プレート導波路の伝達問題をFEMで定式化し、亀裂の反射特性を計算しました。

PML(Perfectly Matched Layers,完全整合層)を用いたLamb wave弾性波導波路の伝達問題のFEM計算

1. はじめに

前に弾性波導波路(elastic waveguides)のLamb wave伝達問題(transmission problem, 不連続問題 discontinuity problem、散乱問題 scattering problem)をFEMで解きましたが、ポート境界は固有モード展開を用いて表現しました。

本記事では、ポート境界に完全整合層(PML, Perfectly Matched Layers)を用いた周波数領域FEM定式化を行い同じ亀裂の問題を解きます。

 

2. 問題設定と定式化

x軸に関して上下対称な導波路不連続部を考えます。

不連続部は2 ポートの半無限導波路に接続されている問題を解くのですが、今回はこの半無限導波路の代わりにポートをPML て終端する方法をとります。
入射波の励振面 (入射面、incident-plane)が別に必要となりますが、これを内部領域に設けます。

f:id:ryujimiya:20200812090834j:plain

定式化をpdfにまとめました。

FEM Formulations for Lamb Wave Elastic Plate Waveguide
Discontinuity Problem Using Perfectly Matched Layers(PML)

http://starlightparade.usamimi.info/ivyfem/doc/ElasticLambWaveguidePML.pdf?p=0

 

3.計算結果

ポアソン比ν=0.31とします。ヤング率Eや密度ρは反射、透過特性には無関係なので、E=1、ρ=1としておきます。

導波路に亀裂が入ったときの反射、透過特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

導波路幅をWとすると、図の幅d=W/2になります。

亀裂の深さhは、h/d = 0.8としています。

また亀裂の角度は2θ=10°としています。

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

不連続部の変位(displacement)ベクトル分布 Re(ux, -j uy)

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

入力導波路の固有モードの分散特性、分布

の計算結果は次のようになりました。

f:id:ryujimiya:20200811214239j:plain

f:id:ryujimiya:20200811214258j:plain

f:id:ryujimiya:20200811214316j:plain

f:id:ryujimiya:20200811214328j:plain

f:id:ryujimiya:20200811214344j:plain

f:id:ryujimiya:20200811214358j:plain

f:id:ryujimiya:20200811214410j:plain

f:id:ryujimiya:20200811214421j:plain

f:id:ryujimiya:20200811214432j:plain

 

4. まとめ

PMLを用いたLamb wave 弾性波導波路の伝達問題の周波数領域FEM定式化を行い、亀裂の反射、透過の周波数特性を計算しました。

弾性波プレート導波路(elastic plate waveguides)、SH wave(Shear Horizontal wave)の伝達問題(不連続問題、散乱問題)

1. はじめに

前記事で弾性波プレート導波路を伝搬するLamb waveの伝達問題を定式化、計算しました。

本記事ではSH wave(Shear Horizontal wave)の伝達問題を定式化、計算します。

 

2. SH waveの伝達問題の定式化

pdfにまとめました。

Frequency Domain FEM Formulations for SH Wave (Shear
Horizontal Wave) Elastic Plate Waveguides – Eigenvalue
Problem and Discontinuity Problem –

http://starlightparade.usamimi.info/ivyfem/doc/ElasticSHWaveguide.pdf?p=0

 

3.計算結果

材料

ポアソン比ν=0.31

ヤング率E=1

密度ρ=1

入力導波路の幅の50%の出力導波路をdown stepで接続したときの反射特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

入力導波路幅をW=W1とすると、図の幅d=W(1/2)になります。

図は導波路をmid-planeで半分にしたものです。

down step 

出力導波路幅W2=(1/2)W1

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

不連続部の変位(displacement)の分布(サーモグラフィー) Re(uz)

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

入力導波路の固有モードの分散特性、分布

の計算結果は次のようになりました。

f:id:ryujimiya:20200816072322j:plain

f:id:ryujimiya:20200816072335j:plain

f:id:ryujimiya:20200816072347j:plain

f:id:ryujimiya:20200816072401j:plain

f:id:ryujimiya:20200816072413j:plain

f:id:ryujimiya:20200816072425j:plain

f:id:ryujimiya:20200816072437j:plain

f:id:ryujimiya:20200816072448j:plain

 

下記文献に同じ50% down stepの計算結果が載っています。

Alessandro emma,
"The Interaction of Guided Waves with Discontinuities in Structures",
thesis submitted to the University of London for the degree of
Doctor of Philosophy,
Department of Mechanical Engineering
Imperial College of Science, Technology and Medicine,
London SW7 2BX,
January 2003
Figure 3.4

文献によると0周波数では|S11|=0.333です。

今回の計算結果では規格化周波数fn=0のとき|S11|=0.330なので文献と一致した結果が得られています。

また|S11|の周波数特性の傾向も文献と同じ様です。(文献は周波数*導波路幅を横軸に示していますが、材料定数の記載がないため具体的な対応はとれませんでした、ただ0周波数はfn=0に対応するので上記に示しました)

 

4.まとめ

SH wave弾性波プレート導波路の伝達問題の定式化を行い、down stepの反射特性を計算しました。

 

弾性波プレート導波路(elastic plate waveguides)、Lamb waveの伝達問題(不連続問題、散乱問題)

1. はじめに

弾性波プレート導波路を伝搬するLamb waveの基本モードをポート1から入射したとき、不連続部により反射する、またはポート2に透過するのをFEMで計算しています。

基本モードは対称な分布のモードのため、中央面(mid-plane)を境に半分にした領域で計算できます。

f:id:ryujimiya:20200730072837j:plain

f:id:ryujimiya:20200730072839j:plain

2. 弾性波プレート導波路、Lamb waveの伝達問題の定式化

pdfにまとめました。

Frequency Domain FEM Formulations for Lamb Wave Elastic
Plate Waveguides – Eigenvalue Problem and Discontinuity
Problem

http://starlightparade.usamimi.info/ivyfem/doc/ElasticLambWaveguide.pdf?p=0

 

3. 計算結果

ポアソン比ν=0.31とします。ヤング率Eや密度ρは反射、透過特性には無関係なので、E=1、ρ=1としておきます。

導波路に亀裂が入ったときの反射、透過特性を計算しました。

入出力導波路の固有モードも(解析解でなく)FEMで計算してます。

図は導波路をmid-planeで半分にしたものです。

導波路幅をWとすると、図の幅d=W/2になります。

亀裂の深さhは、h/d = 0.8としています。

また亀裂の角度は2θ=10°としています。

下記文献を参照してください。

Masanori Koshiba, Shoji Karakida, and Michio Suzuki, ”Finite-Element Analysis of Lamb Wave Scattering in an Elastic Plate Waveguide”, IEEE Transactions on Sonics and Ultrasonics, vol. su-31, no. 1,
January 1984

規格化周波数

fn = ksW/π

波数

ks=√(ρ/μ) 

不連続部の変位(displacement)ベクトル分布 Re(ux, -j uy)

入力導波路ポート1の反射係数|S11|、出力導波路ポート2の透過係数|S21|

入力導波路の固有モードの分散特性、分布

の計算結果は次のようになりました。

f:id:ryujimiya:20200802155937j:plain

f:id:ryujimiya:20200802155949j:plain

f:id:ryujimiya:20200802160003j:plain

f:id:ryujimiya:20200802160017j:plain

f:id:ryujimiya:20200802160029j:plain

f:id:ryujimiya:20200802160038j:plain

f:id:ryujimiya:20200802160049j:plain









 

 

 

 

 

 

 

h/W=0.8のとき、fn=ksW/π=1で、反射係数|S11|=0.7になることが文献に載っています。今回の計算結果と一致しています。

また、h/W=0.5のときfn=1で|S11|=0.4になることが文献に載っていますがこれも計算したところ一致することを確認しました。

 

4. まとめ

弾性波プレート導波路の伝達問題をFEMで定式化し、導波路に亀裂が入ったときの反射、透過特性を計算しました。

FEMによる渦電流(eddy current)場の過渡応答の計算

1. はじめに

変位電流が無視できるような周波数帯は準定常場とよばれ、渦電流場問題で使われます。

支配方程式は、ベクトルポテンシャルAとスカラーポテンシャルφで表したものが使われることが多いようです(A-φ法)。ポテンシャルにはなんらかのゲージを課す必要がありますが、ここではクーロンゲージ(∇・A= 0)を用いています。

 

2. 渦電流場の過渡解析の定式化

pdfにまとめました。

Two Dimensional Transient Eddy Current Problem

http://starlightparade.usamimi.info/ivyfem/doc/EddyCurrent2D.pdf?p=0

 

3, 計算結果

3.1. 鋼鉄の角材の両端に直流電圧をかける

f:id:ryujimiya:20200701201234j:plain

図の色を付けた鋼鉄の断面を解析領域とします。

鋼鉄
 断面の幅 w = 1.0;
 断面の高さ h = w / 2
 比透磁率μr = 300

 導電率σ = 8.33 x 10^5 S/m

印加電圧の勾配 ∇φ = -1 V/m

空気と鋼鉄の境界ではベクトルポテンシャルのz成分A=0としています。

時間の経過順に分布を示します。

サーモグラフィーはベクトルポテンシャルのz成分A、矢印は断面内磁束密度ベクトルBです。

f:id:ryujimiya:20200701201432j:plain

f:id:ryujimiya:20200701201502j:plain

f:id:ryujimiya:20200701201514j:plain

f:id:ryujimiya:20200701201526j:plain

時間が経過するにつれて磁束が鋼鉄内部に浸透していく様子が見られます。

3.2. スイッチトリラクタンスモータ(Switched reluctance motor, SRM)

f:id:ryujimiya:20200701201859j:plain

ステーター
 比透磁率 = 1000

 導電率 = 1.39 x 10^6 S/m

ローター
 比透磁率 = 1000

 導電率 = 1.39 x 10^6 S/m

コイル

 導電率 = 59.0 x 10^6 S/m

 100Vの三相交流を印加

f:id:ryujimiya:20200701202349j:plain

f:id:ryujimiya:20200701202400j:plain

f:id:ryujimiya:20200701202409j:plain

 

4. まとめ

2次元渦電流場をA-φ法で計算しました。