IGMBaseLib 1.0
|
00001 00020 module Derivate_Field_IcGrid2D_Manager 00021 00022 ! モジュール引用 ; Use statements 00023 ! 00024 00025 ! 種類型パラメータ 00026 ! Kind type parameter 00027 ! 00028 use dc_types, only: DP ! 倍精度実数型. Double precision 00029 00030 ! 球面三角法 00031 ! Spherical trigonometry 00032 ! 00033 use igmcore_spherical_trigonometry, only: & 00034 & spherical_tri_area 00035 00036 00037 ! 一つの格子点に対する微分演算子の提供 00038 ! Providing with the differential operation for a single grid-point 00039 ! 00040 use FVM_HDiff_operator, only: & 00041 & eval_div_operator, eval_rot_operator, eval_grad_operator 00042 00043 ! 正二十面体格子データの管理 00044 ! Manager of icosahedral grid data 00045 ! 00046 use IcGrid2D_FVM_Manager, only: & 00047 & IcGrid2D_FVM, & 00048 & paste_margin_width, generate_CV5_GPindex, generate_CV6_GPindex, & 00049 & get_IdMin, get_EffSize_Min, get_EffSize_Max, get_IcRadius, & 00050 & RC_REGIONS_NUM 00051 00052 ! 正二十面体格子上の物理場データの管理 00053 ! Manager of the physical field data on icosahedral grid 00054 ! 00055 use Field_IcGrid2D_Manager, only: & 00056 & Field_IcGrid2D 00057 00058 ! OpenMP 00059 ! 00060 ! 00061 use omp_lib 00062 00063 ! 宣言文 ; Declaration statements 00064 ! 00065 implicit none 00066 private 00067 00068 ! 公開変数 00069 ! Public variables 00070 ! 00071 00072 ! 00079 type, public :: Derivate_Field_IcGrid2D 00080 00083 type(IcGrid2D_FVM), pointer :: icgrid 00084 ! f2003: type(IcGrid2D_FVM), pointer, private :: icgrid 00085 00086 00093 real(DP), pointer :: rcs_areaGPP(:,:,:,:,:) 00094 ! f2003: real(DP), allocatable, private :: rcs_areaGPP(:,:,:,:,:) 00095 00098 integer, pointer :: rcU_CVGPindex(:,:,:,:,:) 00099 ! f2003: integer, allocatable, private :: rcU_CVGPindex(:,:,:,:,:) 00100 00103 integer, pointer :: rcL_CVGPindex(:,:,:,:,:) 00104 ! integer, allocatable, private :: rcL_CVGPindex(:,:,:,:,:) 00105 00106 end type Derivate_Field_IcGrid2D 00107 00108 ! 公開メンバ関数 00109 ! Public memeber function 00110 ! 00111 public :: Derivate_Field_IcGrid2D_Init, Derivate_Field_IcGrid2D_Final 00112 public :: divergence_op, vertical_curl_op, gradient_op 00113 00114 contains 00115 00116 ! 00129 subroutine Derivate_Field_IcGrid2D_Init( & 00130 & self, & ! (inout) 00131 & icgrid & ! (in) 00132 & ) 00133 00134 ! 宣言文 ; Declaration statements 00135 ! 00136 type(Derivate_Field_IcGrid2D), intent(inout) :: self 00137 type(IcGrid2D_FVM), intent(in), target :: icgrid 00138 00139 ! 作業変数 00140 ! Work variables 00141 ! 00142 integer :: EMin, EMax 00143 00144 ! 実行文 ; Declaration statements 00145 ! 00146 00147 self%icgrid => icgrid 00148 EMin = get_EffSize_Min(self%icgrid) 00149 EMax = get_EffSize_Max(self%icgrid) 00150 00151 ! Derivate_Field_IcGrid2D クラスが作業変数として保持する配列を確保する. 00152 allocate(self%rcs_areaGPP(RC_REGIONS_NUM, EMin:EMax, EMin:EMax, 6, 3)) 00153 allocate(self%rcU_CVGPindex(EMin:EMax, EMin:EMax, 6,3,2 )) 00154 allocate(self%rcL_CVGPindex(EMin:EMax, EMin:EMax, 6,3,2 )) 00155 00156 ! コントールボリュームの頂点座標を計算するために必要な格子点のインデックスを, あらかじめ計算しておく. 00157 call calc_CVGPindex(self) 00158 ! コントールボリュームの頂点上での物理場の値を計算する際に必要となる, 面積データをあらかじめ計算しておく. 00159 call calc_GPP_area_weight(self) 00160 00161 end subroutine Derivate_Field_IcGrid2D_Init 00162 00163 ! 00174 subroutine Derivate_Field_IcGrid2D_Final( & 00175 & self & ! (inout) 00176 & ) 00177 00178 ! 宣言文 ; Declaration statements 00179 ! 00180 type(Derivate_Field_IcGrid2D), intent(inout) :: self 00181 00182 ! 実行文 ; Declaration statements 00183 ! 00184 00185 deallocate(self%rcs_areaGPP) 00186 deallocate(self%rcU_CVGPindex) 00187 deallocate(self%rcL_CVGPindex) 00188 00189 end subroutine Derivate_Field_IcGrid2D_Final 00190 00191 ! 00206 subroutine divergence_op( & 00207 & self, & ! (inout) 00208 & vector_field, & ! (in) 00209 & ret_scalar_field & ! (inout) 00210 & ) 00211 00212 ! 宣言文 ; Declaration statements 00213 ! 00214 type(Derivate_Field_IcGrid2D), intent(inout), target :: self 00215 type(Field_IcGrid2D), intent(in) :: vector_field 00216 type(Field_IcGrid2D), intent(inout) :: ret_scalar_field 00217 00218 ! 作業変数 00219 ! Work variables 00220 ! 00221 integer :: idMin 00222 00223 ! 実行文 ; Executable statements 00224 ! 00225 00226 idMin = get_IdMin(self%icgrid) 00227 00228 call eval_diffOptr( & 00229 & self=self, & 00230 & diffOptr_fnc=eval_div_operator, rcs_GP_val=vector_field%val, val_dim=3, & ! (in) 00231 & ret_rcs_val=ret_scalar_field%val, & ! (inout) 00232 & ret_dim=1, idMin=idMin ) ! (in) 00233 00234 end subroutine divergence_op 00235 00236 ! 00251 subroutine vertical_curl_op( & 00252 & self, & ! (inout) 00253 & vector_field, & ! (in) 00254 & ret_scalar_field & ! (inout) 00255 & ) 00256 00257 ! 宣言文 ; Declaration statements 00258 ! 00259 type(Derivate_Field_IcGrid2D), intent(inout), target :: self 00260 type(Field_IcGrid2D), intent(in) :: vector_field 00261 type(Field_IcGrid2D), intent(inout) :: ret_scalar_field 00262 00263 ! 作業変数 00264 ! Work variables 00265 ! 00266 integer :: idMin 00267 00268 ! 実行文 ; Executable statements 00269 ! 00270 00271 idMin = get_IdMin(self%icgrid) 00272 00273 call eval_diffOptr( & 00274 & self=self, & 00275 & diffOptr_fnc=eval_rot_operator, rcs_GP_val=vector_field%val, val_dim=3, & ! (in) 00276 & ret_rcs_val=ret_scalar_field%val, & ! (inout) 00277 & ret_dim=1, idMin=idMin ) ! (in) 00278 00279 end subroutine vertical_curl_op 00280 00281 ! 00296 subroutine gradient_op( & 00297 & self, & ! (inout) 00298 & scalar_field, & ! (in) 00299 & ret_vector_field & ! (inout) 00300 & ) 00301 00302 ! 宣言文 ; Declaration statements 00303 ! 00304 type(Derivate_Field_IcGrid2D), intent(inout), target :: self 00305 type(Field_IcGrid2D), intent(in) :: scalar_field 00306 type(Field_IcGrid2D), intent(inout) :: ret_vector_field 00307 00308 ! 作業変数 00309 ! Work variables 00310 ! 00311 integer :: idMin 00312 00313 ! 実行文 ; Executable statements 00314 ! 00315 00316 idMin = get_IdMin(self%icgrid) 00317 00318 call eval_diffOptr( & 00319 & self=self, & 00320 & diffOptr_fnc=eval_grad_operator, rcs_GP_val=scalar_field%val, val_dim=1, & ! (in) 00321 & ret_rcs_val=ret_vector_field%val, & ! (inout) 00322 & ret_dim=3, idMin=idMin ) ! (in) 00323 00324 end subroutine gradient_op 00325 00326 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!1 00327 ! 00328 ! モジュール外非公開手続き 00329 ! 00330 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00331 00332 ! 00360 subroutine eval_diffOptr( & 00361 & self, & ! (inout) 00362 & diffOptr_fnc, rcs_GP_val, val_dim, & ! (in) 00363 & ret_rcs_val, & ! (inout) 00364 & ret_dim, idMin & ! (in) 00365 & ) 00366 00367 ! 宣言文 ; Declaration statements 00368 ! 00369 type(Derivate_Field_IcGrid2D), intent(inout), target :: self 00370 integer, intent(in) :: idMin 00371 integer, intent(in) :: val_dim 00372 integer, intent(in) :: ret_dim 00373 00374 interface 00375 !! 00376 !import DP 00377 00378 ! 00379 function diffOptr_fnc( & 00380 & val_GP0, val_CV, GP0, CV, & ! (in) 00381 & CVArea, ic_radius, CV_num, val_dim, ret_dim & ! (in) 00382 & ) result(val) 00383 00384 use dc_types, only: DP 00385 implicit none 00386 00387 ! 宣言文 ; Declaration statements 00388 ! 00389 00390 integer, intent(in) :: val_dim 00391 integer, intent(in) :: ret_dim 00392 integer, intent(in) :: CV_num 00393 real(DP), intent(in) :: val_GP0(val_dim) 00394 real(DP), intent(in) :: GP0(3) 00395 real(DP), intent(in) :: val_CV(CV_num, val_dim) 00396 real(DP), intent(in) :: CV(CV_num, val_dim) 00397 real(DP), intent(in) :: CVArea 00398 real(DP), intent(in) :: ic_radius 00399 real(DP) val(ret_dim) 00400 00401 end function diffOptr_fnc 00402 end interface 00403 00404 real(DP), intent(in) :: rcs_GP_val(:,idMin:,idMin:,:) 00405 real(DP), intent(inout) :: ret_rcs_val(:,idMin:,idMin:,:) 00406 00407 ! 作業変数 00408 ! Work variables 00409 ! 00410 real(DP) :: ic_radius 00411 real(DP), allocatable :: CV_val(:,:,:,:,:) 00412 integer :: rcID, i,j 00413 integer :: EMin, EMax 00414 integer :: CV_num 00415 00416 ! 実行文 ; Executable statements 00417 ! 00418 00419 EMin = get_EffSize_Min(self%icgrid) 00420 EMax = get_EffSize_Max(self%icgrid) 00421 ic_radius = get_IcRadius(self%icgrid) 00422 00423 ! コントロールボリュームの頂点上の物理場の値を保持するための配列のメモリを確保する. 00424 allocate(CV_Val(RC_REGIONS_NUM, EMin:EMax, EMin:EMax, 6, 3)) 00425 00426 ! 00427 do rcID=1, RC_REGIONS_NUM 00428 !$omp parallel do private(CV_num) 00429 do j=EMin, EMax 00430 do i=EMin,EMax 00431 00432 if ( ( i==EMin .and. j==EMin ) .or. ( i==EMax .and. j==EMin ) .or. & 00433 ( i==EMin .and. j==EMax ) .or. ( i==EMax .and. j==EMax ) ) then 00434 ! 特異格子点のコントールボリュームの頂点数は 5 個 00435 CV_num = 5 00436 else 00437 ! 通常ののコントールボリュームの頂点数は 6 個 00438 CV_num = 6 00439 end if 00440 00441 ! 格子点(rcID,i,j) のコントールボリュームの頂点上における物理場の値を計算する. 00442 call calc_CV_val( & 00443 & self, & 00444 & CV_val(rcID, i, j, 1:CV_num, 1:val_dim), & ! (inout) 00445 & rcs_GP_val(:,:,:,1:val_dim), & ! (in) 00446 & rcID, i,j, val_dim, CV_num, idMin, EMin, EMax ) ! (in) 00447 00448 ! 動的に接続された微分演算オブジェクトに対して微分演算を実行し, 格子点(rcID,i,j)上の微分値を計算する. 00449 ! また, その結果を保存する. 00450 ret_rcs_val(rcID, i,j, 1:ret_dim) = diffOptr_fnc( & 00451 & val_GP0=rcs_GP_val(rcID,i,j,1:val_dim), val_CV=CV_val(rcID,i,j, 1:CV_num,1:val_dim), & ! (in) 00452 & GP0=self%icgrid%rcs_AGrid(rcID,i,j,:), CV=self%icgrid%rcs_CV(rcID,i,j,1:CV_num,:), & ! (in) 00453 & CVArea=self%icgrid%rcs_CVArea(rcID,i,j), & ! (in) 00454 & ic_radius=ic_radius, CV_num=CV_num, val_dim=val_dim, ret_dim=ret_dim ) ! (in) 00455 00456 end do 00457 end do 00458 end do 00459 00460 end subroutine eval_diffOptr 00461 00462 ! 00493 subroutine calc_CV_val( self, CV_val, rcs_GP_val, rcID, gp_i, gp_j, val_dim, CV_num, idMin, EMin, EMax) 00494 00495 ! 宣言文 ; Declaration statement 00496 ! 00497 type(Derivate_Field_IcGrid2D), intent(inout) :: self 00498 integer, intent(in) :: val_dim, CV_num, idMin, EMin, EMax 00499 real(DP), intent(inout) :: CV_val(CV_num, val_dim) 00500 real(DP), intent(in) :: rcs_GP_val(:,idMin:,idMin:,:) 00501 integer, intent(in) :: rcID, gp_i, gp_j 00502 00503 ! 作業変数 00504 ! Work variable 00505 ! 00506 integer :: cvID, k 00507 real(DP) :: GP_val3(3,val_dim) 00508 real(DP) :: area3(3) 00509 00510 ! 実行文 ; Executable statement 00511 ! 00512 00513 do cvID=1, CV_num 00514 00515 ! 矩形領域 rcID 内の格子点 (gp_i, gp_j) に付随するコントロールボリュームの一頂点 cvID と近傍の格子点が作る 3 つの三角形の 00516 ! 面積を取得する. 00517 area3(:) = self%rcs_areaGPP(rcID,gp_i,gp_j,cvID, :) 00518 00519 ! 00520 if( rcID <= 5 ) then 00521 if ( gp_i==EMin .and. gp_j==EMin )then 00522 GP_val3(1,:) = rcs_GP_val(cvID, EMin, EMin, :); GP_val3(2,:) = rcs_GP_val(cvID, EMin+1, EMin, :); 00523 GP_val3(3,:) = rcs_GP_val(cvID, EMin, EMin+1, :); 00524 else 00525 do k=1,3 00526 GP_val3(k,:) = rcs_GP_val(rcID, self%rcU_CVGPindex(gp_i,gp_j,cvID,k,1), self%rcU_CVGPindex(gp_i,gp_j,cvID,k,2), :) 00527 end do 00528 end if 00529 00530 else 00531 if ( gp_i==EMax .and. gp_j==EMax )then 00532 GP_val3(1,:) = rcs_GP_val(11-cvID, EMax, EMax, :); GP_val3(2,:) = rcs_GP_val(11-cvID, EMax-1, EMax, :); 00533 GP_val3(3,:) = rcs_GP_val(11-cvID, EMax, EMax-1, :); 00534 else 00535 do k=1,3 00536 GP_val3(k,:) = rcs_GP_val(rcID, self%rcL_CVGPindex(gp_i,gp_j,cvID,k,1), self%rcL_CVGPindex(gp_i,gp_j,cvID,k,2), :) 00537 end do 00538 end if 00539 00540 end if 00541 00542 ! コントロールボリュームの各頂点(点 G)上の物理場の値を, 近傍の格子点(点 P0,P1,P2)上の物理場の値から補間する. 00543 ! 各格子点の物理場の値は, コントロールボリュームの頂点と格子点が作る三角形(GP0P1, GP1P2, GP2P0)の面積によって重みを付けられる. 00544 CV_val(cvID,1:val_dim) = ( & 00545 & area3(1) * GP_val3(1,:) + area3(2) * GP_val3(2,:) + area3(3) * GP_val3(3,:) & 00546 & ) / ( area3(1) + area3(2) + area3(3) ) 00547 00548 end do 00549 00550 end subroutine calc_CV_val 00551 00552 ! 00563 subroutine calc_CVGPindex(self) 00564 00565 ! 宣言文 ; Declaration statement 00566 ! 00567 type(Derivate_Field_IcGrid2D), intent(inout) :: self 00568 00569 ! 作業変数 00570 ! Work variable 00571 ! 00572 integer :: i, j 00573 integer :: EMin, EMax, idMin 00574 00575 ! 実行文 ; Executable statement 00576 ! 00577 00578 EMin = get_EffSize_Min(self%icgrid) 00579 EMax = get_EffSize_Max(self%icgrid) 00580 idMin = get_IdMin(self%icgrid) 00581 00582 do j=EMin,EMax 00583 do i=EMin,EMax 00584 if ( ( i==EMax .and. j==EMin ) .or. ( i==EMin .and. j==EMax ) ) then 00585 ! 00586 ! 00587 self%rcU_CVGPindex(i,j,1:5,:,:) = & 00588 & generate_CV5_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=1 ) 00589 self%rcL_CVGPindex(i,j,1:5,:,:) = & 00590 & generate_CV5_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=6 ) 00591 else 00592 00593 ! 00594 self%rcU_CVGPindex(i,j,1:6,:,:) = & 00595 & generate_CV6_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=1 ) 00596 ! 00597 self%rcL_CVGPindex(i,j,1:6,:,:) = & 00598 & generate_CV6_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=6 ) 00599 00600 end if 00601 end do 00602 end do 00603 00604 self%rcU_CVGPindex(EMax,EMax,1:5,:,:) = & 00605 & generate_CV5_GPindex(self%icgrid, GP_i=EMax, GP_j=EMax, rcID=1 ) 00606 self%rcL_CVGPindex(EMin,EMin,1:5,:,:) = & 00607 & generate_CV5_GPindex(self%icgrid, GP_i=EMin, GP_j=EMin, rcID=6 ) 00608 00609 end subroutine calc_CVGPindex 00610 00611 ! 00622 subroutine calc_GPP_area_weight( & 00623 & self & ! (inout) 00624 & ) 00625 00626 ! 宣言文 ; Declaration statement 00627 ! 00628 type(Derivate_Field_IcGrid2D), intent(inout) :: self 00629 00630 ! 作業変数 00631 ! Work variables 00632 ! 00633 integer :: rcID, i, j, cvID 00634 integer :: EMin, EMax, idMin 00635 real(DP) :: ic_radius 00636 integer :: CV5_GPindx(5,3,2), CV6_GPindx(6,3,2) 00637 00638 ! 実行文 ; Executable statement 00639 ! 00640 00641 EMin = get_EffSize_Min(self%icgrid) 00642 EMax = get_EffSize_Max(self%icgrid) 00643 idMin = get_IdMin(self%icgrid) 00644 ic_radius = get_IcRadius(self%icgrid) 00645 00646 do rcID=1, RC_REGIONS_NUM 00647 do j=EMin,EMax 00648 do i=EMin,EMax 00649 if ( ( rcID > 5 .and. i==EMin .and. j==EMin ) .or. ( i==EMax .and. j==EMin ) .or. & 00650 ( i==EMin .and. j==EMax ) .or. ( rcID <= 5 .and. i==EMax .and. j==EMax ) ) then 00651 ! 00652 CV5_GPindx(:,:,:) = generate_CV5_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=rcID ) 00653 00654 do cvID=1,5 00655 self%rcs_areaGPP(rcID,i,j,cvID, 1:3) = calc_usual_area3( CV_GPindx=CV5_GPindx(cvID, :,:), & 00656 & g=self%icgrid%rcs_CV(rcID,i,j,cvID,:), & 00657 & rc_AGrid=self%icgrid%rcs_AGrid(rcID,:,:,:), ic_radius=ic_radius, idMin=idMin ) 00658 end do 00659 else if( rcID <= 5 .and. i==EMin .and. j==EMin ) then 00660 ! 00661 do cvID=1,5 00662 self%rcs_areaGPP(rcID,i,j,cvID, 1:3) = calc_area3_from_4pts( ic_radius=ic_radius, & 00663 & p1=self%icgrid%rcs_AGrid(cvID,EMin,EMin,:), & 00664 & p2=self%icgrid%rcs_AGrid(cvID,EMin+1,EMin,:), & 00665 & p3=self%icgrid%rcs_AGrid(cvID,EMin,EMin+1,:), & 00666 & g=self%icgrid%rcs_CV(rcID,EMin,EMin,cvID,:) ) 00667 end do 00668 00669 else if( rcID > 5 .and. i==EMax .and. j==EMax ) then 00670 ! 00671 do cvID=1,5 00672 self%rcs_areaGPP(rcID,i,j,cvID, 1:3) = calc_area3_from_4pts( ic_radius=ic_radius, & 00673 & p1=self%icgrid%rcs_AGrid(11-cvID,EMax,EMax,:), & 00674 & p2=self%icgrid%rcs_AGrid(11-cvID,EMax-1,EMax,:), & 00675 & p3=self%icgrid%rcs_AGrid(11-cvID,EMax,EMax-1,:), & 00676 & g=self%icgrid%rcs_CV(rcID,EMax,EMax,cvID,:) ) 00677 end do 00678 00679 else 00680 ! 00681 CV6_GPindx(:,:,:) = generate_CV6_GPindex(self%icgrid, GP_i=i, GP_j=j, rcID=rcID ) 00682 00683 do cvID=1,6 00684 self%rcs_areaGPP(rcID,i,j,cvID, 1:3) = calc_usual_area3( CV_GPindx=CV6_GPindx(cvID, :,:), & 00685 & g=self%icgrid%rcs_CV(rcID,i,j,cvID,:), & 00686 & rc_AGrid=self%icgrid%rcs_AGrid(rcID,:,:,:), ic_radius=ic_radius, idMin=idMin ) 00687 end do 00688 00689 end if 00690 end do 00691 end do 00692 end do 00693 00694 end subroutine calc_GPP_area_weight 00695 00696 ! 00717 function calc_usual_area3( & 00718 & CV_GPindx, g, rc_AGrid, ic_radius, idMin & ! (in) 00719 & ) result(area3) 00720 00721 ! 宣言文 ; Declaration statement 00722 ! 00723 integer, intent(in) :: idMin 00724 integer, intent(in) :: CV_GPindx(3,2) 00725 real(DP), intent(in) :: g(3) 00726 real(DP), intent(in) :: rc_AGrid(idMin:,idMin:,:) 00727 real(DP), intent(in) :: ic_radius 00728 00729 ! 作業変数 00730 ! Work variable 00731 ! 00732 real(DP) area3(3) 00733 00734 ! 実行文 ; Executable statement 00735 ! 00736 00737 area3(:) = calc_area3_from_4pts( & 00738 & p1=rc_AGrid(CV_GPindx(1,1), CV_GPindx(1,2), :), & 00739 & p2=rc_AGrid(CV_GPindx(2,1), CV_GPindx(2,2), :), & 00740 & p3=rc_AGrid(CV_GPindx(3,1), CV_GPindx(3,2), :), & 00741 & g=g(:), & 00742 & ic_radius=ic_radius ) 00743 00744 end function calc_usual_area3 00745 00746 ! 00767 function calc_area3_from_4pts(p1, p2, p3, g, ic_radius) result(area3) 00768 00769 ! 宣言文 ; Declaration statement 00770 ! 00771 real(DP), intent(in) :: p1(3), p2(3), p3(3) 00772 real(DP), intent(in) :: g(3) 00773 real(DP), intent(in) :: ic_radius 00774 real(DP) area3(3) 00775 00776 ! 実行文 ; Executable statement 00777 ! 00778 00779 area3(1) = spherical_tri_area(ic_radius, g, p2, p3) 00780 area3(2) = spherical_tri_area(ic_radius, g, p3, p1) 00781 area3(3) = spherical_tri_area(ic_radius, g, p1, p2) 00782 00783 end function calc_area3_from_4pts 00784 00785 end module Derivate_Field_IcGrid2D_Manager