!==================================================================== ! subroutine Agrif_u2dbc_interp_tile !==================================================================== ! #include "cppdefs.h" #ifdef AGRIF subroutine u2dbc_interp_tile(Istr,Iend,Jstr,Jend & ,DU_west4,DU_east4,DU_west6,DU_east6 & ,DU_south4,DU_north4,DU_south6,DU_north6) use AGRIF_Util ! implicit none # include "param.h" # include "grid.h" # include "ocean2d.h" # include "scalars.h" # include "climat.h" # include "boundary.h" # include "zoom.h" # include "coupling.h" integer Istr,Iend,Jstr,Jend, i,j real t1,t2,t3,t7,t8,t9,t10,t11,t12,t13,t14,t15,t16,t17,t18,t19 logical ptinterp real t4,t5,t6,tfin,c1,c2,c3,cff1,cff2,cff3 external :: u2dinterp integer :: irhox, irhoy, irhot real :: rrhox, rrhoy, rrhot real :: tinterp, onemtinterp integer :: iter integer :: ipu, jpu, parentnnew integer :: parentnbstep real :: cffx, cffy real,dimension(PRIVATE_1DETA_SCRATCH_ARRAY) :: DU_west, DU_east, & DU_west6,DU_east6 real,dimension(PRIVATE_1DETA_SCRATCH_ARRAY,0:NWEIGHT) :: & DU_west4, DU_east4 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY) :: DU_south,DU_north, & DU_south6,DU_north6 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY,0:NWEIGHT) :: & DU_south4,DU_north4 real :: lastcff !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT # ifdef MPI include 'mpif.h' # endif # ifdef MPI # define LOCALLM Lmmpi # define LOCALMM Mmmpi # else # define LOCALLM Lm # define LOCALMM Mm # endif ! # include "compute_auxiliary_bounds.h" ! return irhox=Agrif_Irhox() irhoy=Agrif_Irhoy() irhot=Agrif_Irhot() rrhox = real(irhox) rrhoy = real(irhoy) rrhot = real(irhot) parentnbstep=Agrif_Parent_Nb_Step() if (U2DTimeindex .NE. parentnbstep) then do J=JstrR,JendR do i=IstrR,IendR dinterp(i,j) = 0. enddo enddo C$OMP BARRIER C$OMP MASTER tinterp=1. # ifdef MASKING Agrif_UseSpecialValue = .true. # endif Agrif_SpecialValue = 0. nbgrid = Agrif_Fixed() C indinterp = 0 indinterp is set to zero in zetabc Call Agrif_Bc_variable(dinterp,DU_avg2,calledweight=tinterp, & procname = u2dinterp) Agrif_UseSpecialvalue=.false. C$OMP END MASTER C$OMP BARRIER lastcff=0. do iter=iif,iif+irhot-1 lastcff=lastcff+((iter+1-iif)/rrhot)* & weight2(iif+irhot-1,iter) enddo lastcff=1./lastcff # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then i = Istr do j=Jstr,Jend DU_west4(j,iif-1) = DU_west3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DU_west6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend DU_west6(j) = DU_west6(j) +cff1*DU_west4(j,iter) enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend DU_west6(j)=DU_west6(j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DU_west4(j,iif-1) enddo enddo do j=Jstr,Jend DU_west6(j)=lastcff*((dinterp(i,j)/rrhoy)-DU_west6(j)) # ifdef MASKING & * umask(i,j) # endif enddo do j=Jstr,Jend DU_west6(j)=(DU_west6(j)-DU_west4(j,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then i=Iend+1 do j=Jstr,Jend DU_east4(j,iif-1) = DU_east3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DU_east6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend DU_east6(j) = DU_east6(j) +cff1*DU_east4(j,iter) enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend DU_east6(j)=DU_east6(j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DU_east4(j,iif-1) enddo enddo do j=Jstr,Jend DU_east6(j)=lastcff*((dinterp(i,j)/rrhoy)-DU_east6(j)) # ifdef MASKING & * umask(i,j) # endif enddo do j=Jstr,Jend DU_east6(j)=(DU_east6(j)-DU_east4(j,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then j=Jstr-1 do i=Istr,Iend DU_south4(i,iif-1) = DU_south3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DU_south6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do i=Istr,Iend DU_south6(i) = DU_south6(i) +cff1*DU_south4(i,iter) enddo enddo do iter=iif,iif+irhot-2 do i=Istr,Iend DU_south6(i)=DU_south6(i)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DU_south4(i,iif-1) enddo enddo do i=Istr,Iend DU_south6(i)=lastcff*((dinterp(i,j)/rrhox)-DU_south6(i)) # ifdef MASKING & * umask(i,j) # endif enddo do i=Istr,Iend DU_south6(i)=(DU_south6(i)-DU_south4(i,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then j=Jend+1 do i=Istr,Iend DU_north4(i,iif-1) = DU_north3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DU_north6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do i=Istr,Iend DU_north6(i) = DU_north6(i) +cff1*DU_north4(i,iter) enddo enddo do iter=iif,iif+irhot-2 do i=Istr,Iend DU_north6(i)=DU_north6(i)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DU_north4(i,iif-1) enddo enddo do i=Istr,Iend DU_north6(i)=lastcff*((dinterp(i,j)/rrhox)-DU_north6(i)) # ifdef MASKING & * umask(i,j) # endif enddo do i=Istr,Iend DU_north6(i)=(DU_north6(i)-DU_north4(i,iif-1))/rrhot enddo endif endif # endif U2DTimeindex = parentnbstep endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then i = Istr do j=Jstr,Jend DU_west4(j,iif-1) = DU_west3(i,j,iif-1) enddo do j=Jstr,Jend DU_west(j)=DU_west4(j,iif-1)+DU_west6(j) enddo endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then i=Iend+1 do j=Jstr,Jend DU_east4(j,iif-1) = DU_east3(i,j,iif-1) enddo do j=Jstr,Jend DU_east(j)=DU_east4(j,iif-1)+DU_east6(j) enddo endif # endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then j=Jstr-1 do i=Istr,Iend DU_south4(i,iif-1) = DU_south3(i,j,iif-1) enddo do i=Istr,Iend DU_south(i)=DU_south4(i,iif-1)+DU_south6(i) enddo endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then j=Jend+1 do i=Istr,Iend DU_north4(i,iif-1) = DU_north3(i,j,iif-1) enddo do i=Istr,Iend DU_north(i)=DU_north4(i,iif-1)+DU_north6(i) enddo endif # endif # if defined M2_FRC_BRY || defined M2NUDGING ! ! Apply the value to ubclm or ubarbry ! cffx = g*dtfast*2./(1.+rrhox) # ifdef AGRIF_2WAY cffx = 0. # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then do j=Jstr,Jend # ifdef M2_FRC_BRY ubarbry_west(j)= # else ubclm(Istr,j)= (cffx/om_u(Istr,j))* & (SSH(Istr,j)-zeta(Istr,j,knew)) + # endif # ifdef AGRIF_FLUX_BC & (2.*DU_west(j)/((h(Istr-1,j)+zeta(Istr-1,j,knew) & +h(Istr,j)+zeta(Istr,j,knew)) & *on_u(Istr,j))) # else & DU_west(j) # endif # ifdef MASKING & *umask(Istr,j) # endif enddo endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then do j=Jstr,Jend # ifdef M2_FRC_BRY ubarbry_east(j)= # else ubclm(Iend+1,j)=-(cffx/om_u(Iend+1,j))* & (SSH(Iend,j)-zeta(Iend,j,knew)) + # endif # ifdef AGRIF_FLUX_BC & (2.*DU_east(j)/(( h(Iend,j)+zeta(Iend,j,knew) & +h(Iend+1,j)+zeta(Iend+1,j,knew)) & *on_u(Iend+1,j))) # else & DU_east(j) # endif # ifdef MASKING & *umask(Iend+1,j) # endif enddo endif # endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then do i=Istr,Iend # ifdef M2_FRC_BRY ubarbry_south(i)= # else ubclm(i,Jstr-1)= # endif # ifdef AGRIF_FLUX_BC & (2.*DU_south(i)/(( h(i,Jstr-1)+zeta(i,Jstr-1,knew) & +h(i-1,Jstr-1)+zeta(i-1,Jstr-1,knew)) & *on_u(i,Jstr-1))) # else & DU_south(i) # endif # ifdef MASKING & *umask(i,Jstr-1) # endif enddo endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then do i=Istr,Iend # ifdef M2_FRC_BRY ubarbry_north(i)= # else ubclm(i,Jend+1)= # endif # ifdef AGRIF_FLUX_BC & (2.*DU_north(i)/(( h(i,Jend+1)+zeta(i,Jend+1,knew) & +h(i-1,Jend+1)+zeta(i-1,Jend+1,knew)) & *on_u(i,Jend+1))) # else & DU_north(i) # endif # ifdef MASKING & *umask(i,Jend+1) # endif enddo endif # endif # endif /* M2_FRC_BRY || M2NUDGING */ return end subroutine u2Dinterp(tabres,i1,i2,j1,j2) implicit none # include "param.h" # include "grid.h" # include "ocean2d.h" # include "scalars.h" integer i1,i2,j1,j2 real tabres(i1:i2,j1:j2) integer i,j,iter, isize real :: cff1 !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT integer :: oldindinterp isize = (j2-j1+1)*(i2-i1+1) IF (iif == 1) THEN IF (.NOT.Associated(coarsevalues)) THEN Allocate(coarsevalues(isize,0:nfast,1)) ELSE CALL checksizeinterp(indinterp+isize,nfast) ENDIF oldindinterp = indinterp do j=j1,j2 do i=max(i1,lbound(h,1)+1),i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,0,nbgrid) = & 0.5*(h(i-1,j)+zeta(i-1,j,kstp)+h(i,j)+ & zeta(i,j,kstp))*on_u(i,j)*ubar(i,j,kstp) enddo enddo ENDIF oldindinterp = indinterp do j=j1,j2 do i=max(i1,lbound(h,1)+1),i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,iif,nbgrid) = & 0.5*(h(i-1,j)+zeta(i-1,j,knew)+h(i,j)+ & zeta(i,j,knew))*on_u(i,j)*ubar(i,j,knew) enddo enddo do iter=0,iif cff1 = weight2(iif,iter) call copy1d(tabres,coarsevalues(indinterp+1,iter,nbgrid), & cff1,isize) enddo indinterp = oldindinterp return end ! !==================================================================== ! subroutine Agrif_v2dbc_interp_tile !==================================================================== ! subroutine v2dbc_interp_tile(Istr,Iend,Jstr,Jend & ,DV_west4,DV_east4,DV_west6,DV_east6 & ,DV_south4,DV_north4,DV_south6,DV_north6) use AGRIF_Util ! implicit none # include "param.h" # include "grid.h" # include "ocean2d.h" # include "scalars.h" # include "climat.h" # include "boundary.h" # include "zoom.h" # include "coupling.h" integer Istr,Iend,Jstr,Jend, i,j real t1,t2,t3,t4,t5,t6,dv1,dv1np1,dv2,tfin,c1,c2,c3 real t7,t8,t9,t10,t11 external :: v2dinterp integer :: irhox, irhoy, irhot real :: rrhox, rrhoy, rrhot real :: tinterp, onemtinterp integer :: iter integer :: parentnbstep real :: cffy real,dimension(PRIVATE_1DETA_SCRATCH_ARRAY) :: DV_west, DV_east, & DV_west6,DV_east6 real,dimension(PRIVATE_1DETA_SCRATCH_ARRAY,0:NWEIGHT) :: & DV_west4, DV_east4 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY) :: DV_south,DV_north, & DV_south6,DV_north6 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY,0:NWEIGHT) :: & DV_south4,DV_north4 real :: lastcff !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT #ifdef MPI include 'mpif.h' #endif # ifdef MPI # define LOCALLM Lmmpi # define LOCALMM Mmmpi # else # define LOCALLM Lm # define LOCALMM Mm # endif ! # include "compute_auxiliary_bounds.h" ! ! return irhox=Agrif_Irhox() irhoy=Agrif_Irhoy() irhot=Agrif_Irhot() rrhox = real(irhox) rrhoy = real(irhoy) rrhot = real(irhot) parentnbstep=Agrif_Parent_Nb_Step() if (V2DTimeindex .NE. parentnbstep) then do J=JstrR,JendR do i=IstrR,IendR dinterp(i,j) = 0. enddo enddo C$OMP BARRIER C$OMP MASTER tinterp=1. #ifdef MASKING Agrif_UseSpecialValue = .true. #endif Agrif_SpecialValue = 0. nbgrid = Agrif_Fixed() C indinterp = 0 indinterp is set to zero in zetabc Call Agrif_Bc_variable(dinterp,DV_avg2,calledweight=tinterp, & procname = v2dinterp) Agrif_UseSpecialvalue=.false. C$OMP END MASTER C$OMP BARRIER lastcff=0. do iter=iif,iif+irhot-1 lastcff=lastcff+((iter+1-iif)/rrhot)* & weight2(iif+irhot-1,iter) enddo lastcff=1./lastcff # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then j=Jstr do i=Istr,Iend DV_south4(i,iif-1) = DV_south3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DV_south6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do i=Istr,Iend DV_south6(i) = DV_south6(i) +cff1*DV_south4(i,iter) enddo enddo do iter=iif,iif+irhot-2 do i=Istr,Iend DV_south6(i)=DV_south6(i)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DV_south4(i,iif-1) enddo enddo do i=Istr,Iend DV_south6(i)=lastcff*((dinterp(i,j)/rrhox)-DV_south6(i)) # ifdef MASKING & * vmask(i,j) # endif enddo do i=Istr,Iend DV_south6(i)=(DV_south6(i)-DV_south4(i,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then j=Jend+1 do i=Istr,Iend DV_north4(i,iif-1) = DV_north3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DV_north6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do i=Istr,Iend DV_north6(i) = DV_north6(i) +cff1*DV_north4(i,iter) enddo enddo do iter=iif,iif+irhot-2 do i=Istr,Iend DV_north6(i)=DV_north6(i)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DV_north4(i,iif-1) enddo enddo do i=Istr,Iend DV_north6(i)=lastcff*((dinterp(i,j)/rrhox)-DV_north6(i)) # ifdef MASKING & * vmask(i,j) # endif enddo do i=Istr,Iend DV_north6(i)=(DV_north6(i)-DV_north4(i,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then i = Istr-1 do j=Jstr,Jend DV_west4(j,iif-1) = DV_west3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DV_west6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend DV_west6(j) = DV_west6(j) +cff1*DV_west4(j,iter) enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend DV_west6(j)=DV_west6(j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DV_west4(j,iif-1) enddo enddo do j=Jstr,Jend DV_west6(j)=lastcff*((dinterp(i,j)/rrhoy)-DV_west6(j)) # ifdef MASKING & * vmask(i,j) # endif enddo do j=Jstr,Jend DV_west6(j)=(DV_west6(j)-DV_west4(j,iif-1))/rrhot enddo endif endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then i=Iend+1 do j=Jstr,Jend DV_east4(j,iif-1) = DV_east3(i,j,iif-1) enddo if (iif+irhot-1<=nfast) then DV_east6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend DV_east6(j) = DV_east6(j) +cff1*DV_east4(j,iter) enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend DV_east6(j)=DV_east6(j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*DV_east4(j,iif-1) enddo enddo do j=Jstr,Jend DV_east6(j)=lastcff*((dinterp(i,j)/rrhoy)-DV_east6(j)) # ifdef MASKING & * vmask(i,j) # endif enddo do j=Jstr,Jend DV_east6(j)=(DV_east6(j)-DV_east4(j,iif-1))/rrhot enddo endif endif # endif V2DTimeindex = parentnbstep endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then j=Jstr do i=Istr,Iend DV_south4(i,iif-1) = DV_south3(i,j,iif-1) enddo do i=Istr,Iend DV_south(i)=DV_south4(i,iif-1)+DV_south6(i) enddo endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then j=Jend+1 do i=Istr,Iend DV_north4(i,iif-1) = DV_north3(i,j,iif-1) enddo do i=Istr,Iend DV_north(i)=DV_north4(i,iif-1)+DV_north6(i) enddo endif # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then i = Istr-1 do j=Jstr,Jend DV_west4(j,iif-1) = DV_west3(i,j,iif-1) enddo do j=Jstr,Jend DV_west(j)=DV_west4(j,iif-1)+DV_west6(j) enddo endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then i=Iend+1 do j=Jstr,Jend DV_east4(j,iif-1) = DV_east3(i,j,iif-1) enddo do j=Jstr,Jend DV_east(j)=DV_east4(j,iif-1)+DV_east6(j) enddo endif # endif # if defined M2_FRC_BRY || defined M2NUDGING ! ! Apply the value to vbclm or vbarbry ! cffy = g*dtfast*2./(1.+rrhoy) # ifdef AGRIF_2WAY cffy = 0. # endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then do i=Istr,Iend # ifdef M2_FRC_BRY vbarbry_south(i)= # else vbclm(i,Jstr)=(cffy/on_v(i,Jstr))* & (SSH(i,Jstr)-zeta(i,Jstr,knew)) + # endif # ifdef AGRIF_FLUX_BC & (2.*DV_south(i)/((h(i,Jstr-1)+zeta(i,Jstr-1,knew) & +h(i,Jstr)+zeta(i,Jstr,knew)) & *om_v(i,Jstr))) # else & DV_south(i) # endif # ifdef MASKING & *vmask(i,Jstr) # endif enddo endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then do i=Istr,Iend # ifdef M2_FRC_BRY vbarbry_north(i)= # else vbclm(i,Jend+1)=-(cffy/on_v(i,Jend+1))* & (SSH(i,Jend)-zeta(i,Jend,knew)) + # endif # ifdef AGRIF_FLUX_BC & (2.*DV_north(i)/(( h(i,Jend)+zeta(i,Jend,knew) & +h(i,Jend+1)+zeta(i,Jend+1,knew)) & *om_v(i,Jend+1))) # else & DV_north(i) # endif # ifdef MASKING & *vmask(i,Jend+1) # endif enddo endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then do j=Jstr,Jend # ifdef M2_FRC_BRY vbarbry_east(j)= # else vbclm(Iend+1,j)= # endif # ifdef AGRIF_FLUX_BC & (2.*DV_east(j)/(( h(Iend+1,j)+zeta(Iend+1,j,knew) & +h(Iend+1,j-1)+zeta(Iend+1,j-1,knew)) & *om_v(Iend+1,j))) # else & DV_east(j) # endif # ifdef MASKING & *vmask(Iend+1,j) # endif enddo endif # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then do j=Jstr,Jend # ifdef M2_FRC_BRY vbarbry_west(j)= # else vbclm(Istr-1,j)= # endif # ifdef AGRIF_FLUX_BC & (2.*DV_west(j)/((h(Istr-1,j)+zeta(Istr-1,j,knew) & +h(Istr-1,j-1)+zeta(Istr-1,j-1,knew)) & *om_v(Istr-1,j))) # else & DV_west(j) # endif # ifdef MASKING & *vmask(Istr-1,j) # endif enddo endif # endif # endif /* M2_FRC_BRY || M2NUDGING */ return end subroutine v2Dinterp(tabres,i1,i2,j1,j2) implicit none # include "param.h" # include "grid.h" # include "ocean2d.h" # include "scalars.h" integer i1,i2,j1,j2 real tabres(i1:i2,j1:j2) integer i,j,iter, isize real :: cff1 !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT integer :: oldindinterp isize = (j2-j1+1)*(i2-i1+1) IF (iif == 1) THEN IF (.NOT.Associated(coarsevalues)) THEN Allocate(coarsevalues(isize,0:nfast,1)) ELSE CALL checksizeinterp(indinterp+isize,nfast) ENDIF oldindinterp = indinterp do j=max(j1,lbound(h,2)+1),j2 do i=i1,i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,0,nbgrid) = & 0.5*(h(i,j-1)+zeta(i,j-1,kstp)+h(i,j)+ & zeta(i,j,kstp))*om_v(i,j)*vbar(i,j,kstp) enddo enddo ENDIF oldindinterp = indinterp do j=max(j1,lbound(h,2)+1),j2 do i=i1,i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,iif,nbgrid) = & 0.5*(h(i,j-1)+zeta(i,j-1,knew)+h(i,j)+ & zeta(i,j,knew))*om_v(i,j)*vbar(i,j,knew) enddo enddo do iter=0,iif cff1 = weight2(iif,iter) call copy1d(tabres,coarsevalues(indinterp+1,iter,nbgrid), & cff1,isize) enddo indinterp = oldindinterp return end subroutine copy1d(tab1,tab2,cff,isize) real,dimension(isize) :: tab1,tab2 integer :: i do i=1,isize tab1(i)=tab1(i)+cff*tab2(i) enddo return end subroutine copy1d ! subroutine checksizeinterp(isize,nfast) integer :: isize,nfast real,dimension(:,:,:),allocatable :: tempvalues integer :: n1,n3 !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT IF (size(coarsevalues,1).LT.(isize)) THEN n1 = size(coarsevalues,1) n3 = size(coarsevalues,3) allocate(tempvalues(n1,0:nfast,n3)) tempvalues=coarsevalues(1:n1,0:nfast,1:n3) deallocate(coarsevalues) allocate(coarsevalues(isize,0:nfast,n3)) coarsevalues(1:n1,0:nfast,1:n3) = tempvalues deallocate(tempvalues) ELSE IF (nbgrid.GT.size(coarsevalues,3)) THEN n1 = size(coarsevalues,1) n3 = size(coarsevalues,3) allocate(tempvalues(n1,0:nfast,n3)) tempvalues=coarsevalues(1:n1,0:nfast,1:n3) deallocate(coarsevalues) allocate(coarsevalues(n1,0:nfast,nbgrid)) coarsevalues(1:n1,0:nfast,1:n3) = tempvalues deallocate(tempvalues) ENDIF return end ! !==================================================================== ! subroutine Agrif_zetabc_interp_tile !==================================================================== ! subroutine zetabc_interp_tile(Istr,Iend,Jstr,Jend & ,Zeta_west4,Zeta_east4,Zeta_west6,Zeta_east6 & ,Zeta_south4,Zeta_north4,Zeta_south6,Zeta_north6) use AGRIF_Util ! implicit none # include "param.h" # include "boundary.h" # include "climat.h" # include "grid.h" # include "ocean2d.h" # include "scalars.h" # include "zoom.h" integer Istr,Iend,Jstr,Jend, i,j, i1, j1 real tinterp Integer itrcind INTEGER :: parentknew, parentkstp,nbstep3dparent real t1,t2,t3,t4,t5,t6,t9,t10,t11,t7 real tin(2), tout(2) real cff1 real zeta1,zeta2 external zetainterp integer :: irhot, irhox, irhoy real :: rrhot integer :: iter real :: onemtinterp integer :: parentnbstep real,dimension(Istr-1:Istr,PRIVATE_1DETA_SCRATCH_ARRAY) :: & Zeta_west6 real,dimension(Iend:Iend+1,PRIVATE_1DETA_SCRATCH_ARRAY) :: & Zeta_east6 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY,Jstr-1:Jstr) :: & Zeta_south6 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY,Jend:Jend+1) :: & Zeta_north6 real,dimension(Istr-1:Istr, & PRIVATE_1DETA_SCRATCH_ARRAY,0:NWEIGHT) :: Zeta_west4 real,dimension(Iend:Iend+1, & PRIVATE_1DETA_SCRATCH_ARRAY,0:NWEIGHT) :: Zeta_east4 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY, & Jstr-1:Jstr,0:NWEIGHT) :: Zeta_south4 real,dimension(PRIVATE_1DXI_SCRATCH_ARRAY, & Jend:Jend+1,0:NWEIGHT) :: Zeta_north4 real lastcff !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT # ifdef MPI # define LOCALLM Lmmpi # define LOCALMM Mmmpi # else # define LOCALLM Lm # define LOCALMM Mm # endif ! # include "compute_auxiliary_bounds.h" ! ! return irhot = Agrif_Irhot() irhox = Agrif_Irhox() irhoy = Agrif_Irhoy() rrhot = real(irhot) parentnbstep=Agrif_Parent_Nb_Step() if (ZetaTimeindex .NE. parentnbstep) then do J=JstrR,JendR do i=IstrR,IendR dinterp(i,j) = 0. enddo enddo C$OMP BARRIER C$OMP MASTER tinterp=1. #ifdef MASKING Agrif_UseSpecialvalue=.true. #endif Agrif_Specialvalue=0. nbgrid = Agrif_Fixed() indinterp = 0 Call Agrif_Bc_variable(dinterp,Zt_avg1, & calledweight=tinterp, & procname=zetainterp) Agrif_UseSpecialvalue=.false. C$OMP END MASTER C$OMP BARRIER lastcff=0. do iter=iif,iif+irhot-1 lastcff=lastcff+((iter+1-iif)/rrhot)* & weight2(iif+irhot-1,iter) enddo lastcff=1./lastcff # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south4(i,j,iif-1) = Zeta_south3(i,j,iif-1) enddo enddo if (iif+irhot-1<=nfast) then Zeta_south6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south6(i,j) = Zeta_south6(i,j) & +cff1*Zeta_south4(i,j,iter) enddo enddo enddo do iter=iif,iif+irhot-2 do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south6(i,j)=Zeta_south6(i,j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*Zeta_south4(i,j,iif-1) enddo enddo enddo do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south6(i,j)=lastcff*(dinterp(i,j)-Zeta_south6(i,j)) enddo enddo do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south6(i,j)=(Zeta_south6(i,j)-Zeta_south4(i,j,iif-1)) & /rrhot enddo enddo endif endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then do j=Jend,Jend+1 do i=Istr,Iend Zeta_north4(i,j,iif-1) = Zeta_north3(i,j,iif-1) enddo enddo if (iif+irhot-1<=nfast) then Zeta_north6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jend,Jend+1 do i=Istr,Iend Zeta_north6(i,j) = Zeta_north6(i,j) & +cff1*Zeta_north4(i,j,iter) enddo enddo enddo do iter=iif,iif+irhot-2 do j=Jend,Jend+1 do i=Istr,Iend Zeta_north6(i,j)=Zeta_north6(i,j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*Zeta_north4(i,j,iif-1) enddo enddo enddo do j=Jend,Jend+1 do i=Istr,Iend Zeta_north6(i,j)=lastcff*(dinterp(i,j)-Zeta_north6(i,j)) enddo enddo do j=Jend,Jend+1 do i=Istr,Iend Zeta_north6(i,j)=(Zeta_north6(i,j)-Zeta_north4(i,j,iif-1)) & /rrhot enddo enddo endif endif # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then do j=Jstr,Jend do i=Istr-1,Istr Zeta_west4(i,j,iif-1) = Zeta_west3(i,j,iif-1) enddo enddo if (iif+irhot-1<=nfast) then Zeta_west6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend do i=Istr-1,Istr Zeta_west6(i,j) = Zeta_west6(i,j) & +cff1*Zeta_west4(i,j,iter) enddo enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend do i=Istr-1,Istr Zeta_west6(i,j)=Zeta_west6(i,j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*Zeta_west4(i,j,iif-1) enddo enddo enddo do j=Jstr,Jend do i=Istr-1,Istr Zeta_west6(i,j)=lastcff*(dinterp(i,j)-Zeta_west6(i,j)) enddo enddo do j=Jstr,Jend do i=Istr-1,Istr Zeta_west6(i,j)=(Zeta_west6(i,j)-Zeta_west4(i,j,iif-1)) & /rrhot enddo enddo endif endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east4(i,j,iif-1) = Zeta_east3(i,j,iif-1) enddo enddo if (iif+irhot-1<=nfast) then Zeta_east6=0. do iter=0,iif-1 cff1=weight2(iif+irhot-1,iter) do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east6(i,j) = Zeta_east6(i,j) & +cff1*Zeta_east4(i,j,iter) enddo enddo enddo do iter=iif,iif+irhot-2 do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east6(i,j)=Zeta_east6(i,j)+((iif+irhot-1-iter)/rrhot)* & weight2(iif+irhot-1,iter)*Zeta_east4(i,j,iif-1) enddo enddo enddo do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east6(i,j)=lastcff*(dinterp(i,j)-Zeta_east6(i,j)) enddo enddo do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east6(i,j)=(Zeta_east6(i,j)-Zeta_east4(i,j,iif-1)) & /rrhot enddo enddo endif endif # endif ZetaTimeindex = parentnbstep endif # ifdef AGRIF_OBC_SOUTH if (SOUTHERN_EDGE) then do j=Jstr-1,Jstr do i=Istr,Iend Zeta_south4(i,j,iif-1) = Zeta_south3(i,j,iif-1) enddo enddo # ifdef ZNUDGING do j=Jstr-1,Jstr do i=Istr,Iend SSH(i,j)=Zeta_south4(i,j,iif-1)+Zeta_south6(i,j) enddo enddo # endif endif # endif # ifdef AGRIF_OBC_NORTH if (NORTHERN_EDGE) then do j=Jend,Jend+1 do i=Istr,Iend Zeta_north4(i,j,iif-1) = Zeta_north3(i,j,iif-1) enddo enddo # ifdef ZNUDGING do j=Jend,Jend+1 do i=Istr,Iend SSH(i,j)=Zeta_north4(i,j,iif-1)+Zeta_north6(i,j) enddo enddo # endif endif # endif # ifdef AGRIF_OBC_WEST if (WESTERN_EDGE) then do j=Jstr,Jend do i=Istr-1,Istr Zeta_west4(i,j,iif-1) = Zeta_west3(i,j,iif-1) enddo enddo # ifdef ZNUDGING do j=Jstr,Jend do i=Istr-1,Istr SSH(i,j)=Zeta_west4(i,j,iif-1)+Zeta_west6(i,j) enddo enddo # endif endif # endif # ifdef AGRIF_OBC_EAST if (EASTERN_EDGE) then do j=Jstr,Jend do i=Iend,Iend+1 Zeta_east4(i,j,iif-1) = Zeta_east3(i,j,iif-1) enddo enddo # ifdef ZNUDGING do j=Jstr,Jend do i=Iend,Iend+1 SSH(i,j)=Zeta_east4(i,j,iif-1)+Zeta_east6(i,j) enddo enddo # endif endif # endif return end subroutine zetainterp(tabres,i1,i2,j1,j2) implicit none # include "param.h" # include "ocean2d.h" # include "scalars.h" integer i1,i2,j1,j2 real tabres(i1:i2,j1:j2) integer i,j,iter, isize real :: cff1 !$AGRIF_DO_NOT_TREAT real,dimension(:,:,:),pointer :: coarsevalues integer :: nbgrid, indinterp common/interp2d/coarsevalues, nbgrid, indinterp !$AGRIF_END_DO_NOT_TREAT integer :: oldindinterp isize = (j2-j1+1)*(i2-i1+1) IF (iif == 1) THEN IF (.NOT.Associated(coarsevalues)) THEN Allocate(coarsevalues(isize,0:nfast,1)) ELSE CALL checksizeinterp(indinterp+isize,nfast) ENDIF oldindinterp = indinterp do j=j1,j2 do i=i1,i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,0,nbgrid) = & zeta(i,j,kstp) enddo enddo ENDIF oldindinterp = indinterp do j=j1,j2 do i=i1,i2 oldindinterp=oldindinterp+1 coarsevalues(oldindinterp,iif,nbgrid) = & zeta(i,j,knew) enddo enddo do iter=0,iif cff1 = weight2(iif,iter) call copy1d(tabres,coarsevalues(indinterp+1,iter,nbgrid), & cff1,isize) enddo indinterp = oldindinterp return end #else subroutine zommbc_2D_empty() return end #endif