- Timestamp:
- 10/05/16 14:25:40 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_v3.20_2016_gravity_drainage/SOURCES/source_3.20/ice_bio_sms.f
r4 r20 82 82 !-----------------------------------------------------------------------------! 83 83 ! 84 consN = 0 85 86 DO layer = layer_00, nlay_bio 87 88 trucN = cbu_i_bio(jn_din,layer) + cbu_i_bio(jn_eon,layer) 89 & + cbu_i_bio(jn_aon,layer) + consN 90 consN = trucN 91 92 END DO 93 94 WRITE(numout,*) ' Conservation BIO N :', consN 95 84 96 85 97 ! Reference of the upper layer … … 246 258 END DO 247 259 248 ! Remineralization 249 IF ( ln_rem ) THEN 250 251 SELECT CASE ( nn_bio_opt ) 252 253 CASE(0) ! NP 254 255 rem_bio(layer_00:nlay_bio) = frem_bio * 256 & lys_bio(layer_00:nlay_bio) 257 258 CASE(1) ! NPDr 259 260 IF ( nn_rem .GE. 1 ) ! no T-dependence 261 & rem_bio(layer_00:nlay_bio) = krem_bio * 262 & cbu_i_bio(jn_eoc,layer_00:nlay_bio) 263 264 265 IF ( nn_rem .EQ. 2 ) ! T-dependence 266 & rem_bio(layer_00:nlay_bio) = rem_bio(layer_00:nlay_bio) * 267 & EXP( rg_bac * tc_bio(layer_00:nlay_bio) ) 268 269 END SELECT 270 271 ENDIF 272 260 DO layer = layer_00, nlay_bio 261 262 IF ( ln_rem ) THEN 263 264 SELECT CASE ( nn_bio_opt ) 265 266 CASE(0) ! NP 267 268 rem_bio(layer) = frem_bio * 269 & lys_bio(layer) 270 271 CASE(1) ! NPD 272 273 IF ( nn_rem .EQ. 1 ) ! no T-dependence 274 & rem_bio(layer) = krem_bio * 275 & cbu_i_bio(jn_eoc,layer) 276 277 IF ( nn_rem .EQ. 2 ) ! T-dependence 278 & rem_bio(layer) = krem_bio * cbu_i_bio(jn_eoc,layer) 279 & * EXP( rg_bac * tc_bio(layer) ) 280 281 END SELECT 282 283 ENDIF 284 285 END DO 273 286 274 287 !------------------------------------------------------------------------ 275 288 ! 2.3 Diatom Synthesis 276 289 !------------------------------------------------------------------------ 290 zeps = 0.000000001 277 291 278 292 DO layer = layer_00, nlay_bio 279 293 280 294 IF ( ln_syn ) THEN 281 295 ! raw synthesis (in s-1) … … 287 301 zsyn1 = z_syn * cbu_i_bio(jn_aoc,layer) ! normal carbon production rate (mmolC/m3/s) 288 302 zsyn2 = zsyn1; zsyn3 = zsyn1; zsyn4 = zsyn1 303 289 304 290 305 IF ( ln_lim_dsi ) ! prevent exhaustion of Si 291 306 & zsyn2 = cbu_i_bio(jn_dsi,layer) / ddtb / si_c - zeps 292 307 293 IF ( ln_lim_no3 ) ! prevent exhaustion of N 294 & zsyn3 = cbu_i_bio(jn_din,layer) / ddtb / no3_c - zeps 295 296 IF ( ln_lim_po4 ) ! prevent exhaustion of P 297 & zsyn4 = cbu_i_bio(jn_dip,layer) / ddtb / po4_c - zeps 298 299 ! nutrient stock-limited PP rate 308 IF (ln_decoupNC) THEN ! C aussi limité en stock par DIN 309 zC = MAX (cbu_i_bio(jn_aoc,layer) , zeps ) 310 zN = MAX (cbu_i_bio(jn_aon,layer) , zeps ) 311 N_C_alg(layer) = zN / zC 312 ! N_C_alg(layer) = 0.14 313 zsyn3 = cbu_i_bio(jn_din,layer) / (ddtb*N_C_alg(layer)) 314 & - zeps 315 ELSE 316 zsyn3 = cbu_i_bio(jn_din,layer) / (ddtb * no3_c) - zeps 317 ENDIF 318 319 320 IF (ln_decoupPC) THEN ! C aussi limité en stock par DIP 321 zC = MAX (cbu_i_bio(jn_aoc,layer) , zeps ) 322 zP = MAX (cbu_i_bio(jn_aop,layer) , zeps ) 323 P_C_alg(layer) = zP / zC 324 ! P_C_alg(layer) = 0.0094 325 zsyn4 = cbu_i_bio(jn_dip,layer)/(ddtb * P_C_alg(layer)) 326 & - zeps 327 ELSE 328 zsyn4 = cbu_i_bio(jn_dip,layer) / ddtb / po4_c - zeps 329 ENDIF 330 331 332 ! nutrient stock-limited PP rate 300 333 syn_bio(layer) = MIN( zsyn1, zsyn2, zsyn3, zsyn4 ) 301 334 … … 312 345 313 346 314 END DO ! layer 315 316 !------------------------------------------------------------------------ 317 ! 2.4 Update reservoirs 318 !------------------------------------------------------------------------ 319 320 DO layer = layer_00, nlay_bio 347 ! ------------------------- 348 ! -- Limitations en stocks 349 ! ------------------------ 350 351 IF (ln_decoupNC) THEN 352 zzdin = syn_bio(layer) / (((1/N_C_alg(layer)) * 353 & (cbu_i_bio(jn_din,layer)/ddtb)) - zeps) 354 355 ELSE 356 zzdin = syn_bio(layer) / (((1/no3_c) * 357 & (cbu_i_bio(jn_din,layer)/ddtb)) - zeps) 358 ENDIF 359 360 361 IF (ln_decoupPC) THEN 362 zzpo4 = syn_bio(layer) / (((1/P_C_alg(layer)) * 363 & (cbu_i_bio(jn_dip,layer)/ddtb)) - zeps) 364 365 ELSE 366 zzpo4 = syn_bio(layer) / (((1/po4_c) * 367 & (cbu_i_bio(jn_dip,layer)/ddtb)) - zeps) 368 ENDIF 369 370 371 372 zzsi = syn_bio(layer) / (((1/si_c) * 373 & (cbu_i_bio(jn_dsi,layer)/ddtb)) - zeps) 374 375 376 IF ( zzdin .EQ. 1 ) THEN 377 lim_din_stock(layer) = 1 !DIN est limitant en stocks 378 ELSE 379 lim_din_stock(layer) = 0 !DIN PAS limitant en stocks 380 ENDIF 381 382 IF ( zzpo4 .EQ. 1 ) THEN 383 lim_dip_stock(layer) = 1 384 ELSE 385 lim_dip_stock(layer) = 0 386 ENDIF 387 388 IF ( zzsi .EQ. 1) THEN 389 lim_dsi_stock(layer) = 1 390 ELSE 391 lim_dsi_stock(layer) = 0 392 ENDIF 393 394 395 !------------------------------------------------------------------------ 396 ! N cycle 397 !------------------------------------------------------------------------ 398 399 400 IF ( ln_decoupNC ) THEN 401 402 !--- Nitrogen carbon algal ratio 403 zC = MAX( cbu_i_bio(jn_aoc,layer) , zeps ) 404 zN = MAX( cbu_i_bio(jn_aon,layer) , zeps ) 405 N_C_alg(layer) = zN / zC 406 ! zN_C = MIN( cbu_i_bio(jn_aon,layer) / zC , N_C_max ) 407 ! zN_C = MAX( N_C_min , zN_C ) 408 ! N_C_alg(layer) = zN_C ! => marche pas 409 ! N_C_alg(layer) = 0.14 410 411 !--- Nitrogen carbon detritic ratio 412 zC = MAX( cbu_i_bio(jn_eoc,layer), zeps ) 413 zNd = MAX( cbu_i_bio(jn_eon,layer), zeps ) 414 N_C_det(layer) = zNd / zC 415 ! N_C_det(layer) = 0.14 416 417 !--- Process rates 418 syn_N(layer) = N_C_alg(layer) * ksyn_N * syn_bio(layer) 419 rsp_N(layer) = N_C_alg(layer) * krsp_N * rsp_bio(layer) 420 lys_N(layer) = N_C_alg(layer) * klys_N * lys_bio(layer) 421 rem_N(layer) = N_C_det(layer) * krem_N * rem_bio(layer) 422 423 424 ENDIF ! ln_decoupNC 425 426 !------------------------------------------------------------------------ 427 ! P cycle 428 !------------------------------------------------------------------------ 429 430 431 IF ( ln_decoupPC ) THEN 432 433 !--- P carbon algal ratio 434 zC = MAX( cbu_i_bio(jn_aoc,layer) , zeps ) 435 zP = MAX( cbu_i_bio(jn_aop,layer) , zeps ) 436 P_C_alg(layer) = zP / zC 437 438 !--- P carbon detritic ratio 439 zC = MAX( cbu_i_bio(jn_eoc,layer), zeps ) 440 zPd = MAX( cbu_i_bio(jn_eop,layer), zeps ) 441 P_C_det(layer) = zPd / zC 442 443 !--- Process rates 444 syn_P(layer) = P_C_alg(layer) * ksyn_P * syn_bio(layer) 445 rsp_P(layer) = P_C_alg(layer) * krsp_P * rsp_bio(layer) 446 lys_P(layer) = P_C_alg(layer) * klys_P * lys_bio(layer) 447 rem_P(layer) = P_C_det(layer) * krem_P * rem_bio(layer) 448 449 ENDIF ! ln_decoupPC 450 451 452 zPa = MAX( cbu_i_bio(jn_aop,layer) , zeps ) 453 zNt = MAX( cbu_i_bio(jn_aon,layer) , zeps ) 454 N_P_alg(layer) = zNt / zPa 455 zPe = MAX( cbu_i_bio(jn_eop,layer) , zeps ) 456 zNtd = MAX( cbu_i_bio(jn_eon,layer) , zeps ) 457 N_P_det(layer) = zNtd / zPe 458 459 460 461 !------------------------------------------------------------------------ 462 ! 2.4 Update reservoirs (including N cycle reservoirs) 463 !------------------------------------------------------------------------ 464 465 ! ------------------------------------------------------- 466 ! -- 1) Increments 467 ! ------------------------------------------------------ 468 321 469 322 470 ! DSi uptake … … 325 473 zd_dsi2 = - cbu_i_bio(jn_dsi,layer) ! Maximum DSi uptake 326 474 zd_dsi = MAX ( zd_dsi1, zd_dsi2 ) 327 328 ! NO3 uptake329 zd_no31 = no3_c * ( rem_bio(layer) +330 & rsp_bio(layer) - syn_bio(layer) ) * ddtb331 zd_no32 = - cbu_i_bio(jn_din,layer) ! Max NO3 uptake332 zd_no3 = MAX ( zd_no31, zd_no32 )333 334 ! update PO4335 zd_po41 = po4_c * ( rem_bio(layer) +336 & rsp_bio(layer) - syn_bio(layer) ) * ddtb337 zd_po42 = - cbu_i_bio(jn_dip,layer) ! Max PO4 uptake338 zd_po4 = MAX ( zd_po41, zd_po42 )339 475 340 476 ! update algae … … 350 486 zd_eoc2 = - cbu_i_bio(jn_eoc,layer) ! Max detrital carbon loss 351 487 zd_eoc = MAX( zd_eoc1, zd_eoc2 ) 352 ENDIF 488 ENDIF ! nn_bio_opt 489 490 491 IF (ln_decoupNC) THEN 492 493 !--- Increments 494 zd_din1 = ( rsp_N(layer) - syn_N(layer) + 495 & rem_N(layer) ) * ddtb 496 zd_din2 = - cbu_i_bio(jn_din,layer) !PEU PAS EPUISER PLUS 497 !QUE CE QUIL Y A !!!! 498 zd_din = MAX ( zd_din1, zd_din2 ) ! interessant dans le cas 499 ! ou zd_din1 est négatif 500 ! zd_din1 peut pas être 501 ! plus petit que zd_din2 502 503 zd_aon1 = ( syn_N(layer) - rsp_N(layer) - lys_N(layer) ) * 504 & ddtb 505 zd_aon2 = - cbu_i_bio(jn_aon,layer) 506 zd_aon = MAX ( zd_aon1, zd_aon2 ) 507 508 zd_eon1 = ( lys_N(layer) - rem_N(layer) ) * ddtb 509 zd_eon2 = - cbu_i_bio(jn_eon,layer) 510 zd_eon = MAX ( zd_eon1, zd_eon2 ) 511 512 !--- Correct rates to prevent leaks 513 IF ( zd_din .EQ. zd_din2 ) 514 & syn_N(layer) = MIN( syn_N(layer), -zd_din2/ddtb) 515 516 IF ( zd_aon .EQ. zd_aon2 ) 517 & lys_N(layer) = MIN( lys_N(layer), -zd_aon2/ddtb) 518 519 IF ( zd_eon .EQ. zd_eon2 ) 520 & rem_N(layer) = MIN( rem_N(layer), -zd_eon2/ddtb) 521 522 523 ELSE 524 zd_din1 = no3_c * ( rem_bio(layer) + 525 & rsp_bio(layer) - syn_bio(layer) ) * ddtb 526 zd_din2 = - cbu_i_bio(jn_din,layer) ! Max NO3 uptake 527 zd_din = MAX ( zd_din1, zd_din2 ) 528 529 ENDIF ! ln_decoupNC 530 531 !!!!!!!!!!!!!!!!!! 532 533 IF (ln_decoupPC) THEN 534 535 !--- Increments 536 zd_dip1 = ( rsp_P(layer) - syn_P(layer) + 537 & rem_P(layer) ) * ddtb 538 zd_dip2 = - cbu_i_bio(jn_dip,layer) !PEU PAS EPUISER PLUS 539 !QUE CE QUIL Y A !!!! 540 zd_dip = MAX ( zd_dip1, zd_dip2 ) ! interessant dans le 541 ! ou zd_din1 est négatif 542 ! zd_din1 peut pas être 543 ! plus petit que zd_din2 544 545 zd_aop1 = ( syn_P(layer) - rsp_P(layer) - lys_P(layer) ) * 546 & ddtb 547 zd_aop2 = - cbu_i_bio(jn_aop,layer) 548 zd_aop = MAX ( zd_aop1, zd_aop2 ) 549 550 zd_eop1 = ( lys_P(layer) - rem_P(layer) ) * ddtb 551 zd_eop2 = - cbu_i_bio(jn_eop,layer) 552 zd_eop = MAX ( zd_eop1, zd_eop2 ) 553 554 !--- Correct rates to prevent leaks 555 IF ( zd_dip .EQ. zd_dip2 ) 556 & syn_P(layer) = MIN( syn_P(layer), -zd_dip2/ddtb) 557 558 IF ( zd_aop .EQ. zd_aop2 ) 559 & lys_P(layer) = MIN( lys_P(layer), -zd_aop2/ddtb) 560 561 IF ( zd_eop .EQ. zd_eop2 ) 562 & rem_P(layer) = MIN( rem_P(layer), -zd_eop2/ddtb) 563 564 565 ELSE 566 zd_dip1 = po4_c * ( rem_bio(layer) + 567 & rsp_bio(layer) - syn_bio(layer) ) * ddtb 568 zd_dip2 = - cbu_i_bio(jn_dip,layer) ! Max DIP uptake 569 zd_dip = MAX ( zd_dip1, zd_dip2 ) 570 571 ENDIF ! ln_decoupPC 572 573 574 ! ------------------------------------------------------- 575 ! -- 2) Add Increments 576 ! ------------------------------------------------------ 577 578 IF (ln_decoupNC) THEN 579 cbu_i_bio(jn_din,layer) = cbu_i_bio(jn_din,layer) + zd_din 580 cbu_i_bio(jn_aon,layer) = cbu_i_bio(jn_aon,layer) + zd_aon 581 cbu_i_bio(jn_eon,layer) = cbu_i_bio(jn_eon,layer) + zd_eon 582 ELSE 583 cbu_i_bio(jn_din,layer) = cbu_i_bio(jn_din,layer) + zd_din 584 ENDIF 585 586 IF (ln_decoupPC) THEN 587 cbu_i_bio(jn_dip,layer) = cbu_i_bio(jn_dip,layer) + zd_dip 588 cbu_i_bio(jn_aop,layer) = cbu_i_bio(jn_aop,layer) + zd_aop 589 cbu_i_bio(jn_eop,layer) = cbu_i_bio(jn_eop,layer) + zd_eop 590 ELSE 591 cbu_i_bio(jn_dip,layer) = cbu_i_bio(jn_dip,layer) + zd_dip 592 ENDIF 593 353 594 354 595 cbu_i_bio(jn_dsi,layer) = cbu_i_bio(jn_dsi,layer) + zd_dsi 355 cbu_i_bio(jn_din,layer) = cbu_i_bio(jn_din,layer) + zd_no3356 cbu_i_bio(jn_dip,layer) = cbu_i_bio(jn_dip,layer) + zd_po4357 596 cbu_i_bio(jn_aoc,layer) = cbu_i_bio(jn_aoc,layer) + zd_daf 358 597 359 598 IF ( nn_bio_opt .GE. 1 ) 360 599 & cbu_i_bio(jn_eoc,layer) = cbu_i_bio(jn_eoc,layer) + zd_eoc 600 361 601 362 602 IF ( ln_write_bio ) THEN … … 375 615 WRITE(numout,*) ' lim_tem : ', lim_tem(layer) 376 616 WRITE(numout,*) ' lim_sal : ', lim_sal(layer) 617 WRITE(numout,*) ' lim_din_stock : ', lim_din_stock(layer) 618 WRITE(numout,*) ' lim_dip_stock : ', lim_dip_stock(layer) 619 WRITE(numout,*) ' lim_dsi_stock : ', lim_dsi_stock(layer) 377 620 WRITE(numout,*) 378 621 WRITE(numout,*) ' tc_bio : ', tc_bio(layer) … … 381 624 WRITE(numout,*) 382 625 WRITE(numout,*) ' zd_dsi : ', zd_dsi 383 WRITE(numout,*) ' zd_no3 : ', zd_no3384 626 WRITE(numout,*) ' zd_po4 : ', zd_po4 385 627 WRITE(numout,*) ' zd_daf : ', zd_daf 386 628 WRITE(numout,*) ' zd_eoc : ', zd_eoc 629 WRITE(numout,*) ' zd_din : ', zd_din 387 630 WRITE(numout,*) 388 631 WRITE(numout,*) ' dSi : ', cbu_i_bio(jn_dsi,layer) 389 WRITE(numout,*) ' dIN : ', cbu_i_bio(jn_din,layer) 632 WRITE(numout,*) ' dIN : ', cbu_i_bio(jn_din,layer) 390 633 WRITE(numout,*) ' dIP : ', cbu_i_bio(jn_dip,layer) 391 634 WRITE(numout,*) ' AoC : ', cbu_i_bio(jn_aoc,layer) 392 635 WRITE(numout,*) ' eoC : ', cbu_i_bio(jn_eoc,layer) 393 WRITE(numout,*) 636 WRITE(numout,*) ' krem_N : ', krem_N 637 638 639 IF (ln_decoupNC) THEN 640 WRITE(numout,*) ' syn_N : ', syn_N(layer) 641 WRITE(numout,*) ' rsp_N : ', rsp_N(layer) 642 WRITE(numout,*) ' lys_N : ', lys_N(layer) 643 WRITE(numout,*) ' rem_N : ', rem_N(layer) 644 WRITE(numout,*) 645 WRITE(numout,*) ' zd_aon : ', zd_aon 646 WRITE(numout,*) ' zd_eon : ', zd_eon 647 WRITE(numout,*) 648 WRITE(numout,*) ' AoN : ', cbu_i_bio(jn_aon,layer) 649 WRITE(numout,*) ' eoN : ', cbu_i_bio(jn_eon,layer) 650 WRITE(numout,*) 651 WRITE(numout,*) ' N_C_alg : ', N_C_alg(layer) 652 WRITE(numout,*) ' N_C_det : ', N_C_det(layer) 653 ENDIF ! ln_decoupNC 394 654 395 655 ENDIF ! ln_write_bio 396 656 397 657 END DO ! layer 658 659 398 660 399 661 !------------------------------------------------------------------------ … … 408 670 & .AND. flag_active(jn_alk) ) THEN 409 671 410 zd_dic1 = - zncp * ddtb 672 zd_dic1 = - zncp * ddtb !! Bizarre 411 673 zd_dic2 = - cbu_i_bio(jn_dic,layer) 412 674 zd_dic = MAX( zd_dic1 , zd_dic2 ) 413 675 414 zd_alk1 = no3_c * zncp * ddtb 676 zd_alk1 = no3_c * zncp * ddtb !! Bizarre 415 677 zd_alk2 = - cbu_i_bio(jn_alk,layer) 416 678 zd_alk = MAX( zd_alk1 , zd_alk2 )
Note: See TracChangeset
for help on using the changeset viewer.