MODULE obcdta !!============================================================================== !! *** MODULE obcdta *** !! Open boundary data : read the data for the open boundaries. !!============================================================================== #if defined key_obc !!------------------------------------------------------------------------------ !! 'key_obc' : Open Boundary Conditions !!------------------------------------------------------------------------------ !! obc_dta : read u, v, t, s data along each open boundary !! obc_dta_psi : read psi data along each open boundary (rigid lid only) !!------------------------------------------------------------------------------ !! * Modules used USE oce ! ocean dynamics and tracers USE dom_oce ! ocean space and time domain USE lbclnk ! ocean lateral boundary conditions (or mpp link) USE phycst ! physical constants USE obc_oce ! ocean open boundary conditions USE daymod ! calendar USE in_out_manager ! I/O logical units USE lib_mpp ! distribued memory computing USE dynspg_rl ! # if ! defined key_dynspg_fsc USE obccli # endif IMPLICIT NONE PRIVATE !! * Accessibility PUBLIC obc_dta ! routines called by step.F90 !! * Substitutions # include "obc_vectopt_loop_substitute.h90" !!--------------------------------------------------------------------------------- !! OPA 9.0 , LOCEAN-IPSL (2005) !! $Header$ !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt !!--------------------------------------------------------------------------------- CONTAINS SUBROUTINE obc_dta( kt ) !!--------------------------------------------------------------------------- !! *** SUBROUTINE obc_dta *** !! !! ** Purpose : Find the climatological boundary arrays for the specified date, !! Originally this routine interpolated between monthly fields !! of a climatology. !! When using hydrological section data, we have only on snapshot !! and do not need to interpolate. !! !! ** Method : Determine the current month from kdat, and interpolate for !! the current day. !! !! History : !! ! 98-05 (J.M. Molines) Original code !! 8.5 ! 02-10 (C. Talandier, A-M. Treguier) Free surface, F90 !!--------------------------------------------------------------------------- !! * Arguments INTEGER, INTENT( in ) :: kt ! ocean time-step index !! * Local declarations INTEGER :: ji, jj, jk, ii, ij ! dummy loop indices INTEGER :: inum = 11 ! temporary logical unit INTEGER :: imois, iman, imoisp1 INTEGER :: i15 INTEGER :: irecl, irec, ios INTEGER, PARAMETER :: & jpmois = 1, & jpf = 3 REAL(wp) :: zxy CHARACTER (len=4 ) :: clversion CHARACTER (len=80) :: clcom !!--------------------------------------------------------------------- IF( lk_dynspg_rl ) CALL obc_dta_psi( kt ) ! update bsf data at open boundaries ! 0. Initialization of date ! imois is the index (1 to 12) of the first month to be used in the ! interpolation. ndastp is the date (integer format yymmdd) calculated ! in routine day.F starting from ndate0 given in namelist. ! ----------------------------------------------------------------------- iman = jpmois i15 = nday/16 imois = nmonth + i15 - 1 IF( imois == 0 ) imois = iman imoisp1 = MOD( imois + 1, iman ) IF( imoisp1 == 0 ) imoisp1 = iman ! ... zxy is the weight of the after field zxy = FLOAT( nday + 15 - 30 * i15 ) / 30. IF( zxy > 1.01 .OR. zxy < 0. ) THEN IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)'obc_dta: Pbm with the the weight of the after field zxy ' IF(lwp) WRITE(numout,*)'~~~~~~~~~~~' nstop = nstop + 1 END IF ! 1. First call: ! ---------------- IF( kt == nit000 ) THEN IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)'obcdta: initial step in obc_dta' IF(lwp) WRITE(numout,*)'~~~~~~ months ',imois,' and', imoisp1,' read' ! 1.1 Tangential velocities set to zero ! -------------------------------------- IF( lp_obc_east ) vfoe(:,1:jpkm1) = 0.e0 IF( lp_obc_west ) vfow(:,1:jpkm1) = 0.e0 IF( lp_obc_south ) ufos(:,1:jpkm1) = 0.e0 IF( lp_obc_north ) ufon(:,1:jpkm1) = 0.e0 ! 1.2 Set the Normal velocity and tracers data for the EAST OBC ! ------------------------------------------------------------- IF( lp_obc_east ) THEN ! initialisation to zero sedta(:,:,1) = 0.e0 tedta(:,:,1) = 0.e0 uedta(:,:,1) = 0.e0 ! ! ================== ! IF( nobc_dta == 0 ) THEN ! initial state used ! ! ================== ! DO ji = fs_nie0, fs_nie1 ! Vector opt. DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp sedta(ij,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk) tedta(ij,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk) uedta(ij,jk,1) = un(ji,jj,jk)*umask(ji,jj,jk) END DO END DO END DO DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp sfoe(jj,jk) = sedta(ij,jk,1)*temsk(jj,jk) tfoe(jj,jk) = tedta(ij,jk,1)*temsk(jj,jk) ufoe(jj,jk) = uedta(ij,jk,1)*uemsk(jj,jk) END DO END DO ! ! =================== ! ELSE ! read in obceast.dta ! ! =================== ! OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obceast.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = 4096 ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obceast.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF READ(inum,REC=1) clversion,clcom,irecl CLOSE(inum) IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)' opening obceast.dta with irecl=',irecl OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obceast.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = irecl ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obceast.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF ! ... Read datafile and set temperature, salinity and normal velocity ! ... initialise the sedta, tedta arrays ! ... index 1 refer to before, 2 to after DO jk = 1, jpkm1 irec = 2 + (jk -1)* jpf READ(inum,REC=irec )((sedta(jj,jk,1),ji=1,1),jj=jpjed, jpjef) READ(inum,REC=irec+1)((tedta(jj,jk,1),ji=1,1),jj=jpjed, jpjef) READ(inum,REC=irec+2)((uedta(jj,jk,1),ji=1,1),jj=jpjed, jpjef) ! ... set climatological temperature and salinity at the boundary ! ... do not interpolate between months : impose values of climatology DO jj = 1, jpj ij = jj -1 + njmpp sfoe(jj,jk) = sedta(ij,jk,1)*temsk(jj,jk) tfoe(jj,jk) = tedta(ij,jk,1)*temsk(jj,jk) END DO END DO CLOSE(inum) # if ! defined key_dynspg_fsc ! ... Rigid lid case: make sure uedta is baroclinic velocity ! ... In rigid lid case uedta needs to be the baroclinic component. CALL obc_cli( uedta, uclie, fs_nie0, fs_nie1, 0, jpj, njmpp ) # endif ! ... Set normal velocity (on nie0, nie1 <=> jpieob) DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp ufoe(jj,jk) = uedta(ij,jk,1)*uemsk(jj,jk) END DO END DO ENDIF ENDIF ! 1.3 Set the Normal velocity and tracers data for the WEST OBC ! ------------------------------------------------------------- IF( lp_obc_west ) THEN ! initialisation to zero swdta(:,:,1) = 0.e0 twdta(:,:,1) = 0.e0 uwdta(:,:,1) = 0.e0 ! ! ================== ! IF( nobc_dta == 0 ) THEN ! initial state used ! ! ================== ! DO ji = fs_niw0, fs_niw1 ! Vector opt. DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp swdta(ij,jk,1) = sn(ji,jj,jk)*tmask(ji,jj,jk) twdta(ij,jk,1) = tn(ji,jj,jk)*tmask(ji,jj,jk) uwdta(ij,jk,1) = un(ji,jj,jk)*umask(ji,jj,jk) END DO END DO END DO DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp sfow(jj,jk) = swdta(ij,jk,1)*twmsk(jj,jk) tfow(jj,jk) = twdta(ij,jk,1)*twmsk(jj,jk) ufow(jj,jk) = uwdta(ij,jk,1)*uwmsk(jj,jk) END DO END DO ! ! =================== ! ELSE ! read in obceast.dta ! ! =================== ! OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcwest.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = 4096 ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcwest.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF READ(inum,REC=1) clversion, clcom,irecl CLOSE(inum) IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)' opening obcwest.dta with irecl=',irecl OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcwest.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = irecl ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcwest.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF ! ... Read datafile and set temperature, salinity and normal velocity ! ... initialise the swdta, twdta arrays ! ... index 1 refer to before, 2 to after DO jk = 1, jpkm1 irec = 2 + (jk -1)* jpf READ(inum,REC=irec )((swdta(jj,jk,1),ji=1,1),jj=jpjwd, jpjwf) READ(inum,REC=irec+1)((twdta(jj,jk,1),ji=1,1),jj=jpjwd, jpjwf) READ(inum,REC=irec+2)((uwdta(jj,jk,1),ji=1,1),jj=jpjwd, jpjwf) DO jj = 1, jpj ij = jj -1 + njmpp sfow(jj,jk) = swdta(ij,jk,1)*twmsk(jj,jk) tfow(jj,jk) = twdta(ij,jk,1)*twmsk(jj,jk) END DO END DO CLOSE(inum) #if ! defined key_dynspg_fsc ! ... Rigid lid case: make sure uwdta is baroclinic velocity ! ... In rigid lid case uwdta needs to be the baroclinic component. CALL obc_cli( uwdta, ucliw, fs_niw0, fs_niw1, 0, jpj, njmpp ) # endif ! ... Set normal velocity (on niw0, niw1 <=> jpiwob) DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp ufow(jj,jk) = uwdta(ij,jk,1)*uwmsk(jj,jk) END DO END DO ENDIF ENDIF ! 1.4 Set the Normal velocity and tracers data for the NORTH OBC ! --------------------------------------------------------------- IF( lp_obc_north ) THEN ! initialisation to zero sndta(:,:,1) = 0.e0 tndta(:,:,1) = 0.e0 vndta(:,:,1) = 0.e0 OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcnorth.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = 4096 ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcnorth.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF READ(inum,REC=1) clversion,clcom,irecl CLOSE(inum) IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)' opening obcnorth.dta with irecl=',irecl OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcnorth.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = irecl ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcnorth.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF ! ... Read datafile and set temperature and salinity ! ... initialise the sndta, tndta arrays ! ... index 1 refer to before, 2 to after DO jk = 1, jpkm1 irec = 2 + (jk -1)* jpf READ(inum,REC=irec )((sndta(jj,jk,1),ji=1,1),jj=jpind, jpinf) READ(inum,REC=irec+1)((tndta(jj,jk,1),ji=1,1),jj=jpind, jpinf) READ(inum,REC=irec+2)((vndta(jj,jk,1),ji=1,1),jj=jpind, jpinf) ! ... do not interpolate: impose values of climatology DO ji = 1, jpi ii = ji -1 + nimpp sfon(ji,jk) = sndta(ii,jk,1)*tnmsk(ji,jk) tfon(ji,jk) = tndta(ii,jk,1)*tnmsk(ji,jk) END DO END DO CLOSE(inum) # if ! defined key_dynspg_fsc ! ... Rigid lid case: make sure vndta is baroclinic velocity ! In rigid lid case vndta needs to be the baroclinic component. ! substract the barotropic velocity component (zvnbtpe) ! from the total one (vndta). CALL obc_cli( vndta, vclin, fs_njn0, fs_njn1, 1, jpi, nimpp ) # endif ! ... Set normal velocity (on njn0, njn1 <=> jpjnob) DO jk = 1, jpkm1 DO ji = 1, jpi ii = ji -1 + nimpp vfon(ji,jk) = vndta(ii,jk,1)*vnmsk(ji,jk) END DO END DO END IF ! 1.5 Set the Normal velocity and tracers data for the SOUTH OBC ! --------------------------------------------------------------- IF( lp_obc_south ) THEN ! initialisation to zero ssdta(:,:,1) = 0.e0 tsdta(:,:,1) = 0.e0 vsdta(:,:,1) = 0.e0 OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcsouth.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = 4096 ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcsouth.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF READ(inum,REC=1) clversion,clcom,irecl CLOSE(inum) IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)' opening obcsouth.dta with irecl=',irecl OPEN(UNIT = inum, & IOSTAT = ios, & FILE ='obcsouth.dta', & FORM ='UNFORMATTED', & ACCESS ='DIRECT', & RECL = irecl ) IF( ios > 0 ) THEN IF(lwp) WRITE(numout,*) 'obc_dta: Pbm to OPEN the obcsouth.dta file ' IF(lwp) WRITE(numout,*) '~~~~~~~' nstop = nstop + 1 END IF ! ... Read datafile and set temperature, salinity and normal velocity ! ... initialise the ssdta, tsdta arrays ! ... index 1 refer to before, 2 to after DO jk = 1, jpkm1 irec = 2 + (jk -1)* jpf READ(inum,REC=irec )((ssdta(jj,jk,1),ji=1,1),jj=jpisd, jpisf) READ(inum,REC=irec+1)((tsdta(jj,jk,1),ji=1,1),jj=jpisd, jpisf) READ(inum,REC=irec+2)((vsdta(jj,jk,1),ji=1,1),jj=jpisd, jpisf) ! ... do not interpolate: impose values of climatology DO ji = 1, jpi ii = ji -1 + nimpp sfos(ji,jk) = ssdta(ii,jk,1)*tsmsk(ji,jk) tfos(ji,jk) = tsdta(ii,jk,1)*tsmsk(ji,jk) END DO END DO CLOSE(inum) # if ! defined key_dynspg_fsc ! ... Rigid lid case: make sure vsdta is baroclinic velocity ! ... In rigid lid case vsdta needs to be the baroclinic component. ! ... substract the barotropic velocity component (zvsbtpe) ! ... from the total one (vsdta). CALL obc_cli( vsdta, vclis, fs_njs0, fs_njs1, 1, jpi, nimpp ) # endif ! ... Set normal velocity (on njs0, njs1 <=> jpjsob) DO jk = 1, jpkm1 DO ji = 1, jpi ii = ji -1 + nimpp vfos(ji,jk) = vsdta(ii,jk,1)*vsmsk(ji,jk) END DO END DO END IF ELSE ! 2. CALL not the first step ... ! do not read but interpolate between months. ! here no interpolation is done but the code is kept as a reminder ! ---------------------------------------------------------------------- IF( lp_obc_east ) THEN DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp sfoe(jj,jk) = sedta(ij,jk,1)*temsk(jj,jk) tfoe(jj,jk) = tedta(ij,jk,1)*temsk(jj,jk) ufoe(jj,jk) = uedta(ij,jk,1)*uemsk(jj,jk) END DO END DO END IF IF( lp_obc_west ) THEN DO jk = 1, jpkm1 DO jj = 1, jpj ij = jj -1 + njmpp sfow(jj,jk) = swdta(ij,jk,1)*twmsk(jj,jk) tfow(jj,jk) = twdta(ij,jk,1)*twmsk(jj,jk) ufow(jj,jk) = uwdta(ij,jk,1)*uwmsk(jj,jk) END DO END DO END IF IF( lp_obc_north ) THEN DO jk = 1, jpkm1 DO ji = 1, jpi ii = ji -1 + nimpp sfon(ji,jk) = sndta(ii,jk,1)*tnmsk(ji,jk) tfon(ji,jk) = tndta(ii,jk,1)*tnmsk(ji,jk) vfon(ji,jk) = vndta(ii,jk,1)*vnmsk(ji,jk) END DO END DO END IF IF( lp_obc_south ) THEN DO jk = 1, jpkm1 DO ji = 1, jpi ii = ji -1 + nimpp sfos(ji,jk) = ssdta(ii,jk,1)*tsmsk(ji,jk) tfos(ji,jk) = tsdta(ii,jk,1)*tsmsk(ji,jk) vfos(ji,jk) = vsdta(ii,jk,1)*vsmsk(ji,jk) END DO END DO END IF END IF END SUBROUTINE obc_dta # if defined key_dynspg_fsc !!----------------------------------------------------------------------------- !! 'key_dynspg_fsc' free surface with constant volume !!----------------------------------------------------------------------------- SUBROUTINE obc_dta_psi ( kt ) ! Empty routine !! * Arguments INTEGER,INTENT(in) :: kt WRITE(*,*) 'obc_dta_psi: You should not have seen this print! error?', kt END SUBROUTINE obc_dta_psi #else !!----------------------------------------------------------------------------- !! Default option Rigid-lid !!----------------------------------------------------------------------------- SUBROUTINE obc_dta_psi ( kt ) !!----------------------------------------------------------------------------- !! *** SUBROUTINE obc_dta_psi *** !! !! ** Purpose : !! Update the climatological streamfunction OBC at each time step. !! Depends on the user's configuration. Here data are read only once !! at the beginning of the run. !! !! ** Method : !! 1. initialization !! kbsfstart: number of time steps over which increase bsf !! during initialization. This is provided for a smooth start !! in cases where the transport is large (like on the Antarctic !! continent). also note that when kbfstart=1, the transport !! increases a lot in one time step and the precision usually !! required for the solver may not be enough. !! 2. set the time evolution of the climatological barotropic streamfunction !! along the isolated coastlines ( gcbic(jnic) ). !! 3. set the climatological barotropic streamfunction at the boundary. !! !! The last two steps are done only at first step (nit000) or if kt <= kbfstart !! !! History : !! ! 97-08 (G. Madec, J.M. Molines) !! 8.5 ! 02-10 (C. Talandier, A-M. Treguier) Free surface, F90 !!---------------------------------------------------------------------------- !! * Arguments INTEGER, INTENT( in ) :: kt ! ocean time-step index !! * Local declarations INTEGER :: ji, jj, jnic, jip ! dummy loop indices INTEGER :: inum = 11 ! temporary logical unit INTEGER :: ip, ii, ij, iii, ijj INTEGER :: kbsfstart REAL(wp) :: zsver1, zsver2, zsver3, z2dtr, zcoef !!---------------------------------------------------------------------------- ! 1. initialisation ! ----------------- kbsfstart = 1 zsver1 = bsfic0(1) zsver2 = zsver1 IF( kt <= kbsfstart ) THEN zcoef = float(kt)/float(kbsfstart) ELSE zcoef = 1.e0 END IF bsfic(1) = zsver1*zcoef IF( lwp .AND. ( kt <= kbsfstart ) ) THEN IF(lwp) WRITE(numout,*)' ' IF(lwp) WRITE(numout,*)'obcdta: spinup phase in obc_dta_psi routine' IF(lwp) WRITE(numout,*)'~~~~~~ it=',kt,' OBC: spinup coef: ', & zcoef, ' and transport: ',bsfic(1) END IF zsver2 = bsfic(1)-bsfic(2) zsver3 = bsfic(2) ! 2. Right hand side of the barotropic elliptic equation (isolated coastlines) ! ---------------------------------------------------------------------------- IF( ( neuler == 0 ) .AND. ( kt == nit000 ) ) THEN z2dtr = 1./rdt ELSE z2dtr = 1./2./rdt END IF ! ... bsfb(ii,ij) should be constant but due to the Asselin filter it ! ... converges asymptotically towards bsfic(jnic) ! ... However, bsfb(ii,ij) is constant along the same coastlines ! ... ---> can be improved using an extra array for storing bsficb (before) IF( nbobc > 1 ) THEN DO jnic = 1,nbobc - 1 gcbic(jnic) = 0.e0 ip=mnic(0,jnic) DO jip = 1,ip ii = miic(jip,0,jnic) ij = mjic(jip,0,jnic) IF( ii >= nldi+ nimpp - 1 .AND. ii <= nlei+ nimpp - 1 .AND. & ij >= nldj+ njmpp - 1 .AND. ij <= nlej+ njmpp - 1 ) THEN iii=ii-nimpp+1 ijj=ij-njmpp+1 gcbic(jnic) = ( bsfic(jnic) - bsfb(iii,ijj) ) * z2dtr END IF END DO END DO END IF IF( lk_mpp ) CALL mpp_isl( gcbic, 3 ) ! 3. Update the climatological barotropic function at the boundary ! ---------------------------------------------------------------- IF( lp_obc_east ) THEN IF( kt == nit000 .OR. kt <= kbsfstart ) THEN OPEN(inum,file='obceastbsf.dta') READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) (bfoe(jj),jj=jpjed, jpjef) CLOSE(inum) END IF DO jj=jpjed, jpjefm1 bfoe(jj)=bfoe(jj)*zcoef END DO END IF IF( lp_obc_west ) THEN IF( kt == nit000 .OR. kt <= kbsfstart ) then OPEN(inum,file='obcwestbsf.dta') READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) (bfow(jj),jj=jpjwd, jpjwf) CLOSE(inum) END IF DO jj=jpjwd, jpjwfm1 bfow(jj) = bfow(jj) * zcoef END DO END IF IF( lp_obc_south ) THEN IF( kt == nit000 .OR. kt <= kbsfstart ) THEN OPEN(inum,file='obcsouthbsf.dta') READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) (bfos(jj),jj=jpisd, jpisf) CLOSE(inum) END IF DO ji=jpisd, jpisfm1 bfos(ji)=bfos(ji)*zcoef END DO END IF IF( lp_obc_north ) THEN IF( kt == nit000 .OR. kt <= kbsfstart ) THEN OPEN(inum,file='obcnorthbsf.dta') READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) READ(inum,*) (bfon(jj),jj=jpind, jpinf) CLOSE(inum) END IF DO ji=jpind, jpinfm1 bfon(ji)=bfon(ji)*zcoef END DO END IF END SUBROUTINE obc_dta_psi # endif #else !!------------------------------------------------------------------------------ !! default option: Dummy module NO Open Boundary Conditions !!------------------------------------------------------------------------------ CONTAINS SUBROUTINE obc_dta( kt ) ! Dummy routine INTEGER, INTENT (in) :: kt WRITE(*,*) 'obc_dta: You should not have seen this print! error?', kt END SUBROUTINE obc_dta #endif !!============================================================================== END MODULE obcdta