source: trunk/NEMO/TOP_SRC/PISCES/p4zsink2.F90 @ 935

Last change on this file since 935 was 935, checked in by cetlod, 13 years ago

adding modules for PISCES SMS model, see ticket 141

File size: 5.8 KB
Line 
1MODULE p4zsink2
2   !!======================================================================
3   !!                         ***  MODULE p4zsink2  ***
4   !! TOP :   PISCES Compute vertical flux of particulate matter due to gravitational sinking
5   !!======================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!----------------------------------------------------------------------
9#if defined key_pisces
10   !!----------------------------------------------------------------------
11   !!   'key_pisces'                                       PISCES bio-model
12   !!----------------------------------------------------------------------
13   !!   p4z_sink2       :  Compute vertical flux of particulate matter due to gravitational sinking
14   !!----------------------------------------------------------------------
15   USE oce_trc         !
16   USE trp_trc
17   USE sms
18
19   IMPLICIT NONE
20   PRIVATE
21
22   PUBLIC   p4z_sink2    ! called in p4zbio.F90
23
24   !!* Substitution
25#  include "domzgr_substitute.h90"
26   !!----------------------------------------------------------------------
27   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
28   !! $Header:$
29   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
30   !!----------------------------------------------------------------------
31
32CONTAINS
33
34   SUBROUTINE p4z_sink2( wstmp, sinktemp, jn )
35      !!---------------------------------------------------------------------
36      !!                     ***  ROUTINE p4z_sink2  ***
37      !!
38      !! ** Purpose :   Compute the sedimentation terms for the various sinking
39      !!     particles. The scheme used to compute the trends is based
40      !!     on MUSCL.
41      !!
42      !! ** Method  : - this ROUTINE compute not exactly the advection but the
43      !!      transport term, i.e.  div(u*tra).
44      !!---------------------------------------------------------------------
45      INTEGER , INTENT(in   )                         ::   jn         ! tracer index index     
46      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   wstmp      ! ???
47      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   sinktemp   ! ???
48      !!
49      INTEGER  ::   ji, jj, jk, jnt
50      REAL(wp) ::   zigma,zew,zstep,zign
51      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztraz, zakz
52      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zkz  , zwstmp2
53      !!---------------------------------------------------------------------
54
55      zstep  = rfact2 / 2.
56
57      ztraz(:,:,:) = 0.e0
58      zkz  (:,:,:) = 0.e0
59      zakz (:,:,:) = 0.e0
60
61      DO jk = 1, jpkm1
62# if defined key_off_degrad
63         zwstmp2(:,:,jk+1)=-wstmp(:,:,jk)/rjjss*tmask(:,:,jk+1)*facvol(:,:,jk)
64# else
65         zwstmp2(:,:,jk+1)=-wstmp(:,:,jk)/rjjss*tmask(:,:,jk+1)
66
67# endif
68      END DO
69 
70      zwstmp2(:,:,1) = 0.e0
71!
72! Vertical advective flux
73!-------------------------------
74
75      DO jnt = 1, 2
76
77! ... first guess of the slopes interior values
78
79         DO jk = 2, jpkm1
80            ztraz(:,:,jk) = (trn(:,:,jk-1,jn) - trn(:,:,jk,jn)) *tmask(:,:,jk)
81         END DO
82
83         ztraz(:,:,1  ) = 0.0
84         ztraz(:,:,jpk) = 0.0
85!
86! slopes
87         DO jk=2,jpkm1
88            DO jj = 1,jpj
89               DO ji = 1, jpi
90                  zign = 0.25 + SIGN( 0.25,ztraz(ji,jj,jk)*ztraz(ji,jj,jk+1) )
91                  zakz(ji,jj,jk) = (ztraz(ji,jj,jk) + ztraz(ji,jj,jk+1) ) * zign
92               END DO
93            END DO
94         END DO
95         !
96         ! Slopes limitation
97         DO jk = 2, jpkm1
98            DO jj = 1, jpj
99               DO ji = 1, jpi
100                  zakz(ji,jj,jk) = SIGN(1.,zakz(ji,jj,jk)) *        &
101                     &             MIN(ABS(zakz(ji,jj,jk)),         &
102                     &             2.*ABS(ztraz(ji,jj,jk+1)),       &
103                     &             2.*ABS(ztraz(ji,jj,jk)))
104               END DO
105            END DO
106         END DO
107         
108         ! vertical advective flux
109         DO jk = 1, jpkm1
110            DO jj = 1, jpj     
111               DO ji = 1, jpi   
112                  zigma = zwstmp2(ji,jj,jk+1)*zstep/fse3w(ji,jj,jk+1)
113                  zew   = zwstmp2(ji,jj,jk+1)
114                  sinktemp(ji,jj,jk+1) = -zew*(trn(ji,jj,jk,jn)                &
115                     &                 -0.5*(1+zigma)*zakz(ji,jj,jk))*zstep
116               END DO
117            END DO
118         END DO
119         !
120         ! Boundary conditions
121         sinktemp(:,:,1  ) = 0.e0
122         sinktemp(:,:,jpk) = 0.e0
123         
124         DO jk=1,jpkm1
125            DO jj = 1,jpj
126               DO ji = 1, jpi
127                  trn(ji,jj,jk,jn) = trn(ji,jj,jk,jn)               &
128                     &        + (sinktemp(ji,jj,jk)-sinktemp(ji,jj,jk+1))     &
129                     &        /fse3t(ji,jj,jk)
130               END DO
131            END DO
132         END DO
133
134      ENDDO
135
136      DO jk=1,jpkm1
137         DO jj = 1,jpj
138            DO ji = 1, jpi
139               !
140               trb(ji,jj,jk,jn) = trb(ji,jj,jk,jn)      &
141               &        + 2.*(sinktemp(ji,jj,jk)-sinktemp(ji,jj,jk+1)) &
142               &        /fse3t(ji,jj,jk)
143!
144          ENDDO
145        ENDDO
146      ENDDO
147!
148        trn(:,:,:,jn)   = trb(:,:,:,jn)
149        sinktemp(:,:,:) = 2. * sinktemp(:,:,:)
150
151      !
152   END SUBROUTINE p4z_sink2
153
154#else
155   !!======================================================================
156   !!  Dummy module :                                   No PISCES bio-model
157   !!======================================================================
158CONTAINS
159   SUBROUTINE p4z_sink2( wstmp, sinktemp, jn )         ! Empty routine
160      INTEGER, INTENT( in ) ::   jn
161      REAL   , INTENT( in ) ::   wstmp,sinktemp   
162      WRITE(*,*) 'p4z_sink2: You should not have seen this print! error?', jn, wstmp, sinktemp
163   END SUBROUTINE p4z_sink2
164#endif 
165
166   !!======================================================================
167END MODULE  p4zsink2
Note: See TracBrowser for help on using the repository browser.