pkg/dic/bio_export.F

C $Header: /u/gcmpack/MITgcm/pkg/dic/bio_export.F,v 1.27 2010/07/30 15:02:43 mlosch Exp $
C $Name:  $

#include "DIC_OPTIONS.h"

CBOP
C !ROUTINE: BIO_EXPORT

C !INTERFACE: ==========================================================
      SUBROUTINE BIO_EXPORT( PTR_PO4 ,
#ifdef ALLOW_FE
     I           PTR_FE,
#endif
     I           bioac,
     I           bi,bj,imin,imax,jmin,jmax,
     I           myIter,myTime,myThid)

c !DESCRIPTION:
C  Calculate biological activity and export

C !USES: ===============================================================
      IMPLICIT NONE
#include "SIZE.h"
#include "DYNVARS.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "DIC_VARS.h"
#ifdef USE_QSW
#include "FFIELDS.h"
#endif

C !INPUT PARAMETERS: ===================================================
C  myThid               :: thread number
C  myIter               :: current timestep
C  myTime               :: current time
C  PTR_PO4              :: phosphate tracer field
C  PTR_FE               :: iron tracer field
      INTEGER myIter
      _RL myTime
      INTEGER myThid
      _RL  PTR_PO4(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
#ifdef ALLOW_FE
      _RL  PTR_FE(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
#endif
      INTEGER imin, imax, jmin, jmax, bi, bj


C !OUTPUT PARAMETERS: ==================================================
C  bioac               :: biological productivity (will be split
C                         between export and dissolved pool)
      _RL bioac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR)

#if (defined ALLOW_PTRACERS && defined DIC_BIOTIC)

C !LOCAL VARIABLES: ====================================================
C  i,j,k                  :: loop indices
       INTEGER I,J,k
      _RL sfac    (1-OLy:sNy+OLy)
      _RL lit     (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
      _RL atten
      _RL nutlimit
      _RL tmppo4
#ifdef ALLOW_FE      
      _RL tmpfe
#endif
#ifdef AD_SAFE
      _RL thx, thy, theps, thaux
#endif
CEOP


#ifndef READ_PAR
#ifndef USE_QSW
      CALL INSOL(myTime,sfac,bi,bj,myThid)
#endif
#endif
cph: following init assumes nlev=nr
cph: set in dic_init_fixed.F (ALLOW_AUTODIFF)
C$TAF INIT bio_export = static, nsx*nsy*nr
C FORTRAN-77 dynamic memory uses adstore adresto
CRG C$TAF INIT bio_export = memory
      DO j=jmin,jmax
       DO i=imin,imax
#ifdef READ_PAR
        lit(i,j)=PAR(i,j,bi,bj)
#elif (defined USE_QSW)
        lit(i,j)=-parfrac*Qsw(i,j,bi,bj)*maskC(i,j,1,bi,bj)
#else
        lit(i,j)=sfac(j)
#endif

        IF ( .NOT. QSW_underice ) THEN
c if using Qsw but not seaice/thsice or coupled, then
c ice fraction needs to be taken into account
         lit(i,j)=lit(i,j)*(1. _d 0 - FIce(i,j,bi,bj))
        ENDIF
       ENDDO
      ENDDO
c
      DO k=1,nlev
C$TAF STORE lit = bio_export
       DO j=jmin,jmax
        DO i=imin,imax
         atten=(k0*drF(k)*hFacC(i,j,k,bi,bj)*.5 _d 0)
         if (k.gt.1) atten=atten+(k0*drF(k-1)
     &                         *hFacC(i,j,k-1,bi,bj)*.5 _d 0)
         lit(i,j)=lit(i,j)*exp(-atten)
#ifndef TARGET_NEC_SX
C     this statement breaks vectorization and causes a dramatic
C     performance drop on vector computers
         IF (lit(i,j).LT.0. _d 0.OR.lit(i,j).GT.350. _d 0) THEN
          print*,'QQ lit',i,j,lit(i,j)
         ENDIF
#endif 

#ifdef DIC_NO_NEG
         tmppo4=max(0. _d 0, PTR_PO4(i,j,k))
         lit(i,j)=max(0. _d 0,lit(i,j))
#else
         tmppo4=PTR_PO4(i,j,k)
#endif

#ifdef ALLOW_FE
#ifdef DIC_NO_NEG
         tmpfe=max(0. _d 0,PTR_FE(i,j,k))
#else
         tmpfe=PTR_FE(i,j,k)
#endif
#ifdef AD_SAFE
         thx = tmppo4/(tmppo4+KPO4)
         thy = tmpfe/(tmpfe+KFE)
c        thx = PTR_PO4(i,j,k)/(PTR_PO4(i,j,k)+KPO4)
c        thy = PTR_FE(i,j,k)/(PTR_FE(i,j,k)+KFE)
c        theps = 1. _d -6
c        thaux = tanh( (thx-thy)/theps )
         thaux = tanh( (thx-thy)*1. _d 6 )
         nutlimit= ( 1. _d 0 - thaux ) * thx * 0.5 _d 0
     &        +    ( 1. _d 0 + thaux ) * thy * 0.5 _d 0
#else
         nutlimit=min( tmppo4/(tmppo4+KPO4),tmpfe/(tmpfe+KFE) )
#endif
#else
         nutlimit=     tmppo4/(tmppo4+KPO4)
#endif

         bioac(i,j,k)=alpha(i,j,bi,bj)*
     &        lit(i,j)/(lit(i,j)+lit0)*maskC(i,j,k,bi,bj)*
     &        nutlimit
        ENDDO
       ENDDO
      ENDDO
c
#endif /* definded ALLOW_PTRACERS && defined DIC_BIOTIC */
      RETURN
      END
1	mlosch	1.28	C $Header: /u/gcmpack/MITgcm/pkg/dic/bio_export.F,v 1.27 2010/07/30 15:02:43 mlosch Exp $
2	jmc	1.15	C $Name: $
3
4	edhill	1.3	#include "DIC_OPTIONS.h"
5	stephd	1.1
6	stephd	1.5	CBOP
7			C !ROUTINE: BIO_EXPORT
8
9			C !INTERFACE: ==========================================================
10	jmc	1.15	SUBROUTINE BIO_EXPORT( PTR_PO4 ,
11	stephd	1.1	#ifdef ALLOW_FE
12	jmc	1.15	I PTR_FE,
13			#endif
14	stephd	1.1	I bioac,
15			I bi,bj,imin,imax,jmin,jmax,
16			I myIter,myTime,myThid)
17
18	stephd	1.5	c !DESCRIPTION:
19	jmc	1.15	C Calculate biological activity and export
20	stephd	1.5
21			C !USES: ===============================================================
22	stephd	1.1	IMPLICIT NONE
23			#include "SIZE.h"
24			#include "DYNVARS.h"
25			#include "EEPARAMS.h"
26			#include "PARAMS.h"
27			#include "GRID.h"
28	dfer	1.22	#include "DIC_VARS.h"
29	stephd	1.12	#ifdef USE_QSW
30	stephd	1.11	#include "FFIELDS.h"
31			#endif
32	stephd	1.1
33	stephd	1.5	C !INPUT PARAMETERS: ===================================================
34			C myThid :: thread number
35			C myIter :: current timestep
36			C myTime :: current time
37			C PTR_PO4 :: phosphate tracer field
38			C PTR_FE :: iron tracer field
39	stephd	1.1	INTEGER myIter
40			_RL myTime
41			INTEGER myThid
42			_RL PTR_PO4(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
43			#ifdef ALLOW_FE
44			_RL PTR_FE(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
45			#endif
46	stephd	1.5	INTEGER imin, imax, jmin, jmax, bi, bj
47
48
49			C !OUTPUT PARAMETERS: ==================================================
50			C bioac :: biological productivity (will be split
51			C between export and dissolved pool)
52	stephd	1.1	_RL bioac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nR)
53
54	mlosch	1.27	#if (defined ALLOW_PTRACERS && defined DIC_BIOTIC)
55	stephd	1.5
56			C !LOCAL VARIABLES: ====================================================
57			C i,j,k :: loop indices
58	stephd	1.1	INTEGER I,J,k
59	mlosch	1.27	_RL sfac (1-OLy:sNy+OLy)
60			_RL lit (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
61			_RL atten
62	stephd	1.10	_RL nutlimit
63	jmc	1.26	_RL tmppo4
64	mlosch	1.27	#ifdef ALLOW_FE
65	jmc	1.26	_RL tmpfe
66			#endif
67	stephd	1.10	#ifdef AD_SAFE
68	mlosch	1.27	_RL thx, thy, theps, thaux
69	stephd	1.10	#endif
70	stephd	1.5	CEOP
71	stephd	1.10
72
73	stephd	1.8	#ifndef READ_PAR
74	stephd	1.12	#ifndef USE_QSW
75	mlosch	1.27	CALL INSOL(myTime,sfac,bi,bj,myThid)
76	stephd	1.8	#endif
77	stephd	1.11	#endif
78	heimbach	1.25	cph: following init assumes nlev=nr
79			cph: set in dic_init_fixed.F (ALLOW_AUTODIFF)
80	mlosch	1.28	C$TAF INIT bio_export = static, nsxnsynr
81	stephd	1.2	C FORTRAN-77 dynamic memory uses adstore adresto
82			CRG C$TAF INIT bio_export = memory
83	mlosch	1.28	DO j=jmin,jmax
84			DO i=imin,imax
85	stephd	1.8	#ifdef READ_PAR
86	mlosch	1.27	lit(i,j)=PAR(i,j,bi,bj)
87	dfer	1.20	#elif (defined USE_QSW)
88	mlosch	1.27	lit(i,j)=-parfracQsw(i,j,bi,bj)maskC(i,j,1,bi,bj)
89	stephd	1.11	#else
90	mlosch	1.27	lit(i,j)=sfac(j)
91	stephd	1.8	#endif
92	stephd	1.18
93	mlosch	1.27	IF ( .NOT. QSW_underice ) THEN
94	dfer	1.20	c if using Qsw but not seaice/thsice or coupled, then
95			c ice fraction needs to be taken into account
96	mlosch	1.27	lit(i,j)=lit(i,j)*(1. _d 0 - FIce(i,j,bi,bj))
97			ENDIF
98			ENDDO
99			ENDDO
100	stephd	1.19	c
101	mlosch	1.27	DO k=1,nlev
102	mlosch	1.28	C$TAF STORE lit = bio_export
103	mlosch	1.27	DO j=jmin,jmax
104			DO i=imin,imax
105			atten=(k0drF(k)hFacC(i,j,k,bi,bj)*.5 _d 0)
106			if (k.gt.1) atten=atten+(k0*drF(k-1)
107	dfer	1.17	& hFacC(i,j,k-1,bi,bj).5 _d 0)
108	mlosch	1.27	lit(i,j)=lit(i,j)*exp(-atten)
109			#ifndef TARGET_NEC_SX
110			C this statement breaks vectorization and causes a dramatic
111			C performance drop on vector computers
112			IF (lit(i,j).LT.0. _d 0.OR.lit(i,j).GT.350. _d 0) THEN
113			print*,'QQ lit',i,j,lit(i,j)
114			ENDIF
115			#endif
116	stephd	1.10
117	stephd	1.21	#ifdef DIC_NO_NEG
118	mlosch	1.27	tmppo4=max(0. _d 0, PTR_PO4(i,j,k))
119			lit(i,j)=max(0. _d 0,lit(i,j))
120	stephd	1.21	#else
121	mlosch	1.27	tmppo4=PTR_PO4(i,j,k)
122	stephd	1.21	#endif
123
124	stephd	1.1	#ifdef ALLOW_FE
125	stephd	1.21	#ifdef DIC_NO_NEG
126	mlosch	1.27	tmpfe=max(0. _d 0,PTR_FE(i,j,k))
127	stephd	1.21	#else
128	mlosch	1.27	tmpfe=PTR_FE(i,j,k)
129	stephd	1.21	#endif
130	stephd	1.10	#ifdef AD_SAFE
131	mlosch	1.27	thx = tmppo4/(tmppo4+KPO4)
132			thy = tmpfe/(tmpfe+KFE)
133			c thx = PTR_PO4(i,j,k)/(PTR_PO4(i,j,k)+KPO4)
134			c thy = PTR_FE(i,j,k)/(PTR_FE(i,j,k)+KFE)
135			c theps = 1. _d -6
136			c thaux = tanh( (thx-thy)/theps )
137			thaux = tanh( (thx-thy)*1. _d 6 )
138			nutlimit= ( 1. _d 0 - thaux ) * thx * 0.5 _d 0
139			& + ( 1. _d 0 + thaux ) * thy * 0.5 _d 0
140	stephd	1.5	#else
141	mlosch	1.27	nutlimit=min( tmppo4/(tmppo4+KPO4),tmpfe/(tmpfe+KFE) )
142	stephd	1.1	#endif
143	stephd	1.10	#else
144	mlosch	1.27	nutlimit= tmppo4/(tmppo4+KPO4)
145	stephd	1.10	#endif
146
147	mlosch	1.27	bioac(i,j,k)=alpha(i,j,bi,bj)*
148			& lit(i,j)/(lit(i,j)+lit0)maskC(i,j,k,bi,bj)
149			& nutlimit
150	stephd	1.1	ENDDO
151	mlosch	1.27	ENDDO
152			ENDDO
153	stephd	1.1	c
154	mlosch	1.27	#endif /* definded ALLOW_PTRACERS && defined DIC_BIOTIC */
155			RETURN
156			END