Index: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90
===================================================================
--- NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90	(revision 10954)
+++ NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tranpc.F90	(revision 10985)
@@ -42,5 +42,5 @@
 CONTAINS
 
-   SUBROUTINE tra_npc( kt, Kmm, Krhs )
+   SUBROUTINE tra_npc( kt, Kmm, Krhs, pts, Kaa )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tranpc  ***
@@ -58,6 +58,7 @@
       !! References :     Madec, et al., 1991, JPO, 21, 9, 1349-1371.
       !!----------------------------------------------------------------------
-      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
-      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
+      INTEGER,                                   INTENT(in   ) :: kt              ! ocean time-step index
+      INTEGER,                                   INTENT(in   ) :: Kmm, Krhs, Kaa  ! time level indices
+      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts             ! active tracers and RHS of tracer equation
       !
       INTEGER  ::   ji, jj, jk   ! dummy loop indices
@@ -67,10 +68,10 @@
       REAL(wp) ::   zta, zalfa, zsum_temp, zsum_alfa, zaw, zdz, zsum_z
       REAL(wp) ::   zsa, zbeta, zsum_sali, zsum_beta, zbw, zrw, z1_r2dt
-      REAL(wp), PARAMETER ::   zn2_zero = 1.e-14_wp      ! acceptance criteria for neutrality (N2==0)
-      REAL(wp), DIMENSION(        jpk     ) ::   zvn2         ! vertical profile of N2 at 1 given point...
-      REAL(wp), DIMENSION(        jpk,jpts) ::   zvts, zvab   ! vertical profile of T & S , and  alpha & betaat 1 given point
-      REAL(wp), DIMENSION(jpi,jpj,jpk     ) ::   zn2          ! N^2 
-      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::   zab          ! alpha and beta
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdt, ztrds   ! 3D workspace
+      REAL(wp), PARAMETER ::   zn2_zero = 1.e-14_wp             ! acceptance criteria for neutrality (N2==0)
+      REAL(wp), DIMENSION(        jpk     )   ::   zvn2         ! vertical profile of N2 at 1 given point...
+      REAL(wp), DIMENSION(        jpk,jpts)   ::   zvts, zvab   ! vertical profile of T & S , and  alpha & betaat 1 given point
+      REAL(wp), DIMENSION(jpi,jpj,jpk     )   ::   zn2          ! N^2 
+      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts)   ::   zab          ! alpha and beta
+      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdt, ztrds ! 3D workspace
       !
       LOGICAL, PARAMETER :: l_LB_debug = .FALSE. ! set to true if you want to follow what is
@@ -85,6 +86,6 @@
          IF( l_trdtra )   THEN                    !* Save initial after fields
             ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) )
-            ztrdt(:,:,:) = ts(:,:,:,jp_tem,Krhs) 
-            ztrds(:,:,:) = ts(:,:,:,jp_sal,Krhs)
+            ztrdt(:,:,:) = pts(:,:,:,jp_tem,Kaa) 
+            ztrds(:,:,:) = pts(:,:,:,jp_sal,Kaa)
          ENDIF
          !
@@ -96,6 +97,6 @@
          ENDIF
          !
-         CALL eos_rab( ts(:,:,:,:,Krhs), zab, Kmm )         ! after alpha and beta (given on T-points)
-         CALL bn2    ( ts(:,:,:,:,Krhs), zab, zn2, Kmm )    ! after Brunt-Vaisala  (given on W-points)
+         CALL eos_rab( pts(:,:,:,:,Kaa), zab, Kmm )         ! after alpha and beta (given on T-points)
+         CALL bn2    ( pts(:,:,:,:,Kaa), zab, zn2, Kmm )    ! after Brunt-Vaisala  (given on W-points)
          !
          inpcc = 0
@@ -106,6 +107,6 @@
                IF( tmask(ji,jj,2) == 1 ) THEN      ! At least 2 ocean points
                   !                                     ! consider one ocean column 
-                  zvts(:,jp_tem) = ts(ji,jj,:,jp_tem,Krhs)      ! temperature
-                  zvts(:,jp_sal) = ts(ji,jj,:,jp_sal,Krhs)      ! salinity
+                  zvts(:,jp_tem) = pts(ji,jj,:,jp_tem,Kaa)      ! temperature
+                  zvts(:,jp_sal) = pts(ji,jj,:,jp_sal,Kaa)      ! salinity
                   !
                   zvab(:,jp_tem)  = zab(ji,jj,:,jp_tem)     ! Alpha 
@@ -287,7 +288,7 @@
                   END DO ! DO WHILE ( .NOT. l_column_treated )
 
-                  !! Updating tsa:
-                  ts(ji,jj,:,jp_tem,Krhs) = zvts(:,jp_tem)
-                  ts(ji,jj,:,jp_sal,Krhs) = zvts(:,jp_sal)
+                  !! Updating pts:
+                  pts(ji,jj,:,jp_tem,Kaa) = zvts(:,jp_tem)
+                  pts(ji,jj,:,jp_sal,Kaa) = zvts(:,jp_sal)
 
                   !! LB:  Potentially some other global variable beside theta and S can be treated here
@@ -303,6 +304,6 @@
          IF( l_trdtra ) THEN         ! send the Non penetrative mixing trends for diagnostic
             z1_r2dt = 1._wp / (2._wp * rdt)
-            ztrdt(:,:,:) = ( ts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:) ) * z1_r2dt
-            ztrds(:,:,:) = ( ts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:) ) * z1_r2dt
+            ztrdt(:,:,:) = ( pts(:,:,:,jp_tem,Kaa) - ztrdt(:,:,:) ) * z1_r2dt
+            ztrds(:,:,:) = ( pts(:,:,:,jp_sal,Kaa) - ztrds(:,:,:) ) * z1_r2dt
             CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_npc, ztrdt )
             CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_npc, ztrds )
@@ -310,5 +311,5 @@
          ENDIF
          !
-         CALL lbc_lnk_multi( 'tranpc', ts(:,:,:,jp_tem,Krhs), 'T', 1., ts(:,:,:,jp_sal,Krhs), 'T', 1. )
+         CALL lbc_lnk_multi( 'tranpc', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. )
          !
          IF( lwp .AND. l_LB_debug ) THEN
