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

Annotation of /MITgcm_contrib/darwin2/pkg/monod/monod_generate_mutants.F

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


Revision 1.1 - (hide annotations) (download)
Wed Apr 13 18:56:25 2011 UTC (14 years, 3 months ago) by jahn
Branch: MAIN
CVS Tags: ctrb_darwin2_ckpt64k_20130723, ctrb_darwin2_ckpt65w_20160512, ctrb_darwin2_ckpt65j_20150225, ctrb_darwin2_ckpt63l_20120405, ctrb_darwin2_ckpt66g_20170424, ctrb_darwin2_ckpt64h_20130528, ctrb_darwin2_ckpt66k_20171025, ctrb_darwin2_ckpt66n_20180118, ctrb_darwin2_ckpt62v_20110413, ctrb_darwin2_ckpt65v_20160409, ctrb_darwin2_ckpt65s_20160114, ctrb_darwin2_ckpt65_20140718, ctrb_darwin2_ckpt64m_20130820, ctrb_darwin2_ckpt66d_20170214, ctrb_darwin2_ckpt64r_20131210, ctrb_darwin2_ckpt65m_20150615, ctrb_darwin2_ckpt65q_20151118, ctrb_darwin2_ckpt65o_20150914, ctrb_darwin2_ckpt64f_20130405, ctrb_darwin2_ckpt63f_20111201, ctrb_darwin2_ckpt64a_20121116, ctrb_darwin2_ckpt65p_20151023, ctrb_darwin2_ckpt64n_20130826, ctrb_darwin2_ckpt65e_20140929, ctrb_darwin2_ckpt64o_20131024, ctrb_darwin2_ckpt64v_20140411, ctrb_darwin2_ckpt64z_20140711, ctrb_darwin2_ckpt65l_20150504, ctrb_darwin2_ckpt65z_20160929, ctrb_darwin2_ckpt65n_20150729, ctrb_darwin2_ckpt62y_20110526, ctrb_darwin2_ckpt64y_20140622, ctrb_darwin2_ckpt65d_20140915, ctrb_darwin2_ckpt64t_20140202, ctrb_darwin2_ckpt66h_20170602, ctrb_darwin2_ckpt64i_20130622, ctrb_darwin2_ckpt62x_20110513, ctrb_darwin2_ckpt64s_20140105, ctrb_darwin2_ckpt62w_20110426, ctrb_darwin2_ckpt64x_20140524, ctrb_darwin2_ckpt63o_20120629, ctrb_darwin2_ckpt64e_20130305, ctrb_darwin2_ckpt65x_20160612, ctrb_darwin2_ckpt66f_20170407, ctrb_darwin2_ckpt63c_20111011, ctrb_darwin2_ckpt63i_20120124, ctrb_darwin2_ckpt65g_20141120, ctrb_darwin2_ckpt63m_20120506, ctrb_darwin2_ckpt63s_20120908, ctrb_darwin2_ckpt65k_20150402, ctrb_darwin2_ckpt63e_20111107, ctrb_darwin2_ckpt64w_20140502, ctrb_darwin2_ckpt63b_20110830, ctrb_darwin2_ckpt63j_20120217, ctrb_darwin2_ckpt66a_20161020, ctrb_darwin2_ckpt63r_20120817, ctrb_darwin2_ckpt64g_20130503, ctrb_darwin2_ckpt64l_20130806, ctrb_darwin2_ckpt63g_20111220, ctrb_darwin2_ckpt65f_20141014, ctrb_darwin2_ckpt64c_20130120, ctrb_darwin2_ckpt63a_20110804, ctrb_darwin2_ckpt66b_20161219, ctrb_darwin2_ckpt64u_20140308, ctrb_darwin2_ckpt64j_20130704, ctrb_darwin2_ckpt65i_20150123, ctrb_darwin2_ckpt66j_20170815, ctrb_darwin2_ckpt65y_20160801, ctrb_darwin2_ckpt63h_20111230, ctrb_darwin2_ckpt63p_20120707, ctrb_darwin2_ckpt66c_20170121, ctrb_darwin2_ckpt65a_20140728, ctrb_darwin2_ckpt65b_20140812, ctrb_darwin2_ckpt65t_20160221, ctrb_darwin2_ckpt64p_20131118, ctrb_darwin2_ckpt63d_20111107, ctrb_darwin2_ckpt63q_20120731, ctrb_darwin2_ckpt63_20110728, ctrb_darwin2_ckpt64b_20121224, ctrb_darwin2_ckpt64d_20130219, ctrb_darwin2_ckpt66o_20180209, ctrb_darwin2_ckpt66e_20170314, ctrb_darwin2_ckpt64_20121012, ctrb_darwin2_ckpt64q_20131118, ctrb_darwin2_baseline, ctrb_darwin2_ckpt64p_20131024, ctrb_darwin2_ckpt65u_20160315, ctrb_darwin2_ckpt65r_20151221, ctrb_darwin2_ckpt66i_20170718, ctrb_darwin2_ckpt63n_20120604, ctrb_darwin2_ckpt63k_20120317, ctrb_darwin2_ckpt65c_20140830, ctrb_darwin2_ckpt62z_20110622, ctrb_darwin2_ckpt66l_20171025, ctrb_darwin2_ckpt65h_20141217, ctrb_darwin2_ckpt66m_20171213, HEAD
darwin2 initial checkin

