source: trunk/SOURCES/source_3.20/ice_bio_interp_phy2bio.f @ 4

Last change on this file since 4 was 4, checked in by vancop, 8 years ago

initial import /Users/ioulianikolskaia/Boulot/CODES/LIM1D/ARCHIVE/TMP/LIM1D_v3.20/
File size: 9.2 KB
Line 
1      SUBROUTINE ice_bio_interp_phy2bio(kideb,kiut,nlay_i,ln_write)
2
3! This routine interpolates salinity, temperature, brine salinity, brine volume
4! on the biological grid
5! (c) Martin Vancoppenolle, May 2007
6 
7      INCLUDE 'type.com'
8      INCLUDE 'para.com'
9      INCLUDE 'const.com'
10      INCLUDE 'ice.com'
11      INCLUDE 'thermo.com'
12      INCLUDE 'bio.com'
13
14      INTEGER :: 
15     &  ji          , ! : index for space
16     &  jk          , ! : index for ice layers
17     &  jn          , ! : index for tracers
18     &  layer1      , ! : relayering index
19     &  layer2        ! : relayering index
20
21      REAL(8), DIMENSION( maxnlay ) ::
22     &  zqs         , ! : scalar content on the physical grid (input)
23     &  zqt           ! : scalar content on the physical grid (input)
24
25      REAL(8), DIMENSION( nlay_bio ) ::
26     &  zq1           ! : scalar content on the biological grid (output)
27
28      REAL(8), DIMENSION( nlay_bio , maxnlay ) ::
29     &  zweight       ! : relayering matrix
30
31      REAL(8) ::   
32     &  zaaa        , ! : dummyfactors for the computation of t_i_bio
33     &  zbbb        ,
34     &  zccc        ,
35     &  zdiscrim    ,
36     &  zsum0         ! : conservation test variable
37     &  zsum1         ! : conservation test variable
38     
39      LOGICAL ::
40     &  ln_write
41
42!=============================================================================!
43
44      IF ( ln_write ) THEN
45         WRITE(numout,*) ' *** ice_bio_interp_phy2bio : '
46         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
47      ENDIF
48
49      DO ji = kideb, kiut
50!
51!-----------------------------------------------------------------------------!
52! 1) Scalar contents
53!-----------------------------------------------------------------------------!
54!
55
56      DO layer = 1, nlay_i
57         zqs(layer)   = s_i_b(ji,layer) * deltaz_i_phy(layer)
58         zqt(layer)   = q_i_b(ji,layer) * deltaz_i_phy(layer)
59      END DO ! layer
60
61      IF ( ln_write ) THEN
62   
63!        WRITE(numout,*) ' s_i_b   : ', ( s_i_b(ji,layer1) ,
64!    &                   layer1 = 1, nlay_i )
65!        WRITE(numout,*) ' q_i_b   : ', ( q_i_b(ji,layer1) ,
66!    &                   layer1 = 1, nlay_i )
67!        WRITE(numout,*) ' t_i_b   : ', ( t_i_b(ji,layer1) ,
68!    &                   layer1 = 1, nlay_i )
69!        WRITE(numout,*) ' zqs     : ', ( zqs(layer1) ,
70!    &                   layer1 = 1, nlay_i )
71!        WRITE(numout,*) ' zqt     : ', ( zqt(layer1) ,
72!    &                   layer1 = 1, nlay_i )
73
74      ENDIF
75 
76!-----------------------------------------------------------------------------!
77! 2) Weights
78!-----------------------------------------------------------------------------!
79!
80      DO layer1 = 1, nlay_bio
81         DO layer0 = 1, nlay_i
82            zweight(layer1,layer0) = MAX ( 0.0 , 
83     &      ( MIN ( zb_i_phy(layer0) ,
84     &      zb_i_bio(layer1) ) 
85     &    - MAX ( zb_i_phy (layer0-1) , zb_i_bio(layer1-1) ) ) / 
86     &      deltaz_i_phy(layer0) )
87         END DO
88      END DO
89!
90!-----------------------------------------------------------------------------!
91! 3) Interpolation
92!-----------------------------------------------------------------------------!
93!
94      !--------------
95      ! Ice salinity
96      !--------------
97      DO layer1 = 1, nlay_bio
98         zq1(layer1) = 0.0
99         DO layer0 = 1, nlay_i
100            zq1(layer1) = zq1(layer1) + zweight(layer1,layer0) *
101     &                    zqs(layer0)
102         END DO
103      END DO
104
105      IF ( ln_write ) THEN
106!        WRITE(numout,*) ' Salt contents '
107!        WRITE(numout,*) ' zq1     : ', ( zq1(layer1) ,
108!    &                   layer1 = 1, nlay_bio )
109      ENDIF
110
111      DO layer1 = 1, nlay_bio
112         s_i_bio(layer1) = zq1(layer1) / deltaz_i_bio(layer1)
113      END DO
114
115      IF ( ln_write ) THEN
116         WRITE(numout,*) ' s_i_bio : ', ( s_i_bio(layer1) , 
117     &                   layer1 = 1, nlay_bio )
118      ENDIF
119
120      !--------------
121      ! Heat content
122      !--------------
123      DO layer1 = 1, nlay_bio
124         zq1(layer1) = 0.0
125         DO layer0 = 1, nlay_i
126            zq1(layer1) = zq1(layer1) + zweight(layer1,layer0) *
127     &                    zqt(layer0)
128         END DO
129      END DO
130
131      IF ( ln_write ) THEN
132!        WRITE(numout,*) ' Heat content '
133!        WRITE(numout,*) ' zq1     : ', ( zq1(layer1) ,
134!    &                   layer1 = 1, nlay_bio )
135      ENDIF
136
137      ! Energy of melting
138      DO layer1 = 1, nlay_bio
139         zq1(layer1) = zq1(layer1) / deltaz_i_bio(layer1)
140      END DO
141
142      ! Invert energy of melting to get temperature back
143      DO layer1 = 1, nlay_bio
144         tmelts = - tmut * s_i_bio(layer1) + tpw
145         zaaa = cpg
146         zbbb = (cpw-cpg)*(tmelts-tpw) + zq1(layer1) / rhog
147     &        - lfus
148         zccc = lfus * (tmelts-tpw)
149         zdiscrim = SQRT( zbbb*zbbb - 4.0*zaaa*zccc )
150         t_i_bio(layer1) = tpw + ( - zbbb - zdiscrim ) / (2.0*zaaa)
151      END DO
152
153      IF ( ln_write ) THEN
154         WRITE(numout,*) ' t_i_bio : ', ( t_i_bio(layer1) , 
155     &                   layer1 = 1, nlay_bio )
156      ENDIF
157
158      !--------------
159      ! Brine volume
160      !--------------
161      DO layer1 = 1, nlay_bio
162         e_i_bio(layer1) = - tmut * s_i_bio(layer1) /
163     &                     ( t_i_bio(layer1) - tpw )
164      END DO ! layer1
165
166      IF ( ln_write ) THEN
167         WRITE(numout,*) ' e_i_bio : ', ( e_i_bio(layer1) , 
168     &                   layer1 = 1, nlay_bio )
169      ENDIF
170
171!-----------------------------------------------------------------------------!
172
173      END DO ! ji
174
175!=============================================================================!
176!-- End of ice_bio_interp_phy2bio --
177 
178      END
179!
180!=============================================================================!
181!=============================================================================!
182!
183
184      SUBROUTINE ice_bio_interp_diffus(kideb,kiut,nlay_i,ln_write)
185
186      INCLUDE 'type.com'
187      INCLUDE 'para.com'
188      INCLUDE 'const.com'
189      INCLUDE 'ice.com'
190      INCLUDE 'thermo.com'
191      INCLUDE 'bio.com'
192
193      INTEGER :: 
194     &  ji          , ! : index for space
195     &  jk          , ! : index for ice layers
196     &  layer_bio   , ! :
197     &  layer_phy   , ! :
198     &  index_mem
199
200      REAL(8), DIMENSION( 0:maxnlay ) :: ! lower interface of the layer
201     &  zz_phy
202
203      REAL(8), DIMENSION( 0:nlay_bio ) :: ! lower interface of the layer
204     &  zz_bio
205
206      LOGICAL ::
207     &  ln_write
208
209!=============================================================================!
210
211      IF ( ln_write ) THEN
212         WRITE(numout,*) ' *** ice_bio_interp_diffus : '
213         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
214      ENDIF
215!
216!-----------------------------------------------------------------------------!
217! 1) Grids
218!-----------------------------------------------------------------------------!
219!
220! compute the coordinates of the interfaces of the layers
221!
222      zz_phy(0) = 0.
223      DO layer_phy = 1, nlay_i
224         zz_phy(layer_phy) = z_i_phy(layer_phy) + 
225     &                       deltaz_i_phy(layer_phy) / 2.
226      END DO
227
228      zz_bio(0) = 0.
229      DO layer_bio = 1, nlay_bio
230         zz_bio(layer_bio) = z_i_bio(layer_bio) +
231     &                       deltaz_i_bio(layer_bio) / 2.
232      END DO
233
234      IF ( ln_write ) THEN
235         WRITE(numout,*) ' zz_phy : ', ( zz_phy(layer_phy), 
236     &                                   layer_phy = 0, nlay_i )
237         WRITE(numout,*) ' zz_bio : ', ( zz_bio(layer_bio), 
238     &                                   layer_bio = 0, nlay_bio )
239      ENDIF
240
241      DO layer_bio = 1, nlay_bio - 1
242         zdist_max = 999.9
243         zdist = zdist_max
244         !WRITE(numout,*) ' '
245         !WRITE(numout,*) ' layer_bio : ', layer_bio
246         DO layer_phy = 1, nlay_i
247            zdist = MIN ( zdist, zz_bio(layer_bio) - zz_phy(layer_phy) )
248            IF ( ( zdist .GE. 0.0 ) .AND. ( zdist .LT. zdist_max ) ) 
249     &      THEN
250               index_mem = layer_phy
251            ENDIF
252         !   WRITE(numout,*) ' layer_phy : ', layer_phy
253         !   WRITE(numout,*) ' zdist : ', zdist
254         !   WRITE(numout,*) ' index_mem ', index_mem
255         END DO ! layer_phy
256         index_mem = MAX ( MIN( index_mem, nlay_i ) , 1 ) ! prevent absurd values sometimes reached in path cases
257         zdummy1 = ( diff_br(index_mem+1) - diff_br(index_mem) ) /
258     &             ( zz_phy(index_mem+1) - zz_phy(index_mem) )
259         zdummy2 = zz_bio(layer_bio) - zz_phy(index_mem)
260 
261         !WRITE(numout,*) ' End of ze loupe '
262         !WRITE(numout,*) ' index_mem : ', index_mem
263 
264         diff_br_bio(layer_bio) = diff_br(index_mem) + zdummy1*zdummy2
265 
266       END DO ! layer_bio
267 
268       diff_br_bio(nlay_bio) = diff_br(nlay_i)
269
270!     DO layer_bio = 1, nlay_bio
271!        diff_br_bio(layer_bio) = diff_br(layer_bio)
272!     END DO
273
274      IF ( ln_write ) THEN
275         WRITE(numout,*)
276         WRITE(numout,*) ' diff_br    : ', ( diff_br(layer_phy), 
277     &                   layer_phy = 1, nlay_i )
278         WRITE(numout,*) ' nlay_i   : ', nlay_i
279         WRITE(numout,*) ' nlay_bio : ', nlay_bio
280         WRITE(numout,*) ' diff_br    : ', ( diff_br(layer_phy), 
281     &                   layer_phy = 1, nlay_i )
282         WRITE(numout,*) ' diff_br_bio : ', ( diff_br_bio(layer_bio),
283     &                   layer_bio = 1, nlay_bio ) 
284      ENDIF
285!     
286!=============================================================================!
287!-- End of ice_bio_interp_diff --
288!
289      END
Note: See TracBrowser for help on using the repository browser.