ocean_inversion_project/code/ptracers_init.F

C $Header: /usr/local/gcmpack/MITgcm_contrib/ocean_inversion_project/code/ptracers_init.F,v 1.3 2003/09/27 07:45:51 dimitri Exp $
C $Name:  $

#include "PTRACERS_OPTIONS.h"

CBOP
C !ROUTINE: PTRACERS_INIT

C !INTERFACE: ==========================================================
      SUBROUTINE PTRACERS_INIT( myThid )

C !DESCRIPTION:
C     Initialize PTRACERS data structures
C     This file is customized to compute CO2 perturbations from
C     30 ocean regions for Gruber's ocean inversion project.

C !USES: ===============================================================
      IMPLICIT NONE
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "PTRACERS.h"


C !INPUT PARAMETERS: ===================================================
C  myThid               :: thread number
      INTEGER myThid

C !OUTPUT PARAMETERS: ==================================================
C  none

#ifdef ALLOW_PTRACERS

C !LOCAL VARIABLES: ====================================================
C  i,j,k,bi,bj,iTracer  :: loop indices
      INTEGER i,j,k,bi,bj,iTracer,iMonth,iUnit,iRec
      CHARACTER*(10) suff
      _RL SumPtracer
CEOP
        
C Loop over tracers
      DO iTracer = 1, PTRACERS_num

C Loop over tiles
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)

C Initialize arrays in common blocks :
         DO k=1,Nr
          DO j=1-Oly,sNy+OLy
           DO i=1-Olx,sNx+Olx
            pTracer(i,j,k,bi,bj,iTracer) = 0. _d 0
            gPtr(i,j,k,bi,bj,iTracer)    = 0. _d 0
            gPtrNM1(i,j,k,bi,bj,iTracer) = 0. _d 0
           ENDDO
          ENDDO
         ENDDO
         DO j=1-Oly,sNy+OLy
          DO i=1-Olx,sNx+Olx
           surfaceTendencyPtr(i,j,bi,bj,iTracer) = 0. _d 0
          ENDDO
         ENDDO

C end bi,bj loops
        ENDDO
       ENDDO

C end of Tracer loop
      ENDDO

C Read from a pickup file if nIter0 is not zero
      IF (nIter0.NE.0) THEN
C--     Suffix for pickup files
       IF (pickupSuff.EQ.' ') THEN
         WRITE(suff,'(I10.10)') nIter0
       ELSE
         WRITE(suff,'(A10)') pickupSuff
       ENDIF
       CALL PTRACERS_READ_CHECKPOINT( nIter0,myThid )
      ENDIF

C     Initialize pTracerMasks
      CALL PTRACERS_READ_MASK ( mythid )

C     Initialize pTracerTakahashi
      CALL PTRACERS_READ_Takahashi ( mythid )

C     Normalize pTracerTakahashi so that 1e18 mol/yr is released
C     from each model region defined in pTracerMasks.
C     It is assumed that each year is 365.25 days (31557600 s)
C     long and that each month is 2629800 s.

      DO iTracer = 1, PTRACERS_num
       SumPtracer = 0.
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)
         DO j=1,sNy
          DO i=1,sNx
           DO iMonth=1,12
            SumPtracer = SumPtracer + 2629800. * rA(I,J,bi,bj) *
     &             pTracerMasks    (I,J,iTracer,bi,bj) *
     &             pTracerTakahashi(I,J,iMonth ,bi,bj)
           ENDDO
          ENDDO
         ENDDO
        ENDDO
       ENDDO
       _GLOBAL_SUM_R8( SumPtracer, myThid )
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)
         DO j=1,sNy
          DO i=1,sNx
           DO iMonth=1,12
            IF ( pTracerMasks(I,J,iTracer,bi,bj) .eq. 1. )
     &             pTracerTakahashi(I,J,iMonth ,bi,bj) = 1.e18 *
     &             pTracerTakahashi(I,J,iMonth ,bi,bj) / SumPtracer
           ENDDO
          ENDDO
         ENDDO
        ENDDO
       ENDDO
      ENDDO

#ifdef OCEAN_INVERSION_PROJECT_TIME_DEPENDENT
C     Initialize atmospheric CO2 array
      call mdsfindunit( iUnit, mythid)
      open(iUnit,file='splco2_cis92a.dat',status='old',form="FORMATTED")
      do iRec=1,pTracerAtm_Nrec
         read(iUnit,*) pTracerAtmYear(iRec), pTracerAtmCO2(iRec)
cdb      print*, '###', iRec, pTracerAtmYear(iRec), pTracerAtmCO2(iRec)
      enddo
      close(iUnit)
#endif /* OCEAN_INVERSION_PROJECT_TIME_DEPENDENT */

#endif /* ALLOW_PTRACERS */

      RETURN
      END