1 jahn 1.1 C $Header$
2     C $Name$
3    
4     #include "CPP_OPTIONS.h"
5     #include "PTRACERS_OPTIONS.h"
6     #include "DARWIN_OPTIONS.h"
7    
8     #ifdef ALLOW_PTRACERS
9     #ifdef ALLOW_MONOD
10    
11     c ==========================================================
12     c SUBROUTINE MONOD_GENERATE_MUTANTS
13     c generate parameters for "functional group" of phyto (index np)
14     c using a "Monte Carlo" approach
15     c Mick Follows, Scott Grant Fall/Winter 2005
16     c Stephanie Dutkiewicz Spring/Summer 2006
17     c modified to set up mutants
18     c Jason Bragg Spring/Summer 2007
19     c ==========================================================
20     SUBROUTINE MONOD_GENERATE_MUTANTS(myThid, np)
21    
22     implicit none
23     #include "MONOD_SIZE.h"
24     #include "MONOD.h"
25    
26     C !INPUT PARAMETERS: ===================================================
27     C myThid :: thread number
28     INTEGER myThid
29     CEOP
30    
31     C === Functions ===
32     _RL DARWIN_RANDOM
33     EXTERNAL DARWIN_RANDOM
34    
35    
36     c local variables
37     _RL RandNo
38     _RL growthdays
39     _RL mortdays
40     _RL pday
41     _RL year
42     _RL month
43     _RL fiveday
44     _RL rtime
45     _RL standin
46     INTEGER np
47     INTEGER nz
48     INTEGER signvar
49    
50     c Mutation --------------------------------------------------------
51     c Scheme to make 'npro' Pro ecotypes, each with 'numtax' sister taxa
52     c some muatation variables [jbmodif]
53    
54     INTEGER nsis
55     INTEGER npro
56     INTEGER numtax
57     INTEGER taxind
58     INTEGER threeNmut
59    
60     #ifdef ALLOW_MUTANTS
61    
62     c initialize mutation variables [jbmodif]
63     npro = 60
64     threeNmut = 1
65     numtax = 4
66     taxind = mod(np,numtax)
67    
68     c end mutation ----------------------------------------------------
69    
70    
71     c
72     standin=0.d0
73    
74     c length of day (seconds)
75     pday = 86400.0d0
76    
77     c each time generate another functional group add one to ngroups
78     ngroups = ngroups + 1
79    
80    
81     c Mutation --------------------------------------------------------
82     c Generate random variables if 1st sis, or non-Pro [jbmodif]
83     if (np .gt. npro .or. taxind .eq. 1.0d0)then
84    
85    
86     c RANDOM NUMBERS
87     c phyto either "small" (physize(np)=0.0) or "big" (physize(np)=1.0)
88     c at this point independent of whether diatom or not
89    
90     c make Pro small, randomize others [jbmodif]
91     if(np.le.npro)physize(np) = 0.0d0
92     if(np.gt.npro)then
93    
94     RandNo = darwin_random(myThid)
95     if(RandNo .gt. 0.000d0)then
96     physize(np) = 1.0d0
97     else
98     physize(np) = 0.0d0
99     end if
100     endif
101     c
102     c phyto either diatoms (diatom=1.0) and use silica or not (diatom=0.0)
103     c if they are large
104     if (physize(np).eq.1.0d0) then
105     RandNo = darwin_random(myThid)
106     if(RandNo .gt. 0.500d0)then
107     diatom(np) = 1.0d0
108     else
109     diatom(np) = 0.0d0
110     end if
111     else
112     diatom(np) = 0.0d0
113     endif
114     c TEST ...........................................
115     c diatom(np) = 0.0d0
116     c write(6,*)'FIXED - no DIATOM '
117     c TEST ...........................................
118    
119    
120    
121     c phyto either diazotrophs (diazotroph=1.0) or not (diazotroph=0.0)
122     RandNo = darwin_random(myThid)
123     if(RandNo .gt. 0.500d0)then
124     diazotroph(np) = 1.0d0
125     else
126     diazotroph(np) = 0.0d0
127     end if
128     c TEST ...........................................
129     #ifndef ALLOW_DIAZ
130     diazotroph(np) = 0.0d0
131     write(6,*)'FIXED - no DIAZO '
132     #endif
133     c TEST ...........................................
134    
135    
136     c growth rates
137     RandNo = darwin_random(myThid)
138     c big/small phyto growth rates..
139     if(physize(np) .eq. 1.0d0)then
140     growthdays = Biggrow +Randno*Biggrowrange
141     else
142     growthdays = Smallgrow +RandNo*Smallgrowrange
143     end if
144     c but diazotrophs always slower due to energetics
145     if(diazotroph(np) .eq. 1.0) then
146     growthdays = growthdays * diaz_growfac
147     endif
148     c now convert to a growth rate
149     mu(np) = 1.0d0/(growthdays*pday)
150    
151     c mortality and export fraction rates
152     RandNo = darwin_random(myThid)
153     c big/small phyto mortality rates..
154     if(physize(np) .eq. 1.0d0)then
155     mortdays = Bigmort +Randno*Bigmortrange
156     ExportFracP(np)=Bigexport
157     else
158     mortdays = Smallmort +RandNo*Smallmortrange
159     ExportFracP(np)=Smallexport
160     end if
161     c now convert to a mortality rate
162     mortphy(np) = 1.0d0/(mortdays*pday)
163    
164    
165    
166     c nutrient source
167     c if using threeNmut=1, all have nsource=3 [jbmodif]
168    
169     if(threeNmut.eq.1)then
170     nsource(np)=3
171     else
172     if(diazotroph(np) .ne. 1.0)then
173     RandNo = darwin_random(myThid)
174     if (physize(np).eq.1.0d0) then
175     nsource(np) = 3
176     else
177     if(RandNo .gt. 0.670d0)then
178     nsource(np) = 1
179     elseif(RandNo .lt. 0.33d0)then
180     nsource(np) = 2
181     else
182     nsource(np) = 3
183     endif
184     endif
185     else
186     nsource(np) = 0
187     end if
188     endif
189    
190     c..........................................................
191     c generate phyto Temperature Function parameters
192     c.......................................................
193     phytoTempCoeff(np) = tempcoeff1
194     phytoTempExp1(np) = tempcoeff3
195     if(physize(np) .eq. 1.0d0)then
196     phytoTempExp2(np) = tempcoeff2_big
197     else
198     phytoTempExp2(np) = tempcoeff2_small
199     endif
200    
201     RandNo = darwin_random(myThid)
202     cswd phytoTempOptimum(np) = 30.0d0 - RandNo*28.0d0
203     phytoTempOptimum(np) = tempmax - RandNo*temprange
204     phytoDecayPower(np) = tempdecay
205    
206     write(6,*)'generate Phyto: np = ',np,' Topt =',
207     & phytoTempOptimum(np)
208    
209    
210     c ...............................................
211     write(6,*)'generate phyto: np = ',np,' growthdays = ',growthdays
212     c ...............................................
213    
214     c stoichiometric ratios for each functional group of phyto
215     c relative to phosphorus - the base currency nutrient
216     c set Si:P
217     if(diatom(np) .eq. 1.0)then
218     R_SiP(np) = val_R_SiP_diatom
219     else
220     R_SiP(np) = 0.0d0
221     end if
222     c set N:P and iron requirement according to diazotroph status
223     if(diazotroph(np) .eq. 1.0)then
224     R_NP(np) = val_R_NP_diaz
225     R_FeP(np) = val_RFeP_diaz
226     else
227     R_NP(np) = val_R_NP
228     R_FeP(np) = val_RFeP
229     end if
230     c set sinking rates according to allometry
231     if(physize(np) .eq. 1.0)then
232     wsink(np) = BigSink
233     else
234     wsink(np) = SmallSink
235     end if
236     c half-saturation coeffs
237    
238     RandNo = darwin_random(myThid)
239     if(physize(np) .eq. 1.0)then
240     ksatPO4(np) = BigPsat + RandNo*BigPsatrange
241     else
242     c ksatPO4(np) = SmallPsat + RandNo*SmallPsatrange
243     c if (nsource(np).lt.3) then
244     c ksatPO4(np) = ksatPO4(np)*prochlPsat
245     c endif
246     c make pro ksat if pro [jbmodif]
247     c if (nsource(np).eq.3) then
248     if (np .gt. npro)then
249     ksatPO4(np) = SmallPsat + RandNo*SmallPsatrange
250     else
251     ksatPO4(np) = ProcPsat + RandNo*ProcPsatrange
252     endif
253     endif
254     ksatNO3(np) = ksatPO4(np)*R_NP(np)
255     ksatNO2(np) = ksatNO3(np)*ksatNO2fac
256     c Made ksatNH4 smaller since it is the preferred source
257     ksatNH4(np) = ksatNO3(np)*ksatNH4fac
258     ksatFeT(np) = ksatPO4(np)*R_FeP(np)
259     ksatSi(np) = val_ksatsi
260    
261     cNEW Light parameters:
262     c ksatPAR {0.1 - 1.3}
263     c 0.35=Av High Light Adapted, 0.8=Av Low Light Adapted
264     c kinhib {0.0 - 3.0}
265     c 0.5 =Av High Light Adapted, 2.0=Av Low Light Adapted
266     c High Light Groups for Large size:
267     if(physize(np) .eq. 1.0d0)then
268     RandNo = darwin_random(myThid)
269     call invnormal(standin,RandNo,Bigksatpar,Bigksatparstd)
270     ksatPAR(np) = abs(standin)
271    
272     RandNo = darwin_random(myThid)
273     CALL invnormal(standin,RandNo,Bigkinhib,Bigkinhibstd)
274     kinhib(np) = abs(standin)
275     else
276     c QQ remove someday
277     RandNo = darwin_random(myThid)
278     c Low Light Groups for Small size:
279     RandNo = darwin_random(myThid)
280     CALL invnormal(standin,RandNo,smallksatpar,
281     & smallksatparstd)
282     ksatPAR(np) = abs(standin)
283    
284     RandNo = darwin_random(myThid)
285     CALL invnormal(standin,RandNo,smallkinhib,
286     & smallkinhibstd)
287     kinhib(np) = abs(standin)
288     endif
289     write(6,*)'generate Phyto: np = ',np,' ksatPAR, kinhib =',
290     & ksatPAR(np), kinhib(np)
291    
292     #ifndef OLD_GRAZE
293     c for zooplankton
294     c assume zoo(1) = small, zoo(2) = big
295     zoosize(1) = 0.0d0
296     zoosize(2) = 1.0d0
297     grazemax(1) = GrazeFast
298     grazemax(2) = GrazeFast
299     ExportFracZ(1)=ZooexfacSmall
300     ExportFracZ(2)=ZooexfacBig
301     mortzoo(1) = ZoomortSmall
302     mortzoo(2) = ZoomortBig
303     ExportFracGraz(1)=ExGrazFracbig
304     ExportFracGraz(2)=ExGrazFracsmall
305     IF ( nzmax.GT.2 ) THEN
306     WRITE(msgBuf,'(2A,I5)') 'MONOD_GENERATE_MUTANTS: ',
307     & 'nzmax = ', nzmax
308     CALL PRINT_ERROR( msgBuf , 1)
309     WRITE(msgBuf,'(2A)') 'MONOD_GENERATE_MUTANTS: ',
310     & 'please provide size info for nz > 2'
311     CALL PRINT_ERROR( msgBuf , 1)
312     STOP 'ABNORMAL END: S/R MONOD_GENERATE_MUTANTS'
313     ENDIF
314     c
315     do nz=1,nzmax
316     c palatibity according to "allometry"
317     c big grazers preferentially eat big phyto etc...
318     if (zoosize(nz).eq.physize(np)) then
319     palat(np,nz)=palathi
320     asseff(np,nz)=GrazeEffmod
321     else
322     palat(np,nz)=palatlo
323     if (physize(np).eq.0.d0) then
324     asseff(np,nz)=GrazeEffhi
325     else
326     asseff(np,nz)=GrazeEfflow
327     endif
328     endif
329     c diatoms even less palatible
330     if (diatom(np).eq.1) then
331     palat(np,nz)= palat(np,nz)*diatomgraz
332     endif
333     enddo
334     #endif
335    
336     c end the if a 1st Pro or non-Pro
337     endif
338    
339     c Mutation --------------------------------------------
340     c make the Pro sister taxa
341     c -----------------------------------------------------
342    
343     if (np .le. npro .and. taxind .ne. 1)then
344    
345     c start by making mutants identical to their sister
346     if (numtax .eq. 2)then
347     if(taxind .eq. 0)nsis = np-1
348     endif
349    
350     if (numtax .eq. 3)then
351     if(taxind .eq. 2)nsis = np-1
352     if(taxind .eq. 0)nsis = np-2
353     endif
354    
355     if (numtax .eq. 4)then
356     if(taxind .eq. 2)nsis = np-1
357     if(taxind .eq. 3)nsis = np-2
358     if(taxind .eq. 0)nsis = np-3
359     endif
360    
361    
362     physize(np) = physize(nsis)
363     diatom(np) = diatom(nsis)
364     diazotroph(np) = diazotroph(nsis)
365     mu(np) = mu(nsis)
366     ExportFracP(np)=ExportFracP(nsis)
367     mortphy(np) = mortphy(nsis)
368     nsource(np) = nsource(nsis)
369     phytoTempCoeff(np) = phytoTempCoeff(nsis)
370     phytoTempExp1(np) = phytoTempExp1(nsis)
371     phytoTempExp2(np) = phytoTempExp2(nsis)
372     phytoTempOptimum(np) = phytoTempOptimum(nsis)
373     phytoDecayPower(np) = phytoDecayPower(nsis)
374     R_SiP(np) = R_SiP(nsis)
375     R_NP(np) = R_NP(nsis)
376     R_FeP(np) = R_FeP(nsis)
377     wsink(np) = wsink(nsis)
378     ksatPO4(np) = ksatPO4(nsis)
379     ksatNO3(np) = ksatNO3(nsis)
380     ksatNO2(np) = ksatNO2(nsis)
381     ksatNH4(np) = ksatNH4(nsis)
382     ksatFeT(np) = ksatFeT(nsis)
383     ksatSi(np) = ksatSi(nsis)
384     ksatPAR(np) = ksatPAR(nsis)
385     kinhib(np) = kinhib(nsis)
386    
387     #ifndef OLD_GRAZE
388     do nz=1,nzmax
389     palat(np,nz)=palat(nsis,nz)
390     asseff(np,nz)=asseff(nsis,nz)
391     enddo
392     #endif
393    
394     c then mutate
395     c here, make nsource mutants
396     if(numtax .eq. 3)then
397     if(taxind .eq. 1.0d0) nsource(np) = 3
398     if(taxind .eq. 2.0d0) nsource(np) = 2
399     if(taxind .eq. 0.0d0) nsource(np) = 1
400     endif
401    
402     if(numtax .eq. 4)then
403     if(taxind .eq. 1.0d0) nsource(np) = 3
404     if(taxind .eq. 2.0d0) nsource(np) = 2
405     if(taxind .eq. 3.0d0) nsource(np) = 1
406     if(taxind .eq. 0.0d0) nsource(np) = 3
407     endif
408    
409    
410     c here, make mu and ksat mutants
411     c if(taxind .eq. 2.0d0) mu(np) = mu(np)*1.1
412     c if(taxind .eq. 0.0d0) ksatPO4(np) = ksatPO4(np)*0.95
413    
414     endif
415    
416     #endif
417    
418    
419     RETURN
420     END
421     #endif /*MONOD*/
422     #endif /*ALLOW_PTRACERS*/
423    
424     c ===========================================================

  ViewVC Help
Powered by ViewVC 1.1.22