/[MITgcm]/MITgcm_contrib/darwin2/pkg/monod/monod_forcing.F
ViewVC logotype

Diff of /MITgcm_contrib/darwin2/pkg/monod/monod_forcing.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.1 by jahn, Wed Apr 13 18:56:25 2011 UTC revision 1.14 by jahn, Thu Oct 25 15:58:24 2012 UTC
# Line 21  C======================================= Line 21  C=======================================
21  #include "PARAMS.h"  #include "PARAMS.h"
22  #include "GRID.h"  #include "GRID.h"
23  #include "DYNVARS.h"  #include "DYNVARS.h"
24  #ifdef USE_QSW  c for Qsw and/or surfaceForcingT
25    c choice which field to take pCO2 from for pCO2limit
26    c this assumes we use Ttendency from offline
27  #include "FFIELDS.h"  #include "FFIELDS.h"
 #endif  
28  #ifdef ALLOW_LONGSTEP  #ifdef ALLOW_LONGSTEP
29  #include "LONGSTEP.h"  #include "LONGSTEP.h"
30  #endif  #endif
# Line 43  c ANNA include wavebands_params.h Line 44  c ANNA include wavebands_params.h
44  #include "WAVEBANDS_PARAMS.h"  #include "WAVEBANDS_PARAMS.h"
45  #endif  #endif
46    
47    
48  C     === Global variables ===  C     === Global variables ===
49  c tracers  c tracers
50         _RL Ptr(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy,nDarwin)         _RL Ptr(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy,nDarwin)
# Line 68  c plankton arrays Line 70  c plankton arrays
70        _RL  Phy_k(npmax,Nr)        _RL  Phy_k(npmax,Nr)
71        _RL  Phyup(npmax)        _RL  Phyup(npmax)
72        _RL  part_k(Nr)        _RL  part_k(Nr)
73    #ifdef ALLOW_CDOM
74          _RL cdom_k(Nr)
75    #endif
76  c iron partitioning  c iron partitioning
77        _RL  freefe(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL  freefe(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
78  c some working variables  c some working variables
# Line 87  c ANNA define variables for wavebands Line 92  c ANNA define variables for wavebands
92         _RL PARw_k(tlam,Nr)         _RL PARw_k(tlam,Nr)
93         _RL PARwup(tlam)         _RL PARwup(tlam)
94         _RL acdom_k(Nr,tlam)         _RL acdom_k(Nr,tlam)
95           _RL Ek_nll(npmax,tlam)
96           _RL EkoverE_nll(npmax,tlam)
97  #ifdef DAR_RADTRANS  #ifdef DAR_RADTRANS
98         integer iday,iyr,imon,isec,lp,wd,mydate(4)         integer iday,iyr,imon,isec,lp,wd,mydate(4)
99         _RL Edwsf(tlam),Eswsf(tlam)         _RL Edwsf(tlam),Eswsf(tlam)
100         _RL Edz(tlam,Nr),Esz(tlam,Nr),Euz(tlam,Nr),Eutop(tlam,Nr)         _RL Edz(tlam,Nr),Esz(tlam,Nr),Euz(tlam,Nr)
101           _RL Estop(tlam,Nr),Eutop(tlam,Nr)
102         _RL tirrq(nr)         _RL tirrq(nr)
103         _RL tirrwq(tlam,nr)         _RL tirrwq(tlam,nr)
104           _RL amp1(tlam,nr), amp2(tlam,nr)
105         _RL solz         _RL solz
106         _RL rmud         _RL rmud
107         _RL actot,bctot,bbctot         _RL actot,bctot,bbctot
108         _RL apart_k(Nr,tlam),bpart_k(Nr,tlam),bbpart_k(Nr,tlam)         _RL apart_k(Nr,tlam),bpart_k(Nr,tlam),bbpart_k(Nr,tlam)
109         _RL bt_k(Nr,tlam), bb_k(Nr,tlam)         _RL bt_k(Nr,tlam), bb_k(Nr,tlam)
110           _RL discEs, discEu
111           INTEGER idiscEs,jdiscEs,kdiscEs,ldiscEs
112           INTEGER idiscEu,jdiscEu,kdiscEu,ldiscEu
113  #else  #else
114         _RL PARwdn(tlam)         _RL PARwdn(tlam)
115  #endif  #endif
# Line 111  C      always need for diagnostics Line 123  C      always need for diagnostics
123        _RL  Diver2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL  Diver2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
124        _RL  Diver3(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL  Diver3(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
125        _RL  Diver4(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)        _RL  Diver4(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
126          _RL Shannon(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
127          _RL Simpson(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
128    
129        _RL  tmpphy(npmax)        _RL  tmpphy(npmax)
130        _RL  totphy, biotot, maxphy, phymax        _RL  totphy, biotot, maxphy, phymax
# Line 119  C      always need for diagnostics Line 133  C      always need for diagnostics
133  #ifdef GEIDER  #ifdef GEIDER
134        _RL phychl(npmax)        _RL phychl(npmax)
135        _RL phychl_k(npmax,Nr)        _RL phychl_k(npmax,Nr)
136          _RL Ekl(npmax)
137          _RL EkoverEl(npmax)
138          _RL chl2cl(npmax)
139  #ifdef DYNAMIC_CHL  #ifdef DYNAMIC_CHL
140        _RL dphychl(npmax)        _RL dphychl(npmax)
141        _RL chlup(npmax)        _RL chlup(npmax)
142          _RL accliml(npmax)
143  #endif  #endif
144  #endif  #endif
145    #ifdef ALLOW_CDOM
146          _RL cdoml  
147          _RL dcdoml
148    #endif
149    
150  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
151  COJ for diagnostics  COJ for diagnostics
# Line 185  c some local variables Line 207  c some local variables
207         _RL PSiupl         _RL PSiupl
208         _RL Tlocal         _RL Tlocal
209         _RL Slocal         _RL Slocal
210           _RL pCO2local
211         _RL Qswlocal         _RL Qswlocal
212         _RL NH4l         _RL NH4l
213         _RL NO2l         _RL NO2l
# Line 248  c tendencies Line 271  c tendencies
271         _RL do2l         _RL do2l
272         _RL dZooCl(nzmax)         _RL dZooCl(nzmax)
273  c air-sea fluxes  c air-sea fluxes
274         _RL flxCO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)         _RL flxCO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
275         _RL flxALK(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)         _RL flxALK(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
276         _RL flxO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)         _RL flxO2(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
277  #endif  #endif
278    
279         _RL tot_Nfix         _RL tot_Nfix
# Line 277  c Line 300  c
300             Diver2(i,j,k)=0. _d 0             Diver2(i,j,k)=0. _d 0
301             Diver3(i,j,k)=0. _d 0             Diver3(i,j,k)=0. _d 0
302             Diver4(i,j,k)=0. _d 0             Diver4(i,j,k)=0. _d 0
303               Shannon(i,j,k)=0. _d 0
304               Simpson(i,j,k)=1. _d 0
305  #endif  #endif
306    
307  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
308  COJ for diagnostics  COJ for diagnostics
309             PParr(i,j,k) = 0. _d 0             PParr(i,j,k) = 0. _d 0
310             Nfixarr(i,j,k) = 0. _d 0             Nfixarr(i,j,k) = 0. _d 0
311    #ifdef DAR_DIAG_CHL
312               GeiderChlarr(i,j,k) = 0. _d 0
313               GeiderChl2Carr(i,j,k) = 0. _d 0
314               DoneyChlarr(i,j,k) = 0. _d 0
315               DoneyChl2Carr(i,j,k) = 0. _d 0
316               CloernChlarr(i,j,k) = 0. _d 0
317               CloernChl2Carr(i,j,k) = 0. _d 0
318    #endif
319  c ANNA_TAVE  c ANNA_TAVE
320  #ifdef WAVES_DIAG_PCHL  #ifdef WAVES_DIAG_PCHL
321             DO np=1,npmax             DO np=1,npmax
# Line 307  COJ Line 340  COJ
340          enddo          enddo
341         ENDDO         ENDDO
342         ENDDO         ENDDO
343    
344    #ifdef DAR_RADTRANS
345           idiscEs = 0
346           jdiscEs = 0
347           kdiscEs = 0
348           ldiscEs = 0
349           idiscEu = 0
350           jdiscEu = 0
351           kdiscEu = 0
352           ldiscEu = 0
353           discEs = 0.
354           discEu = 0.
355    #endif
356  c  c
357  c bio-chemical time loop  c bio-chemical time loop
358  c--------------------------------------------------  c--------------------------------------------------
# Line 405  C ========================== i,j loops = Line 451  C ========================== i,j loops =
451  c ------------ these are convenient ------------------------------------  c ------------ these are convenient ------------------------------------
452           DO k=1,Nr           DO k=1,Nr
453            part_k(k)  = max(Ptr(i,j,k,bi,bj,iPOP),0. _d 0)            part_k(k)  = max(Ptr(i,j,k,bi,bj,iPOP),0. _d 0)
454    #ifdef ALLOW_CDOM
455              cdom_k(k)  = max(Ptr(i,j,k,bi,bj,iCDOM),0. _d 0)
456    #endif
457            DO np = 1,npmax            DO np = 1,npmax
458              Phy_k(np,k) = max(Ptr(i,j,k,bi,bj,iPhy+np-1),0. _d 0)              Phy_k(np,k) = max(Ptr(i,j,k,bi,bj,iPhy+np-1),0. _d 0)
459  #ifdef GEIDER  #ifdef GEIDER
# Line 420  c ------------ these are convenient ---- Line 469  c ------------ these are convenient ----
469  c ------------ GET CDOM_k FOR WAVEBANDS_3D and RADTRANS ----------------  c ------------ GET CDOM_k FOR WAVEBANDS_3D and RADTRANS ----------------
470  #ifdef WAVEBANDS  #ifdef WAVEBANDS
471  #if defined(DAR_CALC_ACDOM) || defined(DAR_RADTRANS)  #if defined(DAR_CALC_ACDOM) || defined(DAR_RADTRANS)
472    #ifdef ALLOW_CDOM
473            call MONOD_ACDOM(cdom_k,
474         O                     acdom_k,
475         I                     myThid)
476    #else
477           call MONOD_ACDOM(phychl_k,aphy_chl,aw,           call MONOD_ACDOM(phychl_k,aphy_chl,aw,
478       O                     acdom_k,       O                     acdom_k,
479       I                     myThid)       I                     myThid)
480    #endif
481  #else  #else
482           DO k=1,Nr           DO k=1,Nr
483            DO ilam = 1,tlam            DO ilam = 1,tlam
# Line 450  c ------------ GET INCIDENT NON-SPECTRAL Line 505  c ------------ GET INCIDENT NON-SPECTRAL
505    
506  #else /* not USE_QSW */  #else /* not USE_QSW */
507    
508           lite = sfac(j)*maskC(i,j,1,bi,bj)/86400*1 _d 6  C        convert W/m2 to uEin/s/m2
509             lite = sfac(j)*parconv*maskC(i,j,1,bi,bj)
510    
511  #endif /* not USE_QSW */  #endif /* not USE_QSW */
512  #endif /* not READ_PAR */  #endif /* not READ_PAR */
# Line 593  c           add water and CDOM Line 649  c           add water and CDOM
649              bt_k(k,ilam) = bw(ilam) + bctot + bpart_k(k,ilam)              bt_k(k,ilam) = bw(ilam) + bctot + bpart_k(k,ilam)
650              bb_k(k,ilam) = darwin_bbw*bw(ilam)+bbctot+bbpart_k(k,ilam)              bb_k(k,ilam) = darwin_bbw*bw(ilam)+bbctot+bbpart_k(k,ilam)
651              bb_k(k,ilam) = MAX(darwin_bbmin, bb_k(k,ilam))              bb_k(k,ilam) = MAX(darwin_bbmin, bb_k(k,ilam))
652    c           initialize output variables
653                Edz(ilam,k) = 0.0
654                Esz(ilam,k) = 0.0
655                Euz(ilam,k) = 0.0
656                Estop(ilam,k) = 0.0
657                Eutop(ilam,k) = 0.0
658                amp1(ilam,k) = 0.0
659                amp2(ilam,k) = 0.0
660            ENDDO            ENDDO
661           ENDDO           ENDDO
662    
663  #ifdef DAR_RADTRANS_ITERATIVE           IF (darwin_radtrans_niter.GE.0) THEN
664           call MONOD_RADTRANS_ITER(             call MONOD_RADTRANS_ITER(
665       I                dz_k,rmud,Edwsf,Eswsf,a_k,bt_k,bb_k,       I                dz_k,rmud,Edwsf,Eswsf,a_k,bt_k,bb_k,
666       I                darwin_radtrans_kmax,darwin_radtrans_niter,       I                darwin_radtrans_kmax,darwin_radtrans_niter,
667       O                Edz,Esz,Euz,Eutop,       O                Edz,Esz,Euz,Eutop,
668       O                tirrq,tirrwq,       O                tirrq,tirrwq,
669         O                amp1,amp2,
670       I                myThid)       I                myThid)
671  #else           ELSEIF (darwin_radtrans_niter.EQ.-1) THEN
672  c dzlocal ?????  c dzlocal ?????
673           call MONOD_RADTRANS(             call MONOD_RADTRANS(
674       I                drF,rmud,Edwsf,Eswsf,a_k,bt_k,bb_k,       I                drF,rmud,Edwsf,Eswsf,a_k,bt_k,bb_k,
675       O                Edz,Esz,Euz,Eutop,       O                Edz,Esz,Euz,Eutop,
676       O                tirrq,tirrwq,       O                tirrq,tirrwq,
677       I                myThid)       I                myThid)
678             ELSE
679               call MONOD_RADTRANS_DIRECT(
680         I                dz_k,rmud,Edwsf,Eswsf,a_k,bt_k,bb_k,
681         I                darwin_radtrans_kmax,
682         O                Edz,Esz,Euz,Estop,Eutop,
683         O                tirrq,tirrwq,
684         O                amp1,amp2,
685         I                myThid)
686    #ifdef DAR_CHECK_IRR_CONT
687               IF( dz_k(1) .GT. 0.0 )THEN
688               DO ilam = 1,tlam
689               IF(Eswsf(ilam).GE.darwin_radmodThresh .OR.
690         &        Edwsf(ilam).GE.darwin_radmodThresh ) THEN
691                IF(ABS(Estop(ilam,1)-Eswsf(ilam)) .GT. discEs )THEN
692                  discEs = ABS(Estop(ilam,1)-Eswsf(ilam))
693                  idiscEs = i
694                  jdiscEs = j
695                  kdiscEs = 1
696                  ldiscEs = ilam
697                ENDIF
698                DO k=1,darwin_radtrans_kmax-1
699                 IF(ABS(Estop(ilam,k+1)-Esz(ilam,k)) .GT. discEs)THEN
700                  discEs = ABS(Estop(ilam,k+1)-Esz(ilam,k))
701                  idiscEs = i
702                  jdiscEs = j
703                  kdiscEs = k+1
704                  ldiscEs = ilam
705                 ENDIF
706                 IF(ABS(Eutop(ilam,k+1)-Euz(ilam,k)) .GT. discEu)THEN
707                  discEu = ABS(Eutop(ilam,k+1)-Euz(ilam,k))
708                  idiscEu = i
709                  jdiscEu = j
710                  kdiscEu = k+1
711                  ldiscEu = ilam
712                 ENDIF
713                ENDDO
714               ENDIF
715               ENDDO
716               ENDIF
717  #endif  #endif
718             ENDIF
719  c  c
720  c uses chl from prev timestep (as wavebands does)  c uses chl from prev timestep (as wavebands does)
721  c keep like this in case need to consider upwelling irradiance as affecting the grid box above  c keep like this in case need to consider upwelling irradiance as affecting the grid box above
# Line 657  c brute force... Line 762  c brute force...
762               psil  = max(Ptr(i,j,k,bi,bj,iPOSi ),0. _d 0)               psil  = max(Ptr(i,j,k,bi,bj,iPOSi ),0. _d 0)
763               NH4l  = max(Ptr(i,j,k,bi,bj,iNH4  ),0. _d 0)               NH4l  = max(Ptr(i,j,k,bi,bj,iNH4  ),0. _d 0)
764               NO2l  = max(Ptr(i,j,k,bi,bj,iNO2  ),0. _d 0)               NO2l  = max(Ptr(i,j,k,bi,bj,iNO2  ),0. _d 0)
765    #ifdef ALLOW_CDOM
766                 cdoml = max(Ptr(i,j,k,bi,bj,iCDOM  ),0. _d 0)
767    #endif
768  #ifdef ALLOW_CARBON  #ifdef ALLOW_CARBON
769               dicl  = max(Ptr(i,j,k,bi,bj,iDIC  ),0. _d 0)               dicl  = max(Ptr(i,j,k,bi,bj,iDIC  ),0. _d 0)
770               docl  = max(Ptr(i,j,k,bi,bj,iDOC  ),0. _d 0)               docl  = max(Ptr(i,j,k,bi,bj,iDOC  ),0. _d 0)
# Line 727  c ratio of maximum species Line 835  c ratio of maximum species
835                      Diver4(i,j,k)=Diver4(i,j,k)+1. _d 0                      Diver4(i,j,k)=Diver4(i,j,k)+1. _d 0
836                   endif                   endif
837                 enddo                 enddo
838    c           totphy > thresh0
839                endif
840    c Shannon and Simpson indices
841                Shannon(i,j,k) = 0. _d 0
842    c note: minimal valid value is 1, but we set to zero below threshold
843                Simpson(i,j,k) = 0. _d 0
844                if (totphy.gt.shannon_thresh) then
845                  do np=1,npmax
846                   if (Phy(np) .gt. 0. _d 0) then
847                    tmpphy(np) = Phy(np)/totphy
848                    Shannon(i,j,k)=Shannon(i,j,k)+tmpphy(np)*LOG(tmpphy(np))
849                    Simpson(i,j,k)=Simpson(i,j,k)+tmpphy(np)*tmpphy(np)
850                   endif
851                  enddo
852                  Shannon(i,j,k) = -Shannon(i,j,k)
853                  Simpson(i,j,k) = 1./Simpson(i,j,k)
854              endif              endif
855  #endif  #endif
856    
# Line 763  c for explicit sinking of particulate ma Line 887  c for explicit sinking of particulate ma
887               Slocal = salt(i,j,k,bi,bj)               Slocal = salt(i,j,k,bi,bj)
888  #endif  #endif
889    
890    c choice where to get pCO2 from
891    c taking from igsm dic run - fed through Tflux array
892    c               pCO2local=surfaceForcingT(i,j,bi,bj)
893    c or from darwin carbon module
894    #ifdef ALLOW_CARBON
895                   pCO2local=pCO2(i,j,bi,bj)
896    #else
897                   pCO2local=280. _d -6
898    #endif
899    
900               freefu = max(freefe(i,j,k),0. _d 0)               freefu = max(freefe(i,j,k),0. _d 0)
901               if (k.eq.1) then               if (k.eq.1) then
902                 inputFel = inputFe(i,j,bi,bj)                 inputFel = inputFe(i,j,bi,bj)
# Line 807  c set tendencies to 0 Line 941  c set tendencies to 0
941                 dphychl(np)=0. _d 0                 dphychl(np)=0. _d 0
942               enddo               enddo
943  #endif  #endif
944    #ifdef ALLOW_CDOM
945                dcdoml=0. _d 0
946    #endif
947  #ifdef ALLOW_CARBON  #ifdef ALLOW_CARBON
948               ddicl=0. _d 0               ddicl=0. _d 0
949               ddocl=0. _d 0               ddocl=0. _d 0
# Line 834  c set other arguments to zero Line 971  c set other arguments to zero
971                  NfixPl(np)=0. _d 0                  NfixPl(np)=0. _d 0
972  #endif  #endif
973  #endif  #endif
974    #ifdef DAR_DIAG_PARW
975                    chl2cl(np)=0. _d 0
976    #endif
977    #ifdef DAR_DIAG_EK
978                    Ekl(np)=0. _d 0
979                    EkoverEl(np)=0. _d 0
980                    do ilam=1,tlam
981                      Ek_nll(np,ilam)=0. _d 0
982                      EkoverE_nll(np,ilam)=0. _d 0
983                    enddo
984    #endif
985               enddo               enddo
986    
987    
# Line 876  c ANNA pass extra variables if WAVEBANDS Line 1024  c ANNA pass extra variables if WAVEBANDS
1024       I                       pofeupl, psiupl,       I                       pofeupl, psiupl,
1025       I                       PARl,       I                       PARl,
1026       I                       Tlocal, Slocal,       I                       Tlocal, Slocal,
1027         I                       pCO2local,
1028       I                       freefu, inputFel,       I                       freefu, inputFel,
1029       I                       bottom, dzlocal,       I                       bottom, dzlocal,
1030       O                       Rstarl, RNstarl,       O                       Rstarl, RNstarl,
# Line 902  c ANNA pass extra variables if WAVEBANDS Line 1051  c ANNA pass extra variables if WAVEBANDS
1051  #endif  #endif
1052  #ifdef GEIDER  #ifdef GEIDER
1053       O                       phychl,       O                       phychl,
1054    #ifdef DAR_DIAG_EK
1055         I                       Ekl, EkoverEl,
1056    #endif
1057    #ifdef DAR_DIAG_PARW
1058         I                       chl2cl,
1059    #endif
1060  #ifdef DYNAMIC_CHL  #ifdef DYNAMIC_CHL
1061       I                       dphychl,       I                       dphychl,
1062       I                       chlup,       I                       chlup,
1063    #ifdef DAR_DIAG_EK
1064         O                       accliml,
1065    #endif
1066    #endif
1067    #ifdef ALLOW_CDOM
1068         O                       dcdoml,
1069         I                       cdoml,
1070  #endif  #endif
1071  #ifdef WAVEBANDS  #ifdef WAVEBANDS
1072       I                       PARw_k(1,k),       I                       PARw_k(1,k),
1073    #ifdef DAR_DIAG_EK
1074         I                       Ek_nll, EkoverE_nll,
1075    #endif
1076  #endif  #endif
1077  #endif  #endif
1078  #ifdef ALLOW_PAR_DAY  #ifdef ALLOW_PAR_DAY
# Line 931  c            endif Line 1096  c            endif
1096  c  c
1097  #ifdef IRON_SED_SOURCE  #ifdef IRON_SED_SOURCE
1098  c only above minimum depth (continental shelf)  c only above minimum depth (continental shelf)
1099               if (rF(k).lt.depthfesed) then               if (rF(k).gt.-depthfesed) then
1100  c only if bottom layer  c only if bottom layer
1101                 if (bottom.eq.1.0 _d 0) then                 if (bottom.eq.1.0 _d 0) then
1102  #ifdef IRON_SED_SOURCE_VARIABLE  #ifdef IRON_SED_SOURCE_VARIABLE
# Line 966  c Line 1131  c
1131               picupl = PICl               picupl = PICl
1132  c include surface forcing  c include surface forcing
1133               if (k.eq.1) then               if (k.eq.1) then
1134                ddicl  =  ddicl  + flxCO2(i,j,bi,bj)                ddicl  =  ddicl  + flxCO2(i,j)
1135                dalkl  =  dalkl  + flxALK(i,j,bi,bj)                dalkl  =  dalkl  + flxALK(i,j)
1136                do2l   =  do2l   + flxO2(i,j,bi,bj)                do2l   =  do2l   + flxO2(i,j)
1137               endif               endif
1138  #endif  #endif
1139  c  c
# Line 1112  c in darwin_plankton this is stored for Line 1277  c in darwin_plankton this is stored for
1277  #endif  #endif
1278  #endif  #endif
1279           ENDDO           ENDDO
1280    #ifdef ALLOW_CDOM
1281              Ptr(i,j,k,bi,bj,iCDOM ) = Ptr(i,j,k,bi,bj,iCDOM ) +
1282         &                                  dtplankton*dcdoml
1283    #endif
1284  #ifdef ALLOW_CARBON  #ifdef ALLOW_CARBON
1285            Ptr(i,j,k,bi,bj,iDIC ) = Ptr(i,j,k,bi,bj,iDIC ) +            Ptr(i,j,k,bi,bj,iDIC ) = Ptr(i,j,k,bi,bj,iDIC ) +
1286       &                                  dtplankton*ddicl       &                                  dtplankton*ddicl
# Line 1242  c    &                       deltaTclock Line 1411  c    &                       deltaTclock
1411       &                           phychl(np)*dtplankton       &                           phychl(np)*dtplankton
1412               enddo               enddo
1413  #endif  #endif
1414    #ifdef DAR_DIAG_PARW
1415                do ilam=1,tlam
1416                   PARwave(i,j,k,bi,bj,ilam)=PARwave(i,j,k,bi,bj,ilam)+
1417         &                           PARw_k(ilam,k)*dtplankton
1418                enddo
1419                do np=1,npmax
1420                  chl2cave(i,j,k,bi,bj,np)=chl2cave(i,j,k,bi,bj,np)+
1421         &                          chl2cl(np)*dtplankton
1422                enddo
1423    #endif
1424  #ifdef DAR_DIAG_ACDOM  #ifdef DAR_DIAG_ACDOM
1425  c            print*,'acdom',k,acdom_k(k,darwin_diag_acdom_ilam)  c            print*,'acdom',k,acdom_k(k,darwin_diag_acdom_ilam)
1426               aCDOMave(i,j,k,bi,bj)=aCDOMave(i,j,k,bi,bj)+               aCDOMave(i,j,k,bi,bj)=aCDOMave(i,j,k,bi,bj)+
# Line 1263  Coj            no Eu at surface (yet) Line 1442  Coj            no Eu at surface (yet)
1442                 Euave(i,j,k,bi,bj,ilam)=Euave(i,j,k,bi,bj,ilam)+                 Euave(i,j,k,bi,bj,ilam)=Euave(i,j,k,bi,bj,ilam)+
1443       &                                 Euz(ilam,k-1)*dtplankton       &                                 Euz(ilam,k-1)*dtplankton
1444                endif                endif
1445                  Estave(i,j,k,bi,bj,ilam)=Estave(i,j,k,bi,bj,ilam)+
1446         &                                 Estop(ilam,k)*dtplankton
1447                Eutave(i,j,k,bi,bj,ilam)=Eutave(i,j,k,bi,bj,ilam)+                Eutave(i,j,k,bi,bj,ilam)=Eutave(i,j,k,bi,bj,ilam)+
1448       &                                 Eutop(ilam,k)*dtplankton       &                                 Eutop(ilam,k)*dtplankton
1449               enddo               enddo
1450  #endif  #endif
1451    #ifdef DAR_DIAG_IRR_AMPS
1452                 do ilam = 1,tlam
1453                   amp1ave(i,j,k,bi,bj,ilam)=amp1ave(i,j,k,bi,bj,ilam)+
1454         &                                 amp1(ilam,k)*dtplankton
1455                   amp2ave(i,j,k,bi,bj,ilam)=amp2ave(i,j,k,bi,bj,ilam)+
1456         &                                 amp2(ilam,k)*dtplankton
1457                 enddo
1458    #endif
1459  #ifdef DAR_DIAG_ABSORP  #ifdef DAR_DIAG_ABSORP
1460               do ilam = 1,tlam               do ilam = 1,tlam
1461                 aave(i,j,k,bi,bj,ilam)=aave(i,j,k,bi,bj,ilam)+                 aave(i,j,k,bi,bj,ilam)=aave(i,j,k,bi,bj,ilam)+
# Line 1291  Coj            no Eu at surface (yet) Line 1480  Coj            no Eu at surface (yet)
1480       &                                 bbpart_k(k,ilam)*dtplankton       &                                 bbpart_k(k,ilam)*dtplankton
1481               enddo               enddo
1482  #endif  #endif
1483    #ifdef DAR_RADTRANS
1484                 if (k.eq.1) then
1485                   rmudave(i,j,bi,bj)=rmudave(i,j,bi,bj)+
1486         &                                 rmud*dtplankton
1487                 endif
1488    #endif
1489    #ifdef DAR_DIAG_EK
1490                do np=1,npmax
1491                 Ekave(i,j,k,bi,bj,np)=Ekave(i,j,k,bi,bj,np)+
1492         &                        Ekl(np)*dtplankton
1493                 EkoverEave(i,j,k,bi,bj,np)=EkoverEave(i,j,k,bi,bj,np)+
1494         &                        EkoverEl(np)*dtplankton
1495                 acclimave(i,j,k,bi,bj,np)=acclimave(i,j,k,bi,bj,np)+
1496         &                        accliml(np)*dtplankton
1497                 do ilam=1,tlam
1498                    Ek_nlave(i,j,k,bi,bj,np,ilam)=
1499         &                        Ek_nlave(i,j,k,bi,bj,np,ilam)+
1500         &                        Ek_nll(np,ilam)*dtplankton
1501                    EkoverE_nlave(i,j,k,bi,bj,np,ilam)=
1502         &                        EkoverE_nlave(i,j,k,bi,bj,np,ilam)+
1503         &                        EkoverE_nll(np,ilam)*dtplankton
1504                 enddo
1505                enddo
1506    #endif
1507  #ifdef DAR_DIAG_RSTAR  #ifdef DAR_DIAG_RSTAR
1508               do np=1,npmax               do np=1,npmax
1509                 Rstarave(i,j,k,bi,bj,np)=Rstarave(i,j,k,bi,bj,np)+                 Rstarave(i,j,k,bi,bj,np)=Rstarave(i,j,k,bi,bj,np)+
# Line 1331  Coj            no Eu at surface (yet) Line 1544  Coj            no Eu at surface (yet)
1544  #ifdef ALLOW_CARBON  #ifdef ALLOW_CARBON
1545               if (k.eq.1) then               if (k.eq.1) then
1546                SURave(i,j,bi,bj)    =SURave(i,j,bi,bj)+                SURave(i,j,bi,bj)    =SURave(i,j,bi,bj)+
1547       &                              flxCO2(i,j,bi,bj)*dtplankton       &                              flxCO2(i,j)*dtplankton
1548                SURCave(i,j,bi,bj)    =SURCave(i,j,bi,bj)+                SURCave(i,j,bi,bj)    =SURCave(i,j,bi,bj)+
1549       &                              FluxCO2(i,j,bi,bj)*dtplankton       &                              FluxCO2(i,j,bi,bj)*dtplankton
1550                SUROave(i,j,bi,bj)   =SUROave(i,j,bi,bj)+                SUROave(i,j,bi,bj)   =SUROave(i,j,bi,bj)+
1551       &                              flxO2(i,j,bi,bj)*dtplankton       &                              flxO2(i,j)*dtplankton
1552                pCO2ave(i,j,bi,bj)   =pCO2ave(i,j,bi,bj)+                pCO2ave(i,j,bi,bj)   =pCO2ave(i,j,bi,bj)+
1553       &                              pCO2(i,j,bi,bj)*dtplankton       &                              pCO2(i,j,bi,bj)*dtplankton
1554                pHave(i,j,bi,bj)     =pHave(i,j,bi,bj)+                pHave(i,j,bi,bj)     =pHave(i,j,bi,bj)+
# Line 1396  C       reset the other slot for averagi Line 1609  C       reset the other slot for averagi
1609  C itistime  C itistime
1610  #endif  #endif
1611    
1612    #ifdef DAR_CHECK_IRR_CONT
1613           i = myXGlobalLo-1+(bi-1)*sNx+idiscEs
1614           j = myYGlobalLo-1+(bj-1)*sNy+jdiscEs
1615           write(6,'(I4.4,X,A,4(X,I4),1PE24.16)')myProcId,'max Es disc',
1616         &                                   i,j,kdiscEs,ldiscEs,discEs
1617           i = myXGlobalLo-1+(bi-1)*sNx+idiscEu
1618           j = myYGlobalLo-1+(bj-1)*sNy+jdiscEu
1619           write(6,'(I4.4,X,A,4(X,I4),1PE24.16)')myProcId,'max Eu disc',
1620         &                                   i,j,kdiscEu,ldiscEu,discEu
1621    #endif
1622    
1623  COJ fill diagnostics  COJ fill diagnostics
1624  #ifdef ALLOW_DIAGNOSTICS  #ifdef ALLOW_DIAGNOSTICS
1625         IF ( useDiagnostics ) THEN         IF ( useDiagnostics ) THEN
# Line 1438  c ANNA end TAVE Line 1662  c ANNA end TAVE
1662          WRITE(diagname,'(A8)') 'Diver4  '          WRITE(diagname,'(A8)') 'Diver4  '
1663          CALL DIAGNOSTICS_FILL( Diver4(1-Olx,1-Oly,1), diagname,          CALL DIAGNOSTICS_FILL( Diver4(1-Olx,1-Oly,1), diagname,
1664       &                         0,Nr,2,bi,bj,myThid )       &                         0,Nr,2,bi,bj,myThid )
1665            WRITE(diagname,'(A8)') 'Shannon '
1666            CALL DIAGNOSTICS_FILL( Shannon(1-Olx,1-Oly,1), diagname,
1667         &                         0,Nr,2,bi,bj,myThid )
1668            WRITE(diagname,'(A8)') 'Simpson '
1669            CALL DIAGNOSTICS_FILL( Simpson(1-Olx,1-Oly,1), diagname,
1670         &                         0,Nr,2,bi,bj,myThid )
1671  #endif  #endif
1672  #ifdef ALLOW_DIAZ  #ifdef ALLOW_DIAZ
1673  #ifdef DAR_DIAG_NFIXP  #ifdef DAR_DIAG_NFIXP
# Line 1463  c ANNA end TAVE Line 1693  c ANNA end TAVE
1693       &                         0,Nr,2,bi,bj,myThid )       &                         0,Nr,2,bi,bj,myThid )
1694  #endif  #endif
1695  #ifdef ALLOW_CARBON  #ifdef ALLOW_CARBON
1696          CALL DIAGNOSTICS_FILL( flxCO2(1-Olx,1-Oly,bi,bj), 'DICTFLX ',          CALL DIAGNOSTICS_FILL( flxCO2(1-Olx,1-Oly), 'DICTFLX ',
1697       &                         0,1,2,bi,bj,myThid )       &                         0,1,2,bi,bj,myThid )
1698          CALL DIAGNOSTICS_FILL( FluxCO2(1-Olx,1-Oly,bi,bj), 'DICCFLX ',          CALL DIAGNOSTICS_FILL( FluxCO2(1-Olx,1-Oly,bi,bj), 'DICCFLX ',
1699       &                         0,1,2,bi,bj,myThid )       &                         0,1,2,bi,bj,myThid )
1700          CALL DIAGNOSTICS_FILL( flxO2(1-Olx,1-Oly,bi,bj), 'DICOFLX ',          CALL DIAGNOSTICS_FILL( flxO2(1-Olx,1-Oly), 'DICOFLX ',
1701       &                         0,1,2,bi,bj,myThid )       &                         0,1,2,bi,bj,myThid )
1702          CALL DIAGNOSTICS_FILL( pCO2(1-Olx,1-Oly,bi,bj), 'DICPCO2 ',          CALL DIAGNOSTICS_FILL( pCO2(1-Olx,1-Oly,bi,bj), 'DICPCO2 ',
1703       &                         0,1,2,bi,bj,myThid )       &                         0,1,2,bi,bj,myThid )
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.14
  ViewVC Help
Powered by ViewVC 1.1.22