37 type(mct_smatp),
pointer :: smatp(:)
38 integer(kind=ip_i4_p) :: nwgts
39 character(len=ic_long):: file
40 character(len=ic_med) :: loc
41 character(len=ic_med) :: opt
42 character(len=ic_med) :: optval
44 integer(kind=ip_i4_p) :: spart
45 integer(kind=ip_i4_p) :: dpart
47 type(mct_avect) :: av_ms
48 type(mct_avect) :: av_md
51 integer(kind=ip_i4_p) ,
public :: prism_mmapper
52 integer(kind=ip_i4_p) ,
public :: prism_nmapper = 0
57 integer,
private,
parameter :: r8 = ip_double_p
58 integer,
private,
parameter :: in = ip_i4_p
62 real(R8),
private,
allocatable :: snew(:,:),sold(:,:)
63 integer,
private,
allocatable :: rnew(:),rold(:)
64 integer,
private,
allocatable :: cnew(:),cold(:)
82 integer(ip_i4_p),
intent(in) :: mapid
83 integer(ip_i4_p),
intent(in) :: namid
86 integer(ip_i4_p) :: src_size,src_rank, ncrn_src
87 integer(ip_i4_p) ,
allocatable :: src_dims(:),src_mask(:)
88 real(ip_double_p),
allocatable :: src_lon(:),src_lat(:)
89 real(ip_double_p),
allocatable :: src_corner_lon(:,:),src_corner_lat(:,:)
90 integer(ip_i4_p) :: dst_size,dst_rank, ncrn_dst
91 integer(ip_i4_p) ,
allocatable :: dst_dims(:),dst_mask(:)
92 real(ip_double_p),
allocatable :: dst_lon(:),dst_lat(:)
93 real(ip_double_p),
allocatable :: dst_corner_lon(:,:),dst_corner_lat(:,:)
94 integer(ip_i4_p) ,
allocatable :: ifld2(:,:)
95 real(ip_double_p),
allocatable :: fld2(:,:),fld3(:,:,:)
96 integer(ip_i4_p) :: i,j,k,icnt,nx,ny,nc
97 logical :: lextrapdone
99 character(len=ic_med) :: filename
100 character(len=ic_med) :: fldname
101 character(len=*),
parameter :: subname =
'(oasis_map_genmap)'
105 lextrapdone = .false.
112 if (trim(namscrtyp(namid)) /=
'SCALAR')
then
113 write(nulprt,*) subname,estr,
'only scrip type SCALAR mapping supported'
119 if (trim(namscrmet(namid)) ==
'CONSERV')
then
125 filename =
'grids.nc'
127 fldname = trim(namsrcgrd(namid))//
'.clo'
130 fldname = trim(namsrcgrd(namid))//
'.lon'
133 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',&
134 trim(fldname),nx,ny,nc,do_corners
137 allocate(src_dims(src_rank))
141 allocate(src_mask(src_size))
142 allocate(src_lon(src_size))
143 allocate(src_lat(src_size))
144 allocate(src_corner_lon(ncrn_src,src_size))
145 allocate(src_corner_lat(ncrn_src,src_size))
147 allocate(ifld2(nx,ny))
148 filename =
'masks.nc'
149 fldname = trim(namsrcgrd(namid))//
'.msk'
151 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
152 src_mask(icnt) = ifld2(i,j)
154 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
155 minval(src_mask),maxval(src_mask)
158 allocate(fld2(nx,ny))
159 filename =
'grids.nc'
160 fldname = trim(namsrcgrd(namid))//
'.lon'
162 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
163 src_lon(icnt) = fld2(i,j)
165 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
166 minval(src_lon),maxval(src_lon)
167 fldname = trim(namsrcgrd(namid))//
'.lat'
169 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
170 src_lat(icnt) = fld2(i,j)
172 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
173 minval(src_lat),maxval(src_lat)
177 allocate(fld3(nx,ny,nc))
178 filename =
'grids.nc'
179 fldname = trim(namsrcgrd(namid))//
'.clo'
181 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
183 src_corner_lon(k,icnt) = fld3(i,j,k)
186 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
187 minval(src_corner_lon),maxval(src_corner_lon)
188 fldname = trim(namsrcgrd(namid))//
'.cla'
190 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
192 src_corner_lat(k,icnt) = fld3(i,j,k)
195 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
196 minval(src_corner_lat),maxval(src_corner_lat)
199 src_corner_lon = -9999.
200 src_corner_lat = -9999.
205 filename =
'grids.nc'
207 fldname = trim(namdstgrd(namid))//
'.clo'
210 fldname = trim(namdstgrd(namid))//
'.lon'
213 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname),nx,ny,nc
216 allocate(dst_dims(dst_rank))
220 allocate(dst_mask(dst_size))
221 allocate(dst_lon(dst_size))
222 allocate(dst_lat(dst_size))
223 allocate(dst_corner_lon(ncrn_dst,dst_size))
224 allocate(dst_corner_lat(ncrn_dst,dst_size))
226 allocate(ifld2(nx,ny))
227 filename =
'masks.nc'
228 fldname = trim(namdstgrd(namid))//
'.msk'
230 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
231 dst_mask(icnt) = ifld2(i,j)
233 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
234 minval(dst_mask),maxval(dst_mask)
237 allocate(fld2(nx,ny))
238 filename =
'grids.nc'
239 fldname = trim(namdstgrd(namid))//
'.lon'
241 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
242 dst_lon(icnt) = fld2(i,j)
244 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
245 minval(dst_lon),maxval(dst_lon)
246 fldname = trim(namdstgrd(namid))//
'.lat'
248 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
249 dst_lat(icnt) = fld2(i,j)
251 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
252 minval(dst_lat),maxval(dst_lat)
256 allocate(fld3(nx,ny,nc))
257 filename =
'grids.nc'
258 fldname = trim(namdstgrd(namid))//
'.clo'
260 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
262 dst_corner_lon(k,icnt) = fld3(i,j,k)
265 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
266 minval(dst_corner_lon),maxval(dst_corner_lon)
267 fldname = trim(namdstgrd(namid))//
'.cla'
269 icnt = 0;
do j = 1,ny;
do i = 1,nx; icnt = icnt + 1
271 dst_corner_lat(k,icnt) = fld3(i,j,k)
274 if (oasis_debug >= 15)
write(nulprt,*) subname,
' read ',trim(filename),
' ',trim(fldname), &
275 minval(dst_corner_lat),maxval(dst_corner_lat)
278 dst_corner_lon = -9999.
279 dst_corner_lat = -9999.
282 IF (oasis_debug >= 15)
THEN
283 WRITE(nulprt,*) subname,
' call grid_init '
288 src_mask = 1 - src_mask
289 dst_mask = 1 - dst_mask
290 call grid_init(namscrmet(namid),namscrres(namid),namscrbin(namid), &
291 src_size, dst_size, src_dims, dst_dims, &
292 src_rank, dst_rank, ncrn_src, ncrn_dst, &
293 src_mask, dst_mask, namsrcgrd(namid), namdstgrd(namid), &
294 src_lat, src_lon, dst_lat, dst_lon, &
295 src_corner_lat, src_corner_lon, &
296 dst_corner_lat, dst_corner_lon, &
298 if (oasis_debug >= 15)
then
299 WRITE(nulprt,*) subname,
' done grid_init '
303 IF (oasis_debug >= 15)
THEN
304 WRITE(nulprt,*) subname,
' call scrip '
307 call scrip(prism_mapper(mapid)%file,prism_mapper(mapid)%file,namscrmet(namid), &
308 namscrnor(namid),lextrapdone,namscrvam(namid),namscrnbr(namid))
309 IF (oasis_debug >= 15)
THEN
310 WRITE(nulprt,*) subname,
' done scrip '
314 deallocate(src_dims, dst_dims)
318 deallocate(src_corner_lon)
319 deallocate(src_corner_lat)
323 deallocate(dst_corner_lon)
324 deallocate(dst_corner_lat)
369 filename,mytask,mpicom,nwgts, &
370 areasrc,areadst,ni_i,nj_i,ni_o,nj_o )
375 integer,
parameter :: R8 = ip_double_p
376 integer,
parameter :: IN = ip_i4_p
380 type(mct_smat) ,
intent(out),
pointer :: sMat(:)
381 type(mct_gsmap) ,
intent(in) ,
target :: SgsMap
382 type(mct_gsmap) ,
intent(in) ,
target :: DgsMap
383 character(*) ,
intent(in) :: newdom
385 character(*) ,
intent(in) :: filename
386 integer(IN) ,
intent(in) :: mytask
387 integer(IN) ,
intent(in) :: mpicom
388 integer(IN) ,
intent(out) :: nwgts
389 type(mct_avect) ,
intent(out),
optional :: areasrc
390 type(mct_avect) ,
intent(out),
optional :: areadst
391 integer(IN) ,
intent(out),
optional :: ni_i
392 integer(IN) ,
intent(out),
optional :: nj_i
393 integer(IN) ,
intent(out),
optional :: ni_o
394 integer(IN) ,
intent(out),
optional :: nj_o
421 real(R8) ,
allocatable :: rtemp(:)
422 real(R8) ,
allocatable :: Sbuf(:,:)
423 real(R8) ,
allocatable :: remaps(:,:)
424 integer,
allocatable :: Rbuf(:)
425 integer,
allocatable :: Cbuf(:)
430 integer,
allocatable :: lsstart(:)
431 integer,
allocatable :: lscount(:)
432 type(mct_gsmap),
pointer :: mygsmap
437 real(R8) ,
allocatable :: Snew(:,:),Sold(:,:)
438 integer,
allocatable :: Rnew(:),Rold(:)
439 integer,
allocatable :: Cnew(:),Cold(:)
441 character,
allocatable :: str(:)
442 character(len=ic_long):: attstr
448 integer,
parameter :: rbuf_size = 100000
451 type(mct_avect) :: areasrc0
452 type(mct_avect) :: areadst0
454 character(*),
parameter :: areaAV_field =
'aream'
457 character(*),
parameter :: subName =
'(oasis_map_sMatReaddnc_orig)'
466 if (mytask == 0)
then
467 if (oasis_debug >= 2)
write(nulprt,*) subname,
" reading mapping matrix data decomposed..."
472 if (oasis_debug >=2 )
write(nulprt,*) subname,
" * file name : ",trim(filename)
473 status = nf90_open(trim(filename),nf90_nowrite,fid)
474 if (status /= nf90_noerr)
then
475 write(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
476 WRITE(nulprt,*) subname,estr,
'mapping file not found = ',trim(filename)
482 status = nf90_inq_dimid(fid,
'num_links', did)
483 status = nf90_inquire_dimension(fid, did , len = ns)
485 status = nf90_inq_dimid(fid,
'src_grid_size', did)
486 status = nf90_inquire_dimension(fid, did , len = na)
488 status = nf90_inq_dimid(fid,
'dst_grid_size', did)
489 status = nf90_inquire_dimension(fid, did , len = nb)
490 status = nf90_inq_dimid(fid,
'num_wgts', did)
491 status = nf90_inquire_dimension(fid, did , len = nwgts)
493 if (
present(ni_i) .and.
present(nj_i) .and.
present(ni_o) .and.
present(nj_o))
then
502 status = nf90_inq_varid(fid,
'src_grid_dims', vid)
503 status = nf90_get_var(fid, vid, dims)
506 status = nf90_inq_varid(fid,
'dst_grid_dims', vid)
507 status = nf90_get_var(fid, vid, dims)
512 if (oasis_debug >= 2)
write(nulprt,*) subname,
" * matrix dims src x dst : ",na,
' x',nb
513 if (oasis_debug >= 2)
write(nulprt,*) subname,
" * number of non-zero elements: ",ns
520 if (
present(areasrc))
then
521 if (mytask == 0)
then
522 call mct_avect_init(areasrc0,
' ',areaav_field,na)
524 status = nf90_inq_varid(fid,
'src_grid_area',vid)
525 IF (status /= nf90_noerr)
THEN
526 WRITE(nulprt,*) subname,
' nf90_strerrr = ',trim(nf90_strerror(status))
527 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
530 status = nf90_get_var(fid, vid, areasrc0%rAttr)
531 IF (status /= nf90_noerr)
THEN
532 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
533 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
537 call mct_avect_scatter(areasrc0, areasrc, sgsmap, 0, mpicom, status)
538 if (status /= 0) call mct_die(subname,
"Error on scatter of areasrc0")
539 if (mytask == 0)
then
546 call mct_avect_clean(areasrc0)
553 if (
present(areadst))
then
554 if (mytask == 0)
then
555 call mct_avect_init(areadst0,
' ',areaav_field,nb)
557 status = nf90_inq_varid(fid,
'dst_grid_area',vid)
558 IF (status /= nf90_noerr)
THEN
559 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
560 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
563 status = nf90_get_var(fid, vid, areadst0%rAttr)
564 IF (status /= nf90_noerr)
THEN
565 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
566 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
570 call mct_avect_scatter(areadst0, areadst, dgsmap, 0, mpicom, status)
571 if (status /= 0) call mct_die(subname,
"Error on scatter of areadst0")
572 if (mytask == 0)
then
579 call mct_avect_clean(areadst0)
586 if (
present(ni_i) .and.
present(nj_i) .and.
present(ni_o) .and.
present(nj_o))
then
602 if (newdom ==
'src')
then
604 elseif (newdom ==
'dst')
then
607 write(nulprt,*) subname,estr,
'invalid newdom value, expect src or dst = ',newdom
611 do n = 1,
size(mygsmap%start)
612 if (mygsmap%pe_loc(n) == mytask)
then
616 allocate(lsstart(lsize),lscount(lsize),stat=status)
617 if (status /= 0) call mct_perr_die(subname,
':: allocate Lsstart',status)
623 do n = 1,
size(mygsmap%start)
624 if (mygsmap%pe_loc(n) == mytask)
then
628 do while (.not.found .and. l1 < lsize)
629 if (mygsmap%start(n) < lsstart(l1))
then
630 do l2 = lsize, l1+1, -1
631 lsstart(l2) = lsstart(l2-1)
632 lscount(l2) = lscount(l2-1)
639 lsstart(l1) = mygsmap%start(n)
640 lscount(l1) = mygsmap%length(n)
644 if (lsstart(n) > lsstart(n+1))
then
645 write(nulprt,*) subname,estr,
'lsstart not properly sorted'
653 rsize = min(rbuf_size,ns)
654 bsize = ((ns/commsize) + 1 ) * 1.2
658 nread = (ns-1)/rsize + 1
664 if (mytask == 0)
then
665 allocate(remaps(nwgts,rsize),stat=status)
666 if (status /= 0) call mct_perr_die(subname,
':: allocate remaps',status)
668 allocate(smat(nwgts),stat=status)
669 if (status /= 0) call mct_perr_die(subname,
':: allocate Smat',status)
670 allocate(sbuf(nwgts,rsize),rbuf(rsize),cbuf(rsize),stat=status)
671 if (status /= 0) call mct_perr_die(subname,
':: allocate Sbuf',status)
672 allocate(snew(nwgts,bsize),cnew(bsize),rnew(bsize),stat=status)
673 if (status /= 0) call mct_perr_die(subname,
':: allocate Snew1',status)
680 start(1) = (n-1)*rsize + 1
681 count(1) = min(rsize,ns-start(1)+1)
692 status = nf90_inq_varid(fid,
'remap_matrix' ,vid)
694 status = nf90_get_var(fid,vid,remaps,start2,count2)
695 sbuf(:,:) = remaps(:,:)
696 IF (status /= nf90_noerr)
THEN
697 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
698 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
703 status = nf90_inq_varid(fid,
'dst_address',vid)
704 status = nf90_get_var(fid,vid,rbuf,start,count)
705 IF (status /= nf90_noerr)
THEN
706 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
707 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
712 status = nf90_inq_varid(fid,
'src_address',vid)
713 status = nf90_get_var(fid,vid,cbuf,start,count)
714 IF (status /= nf90_noerr)
THEN
715 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
716 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
734 if (newdom ==
'src')
then
736 elseif (newdom ==
'dst')
then
747 if (cnt > bsize)
then
749 allocate(sold(1:nwgts,cntold),rold(cntold),cold(cntold),stat=status)
750 if (status /= 0) call mct_perr_die(subname,
':: allocate old',status)
751 sold(1:nwgts,1:cntold) = snew(1:nwgts,1:cntold)
752 rold(1:cntold) = rnew(1:cntold)
753 cold(1:cntold) = cnew(1:cntold)
756 deallocate(snew,rnew,cnew,stat=status)
757 if (status /= 0) call mct_perr_die(subname,
':: allocate new',status)
759 if (oasis_debug > 15)
write(nulprt,*) subname,
' reallocate bsize to ',bsize
760 allocate(snew(nwgts,bsize),rnew(bsize),cnew(bsize),stat=status)
761 if (status /= 0) call mct_perr_die(subname,
':: allocate old',status)
764 snew(1:nwgts,1:cntold) = sold(1:nwgts,1:cntold)
765 rnew(1:cntold) = rold(1:cntold)
766 cnew(1:cntold) = cold(1:cntold)
767 deallocate(sold,rold,cold,stat=status)
768 if (status /= 0) call mct_perr_die(subname,
':: deallocate old',status)
771 snew(1:nwgts,cnt) = sbuf(1:nwgts,m)
781 if (mytask == 0)
then
782 deallocate(remaps, stat=status)
783 if (status /= 0) call mct_perr_die(subname,
':: deallocate remaps',status)
785 deallocate(sbuf,rbuf,cbuf, stat=status)
786 if (status /= 0) call mct_perr_die(subname,
':: deallocate Sbuf',status)
795 call mct_smat_init(smat(n), nb, na, cnt)
798 igrow = mct_smat_indexia(smat(1),
'grow')
799 igcol = mct_smat_indexia(smat(1),
'gcol')
800 iwgt = mct_smat_indexra(smat(1),
'weight')
804 smat(n)%data%rAttr(iwgt ,1:cnt) = snew(n,1:cnt)
805 smat(n)%data%iAttr(igrow,1:cnt) = rnew(1:cnt)
806 smat(n)%data%iAttr(igcol,1:cnt) = cnew(1:cnt)
813 deallocate(snew,rnew,cnew, stat=status)
814 deallocate(lsstart,lscount,stat=status)
815 if (status /= 0) call mct_perr_die(subname,
':: deallocate new',status)
817 if (mytask == 0)
then
818 status = nf90_close(fid)
819 IF (oasis_debug >= 2)
THEN
820 WRITE(nulprt,*) subname,
" ... done reading file"
850 integer,
parameter :: R8 = ip_double_p
851 integer,
parameter :: IN = ip_i4_p
856 integer(IN) :: starti(:)
857 integer(IN) :: counti(:)
862 integer(IN) :: nl,nc,nr,ncprev
867 character(*),
parameter :: subName =
'(check_myindex) '
886 do while (.not.stopnow)
887 if (index < starti(nc))
then
889 elseif (index > (starti(nc) + counti(nc) - 1))
then
898 if (nc == ncprev .or. nc < 1 .or. nc > lsize) stopnow = .true.
942 filename,mytask,mpicom,nwgts, &
943 areasrc,areadst,ni_i,nj_i,ni_o,nj_o )
949 type(mct_smat) ,
intent(out),
pointer :: sMat(:)
950 type(mct_gsmap) ,
intent(in) ,
target :: SgsMap
951 type(mct_gsmap) ,
intent(in) ,
target :: DgsMap
952 character(*) ,
intent(in) :: newdom
954 character(*) ,
intent(in) :: filename
955 integer(IN) ,
intent(in) :: mytask
956 integer(IN) ,
intent(in) :: mpicom
957 integer(IN) ,
intent(out) :: nwgts
958 type(mct_avect) ,
intent(out),
optional :: areasrc
959 type(mct_avect) ,
intent(out),
optional :: areadst
960 integer(IN) ,
intent(out),
optional :: ni_i
961 integer(IN) ,
intent(out),
optional :: nj_i
962 integer(IN) ,
intent(out),
optional :: ni_o
963 integer(IN) ,
intent(out),
optional :: nj_o
989 real(R8) ,
allocatable :: rtemp(:)
990 real(R8) ,
allocatable :: Sbuf(:,:)
991 real(R8) ,
allocatable :: remaps(:,:)
992 integer,
allocatable :: Rbuf(:)
993 integer,
allocatable :: Cbuf(:)
994 real(R8),
allocatable :: SReadData(:,:)
995 integer,
allocatable :: RReadData(:)
996 integer,
allocatable :: CReadData(:)
997 integer,
allocatable :: pesave(:)
998 real(R8),
allocatable :: SDistData(:,:)
999 integer,
allocatable :: RDistData(:)
1000 integer,
allocatable :: CDistData(:)
1004 type(mct_gsmap),
pointer :: mygsmap
1005 integer :: l1,l2,lsize
1006 integer,
allocatable :: lsstart(:)
1007 integer,
allocatable :: lscount(:)
1008 integer,
allocatable :: lspeloc(:)
1013 integer,
allocatable :: cntrs(:)
1014 integer :: mpistatus(mpi_status_size)
1016 character,
allocatable :: str(:)
1017 character(len=ic_long):: attstr
1023 integer,
parameter :: rbuf_size = 100000
1024 integer :: dimbuffer(8)
1027 type(mct_avect) :: areasrc0
1028 type(mct_avect) :: areadst0
1030 character(*),
parameter :: areaAV_field =
'aream'
1033 character(*),
parameter :: subName =
'(oasis_map_sMatReaddnc_ceg)'
1034 character*80 :: fname
1043 if (mytask == 0)
then
1044 if (oasis_debug >= 2)
write(nulprt,*) subname,
" reading mapping matrix data decomposed..."
1049 if (oasis_debug >=2 )
write(nulprt,*) subname,
" * file name : ",trim(filename)
1050 status = nf90_open(trim(filename),nf90_nowrite,fid)
1051 if (status /= nf90_noerr)
then
1052 write(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1053 WRITE(nulprt,*) subname,estr,
'mapping file not found = ',trim(filename)
1059 status = nf90_inq_dimid(fid,
'num_links', did)
1060 status = nf90_inquire_dimension(fid, did , len = ns)
1062 status = nf90_inq_dimid(fid,
'src_grid_size', did)
1063 status = nf90_inquire_dimension(fid, did , len = na)
1065 status = nf90_inq_dimid(fid,
'dst_grid_size', did)
1066 status = nf90_inquire_dimension(fid, did , len = nb)
1067 status = nf90_inq_dimid(fid,
'num_wgts', did)
1068 status = nf90_inquire_dimension(fid, did , len = nwgts)
1070 if (
present(ni_i) .and.
present(nj_i) .and.
present(ni_o) .and.
present(nj_o))
then
1079 status = nf90_inq_varid(fid,
'src_grid_dims', vid)
1080 status = nf90_get_var(fid, vid, dims)
1083 status = nf90_inq_varid(fid,
'dst_grid_dims', vid)
1084 status = nf90_get_var(fid, vid, dims)
1089 if (oasis_debug >= 2)
write(nulprt,*) subname,
" * matrix dims src x dst : ",na,
' x',nb
1090 if (oasis_debug >= 2)
write(nulprt,*) subname,
" * number of non-zero elements: ",ns
1097 if (
present(areasrc))
then
1098 if (mytask == 0)
then
1099 call mct_avect_init(areasrc0,
' ',areaav_field,na)
1101 status = nf90_inq_varid(fid,
'src_grid_area',vid)
1102 IF (status /= nf90_noerr)
THEN
1103 WRITE(nulprt,*) subname,
' nf90_strerrr = ',trim(nf90_strerror(status))
1104 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1107 status = nf90_get_var(fid, vid, areasrc0%rAttr)
1108 IF (status /= nf90_noerr)
THEN
1109 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1110 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1114 call mct_avect_scatter(areasrc0, areasrc, sgsmap, 0, mpicom, status)
1115 if (status /= 0) call mct_die(subname,
"Error on scatter of areasrc0")
1116 if (mytask == 0)
then
1122 call mct_avect_clean(areasrc0)
1129 if (
present(areadst))
then
1130 if (mytask == 0)
then
1131 call mct_avect_init(areadst0,
' ',areaav_field,nb)
1133 status = nf90_inq_varid(fid,
'dst_grid_area',vid)
1134 IF (status /= nf90_noerr)
THEN
1135 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1136 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1139 status = nf90_get_var(fid, vid, areadst0%rAttr)
1140 IF (status /= nf90_noerr)
THEN
1141 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1142 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1146 call mct_avect_scatter(areadst0, areadst, dgsmap, 0, mpicom, status)
1147 if (status /= 0) call mct_die(subname,
"Error on scatter of areadst0")
1148 if (mytask == 0)
then
1154 call mct_avect_clean(areadst0)
1163 if (mpi_rank_local.eq.0)
then
1167 dimbuffer(4) = nwgts
1168 if (
present(ni_i) .and.
present(nj_i) .and.
present(ni_o) .and.
present(nj_o))
then
1176 if (mpi_rank_local.ne.0)
then
1180 nwgts = dimbuffer(4)
1181 if (
present(ni_i) .and.
present(nj_i) .and.
present(ni_o) .and.
present(nj_o))
then
1195 if (newdom ==
'src')
then
1197 elseif (newdom ==
'dst')
then
1200 write(nulprt,*) subname,estr,
'invalid newdom value, expect src or dst = ',newdom
1207 rsize = min(rbuf_size,ns)
1208 bsize = ((ns/commsize) + 1 ) * 1.2
1212 nread = (ns-1)/rsize + 1
1218 allocate(smat(nwgts),stat=status)
1219 if (status /= 0) call mct_perr_die(subname,
':: allocate Smat',status)
1220 allocate(sbuf(nwgts,rsize),rbuf(rsize),cbuf(rsize),stat=status)
1221 if (status /= 0) call mct_perr_die(subname,
':: allocate Sbuf',status)
1222 allocate(snew(nwgts,bsize),cnew(bsize),rnew(bsize),stat=status)
1223 if (status /= 0) call mct_perr_die(subname,
':: allocate Snew1',status)
1231 if (mytask == 0)
then
1236 lsize =
size(mygsmap%start)
1237 allocate(lsstart(lsize),lscount(lsize),lspeloc(lsize),stat=status)
1238 if (status /= 0) call mct_perr_die(subname,
':: allocate Lsstart',status)
1240 do n = 1,
size(mygsmap%start)
1243 do while (.not.found .and. l1 < n)
1244 if (mygsmap%start(n) < lsstart(l1))
then
1246 lsstart(l2) = lsstart(l2-1)
1247 lscount(l2) = lscount(l2-1)
1248 lspeloc(l2) = lspeloc(l2-1)
1255 lsstart(l1) = mygsmap%start(n)
1256 lscount(l1) = mygsmap%length(n)
1257 lspeloc(l1) = mygsmap%pe_loc(n)
1259 do n = 1,
size(mygsmap%start)-1
1260 if (lsstart(n) > lsstart(n+1))
then
1261 write(nulprt,*) subname,estr,
'lsstart not properly sorted'
1269 allocate(remaps(nwgts,rsize),stat=status)
1270 if (status /= 0) call mct_perr_die(subname,
':: allocate remaps',status)
1271 allocate(sreaddata(nwgts,rsize),rreaddata(rsize),creaddata(rsize),stat=status)
1272 if (status /= 0) call mct_perr_die(subname,
':: allocate SReadData',status)
1273 allocate(pesave(rsize),stat=status)
1274 if (status /= 0) call mct_perr_die(subname,
':: allocate pesave',status)
1275 allocate(sdistdata(nwgts,rsize),rdistdata(rsize),cdistdata(rsize), stat=status)
1276 if (status /= 0) call mct_perr_die(subname,
':: allocate SDistData',status)
1277 allocate(cntrs(0:mpi_size_local), stat=status)
1278 if (status /= 0) call mct_perr_die(subname,
':: allocate cntrs',status)
1284 start(1) = (n-1)*rsize + 1
1285 count(1) = min(rsize,ns-start(1)+1)
1288 start2(2) = start(1)
1289 count2(2) = count(1)
1295 status = nf90_inq_varid(fid,
'remap_matrix' ,vid)
1297 status = nf90_get_var(fid,vid,remaps,start2,count2)
1298 sreaddata(:,:) = remaps(:,:)
1299 IF (status /= nf90_noerr)
THEN
1300 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1301 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1305 status = nf90_inq_varid(fid,
'dst_address',vid)
1306 status = nf90_get_var(fid,vid,rreaddata,start,count)
1307 IF (status /= nf90_noerr)
THEN
1308 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1309 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1314 status = nf90_inq_varid(fid,
'src_address',vid)
1315 status = nf90_get_var(fid,vid,creaddata,start,count)
1316 IF (status /= nf90_noerr)
THEN
1317 WRITE(nulprt,*) subname,
' nf90_strerror = ',trim(nf90_strerror(status))
1318 WRITE(nulprt,*) subname,
'model :',compid,
' proc :',mpi_rank_local
1336 if (newdom ==
'src')
then
1337 pe =
get_cegindex(rreaddata(m),lsstart,lscount,lspeloc)
1338 else if (newdom ==
'dst')
then
1339 pe =
get_cegindex(creaddata(m),lsstart,lscount,lspeloc)
1344 if (pe+1 < 1 .or. pe+1 > mpi_size_local)
then
1345 write(nulprt,*) subname,wstr,
'get_cegindex search error', m,count(1),creaddata(m),rreaddata(m),sreaddata(:,m)
1347 cntrs(pe+1) = cntrs(pe+1) + 1
1353 do pe = 1, mpi_size_local
1354 cntrs(pe) = cntrs(pe-1) + cntrs(pe)
1362 if (newdom ==
'src')
then
1363 pe =
get_cegindex(rreaddata(m),lsstart,lscount,lspeloc)
1364 else if (newdom ==
'dst')
then
1365 pe =
get_cegindex(creaddata(m),lsstart,lscount,lspeloc)
1371 if (pe+1 < 1 .or. pe+1 > mpi_size_local)
then
1376 sdistdata(:,cntrs(pe)) = sreaddata(:,m)
1377 rdistdata(cntrs(pe)) = rreaddata(m)
1378 cdistdata(cntrs(pe)) = creaddata(m)
1380 cntrs(pe) = cntrs(pe) + 1
1390 if (cntrs(0).gt.1)
then
1396 snew(1:nwgts,cnt:cnt+reclen-1) = sdistdata(1:nwgts,1:reclen)
1397 rnew(cnt:cnt+reclen-1) = rdistdata(1:reclen)
1398 cnew(cnt:cnt+reclen-1) = cdistdata(1:reclen)
1403 do pe = 1, mpi_size_local
1404 if (cntrs(pe).gt.cntrs(pe-1))
then
1406 m = cntrs(pe)-cntrs(pe-1)
1408 call mpi_send(m, 1, mpi_integer, pe, 4000, mpicom, ierr)
1409 call mpi_send(sdistdata(1,cntrs(pe-1)), nwgts*m, mpi_real8, pe, 1000, mpicom, ierr)
1410 call mpi_send(rdistdata(cntrs(pe-1)), m, mpi_integer, pe, 2000, mpicom, ierr)
1411 call mpi_send(cdistdata(cntrs(pe-1)), m, mpi_integer, pe, 3000, mpicom, ierr)
1421 do pe = 1, mpi_size_local-1
1423 call mpi_send(m, 1, mpi_integer, pe, 4000, mpicom, ierr)
1427 deallocate(lsstart,lscount,lspeloc, stat=status)
1428 if (status /= 0) call mct_perr_die(subname,
':: deallocate lsstart',status)
1429 deallocate(pesave, stat=status)
1430 if (status /= 0) call mct_perr_die(subname,
':: deallocate pesave',status)
1431 deallocate(sreaddata,rreaddata,creaddata, stat=status)
1432 if (status /= 0) call mct_perr_die(subname,
':: deallocate SReadData',status)
1433 deallocate(sdistdata,rdistdata,cdistdata, stat=status)
1434 if (status /= 0) call mct_perr_die(subname,
':: deallocate SDistData',status)
1435 deallocate(cntrs, stat=status)
1436 if (status /= 0) call mct_perr_die(subname,
':: deallocate cntrs',status)
1437 deallocate(remaps, stat=status)
1438 if (status /= 0) call mct_perr_die(subname,
':: deallocate remaps',status)
1449 call
oasis_mpi_recv(reclen, 0, 4000, mpicom, subname//
" MPI in reclen recv")
1451 do while (reclen.ne.-1)
1455 call mpi_recv(sbuf, reclen*nwgts, mpi_real8, 0, 1000, mpicom, mpi_status_ignore, ierr)
1456 call mpi_recv(rbuf, reclen, mpi_integer, 0, 2000, mpicom, mpi_status_ignore, ierr)
1457 call mpi_recv(cbuf, reclen, mpi_integer, 0, mpi_any_tag, mpicom, mpistatus, ierr)
1465 snew(1:nwgts,cnt:cnt+reclen-1) = sbuf(1:nwgts,1:reclen)
1466 rnew(cnt:cnt+reclen-1) = rbuf(1:reclen)
1467 cnew(cnt:cnt+reclen-1) = cbuf(1:reclen)
1470 call
oasis_mpi_recv(reclen, 0, 4000, mpicom, subname//
" MPI in reclen recv")
1481 deallocate(sbuf,rbuf,cbuf, stat=status)
1482 if (status /= 0) call mct_perr_die(subname,
':: deallocate Sbuf',status)
1491 call mct_smat_init(smat(n), nb, na, cnt)
1494 igrow = mct_smat_indexia(smat(1),
'grow')
1495 igcol = mct_smat_indexia(smat(1),
'gcol')
1496 iwgt = mct_smat_indexra(smat(1),
'weight')
1500 smat(n)%data%rAttr(iwgt ,1:cnt) = snew(n,1:cnt)
1501 smat(n)%data%iAttr(igrow,1:cnt) = rnew(1:cnt)
1502 smat(n)%data%iAttr(igcol,1:cnt) = cnew(1:cnt)
1509 deallocate(snew,rnew,cnew, stat=status)
1510 if (status /= 0) call mct_perr_die(subname,
':: deallocate new',status)
1512 if (mytask == 0)
then
1513 status = nf90_close(fid)
1514 IF (oasis_debug >= 2)
THEN
1515 WRITE(nulprt,*) subname,
" ... done reading file"
1532 integer,
intent(inout) :: cnt
1533 integer,
intent(in) :: reclen
1534 integer,
intent(inout) :: bsize
1535 integer,
intent(in) :: nwgts
1539 character(*),
parameter :: subName =
'(ceg_coupler_augment_arrays)'
1543 if (cnt+reclen > bsize)
then
1548 allocate(sold(1:nwgts,cnt),rold(cnt),cold(cnt),stat=status)
1549 if (status /= 0) call mct_perr_die(subname,
':: allocate old',status)
1550 sold(1:nwgts,1:cnt-1) = snew(1:nwgts,1:cnt-1)
1551 rold(1:cnt-1) = rnew(1:cnt-1)
1552 cold(1:cnt-1) = cnew(1:cnt-1)
1555 deallocate(snew,rnew,cnew,stat=status)
1556 if (status /= 0) call mct_perr_die(subname,
':: allocate new',status)
1557 bsize = 1.5 * (cnt+reclen)
1558 if (oasis_debug > 15)
write(nulprt,*) subname,
' reallocate bsize to ',bsize
1559 allocate(snew(nwgts,bsize),rnew(bsize),cnew(bsize),stat=status)
1560 if (status /= 0) call mct_perr_die(subname,
':: allocate old',status)
1563 snew(1:nwgts,1:cnt-1) = sold(1:nwgts,1:cnt-1)
1564 rnew(1:cnt-1) = rold(1:cnt-1)
1565 cnew(1:cnt-1) = cold(1:cnt-1)
1566 deallocate(sold,rold,cold,stat=status)
1567 if (status /= 0) call mct_perr_die(subname,
':: deallocate old',status)
1591 integer,
parameter :: R8 = ip_double_p
1592 integer,
parameter :: IN = ip_i4_p
1596 integer(IN) :: index
1597 integer(IN) :: starti(:)
1598 integer(IN) :: counti(:)
1599 integer(IN) :: peloci(:)
1605 integer(IN) :: nl,nc,nr,ncprev
1606 integer(IN) :: lsize
1610 character(*),
parameter :: subName =
'(get_cegindex) '
1618 lsize =
size(starti)
1628 do while (.not.stopnow)
1629 if (index < starti(nc))
then
1631 elseif (index > (starti(nc) + counti(nc) - 1))
then
1640 if (nc == ncprev .or. nc < 1 .or. nc > lsize) stopnow = .true.
subroutine, public oasis_debug_exit(string)
Used when a subroutine is exited, write info to log file at some debug level.
OASIS partition data and methods.
Reads the namcouple file for use in OASIS.
Provides a generic and simpler interface into MPI calls for OASIS.
Generic overloaded interface into MPI broadcast.
subroutine, public oasis_map_genmap(mapid, namid)
Routine to generate mapping weights data via a direct SCRIP call.
subroutine, public oasis_mpi_commsize(comm, size, string)
Get the total number of tasks associated with a communicator.
Provides a common location for several OASIS variables.
OASIS map (interpolation) data and methods.
subroutine, public oasis_debug_enter(string)
Used when a subroutine is entered, write info to log file at some debug level.
subroutine, public oasis_flush(nu)
Flushes output to file.
subroutine augment_arrays(cnt, reclen, bsize, nwgts)
Character string manipulation methods.
Performance timer methods.
subroutine, public oasis_abort(id_compid, cd_routine, cd_message)
OASIS abort method, publically available to users.
logical function check_myindex(index, starti, counti)
Function that checks whether an index is part of a start and count list.
subroutine, public oasis_map_smatreaddnc_ceg(sMat, SgsMap, DgsMap, newdom, fileName, mytask, mpicom, nwgts, areasrc, areadst, ni_i, nj_i, ni_o, nj_o)
Read in mapping matrix data from a SCRIP netCDF file using smart scatter (ceg)
subroutine, public oasis_map_smatreaddnc_orig(sMat, SgsMap, DgsMap, newdom, fileName, mytask, mpicom, nwgts, areasrc, areadst, ni_i, nj_i, ni_o, nj_o)
Read in mapping matrix data from a SCRIP netCDF weights file.
subroutine, public oasis_io_read_field_fromroot(filename, fldname, ifld2, fld2, fld3, nx, ny, nz)
Read a field on the root task from a file into an array.
Provides reusable IO routines for OASIS.
integer function get_cegindex(index, starti, counti, peloci)
Mapper data for interpolating data between grids.
OASIS variable data and methods.
Generic overloaded interface into MPI receive.
Defines parameters for OASIS.