102 mpp_clock_begin, mpp_clock_end, mpp_clock_sync, &
103 comm_tag_1, comm_tag_2, comm_tag_3, comm_tag_4, &
104 comm_tag_5, comm_tag_6, comm_tag_7, comm_tag_8, &
105 comm_tag_9, comm_tag_10
117 domainug, mpp_get_ug_compute_domains, &
118 mpp_get_ug_domains_index, mpp_get_ug_domain_grid_index, &
120 use constants_mod,
only: pi, radius
125 use stock_constants_mod,
only: istock_top, istock_bottom, istock_side, stock_names, &
127 use gradient_mod,
only: gradient_cubic
129 use fms2_io_mod,
only: fmsnetcdfdomainfile_t,
read_data, get_dimension_size
130 use fms2_io_mod,
only: get_variable_units, dimension_exists
131 use platform_mod,
only: r8_kind, i8_kind, fms_file_len
145 integer,
parameter :: FIRST_ORDER = 1
146 integer,
parameter :: SECOND_ORDER = 2
149 integer,
parameter :: max_fields = 100
156 logical :: debug_stocks = .false.
157 logical :: xgrid_clocks_on = .false.
158 logical :: monotonic_exchange = .false.
166 logical :: do_alltoall = .true.
167 logical :: do_alltoallv = .false.
172 monotonic_exchange,
nsubset, do_alltoall, do_alltoallv, &
178 real(r8_kind),
allocatable,
dimension(:,:) ::
area_atm_model, area_lnd_model, area_ocn_model
180 real(r8_kind),
allocatable,
dimension(:,:) ::
area_atm_sphere, area_lnd_sphere, area_ocn_sphere
215 module procedure put_side2_to_xgrid_ug
223 module procedure get_side2_from_xgrid_ug
224 module procedure get_side1_from_xgrid_ug
277 real(r8_kind) :: area
282 real(r8_kind) :: scale
288 real(r8_kind),
dimension(:,:),
pointer :: dx => null()
289 real(r8_kind),
dimension(:,:),
pointer :: dy => null()
290 real(r8_kind),
dimension(:,:),
pointer :: area => null()
291 real(r8_kind),
dimension(:),
pointer :: edge_w => null()
292 real(r8_kind),
dimension(:),
pointer :: edge_e => null()
293 real(r8_kind),
dimension(:),
pointer :: edge_s => null()
294 real(r8_kind),
dimension(:),
pointer :: edge_n => null()
295 real(r8_kind),
dimension(:,:,:),
pointer :: en1 => null()
296 real(r8_kind),
dimension(:,:,:),
pointer :: en2 => null()
297 real(r8_kind),
dimension(:,:,:),
pointer :: vlon => null()
298 real(r8_kind),
dimension(:,:,:),
pointer :: vlat => null()
304 character(len=3) :: id
306 logical :: on_this_pe
308 integer,
pointer,
dimension(:) :: pelist
312 integer,
pointer,
dimension(:) :: tile =>null()
313 integer,
pointer,
dimension(:) :: is =>null()
314 integer,
pointer,
dimension(:) :: ie =>null()
315 integer,
pointer,
dimension(:) :: js =>null()
316 integer,
pointer,
dimension(:) :: je =>null()
317 integer,
pointer :: is_me =>null()
318 integer,
pointer :: ie_me =>null()
319 integer,
pointer :: js_me =>null()
320 integer,
pointer :: je_me =>null()
329 integer,
pointer :: tile_me
333 real(r8_kind),
pointer,
dimension(:) :: lon =>null()
334 real(r8_kind),
pointer,
dimension(:) :: lat =>null()
335 real(r8_kind),
pointer,
dimension(:,:) :: geolon=>null()
336 real(r8_kind),
pointer,
dimension(:,:) :: geolat=>null()
337 real(r8_kind),
pointer,
dimension(:,:,:) :: frac_area =>null()
338 real(r8_kind),
pointer,
dimension(:,:) :: area =>null()
339 real(r8_kind),
pointer,
dimension(:,:) :: area_inv =>null()
346 integer :: size_repro
355 integer,
pointer :: ls_me =>null()
356 integer,
pointer :: le_me =>null()
357 integer,
pointer,
dimension(:) :: ls =>null(), le =>null()
358 integer,
pointer :: gs_me =>null(), ge_me =>null()
359 integer,
pointer,
dimension(:) :: gs =>null(), ge =>null()
360 integer,
pointer,
dimension(:) :: l_index =>null()
369 real(r8_kind) :: area
380 integer :: i, j, l, k, pos
381 real(r8_kind) :: area
390 integer :: buffer_pos
391 integer,
allocatable :: i(:)
392 integer,
allocatable :: j(:)
393 integer,
allocatable :: g(:)
394 integer,
allocatable :: xLoc(:)
395 integer,
allocatable :: tile(:)
396 real(r8_kind),
allocatable :: di(:)
397 real(r8_kind),
allocatable :: dj(:)
403 integer :: nsend, nrecv
404 integer :: sendsize, recvsize
405 integer,
pointer,
dimension(:) :: unpack_ind=>null()
406 type(overlap_type),
pointer,
dimension(:) :: send=>null()
407 type(overlap_type),
pointer,
dimension(:) :: recv=>null()
417 integer :: me, npes, root_pe
418 logical,
pointer,
dimension(:) :: your1my2 =>null()
421 logical,
pointer,
dimension(:) :: your2my1 =>null()
424 integer,
pointer,
dimension(:) :: your2my1_size=>null()
429 type (
grid_type),
pointer,
dimension(:) :: grids =>null()
434 type(
x1_type),
pointer,
dimension(:) :: x1 =>null()
435 type(
x1_type),
pointer,
dimension(:) :: x1_put =>null()
436 type(
x2_type),
pointer,
dimension(:) :: x2 =>null()
437 type(
x2_type),
pointer,
dimension(:) :: x2_get =>null()
439 integer,
pointer,
dimension(:) :: send_count_repro =>null()
440 integer,
pointer,
dimension(:) :: recv_count_repro =>null()
441 integer :: send_count_repro_tot
442 integer :: recv_count_repro_tot
445 integer,
pointer,
dimension(:) :: ind_get1 =>null()
446 integer,
pointer,
dimension(:) :: ind_put1 =>null()
456 #include<file_version.h>
458 real(r8_kind),
parameter :: eps = 1.0e-10_r8_kind
459 real(r8_kind),
parameter :: large_number = 1.e20_r8_kind
460 logical :: module_is_initialized = .false.
461 integer :: id_put_1_to_xgrid_order_1 = 0
462 integer :: id_put_1_to_xgrid_order_2 = 0
463 integer :: id_get_1_from_xgrid = 0
464 integer :: id_get_1_from_xgrid_repro = 0
465 integer :: id_get_2_from_xgrid = 0
466 integer :: id_put_2_to_xgrid = 0
467 integer :: id_setup_xmap = 0
468 integer :: id_load_xgrid1, id_load_xgrid2, id_load_xgrid3
469 integer :: id_load_xgrid4, id_load_xgrid5
470 integer :: id_load_xgrid, id_set_comm, id_regen, id_conservation_check
474 integer :: nnest=0, tile_nest, tile_parent
475 integer :: is_nest=0, ie_nest=0, js_nest=0, je_nest=0
476 integer :: is_parent=0, ie_parent=0, js_parent=0, je_parent=0
488 module procedure stock_move_ug_3d
510 logical function in_box(i, j, is, ie, js, je)
511 integer,
intent(in) :: i, j, is, ie, js, je
513 in_box = (i>=is) .and. (i<=ie) .and. (j>=js) .and. (j<=je)
522 integer,
intent(out) :: remap_method
525 integer :: iunit, ierr, io, out_unit
527 if (module_is_initialized)
return
528 module_is_initialized = .true.
530 read (input_nml_file, xgrid_nml, iostat=io)
538 if (
mpp_pe() == mpp_root_pe() )
write (iunit,nml=xgrid_nml)
543 'MPP_IO is no longer supported. Please remove use_mpp_io from namelists',&
551 remap_method = first_order
552 if( monotonic_exchange )
call error_mesg(
'xgrid_mod', &
553 'xgrid_nml monotonic_exchange must be .false. when interp_method = first_order', fatal)
554 write(out_unit,*)
"NOTE from xgrid_mod: use first_order conservative exchange"
556 if(monotonic_exchange)
then
557 write(out_unit,*)
"NOTE from xgrid_mod: use monotonic second_order conservative exchange"
559 write(out_unit,*)
"NOTE from xgrid_mod: use second_order conservative exchange"
561 remap_method = second_order
564 ' is not a valid namelist option', fatal)
567 if(xgrid_clocks_on)
then
568 id_put_1_to_xgrid_order_1 =
mpp_clock_id(
"put_1_to_xgrid_order_1", flags=mpp_clock_sync)
569 id_put_1_to_xgrid_order_2 =
mpp_clock_id(
"put_1_to_xgrid_order_2", flags=mpp_clock_sync)
570 id_get_1_from_xgrid =
mpp_clock_id(
"get_1_from_xgrid", flags=mpp_clock_sync)
571 id_get_1_from_xgrid_repro =
mpp_clock_id(
"get_1_from_xgrid_repro", flags=mpp_clock_sync)
572 id_get_2_from_xgrid =
mpp_clock_id(
"get_2_from_xgrid", flags=mpp_clock_sync)
573 id_put_2_to_xgrid =
mpp_clock_id(
"put_2_to_xgrid", flags=mpp_clock_sync)
574 id_setup_xmap =
mpp_clock_id(
"setup_xmap", flags=mpp_clock_sync)
577 id_conservation_check =
mpp_clock_id(
"conservation_check")
592 subroutine load_xgrid (xmap, grid, grid_file, grid1_id, grid_id, tile1, tile2, use_higher_order)
595 character(len=*),
intent(in) :: grid_file
596 character(len=3),
intent(in) :: grid1_id, grid_id
597 integer,
intent(in) :: tile1, tile2
598 logical,
intent(in) :: use_higher_order
600 integer,
pointer,
dimension(:) :: i1=>null(), j1=>null()
601 integer,
pointer,
dimension(:) :: i2=>null(), j2=>null()
602 real(r8_kind),
pointer,
dimension(:) :: di=>null(), dj=>null()
603 real(r8_kind),
pointer,
dimension(:) :: area =>null()
604 integer,
pointer,
dimension(:) :: i1_tmp=>null(), j1_tmp=>null()
605 integer,
pointer,
dimension(:) :: i2_tmp=>null(), j2_tmp=>null()
606 real(r8_kind),
pointer,
dimension(:) :: di_tmp=>null(), dj_tmp=>null()
607 real(r8_kind),
pointer,
dimension(:) :: area_tmp =>null()
608 integer,
pointer,
dimension(:) :: i1_side1=>null(), j1_side1=>null()
609 integer,
pointer,
dimension(:) :: i2_side1=>null(), j2_side1=>null()
610 real(r8_kind),
pointer,
dimension(:) :: di_side1=>null(), dj_side1=>null()
611 real(r8_kind),
pointer,
dimension(:) :: area_side1 =>null()
613 real(r8_kind),
allocatable,
dimension(:,:) :: tmp
614 real(r8_kind),
allocatable,
dimension(:) :: send_buffer, recv_buffer
615 type (
grid_type),
pointer,
save :: grid1 =>null()
616 integer :: l, ll, ll_repro, p, nxgrid, size_prev
618 integer :: size_repro, out_unit
619 logical :: scale_exist = .false.
620 logical :: is_distribute = .false.
621 real(r8_kind),
allocatable,
dimension(:) :: scale
622 real(r8_kind) :: garea
623 integer :: npes, isc, iec, nxgrid_local, pe, nxgrid_local_orig
624 integer :: nxgrid1, nxgrid2, nset1, nset2, ndivs, cur_ind
625 integer :: pos, nsend, nrecv, l1, l2, n, mypos
626 integer :: start(4), nread(4)
628 character(len=128) :: attvalue
629 integer,
dimension(0:xmap%npes-1) :: pelist
630 logical,
dimension(0:xmap%npes-1) :: subset_rootpe
631 integer,
dimension(0:xmap%npes-1) :: nsend1, nsend2, nrecv1, nrecv2
632 integer,
dimension(0:xmap%npes-1) :: send_cnt, recv_cnt
633 integer,
dimension(0:xmap%npes-1) :: send_buffer_pos, recv_buffer_pos
634 integer,
dimension(0:xmap%npes-1) :: ibegin, iend, pebegin, peend
635 integer,
dimension(2*xmap%npes) :: ibuf1, ibuf2
636 integer,
dimension(0:xmap%npes-1) :: pos_x, y2m1_size
637 integer,
allocatable,
dimension(:) :: y2m1_pe
638 integer,
pointer,
save :: iarray(:), jarray(:)
639 integer,
allocatable,
save :: pos_s(:)
640 integer,
pointer,
dimension(:) :: iarray2(:)=>null(), jarray2(:)=>null()
642 integer :: nxgrid1_old
644 type(fmsnetcdffile_t) :: fileobj
646 if(.not.
open_file(fileobj, grid_file,
'read' ))
then
647 call error_mesg(
'xgrid_mod(load_xgrid)',
'Error in opening file '//trim(grid_file), fatal)
650 scale_exist = .false.
651 grid1 => xmap%grids(1)
655 mypos =
mpp_pe()-mpp_root_pe()
657 call mpp_get_current_pelist(pelist)
659 if( npes .NE. pelist(npes-1) - pelist(0) + 1 )
then
660 print*,
"npes =", npes,
", pelist(npes-1)=", pelist(npes-1),
", pelist(0)=", pelist(0)
661 call error_mesg(
'xgrid_mod', .NE.
'npes pelist(npes-1) - pelist(0)', fatal)
664 select case(xmap%version)
667 if (dimension_exists(fileobj,
'i_'//lowercase(grid1_id)//
'X'//lowercase(grid_id)))
then
668 call get_dimension_size(fileobj,
'i_'//lowercase(grid1_id)//
'X'//lowercase(grid_id), nxgrid)
670 if(nxgrid .LE. 0)
return
674 if(nxgrid .LE. 0)
return
678 if(nxgrid > npes)
then
682 if(npes == ndivs)
then
683 p =
mpp_pe()-mpp_root_pe()
686 subset_rootpe(:) = .true.
691 if(pe == pebegin(n))
then
698 subset_rootpe(:) = .false.
701 if(pelist(n) == pebegin(cur_ind))
then
702 subset_rootpe(n) = .true.
704 if(cur_ind == ndivs)
exit
708 is_distribute = .true.
710 is_distribute = .false.
711 isc = 1; iec = nxgrid
716 if(use_higher_order)
then
720 if(scale_exist) nset2 = nset1 + 1
722 call mpp_clock_begin(id_load_xgrid1)
723 if(iec .GE. isc)
then
724 nxgrid_local = iec - isc + 1
725 allocate(i1_tmp(isc:iec), j1_tmp(isc:iec), i2_tmp(isc:iec), j2_tmp(isc:iec), area_tmp(isc:iec) )
726 if(use_higher_order)
allocate(di_tmp(isc:iec), dj_tmp(isc:iec))
730 select case(xmap%version)
732 start(1) = isc; nread(1) = nxgrid_local
733 allocate(tmp(nxgrid_local,1))
734 call read_data(fileobj,
'I_'//grid1_id//
'_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
735 i1_tmp = int(tmp(:,1))
736 call read_data(fileobj,
'J_'//grid1_id//
'_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
737 j1_tmp = int(tmp(:,1))
738 call read_data(fileobj,
'I_'//grid_id//
'_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
739 i2_tmp = int(tmp(:,1))
740 call read_data(fileobj,
'J_'//grid_id//
'_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
741 j2_tmp = int(tmp(:,1))
742 call read_data(fileobj,
'AREA_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
744 if(use_higher_order)
then
745 call read_data(fileobj,
'DI_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
747 call read_data(fileobj,
'DJ_'//grid1_id//
'x'//grid_id, tmp, corner=start, edge_lengths=nread)
752 nread(1) = 2; start(2) = isc; nread(2) = nxgrid_local
753 allocate(tmp(2, isc:iec))
754 call read_data(fileobj,
"tile1_cell", tmp, corner=start, edge_lengths=nread)
755 i1_tmp(isc:iec) = int(tmp(1, isc:iec))
756 j1_tmp(isc:iec) = int(tmp(2, isc:iec))
757 call read_data(fileobj,
"tile2_cell", tmp, corner=start, edge_lengths=nread)
758 i2_tmp(isc:iec) = int(tmp(1, isc:iec))
759 j2_tmp(isc:iec) = int(tmp(2, isc:iec))
760 if(use_higher_order)
then
761 call read_data(fileobj,
"tile1_distance", tmp, corner=start, edge_lengths=nread)
762 di_tmp(isc:iec) = tmp(1, isc:iec)
763 dj_tmp(isc:iec) = tmp(2, isc:iec)
766 start(1) = isc; nread(1) = nxgrid_local
768 allocate(tmp(isc:iec,1) )
770 call read_data(fileobj,
"xgrid_area", tmp(:,1:1), corner=start, edge_lengths=nread)
772 call get_variable_units(fileobj,
"xgrid_area", attvalue)
774 if( trim(attvalue) ==
'm2' )
then
775 garea = 4.0_r8_kind * pi * radius * radius;
776 area_tmp = tmp(:,1)/garea
777 else if( trim(attvalue) ==
'none' )
then
780 call error_mesg(
'xgrid_mod',
'In file '//trim(grid_file)//
', xgrid_area units = '// &
781 trim(attvalue)//
' should be "m2" or "none"', fatal)
786 if(grid1_id ==
'LND' .AND. grid_id ==
'OCN')
then
787 if(variable_exists(fileobj,
"scale"))
then
788 allocate(scale(isc:iec))
789 write(out_unit, *)
"NOTE from load_xgrid(xgrid_mod): field 'scale' exist in the file "// &
790 & trim(grid_file)//
", this field will be read and the exchange grid cell area will be"// &
791 &
" multiplied by scale"
792 call read_data(fileobj,
"scale", tmp, corner=start, edge_lengths=nread)
802 nxgrid_local_orig = nxgrid_local
803 allocate(i1(isc:iec), j1(isc:iec), i2(isc:iec), j2(isc:iec), area(isc:iec) )
804 if(use_higher_order)
allocate(di(isc:iec), dj(isc:iec))
810 if(grid1%tile(p) == tile1)
then
811 if(
in_box_nbr(i1_tmp(l), j1_tmp(l), grid1, p))
then
820 if(grid%tile(p) == tile2)
then
821 if (
in_box_nbr(i2_tmp(l), j2_tmp(l), grid, p))
then
827 area(pos) = area_tmp(l)
828 if(use_higher_order)
then
839 deallocate(i1_tmp, i2_tmp, j1_tmp, j2_tmp, area_tmp)
840 if(use_higher_order)
deallocate( di_tmp, dj_tmp)
842 if(iec .GE. isc)
then
843 nxgrid_local = iec - isc + 1
849 nxgrid_local_orig = 0
854 call mpp_clock_end(id_load_xgrid1)
856 if(is_distribute)
then
860 nsend1(:) = 0; nrecv1(:) = 0
861 nsend2(:) = 0; nrecv2(:) = 0
862 ibuf1(:)= 0; ibuf2(:)= 0
864 call mpp_clock_begin(id_load_xgrid2)
865 if(nxgrid_local>0)
then
866 allocate( send_buffer(nxgrid_local * (nset1+nset2)) )
869 send_buffer_pos(p) = pos
870 if(grid%tile(p) == tile2)
then
873 nsend2(p) = nsend2(p) + 1
874 send_buffer(pos+1) = real(i1(l), r8_kind)
875 send_buffer(pos+2) = real(j1(l), r8_kind)
876 send_buffer(pos+3) = real(i2(l), r8_kind)
877 send_buffer(pos+4) = real(j2(l), r8_kind)
878 send_buffer(pos+5) = area(l)
879 if(use_higher_order)
then
880 send_buffer(pos+6) = di(l)
881 send_buffer(pos+7) = dj(l)
883 if(scale_exist) send_buffer(pos+nset2) = scale(l)
888 if(grid1%tile(p) == tile1)
then
891 nsend1(p) = nsend1(p) + 1
892 send_buffer(pos+1) = real(i1(l), r8_kind)
893 send_buffer(pos+2) = real(j1(l), r8_kind)
894 send_buffer(pos+3) = real(i2(l), r8_kind)
895 send_buffer(pos+4) = real(j2(l), r8_kind)
896 send_buffer(pos+5) = area(l)
897 if(use_higher_order)
then
898 send_buffer(pos+6) = di(l)
899 send_buffer(pos+7) = dj(l)
907 call mpp_clock_end(id_load_xgrid2)
910 call mpp_clock_begin(id_load_xgrid3)
912 if (do_alltoall)
then
914 ibuf1(2*p+1) = nsend1(p)
915 ibuf1(2*p+2) = nsend2(p)
920 p = mod(mypos+npes-n, npes)
921 if(.not. subset_rootpe(p)) cycle
922 call mpp_recv( ibuf2(2*p+1), glen=2, from_pe=pelist(p), block=.false., tag=comm_tag_1)
925 if(nxgrid_local_orig>0)
then
927 p = mod(mypos+n, npes)
928 ibuf1(2*p+1) = nsend1(p)
929 ibuf1(2*p+2) = nsend2(p)
930 call mpp_send( ibuf1(2*p+1), plen=2, to_pe=pelist(p), tag=comm_tag_1)
936 nrecv1(p) = ibuf2(2*p+1)
937 nrecv2(p) = ibuf2(2*p+2)
941 call mpp_clock_end(id_load_xgrid3)
942 call mpp_clock_begin(id_load_xgrid4)
945 recv_buffer_pos(p) = pos
946 pos = pos + nrecv1(p) * nset1 + nrecv2(p) * nset2
950 nxgrid1 = sum(nrecv1)
951 nxgrid2 = sum(nrecv2)
952 if(nxgrid1>0 .OR. nxgrid2>0)
allocate(recv_buffer(nxgrid1*nset1+nxgrid2*nset2))
954 if (do_alltoallv)
then
956 send_cnt(:) = nset1 * nsend1(:) + nset2 * nsend2(:)
957 recv_cnt(:) = nset1 * nrecv1(:) + nset2 * nrecv2(:)
959 call mpp_alltoall(send_buffer, send_cnt, send_buffer_pos, &
960 recv_buffer, recv_cnt, recv_buffer_pos)
963 p = mod(mypos+npes-n, npes)
964 nrecv = nrecv1(p)*nset1+nrecv2(p)*nset2
966 pos = recv_buffer_pos(p)
967 call mpp_recv(recv_buffer(pos+1), glen=nrecv, from_pe=pelist(p), &
968 block=.false., tag=comm_tag_2)
972 p = mod(mypos+n, npes)
973 nsend = nsend1(p)*nset1 + nsend2(p)*nset2
975 pos = send_buffer_pos(p)
976 call mpp_send(send_buffer(pos+1), plen=nsend, to_pe=pelist(p), &
981 call mpp_clock_end(id_load_xgrid4)
983 if( nxgrid_local>0)
then
984 deallocate(i1,j1,i2,j2,area)
987 allocate(i1(nxgrid2), j1(nxgrid2))
988 allocate(i2(nxgrid2), j2(nxgrid2))
989 allocate(area(nxgrid2))
990 allocate(i1_side1(nxgrid1), j1_side1(nxgrid1))
991 allocate(i2_side1(nxgrid1), j2_side1(nxgrid1))
992 allocate(area_side1(nxgrid1))
993 if(use_higher_order)
then
994 if(nxgrid_local>0)
deallocate(di,dj)
995 allocate(di(nxgrid2), dj(nxgrid2))
996 allocate(di_side1(nxgrid1), dj_side1(nxgrid1))
999 if(nxgrid_local>0)
deallocate(scale)
1000 allocate(scale(nxgrid2))
1007 i1(l2) = int(recv_buffer(pos+1))
1008 j1(l2) = int(recv_buffer(pos+2))
1009 i2(l2) = int(recv_buffer(pos+3))
1010 j2(l2) = int(recv_buffer(pos+4))
1011 area(l2) = recv_buffer(pos+5)
1012 if(use_higher_order)
then
1013 di(l2) = recv_buffer(pos+6)
1014 dj(l2) = recv_buffer(pos+7)
1016 if(scale_exist)scale(l2) = recv_buffer(pos+nset2)
1021 i1_side1(l1) = int(recv_buffer(pos+1))
1022 j1_side1(l1) = int(recv_buffer(pos+2))
1023 i2_side1(l1) = int(recv_buffer(pos+3))
1024 j2_side1(l1) = int(recv_buffer(pos+4))
1025 area_side1(l1) = recv_buffer(pos+5)
1026 if(use_higher_order)
then
1027 di_side1(l1) = recv_buffer(pos+6)
1028 dj_side1(l1) = recv_buffer(pos+7)
1034 if(
allocated(send_buffer))
deallocate(send_buffer)
1035 if(
allocated(recv_buffer))
deallocate(recv_buffer)
1040 i1_side1 => i1; j1_side1 => j1
1041 i2_side1 => i2; j2_side1 => j2
1043 if(use_higher_order)
then
1049 call mpp_clock_begin(id_load_xgrid5)
1052 size_prev = grid%size
1054 if(grid%tile_me == tile2)
then
1056 if (
in_box_me(i2(l), j2(l), grid) )
then
1057 grid%size = grid%size + 1
1059 if( grid1_id .NE.
"ATM" .OR. tile1 .NE. tile_parent .OR. &
1060 .NOT.
in_box(i1(l), j1(l), is_parent, ie_parent, js_parent, je_parent) )
then
1062 lll = grid%l_index((j2(l)-1)*grid%im+i2(l))
1063 grid%area(lll,1) = grid%area(lll,1)+area(l)
1065 grid%area(i2(l),j2(l)) = grid%area(i2(l),j2(l))+area(l)
1069 if(grid1%tile(p) == tile1)
then
1071 xmap%your1my2(p) = .true.
1079 if(grid%size > size_prev)
then
1080 if(size_prev > 0)
then
1081 allocate(x_local(size_prev))
1083 if(
ASSOCIATED(grid%x))
deallocate(grid%x)
1084 allocate( grid%x( grid%size ) )
1085 grid%x(1:size_prev) = x_local
1088 if(
ASSOCIATED(grid%x))
deallocate(grid%x)
1089 allocate( grid%x( grid%size ) )
1090 grid%x%di = 0.0_r8_kind; grid%x%dj = 0.0_r8_kind
1095 if( grid%tile_me == tile2 )
then
1100 grid%x(ll)%i1 = i1(l); grid%x(ll)%i2 = i2(l)
1101 grid%x(ll)%j1 = j1(l); grid%x(ll)%j2 = j2(l)
1103 grid%x(ll)%l2 = grid%l_index((j2(l)-1)*grid%im + i2(l))
1108 grid%x(ll)%tile = tile1
1109 grid%x(ll)%area = area(l)
1110 if(scale_exist)
then
1111 grid%x(ll)%scale = scale(l)
1113 grid%x(ll)%scale = 1.0_r8_kind
1115 if(use_higher_order)
then
1116 grid%x(ll)%di = di(l)
1117 grid%x(ll)%dj = dj(l)
1122 if(grid1%tile(p) == tile1)
then
1124 grid%x(ll)%pe = p + xmap%root_pe
1133 if(grid%id == xmap%grids(
size(xmap%grids(:)))%id)
then
1140 if(grid1%tile_me == tile1)
then
1141 if(
associated(iarray))
then
1142 nxgrid1_old =
size(iarray(:))
1147 allocate(y2m1_pe(nxgrid1))
1148 if(.not. last_grid )
allocate(pos_s(0:xmap%npes-1))
1150 if(nxgrid1_old > 0)
then
1152 y2m1_size(p) = xmap%your2my1_size(p)
1159 if (
in_box_me(i1_side1(l), j1_side1(l), grid1) )
then
1160 if(grid1%is_ug)
then
1161 lll = grid1%l_index((j1_side1(l)-1)*grid1%im+i1_side1(l))
1162 grid1%area(lll,1) = grid1%area(lll,1) + area_side1(l)
1164 grid1%area(i1_side1(l),j1_side1(l)) = grid1%area(i1_side1(l),j1_side1(l))+area_side1(l)
1167 if (grid%tile(p) == tile2)
then
1168 if (
in_box_nbr(i2_side1(l), j2_side1(l), grid, p))
then
1169 xmap%your2my1(p) = .true.
1171 y2m1_size(p) = y2m1_size(p) + 1
1175 size_repro = size_repro + 1
1180 pos_x(p) = pos_x(p-1) + y2m1_size(p-1)
1183 if(.not. last_grid) pos_s(:) = pos_x(:)
1185 if(nxgrid1_old > 0)
then
1186 y2m1_size(:) = xmap%your2my1_size(:)
1189 allocate(iarray(nxgrid1+nxgrid1_old), jarray(nxgrid1+nxgrid1_old))
1192 do n = 1, xmap%your2my1_size(p)
1193 iarray(pos_x(p)+n) = iarray2(pos_s(p)+n)
1194 jarray(pos_x(p)+n) = jarray2(pos_s(p)+n)
1197 deallocate(iarray2, jarray2)
1199 allocate(iarray(nxgrid1), jarray(nxgrid1))
1209 if(y2m1_size(p) > 0)
then
1210 pos = pos_x(p)+y2m1_size(p)
1211 if( i1_side1(l) == iarray(pos) .AND. j1_side1(l) == jarray(pos) )
then
1215 do n = 1, y2m1_size(p)
1217 if(i1_side1(l) == iarray(pos) .AND. j1_side1(l) == jarray(pos))
then
1224 if( (.NOT. found) .OR. monotonic_exchange )
then
1225 y2m1_size(p) = y2m1_size(p)+1
1226 pos = pos_x(p)+y2m1_size(p)
1227 iarray(pos) = i1_side1(l)
1228 jarray(pos) = j1_side1(l)
1231 xmap%your2my1_size(:) = y2m1_size(:)
1234 deallocate(iarray, jarray)
1235 if(
allocated(pos_s))
deallocate(pos_s)
1239 if (grid1%tile_me == tile1 .and. size_repro > 0)
then
1240 ll_repro = grid%size_repro
1241 grid%size_repro = ll_repro + size_repro
1242 if(ll_repro > 0)
then
1243 allocate(x_local(ll_repro))
1244 x_local = grid%x_repro
1245 if(
ASSOCIATED(grid%x_repro))
deallocate(grid%x_repro)
1246 allocate( grid%x_repro(grid%size_repro ) )
1247 grid%x_repro(1:ll_repro) = x_local
1250 if(
ASSOCIATED(grid%x_repro))
deallocate(grid%x_repro)
1251 allocate( grid%x_repro( grid%size_repro ) )
1252 grid%x_repro%di = 0.0_r8_kind; grid%x_repro%dj = 0.0_r8_kind
1255 if (
in_box_me(i1_side1(l),j1_side1(l), grid1) )
then
1256 ll_repro = ll_repro + 1
1257 grid%x_repro(ll_repro)%i1 = i1_side1(l); grid%x_repro(ll_repro)%i2 = i2_side1(l)
1258 grid%x_repro(ll_repro)%j1 = j1_side1(l); grid%x_repro(ll_repro)%j2 = j2_side1(l)
1259 if(grid1%is_ug)
then
1260 grid%x_repro(ll_repro)%l1 = grid1%l_index((j1_side1(l)-1)*grid1%im+i1_side1(l))
1265 grid%x_repro(ll_repro)%tile = tile1
1266 grid%x_repro(ll_repro)%area = area_side1(l)
1267 if(use_higher_order)
then
1268 grid%x_repro(ll_repro)%di = di_side1(l)
1269 grid%x_repro(ll_repro)%dj = dj_side1(l)
1273 if(grid%tile(p) == tile2)
then
1274 if (
in_box_nbr(i2_side1(l), j2_side1(l), grid, p))
then
1275 grid%x_repro(ll_repro)%pe = p + xmap%root_pe
1283 deallocate(i1, j1, i2, j2, area)
1284 if(use_higher_order)
deallocate(di, dj)
1285 if(scale_exist)
deallocate(scale)
1286 if(is_distribute)
then
1287 deallocate(i1_side1, j1_side1, i2_side1, j2_side1, area_side1)
1288 if(use_higher_order)
deallocate(di_side1, dj_side1)
1291 i1=>null(); j1=>null(); i2=>null(); j2=>null()
1292 call mpp_clock_end(id_load_xgrid5)
1296 end subroutine load_xgrid
1305 character(len=3),
intent(in) :: grid_id
1306 character(len=*),
intent(in) :: grid_file
1308 real(r8_kind),
dimension(grid%im) :: lonb
1309 real(r8_kind),
dimension(grid%jm) :: latb
1310 real(r8_kind) :: d2r
1311 integer :: is, ie, js, je
1312 type(fmsnetcdfdomainfile_t) :: fileobj
1314 d2r = pi / 180.0_r8_kind
1316 if(.not.
open_file(fileobj, grid_file,
'read', grid%domain) )
then
1317 call error_mesg(
'xgrid_mod(get_grid_version1)',
'Error in opening file '//trim(grid_file), fatal)
1321 if (
associated(grid%lon))
deallocate(grid%lon)
1322 if (
associated(grid%lat))
deallocate(grid%lat)
1323 allocate(grid%lon(grid%im), grid%lat(grid%jm))
1324 if(grid_id ==
'ATM')
then
1336 else if(grid_id ==
'LND')
then
1339 if(.not.
allocated(area_lnd_model))
then
1340 allocate(area_lnd_model(is:ie, js:je))
1343 if(.not.
allocated(area_lnd_sphere))
then
1344 allocate(area_lnd_sphere(is:ie, js:je))
1347 else if(grid_id ==
'OCN' )
then
1348 if(.not.
allocated(area_ocn_sphere))
then
1349 allocate(area_ocn_sphere(is:ie, js:je))
1354 if(grid_id ==
'LND' .or. grid_id ==
'ATM')
then
1355 grid%lon = lonb * d2r
1356 grid%lat = latb * d2r
1358 grid%is_latlon = .true.
1373 character(len=3),
intent(in) :: grid_id
1374 character(len=*),
intent(in) :: grid_file
1376 real(r8_kind),
allocatable :: tmpx(:,:), tmpy(:,:)
1377 real(r8_kind) :: d2r
1378 integer :: is, ie, js, je, nlon, nlat, i, j
1379 integer :: start(4), nread(4), isc2, iec2, jsc2, jec2
1380 type(fmsnetcdffile_t) :: fileobj
1382 if(.not.
open_file(fileobj, grid_file,
'read') )
then
1383 call error_mesg(
'xgrid_mod(get_grid_version2)',
'Error in opening file '//trim(grid_file), fatal)
1386 d2r = pi / 180.0_r8_kind
1390 call get_dimension_size(fileobj,
"nx", nlon)
1391 call get_dimension_size(fileobj,
"ny", nlat)
1392 if( mod(nlon,2) .NE. 0)
call error_mesg(
'xgrid_mod', &
1393 'flux_exchange_mod: atmos supergrid longitude size can not be divided by 2', fatal)
1394 if( mod(nlat,2) .NE. 0)
call error_mesg(
'xgrid_mod', &
1395 'flux_exchange_mod: atmos supergrid latitude size can not be divided by 2', fatal)
1398 if(nlon .NE. grid%im .OR. nlat .NE. grid%jm)
call error_mesg(
'xgrid_mod', &
1399 'grid size in tile_file does not match the global grid size', fatal)
1401 if( grid_id ==
'LND' .or. grid_id ==
'ATM' .or. grid_id ==
'WAV' )
then
1402 isc2 = 2*grid%is_me-1; iec2 = 2*grid%ie_me+1
1403 jsc2 = 2*grid%js_me-1; jec2 = 2*grid%je_me+1
1404 allocate(tmpx(isc2:iec2, jsc2:jec2) )
1405 allocate(tmpy(isc2:iec2, jsc2:jec2) )
1406 start = 1; nread = 1
1407 start(1) = isc2; nread(1) = iec2 - isc2 + 1
1408 start(2) = jsc2; nread(2) = jec2 - jsc2 + 1
1409 call read_data(fileobj,
'x', tmpx, corner=start, edge_lengths=nread)
1410 call read_data(fileobj,
'y', tmpy, corner=start, edge_lengths=nread)
1411 if(is_lat_lon(tmpx, tmpy) )
then
1412 deallocate(tmpx, tmpy)
1413 start = 1; nread = 1
1414 start(2) = 2; nread(1) = nlon*2+1
1415 allocate(tmpx(nlon*2+1, 1), tmpy(1, nlat*2+1))
1416 call read_data(fileobj,
"x", tmpx, corner=start, edge_lengths=nread)
1417 if (
associated(grid%lon))
deallocate(grid%lon)
1418 if (
associated(grid%lat))
deallocate(grid%lat)
1419 allocate(grid%lon(grid%im), grid%lat(grid%jm))
1421 grid%lon(i) = tmpx(2*i,1) * d2r
1423 start = 1; nread = 1
1424 start(1) = 2; nread(2) = nlat*2+1
1425 call read_data(fileobj,
"y", tmpy, corner=start, edge_lengths=nread)
1427 grid%lat(j) = tmpy(1, 2*j) * d2r
1429 grid%is_latlon = .true.
1431 if (
associated(grid%geolon))
deallocate(grid%geolon)
1432 if (
associated(grid%geolat))
deallocate(grid%geolat)
1433 allocate(grid%geolon(grid%isd_me:grid%ied_me, grid%jsd_me:grid%jed_me))
1434 allocate(grid%geolat(grid%isd_me:grid%ied_me, grid%jsd_me:grid%jed_me))
1435 grid%geolon = 1.0e10_r8_kind
1436 grid%geolat = 1.0e10_r8_kind
1438 do j = grid%js_me,grid%je_me
1439 do i = grid%is_me,grid%ie_me
1440 grid%geolon(i, j) = tmpx(i*2,j*2)*d2r
1441 grid%geolat(i, j) = tmpy(i*2,j*2)*d2r
1446 grid%is_latlon = .false.
1448 deallocate(tmpx, tmpy)
1460 type(fmsnetcdfdomainfile_t),
intent(in) :: fileobj
1461 character(len=*),
intent(in) :: name
1462 real(r8_kind),
intent(out) :: get_area_data(:,:)
1464 if(variable_exists(fileobj, name))
then
1465 call read_data(fileobj, name, get_area_data)
1467 call error_mesg(
'xgrid_mod',
'no field named '//trim(name)//
' in grid file '//trim(fileobj%path)// &
1468 ' Will set data to negative values...', note)
1470 get_area_data = -1.0_r8_kind
1484 type(
domain2d),
intent(in) :: domain
1485 character(len=*),
intent(in) :: grid_file
1486 integer :: is, ie, js, je
1487 type(fmsnetcdffile_t) :: fileobj
1489 if(
allocated(area_ocn_model))
return
1494 allocate(area_ocn_model(is:ie, js:je))
1495 if(ie < is .or. je < js )
return
1497 if(.not.
open_file(fileobj, grid_file,
'read') )
then
1498 call error_mesg(
'xgrid_mod(get_ocean_model_area_elements)',
'Error in opening file '//trim(grid_file), fatal)
1501 if(variable_exists(fileobj,
'AREA_OCN_MODEL') )
then
1502 call read_data(fileobj,
'AREA_OCN_MODEL', area_ocn_model)
1504 deallocate(area_ocn_model)
1515 subroutine setup_xmap(xmap, grid_ids, grid_domains, grid_file, atm_grid, lnd_ug_domain)
1517 character(len=3),
dimension(:),
intent(in ) :: grid_ids
1518 type(
domain2d),
dimension(:),
intent(in ) :: grid_domains
1519 character(len=*),
intent(in ) :: grid_file
1521 type(
domainug),
optional,
intent(in ) :: lnd_ug_domain
1524 integer :: nxgrid_file, i1, i2, i3, tile1, tile2, j
1525 integer :: nxc, nyc, out_unit
1526 type(
grid_type),
pointer :: grid => null()
1527 type(
grid_type),
pointer,
save :: grid1 => null()
1528 real(r8_kind),
dimension(3) :: xxx
1529 real(r8_kind),
dimension(:,:),
allocatable :: check_data
1530 real(r8_kind),
dimension(:,:,:),
allocatable :: check_data_3d
1531 real(r8_kind),
allocatable :: tmp_2d(:,:), tmp_3d(:,:,:)
1532 character(len=FMS_FILE_LEN) :: xgrid_file, xgrid_name
1533 character(len=FMS_FILE_LEN) :: tile_file, mosaic_file
1534 character(len=256) :: mosaic1, mosaic2, contact, xgrid_dimname
1535 character(len=256) :: tile1_name, tile2_name
1536 character(len=256),
allocatable :: tile1_list(:), tile2_list(:)
1537 character(len=FMS_FILE_LEN),
allocatable :: xgrid_filelist(:)
1538 integer :: npes, npes2
1539 integer,
allocatable :: pelist(:)
1541 logical :: use_higher_order = .false.
1542 integer :: lnd_ug_id, l
1543 integer,
allocatable :: grid_index(:)
1544 type(fmsnetcdffile_t) :: gridfileobj, mosaicfileobj, fileobj
1545 type(
grid_type),
allocatable,
target :: grids_tmp(:)
1548 call mpp_clock_begin(id_setup_xmap)
1550 if(
interp_method .ne.
'first_order') use_higher_order = .true.
1555 xmap%root_pe = mpp_root_pe()
1557 if (
associated(xmap%grids))
deallocate(xmap%grids)
1558 allocate( xmap%grids(1:
size(grid_ids(:))) )
1560 if (
associated(xmap%your1my2))
deallocate(xmap%your1my2)
1561 if (
associated(xmap%your2my1))
deallocate(xmap%your2my1)
1562 if (
associated(xmap%your2my1_size))
deallocate(xmap%your2my1_size)
1563 allocate ( xmap%your1my2(0:xmap%npes-1), xmap%your2my1(0:xmap%npes-1) )
1564 allocate ( xmap%your2my1_size(0:xmap%npes-1) )
1566 xmap%your1my2 = .false.; xmap%your2my1 = .false.;
1567 xmap%your2my1_size = 0
1569 if(.not.
open_file(gridfileobj,trim(grid_file),
"read"))
then
1570 call error_mesg(
'xgrid_mod',
'Error when opening file'//trim(grid_file), fatal)
1574 if(variable_exists(gridfileobj,
"AREA_ATMxOCN" ) )
then
1577 else if(variable_exists(gridfileobj,
"ocn_mosaic_file" ) )
then
1580 call error_mesg(
'xgrid_mod',
'both AREA_ATMxOCN and ocn_mosaic_file does not exist in '//trim(grid_file), fatal)
1585 call error_mesg(
'xgrid_mod',
'reading exchange grid information from grid spec file', note)
1587 call error_mesg(
'xgrid_mod',
'reading exchange grid information from mosaic grid file', note)
1592 if(
present(lnd_ug_domain))
then
1593 do g=1,
size(grid_ids(:))
1594 if(grid_ids(g) ==
'LND') lnd_ug_id = g
1598 call mpp_clock_begin(id_load_xgrid)
1602 grids_tmp = xmap%grids
1604 grid1 => xmap%grids(1)
1606 do g=1,
size(grid_ids(:))
1608 grid => grids_tmp(g)
1610 grid%id = grid_ids(g)
1611 grid%domain = grid_domains(g)
1613 if (
associated(grid%is))
deallocate(grid%is)
1614 if (
associated(grid%ie))
deallocate(grid%ie)
1615 if (
associated(grid%js))
deallocate(grid%js)
1616 if (
associated(grid%je))
deallocate(grid%je)
1617 if (
associated(grid%tile))
deallocate(grid%tile)
1618 allocate ( grid%is(0:xmap%npes-1), grid%ie(0:xmap%npes-1) )
1619 allocate ( grid%js(0:xmap%npes-1), grid%je(0:xmap%npes-1) )
1620 allocate ( grid%tile(0:xmap%npes-1) )
1630 select case(xmap%version)
1634 call read_data(gridfileobj, lowercase(grid_ids(g))//
'_mosaic_file', mosaic_file)
1635 if(.not.
open_file(mosaicfileobj,
'INPUT/'//trim(mosaic_file),
"read"))
then
1636 call error_mesg(
'xgrid_mod',
'Error when opening solo mosaic file INPUT/'//trim(mosaic_file), fatal)
1638 call get_dimension_size(mosaicfileobj,
'ntiles', grid%ntile)
1641 if( g == 1 .AND. grid_ids(1) ==
'ATM' )
then
1642 if( .NOT. grid%on_this_pe )
call error_mesg(
'xgrid_mod',
'ATM domain is not defined on some processor' ,fatal)
1645 if( xmap%npes > grid%npes .AND. g == 1 .AND. grid_ids(1) ==
'ATM' )
then
1647 else if(xmap%npes > grid%npes)
then
1653 allocate(grid%pelist(0:npes-1))
1657 call mpp_get_data_domain(grid%domain, grid%isd_me, grid%ied_me, grid%jsd_me, grid%jed_me, &
1658 xsize=grid%nxd_me, ysize=grid%nyd_me)
1661 if( grid%root_pe == xmap%root_pe )
then
1663 ybegin=grid%js(0:npes-1), yend=grid%je(0:npes-1) )
1665 if( xmap%npes > npes .AND. g == 1 .AND. grid_ids(1) ==
'ATM' )
then
1667 ybegin=grid%js(npes:xmap%npes-1), yend=grid%je(npes:xmap%npes-1) )
1671 npes2 = xmap%npes-npes
1673 ybegin=grid%js(0:npes2-1), yend=grid%je(0:npes2-1) )
1674 call mpp_get_compute_domains(grid%domain, xbegin=grid%is(npes2:xmap%npes-1), xend=grid%ie(npes2:xmap%npes-1), &
1675 ybegin=grid%js(npes2:xmap%npes-1), yend=grid%je(npes2:xmap%npes-1) )
1679 if( xmap%npes > grid%npes .AND. g == 1 .AND. grid_ids(1) ==
'ATM' )
then
1683 if( g == 1 .AND. grid_ids(1) ==
'ATM' ) npes = xmap%npes
1685 if(grid%tile(p) > grid%ntile .or. grid%tile(p) < 1)
call error_mesg(
'xgrid_mod', &
1686 'tile id should between 1 and ntile', fatal)
1694 grid%is_me => grid%is(xmap%me-xmap%root_pe); grid%ie_me => grid%ie(xmap%me-xmap%root_pe)
1695 grid%js_me => grid%js(xmap%me-xmap%root_pe); grid%je_me => grid%je(xmap%me-xmap%root_pe)
1696 grid%nxc_me = grid%ie_me - grid%is_me + 1
1697 grid%nyc_me = grid%je_me - grid%js_me + 1
1698 grid%tile_me => grid%tile(xmap%me-xmap%root_pe)
1701 grid%is_ug = .false.
1703 if( g == lnd_ug_id )
then
1705 'does not support unstructured grid for VERSION1 grid' ,fatal)
1707 grid%ug_domain = lnd_ug_domain
1708 if (
associated(grid%ls))
deallocate(grid%ls)
1709 if (
associated(grid%le))
deallocate(grid%le)
1710 if (
associated(grid%gs))
deallocate(grid%gs)
1711 if (
associated(grid%ge))
deallocate(grid%ge)
1712 allocate ( grid%ls(0:xmap%npes-1), grid%le(0:xmap%npes-1) )
1713 allocate ( grid%gs(0:xmap%npes-1), grid%ge(0:xmap%npes-1) )
1718 if(xmap%npes > grid%npes)
then
1721 call mpp_get_ug_compute_domains(grid%ug_domain, begin=grid%ls(0:npes-1),
end=grid%le(0:npes-1) )
1722 call mpp_get_ug_domains_index(grid%ug_domain, grid%gs(0:npes-1), grid%ge(0:npes-1) )
1723 call mpp_get_ug_domain_tile_list(grid%ug_domain, grid%tile(0:npes-1))
1724 grid%ls_me => grid%ls(xmap%me-xmap%root_pe); grid%le_me => grid%le(xmap%me-xmap%root_pe)
1725 grid%gs_me => grid%gs(xmap%me-xmap%root_pe); grid%ge_me => grid%ge(xmap%me-xmap%root_pe)
1726 grid%tile_me => grid%tile(xmap%me-xmap%root_pe)
1727 grid%nxl_me = grid%le_me - grid%ls_me + 1
1728 if (
associated(grid%l_index))
deallocate(grid%l_index)
1729 allocate(grid%l_index(grid%gs_me:grid%ge_me))
1730 allocate(grid_index(grid%ls_me:grid%le_me))
1731 call mpp_get_ug_domain_grid_index(grid%ug_domain, grid_index)
1734 do l = grid%ls_me,grid%le_me
1735 grid%l_index(grid_index(l)) = l
1738 if( grid%on_this_pe )
then
1739 if (
associated(grid%area))
deallocate(grid%area)
1740 if (
associated(grid%area_inv))
deallocate(grid%area_inv)
1741 allocate( grid%area (grid%ls_me:grid%le_me,1) )
1742 allocate( grid%area_inv(grid%ls_me:grid%le_me,1) )
1743 grid%area = 0.0_r8_kind
1747 else if( grid%on_this_pe )
then
1748 if (
associated(grid%area))
deallocate(grid%area)
1749 if (
associated(grid%area_inv))
deallocate(grid%area_inv)
1750 allocate( grid%area (grid%is_me:grid%ie_me, grid%js_me:grid%je_me) )
1751 allocate( grid%area_inv(grid%is_me:grid%ie_me, grid%js_me:grid%je_me) )
1752 grid%area = 0.0_r8_kind
1758 if(.not. grid%is_ug)
then
1759 select case(xmap%version)
1761 if( grid%npes .NE. xmap%npes )
then
1762 call error_mesg(
'xgrid_mod', .NE.
' grid%npes xmap%npes ', fatal)
1766 allocate(pelist(0:xmap%npes-1))
1767 call mpp_get_current_pelist(pelist)
1768 if( grid%on_this_pe )
then
1776 if( g == 1 .AND. grid_ids(1) ==
'ATM' )
then
1778 ie_nest, js_nest, je_nest, is_parent, ie_parent, js_parent, je_parent)
1783 if( use_higher_order .AND. grid%id ==
'ATM')
then
1784 if( nnest > 0 )
call error_mesg(
'xgrid_mod',
'second_order is not supported for nested coupler', fatal)
1785 if( grid%is_latlon )
then
1786 call mpp_modify_domain(grid%domain, grid%domain_with_halo, whalo=1, ehalo=1, shalo=1, nhalo=1)
1787 call mpp_get_data_domain(grid%domain_with_halo, grid%isd_me, grid%ied_me, grid%jsd_me, grid%jed_me, &
1788 xsize=grid%nxd_me, ysize=grid%nyd_me)
1790 if(.NOT.
present(atm_grid))
call error_mesg(
'xgrid_mod', &
1791 'when first grid is "ATM", atm_grid should be present', fatal)
1792 if(grid%is_me-grid%isd_me .NE. 1 .or. grid%ied_me-grid%ie_me .NE. 1 .or. &
1793 grid%js_me-grid%jsd_me .NE. 1 .or. grid%jed_me-grid%je_me .NE. 1 ) &
1794 &
call error_mesg(
'xgrid_mod',
'for non-latlon grid (cubic grid), '//&
1795 &
'the halo size should be 1 in all four direction', fatal)
1796 if(.NOT.(
ASSOCIATED(atm_grid%dx) .AND.
ASSOCIATED(atm_grid%dy) .AND.
ASSOCIATED(atm_grid%edge_w) .AND. &
1797 ASSOCIATED(atm_grid%edge_e) .AND.
ASSOCIATED(atm_grid%edge_s) .AND.
ASSOCIATED(atm_grid%edge_n).AND.&
1798 ASSOCIATED(atm_grid%en1) .AND.
ASSOCIATED(atm_grid%en2) .AND.
ASSOCIATED(atm_grid%vlon) .AND. &
1799 ASSOCIATED(atm_grid%vlat) ) )
call error_mesg(
'xgrid_mod', &
1800 'for non-latlon grid (cubic grid), all the fields in atm_grid data type should be allocated', fatal)
1801 nxc = grid%ie_me - grid%is_me + 1
1802 nyc = grid%je_me - grid%js_me + 1
1803 if(
size(atm_grid%dx,1) .NE. nxc .OR.
size(atm_grid%dx,2) .NE. nyc+1) &
1804 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%dx', fatal)
1805 if(
size(atm_grid%dy,1) .NE. nxc+1 .OR.
size(atm_grid%dy,2) .NE. nyc) &
1806 call error_mesg(
'xgrid_mod',
'incorrect dimension sizeof atm_grid%dy', fatal)
1807 if(
size(atm_grid%area,1) .NE. nxc .OR.
size(atm_grid%area,2) .NE. nyc) &
1808 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%area', fatal)
1809 if(
size(atm_grid%edge_w(:)) .NE. nyc+1 .OR.
size(atm_grid%edge_e(:)) .NE. nyc+1) &
1810 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%edge_w/edge_e', fatal)
1811 if(
size(atm_grid%edge_s(:)) .NE. nxc+1 .OR.
size(atm_grid%edge_n(:)) .NE. nxc+1) &
1812 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%edge_s/edge_n', fatal)
1813 if(
size(atm_grid%en1,1) .NE. 3 .OR.
size(atm_grid%en1,2) .NE. nxc .OR.
size(atm_grid%en1,3) .NE. nyc+1) &
1814 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%en1', fatal)
1815 if(
size(atm_grid%en2,1) .NE. 3 .OR.
size(atm_grid%en2,2) .NE. nxc+1 .OR.
size(atm_grid%en2,3) .NE. nyc) &
1816 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%en2', fatal)
1817 if(
size(atm_grid%vlon,1) .NE. 3 .OR.
size(atm_grid%vlon,2) .NE. nxc .OR.
size(atm_grid%vlon,3) .NE. nyc)&
1818 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%vlon', fatal)
1819 if(
size(atm_grid%vlat,1) .NE. 3 .OR.
size(atm_grid%vlat,2) .NE. nxc .OR.
size(atm_grid%vlat,3) .NE. nyc)&
1820 call error_mesg(
'xgrid_mod',
'incorrect dimension size of atm_grid%vlat', fatal)
1821 if (
associated(grid%box%dx))
deallocate(grid%box%dx)
1822 if (
associated(grid%box%dy))
deallocate(grid%box%dy)
1823 if (
associated(grid%box%area))
deallocate(grid%box%area)
1824 if (
associated(grid%box%edge_w))
deallocate(grid%box%edge_w)
1825 if (
associated(grid%box%edge_e))
deallocate(grid%box%edge_e)
1826 if (
associated(grid%box%edge_s))
deallocate(grid%box%edge_s)
1827 if (
associated(grid%box%edge_n))
deallocate(grid%box%edge_n)
1828 if (
associated(grid%box%en1))
deallocate(grid%box%en1)
1829 if (
associated(grid%box%en2))
deallocate(grid%box%en2)
1830 if (
associated(grid%box%vlon))
deallocate(grid%box%vlon)
1831 if (
associated(grid%box%vlat))
deallocate(grid%box%vlat)
1832 allocate(grid%box%dx (grid%is_me:grid%ie_me, grid%js_me:grid%je_me+1 ))
1833 allocate(grid%box%dy (grid%is_me:grid%ie_me+1, grid%js_me:grid%je_me ))
1834 allocate(grid%box%area (grid%is_me:grid%ie_me, grid%js_me:grid%je_me ))
1835 allocate(grid%box%edge_w(grid%js_me:grid%je_me+1))
1836 allocate(grid%box%edge_e(grid%js_me:grid%je_me+1))
1837 allocate(grid%box%edge_s(grid%is_me:grid%ie_me+1))
1838 allocate(grid%box%edge_n(grid%is_me:grid%ie_me+1))
1839 allocate(grid%box%en1 (3, grid%is_me:grid%ie_me, grid%js_me:grid%je_me+1 ))
1840 allocate(grid%box%en2 (3, grid%is_me:grid%ie_me+1, grid%js_me:grid%je_me ))
1841 allocate(grid%box%vlon (3, grid%is_me:grid%ie_me, grid%js_me:grid%je_me ))
1842 allocate(grid%box%vlat (3, grid%is_me:grid%ie_me, grid%js_me:grid%je_me ))
1843 grid%box%dx = atm_grid%dx
1844 grid%box%dy = atm_grid%dy
1845 grid%box%area = atm_grid%area
1846 grid%box%edge_w = atm_grid%edge_w
1847 grid%box%edge_e = atm_grid%edge_e
1848 grid%box%edge_s = atm_grid%edge_s
1849 grid%box%edge_n = atm_grid%edge_n
1850 grid%box%en1 = atm_grid%en1
1851 grid%box%en2 = atm_grid%en2
1852 grid%box%vlon = atm_grid%vlon
1853 grid%box%vlat = atm_grid%vlat
1859 if(grid%on_this_pe)
then
1860 if (
associated(grid%frac_area))
deallocate(grid%frac_area)
1862 allocate( grid%frac_area(grid%ls_me:grid%le_me, 1, grid%km) )
1864 allocate( grid%frac_area(grid%is_me:grid%ie_me, grid%js_me:grid%je_me, grid%km) )
1866 grid%frac_area = 1.0_r8_kind
1870 xmap%grids(g) = grid
1873 select case(xmap%version)
1875 call load_xgrid (xmap, grid, grid_file, grid_ids(1), grid_ids(g), 1, 1, use_higher_order)
1877 select case(grid_ids(1))
1885 call error_mesg(
'xgrid_mod',
'grid_ids(1) should be ATM, LND or WAV', fatal)
1887 select case(grid_ids(g))
1889 xgrid_dimname =
'nfile_'//trim(xgrid_name)//
'Xl'
1890 xgrid_name = trim(xgrid_name)//
'Xl_file'
1892 xgrid_dimname =
'nfile_'//trim(xgrid_name)//
'Xo'
1893 xgrid_name = trim(xgrid_name)//
'Xo_file'
1895 xgrid_dimname =
'nfile_'//trim(xgrid_name)//
'Xw'
1896 xgrid_name = trim(xgrid_name)//
'Xw_file'
1898 call error_mesg(
'xgrid_mod',
'grid_ids(g) should be LND, OCN or WAV', fatal)
1902 call read_data(gridfileobj, lowercase(grid_ids(1))//
'_mosaic_file', mosaic1)
1903 call read_data(gridfileobj, lowercase(grid_ids(g))//
'_mosaic_file', mosaic2)
1905 mosaic1 =
'INPUT/'//trim(mosaic1)
1906 mosaic2 =
'INPUT/'//trim(mosaic2)
1908 allocate(tile1_list(grid1%ntile), tile2_list(grid%ntile) )
1909 if(.not.
open_file(fileobj,mosaic1,
"read"))
then
1910 call error_mesg(
'xgrid_mod(setup_xmap)',
'Error when opening mosaic1 file '//trim(mosaic1), fatal)
1912 call read_data(fileobj,
'gridtiles', tile1_list)
1915 if(.not.
open_file(fileobj,mosaic2,
"read"))
then
1916 call error_mesg(
'xgrid_mod(setup_xmap)',
'Error when opening mosaic2 file '//trim(mosaic2), fatal)
1918 call read_data(fileobj,
'gridtiles', tile2_list)
1921 if(variable_exists(gridfileobj, xgrid_name))
then
1922 call get_dimension_size(gridfileobj, xgrid_dimname, nxgrid_file)
1923 if(nxgrid_file>0)
then
1924 allocate(xgrid_filelist(nxgrid_file))
1925 call read_data(gridfileobj, xgrid_name, xgrid_filelist)
1928 do i = 1, nxgrid_file
1929 xgrid_file =
'INPUT/'//trim(xgrid_filelist(i))
1930 if(.not.
open_file(fileobj,xgrid_file,
"read"))
then
1931 call error_mesg(
'xgrid_mod(setup_xmap)',
'Error when opening xgrid file '// &
1932 & trim(xgrid_file), fatal)
1936 call read_data(fileobj,
"contact", contact)
1937 i1 = index(contact,
":")
1938 i2 = index(contact,
"::")
1939 i3 = index(contact,
":", back=.true. )
1940 if(i1 == 0 .OR. i2 == 0)
call error_mesg(
'xgrid_mod', &
1941 'field contact in file '//trim(xgrid_file)//
' should contains ":" and "::" ', fatal)
1943 'field contact in file '//trim(xgrid_file)//
' should contains two ":"', fatal)
1944 tile1_name = contact(i1+1:i2-1)
1945 tile2_name = contact(i3+1:len_trim(contact))
1946 tile1 = 0; tile2 = 0
1947 do j = 1, grid1%ntile
1948 if( trim(tile1_name) == trim(tile1_list(j)) )
then
1953 do j = 1, grid%ntile
1954 if( tile2_name == tile2_list(j) )
then
1960 if(tile1 == 0)
call error_mesg(
'xgrid_mod', &
1961 trim(tile1_name)//
' is not a tile of mosaic '//trim(mosaic1), fatal)
1962 if(tile2 == 0)
call error_mesg(
'xgrid_mod', &
1963 trim(tile2_name)//
' is not a tile of mosaic '//trim(mosaic2), fatal)
1965 call load_xgrid (xmap, grid, xgrid_file, grid_ids(1), grid_ids(g), tile1, tile2, &
1968 deallocate(xgrid_filelist)
1970 deallocate(tile1_list, tile2_list)
1972 if(grid%on_this_pe)
then
1973 grid%area_inv = 0.0_r8_kind;
1974 where (grid%area>0.0_r8_kind) grid%area_inv = 1.0_r8_kind/grid%area
1979 xmap%grids(g) = grid
1984 call mpp_clock_end(id_load_xgrid)
1986 grid1%area_inv = 0.0_r8_kind;
1987 where (grid1%area>0.0_r8_kind)
1988 grid1%area_inv = 1.0_r8_kind/grid1%area
1991 xmap%your1my2(xmap%me-xmap%root_pe) = .false.
1992 xmap%your2my1(xmap%me-xmap%root_pe) = .false.
1995 if (
associated(xmap%send_count_repro))
deallocate(xmap%send_count_repro)
1996 if (
associated(xmap%recv_count_repro))
deallocate(xmap%recv_count_repro)
1997 allocate( xmap%send_count_repro(0:xmap%npes-1) )
1998 allocate( xmap%recv_count_repro(0:xmap%npes-1) )
1999 xmap%send_count_repro = 0
2000 xmap%recv_count_repro = 0
2001 do g=2,
size(xmap%grids(:))
2003 if(xmap%grids(g)%size >0) &
2004 xmap%send_count_repro(p) = xmap%send_count_repro(p) &
2005 +count(xmap%grids(g)%x (:)%pe==p+xmap%root_pe)
2006 if(xmap%grids(g)%size_repro >0) &
2007 xmap%recv_count_repro(p) = xmap%recv_count_repro(p) &
2008 +count(xmap%grids(g)%x_repro(:)%pe==p+xmap%root_pe)
2011 xmap%send_count_repro_tot = sum(xmap%send_count_repro)
2012 xmap%recv_count_repro_tot = sum(xmap%recv_count_repro)
2014 xmap%send_count_repro_tot = 0
2015 xmap%recv_count_repro_tot = 0
2018 if (
associated(xmap%x1))
deallocate(xmap%x1)
2019 if (
associated(xmap%x2))
deallocate(xmap%x2)
2020 if (
associated(xmap%x1_put))
deallocate(xmap%x1_put)
2021 if (
associated(xmap%x2_get))
deallocate(xmap%x2_get)
2022 allocate( xmap%x1(1:sum(xmap%grids(2:
size(xmap%grids(:)))%size)) )
2023 allocate( xmap%x2(1:sum(xmap%grids(2:
size(xmap%grids(:)))%size)) )
2024 allocate( xmap%x1_put(1:sum(xmap%grids(2:
size(xmap%grids(:)))%size)) )
2025 allocate( xmap%x2_get(1:sum(xmap%grids(2:
size(xmap%grids(:)))%size)) )
2028 if (
associated(xmap%get1))
deallocate(xmap%get1)
2029 if (
associated(xmap%put1))
deallocate(xmap%put1)
2030 allocate(xmap%get1, xmap%put1)
2031 call mpp_clock_begin(id_set_comm)
2033 call set_comm_get1(xmap)
2035 call set_comm_put1(xmap)
2038 if (
associated(xmap%get1_repro))
deallocate(xmap%get1_repro)
2039 allocate(xmap%get1_repro)
2040 call set_comm_get1_repro(xmap)
2043 call mpp_clock_end(id_set_comm)
2045 call mpp_clock_begin(id_regen)
2047 call mpp_clock_end(id_regen)
2049 call mpp_clock_begin(id_conservation_check)
2051 if(lnd_ug_id ==0)
then
2054 allocate(tmp_2d(grid1%is_me:grid1%ie_me, grid1%js_me:grid1%je_me))
2055 tmp_2d = 1.0_r8_kind
2059 write(out_unit,* )
"Checked data is array of constant 1"
2060 write(out_unit,* )grid1%id,
'(',xmap%grids(:)%id,
')=', xxx
2062 if(lnd_ug_id == 0)
then
2063 do g=2,
size(xmap%grids(:))
2064 xxx =
conservation_check(xmap%grids(g)%frac_area*0.0_r8_kind+1.0_r8_kind, xmap%grids(g)%id, xmap )
2065 write( out_unit,* )xmap%grids(g)%id,
'(',xmap%grids(:)%id,
')=', xxx
2068 do g=2,
size(xmap%grids(:))
2069 grid => xmap%grids(g)
2070 allocate(tmp_3d(grid%is_me:grid%ie_me, grid%js_me:grid%je_me,grid%km))
2071 tmp_3d = 1.0_r8_kind
2073 write( out_unit,* )xmap%grids(g)%id,
'(',xmap%grids(:)%id,
')=', xxx
2078 if(grid1%id ==
"ATM")
then
2079 allocate(check_data(
size(grid1%area,1),
size(grid1%area,2)))
2080 call random_number(check_data)
2083 if(lnd_ug_id ==0)
then
2088 write( out_unit,* ) &
2089 "Checked data is array of random number between 0 and 1 using "//trim(
interp_method)
2090 write( out_unit,* )grid1%id,
'(',xmap%grids(:)%id,
')=', xxx
2092 deallocate(check_data)
2093 do g=2,
size(xmap%grids(:))
2094 allocate(check_data_3d(xmap%grids(g)%is_me:xmap%grids(g)%ie_me, &
2095 xmap%grids(g)%js_me:xmap%grids(g)%je_me, grid1%km))
2096 call random_number(check_data_3d)
2097 if(lnd_ug_id ==0)
then
2102 write( out_unit,* )xmap%grids(g)%id,
'(',xmap%grids(:)%id,
')=', xxx
2103 deallocate( check_data_3d)
2106 call mpp_clock_end(id_conservation_check)
2108 call mpp_clock_end(id_setup_xmap)
2117 ie_nest_out, js_nest_out, je_nest_out, is_parent_out, &
2118 ie_parent_out, js_parent_out, je_parent_out) &
2120 type(fmsnetcdffile_t),
intent(in) :: fileobj
2121 integer,
intent(out) :: tile_nest_out, tile_parent_out
2122 integer,
intent(out) :: is_nest_out, ie_nest_out
2123 integer,
intent(out) :: js_nest_out, je_nest_out
2124 integer,
intent(out) :: is_parent_out, ie_parent_out
2125 integer,
intent(out) :: js_parent_out, je_parent_out
2128 integer :: ntiles, ncontacts, n, t1, t2
2129 integer :: nx1_contact, ny1_contact
2130 integer :: nx2_contact, ny2_contact
2131 integer,
allocatable,
dimension(:) :: nx, ny
2132 integer,
allocatable,
dimension(:) :: tile1, tile2
2133 integer,
allocatable,
dimension(:) :: istart1, iend1, jstart1, jend1
2134 integer,
allocatable,
dimension(:) :: istart2, iend2, jstart2, jend2
2136 tile_nest_out = 0; tile_parent_out = 0
2137 is_nest_out = 0; ie_nest_out = 0
2138 js_nest_out = 0; je_nest_out = 0
2139 is_parent_out = 0; ie_parent_out = 0
2140 js_parent_out = 0; je_parent_out = 0
2144 ntiles = get_mosaic_ntiles(fileobj)
2145 if( ntiles == 1 )
return
2147 allocate(nx(ntiles), ny(ntiles))
2148 call get_mosaic_grid_sizes(fileobj, nx, ny)
2150 ncontacts = get_mosaic_ncontacts(fileobj)
2151 if(ncontacts == 0)
return
2152 allocate(tile1(ncontacts), tile2(ncontacts))
2153 allocate(istart1(ncontacts), iend1(ncontacts))
2154 allocate(jstart1(ncontacts), jend1(ncontacts))
2155 allocate(istart2(ncontacts), iend2(ncontacts))
2156 allocate(jstart2(ncontacts), jend2(ncontacts))
2158 call get_mosaic_contact( fileobj, tile1, tile2, istart1, iend1, jstart1, jend1, &
2159 istart2, iend2, jstart2, jend2)
2162 if( tile1(n) == tile2(n) ) cycle
2164 nx1_contact = iend1(n)-istart1(n)+1
2165 ny1_contact = jend1(n)-jstart1(n)+1
2166 nx2_contact = iend2(n)-istart2(n)+1
2167 ny2_contact = jend2(n)-jstart2(n)+1
2171 if( (nx(t1) .NE. nx1_contact .OR. ny(t1) .NE. ny1_contact ) .AND. &
2172 (nx(t2) .NE. nx2_contact .OR. ny(t2) .NE. ny2_contact ) ) cycle
2173 if(nx1_contact == nx2_contact .AND. ny1_contact == ny2_contact)
then
2174 call error_mesg(
'xgrid_mod',
'There is no refinement for the overlapping region', fatal)
2179 call error_mesg(
'xgrid_mod',
'only support one nest region, contact developer' ,fatal)
2181 if(nx2_contact*ny2_contact > nx1_contact*ny1_contact)
then
2182 is_nest_out = istart2(n);
2183 ie_nest_out = iend2(n);
2184 js_nest_out = jstart2(n);
2185 je_nest_out = jend2(n);
2186 tile_nest_out = tile2(n);
2187 is_parent_out = istart1(n);
2188 ie_parent_out = iend1(n);
2189 js_parent_out = jstart1(n);
2190 je_parent_out = jend1(n);
2191 tile_parent_out = tile1(n);
2193 is_nest_out = istart1(n);
2194 ie_nest_out = iend1(n);
2195 js_nest_out = jstart1(n);
2196 je_nest_out = jend1(n);
2197 tile_nest_out = tile1(n);
2198 is_parent_out = istart2(n);
2199 ie_parent_out = iend2(n);
2200 js_parent_out = jstart2(n);
2201 je_parent_out = jend2(n);
2202 tile_parent_out = tile2(n);
2206 deallocate(nx, ny, tile1, tile2)
2207 deallocate(istart1, iend1, jstart1, jend1)
2208 deallocate(istart2, iend2, jstart2, jend2)
2216 subroutine set_comm_get1_repro(xmap)
2218 integer,
dimension(xmap%npes) :: pe_ind, cnt
2219 integer,
dimension(0:xmap%npes-1) :: send_ind, pl
2220 integer :: npes, nsend, nrecv, mypos
2221 integer :: m, p, pos, n, g, l, im, i, j
2222 type(
comm_type),
pointer,
save :: comm => null()
2224 comm => xmap%get1_repro
2228 mypos =
mpp_pe() - mpp_root_pe()
2230 p = mod(mypos+npes-m, npes)
2231 if( xmap%recv_count_repro(p) > 0 )
then
2238 if( nrecv > 0 )
then
2239 if (
associated(comm%recv))
deallocate(comm%recv)
2240 allocate(comm%recv(nrecv))
2244 comm%recv(n)%count = xmap%recv_count_repro(p)
2245 comm%recv(n)%pe = p + xmap%root_pe
2246 comm%recv(n)%buffer_pos = pos
2247 pos = pos + comm%recv(n)%count
2254 mypos =
mpp_pe() - mpp_root_pe()
2256 p = mod(mypos+m, npes)
2257 if( xmap%send_count_repro(p) > 0 )
then
2265 if( nsend > 0 )
then
2266 if (
associated(comm%send))
deallocate(comm%send)
2267 allocate(comm%send(nsend))
2272 comm%send(n)%count = xmap%send_count_repro(p)
2273 comm%send(n)%pe = p + xmap%root_pe
2274 comm%send(n)%buffer_pos = pos
2275 pos = pos + comm%send(n)%count
2276 allocate(comm%send(n)%i(comm%send(n)%count))
2277 allocate(comm%send(n)%j(comm%send(n)%count))
2278 allocate(comm%send(n)%g(comm%send(n)%count))
2279 allocate(comm%send(n)%xLoc(comm%send(n)%count))
2282 do g=2,
size(xmap%grids(:))
2283 im = xmap%grids(g)%im
2284 do l=1,xmap%grids(g)%size
2285 p = xmap%grids(g)%x(l)%pe-xmap%root_pe
2289 i = xmap%grids(g)%x(l)%i2
2290 j = xmap%grids(g)%x(l)%j2
2291 if(xmap%grids(g)%is_ug)
then
2292 comm%send(n)%i(pos) = xmap%grids(g)%l_index((j-1)*im+i)
2293 comm%send(n)%j(pos) = 1
2295 comm%send(n)%i(pos) = xmap%grids(g)%x(l)%i2
2296 comm%send(n)%j(pos) = xmap%grids(g)%x(l)%j2
2298 comm%send(n)%g(pos) = g
2303 if( comm%send(n)%count .NE. cnt(n) )
call error_mesg(
'xgrid_mod', &
2304 .NE.
'comm%send(n)%count cnt(n)', fatal)
2310 do g=2,
size(xmap%grids(:))
2311 do l=1,xmap%grids(g)%size_repro
2312 p = xmap%grids(g)%x_repro(l)%pe-xmap%root_pe
2313 xmap%grids(g)%x_repro(l)%recv_pos = pl(p)
2320 end subroutine set_comm_get1_repro
2323 subroutine set_comm_get1(xmap)
2325 type (
grid_type),
pointer,
save :: grid1 =>null()
2326 integer,
allocatable :: send_size(:)
2327 integer,
allocatable :: recv_size(:)
2328 integer :: max_size, g, npes, l, ll, nset, m
2329 integer :: i1, j1, tile1, p, n, pos, buffer_pos, mypos
2330 integer :: nsend, nrecv, rbuf_size, sbuf_size, msgsize
2332 real(r8_kind),
allocatable :: recv_buf(:), send_buf(:)
2333 real(r8_kind),
allocatable :: diarray(:), djarray(:)
2334 integer,
allocatable :: iarray(:), jarray(:), tarray(:)
2335 integer,
allocatable :: pos_x(:), pelist(:), size_pe(:), pe_side1(:)
2336 integer :: recv_buffer_pos(0:xmap%npes)
2337 integer :: send_buffer_pos(0:xmap%npes)
2338 type(
comm_type),
pointer,
save :: comm => null()
2342 do g=2,
size(xmap%grids(:))
2343 max_size = max_size + xmap%grids(g)%size
2346 grid1 => xmap%grids(1)
2351 allocate(pelist(0:npes-1))
2352 call mpp_get_current_pelist(pelist)
2353 allocate(send_size(0:npes-1))
2354 allocate(recv_size(0:npes-1))
2355 allocate(size_pe(0:npes-1))
2356 allocate(pos_x(0:npes-1))
2361 if(max_size > 0)
then
2362 allocate(pe_side1(max_size))
2363 if (
associated(xmap%ind_get1))
deallocate(xmap%ind_get1)
2364 allocate(xmap%ind_get1(max_size))
2368 do g=2,
size(xmap%grids(:))
2369 do l=1,xmap%grids(g)%size
2370 i1 = xmap%grids(g)%x(l)%i1
2371 j1 = xmap%grids(g)%x(l)%j1
2372 tile1 = xmap%grids(g)%x(l)%tile
2374 if(grid1%tile(p) == tile1)
then
2376 size_pe(p) = size_pe(p) + 1
2381 if( p == npes )
then
2382 call error_mesg(
'xgrid_mod',
'tile is not in grid1%tile(:)', fatal)
2391 pos_x(p) = pos_x(p-1) + size_pe(p-1)
2395 allocate(iarray(max_size))
2396 allocate(jarray(max_size))
2397 allocate(tarray(max_size))
2398 if(monotonic_exchange)
then
2399 allocate(diarray(max_size))
2400 allocate(djarray(max_size))
2405 do g=2,
size(xmap%grids(:))
2406 do l=1,xmap%grids(g)%size
2407 i1 = xmap%grids(g)%x(l)%i1
2408 j1 = xmap%grids(g)%x(l)%j1
2409 tile1 = xmap%grids(g)%x(l)%tile
2414 if(send_size(p) > 0)
then
2415 if( i1 == iarray(pos_x(p)+send_size(p)) .AND. j1 == jarray(pos_x(p)+send_size(p)) &
2416 .AND. tile1 == tarray(pos_x(p)+send_size(p)))
then
2421 do n = 1, send_size(p)
2422 if(i1 == iarray(pos_x(p)+n) .AND. j1 == jarray(pos_x(p)+n) .AND. tile1 == tarray(pos_x(p)+n))
then
2429 if( (.NOT. found) .OR. monotonic_exchange )
then
2430 send_size(p) = send_size(p)+1
2431 pos = pos_x(p)+send_size(p)
2435 if(monotonic_exchange)
then
2436 diarray(pos) = xmap%grids(g)%x(l)%di
2437 djarray(pos) = xmap%grids(g)%x(l)%dj
2441 xmap%ind_get1(ll) = n
2447 pos_x(p) = pos_x(p-1) + send_size(p-1)
2451 do g=2,
size(xmap%grids(:))
2452 do l=1,xmap%grids(g)%size
2455 xmap%ind_get1(ll) = pos_x(p) + xmap%ind_get1(ll)
2460 mypos =
mpp_pe()-mpp_root_pe()
2463 recv_size(:) = xmap%your2my1_size(:)
2464 nsend = count( send_size> 0)
2467 if (
associated(comm%send))
deallocate(comm%send)
2468 allocate(comm%send(nsend))
2469 comm%send(:)%count = 0
2474 send_buffer_pos(p) = pos
2475 pos = pos + send_size(p)
2481 p = mod(mypos+n, npes)
2482 if(send_size(p)>0)
then
2484 allocate(comm%send(pos)%i(send_size(p)))
2485 comm%send(pos)%buffer_pos = send_buffer_pos(p)
2486 comm%send(pos)%count = send_size(p)
2487 comm%send(pos)%pe = pelist(p)
2488 comm%sendsize = comm%sendsize + send_size(p)
2493 if(monotonic_exchange) nset = 5
2494 rbuf_size = sum(recv_size)*nset
2495 sbuf_size = sum(send_size)*nset
2496 if(rbuf_size>0)
allocate(recv_buf(rbuf_size))
2497 if(sbuf_size>0)
allocate(send_buf(sbuf_size))
2501 p = mod(mypos+npes-n, npes)
2502 if(recv_size(p) ==0) cycle
2503 msgsize = recv_size(p)*nset
2504 call mpp_recv(recv_buf(pos+1), glen=msgsize, from_pe=pelist(p), block=.false., tag=comm_tag_4)
2510 pos_x(p) = pos_x(p-1) + size_pe(p-1)
2515 p = mod(mypos+n, npes)
2516 do l = 1, send_size(p)
2517 send_buf(pos+1) = real(iarray(pos_x(p)+l), r8_kind)
2518 send_buf(pos+2) = real(jarray(pos_x(p)+l), r8_kind)
2519 send_buf(pos+3) = real(tarray(pos_x(p)+l), r8_kind)
2520 if(monotonic_exchange)
then
2521 send_buf(pos+4) = diarray(pos_x(p)+l)
2522 send_buf(pos+5) = djarray(pos_x(p)+l)
2530 p = mod(mypos+n, npes)
2531 if(send_size(p) ==0) cycle
2532 msgsize = send_size(p)*nset
2533 call mpp_send(send_buf(pos+1), plen=msgsize, to_pe=pelist(p), tag=comm_tag_4 )
2538 nrecv = count(recv_size>0)
2543 if (
associated(comm%recv))
deallocate(comm%recv)
2544 allocate(comm%recv(nrecv))
2545 comm%recv(:)%count = 0
2549 recv_buffer_pos(p) = buffer_pos
2550 buffer_pos = buffer_pos + recv_size(p)
2555 p = mod(mypos+npes-m, npes)
2556 if(recv_size(p)>0)
then
2558 allocate(comm%recv(pos)%i(recv_size(p)))
2559 allocate(comm%recv(pos)%j(recv_size(p)))
2560 allocate(comm%recv(pos)%tile(recv_size(p)))
2561 comm%recv(pos)%buffer_pos = recv_buffer_pos(p)
2562 comm%recv(pos)%pe = pelist(p)
2563 comm%recv(pos)%count = recv_size(p)
2564 comm%recvsize = comm%recvsize + recv_size(p)
2565 if(monotonic_exchange)
then
2566 allocate(comm%recv(pos)%di(recv_size(p)))
2567 allocate(comm%recv(pos)%dj(recv_size(p)))
2569 if(grid1%is_ug)
then
2570 do n = 1, recv_size(p)
2571 i = int(recv_buf(buffer_pos+1))
2572 j = int(recv_buf(buffer_pos+2))
2573 comm%recv(pos)%i(n) = grid1%l_index((j-1)*grid1%im+i)
2574 comm%recv(pos)%j(n) = 1
2575 comm%recv(pos)%tile(n) = int(recv_buf(buffer_pos+3))
2576 if(monotonic_exchange)
then
2577 comm%recv(pos)%di(n) = recv_buf(buffer_pos+4)
2578 comm%recv(pos)%dj(n) = recv_buf(buffer_pos+5)
2580 buffer_pos = buffer_pos + nset
2583 do n = 1, recv_size(p)
2584 comm%recv(pos)%i(n) = int(recv_buf(buffer_pos+1) )- grid1%is_me + 1
2585 comm%recv(pos)%j(n) = int(recv_buf(buffer_pos+2) )- grid1%js_me + 1
2586 comm%recv(pos)%tile(n) = int(recv_buf(buffer_pos+3))
2587 if(monotonic_exchange)
then
2588 comm%recv(pos)%di(n) = recv_buf(buffer_pos+4)
2589 comm%recv(pos)%dj(n) = recv_buf(buffer_pos+5)
2591 buffer_pos = buffer_pos + nset
2596 if (
associated(comm%unpack_ind))
deallocate(comm%unpack_ind)
2597 allocate(comm%unpack_ind(nrecv))
2600 if(recv_size(p)>0)
then
2603 if(comm%recv(m)%pe == pelist(p))
then
2604 comm%unpack_ind(pos) = m
2613 if(
allocated(send_buf) )
deallocate(send_buf)
2614 if(
allocated(recv_buf) )
deallocate(recv_buf)
2615 if(
allocated(pelist) )
deallocate(pelist)
2616 if(
allocated(pos_x) )
deallocate(pos_x)
2617 if(
allocated(pelist) )
deallocate(pelist)
2618 if(
allocated(iarray) )
deallocate(iarray)
2619 if(
allocated(jarray) )
deallocate(jarray)
2620 if(
allocated(tarray) )
deallocate(tarray)
2621 if(
allocated(size_pe) )
deallocate(size_pe)
2623 end subroutine set_comm_get1
2626 subroutine set_comm_put1(xmap)
2628 type (
grid_type),
pointer,
save :: grid1 =>null()
2629 integer,
allocatable :: send_size(:)
2630 integer,
allocatable :: recv_size(:)
2631 integer :: max_size, g, npes, l, ll, m, mypos
2632 integer :: i1, j1, tile1, p, n, pos, buffer_pos
2633 integer :: nsend, nrecv, msgsize, nset, rbuf_size, sbuf_size
2635 real(r8_kind),
allocatable :: recv_buf(:), send_buf(:)
2636 real(r8_kind),
allocatable :: diarray(:), djarray(:)
2637 integer,
allocatable :: iarray(:), jarray(:), tarray(:)
2638 integer,
allocatable :: pos_x(:), pelist(:), size_pe(:), pe_put1(:)
2639 integer :: recv_buffer_pos(0:xmap%npes)
2640 type(
comm_type),
pointer,
save :: comm => null()
2644 if(nnest == 0 .OR. xmap%grids(1)%id .NE.
'ATM' )
then
2645 comm%nsend = xmap%get1%nrecv
2646 comm%nrecv = xmap%get1%nsend
2647 comm%sendsize = xmap%get1%recvsize
2648 comm%recvsize = xmap%get1%sendsize
2649 comm%send => xmap%get1%recv
2650 comm%recv => xmap%get1%send
2651 xmap%ind_put1 => xmap%ind_get1
2656 do g=2,
size(xmap%grids(:))
2657 max_size = max_size + xmap%grids(g)%size
2659 grid1 => xmap%grids(1)
2663 allocate(pelist(0:npes-1))
2664 call mpp_get_current_pelist(pelist)
2665 allocate(send_size(0:npes-1))
2666 allocate(recv_size(0:npes-1))
2667 allocate(size_pe(0:npes-1))
2668 allocate(pos_x(0:npes-1))
2673 if(max_size > 0)
then
2674 allocate(pe_put1(max_size))
2675 if (
associated(xmap%ind_put1))
deallocate(xmap%ind_put1)
2676 allocate(xmap%ind_put1(max_size))
2680 do g=2,
size(xmap%grids(:))
2681 do l=1,xmap%grids(g)%size
2682 i1 = xmap%grids(g)%x(l)%i1
2683 j1 = xmap%grids(g)%x(l)%j1
2684 tile1 = xmap%grids(g)%x(l)%tile
2686 if(grid1%tile(p) == tile1)
then
2687 if(
in_box(i1, j1, grid1%is(p), grid1%ie(p), grid1%js(p), grid1%je(p)))
then
2688 size_pe(p) = size_pe(p) + 1
2700 pos_x(p) = pos_x(p-1) + size_pe(p-1)
2704 allocate(iarray(max_size))
2705 allocate(jarray(max_size))
2706 allocate(tarray(max_size))
2707 if(monotonic_exchange)
then
2708 allocate(diarray(max_size))
2709 allocate(djarray(max_size))
2714 do g=2,
size(xmap%grids(:))
2715 do l=1,xmap%grids(g)%size
2716 i1 = xmap%grids(g)%x(l)%i1
2717 j1 = xmap%grids(g)%x(l)%j1
2718 tile1 = xmap%grids(g)%x(l)%tile
2723 if(send_size(p) > 0)
then
2724 if( i1 == iarray(pos_x(p)+send_size(p)) .AND. j1 == jarray(pos_x(p)+send_size(p)) &
2725 .AND. tile1 == tarray(pos_x(p)+send_size(p)))
then
2730 do n = 1, send_size(p)
2731 if(i1 == iarray(pos_x(p)+n) .AND. j1 == jarray(pos_x(p)+n) .AND. tile1 == tarray(pos_x(p)+n))
then
2738 if( (.NOT. found) .OR. monotonic_exchange )
then
2739 send_size(p) = send_size(p)+1
2740 pos = pos_x(p)+send_size(p)
2744 if(monotonic_exchange)
then
2745 diarray(pos) = xmap%grids(g)%x(l)%di
2746 djarray(pos) = xmap%grids(g)%x(l)%dj
2750 xmap%ind_put1(ll) = n
2756 pos_x(p) = pos_x(p-1) + send_size(p-1)
2760 do g=2,
size(xmap%grids(:))
2761 do l=1,xmap%grids(g)%size
2762 i1 = xmap%grids(g)%x(l)%i1
2763 j1 = xmap%grids(g)%x(l)%j1
2764 tile1 = xmap%grids(g)%x(l)%tile
2767 xmap%ind_put1(ll) = pos_x(p) + xmap%ind_put1(ll)
2773 mypos =
mpp_pe()-mpp_root_pe()
2775 if (do_alltoall)
then
2776 call mpp_alltoall(send_size, 1, recv_size, 1)
2779 p = mod(mypos+npes-n, npes)
2780 call mpp_recv(recv_size(p), glen=1, from_pe=pelist(p), block=.false., tag=comm_tag_5)
2785 p = mod(mypos+n, npes)
2786 call mpp_send(send_size(p), plen=1, to_pe=pelist(p), tag=comm_tag_5)
2793 nrecv = count( send_size> 0)
2796 if (
associated(comm%recv))
deallocate(comm%recv)
2797 allocate(comm%recv(nrecv))
2798 comm%recv(:)%count = 0
2803 recv_buffer_pos(p) = pos
2804 pos = pos + send_size(p)
2809 p = mod(mypos+npes-n, npes)
2810 if(send_size(p)>0)
then
2812 allocate(comm%recv(pos)%i(send_size(p)))
2813 comm%recv(pos)%buffer_pos = recv_buffer_pos(p)
2814 comm%recv(pos)%count = send_size(p)
2815 comm%recv(pos)%pe = pelist(p)
2816 comm%recvsize = comm%recvsize + send_size(p)
2821 if(monotonic_exchange) nset = 5
2822 rbuf_size = sum(recv_size)*nset
2823 sbuf_size = sum(send_size)*nset
2824 if(rbuf_size>0)
allocate(recv_buf(rbuf_size))
2825 if(sbuf_size>0)
allocate(send_buf(sbuf_size))
2829 p = mod(mypos+npes-n, npes)
2830 if(recv_size(p) ==0) cycle
2831 msgsize = recv_size(p)*nset
2832 call mpp_recv(recv_buf(pos+1), glen=msgsize, from_pe=pelist(p), block=.false., tag=comm_tag_6)
2838 pos_x(p) = pos_x(p-1) + size_pe(p-1)
2843 p = mod(mypos+n, npes)
2844 do l = 1, send_size(p)
2845 send_buf(pos+1) = real(iarray(pos_x(p)+l), r8_kind)
2846 send_buf(pos+2) = real(jarray(pos_x(p)+l), r8_kind)
2847 send_buf(pos+3) = real(tarray(pos_x(p)+l), r8_kind)
2848 if(monotonic_exchange)
then
2849 send_buf(pos+4) = diarray(pos_x(p)+l)
2850 send_buf(pos+5) = djarray(pos_x(p)+l)
2858 p = mod(mypos+n, npes)
2859 if(send_size(p) ==0) cycle
2860 msgsize = send_size(p)*nset
2861 call mpp_send(send_buf(pos+1), plen=msgsize, to_pe=pelist(p), tag=comm_tag_6 )
2866 nsend = count(recv_size>0)
2871 if (
associated(comm%send))
deallocate(comm%send)
2872 allocate(comm%send(nsend))
2873 comm%send(:)%count = 0
2877 p = mod(mypos+npes-m, npes)
2878 if(recv_size(p)>0)
then
2880 allocate(comm%send(pos)%i(recv_size(p)))
2881 allocate(comm%send(pos)%j(recv_size(p)))
2882 allocate(comm%send(pos)%tile(recv_size(p)))
2883 comm%send(pos)%pe = pelist(p)
2884 comm%send(pos)%count = recv_size(p)
2885 comm%sendsize = comm%sendsize + recv_size(p)
2886 if(monotonic_exchange)
then
2887 allocate(comm%send(pos)%di(recv_size(p)))
2888 allocate(comm%send(pos)%dj(recv_size(p)))
2890 do n = 1, recv_size(p)
2891 comm%send(pos)%i(n) = int(recv_buf(buffer_pos+1) )- grid1%is_me + 1
2892 comm%send(pos)%j(n) = int(recv_buf(buffer_pos+2) )- grid1%js_me + 1
2893 comm%send(pos)%tile(n) = int(recv_buf(buffer_pos+3))
2894 if(monotonic_exchange)
then
2895 comm%send(pos)%di(n) = recv_buf(buffer_pos+4)
2896 comm%send(pos)%dj(n) = recv_buf(buffer_pos+5)
2898 buffer_pos = buffer_pos + nset
2905 if(
allocated(send_buf) )
deallocate(send_buf)
2906 if(
allocated(recv_buf) )
deallocate(recv_buf)
2907 if(
allocated(pelist) )
deallocate(pelist)
2908 if(
allocated(pos_x) )
deallocate(pos_x)
2909 if(
allocated(pelist) )
deallocate(pelist)
2910 if(
allocated(iarray) )
deallocate(iarray)
2911 if(
allocated(jarray) )
deallocate(jarray)
2912 if(
allocated(tarray) )
deallocate(tarray)
2913 if(
allocated(size_pe) )
deallocate(size_pe)
2915 end subroutine set_comm_put1
2939 type (xmap_type),
intent(inout) :: xmap
2941 integer :: g, l, k, max_size
2942 integer :: i1, j1, i2, j2, p
2945 logical :: overlap_with_nest
2946 integer :: cnt(xmap%get1%nsend)
2947 integer :: i,j,n,xloc,pos,nsend,m,npes, mypos
2948 integer :: send_ind(0:xmap%npes-1)
2952 do g=2,
size(xmap%grids(:))
2953 max_size = max_size + xmap%grids(g)%size * xmap%grids(g)%km
2956 if (max_size>
size(xmap%x1(:)))
then
2957 if (
associated(xmap%x1))
deallocate(xmap%x1)
2958 if (
associated(xmap%x2))
deallocate(xmap%x2)
2959 allocate( xmap%x1(1:max_size) )
2960 allocate( xmap%x2(1:max_size) )
2964 do g=2,
size(xmap%grids(:))
2965 xmap%grids(g)%first = 1
2966 xmap%grids(g)%last = 0
2971 do g=2,
size(xmap%grids(:))
2972 xmap%grids(g)%first = xmap%size + 1;
2974 do l=1,xmap%grids(g)%size
2975 i1 = xmap%grids(g)%x(l)%i1
2976 j1 = xmap%grids(g)%x(l)%j1
2977 i2 = xmap%grids(g)%x(l)%i2
2978 j2 = xmap%grids(g)%x(l)%j2
2979 tile1 = xmap%grids(g)%x(l)%tile
2981 if(xmap%grids(g)%is_ug)
then
2982 do k=1,xmap%grids(g)%km
2983 lll = xmap%grids(g)%l_index((j2-1)*xmap%grids(g)%im+i2)
2984 if (xmap%grids(g)%frac_area(lll,1,k)/=0.0_r8_kind)
then
2985 xmap%size = xmap%size+1
2986 xmap%x1(xmap%size)%pos = xmap%ind_get1(ll)
2987 xmap%x1(xmap%size)%i = xmap%grids(g)%x(l)%i1
2988 xmap%x1(xmap%size)%j = xmap%grids(g)%x(l)%j1
2989 xmap%x1(xmap%size)%tile = xmap%grids(g)%x(l)%tile
2990 xmap%x1(xmap%size)%area = xmap%grids(g)%x(l)%area &
2991 *xmap%grids(g)%frac_area(lll,1,k)
2992 xmap%x1(xmap%size)%di = xmap%grids(g)%x(l)%di
2993 xmap%x1(xmap%size)%dj = xmap%grids(g)%x(l)%dj
2994 xmap%x2(xmap%size)%i = xmap%grids(g)%x(l)%i2
2995 xmap%x2(xmap%size)%j = xmap%grids(g)%x(l)%j2
2996 xmap%x2(xmap%size)%l = lll
2997 xmap%x2(xmap%size)%k = k
2998 xmap%x2(xmap%size)%area = xmap%grids(g)%x(l)%area * xmap%grids(g)%x(l)%scale
3002 do k=1,xmap%grids(g)%km
3003 if (xmap%grids(g)%frac_area(i2,j2,k)/=0.0_r8_kind)
then
3004 xmap%size = xmap%size+1
3005 xmap%x1(xmap%size)%pos = xmap%ind_get1(ll)
3006 xmap%x1(xmap%size)%i = xmap%grids(g)%x(l)%i1
3007 xmap%x1(xmap%size)%j = xmap%grids(g)%x(l)%j1
3008 xmap%x1(xmap%size)%tile = xmap%grids(g)%x(l)%tile
3009 xmap%x1(xmap%size)%area = xmap%grids(g)%x(l)%area &
3010 *xmap%grids(g)%frac_area(i2,j2,k)
3011 xmap%x1(xmap%size)%di = xmap%grids(g)%x(l)%di
3012 xmap%x1(xmap%size)%dj = xmap%grids(g)%x(l)%dj
3013 xmap%x2(xmap%size)%i = xmap%grids(g)%x(l)%i2
3014 xmap%x2(xmap%size)%j = xmap%grids(g)%x(l)%j2
3015 xmap%x2(xmap%size)%k = k
3016 xmap%x2(xmap%size)%area = xmap%grids(g)%x(l)%area * xmap%grids(g)%x(l)%scale
3021 xmap%grids(g)%last = xmap%size
3025 if (max_size>
size(xmap%x1_put(:)))
then
3026 if (
associated(xmap%x1_put))
deallocate(xmap%x1_put)
3027 allocate( xmap%x1_put(1:max_size) )
3029 if (max_size>
size(xmap%x2_get(:)))
then
3030 if (
associated(xmap%x2_get))
deallocate(xmap%x2_get)
3031 allocate( xmap%x2_get(1:max_size) )
3034 do g=2,
size(xmap%grids(:))
3035 xmap%grids(g)%first_get = 1
3036 xmap%grids(g)%last_get = 0
3042 do g=2,
size(xmap%grids(:))
3043 xmap%grids(g)%first_get = xmap%size_get2 + 1;
3045 do l=1,xmap%grids(g)%size
3046 i1 = xmap%grids(g)%x(l)%i1
3047 j1 = xmap%grids(g)%x(l)%j1
3048 i2 = xmap%grids(g)%x(l)%i2
3049 j2 = xmap%grids(g)%x(l)%j2
3050 tile1 = xmap%grids(g)%x(l)%tile
3052 overlap_with_nest = .false.
3053 if( xmap%grids(1)%id ==
"ATM" .AND. tile1 == tile_parent .AND. &
3054 in_box(i1, j1, is_parent, ie_parent, js_parent, je_parent) ) overlap_with_nest = .true.
3055 if(xmap%grids(g)%is_ug)
then
3056 do k=1,xmap%grids(g)%km
3057 lll = xmap%grids(g)%l_index((j2-1)*xmap%grids(g)%im+i2)
3058 if (xmap%grids(g)%frac_area(lll,1,k)/=0.0_r8_kind)
then
3059 xmap%size_put1 = xmap%size_put1+1
3060 xmap%x1_put(xmap%size_put1)%pos = xmap%ind_put1(ll)
3061 xmap%x1_put(xmap%size_put1)%i = xmap%grids(g)%x(l)%i1
3062 xmap%x1_put(xmap%size_put1)%j = xmap%grids(g)%x(l)%j1
3063 xmap%x1_put(xmap%size_put1)%tile = xmap%grids(g)%x(l)%tile
3064 xmap%x1_put(xmap%size_put1)%area = xmap%grids(g)%x(l)%area &
3065 *xmap%grids(g)%frac_area(lll,1,k)
3066 xmap%x1_put(xmap%size_put1)%di = xmap%grids(g)%x(l)%di
3067 xmap%x1_put(xmap%size_put1)%dj = xmap%grids(g)%x(l)%dj
3068 if( .not. overlap_with_nest)
then
3069 xmap%size_get2 = xmap%size_get2+1
3070 xmap%x2_get(xmap%size_get2)%i = xmap%grids(g)%x(l)%i2
3071 xmap%x2_get(xmap%size_get2)%j = xmap%grids(g)%x(l)%j2
3072 xmap%x2_get(xmap%size_get2)%l = lll
3073 xmap%x2_get(xmap%size_get2)%k = k
3074 xmap%x2_get(xmap%size_get2)%area = xmap%grids(g)%x(l)%area * xmap%grids(g)%x(l)%scale
3075 xmap%x2_get(xmap%size_get2)%pos = xmap%size_put1
3080 do k=1,xmap%grids(g)%km
3081 if (xmap%grids(g)%frac_area(i2,j2,k)/=0.0_r8_kind)
then
3082 xmap%size_put1 = xmap%size_put1+1
3083 xmap%x1_put(xmap%size_put1)%pos = xmap%ind_put1(ll)
3084 xmap%x1_put(xmap%size_put1)%i = xmap%grids(g)%x(l)%i1
3085 xmap%x1_put(xmap%size_put1)%j = xmap%grids(g)%x(l)%j1
3086 xmap%x1_put(xmap%size_put1)%tile = xmap%grids(g)%x(l)%tile
3087 xmap%x1_put(xmap%size_put1)%area = xmap%grids(g)%x(l)%area &
3088 *xmap%grids(g)%frac_area(i2,j2,k)
3089 xmap%x1_put(xmap%size_put1)%di = xmap%grids(g)%x(l)%di
3090 xmap%x1_put(xmap%size_put1)%dj = xmap%grids(g)%x(l)%dj
3091 if( .not. overlap_with_nest)
then
3092 xmap%size_get2 = xmap%size_get2+1
3093 xmap%x2_get(xmap%size_get2)%i = xmap%grids(g)%x(l)%i2
3094 xmap%x2_get(xmap%size_get2)%j = xmap%grids(g)%x(l)%j2
3095 xmap%x2_get(xmap%size_get2)%k = k
3096 xmap%x2_get(xmap%size_get2)%area = xmap%grids(g)%x(l)%area * xmap%grids(g)%x(l)%scale
3097 xmap%x2_get(xmap%size_get2)%pos = xmap%size_put1
3103 xmap%grids(g)%last_get = xmap%size_get2
3108 if (xmap%get1_repro%nsend > 0)
then
3112 mypos =
mpp_pe() - mpp_root_pe()
3115 p = mod(mypos+m, npes)
3116 if( xmap%send_count_repro(p) > 0 )
then
3121 do g=2,
size(xmap%grids(:))
3122 do l=1,xmap%grids(g)%size
3123 p = xmap%grids(g)%x(l)%pe-xmap%root_pe
3127 xmap%get1_repro%send(n)%xLoc(pos) = xloc
3128 if( xmap%grids(g)%is_ug )
then
3129 i = xmap%grids(g)%x(l)%l2
3130 xloc = xloc + count(xmap%grids(g)%frac_area(i,1,:)/=0.0_r8_kind)
3132 i = xmap%grids(g)%x(l)%i2
3133 j = xmap%grids(g)%x(l)%j2
3134 xloc = xloc + count(xmap%grids(g)%frac_area(i,j,:)/=0.0_r8_kind)
3141 end subroutine regen
3172 real(r8_kind),
dimension(:,:,:),
intent(in) :: f
3173 character(len=3),
intent(in) :: grid_id
3174 type (xmap_type),
intent(inout) :: xmap
3177 type(
grid_type),
pointer,
save :: grid =>null()
3179 if (grid_id==xmap%grids(1)%id)
call error_mesg (
'xgrid_mod', &
3180 'set_frac_area called on side 1 grid', fatal)
3181 do g=2,
size(xmap%grids(:))
3182 grid => xmap%grids(g)
3183 if (grid_id==grid%id)
then
3184 if (
size(f,3)/=
size(grid%frac_area,3))
then
3185 if (
associated(grid%frac_area))
deallocate (grid%frac_area)
3186 grid%km =
size(f,3);
3187 allocate( grid%frac_area(grid%is_me:grid%ie_me, grid%js_me:grid%je_me, &
3196 call error_mesg (
'xgrid_mod',
'set_frac_area: could not find grid id', fatal)
3204 real(r8_kind),
dimension(:,:),
intent(in) :: f
3205 character(len=3),
intent(in) :: grid_id
3206 type (xmap_type),
intent(inout) :: xmap
3209 type(
grid_type),
pointer,
save :: grid =>null()
3211 if (grid_id==xmap%grids(1)%id)
call error_mesg (
'xgrid_mod', &
3212 'set_frac_area_ug called on side 1 grid', fatal)
3213 if (grid_id .NE.
'LND' )
call error_mesg (
'xgrid_mod', &
3214 .NE.
'set_frac_area_ug called for grid_id LND', fatal)
3215 do g=2,
size(xmap%grids(:))
3216 grid => xmap%grids(g)
3217 if (grid_id==grid%id)
then
3218 if (
size(f,2)/=
size(grid%frac_area,3))
then
3219 if (
associated(grid%frac_area))
deallocate (grid%frac_area)
3220 grid%km =
size(f,2);
3221 allocate( grid%frac_area(grid%ls_me:grid%le_me, 1, grid%km) )
3223 grid%frac_area(:,1,:) = f(:,:);
3229 call error_mesg (
'xgrid_mod',
'set_frac_area_ug: could not find grid id', fatal)
3247 real(r8_kind),
dimension(:,:),
intent(in) :: d
3248 character(len=3),
intent(in) :: grid_id
3249 real(r8_kind),
dimension(:),
intent(inout) :: x
3250 type (xmap_type),
intent(inout) :: xmap
3251 integer,
intent(in),
optional :: remap_method
3253 logical,
intent(in),
optional :: complete
3255 logical :: is_complete, set_mismatch
3256 integer :: g, method
3257 character(len=2) :: text
3258 integer,
save :: isize=0
3259 integer,
save :: jsize=0
3260 integer,
save :: lsize=0
3261 integer,
save :: xsize=0
3262 integer,
save :: method_saved=0
3263 character(len=3),
save :: grid_id_saved=
""
3264 integer(i8_kind),
dimension(MAX_FIELDS),
save :: d_addrs = -9999_i8_kind
3265 integer(i8_kind),
dimension(MAX_FIELDS),
save :: x_addrs = -9999_i8_kind
3267 if (grid_id==xmap%grids(1)%id)
then
3268 method = first_order
3269 if(
present(remap_method)) method = remap_method
3270 is_complete = .true.
3271 if(
present(complete)) is_complete=complete
3273 if( lsize > max_fields )
then
3274 write( text,
'(i2)' ) max_fields
3275 call error_mesg (
'xgrid_mod',
'MAX_FIELDS='//trim(text)//
' exceeded for group put_side1_to_xgrid', fatal)
3277 d_addrs(lsize) = loc(d)
3278 x_addrs(lsize) = loc(x)
3284 method_saved = method
3285 grid_id_saved = grid_id
3287 set_mismatch = .false.
3288 set_mismatch = set_mismatch .OR. (isize /=
size(d,1))
3289 set_mismatch = set_mismatch .OR. (jsize /=
size(d,2))
3290 set_mismatch = set_mismatch .OR. (xsize /=
size(x(:)))
3291 set_mismatch = set_mismatch .OR. (method_saved /= method)
3292 set_mismatch = set_mismatch .OR. (grid_id_saved /= grid_id)
3293 if(set_mismatch)
then
3294 write( text,
'(i2)' ) lsize
3295 call error_mesg (
'xgrid_mod',
'Incompatible field at count '//text//
' for group put_side1_to_xgrid', fatal )
3299 if(is_complete)
then
3302 if(monotonic_exchange .AND. grid_id ==
'ATM')
then
3303 call put_1_to_xgrid_order_2(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3304 else if(method == first_order)
then
3305 call put_1_to_xgrid_order_1(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3307 if(grid_id .NE.
'ATM')
call error_mesg (
'xgrid_mod', &
3308 "second order put_to_xgrid should only be applied to 'ATM' model, "//&
3309 "contact developer", fatal)
3310 call put_1_to_xgrid_order_2(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3313 d_addrs = -9999_i8_kind
3314 x_addrs = -9999_i8_kind
3325 do g=2,
size(xmap%grids(:))
3326 if (grid_id==xmap%grids(g)%id) &
3327 call error_mesg (
'xgrid_mod', &
3328 'put_to_xgrid expects a 3D side 2 grid', fatal)
3331 call error_mesg (
'xgrid_mod',
'put_to_xgrid: could not find grid id', fatal)
3339 real(r8_kind),
dimension(:,:,:),
intent(in) :: d
3340 character(len=3),
intent(in) :: grid_id
3341 real(r8_kind),
dimension(:),
intent(inout) :: x
3342 type (xmap_type),
intent(inout) :: xmap
3346 if (grid_id==xmap%grids(1)%id) &
3347 call error_mesg (
'xgrid_mod', &
3348 'put_side2_to_xgrid expects a 3D side 2 grid', fatal)
3350 do g=2,
size(xmap%grids(:))
3351 if (grid_id==xmap%grids(g)%id)
then
3352 call put_2_to_xgrid(d, xmap%grids(g), x, xmap)
3357 call error_mesg (
'xgrid_mod',
'put_to_xgrid: could not find grid id', fatal)
3364 real(r8_kind),
dimension(:,:),
intent(out) :: d
3365 character(len=3),
intent(in) :: grid_id
3366 real(r8_kind),
dimension(:),
intent(in) :: x
3367 type (xmap_type),
intent(inout) :: xmap
3368 logical,
intent(in),
optional :: complete
3370 logical :: is_complete, set_mismatch
3372 character(len=2) :: text
3373 integer,
save :: isize=0
3374 integer,
save :: jsize=0
3375 integer,
save :: lsize=0
3376 integer,
save :: xsize=0
3377 character(len=3),
save :: grid_id_saved=
""
3378 integer(i8_kind),
dimension(MAX_FIELDS),
save :: d_addrs = -9999_i8_kind
3379 integer(i8_kind),
dimension(MAX_FIELDS),
save :: x_addrs = -9999_i8_kind
3382 if (grid_id==xmap%grids(1)%id)
then
3383 is_complete = .true.
3384 if(
present(complete)) is_complete=complete
3386 if( lsize > max_fields )
then
3387 write( text,
'(i2)' ) max_fields
3388 call error_mesg (
'xgrid_mod',
'MAX_FIELDS='//trim(text)//
' exceeded for group get_side1_from_xgrid', fatal)
3390 d_addrs(lsize) = loc(d)
3391 x_addrs(lsize) = loc(x)
3397 grid_id_saved = grid_id
3399 set_mismatch = .false.
3400 set_mismatch = set_mismatch .OR. (isize /=
size(d,1))
3401 set_mismatch = set_mismatch .OR. (jsize /=
size(d,2))
3402 set_mismatch = set_mismatch .OR. (xsize /=
size(x(:)))
3403 set_mismatch = set_mismatch .OR. (grid_id_saved /= grid_id)
3404 if(set_mismatch)
then
3405 write( text,
'(i2)' ) lsize
3406 call error_mesg (
'xgrid_mod',
'Incompatible field at count '//text// &
3407 &
' for group get_side1_from_xgrid', fatal )
3411 if(is_complete)
then
3413 call get_1_from_xgrid_repro(d_addrs, x_addrs, xmap, xsize, lsize)
3415 call get_1_from_xgrid(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3417 d_addrs(1:lsize) = -9999
3418 x_addrs(1:lsize) = -9999
3428 do g=2,
size(xmap%grids(:))
3429 if (grid_id==xmap%grids(g)%id) &
3430 call error_mesg (
'xgrid_mod', &
3431 'get_from_xgrid expects a 3D side 2 grid', fatal)
3434 call error_mesg (
'xgrid_mod',
'get_from_xgrid: could not find grid id', fatal)
3441 real(r8_kind),
dimension(:,:,:),
intent(out) :: d
3442 character(len=3),
intent(in) :: grid_id
3443 real(r8_kind),
dimension(:),
intent(in) :: x
3444 type (xmap_type),
intent(in) :: xmap
3448 if (grid_id==xmap%grids(1)%id) &
3449 call error_mesg (
'xgrid_mod', &
3450 'get_from_xgrid expects a 2D side 1 grid', fatal)
3452 do g=2,
size(xmap%grids(:))
3453 if (grid_id==xmap%grids(g)%id)
then
3454 call get_2_from_xgrid(d, xmap%grids(g), x, xmap)
3459 call error_mesg (
'xgrid_mod',
'get_from_xgrid: could not find grid id', fatal)
3466 subroutine some(xmap, some_arr, grid_id)
3468 character(len=3),
optional,
intent(in) :: grid_id
3469 logical,
dimension(:),
intent(out) :: some_arr
3473 if (.not.
present(grid_id))
then
3475 if(xmap%size > 0)
then
3483 if (grid_id==xmap%grids(1)%id) &
3484 call error_mesg (
'xgrid_mod',
'some expects a side 2 grid id', fatal)
3486 do g=2,
size(xmap%grids(:))
3487 if (grid_id==xmap%grids(g)%id)
then
3489 some_arr(xmap%grids(g)%first:xmap%grids(g)%last) = .true.;
3494 call error_mesg (
'xgrid_mod',
'some could not find grid id', fatal)
3500 subroutine put_2_to_xgrid(d, grid, x, xmap)
3502 real(r8_kind),
dimension(grid%is_me:grid%ie_me, &
grid%js_me:grid%je_me, grid%km),
intent(in) :: d
3503 real(r8_kind),
dimension(:),
intent(inout) :: x
3507 call mpp_clock_begin(id_put_2_to_xgrid)
3509 do l=grid%first,grid%last
3510 x(l) = d(xmap%x2(l)%i,xmap%x2(l)%j,xmap%x2(l)%k)
3513 call mpp_clock_end(id_put_2_to_xgrid)
3514 end subroutine put_2_to_xgrid
3518 subroutine get_2_from_xgrid(d, grid, x, xmap)
3520 real(r8_kind),
dimension(grid%is_me:grid%ie_me, &
grid%js_me:grid%je_me, grid%km),
intent(out) :: d
3521 real(r8_kind),
dimension(:),
intent(in) :: x
3526 call mpp_clock_begin(id_get_2_from_xgrid)
3529 do l=grid%first_get,grid%last_get
3530 d(xmap%x2_get(l)%i,xmap%x2_get(l)%j,xmap%x2_get(l)%k) = &
3531 d(xmap%x2_get(l)%i,xmap%x2_get(l)%j,xmap%x2_get(l)%k) + xmap%x2_get(l)%area*x(xmap%x2_get(l)%pos)
3537 d(:,:,k) = d(:,:,k) * grid%area_inv
3540 call mpp_clock_end(id_get_2_from_xgrid)
3542 end subroutine get_2_from_xgrid
3546 subroutine put_1_to_xgrid_order_1(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3547 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
3548 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
3550 integer,
intent(in) :: isize, jsize, xsize, lsize
3552 integer :: i, j, p, buffer_pos, msgsize
3553 integer :: from_pe, to_pe, pos, n, l, count
3554 integer :: ibegin, istart, iend, start_pos
3555 type (
comm_type),
pointer,
save :: comm =>null()
3556 real(r8_kind) :: recv_buffer(xmap%put1%recvsize*lsize)
3557 real(r8_kind) :: send_buffer(xmap%put1%sendsize*lsize)
3558 real(r8_kind) :: unpack_buffer(xmap%put1%recvsize)
3560 real(r8_kind),
dimension(isize, jsize) :: d
3561 real(r8_kind),
dimension(xsize) :: x
3565 call mpp_clock_begin(id_put_1_to_xgrid_order_1)
3569 do p = 1, comm%nrecv
3570 msgsize = comm%recv(p)%count*lsize
3571 from_pe = comm%recv(p)%pe
3572 buffer_pos = comm%recv(p)%buffer_pos*lsize
3573 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = from_pe, block=.false., tag=comm_tag_7)
3578 do p = 1, comm%nsend
3579 msgsize = comm%send(p)%count*lsize
3580 to_pe = comm%send(p)%pe
3584 do n = 1, comm%send(p)%count
3586 i = comm%send(p)%i(n)
3587 j = comm%send(p)%j(n)
3588 send_buffer(pos) = d(i,j)
3591 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = to_pe, tag=comm_tag_7 )
3592 buffer_pos = buffer_pos + msgsize
3598 if( lsize == 1)
then
3600 do l=1,xmap%size_put1
3601 x(l) = recv_buffer(xmap%x1_put(l)%pos)
3609 do p = 1, comm%nrecv
3610 count = comm%recv(p)%count
3611 ibegin = comm%recv(p)%buffer_pos*lsize + 1
3612 istart = ibegin + (l-1)*count
3613 iend = istart + count - 1
3614 pos = comm%recv(p)%buffer_pos
3617 unpack_buffer(pos) = recv_buffer(n)
3620 do i=1,xmap%size_put1
3621 x(i) = unpack_buffer(xmap%x1_put(i)%pos)
3628 call mpp_clock_end(id_put_1_to_xgrid_order_1)
3630 end subroutine put_1_to_xgrid_order_1
3635 subroutine put_1_to_xgrid_order_2(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3636 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
3637 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
3639 integer,
intent(in) :: isize, jsize, xsize, lsize
3642 real(r8_kind),
dimension(0:isize+1, 0:jsize+1, lsize) :: tmp
3643 real(r8_kind),
dimension(isize, jsize, lsize) :: tmpx, tmpy
3644 real(r8_kind),
dimension(isize, jsize, lsize) :: d_bar_max, d_bar_min
3645 real(r8_kind),
dimension(isize, jsize, lsize) :: d_max, d_min
3646 real(r8_kind) :: d_bar
3647 integer :: i, is, ie, j, js, je, ii, jj
3648 integer :: p, l, isd, jsd
3649 type (
grid_type),
pointer,
save :: grid1 =>null()
3650 type (
comm_type),
pointer,
save :: comm =>null()
3651 integer :: buffer_pos, msgsize, from_pe, to_pe, pos, n
3652 integer :: ibegin, count, istart, iend
3653 real(r8_kind) :: recv_buffer(xmap%put1%recvsize*lsize*3)
3654 real(r8_kind) :: send_buffer(xmap%put1%sendsize*lsize*3)
3655 real(r8_kind) :: unpack_buffer(xmap%put1%recvsize*3)
3656 logical :: on_west_edge, on_east_edge, on_south_edge, on_north_edge
3657 real(r8_kind),
dimension(isize, jsize) :: d
3658 real(r8_kind),
dimension(xsize) :: x
3662 call mpp_clock_begin(id_put_1_to_xgrid_order_2)
3663 grid1 => xmap%grids(1)
3665 is = grid1%is_me; ie = grid1%ie_me
3666 js = grid1%js_me; je = grid1%je_me
3672 tmp(:,:,l) = large_number
3674 tmp(1:isize,1:jsize,l) = d(:,:)
3677 if(grid1%is_latlon)
then
3678 call mpp_update_domains(tmp,grid1%domain_with_halo)
3681 tmpy(:,:,l) =
grad_merid_latlon(tmp(:,:,l), grid1%lat, is, ie, js, je, isd, jsd)
3682 tmpx(:,:,l) =
grad_zonal_latlon(tmp(:,:,l), grid1%lon, grid1%lat, is, ie, js, je, isd, jsd)
3685 call mpp_update_domains(tmp,grid1%domain)
3686 on_west_edge = (is==1)
3687 on_east_edge = (ie==grid1%im)
3688 on_south_edge = (js==1)
3689 on_north_edge = (je==grid1%jm)
3693 call gradient_cubic(tmp(:,:,l), grid1%box%dx, grid1%box%dy, grid1%box%area, &
3694 grid1%box%edge_w, grid1%box%edge_e, grid1%box%edge_s, &
3695 grid1%box%edge_n, grid1%box%en1, grid1%box%en2, &
3696 grid1%box%vlon, grid1%box%vlat, tmpx(:,:,l), tmpy(:,:,l), &
3697 on_west_edge, on_east_edge, on_south_edge, on_north_edge)
3704 do p = 1, comm%nrecv
3705 msgsize = comm%recv(p)%count*lsize
3706 buffer_pos = comm%recv(p)%buffer_pos*lsize
3707 if(.NOT. monotonic_exchange)
then
3709 buffer_pos = buffer_pos*3
3711 from_pe = comm%recv(p)%pe
3712 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = from_pe, block=.false., tag=comm_tag_8)
3716 if(monotonic_exchange)
then
3721 d_bar_max(i,j,l) = -large_number
3722 d_bar_min(i,j,l) = large_number
3723 d_max(i,j,l) = -large_number
3724 d_min(i,j,l) = large_number
3727 if(tmp(i,j,l) .NE. large_number)
then
3728 if(tmp(i,j,l) > d_bar_max(i,j,l)) d_bar_max(i,j,l) = tmp(i,j,l)
3729 if(tmp(i,j,l) < d_bar_min(i,j,l)) d_bar_min(i,j,l) = tmp(i,j,l)
3740 if(monotonic_exchange)
then
3742 do p = 1, comm%nsend
3743 msgsize = comm%send(p)%count*lsize
3744 to_pe = comm%send(p)%pe
3747 do n = 1, comm%send(p)%count
3749 i = comm%send(p)%i(n)
3750 j = comm%send(p)%j(n)
3751 send_buffer(pos) = d(i,j) + tmpy(i,j,l)*comm%send(p)%dj(n) + tmpx(i,j,l)*comm%send(p)%di(n)
3752 if(send_buffer(pos) > d_max(i,j,l)) d_max(i,j,l) = send_buffer(pos)
3753 if(send_buffer(pos) < d_min(i,j,l)) d_min(i,j,l) = send_buffer(pos)
3758 do p = 1, comm%nsend
3759 msgsize = comm%send(p)%count*lsize
3760 to_pe = comm%send(p)%pe
3764 do n = 1, comm%send(p)%count
3766 i = comm%send(p)%i(n)
3767 j = comm%send(p)%j(n)
3769 if( d_max(i,j,l) > d_bar_max(i,j,l) )
then
3770 send_buffer(pos) = d_bar + ((send_buffer(pos)-d_bar)/(d_max(i,j,l)-d_bar)) * (d_bar_max(i,j,l)-d_bar)
3771 else if( d_min(i,j,l) < d_bar_min(i,j,l) )
then
3772 send_buffer(pos) = d_bar + ((send_buffer(pos)-d_bar)/(d_min(i,j,l)-d_bar)) * (d_bar_min(i,j,l)-d_bar)
3776 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = to_pe, tag=comm_tag_8 )
3777 buffer_pos = buffer_pos + msgsize
3780 do p = 1, comm%nsend
3781 msgsize = comm%send(p)%count*lsize*3
3782 to_pe = comm%send(p)%pe
3786 do n = 1, comm%send(p)%count
3788 i = comm%send(p)%i(n)
3789 j = comm%send(p)%j(n)
3790 send_buffer(pos-2) = d(i,j)
3791 send_buffer(pos-1) = tmpy(i,j,l)
3792 send_buffer(pos ) = tmpx(i,j,l)
3795 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = to_pe, tag=comm_tag_8 )
3796 buffer_pos = buffer_pos + msgsize
3803 if(monotonic_exchange)
then
3804 if( lsize == 1)
then
3806 do l=1,xmap%size_put1
3807 pos = xmap%x1_put(l)%pos
3808 x(l) = recv_buffer(pos)
3814 do p = 1, comm%nsend
3815 count = comm%send(p)%count
3816 ibegin = comm%recv(p)%buffer_pos*lsize + 1
3817 istart = ibegin + (l-1)*count
3818 iend = istart + count - 1
3819 pos = comm%recv(p)%buffer_pos
3822 unpack_buffer(pos) = recv_buffer(n)
3825 do i=1,xmap%size_put1
3826 pos = xmap%x1_put(i)%pos
3827 x(i) = unpack_buffer(pos)
3832 if( lsize == 1)
then
3835 do l=1,xmap%size_put1
3836 pos = xmap%x1_put(l)%pos
3837 x(l) = recv_buffer(3*pos-2) + recv_buffer(3*pos-1)*xmap%x1_put(l)%dj + recv_buffer(3*pos)*xmap%x1_put(l)%di
3846 do p = 1, comm%nrecv
3847 count = comm%recv(p)%count*3
3848 ibegin = comm%recv(p)%buffer_pos*lsize*3 + 1
3849 istart = ibegin + (l-1)*count
3850 iend = istart + count - 1
3851 pos = comm%recv(p)%buffer_pos*3
3854 unpack_buffer(pos) = recv_buffer(n)
3857 do i=1,xmap%size_put1
3858 pos = xmap%x1_put(i)%pos
3859 x(i) = unpack_buffer(3*pos-2) + unpack_buffer(3*pos-1)*xmap%x1_put(i)%dj + unpack_buffer(3*pos) &
3860 & * xmap%x1_put(i)%di
3867 call mpp_clock_end(id_put_1_to_xgrid_order_2)
3869 end subroutine put_1_to_xgrid_order_2
3873 subroutine get_1_from_xgrid(d_addrs, x_addrs, xmap, isize, jsize, xsize, lsize)
3874 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
3875 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
3877 integer,
intent(in) :: isize, jsize, xsize, lsize
3879 real(r8_kind),
dimension(xmap%size),
target :: dg(xmap%size, lsize)
3880 integer :: i, j, l, p, n, m
3881 integer :: msgsize, buffer_pos, pos
3882 integer :: istart, iend, count
3883 real(r8_kind) ,
pointer,
save :: dgp =>null()
3884 type(
grid_type) ,
pointer,
save :: grid1 =>null()
3885 type(
comm_type) ,
pointer,
save :: comm =>null()
3888 real(r8_kind) :: recv_buffer(xmap%get1%recvsize*lsize*3)
3889 real(r8_kind) :: send_buffer(xmap%get1%sendsize*lsize*3)
3890 real(r8_kind) :: d(isize,jsize)
3891 real(r8_kind),
dimension(xsize) :: x
3895 call mpp_clock_begin(id_get_1_from_xgrid)
3898 grid1 => xmap%grids(1)
3900 do p = 1, comm%nrecv
3901 recv => comm%recv(p)
3902 msgsize = recv%count*lsize
3903 buffer_pos = recv%buffer_pos*lsize
3904 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = recv%pe, block=.false., tag=comm_tag_9)
3912 dgp => dg(xmap%x1(i)%pos,l)
3913 dgp = dgp + xmap%x1(i)%area*x(i)
3921 do p = 1, comm%nsend
3922 send => comm%send(p)
3923 msgsize = send%count*lsize
3925 istart = send%buffer_pos+1
3926 iend = istart + send%count - 1
3930 send_buffer(pos) = dg(n,l)
3933 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = send%pe, tag=comm_tag_9 )
3934 buffer_pos = buffer_pos + msgsize
3947 do p = 1, comm%nrecv
3948 recv => comm%recv(p)
3950 buffer_pos = recv%buffer_pos*lsize
3951 if( recv%pe == xmap%me )
then
3955 pos = buffer_pos + (l-1)*count
3961 d(i,j) = recv_buffer(pos)
3969 do m = 1, comm%nrecv
3970 p = comm%unpack_ind(m)
3971 recv => comm%recv(p)
3972 if( recv%pe == xmap%me )
then
3975 buffer_pos = recv%buffer_pos*lsize
3979 pos = buffer_pos + (l-1)*recv%count
3981 do n = 1, recv%count
3985 d(i,j) = d(i,j) + recv_buffer(pos)
3996 d = d * grid1%area_inv
3999 call mpp_clock_end(id_get_1_from_xgrid)
4001 end subroutine get_1_from_xgrid
4005 subroutine get_1_from_xgrid_repro(d_addrs, x_addrs, xmap, xsize, lsize)
4006 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
4007 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
4009 integer,
intent(in) :: xsize, lsize
4011 integer :: g, i, j, k, p, l, n, l2, l3
4012 integer :: msgsize, buffer_pos, pos
4013 type (
grid_type),
pointer,
save :: grid =>null()
4014 type(
comm_type),
pointer,
save :: comm => null()
4017 integer,
dimension(0:xmap%npes-1) :: pl, ml
4018 real(r8_kind) :: recv_buffer(xmap%recv_count_repro_tot*lsize)
4019 real(r8_kind) :: send_buffer(xmap%send_count_repro_tot*lsize)
4020 real(r8_kind) :: d(xmap%grids(1)%is_me:xmap%grids(1)%ie_me, &
4021 xmap%grids(1)%js_me:xmap%grids(1)%je_me)
4022 real(r8_kind),
dimension(xsize) :: x
4026 call mpp_clock_begin(id_get_1_from_xgrid_repro)
4027 comm => xmap%get1_repro
4029 do p = 1, comm%nrecv
4030 recv => comm%recv(p)
4031 msgsize = recv%count*lsize
4032 buffer_pos = recv%buffer_pos*lsize
4033 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = recv%pe, block=.false., tag=comm_tag_10)
4034 n = recv%pe -xmap%root_pe
4040 send_buffer(:) = 0.0_r8_kind
4043 do p = 1, comm%nsend
4044 pos = comm%send(p)%buffer_pos*lsize
4045 send => comm%send(p)
4048 do n = 1, send%count
4054 do k =1, xmap%grids(g)%km
4055 if(xmap%grids(g)%frac_area(i,j,k)/=0.0_r8_kind)
then
4057 send_buffer(pos) = send_buffer(pos) + xmap%x1(l2)%area *x(l2)
4065 buffer_pos = comm%send(p)%buffer_pos*lsize
4066 msgsize = comm%send(p)%count*lsize
4067 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe=comm%send(p)%pe, tag=comm_tag_10)
4081 do g=2,
size(xmap%grids(:))
4082 grid => xmap%grids(g)
4083 do l3=1,grid%size_repro
4084 i = grid%x_repro(l3)%i1
4085 j = grid%x_repro(l3)%j1
4086 p = grid%x_repro(l3)%pe-xmap%root_pe
4087 pos = pl(p) + (l-1)*ml(p) + grid%x_repro(l3)%recv_pos
4088 d(i,j) = d(i,j) + recv_buffer(pos)
4092 d = d * xmap%grids(1)%area_inv
4097 call mpp_clock_end(id_get_1_from_xgrid_repro)
4099 end subroutine get_1_from_xgrid_repro
4108 real(r8_kind),
dimension(:,:),
intent(in) :: d
4109 character(len=3),
intent(in) :: grid_id
4112 integer,
intent(in),
optional :: remap_method
4115 real(r8_kind),
dimension(xmap%size) :: x_over, x_back
4116 real(r8_kind),
dimension(size(d,1),size(d,2)) :: d1
4117 real(r8_kind),
dimension(:,:,:),
allocatable :: d2
4119 type (
grid_type),
pointer,
save :: grid1 =>null(), grid2 =>null()
4121 grid1 => xmap%grids(1)
4127 call put_to_xgrid (d, grid1%id, x_over, xmap, remap_method)
4128 do g=2,
size(xmap%grids(:))
4129 grid2 => xmap%grids(g)
4130 if(grid2%on_this_pe)
then
4131 allocate (d2(grid2%is_me:grid2%ie_me, grid2%js_me:grid2%je_me, grid2%km) )
4134 if(grid2%on_this_pe)
then
4138 if(
allocated(d2))
deallocate (d2)
4155 real(r8_kind),
dimension(:,:,:),
intent(in) :: d
4156 character(len=3),
intent(in) :: grid_id
4159 integer,
intent(in),
optional :: remap_method
4162 real(r8_kind),
dimension(xmap%size) :: x_over, x_back
4163 real(r8_kind),
dimension(:,: ),
allocatable :: d1
4164 real(r8_kind),
dimension(:,:,:),
allocatable :: d2
4166 type (
grid_type),
pointer,
save :: grid1 =>null(), grid2 =>null()
4168 grid1 => xmap%grids(1)
4170 do g = 2,
size(xmap%grids(:))
4171 grid2 => xmap%grids(g)
4172 if (grid_id==grid2%id)
then
4173 if(grid2%on_this_pe)
then
4178 call put_to_xgrid(0.0_r8_kind * grid2%frac_area, grid2%id, x_over, xmap)
4182 allocate ( d1(
size(grid1%area,1),
size(grid1%area,2)) )
4185 call put_to_xgrid(d1, grid1%id, x_back, xmap,remap_method)
4189 do g = 2,
size(xmap%grids(:))
4190 grid2 => xmap%grids(g)
4191 if(grid2%on_this_pe)
then
4192 allocate ( d2(
size(grid2%frac_area, 1),
size(grid2%frac_area, 2), &
4193 size(grid2%frac_area, 3) ) )
4197 if(
allocated(d2) )
deallocate ( d2 )
4211 real(r8_kind),
dimension(:,:),
intent(in) :: d
4212 character(len=3),
intent(in) :: grid_id
4215 integer,
intent(in),
optional :: remap_method
4217 real(r8_kind),
dimension(xmap%size) :: x_over, x_back
4218 real(r8_kind),
dimension(size(d,1),size(d,2)) :: d1
4219 real(r8_kind),
dimension(:,:,:),
allocatable :: d2
4220 real(r8_kind),
dimension(: ),
allocatable :: d_ug
4221 real(r8_kind),
dimension(:,:),
allocatable :: d2_ug
4223 type (
grid_type),
pointer,
save :: grid1 =>null(), grid2 =>null()
4225 grid1 => xmap%grids(1)
4229 if(grid1%is_ug)
then
4230 allocate(d_ug(grid1%ls_me:grid1%le_me))
4231 call mpp_pass_sg_to_ug(grid1%ug_domain, d, d_ug)
4236 call put_to_xgrid (d, grid1%id, x_over, xmap, remap_method)
4238 do g=2,
size(xmap%grids(:))
4239 grid2 => xmap%grids(g)
4240 if(grid2%is_ug)
then
4241 if(grid2%on_this_pe)
then
4242 allocate (d2_ug(grid2%ls_me:grid2%le_me, grid2%km) )
4246 if(grid2%on_this_pe)
then
4248 sum( grid2%area(:,1) * sum(grid2%frac_area(:,1,:)*d2_ug,dim=2) )
4251 if(
allocated(d2_ug))
deallocate (d2_ug)
4253 if(grid2%on_this_pe)
then
4254 allocate (d2(grid2%is_me:grid2%ie_me, grid2%js_me:grid2%je_me, grid2%km) )
4257 if(grid2%on_this_pe)
then
4259 & + sum( grid2%area * sum(grid2%frac_area*d2,dim=3) )
4262 if(
allocated(d2))
deallocate (d2)
4265 if(grid1%is_ug)
then
4272 if(
allocated(d_ug))
deallocate(d_ug)
4284 real(r8_kind),
dimension(:,:,:),
intent(in) :: d
4285 character(len=3),
intent(in) :: grid_id
4288 integer,
intent(in),
optional :: remap_method
4291 real(r8_kind),
dimension(xmap%size) :: x_over, x_back
4292 real(r8_kind),
dimension(:,: ),
allocatable :: d1, d_ug
4293 real(r8_kind),
dimension(:,:,:),
allocatable :: d2
4295 type (
grid_type),
pointer,
save :: grid1 =>null(), grid2 =>null()
4297 grid1 => xmap%grids(1)
4299 do g = 2,
size(xmap%grids(:))
4300 grid2 => xmap%grids(g)
4301 if (grid_id==grid2%id)
then
4302 if(grid2%on_this_pe)
then
4303 if(grid2%is_ug)
then
4304 allocate(d_ug(grid2%ls_me:grid2%le_me,grid2%km))
4305 call mpp_pass_sg_to_ug(grid2%ug_domain, d, d_ug)
4311 if(grid2%is_ug)
then
4316 if(
allocated(d_ug))
deallocate(d_ug)
4318 if(grid2%is_ug)
then
4319 call put_to_xgrid_ug(0.0_r8_kind * grid2%frac_area(:,1,:), grid2%id, x_over, xmap)
4321 call put_to_xgrid(0.0_r8_kind * grid2%frac_area, grid2%id, x_over, xmap)
4326 allocate ( d1(
size(grid1%area,1),
size(grid1%area,2)) )
4327 if(grid1%is_ug)
then
4333 if(grid1%is_ug)
then
4336 call put_to_xgrid(d1, grid1%id, x_back, xmap,remap_method)
4341 do g = 2,
size(xmap%grids(:))
4342 grid2 => xmap%grids(g)
4343 if(grid2%on_this_pe)
then
4344 allocate ( d2(
size(grid2%frac_area, 1),
size(grid2%frac_area, 2), &
4345 size(grid2%frac_area, 3) ) )
4347 if(grid2%is_ug)
then
4353 if(
allocated(d2) )
deallocate ( d2 )
4364 character(len=3),
intent(in) :: id
4366 real(r8_kind),
dimension(:,:),
intent(out) :: area
4371 do g = 1,
size(xmap%grids(:))
4372 if (id==xmap%grids(g)%id )
then
4373 if(
size(area,1) .NE.
size(xmap%grids(g)%area,1) .OR.
size(area,2) .NE.
size(xmap%grids(g)%area,2) ) &
4374 call error_mesg(
"xgrid_mod",
"size mismatch between area and xmap%grids(g)%area", fatal)
4375 area = xmap%grids(g)%area
4381 if(.not. found)
call error_mesg(
"xgrid_mod", id//
" is not found in xmap%grids id", fatal)
4393 integer,
intent(in) :: isd, jsd
4394 real(r8_kind),
dimension(isd:,jsd:),
intent(in) :: d
4395 real(r8_kind),
dimension(:),
intent(in) :: lon
4396 real(r8_kind),
dimension(:),
intent(in) :: lat
4397 integer,
intent(in) :: is, ie, js, je
4399 real(r8_kind) :: dx, costheta
4400 integer :: i, j, ip1, im1
4406 else if(i==
size(lon(:)) )
then
4409 ip1 = i+1; im1 = i-1
4411 dx = lon(ip1) - lon(im1)
4412 if(abs(dx).lt.eps )
call error_mesg(
'xgrids_mod(grad_zonal_latlon)',
'Improper grid size in lontitude', fatal)
4413 if(dx .gt. pi) dx = dx - 2.0_r8_kind* pi
4414 if(dx .lt. -pi) dx = dx + 2.0_r8_kind* pi
4416 costheta = cos(lat(j))
4417 if(abs(costheta) .lt. eps)
call error_mesg(
'xgrids_mod(grad_zonal_latlon)',
'Improper latitude grid', fatal)
4432 integer,
intent(in) :: isd, jsd
4433 real(r8_kind),
dimension(isd:,jsd:),
intent(in) :: d
4434 real(r8_kind),
dimension(:),
intent(in) :: lat
4435 integer,
intent(in) :: is, ie, js, je
4438 integer :: i, j, jp1, jm1
4444 else if(j ==
size(lat(:)) )
then
4447 jp1 = j+1; jm1 = j-1
4449 dy = lat(jp1) - lat(jm1)
4450 if(abs(dy).lt.eps)
call error_mesg(
'xgrids_mod(grad_merid_latlon)',
'Improper grid size in latitude', fatal)
4461 subroutine get_index_range(xmap, grid_index, is, ie, js, je, km)
4464 integer,
intent(in) :: grid_index
4465 integer,
intent(out) :: is, ie, js, je, km
4467 is = xmap % grids(grid_index) % is_me
4468 ie = xmap % grids(grid_index) % ie_me
4469 js = xmap % grids(grid_index) % js_me
4470 je = xmap % grids(grid_index) % je_me
4471 km = xmap % grids(grid_index) % km
4473 end subroutine get_index_range
4480 subroutine stock_move_3d(from, to, grid_index, stock_data3d, xmap, &
4481 & delta_t, from_side, to_side, radius, verbose, ier)
4491 type(stock_type),
intent(inout),
optional :: from, to
4492 integer,
intent(in) :: grid_index
4493 real(r8_kind),
intent(in) :: stock_data3d(:,:,:)
4495 real(r8_kind),
intent(in) :: delta_t
4496 integer,
intent(in) :: from_side
4497 integer,
intent(in) :: to_side
4498 real(r8_kind),
intent(in) :: radius
4499 character(len=*),
intent(in),
optional :: verbose
4500 integer,
intent(out) :: ier
4502 real(r8_kind) :: from_dq, to_dq
4505 if(grid_index == 1)
then
4511 if(.not.
associated(xmap%grids) )
then
4516 from_dq = delta_t * 4.0_r8_kind * pi * radius**2 * sum( sum(xmap%grids(grid_index)%area * &
4517 & sum(xmap%grids(grid_index)%frac_area * stock_data3d, dim=3), dim=1))
4521 if(
present(to )) to % dq( to_side) = to % dq( to_side) + to_dq
4522 if(
present(from)) from % dq(from_side) = from % dq(from_side) - from_dq
4524 if(
present(verbose).and.debug_stocks)
then
4527 from_dq = from_dq/(4.0_r8_kind*pi*radius**2)
4528 to_dq = to_dq /(4.0_r8_kind*pi*radius**2)
4529 if(
mpp_pe()==mpp_root_pe())
then
4530 write(stocks_file,
'(a,es19.12,a,es19.12,a)') verbose, from_dq,
' [*/m^2]'
4540 subroutine stock_move_2d(from, to, grid_index, stock_data2d, xmap, &
4541 & delta_t, from_side, to_side, radius, verbose, ier)
4550 type(stock_type),
intent(inout),
optional :: from, to
4551 integer,
optional,
intent(in) :: grid_index
4552 real(r8_kind),
intent(in) :: stock_data2d(:,:)
4554 real(r8_kind),
intent(in) :: delta_t
4555 integer,
intent(in) :: from_side
4556 integer,
intent(in) :: to_side
4557 real(r8_kind),
intent(in) :: radius
4558 character(len=*),
intent(in) :: verbose
4559 integer,
intent(out) :: ier
4561 real(r8_kind) :: to_dq, from_dq
4565 if(.not.
associated(xmap%grids) )
then
4570 if( .not.
present(grid_index) .or. grid_index==1 )
then
4573 from_dq = delta_t * 4.0_r8_kind*pi*radius**2 * sum(sum(xmap%grids(1)%area * stock_data2d, dim=1))
4584 if(
present(to )) to % dq( to_side) = to % dq( to_side) + to_dq
4585 if(
present(from)) from % dq(from_side) = from % dq(from_side) - from_dq
4587 if(debug_stocks)
then
4590 from_dq = from_dq/(4.0_r8_kind*pi*radius**2)
4591 to_dq = to_dq /(4.0_r8_kind*pi*radius**2)
4592 if(
mpp_pe()==mpp_root_pe())
then
4593 write(stocks_file,
'(a,es19.12,a,es19.12,a)') verbose, from_dq,
' [*/m^2]'
4597 end subroutine stock_move_2d
4604 subroutine stock_move_ug_3d(from, to, grid_index, stock_ug_data3d, xmap, &
4605 & delta_t, from_side, to_side, radius, verbose, ier)
4615 type(stock_type),
intent(inout),
optional :: from, to
4616 integer,
intent(in) :: grid_index
4617 real(r8_kind),
intent(in) :: stock_ug_data3d(:,:)
4619 real(r8_kind),
intent(in) :: delta_t
4620 integer,
intent(in) :: from_side
4621 integer,
intent(in) :: to_side
4622 real(r8_kind),
intent(in) :: radius
4623 character(len=*),
intent(in),
optional :: verbose
4624 integer,
intent(out) :: ier
4625 real(r8_kind),
dimension(size(stock_ug_data3d,1),size(stock_ug_data3d,2)) :: tmp
4627 real(r8_kind) :: from_dq, to_dq
4630 if(grid_index == 1)
then
4636 if(.not.
associated(xmap%grids) )
then
4641 tmp = xmap%grids(grid_index)%frac_area(:,1,:) * stock_ug_data3d
4642 from_dq = delta_t * 4.0_r8_kind * pi * radius**2 * sum( xmap%grids(grid_index)%area(:,1) * &
4647 if(
present(to )) to % dq( to_side) = to % dq( to_side) + to_dq
4648 if(
present(from)) from % dq(from_side) = from % dq(from_side) - from_dq
4650 if(
present(verbose).and.debug_stocks)
then
4653 from_dq = from_dq/(4.0_r8_kind*pi*radius**2)
4654 to_dq = to_dq /(4.0_r8_kind*pi*radius**2)
4655 if(
mpp_pe()==mpp_root_pe())
then
4656 write(stocks_file,
'(a,es19.12,a,es19.12,a)') verbose, from_dq,
' [*/m^2]'
4660 end subroutine stock_move_ug_3d
4666 subroutine stock_integrate_2d(integrate_data2d, xmap, delta_t, radius, res, ier)
4672 real(r8_kind),
intent(in) :: integrate_data2d(:,:)
4674 real(r8_kind),
intent(in) :: delta_t
4675 real(r8_kind),
intent(in) :: radius
4676 real(r8_kind),
intent(out) :: res
4677 integer,
intent(out) :: ier
4682 if(.not.
associated(xmap%grids) )
then
4687 res = delta_t * 4.0_r8_kind * pi * radius**2 * sum(sum(xmap%grids(1)%area * integrate_data2d, dim=1))
4689 end subroutine stock_integrate_2d
4696 subroutine stock_print(stck, Time, comp_name, index, ref_value, radius, pelist)
4702 type(stock_type),
intent(in) :: stck
4704 character(len=*) :: comp_name
4705 integer,
intent(in) :: index
4706 real(r8_kind),
intent(in) :: ref_value
4707 real(r8_kind),
intent(in) :: radius
4708 integer,
intent(in),
optional :: pelist(:)
4710 integer,
parameter :: initID = -2
4714 real(r8_kind) :: f_value, c_value, planet_area
4715 character(len=80) :: formatString
4716 integer :: iday, isec, hours
4717 integer :: diagID, compInd
4718 integer,
dimension(NELEMS,4),
save :: f_valueDiagID = initid
4719 integer,
dimension(NELEMS,4),
save :: c_valueDiagID = initid
4720 integer,
dimension(NELEMS,4),
save :: fmc_valueDiagID = initid
4722 real(r8_kind) :: diagField
4724 character(len=30) :: field_name, units
4726 f_value = sum(stck % dq)
4727 c_value = ref_value - stck % q_start
4728 if(
present(pelist))
then
4729 call mpp_sum(f_value, pelist=pelist)
4730 call mpp_sum(c_value, pelist=pelist)
4736 if(
mpp_pe() == mpp_root_pe())
then
4738 planet_area = 4.0_r8_kind * pi * radius**2
4739 f_value = f_value / planet_area
4740 c_value = c_value / planet_area
4742 if(comp_name ==
'ATM') compind = 1
4743 if(comp_name ==
'LND') compind = 2
4744 if(comp_name ==
'ICE') compind = 3
4745 if(comp_name ==
'OCN') compind = 4
4748 if(f_valuediagid(index,compind) == initid)
then
4749 field_name = trim(comp_name) // trim(stock_names(index))
4750 field_name = trim(field_name) //
'StocksChange_Flux'
4751 units = trim(stock_units(index))
4756 if(c_valuediagid(index,compind) == initid)
then
4757 field_name = trim(comp_name) // trim(stock_names(index))
4758 field_name = trim(field_name) //
'StocksChange_Comp'
4759 units = trim(stock_units(index))
4764 if(fmc_valuediagid(index,compind) == initid)
then
4765 field_name = trim(comp_name) // trim(stock_names(index))
4766 field_name = trim(field_name) //
'StocksChange_Diff'
4767 units = trim(stock_units(index))
4772 diagid=f_valuediagid(index,compind)
4774 if (diagid > 0) used =
send_data(diagid, diagfield, time = time)
4775 diagid=c_valuediagid(index,compind)
4777 if (diagid > 0) used =
send_data(diagid, diagfield, time)
4778 diagid=fmc_valuediagid(index,compind)
4779 diagfield = f_value-c_value
4780 if (diagid > 0) used =
send_data(diagid, diagfield, time=time)
4784 hours = iday*24 + isec/3600
4785 formatstring =
'(a,a,a,i16,2x,es22.15,2x,es22.15,2x,es22.15)'
4786 write(stocks_file,formatstring) trim(comp_name),stock_names(index),stock_units(index) &
4787 ,hours,f_value,c_value,f_value-c_value
4792 end subroutine stock_print
4797 function is_lat_lon(lon, lat)
4798 real(r8_kind),
dimension(:,:),
intent(in) :: lon, lat
4799 logical :: is_lat_lon
4800 integer :: i, j, nlon, nlat, num
4805 loop_lat:
do j = 1, nlat
4807 if(lat(i,j) .NE. lat(1,j))
then
4808 is_lat_lon = .false.
4815 loop_lon:
do i = 1, nlon
4817 if(lon(i,j) .NE. lon(i,1))
then
4818 is_lat_lon = .false.
4826 if(is_lat_lon) num = 1
4831 is_lat_lon = .false.
4835 end function is_lat_lon
4845 subroutine get_side1_from_xgrid_ug(d, grid_id, x, xmap, complete)
4846 real(r8_kind),
dimension(:),
intent(out) :: d
4847 character(len=3),
intent(in) :: grid_id
4848 real(r8_kind),
dimension(:),
intent(in) :: x
4849 type (xmap_type),
intent(inout) :: xmap
4850 logical,
intent(in),
optional :: complete
4852 logical :: is_complete, set_mismatch
4854 character(len=2) :: text
4855 integer,
save :: isize=0
4856 integer,
save :: lsize=0
4857 integer,
save :: xsize=0
4858 character(len=3),
save :: grid_id_saved=
""
4859 integer(i8_kind),
dimension(MAX_FIELDS),
save :: d_addrs = -9999_i8_kind
4860 integer(i8_kind),
dimension(MAX_FIELDS),
save :: x_addrs = -9999_i8_kind
4863 if (grid_id==xmap%grids(1)%id)
then
4864 is_complete = .true.
4865 if(
present(complete)) is_complete=complete
4867 if( lsize > max_fields )
then
4868 write( text,
'(i2)' ) max_fields
4869 call error_mesg (
'xgrid_mod',
'MAX_FIELDS='//trim(text)//
' exceeded for group get_side1_from_xgrid_ug', fatal)
4871 d_addrs(lsize) = loc(d)
4872 x_addrs(lsize) = loc(x)
4877 grid_id_saved = grid_id
4879 set_mismatch = .false.
4880 set_mismatch = set_mismatch .OR. (isize /=
size(d(:)))
4881 set_mismatch = set_mismatch .OR. (xsize /=
size(x(:)))
4882 set_mismatch = set_mismatch .OR. (grid_id_saved /= grid_id)
4883 if(set_mismatch)
then
4884 write( text,
'(i2)' ) lsize
4885 call error_mesg (
'xgrid_mod',
'Incompatible field at count '//text// &
4886 &
' for group get_side1_from_xgrid_ug', fatal )
4890 if(is_complete)
then
4892 call get_1_from_xgrid_ug_repro(d_addrs, x_addrs, xmap, xsize, lsize)
4894 call get_1_from_xgrid_ug(d_addrs, x_addrs, xmap, isize, xsize, lsize)
4896 d_addrs(1:lsize) = -9999_i8_kind
4897 x_addrs(1:lsize) = -9999_i8_kind
4906 do g=2,
size(xmap%grids(:))
4907 if (grid_id==xmap%grids(g)%id) &
4908 call error_mesg (
'xgrid_mod', &
4909 'get_from_xgrid_ug expects a 3D side 2 grid', fatal)
4912 call error_mesg (
'xgrid_mod',
'get_from_xgrid_ug: could not find grid id', fatal)
4914 end subroutine get_side1_from_xgrid_ug
4928 real(r8_kind),
dimension(:),
intent(in) :: d
4929 character(len=3),
intent(in) :: grid_id
4930 real(r8_kind),
dimension(:),
intent(inout) :: x
4931 type (xmap_type),
intent(inout) :: xmap
4932 logical,
intent(in),
optional :: complete
4934 logical :: is_complete, set_mismatch
4936 character(len=2) :: text
4937 integer,
save :: dsize=0
4938 integer,
save :: lsize=0
4939 integer,
save :: xsize=0
4940 character(len=3),
save :: grid_id_saved=
""
4941 integer(i8_kind),
dimension(MAX_FIELDS),
save :: d_addrs = -9999_i8_kind
4942 integer(i8_kind),
dimension(MAX_FIELDS),
save :: x_addrs = -9999_i8_kind
4944 if (grid_id==xmap%grids(1)%id)
then
4945 is_complete = .true.
4946 if(
present(complete)) is_complete=complete
4948 if( lsize > max_fields )
then
4949 write( text,
'(i2)' ) max_fields
4950 call error_mesg (
'xgrid_mod',
'MAX_FIELDS='//trim(text)//
' exceeded for group put_side1_to_xgrid_ug', fatal)
4952 d_addrs(lsize) = loc(d)
4953 x_addrs(lsize) = loc(x)
4958 grid_id_saved = grid_id
4960 set_mismatch = .false.
4961 set_mismatch = set_mismatch .OR. (dsize /=
size(d(:)))
4962 set_mismatch = set_mismatch .OR. (xsize /=
size(x(:)))
4963 set_mismatch = set_mismatch .OR. (grid_id_saved /= grid_id)
4964 if(set_mismatch)
then
4965 write( text,
'(i2)' ) lsize
4966 call error_mesg (
'xgrid_mod',
'Incompatible field at count '//text// &
4967 &
' for group put_side1_to_xgrid_ug', fatal )
4971 if(is_complete)
then
4972 call put_1_to_xgrid_ug_order_1(d_addrs, x_addrs, xmap, dsize, xsize, lsize)
4973 d_addrs(1:lsize) = -9999_i8_kind
4974 x_addrs(1:lsize) = -9999_i8_kind
4983 do g=2,
size(xmap%grids(:))
4984 if (grid_id==xmap%grids(g)%id) &
4985 call error_mesg (
'xgrid_mod', &
4986 'put_to_xgrid_ug expects a 2D side 2 grid', fatal)
4989 call error_mesg (
'xgrid_mod',
'put_to_xgrid_ug: could not find grid id', fatal)
5002 subroutine put_side2_to_xgrid_ug(d, grid_id, x, xmap)
5003 real(r8_kind),
dimension(:,:),
intent(in) :: d
5004 character(len=3),
intent(in) :: grid_id
5005 real(r8_kind),
dimension(:),
intent(inout) :: x
5006 type (xmap_type),
intent(inout) :: xmap
5010 if (grid_id==xmap%grids(1)%id) &
5011 call error_mesg (
'xgrid_mod', &
5012 'put_to_xgrid_ug expects a 2D side 1 grid', fatal)
5014 do g=2,
size(xmap%grids(:))
5015 if (grid_id==xmap%grids(g)%id)
then
5016 call put_2_to_xgrid_ug(d, xmap%grids(g), x, xmap)
5021 call error_mesg (
'xgrid_mod',
'put_to_xgrid_ug: could not find grid id', fatal)
5023 end subroutine put_side2_to_xgrid_ug
5034 subroutine get_side2_from_xgrid_ug(d, grid_id, x, xmap)
5035 real(r8_kind),
dimension(:,:),
intent(out) :: d
5036 character(len=3),
intent(in) :: grid_id
5037 real(r8_kind),
dimension(:),
intent(in) :: x
5038 type (xmap_type),
intent(in) :: xmap
5042 if (grid_id==xmap%grids(1)%id) &
5043 call error_mesg (
'xgrid_mod', &
5044 'get_from_xgrid_ug expects a 2D side 1 grid', fatal)
5046 do g=2,
size(xmap%grids(:))
5047 if (grid_id==xmap%grids(g)%id)
then
5048 call get_2_from_xgrid_ug(d, xmap%grids(g), x, xmap)
5053 call error_mesg (
'xgrid_mod',
'get_from_xgrid_ug: could not find grid id', fatal)
5055 end subroutine get_side2_from_xgrid_ug
5061 subroutine put_1_to_xgrid_ug_order_1(d_addrs, x_addrs, xmap, dsize, xsize, lsize)
5062 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
5063 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
5064 type (xmap_type),
intent(inout) :: xmap
5065 integer,
intent(in) :: dsize, xsize, lsize
5067 integer :: i, p, buffer_pos, msgsize
5068 integer :: from_pe, to_pe, pos, n, l, count
5069 integer :: ibegin, istart, iend, start_pos
5070 type (comm_type),
pointer,
save :: comm =>null()
5071 real(r8_kind) :: recv_buffer(xmap%put1%recvsize*lsize)
5072 real(r8_kind) :: send_buffer(xmap%put1%sendsize*lsize)
5073 real(r8_kind) :: unpack_buffer(xmap%put1%recvsize)
5075 real(r8_kind),
dimension(dsize) :: d
5076 real(r8_kind),
dimension(xsize) :: x
5081 call mpp_clock_begin(id_put_1_to_xgrid_order_1)
5085 do p = 1, comm%nrecv
5086 msgsize = comm%recv(p)%count*lsize
5087 from_pe = comm%recv(p)%pe
5088 buffer_pos = comm%recv(p)%buffer_pos*lsize
5089 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = from_pe, block=.false., tag=comm_tag_7)
5094 do p = 1, comm%nsend
5095 msgsize = comm%send(p)%count*lsize
5096 to_pe = comm%send(p)%pe
5100 do n = 1, comm%send(p)%count
5102 lll = comm%send(p)%i(n)
5103 send_buffer(pos) = d(lll)
5106 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = to_pe, tag=comm_tag_7 )
5107 buffer_pos = buffer_pos + msgsize
5113 if( lsize == 1)
then
5115 do l=1,xmap%size_put1
5116 x(l) = recv_buffer(xmap%x1_put(l)%pos)
5124 do p = 1, comm%nrecv
5125 count = comm%recv(p)%count
5126 ibegin = comm%recv(p)%buffer_pos*lsize + 1
5127 istart = ibegin + (l-1)*count
5128 iend = istart + count - 1
5129 pos = comm%recv(p)%buffer_pos
5132 unpack_buffer(pos) = recv_buffer(n)
5135 do i=1,xmap%size_put1
5136 x(i) = unpack_buffer(xmap%x1_put(i)%pos)
5143 call mpp_clock_end(id_put_1_to_xgrid_order_1)
5145 end subroutine put_1_to_xgrid_ug_order_1
5149 subroutine put_2_to_xgrid_ug(d, grid, x, xmap)
5150 type (grid_type),
intent(in) :: grid
5151 real(r8_kind),
dimension(grid%ls_me:grid%le_me, grid%km),
intent(in) :: d
5152 real(r8_kind),
dimension(:),
intent(inout) :: x
5153 type (xmap_type),
intent(in) :: xmap
5156 call mpp_clock_begin(id_put_2_to_xgrid)
5158 do l=grid%first,grid%last
5159 x(l) = d(xmap%x2(l)%l,xmap%x2(l)%k)
5162 call mpp_clock_end(id_put_2_to_xgrid)
5163 end subroutine put_2_to_xgrid_ug
5166 subroutine get_1_from_xgrid_ug(d_addrs, x_addrs, xmap, isize, xsize, lsize)
5167 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
5168 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
5169 type (xmap_type),
intent(inout) :: xmap
5170 integer,
intent(in) :: isize, xsize, lsize
5172 real(r8_kind),
dimension(xmap%size),
target :: dg(xmap%size, lsize)
5173 integer :: i, j, l, p, n, m
5174 integer :: msgsize, buffer_pos, pos
5175 integer :: istart, iend, count
5176 real(r8_kind) ,
pointer,
save :: dgp =>null()
5177 type (grid_type) ,
pointer,
save :: grid1 =>null()
5178 type (comm_type) ,
pointer,
save :: comm =>null()
5181 real(r8_kind) :: recv_buffer(xmap%get1%recvsize*lsize*3)
5182 real(r8_kind) :: send_buffer(xmap%get1%sendsize*lsize*3)
5183 real(r8_kind) :: d(isize)
5184 real(r8_kind),
dimension(xsize) :: x
5188 call mpp_clock_begin(id_get_1_from_xgrid)
5191 grid1 => xmap%grids(1)
5193 do p = 1, comm%nrecv
5194 recv => comm%recv(p)
5195 msgsize = recv%count*lsize
5196 buffer_pos = recv%buffer_pos*lsize
5197 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = recv%pe, block=.false., tag=comm_tag_9)
5205 dgp => dg(xmap%x1(i)%pos,l)
5206 dgp = dgp + xmap%x1(i)%area*x(i)
5214 do p = 1, comm%nsend
5215 send => comm%send(p)
5216 msgsize = send%count*lsize
5218 istart = send%buffer_pos+1
5219 iend = istart + send%count - 1
5223 send_buffer(pos) = dg(n,l)
5226 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe = send%pe, tag=comm_tag_9 )
5227 buffer_pos = buffer_pos + msgsize
5240 do p = 1, comm%nrecv
5241 recv => comm%recv(p)
5243 buffer_pos = recv%buffer_pos*lsize
5244 if( recv%pe == xmap%me )
then
5248 pos = buffer_pos + (l-1)*count
5253 d(i) = recv_buffer(pos)
5261 do m = 1, comm%nrecv
5262 p = comm%unpack_ind(m)
5263 recv => comm%recv(p)
5264 if( recv%pe == xmap%me )
then
5267 buffer_pos = recv%buffer_pos*lsize
5271 pos = buffer_pos + (l-1)*recv%count
5273 do n = 1, recv%count
5276 d(i) = d(i) + recv_buffer(pos)
5287 d = d * grid1%area_inv(:,1)
5290 call mpp_clock_end(id_get_1_from_xgrid)
5292 end subroutine get_1_from_xgrid_ug
5296 subroutine get_1_from_xgrid_ug_repro(d_addrs, x_addrs, xmap, xsize, lsize)
5297 integer(i8_kind),
dimension(:),
intent(in) :: d_addrs
5298 integer(i8_kind),
dimension(:),
intent(in) :: x_addrs
5299 type (xmap_type),
intent(inout) :: xmap
5300 integer,
intent(in) :: xsize, lsize
5302 integer :: g, i, j, k, p, l, n, l2, l3
5303 integer :: msgsize, buffer_pos, pos
5304 type (grid_type),
pointer,
save :: grid =>null()
5305 type(
comm_type),
pointer,
save :: comm => null()
5308 integer,
dimension(0:xmap%npes-1) :: pl, ml
5309 real(r8_kind) :: recv_buffer(xmap%recv_count_repro_tot*lsize)
5310 real(r8_kind) :: send_buffer(xmap%send_count_repro_tot*lsize)
5311 real(r8_kind) :: d(xmap%grids(1)%ls_me:xmap%grids(1)%le_me)
5312 real(r8_kind),
dimension(xsize) :: x
5316 call mpp_clock_begin(id_get_1_from_xgrid_repro)
5317 comm => xmap%get1_repro
5319 do p = 1, comm%nrecv
5320 recv => comm%recv(p)
5321 msgsize = recv%count*lsize
5322 buffer_pos = recv%buffer_pos*lsize
5323 call mpp_recv(recv_buffer(buffer_pos+1), glen=msgsize, from_pe = recv%pe, block=.false., tag=comm_tag_10)
5324 n = recv%pe -xmap%root_pe
5330 send_buffer(:) = 0.0_r8_kind
5333 do p = 1, comm%nsend
5334 pos = comm%send(p)%buffer_pos*lsize
5335 send => comm%send(p)
5338 do n = 1, send%count
5344 do k =1, xmap%grids(g)%km
5345 if(xmap%grids(g)%frac_area(i,j,k)/=0.0_r8_kind)
then
5347 send_buffer(pos) = send_buffer(pos) + xmap%x1(l2)%area *x(l2)
5355 buffer_pos = comm%send(p)%buffer_pos*lsize
5356 msgsize = comm%send(p)%count*lsize
5357 call mpp_send(send_buffer(buffer_pos+1), plen=msgsize, to_pe=comm%send(p)%pe, tag=comm_tag_10)
5371 do g=2,
size(xmap%grids(:))
5372 grid => xmap%grids(g)
5373 do l3=1,grid%size_repro
5374 i = grid%x_repro(l3)%l1
5375 p = grid%x_repro(l3)%pe-xmap%root_pe
5376 pos = pl(p) + (l-1)*ml(p) + grid%x_repro(l3)%recv_pos
5377 d(i) = d(i) + recv_buffer(pos)
5381 d = d * xmap%grids(1)%area_inv(:,1)
5386 call mpp_clock_end(id_get_1_from_xgrid_repro)
5388 end subroutine get_1_from_xgrid_ug_repro
5392 subroutine get_2_from_xgrid_ug(d, grid, x, xmap)
5393 type (grid_type),
intent(in) :: grid
5394 real(r8_kind),
dimension(grid%ls_me:grid%le_me, grid%km),
intent(out) :: d
5395 real(r8_kind),
dimension(:),
intent(in) :: x
5396 type (xmap_type),
intent(in) :: xmap
5400 call mpp_clock_begin(id_get_2_from_xgrid)
5403 do l=grid%first_get,grid%last_get
5404 d(xmap%x2_get(l)%l,xmap%x2_get(l)%k) = &
5405 d(xmap%x2_get(l)%l,xmap%x2_get(l)%k) + xmap%x2_get(l)%area*x(xmap%x2_get(l)%pos)
5411 d(:,k) = d(:,k) * grid%area_inv(:,1)
5414 call mpp_clock_end(id_get_2_from_xgrid)
5416 end subroutine get_2_from_xgrid_ug
5421 integer,
intent(in) :: i, j
5426 g = (j-1)*grid%ni + i
5427 in_box_me = (g>=grid%gs_me) .and. (g<=grid%ge_me)
5429 in_box_me = (i>=grid%is_me) .and. (i<=grid%ie_me) .and. (j>=grid%js_me) .and. (j<=grid%je_me)
5437 integer,
intent(in) :: i, j, p
5442 g = (j-1)*grid%ni + i
5443 in_box_nbr = (g>=grid%gs(p)) .and. (g<=grid%ge(p))
5445 in_box_nbr = (i>=grid%is(p)) .and. (i<=grid%ie(p)) .and. (j>=grid%js(p)) .and. (j<=grid%je(p))
5450 end module xgrid_mod
5453 Register a diagnostic field for a given module.
Send data over to output fields.
Close a netcdf or domain file opened with open_file or open_virtual_file.
Opens a given netcdf or domain file.
Read data from a defined field in a file.
integer function, public check_nml_error(IOSTAT, NML_NAME)
Checks the iostat argument that is returned after reading a namelist and determines if the error code...
subroutine, public write_version_number(version, tag, unit)
Prints to the log file (or a specified unit) the version id string and tag name.
subroutine, public error_mesg(routine, message, level)
Print notes, warnings and error messages; terminates program for warning and error messages....
integer function, public get_mosaic_xgrid_size(fileobj)
return exchange grid size of mosaic xgrid file.
subroutine, public get_mosaic_contact(fileobj, tile1, tile2, istart1, iend1, jstart1, jend1, istart2, iend2, jstart2, jend2)
Get contact information from mosaic_file Example usage: call get_mosaic_contact(mosaic_file,...
integer function, public get_mosaic_ntiles(fileobj)
Get number of tiles in the mosaic_file.
subroutine, public get_mosaic_tile_grid(grid_file, fileobj, domain, tile_count)
Gets the name of a mosaic tile grid file.
integer function, public get_mosaic_ncontacts(fileobj)
Get number of contacts in the mosaic_file.
subroutine, public get_mosaic_grid_sizes(fileobj, nx, ny)
Get grid size of each tile from mosaic_file.
integer function mpp_get_domain_npes(domain)
Set user stack size.
integer function, dimension(size(domain%tile_id(:))) mpp_get_tile_id(domain)
Returns the tile_id on current pe.
integer function mpp_get_current_ntile(domain)
Returns number of tile on current pe.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
subroutine mpp_compute_extent(isg, ieg, ndivs, ibegin, iend, extent)
Computes extents for a grid decomposition with the given indices and divisions.
subroutine mpp_get_tile_list(domain, tiles)
Return the tile_id on current pelist. one-tile-per-pe is assumed.
logical function mpp_domain_is_initialized(domain)
Set user stack size.
subroutine mpp_get_domain_pelist(domain, pelist)
Set user stack size.
integer function mpp_get_domain_root_pe(domain)
Set user stack size.
Broadcasts domain to every pe. Only useful outside the context of it's own pelist.
Deallocate given 1D or 2D domain.
Set up a domain decomposition.
These routines retrieve the axis specifications associated with the compute domains....
Retrieve the entire array of compute domain extents associated with a decomposition.
These routines retrieve the axis specifications associated with the data domains. The domain is a der...
These routines retrieve the axis specifications associated with the global domains....
Global sum of domain-decomposed arrays. mpp_global_sum is used to get the sum of a domain-decomposed...
Modifies the extents (compute, data and global) of a given domain.
Passes data from a structured grid to an unstructured grid Example usage:
Reorganization of distributed global arrays. mpp_redistribute is used to reorganize a distributed ar...
Performs halo updates for a given domain.
One dimensional domain used to manage shared data access between pes.
The domain2D type contains all the necessary information to define the global, compute and data domai...
Domain information for managing data on unstructured grids.
subroutine mpp_sync_self(pelist, check, request, msg_size, msg_type)
This is to check if current PE's outstanding puts are complete but we can't use shmem_fence because w...
integer function stdout()
This function returns the current standard fortran unit numbers for output.
subroutine mpp_set_current_pelist(pelist, no_sync)
Set context pelist.
integer function stdlog()
This function returns the current standard fortran unit numbers for log messages. Log messages,...
integer function mpp_npes()
Returns processor count for current pelist.
integer function mpp_pe()
Returns processor ID.
subroutine mpp_sync(pelist, do_self)
Synchronize PEs in list.
integer function mpp_clock_id(name, flags, grain)
Return an ID for a new or existing clock.
Scatter a vector across all PEs.
Reduction operations. Find the max of scalar a from the PEs in pelist result is also automatically br...
Reduction operations. Find the min of scalar a from the PEs in pelist result is also automatically br...
Receive data from another PE.
Send data to a receiving PE.
Holds stocks amounts per PE values.
subroutine, public get_time(Time, seconds, days, ticks, err_msg)
Returns days and seconds ( < 86400 ) corresponding to a time. err_msg should be checked for any error...
Type to represent amounts of time. Implemented as seconds and days to allow for larger intervals.
subroutine get_area_elements_fms2_io(fileobj, name, get_area_data)
Read the area elements from NetCDF file.
integer function, public xgrid_count(xmap)
Returns current size of exchange grid variables.
subroutine, public some(xmap, some_arr, grid_id)
Returns logical associating exchange grid cells with given side two grid.
subroutine regen(xmap)
Regenerate/Update the xmap.
integer nsubset
Number of processors to read exchange grid information. Those processors that read the exchange grid ...
subroutine put_side1_to_xgrid(d, grid_id, x, xmap, remap_method, complete)
Scatters data to exchange grid.
integer function get_nest_contact_fms2_io(fileobj, tile_nest_out, tile_parent_out, is_nest_out, ie_nest_out, js_nest_out, je_nest_out, is_parent_out, ie_parent_out, js_parent_out, je_parent_out)
currently we are assuming there is only one nest region
subroutine, public get_ocean_model_area_elements(domain, grid_file)
Read Ocean area element data from netCDF file.
subroutine, public set_frac_area_ug(f, grid_id, xmap)
Changes sub-grid portion areas and/or number.
real(r8_kind) function, dimension(3) conservation_check_ug_side1(d, grid_id, xmap, remap_method)
conservation_check_ug - returns three numbers which are the global sum of a variable (1) on its home ...
subroutine get_side1_from_xgrid(d, grid_id, x, xmap, complete)
subroutine stock_move_3d(from, to, grid_index, stock_data3d, xmap, delta_t, from_side, to_side, radius, verbose, ier)
this version takes rank 3 data, it can be used to compute the flux on anything but the first grid,...
real(r8_kind), dimension(:,:), allocatable, public area_atm_sphere
Area elements based on a the spherical model used by the ICE layer.
logical function in_box_me(i, j, grid)
logical function in_box(i, j, is, ie, js, je)
logical use_mpp_io
use_mpp_io Default = .false. When true, uses mpp_io for IO. When false, uses fms2_io for IO.
integer remapping_method
xgrid nml
subroutine get_grid_version2(grid, grid_id, grid_file)
read the center point of the grid from version 1 grid file. only the grid at the side 1 is needed,...
subroutine get_side2_from_xgrid(d, grid_id, x, xmap)
real(r8_kind) function, dimension(3) conservation_check_ug_side2(d, grid_id, xmap, remap_method)
conservation_check_ug - returns three numbers which are the global sum of a variable (1) on its home ...
real(r8_kind) function, dimension(is:ie, js:je) grad_merid_latlon(d, lat, is, ie, js, je, isd, jsd)
This function is used to calculate the gradient along meridinal direction. Maybe need to setup a limi...
logical make_exchange_reproduce
Set to .true. to make xgrid_mod reproduce answers on different numbers of PEs. This option has a cons...
character(len=64) interp_method
Exchange grid interpolation method. It has two options: "first_order", "second_order".
integer, parameter version2
mosaic grid file
subroutine, public setup_xmap(xmap, grid_ids, grid_domains, grid_file, atm_grid, lnd_ug_domain)
Sets up exchange grid connectivity using grid specification file and processor domain decomposition.
subroutine get_grid_version1(grid, grid_id, grid_file)
read the center point of the grid from version 1 grid file. only the grid at the side 1 is needed,...
subroutine put_side1_to_xgrid_ug(d, grid_id, x, xmap, complete)
Currently only support first order.
subroutine put_side2_to_xgrid(d, grid_id, x, xmap)
Scatters data to exchange grid.
subroutine set_frac_area_sg(f, grid_id, xmap)
Changes sub-grid portion areas and/or number.
logical function in_box_nbr(i, j, grid, p)
integer, parameter version1
grid spec file
real(r8_kind) function, dimension(3) conservation_check_side2(d, grid_id, xmap, remap_method)
conservation_check - returns three numbers which are the global sum of a variable (1) on its home mod...
real(r8_kind) function, dimension(3) conservation_check_side1(d, grid_id, xmap, remap_method)
conservation_check - returns three numbers which are the global sum of a variable (1) on its home mod...
subroutine, public xgrid_init(remap_method)
Initialize the xgrid_mod.
real(r8_kind), dimension(:,:), allocatable, public area_atm_model
Area elements used inside each model.
real(r8_kind) function, dimension(is:ie, js:je) grad_zonal_latlon(d, lon, lat, is, ie, js, je, isd, jsd)
This function is used to calculate the gradient along zonal direction. Maybe need to setup a limit fo...
subroutine, public get_xmap_grid_area(id, xmap, area)
This routine is used to get the grid area of component model with id.
Returns three numbers which are the global sum of a variable.
For an unstructured grid, returns three numbers which are the global sum of a variable (1) on its hom...
Sums data from exchange grid to model grid.
get_from_xgrid for unstructured grids.
Scatters data from model grid onto exchange grid.
put_to_xgrid for unstructured grids.
Sets sub-grid area and numbering in the given exchange grid.
Private type used for exchange grid communication.
Type to hold pointers for grid boxes.
Private type to hold all data needed from given grid for an exchange grid.
Private type for overlap exchange grid data.
Private type for exchange grid data.
Private type for exchange grid data.
Private type for cell indices and data in the exchange grid.
Type for an exchange grid, holds pointers to included grids and any necessary data.