| 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 =========================================================== |