1	dimitri	1.4	C $Header: /usr/local/gcmpack/MITgcm_contrib/ocean_inversion_project/code/ptracers_init.F,v 1.3 2003/09/27 07:45:51 dimitri Exp $
2	dimitri	1.1	C $Name: $
3
4			#include "PTRACERS_OPTIONS.h"
5
6			CBOP
7			C !ROUTINE: PTRACERS_INIT
8
9			C !INTERFACE: ==========================================================
10			SUBROUTINE PTRACERS_INIT( myThid )
11
12			C !DESCRIPTION:
13			C Initialize PTRACERS data structures
14	dimitri	1.2	C This file is customized to compute CO2 perturbations from
15			C 30 ocean regions for Gruber's ocean inversion project.
16	dimitri	1.1
17			C !USES: ===============================================================
18			IMPLICIT NONE
19			#include "SIZE.h"
20			#include "EEPARAMS.h"
21			#include "PARAMS.h"
22			#include "GRID.h"
23			#include "PTRACERS.h"
24
25
26			C !INPUT PARAMETERS: ===================================================
27			C myThid :: thread number
28			INTEGER myThid
29
30			C !OUTPUT PARAMETERS: ==================================================
31			C none
32
33			#ifdef ALLOW_PTRACERS
34
35			C !LOCAL VARIABLES: ====================================================
36			C i,j,k,bi,bj,iTracer :: loop indices
37	dimitri	1.4	INTEGER i,j,k,bi,bj,iTracer,iMonth,iUnit,iRec
38	dimitri	1.1	CHARACTER*(10) suff
39	dimitri	1.2	_RL SumPtracer
40	dimitri	1.1	CEOP
41
42			C Loop over tracers
43			DO iTracer = 1, PTRACERS_num
44
45			C Loop over tiles
46	dimitri	1.2	DO bj = myByLo(myThid), myByHi(myThid)
47			DO bi = myBxLo(myThid), myBxHi(myThid)
48	dimitri	1.1
49			C Initialize arrays in common blocks :
50	dimitri	1.2	DO k=1,Nr
51			DO j=1-Oly,sNy+OLy
52			DO i=1-Olx,sNx+Olx
53			pTracer(i,j,k,bi,bj,iTracer) = 0. _d 0
54			gPtr(i,j,k,bi,bj,iTracer) = 0. _d 0
55			gPtrNM1(i,j,k,bi,bj,iTracer) = 0. _d 0
56			ENDDO
57	dimitri	1.3	ENDDO
58			ENDDO
59			DO j=1-Oly,sNy+OLy
60			DO i=1-Olx,sNx+Olx
61			surfaceTendencyPtr(i,j,bi,bj,iTracer) = 0. _d 0
62	dimitri	1.1	ENDDO
63			ENDDO
64
65			C end bi,bj loops
66	dimitri	1.2	ENDDO
67	dimitri	1.1	ENDDO
68
69			C end of Tracer loop
70			ENDDO
71
72	dimitri	1.2	C Read from a pickup file if nIter0 is not zero
73			IF (nIter0.NE.0) THEN
74	dimitri	1.1	C-- Suffix for pickup files
75			IF (pickupSuff.EQ.' ') THEN
76			WRITE(suff,'(I10.10)') nIter0
77			ELSE
78			WRITE(suff,'(A10)') pickupSuff
79			ENDIF
80			CALL PTRACERS_READ_CHECKPOINT( nIter0,myThid )
81			ENDIF
82
83	dimitri	1.2	C Initialize pTracerMasks
84	dimitri	1.1	CALL PTRACERS_READ_MASK ( mythid )
85	dimitri	1.2
86			C Initialize pTracerTakahashi
87			CALL PTRACERS_READ_Takahashi ( mythid )
88
89			C Normalize pTracerTakahashi so that 1e18 mol/yr is released
90			C from each model region defined in pTracerMasks.
91			C It is assumed that each year is 365.25 days (31557600 s)
92			C long and that each month is 2629800 s.
93
94			DO iTracer = 1, PTRACERS_num
95			SumPtracer = 0.
96			DO bj = myByLo(myThid), myByHi(myThid)
97			DO bi = myBxLo(myThid), myBxHi(myThid)
98			DO j=1,sNy
99			DO i=1,sNx
100			DO iMonth=1,12
101			SumPtracer = SumPtracer + 2629800. * rA(I,J,bi,bj) *
102			& pTracerMasks (I,J,iTracer,bi,bj) *
103			& pTracerTakahashi(I,J,iMonth ,bi,bj)
104			ENDDO
105			ENDDO
106			ENDDO
107			ENDDO
108			ENDDO
109			_GLOBAL_SUM_R8( SumPtracer, myThid )
110			DO bj = myByLo(myThid), myByHi(myThid)
111			DO bi = myBxLo(myThid), myBxHi(myThid)
112			DO j=1,sNy
113			DO i=1,sNx
114			DO iMonth=1,12
115			IF ( pTracerMasks(I,J,iTracer,bi,bj) .eq. 1. )
116			& pTracerTakahashi(I,J,iMonth ,bi,bj) = 1.e18 *
117			& pTracerTakahashi(I,J,iMonth ,bi,bj) / SumPtracer
118			ENDDO
119			ENDDO
120			ENDDO
121			ENDDO
122			ENDDO
123			ENDDO
124	dimitri	1.4
125			#ifdef OCEAN_INVERSION_PROJECT_TIME_DEPENDENT
126			C Initialize atmospheric CO2 array
127			call mdsfindunit( iUnit, mythid)
128			open(iUnit,file='splco2_cis92a.dat',status='old',form="FORMATTED")
129			do iRec=1,pTracerAtm_Nrec
130			read(iUnit,*) pTracerAtmYear(iRec), pTracerAtmCO2(iRec)
131			cdb print*, '###', iRec, pTracerAtmYear(iRec), pTracerAtmCO2(iRec)
132			enddo
133			close(iUnit)
134			#endif /* OCEAN_INVERSION_PROJECT_TIME_DEPENDENT */
135	dimitri	1.1
136			#endif /* ALLOW_PTRACERS */
137
138			RETURN
139			END