IGMBaseLib 1.0
|
00001 00021 module IcGrid2D_FVM_Manager 00022 00023 ! モジュール引用 ; Use statements 00024 ! 00025 00026 ! 種類型パラメータ 00027 ! Kind type parameter 00028 ! 00029 use dc_types, only: DP ! 倍精度実数型 ; Double recision 00030 00031 ! 数学定数 00032 ! Mathematic 00033 ! 00034 use igmcore_math, only: PI ! 円周率. PI. 00035 00036 ! 線形代数 00037 ! Linear algebra 00038 ! 00039 use igmcore_linear_algebra, only: vec_length, & 00040 & rotateZ 00041 00042 ! 座標変換 00043 ! Coordinate conversion 00044 ! 00045 use igmcore_coordinate_conversion, only: geo_to_orth_pos 00046 00047 ! 宣言文 ; Declaration statemenets 00048 ! 00049 00050 implicit none 00051 private 00052 00053 ! 公開変数 00054 ! Public variables 00055 ! 00056 00059 integer, parameter, public :: RC_REGIONS_NUM = 10 00060 00063 integer, parameter, public :: ICOSAHEDRON_VERTEX_NUM = 12 00064 00067 integer, parameter, public :: CELL_POINTS_NUM = 6 00068 00071 integer, parameter, public :: NOT_POLE_FLAG = 0 00072 00075 integer, parameter, public :: NORTH_POLE_FLAG = 1 00076 00079 integer, parameter, public :: SOUTH_POLE_FLAG = 2 00080 00081 ! 格子タイプを表す定数 00082 ! 00083 00086 integer, parameter, public :: GTYPE_STDGRID = 1 00087 00090 integer, parameter, public :: GTYPE_STDGCGRID = 2 00091 00094 integer, parameter, public :: GTYPE_SPRGCGRID = 3 00095 00096 00097 ! 00103 type, public :: IcGrid2D_FVM 00104 00105 ! 公開メンバ変数 00106 ! Public memeber variables 00107 ! 00108 00118 real(DP), pointer :: rcs_AGrid(:,:,:,:) => null() 00119 ! f2003: real(DP), public, allocatable :: rcs_AGrid(:,:,:,:) 00120 00131 real(DP), pointer :: rcs_CV(:,:,:,:,:) => null() 00132 ! f2003: real(DP), public, allocatable :: rcs_CV(:,:,:,:,:) 00133 00144 real(DP), pointer :: rcs_CVArea(:,:,:) => null() 00145 ! f2003: real(DP), public, allocatable :: rcs_CVArea(:,:,:) 00146 00147 ! 非公開メンバ変数 00148 ! Private member variables 00149 ! 00150 00153 integer :: h_glevel 00154 ! f2003: integer, private :: h_glevel 00155 00158 real(DP) :: icRadius 00159 ! f2003: real(DP), private :: icRadius 00160 00161 00164 integer :: grid_EffSize 00165 ! f2003: integer, private :: grid_EffSize 00166 00169 integer :: grid_margine = 1 00170 ! f2003: integer, private :: grid_margine = 1 00171 00174 integer :: idMin 00175 ! f2003: integer, private :: idMin 00176 00179 integer :: idMax 00180 ! f2003: integer, private :: idMax 00181 00184 integer :: idEfMin = 1 00185 ! f2003: integer, private :: idEfMin = 1 00186 00189 integer :: idEfMax 00190 ! f2003: integer, private :: idEfMax 00191 00192 end type IcGrid2D_FVM 00193 00194 00195 ! 公開関数 00196 ! Public procedures 00197 ! 00198 public :: IcGrid2D_FVM_Init, IcGrid2D_FVM_Final 00199 00200 public :: malloc_GP_array, malloc_CV_array 00201 public :: generate_CV5_GPindex, generate_CV6_GPindex, check_pole 00202 00203 public :: set_rcregion_bounds, paste_margin_width 00204 public :: get_EffSize_Min, get_EffSize_Max, get_IdMin, get_IdMax 00205 public :: get_glevel, get_IcRadius 00206 00207 public :: calc_icosahedron_vertex 00208 00209 contains 00210 00211 ! 00232 subroutine IcGrid2D_FVM_Init( & 00233 & self, & ! (inout) 00234 & glevel, ic_radius & ! (in) 00235 ) 00236 00237 ! 宣言文 ; Declaration statements 00238 ! 00239 type(IcGrid2D_FVM), intent(inout) :: self 00240 integer, intent(in) :: glevel 00241 real(DP), intent(in) :: ic_radius 00242 00243 ! 実行文 ; Executable statements 00244 ! 00245 self%h_glevel = glevel 00246 self%icRadius = ic_radius 00247 00248 ! 指定された水平格子の分割レベルに対応する要素数を計算する. 00249 self%grid_EffSize = 2**self%h_glevel + 1 00250 00251 ! IcGrid2D_FVM クラスの配列変数のインデックスののりしろを含めた下限・上限, 00252 ! 物理的に意味のある下限・上限を求める. 00253 self%idEfMax = self%idEfMin + self%grid_EffSize - 1 00254 self%idMin = self%idEfMin - self%grid_margine 00255 self%idMax = self%idEfMax + self%grid_margine 00256 00257 end subroutine IcGrid2D_FVM_Init 00258 00259 ! 00270 subroutine malloc_GP_array( & 00271 & self & ! (inout) 00272 & ) 00273 00274 ! 宣言文 ; Declaration statements 00275 ! 00276 type(IcGrid2D_FVM), intent(inout) :: self ! IcGrid2D_FVM クラスのインスタンスの参照. 00277 00278 ! 実行文 ; Executable statements 00279 ! 00280 00281 !write(*,*) 'allocate grid point array' 00282 allocate(self%rcs_AGrid(RC_REGIONS_NUM, self%idMin:self%idMax, self%idMin:self%idMax, 3)) 00283 00284 end subroutine malloc_GP_array 00285 00286 ! 00298 subroutine malloc_CV_array( & 00299 & self & ! (inout) 00300 & ) 00301 00302 ! 宣言文 ; Declaration statements 00303 ! 00304 type(IcGrid2D_FVM), intent(inout) :: self 00305 00306 ! 実行文 ; Executable statements 00307 ! 00308 00309 !write(*,*) 'allocate controll volume array' 00310 allocate( self%rcs_CV(RC_REGIONS_NUM, self%idEfMin:self%idEfMax, self%idEfMin:self%idEfMax, CELL_POINTS_NUM, 3) ) 00311 allocate( self%rcs_CVArea(RC_REGIONS_NUM, self%idEfMin:self%idEfMax, self%idEfMin:self%idEfMax) ) 00312 00313 end subroutine malloc_CV_array 00314 00315 ! 00330 subroutine IcGrid2D_FVM_Final( & 00331 & self & ! (inout) 00332 & ) 00333 00334 ! 宣言文 ; Declaration statements 00335 ! 00336 type(IcGrid2D_FVM), intent(inout) :: self 00337 00338 ! 実行文 ; Executable statements 00339 ! 00340 00341 ! 本モジュールにおいて, 動的に確保した配列を解放する. 00342 deallocate(self%rcs_AGrid) 00343 deallocate(self%rcs_CV) 00344 deallocate(self%rcs_CVArea) 00345 00346 end subroutine IcGrid2D_FVM_Final 00347 00348 ! 00361 function get_EffSize_Min( & 00362 & self & ! (in) 00363 & ) result(gridEffMin) 00364 00365 type(IcGrid2D_FVM), intent(in) :: self 00366 integer gridEffMin 00367 00368 gridEffMin = self%idEfMin 00369 00370 end function get_EffSize_Min 00371 00372 ! 00385 function get_EffSize_Max( & 00386 & self & ! (in) 00387 & ) result(gridEffMax) 00388 00389 type(IcGrid2D_FVM), intent(in) :: self 00390 integer gridEffMax 00391 00392 gridEffMax = self%idEfMax 00393 00394 end function get_EffSize_Max 00395 00396 ! 00410 function get_IdMax( & 00411 & self & ! (in) 00412 & ) result(id_Max) 00413 00414 type(IcGrid2D_FVM), intent(in) :: self 00415 integer id_Max 00416 00417 id_Max = self%idMax 00418 00419 end function get_IdMax 00420 00421 ! 00434 function get_IdMin( & 00435 & self & ! (in) 00436 & ) result(id_Min) 00437 00438 type(IcGrid2D_FVM), intent(in) :: self 00439 integer id_Min 00440 00441 id_Min = self%idMin 00442 00443 end function get_IdMin 00444 00445 ! 00458 function get_IcRadius( & 00459 & self & ! (in) 00460 & ) result(ic_radius) 00461 00462 type(IcGrid2D_FVM), intent(in) :: self 00463 real(DP) ic_radius 00464 00465 ic_radius = self%icRadius 00466 end function get_IcRadius 00467 00468 ! 00481 function get_glevel( & 00482 & self & ! (in) 00483 & ) result(glevel) 00484 00485 type(IcGrid2D_FVM), intent(in) :: self ! IcGrid2D_FVM クラスのインスタンスの参照 00486 integer glevel 00487 00488 glevel = self%h_glevel 00489 end function get_glevel 00490 00491 ! 00512 subroutine set_rcregion_bounds( & 00513 & self, & ! (inout) 00514 & rcID, top, left, bottom, right & ! (in) 00515 & ) 00516 00517 ! 宣言文 ; Declaration statements 00518 ! 00519 type(IcGrid2D_FVM), intent(inout) :: self 00520 integer, intent(in) :: rcID 00521 real(DP), intent(in) :: top(3) 00522 real(DP), intent(in) :: left(3) 00523 real(DP), intent(in) :: bottom(3) 00524 real(DP), intent(in) :: right(3) 00525 00526 self%rcs_AGrid(rcID, self%idEfMin, self%idEfMin, :) = top 00527 self%rcs_AGrid(rcID, self%idEfMax, self%idEfMin, :) = left 00528 self%rcs_AGrid(rcID, self%idEfMax, self%idEfMax, :) = bottom 00529 self%rcs_AGrid(rcID, self%idEfMin, self%idEfMax, :) = right 00530 00531 end subroutine set_rcregion_bounds 00532 00533 ! 00544 subroutine paste_margin_width( & 00545 & self & ! (inout) 00546 & ) 00547 00548 ! 宣言文 ; Declaration statements 00549 ! 00550 type(IcGrid2D_FVM), intent(inout) :: self 00551 00552 ! 作業変数 00553 ! Work variables 00554 ! 00555 integer :: idlist1(7) = (/ 5, 1, 2, 3, 4, 5, 1 /) 00556 integer :: idlist2(7) = (/ 10, 6, 7, 8, 9, 10, 6 /) 00557 integer :: EMin, EMax, marginw 00558 integer i, k 00559 00560 ! 実行文 ; Executable statements 00561 ! 00562 00563 EMin = get_EffSize_Min(self) 00564 EMax = get_EffSize_Max(self) 00565 marginw = EMin - get_IdMin(self) 00566 00567 do i=1, RC_REGIONS_NUM/2 00568 do k=1, marginw 00569 ! 北半球側の矩形領域ののりしろに値をコピーする. 00570 self%rcs_AGrid(i, EMin:EMax, EMin-k, :) = self%rcs_AGrid(idlist1(i), EMin+k, EMin:EMax, :) 00571 self%rcs_AGrid(i, EMax+k, EMin:EMax, :) = self%rcs_AGrid(idlist2(i), EMin+k, EMin:EMax, :) 00572 self%rcs_AGrid(i, EMin:EMax, EMax+k, :) = self%rcs_AGrid(idlist2(i+1), EMin:EMax, EMin+k, :) 00573 self%rcs_AGrid(i, EMin-k, EMin:EMax, :) = self%rcs_AGrid(idlist1(i+2), EMin:EMax, EMin+k, :) 00574 00575 ! 南半球側の矩形領域ののりしろに値をコピーする. 00576 self%rcs_AGrid(i+5, EMin:EMax, EMin-k, :) = self%rcs_AGrid(idlist1(i+1), EMin:EMax, EMax-k, :) 00577 self%rcs_AGrid(i+5, EMax+k, EMin:EMax, :) = self%rcs_AGrid(idlist2(i), EMin:EMax, EMax-k, :) 00578 self%rcs_AGrid(i+5, EMin:EMax, EMax+k, :) = self%rcs_AGrid(idlist2(i+2), EMax-k, EMin:EMax, :) 00579 self%rcs_AGrid(i+5, EMin-k, EMin:EMax, :) = self%rcs_AGrid(idlist1(i+2), EMax-k, EMin:EMax, :) 00580 end do 00581 end do 00582 00583 end subroutine paste_margin_width 00584 00585 ! 00608 function generate_CV5_GPindex( & 00609 & self, & ! (inout) 00610 & GP_i, GP_j, rcID & ! (in) 00611 & ) result(CV_GPId) 00612 00613 ! 宣言文 ; Declaration statements 00614 ! 00615 type(IcGrid2D_FVM), intent(inout) :: self 00616 integer, intent(in) :: GP_i 00617 integer, intent(in) :: GP_j 00618 integer, intent(in) :: rcID 00619 integer CV_GPId(5,3,2) 00620 00621 ! 作業変数 00622 ! Work variables 00623 ! 00624 integer :: EMin, EMax 00625 integer i,j 00626 integer :: CV_GPIdShape(2) = (/ 3, 2 /) 00627 00628 ! 実行文 ; Executable statements 00629 ! 00630 00631 EMin = get_EffSize_Min(self) 00632 EMax = get_EffSize_Max(self) 00633 00634 i = GP_i; j = GP_j 00635 if( rcID <= 5 )then 00636 ! 北半球に対して 00637 if( i==EMax .and. j == EMax ) then 00638 ! 下端特異点 00639 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00640 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00641 CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00642 CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00643 CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00644 else if( j == EMin ) then 00645 ! 左端特異点 00646 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00647 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00648 CV_GPId(3,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape ) 00649 CV_GPId(4,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape ) 00650 CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00651 return 00652 else if( j == EMax ) then 00653 ! 右端特異点 00654 CV_GPId(1,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape ) 00655 CV_GPId(2,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape ) 00656 CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00657 CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00658 CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00659 return 00660 end if 00661 else 00662 ! 南半球に対して 00663 if ( i==EMin .and. j==EMin ) then 00664 ! 上端特異点 00665 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00666 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00667 CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00668 CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00669 CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00670 return 00671 else if( j == EMin ) then 00672 ! 左端特異点 00673 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00674 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00675 CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00676 CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape ) 00677 CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape ) 00678 return 00679 else if( j == EMax ) then 00680 ! 右端特異点 00681 CV_GPId(1,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00682 CV_GPId(2,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00683 CV_GPId(3,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00684 CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape ) 00685 CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape ) 00686 return 00687 end if 00688 end if 00689 00690 end function generate_CV5_GPindex 00691 00692 ! 00711 function generate_CV6_GPindex( & 00712 & self, & ! (inout) 00713 & GP_i, GP_j, rcID & ! (in) 00714 & ) result(CV_GPId) 00715 00716 ! 宣言文 ; Declaration statements 00717 ! 00718 type(IcGrid2D_FVM), intent(inout) :: self 00719 integer, intent(in) :: GP_i 00720 integer, intent(in) :: GP_j 00721 integer, intent(in) :: rcID 00722 integer :: CV_GPId(6,3,2) 00723 00724 ! 作業変数 00725 ! Work variables 00726 ! 00727 integer :: i,j, EMin, EMax 00728 integer :: CV_GPIdShape(2) = (/ 3, 2 /) 00729 00730 ! 実行文 ; Executable statements 00731 ! 00732 00733 EMin = get_EffSize_Min(self) 00734 EMax = get_EffSize_Max(self) 00735 00736 i = GP_i; j = GP_j 00737 00738 ! 先に, インデックス i,j が矩形領域の境界に当たる場合を処理する. 00739 if( rcID <= 5 )then 00740 ! 北半球の矩形領域に対して 00741 if( i == 1 ) then 00742 ! 右上端 00743 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j /), CV_GPIdShape ) 00744 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape ) 00745 CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape ) 00746 CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00747 CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00748 CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00749 return 00750 else if( j == 1 ) then 00751 ! 左上端 00752 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00753 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00754 CV_GPId(3,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape ) 00755 CV_GPId(4,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape ) 00756 CV_GPId(5,:,:) = reshape( (/ i,i,i+1, j,j-1,j /), CV_GPIdShape ) 00757 CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00758 return 00759 end if 00760 00761 else 00762 ! 南半球の矩形領域に対して 00763 if( i == EMax ) then 00764 ! 左下端 00765 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00766 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00767 CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00768 CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j /), CV_GPIdShape ) 00769 CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape ) 00770 CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape ) 00771 return 00772 else if( j == EMax ) then 00773 ! 右下端 00774 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j /), CV_GPIdShape ) 00775 CV_GPId(2,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00776 CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00777 CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00778 CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape ) 00779 CV_GPId(6,:,:) = reshape( (/ i,i,i+1, j,j+1,j+1 /), CV_GPIdShape ) 00780 return 00781 end if 00782 end if 00783 00784 ! インデックス i,j が矩形領域の内部(境界以外)に対応する場合 00785 CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 00786 CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 00787 CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 00788 CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 00789 CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 00790 CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 00791 00792 end function generate_CV6_GPindex 00793 00794 ! 00819 function check_pole(self, rcID, i, j) result (pole_flag) 00820 00821 ! 宣言文 ; Declaration statements 00822 ! 00823 type(IcGrid2D_FVM), intent(inout) :: self 00824 integer, intent(in) :: rcID 00825 integer, intent(in) :: i 00826 integer, intent(in) :: j 00827 integer pole_flag 00828 00829 if( rcID <= 5 .and. i==self%idEfMin .and. j== self%idEfMin ) then 00830 pole_flag = NORTH_POLE_FLAG 00831 else if( rcID > 5 .and. i==self%idEfMax .and. j== self%idEfMax ) then 00832 pole_flag = SOUTH_POLE_FLAG 00833 else 00834 pole_flag = NOT_POLE_FLAG 00835 end if 00836 00837 end function check_pole 00838 00839 ! 00851 subroutine calc_icosahedron_vertex(orth_icvertex) 00852 00853 ! 宣言文 ; Declaration statements 00854 ! 00855 00856 real(DP), intent(inout) :: orth_icvertex(ICOSAHEDRON_VERTEX_NUM, 3) 00857 00858 ! 作業変数 00859 ! Work variable 00860 ! 00861 real(DP) :: unit_radius = 1.0d0 00862 real(DP) :: del_lambda = 2.0d0 * PI / 5.0d0 00863 real(DP) :: pos2_theta 00864 real(DP) :: pos2_lambda 00865 integer i 00866 00867 ! 実行文 ; Executable statements 00868 ! 00869 00870 ! 頂点1, 12 を z 軸上に配置する. 00871 orth_icvertex(1, :) = (/ 0.0d0, 0.0d0, unit_radius /) 00872 orth_icvertex(12, :) = (/ 0.0d0, 0.0d0, - unit_radius /) 00873 00874 ! 頂点 2 の緯度・経度を設定する. 00875 pos2_lambda = 0.0d0 00876 pos2_theta = 2.0d0 * dasin( 1.0d0 / (2.0d0 * dsin(0.2d0 * PI)) ) - PI/2.0d0 00877 00878 ! 地理座標系における頂点 2, 7 の座標をデカルト座標系における座標へ変換する. 00879 ! なお, 頂点 7 は頂点 2 と xy 平面を挟んで反対側に存在する. 00880 orth_icvertex(2, :) = geo_to_orth_pos( (/ pos2_lambda, pos2_theta, unit_radius /) ) 00881 orth_icvertex(7, :) = geo_to_orth_pos( (/ pos2_lambda + del_lambda/2.0d0, - pos2_theta, unit_radius /) ) 00882 00883 ! z 軸を軸にして, 座標が既知である頂点を del_lambda づつ回転させて, 残りの頂点座標を計算する. 00884 do i=3, 6 00885 orth_icvertex(i, :) = rotateZ( orth_icvertex(i-1, :), del_lambda ) 00886 orth_icvertex(i+5, :) = rotateZ( orth_icvertex( (i+5)-1, : ), del_lambda ) 00887 end do 00888 00889 end subroutine calc_icosahedron_vertex 00890 00891 end module IcGrid2D_FVM_Manager