IGMBaseLib 1.0
|
This module provides some subroutines to construct the SPR-GC-grid which is obtained by modifying the STD-Grid using the spring dynamics. More...
Functions/Subroutines | |
subroutine, public | construct_SPRGCGrid (icgrid, tstep_limit) |
Constructs the SPR-GC-grid. | |
subroutine, public | set_iteration_limiter (limiter_flag) |
時間積分に伴う逐次計算ループの反復回数を制限するか否かのフラグを設定する. | |
subroutine | modify_STDGrid_with_spring_dyn (icgrid) |
Generates the SPR-GC-grid obtained by modifying the STD-grid. | |
subroutine | move_GP_with_spring_dyn (icgrid, idMin, idMax, EMin, EMax) |
全球領域の各格子点に対してバネ力学を適用し, 仮想的なつり合い位置を求める. | |
subroutine | check_global_convergence (icgrid, conv_total_ave, conv_max, EMin, EMax) |
バネでつながれた格子点の振動が, 全球で収束したかを調べる. | |
real(8) | check_convergence (r0, pts, pt_num) |
収束度合いを, Tomita etal, 2001 の式(22) を用いてチェックする. | |
integer | generate_near_GP6index (rcID, GP0_i, GP0_j, EMin, EMax) |
格子点(rcID, GP0_i, GP0_j) の近傍にある格子点 6 個に対応する, 配列のインデックスを生成する. | |
integer | generate_near_GP5index (icgrid, rcID, GP0_i, GP0_j) |
格子点(rcID, GP0_i, GP0_j) の近傍にある格子点 5 個に対応する, 配列のインデックスを生成する. | |
subroutine | integ_ode_1step_spring_dynamics (r0, w0, GP, pt_num) |
固定された 6 点とバネでつながれた質点の位置の時間発展を計算する. | |
real(DP), dimension(3) | f1 (r0, w0) |
Calculates the value of right-hand side of the following ordinary differential equation described in Eq(21) in Tomita et al(2001). | |
real(8), dimension(1:3) | f2 (r0, w0, pts, pt_num, alpha, k, eq_d, radius) |
Calculates the value of right-hand side of the following ordinary differential equation described in Eq(20) in Tomita et al(2001). | |
real(8), dimension(3) | calc_ei (r0, ri) |
点 G_0(r0) から 点 Gi(ri) を結ぶベクトルの単位ベクトルを取得する. | |
Variables | |
real(DP) | eq_d |
The length of spring without imposed force. | |
real(DP) | radius |
The radius of a sphere which an icosahedron is embedded into. | |
real(DP) | beta = 0.4d0 |
Tuning parameter. | |
real(DP) | alpha = 2000.0d0 |
Frictional constant. | |
real(DP) | k = 1.0d06 |
Spring constant. | |
real(DP) | dt = 0.001d0 |
Time step. | |
integer | end_tstep = 20000 |
The number of time integration steps. | |
logical | tstep_limiter = .false. |
時間積分に伴う逐次計算ループの反復回数を制限するか否かのフラグ. デフォルトは, ループ回数を制限しない. | |
integer | conv_check_step = 100 |
The step interval for which cheking the convergence is performed. | |
integer, allocatable | GP_IJ |
ある格子点を取り囲む格子点(6 個 or 5 個)のインデックスを保持する配列. |
This module provides some subroutines to construct the SPR-GC-grid which is obtained by modifying the STD-Grid using the spring dynamics.
Copyright (C) GFD Dennou Club, 2011-2012. All rights reserved.
license ??
real(8),dimension(3) SPRGCGrid_Builder::calc_ei | ( | real(8),dimension(3),intent(in) | r0, |
real(8),dimension(3),intent(in) | ri | ||
) | [private] |
点 G_0(r0) から 点 Gi(ri) を結ぶベクトルの単位ベクトルを取得する.
[in] | r0 | The position vector of point r0. |
[in] | ri | The position vector of point ri. |
Definition at line 894 of file SPRGCGrid_Builder.f90.
real(8) SPRGCGrid_Builder::check_convergence | ( | real(8),dimension(3),intent(in) | r0, |
real(8),dimension(pt_num, 3),intent(in) | pts, | ||
integer,intent(in) | pt_num | ||
) | [private] |
収束度合いを, Tomita etal, 2001 の式(22) を用いてチェックする.
戻りがゼロに近ければ近いほど, 質点はつり合いの位置に存在することになる.
[in] | r0 | The position vector of P0. |
[in] | pts | The position vectors of the points surrounding point P0. |
[in] | pt_num | The number of the grid points surrounding point P0. |
Definition at line 560 of file SPRGCGrid_Builder.f90.
subroutine SPRGCGrid_Builder::check_global_convergence | ( | type(IcGrid2D_FVM),intent(in) | icgrid, |
real(DP),intent(out) | conv_total_ave, | ||
real(DP),intent(out) | conv_max, | ||
integer,intent(in) | EMin, | ||
integer,intent(in) | EMax | ||
) | [private] |
バネでつながれた格子点の振動が, 全球で収束したかを調べる.
[in] | conv_toal_ave | |
[in] | conv_max | |
[in] | EMin | 正二十面格子の座標データを保持する配列の袖領域も含めないインデックスの最小値. |
[in] | EMax | 正二十面格子の座標データを保持する配列の袖領域も含めないインデックスの最大値. |
Definition at line 482 of file SPRGCGrid_Builder.f90.
subroutine,public SPRGCGrid_Builder::construct_SPRGCGrid | ( | type(IcGrid2D_FVM),intent(inout) | icgrid, |
integer,intent(in),optional | tstep_limit | ||
) |
Constructs the SPR-GC-grid.
[in,out] | icgrid | The variable of derived type IcGrid2D_FVM. |
Definition at line 156 of file SPRGCGrid_Builder.f90.
real(DP),dimension(3) SPRGCGrid_Builder::f1 | ( | real(DP),dimension(3),intent(in) | r0, |
real(DP),dimension(3),intent(in) | w0 | ||
) | [private] |
Calculates the value of right-hand side of the following ordinary differential equation described in Eq(21) in Tomita et al(2001).
[in] | r0 | The position vector of point r0. |
[in] | w0 | The velocity vector of point r0. |
Definition at line 793 of file SPRGCGrid_Builder.f90.
real(8),dimension(1:3) SPRGCGrid_Builder::f2 | ( | real(DP),dimension(3),intent(in) | r0, |
real(DP),dimension(3),intent(in) | w0, | ||
real(DP),dimension(pt_num,3),intent(in) | pts, | ||
integer,intent(in) | pt_num, | ||
real(DP),intent(in) | alpha, | ||
real(DP),intent(in) | k, | ||
real(DP),intent(in) | eq_d, | ||
real(DP),intent(in) | radius | ||
) | [private] |
Calculates the value of right-hand side of the following ordinary differential equation described in Eq(20) in Tomita et al(2001).
[in] | r0 | The position vector of point r0. |
[in] | w0 | The velocity vector of point r0. |
[in] | pts | |
[in] | pt_num | |
[in] | alpha | |
[in] | k | |
[in] | eq_d | |
[in] | radius | The radius of a sphere which an icosahedron is embedded into. |
Definition at line 842 of file SPRGCGrid_Builder.f90.
integer SPRGCGrid_Builder::generate_near_GP5index | ( | type(IcGrid2D_FVM),intent(inout) | icgrid, |
integer,intent(in) | rcID, | ||
integer,intent(in) | GP0_i, | ||
integer,intent(in) | GP0_j | ||
) | [private] |
格子点(rcID, GP0_i, GP0_j) の近傍にある格子点 5 個に対応する, 配列のインデックスを生成する.
[in,out] | icgrid | A variable of derived type IcGrid2D_FVM. |
[in] | rcID | |
[in] | GP0_i | |
[in] | GP0_j |
Definition at line 667 of file SPRGCGrid_Builder.f90.
integer SPRGCGrid_Builder::generate_near_GP6index | ( | integer,intent(in) | rcID, |
integer,intent(in) | GP0_i, | ||
integer,intent(in) | GP0_j, | ||
integer | EMin, | ||
integer | EMax | ||
) | [private] |
格子点(rcID, GP0_i, GP0_j) の近傍にある格子点 6 個に対応する, 配列のインデックスを生成する.
[in] | rcID | |
[in] | GP0_i | |
[in] | GP0_j | |
[in] | EMin | |
[in] | EMax |
Definition at line 609 of file SPRGCGrid_Builder.f90.
subroutine SPRGCGrid_Builder::integ_ode_1step_spring_dynamics | ( | real(DP),dimension(3),intent(inout) | r0, |
real(DP),dimension(3),intent(inout) | w0, | ||
real(DP),dimension(pt_num,3),intent(in) | GP, | ||
integer,intent(in) | pt_num | ||
) | [private] |
固定された 6 点とバネでつながれた質点の位置の時間発展を計算する.
常微分方程式系の時間積分の方法は, 3 次の Runge=Kutta 法である.
[in,out] | r0 | |
[in,out] | w0 | |
[in] | GP | |
[in] | pt_num | |
[in] | dt | |
[in] | alpha | |
[in] | k | |
[in] | eq_d | |
[in] | ic_radius | The radius of a sphere which an icosahedron is embedded into. |
Definition at line 724 of file SPRGCGrid_Builder.f90.
subroutine SPRGCGrid_Builder::modify_STDGrid_with_spring_dyn | ( | type(IcGrid2D_FVM),intent(inout) | icgrid | ) | [private] |
Generates the SPR-GC-grid obtained by modifying the STD-grid.
Specifically, the STD-grid configuration is modified by employing the spring dynamics. After that, control volumes are defined and grid points are moved to the gravitational centers of the control volumes in the same way as STD-GC-grid.
このバネ力学を用いた修正により, コントロールボリュームの面積・歪みが, 全球領域で単調に変化するようになる. STD-grid の場合, コントロールボリュームの面積・歪みがフラクタル的な分布をしているために, 数値的に微分値を評価する際, 組織的な数値誤差が発生してしまうが, この修正により数値誤差は改善される.
[in,out] | icgrid | A variable of derived type IcGrid2D_FVM. |
Definition at line 268 of file SPRGCGrid_Builder.f90.
subroutine SPRGCGrid_Builder::move_GP_with_spring_dyn | ( | type(IcGrid2D_FVM),intent(inout) | icgrid, |
integer,intent(in) | idMin, | ||
integer,intent(in) | idMax, | ||
integer,intent(in) | EMin, | ||
integer,intent(in) | EMax | ||
) | [private] |
全球領域の各格子点に対してバネ力学を適用し, 仮想的なつり合い位置を求める.
各格子点を仮想的な質点とみなし, それぞれ近傍の質点 6 個とを仮想的なバネでつなぐ. そして, 全質点がつり合いの状態にある位置を求める.
[in] | idMin | 正二十面格子の座標データを保持する配列の袖領域も含めたインデックスの最小値. |
[in] | idMax | 正二十面格子の座標データを保持する配列の袖領域も含めたインデックスの最大値. |
[in] | EMin | 正二十面格子の座標データを保持する配列の袖領域も含めないインデックスの最小値. |
[in] | EMax | 正二十面格子の座標データを保持する配列の袖領域も含めないインデックスの最大値. |
Definition at line 348 of file SPRGCGrid_Builder.f90.
subroutine,public SPRGCGrid_Builder::set_iteration_limiter | ( | logical,intent(in) | limiter_flag | ) |
時間積分に伴う逐次計算ループの反復回数を制限するか否かのフラグを設定する.
デフォルトではフラグ tstep_limiter は真であり, 最大 2 万回で逐次計算ループを抜ける. このフラグを偽に設定した場合には, check_convergence による収束判定が真を返すまで, 逐次計算ループを続行される.
[in] | limiter_flag | 反復回数を制限するか否かのフラグ. |
Definition at line 221 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::alpha = 2000.0d0 |
Frictional constant.
Definition at line 112 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::beta = 0.4d0 |
Tuning parameter.
Definition at line 105 of file SPRGCGrid_Builder.f90.
integer SPRGCGrid_Builder::conv_check_step = 100 |
The step interval for which cheking the convergence is performed.
Definition at line 137 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::dt = 0.001d0 |
Time step.
Definition at line 120 of file SPRGCGrid_Builder.f90.
integer SPRGCGrid_Builder::end_tstep = 20000 |
The number of time integration steps.
Definition at line 124 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::eq_d |
The length of spring without imposed force.
Definition at line 97 of file SPRGCGrid_Builder.f90.
integer,allocatable SPRGCGrid_Builder::GP_IJ |
ある格子点を取り囲む格子点(6 個 or 5 個)のインデックスを保持する配列.
Definition at line 141 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::k = 1.0d06 |
Spring constant.
Definition at line 116 of file SPRGCGrid_Builder.f90.
real(DP) SPRGCGrid_Builder::radius |
The radius of a sphere which an icosahedron is embedded into.
Definition at line 101 of file SPRGCGrid_Builder.f90.
logical SPRGCGrid_Builder::tstep_limiter = .false. |
時間積分に伴う逐次計算ループの反復回数を制限するか否かのフラグ. デフォルトは, ループ回数を制限しない.
Definition at line 130 of file SPRGCGrid_Builder.f90.