Changeset 5841 for branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM
- Timestamp:
- 2015-10-30T12:48:06+01:00 (9 years ago)
- Location:
- branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC
- Files:
-
- 8 added
- 14 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90
r5726 r5841 172 172 !! 173 173 !! UKESM diagnostics 174 INTEGER :: jdms !: include DMS diagnostics ? Jpalm (27-08-2014) 174 INTEGER :: jdms !: include DMS diagnostics ? Jpalm (27-08-2014) 175 INTEGER :: jdms_input !: use instant (0) or diel-average (1) inputs (AXY, 08/07/2015) 176 INTEGER :: jdms_model !: choice of DMS model passed to atmosphere 177 !! 1 = ANDR, 2 = SIMO, 3 = ARAN, 4 = HALL 175 178 !! 176 179 !! … … 217 220 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_sed_ca !: 2D inorganic carbon (now) 218 221 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_sed_ca !: 2D inorganic carbon (after) 222 !! 223 !! 2D fields of temporally averaged properties for DMS calculations (AXY, 07/07/15) 224 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zb_dms_chn !: 2D avg CHN (before) 225 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_dms_chn !: 2D avg CHN (now) 226 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_dms_chn !: 2D avg CHN (after) 227 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zb_dms_chd !: 2D avg CHD (before) 228 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_dms_chd !: 2D avg CHD (now) 229 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_dms_chd !: 2D avg CHD (after) 230 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zb_dms_mld !: 2D avg MLD (before) 231 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_dms_mld !: 2D avg MLD (now) 232 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_dms_mld !: 2D avg MLD (after) 233 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zb_dms_qsr !: 2D avg QSR (before) 234 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_dms_qsr !: 2D avg QSR (now) 235 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_dms_qsr !: 2D avg QSR (after) 236 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zb_dms_din !: 2D avg DIN (before) 237 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zn_dms_din !: 2D avg DIN (now) 238 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: za_dms_din !: 2D avg DIN (after) 219 239 #endif 220 240 … … 415 435 !!---------------------------------------------------------------------- 416 436 USE lib_mpp , ONLY: ctl_warn 417 INTEGER :: ierr( 6) ! Local variables437 INTEGER :: ierr(7) ! Local variables 418 438 !!---------------------------------------------------------------------- 419 439 ierr(:) = 0 … … 439 459 & zb_sed_ca(jpi,jpj) , zn_sed_ca(jpi,jpj) , & 440 460 & za_sed_ca(jpi,jpj) , STAT=ierr(3) ) 461 !* 2D fields of temporally averaged properties for DMS calculations (AXY, 07/07/15) 462 ALLOCATE( zb_dms_chn(jpi,jpj) , zn_dms_chn(jpi,jpj) , & 463 & za_dms_chn(jpi,jpj) , & 464 & zb_dms_chd(jpi,jpj) , zn_dms_chd(jpi,jpj) , & 465 & za_dms_chd(jpi,jpj) , & 466 & zb_dms_mld(jpi,jpj) , zn_dms_mld(jpi,jpj) , & 467 & za_dms_mld(jpi,jpj) , & 468 & zb_dms_qsr(jpi,jpj) , zn_dms_qsr(jpi,jpj) , & 469 & za_dms_qsr(jpi,jpj) , & 470 & zb_dms_din(jpi,jpj) , zn_dms_din(jpi,jpj) , & 471 & za_dms_din(jpi,jpj) , STAT=ierr(4) ) 441 472 # endif 442 473 !* 2D fields of miscellaneous parameters … … 444 475 & dustmo(jpi,jpj,2) , riv_n(jpi,jpj) , & 445 476 & riv_si(jpi,jpj) , riv_c(jpi,jpj) , & 446 & riv_alk(jpi,jpj) , friver_dep(jpk,jpk) , STAT=ierr( 4) )477 & riv_alk(jpi,jpj) , friver_dep(jpk,jpk) , STAT=ierr(5) ) 447 478 !* 2D and 3D fields of light parameters 448 479 ALLOCATE( neln(jpi,jpj) , xze(jpi,jpj) , & 449 & xpar(jpi,jpj,jpk) , STAT=ierr( 5) )480 & xpar(jpi,jpj,jpk) , STAT=ierr(6) ) 450 481 !* 2D and 3D fields of sediment-associated parameters 451 482 ALLOCATE( dminl(jpi,jpj) , dmin3(jpi,jpj,jpk) , & … … 454 485 & fbodf(jpi,jpj) , fbods(jpi,jpj) , & 455 486 & ffln(jpi,jpj,jpk) , fflf(jpi,jpj,jpk) , & 456 & ffls(jpi,jpj,jpk) , cmask(jpi,jpj) , STAT=ierr( 6) )487 & ffls(jpi,jpj,jpk) , cmask(jpi,jpj) , STAT=ierr(7) ) 457 488 #endif 458 489 ! -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90
r5726 r5841 16 16 !! - ! 2013-05 (A. Yool) updated for v3.5 17 17 !! - ! 2014-08 (A. Yool, J. Palm) Add DMS module for UKESM1 model 18 !! - ! 2015-06 (A. Yool) Update to include MOCSY 19 !! - ! 2015-07 (A. Yool) Update for rolling averages 18 20 !!---------------------------------------------------------------------- 19 21 !! … … 36 38 #endif 37 39 !! 40 #if defined key_mocsy 41 !!---------------------------------------------------------------------- 42 !! Updates with the addition of MOCSY include: 43 !! - option to use PML or MOCSY carbonate chemistry (the latter is 44 !! preferred) 45 !! - central calculation of gas transfer velocity, f_kw660; previously 46 !! this was done separately for CO2 and O2 with predictable results 47 !! - distribution of f_kw660 to both PML and MOCSY CO2 air-sea flux 48 !! calculations and to those for O2 air-sea flux 49 !! - extra diagnostics included for MOCSY 50 !!---------------------------------------------------------------------- 51 #endif 52 !! 38 53 #if defined key_medusa 39 54 !!---------------------------------------------------------------------- … … 55 70 # endif 56 71 # if defined key_roam 72 USE gastransfer 73 # if defined key_mocsy 74 USE mocsy_wrapper 75 # else 57 76 USE trcco2_medusa 77 # endif 58 78 USE trcoxy_medusa 59 79 !! Jpalm (08/08/2014) … … 134 154 REAL(wp) :: ztmp, zsal 135 155 # endif 156 # if defined key_mocsy 157 REAL(wp) :: zpho 158 # endif 136 159 !! 137 160 !! integrated source and sink terms … … 266 289 REAL(wp) :: f_kw660, f_o2flux, f_o2sat 267 290 REAL(wp), DIMENSION(jpi,jpj) :: f_omcal, f_omarg 291 !! 292 !! AXY (23/06/15): additional diagnostics for MOCSY and oxygen 293 REAL(wp) :: f_fco2w, f_BetaD, f_rhosw, f_opres, f_insitut, f_pco2atm, f_fco2atm 294 REAL(wp) :: f_schmidtco2, f_kwco2, f_K0, f_co2starair, f_dpco2, f_kwo2 295 !! 268 296 INTEGER :: iters 269 297 REAL(wp) :: f_year … … 305 333 306 334 !!--------------------------------------------------------------------- 335 336 # if defined key_debug_medusa 337 IF (lwp) write (numout,*) 'trc_bio_medusa: variables defined' 338 CALL flush(numout) 339 # endif 307 340 308 341 !! AXY (20/11/14): alter this to report on first MEDUSA call … … 486 519 # endif 487 520 521 # if defined key_debug_medusa 522 IF (lwp) write (numout,*) 'trc_bio_medusa: variables initialised and checked' 523 CALL flush(numout) 524 # endif 525 488 526 # if defined key_roam 489 527 !!---------------------------------------------------------------------- … … 514 552 f_pco2a = fq4 515 553 endif 516 # if defined key_axy_pi_co2554 # if defined key_axy_pi_co2 517 555 f_pco2a = hist_pco2(1) 518 IF(lwp) WRITE(numout,*) ' MEDUSA atm pCO2 = FIXED' 519 # endif 556 # endif 520 557 !! IF(lwp) WRITE(numout,*) ' MEDUSA nyear =', nyear 521 558 !! IF(lwp) WRITE(numout,*) ' MEDUSA nsec_day =', real(nsec_day) … … 529 566 # endif 530 567 568 # if defined key_debug_medusa 569 IF (lwp) write (numout,*) 'trc_bio_medusa: ready for carbonate chemistry' 570 IF (lwp) write (numout,*) 'trc_bio_medusa: kt = ', kt 571 IF (lwp) write (numout,*) 'trc_bio_medusa: nittrc000 = ', nittrc000 572 CALL flush(numout) 573 # endif 574 531 575 # if defined key_roam 532 !! AXY (20/11/14): alter to call on first MEDUSA timestep 576 !! AXY (20/11/14): alter to call on first MEDUSA timestep and then every 577 !! month (this is hardwired as 960 timesteps but should 578 !! be calculated and done properly 533 579 !! IF( kt == nit000 .or. mod(kt,1920) == 0 ) THEN 534 IF( kt == nittrc000 .or. mod(kt, 1920) == 0 ) THEN580 IF( kt == nittrc000 .or. mod(kt,960) == 0 ) THEN 535 581 !!---------------------------------------------------------------------- 536 582 !! Calculate the carbonate chemistry for the whole ocean on the first … … 540 586 !! 541 587 IF(lwp) WRITE(numout,*) ' MEDUSA calculating all carbonate chemistry at kt =', kt 588 CALL flush(numout) 542 589 !! blank flags 543 590 i2_omcal(:,:) = 0 … … 549 596 !! OPEN wet point IF..THEN loop 550 597 if (tmask(ji,jj,jk).eq.1) then 551 !! carbonate chemistry 598 !! do carbonate chemistry 599 !! 552 600 fdep2 = fsdept(ji,jj,jk) !! set up level midpoint 601 !! 553 602 !! set up required state variables 554 603 zdic = max(0.,trn(ji,jj,jk,jpdic)) !! dissolved inorganic carbon … … 556 605 ztmp = tsn(ji,jj,jk,jp_tem) !! temperature 557 606 zsal = tsn(ji,jj,jk,jp_sal) !! salinity 607 # if defined key_mocsy 608 zsil = max(0.,trn(ji,jj,jk,jpsil)) !! silicic acid 609 zpho = max(0.,trn(ji,jj,jk,jpdin)) / 16.0 !! phosphate via DIN and Redfield 610 # endif 558 611 !! 559 612 !! AXY (28/02/14): check input fields … … 571 624 ji, ',', jj, ',', jk, ') at time', kt 572 625 endif 626 !! 627 !! blank input variables not used at this stage (they relate to air-sea flux) 628 f_kw660 = 1.0 629 f_pp0 = 1.0 630 !! 573 631 !! calculate carbonate chemistry at grid cell midpoint 574 CALL trc_co2_medusa( ztmp, zsal, zdic, zalk, fdep2, 5.0, f_pco2a, & ! inputs 575 f_ph, f_pco2w, f_h2co3, f_hco3, f_co3, f_omcal(ji,jj), & ! outputs 632 # if defined key_mocsy 633 !! AXY (22/06/15): use Orr & Epitalon (2015) MOCSY-2 carbonate 634 !! chemistry package 635 CALL mocsy_interface( ztmp, zsal, zalk, zdic, zsil, zpho, & ! inputs 636 f_pp0, fdep2, flatx, f_kw660, f_pco2a, 1, & ! inputs 637 f_ph, f_pco2w, f_fco2w, f_h2co3, f_hco3, f_co3, f_omarg(ji,jj), & ! outputs 638 f_omcal(ji,jj), f_BetaD, f_rhosw, f_opres, f_insitut, & ! outputs 639 f_pco2atm, f_fco2atm, f_schmidtco2, f_kwco2, f_K0, & ! outputs 640 f_co2starair, f_co2flux, f_dpco2 ) ! outputs 641 !! 642 f_TDIC = (zdic / f_rhosw) * 1000. ! mmol / m3 -> umol / kg 643 f_TALK = (zalk / f_rhosw) * 1000. ! meq / m3 -> ueq / kg 644 f_dcf = f_rhosw 645 # else 646 !! AXY (22/06/15): use old PML carbonate chemistry package (the 647 !! MEDUSA-2 default) 648 CALL trc_co2_medusa( ztmp, zsal, zdic, zalk, fdep2, f_kw660, & ! inputs 649 f_pco2a, f_ph, f_pco2w, f_h2co3, f_hco3, f_co3, f_omcal(ji,jj), & ! outputs 576 650 f_omarg(ji,jj), f_co2flux, f_TDIC, f_TALK, f_dcf, f_henry, iters) ! outputs 577 651 !! … … 581 655 iters, ' AT (', ji, ', ', jj, ', ', jk, ') AT ', kt 582 656 endif 657 # endif 658 !! 583 659 !! store 3D outputs 584 660 f3_pH(ji,jj,jk) = f_ph … … 588 664 f3_omcal(ji,jj,jk) = f_omcal(ji,jj) 589 665 f3_omarg(ji,jj,jk) = f_omarg(ji,jj) 666 !! 590 667 !! CCD calculation: calcite 591 668 if (i2_omcal(ji,jj) .eq. 0 .and. f_omcal(ji,jj) .lt. 1.0) then … … 607 684 i2_omcal(ji,jj) = 1 608 685 endif 686 !! 609 687 !! CCD calculation: aragonite 610 688 if (i2_omarg(ji,jj) .eq. 0 .and. f_omarg(ji,jj) .lt. 1.0) then … … 633 711 # endif 634 712 713 # if defined key_debug_medusa 714 IF (lwp) write (numout,*) 'trc_bio_medusa: ready for full domain calculations' 715 CALL flush(numout) 716 # endif 717 635 718 !!---------------------------------------------------------------------- 636 719 !! MEDUSA has unified equation through the water column … … 656 739 !! OPEN wet point IF..THEN loop 657 740 if (tmask(ji,jj,jk).eq.1) then 658 659 741 !!====================================================================== 660 742 !! SETUP LOCAL GRID CELL … … 814 896 !!---------------------------------------------------------------------- 815 897 !! 816 !! a bit of set up ...817 ! f_uwind = zwnd_i(ji,jj)818 ! f_vwind = zwnd_j(ji,jj)819 898 !! AXY (17/07/14): zwind_i and zwind_j do not exist in this 820 899 !! version of NEMO because it does not include … … 827 906 !! revisited when MEDUSA properly interacts 828 907 !! with UKESM1 physics 829 ! f_uwind = zwind_i(ji,jj) 830 ! f_vwind = zwind_j(ji,jj) 831 ! f_wind = ((f_uwind**2.0) + (f_vwind**2.0))**0.5 908 !! 832 909 f_wind = wndm(ji,jj) 833 !! AXY (17/07/14): the current oxygen code takes in separate 834 !! U and V components of the wind; to avoid 835 !! the need to change this, calculate these 836 !! components based on wndm; again, this is 837 !! not ideal, but should suffice as a 838 !! temporary measure; in the long-term all 839 !! of MEDUSA's air-sea gas exchange terms 840 !! will be revisited to ensure that they are 841 !! valid and self-consistent; and the CO2 842 !! code will be wholly replaced with a more 843 !! up-to-date parameterisation 844 f_uwind = ((f_wind**2.0) / 2.0)**0.5 845 f_vwind = f_uwind 910 !! 911 !! AXY (23/06/15): as part of an effort to update the carbonate chemistry 912 !! in MEDUSA, the gas transfer velocity used in the carbon 913 !! and oxygen cycles has been harmonised and is calculated 914 !! by the same function here; this harmonisation includes 915 !! changes to the PML carbonate chemistry scheme so that 916 !! it too makes use of the same gas transfer velocity; the 917 !! preferred parameterisation of this is Wanninkhof (2014), 918 !! option 7 919 !! 920 # if defined key_debug_medusa 921 IF (lwp) write (numout,*) 'trc_bio_medusa: entering gas_transfer' 922 CALL flush(numout) 923 # endif 924 CALL gas_transfer( f_wind, 1, 7, & ! inputs 925 f_kw660 ) ! outputs 926 # if defined key_debug_medusa 927 IF (lwp) write (numout,*) 'trc_bio_medusa: exiting gas_transfer' 928 CALL flush(numout) 929 # endif 930 !! 931 !! air pressure (atm); ultimately this will use air pressure at the base 932 !! of the UKESM1 atmosphere 933 !! 846 934 f_pp0 = 1.0 935 !! 847 936 !! IF(lwp) WRITE(numout,*) ' MEDUSA ztmp =', ztmp 848 937 !! IF(lwp) WRITE(numout,*) ' MEDUSA zwind_i =', zwind_i(ji,jj) … … 852 941 !! 853 942 # if defined key_axy_carbchem 943 # if defined key_mocsy 944 !! 945 !! AXY (22/06/15): use Orr & Epitalon (2015) MOCSY-2 carbonate 946 !! chemistry package; note that depth is set to 947 !! zero in this call 948 CALL mocsy_interface( ztmp, zsal, zalk, zdic, zsil, zpho, & ! inputs 949 f_pp0, 0.0, flatx, f_kw660, f_pco2a, 1, & ! inputs 950 f_ph, f_pco2w, f_fco2w, f_h2co3, f_hco3, f_co3, f_omarg(ji,jj), & ! outputs 951 f_omcal(ji,jj), f_BetaD, f_rhosw, f_opres, f_insitut, & ! outputs 952 f_pco2atm, f_fco2atm, f_schmidtco2, f_kwco2, f_K0, & ! outputs 953 f_co2starair, f_co2flux, f_dpco2 ) ! outputs 954 !! 955 f_TDIC = (zdic / f_rhosw) * 1000. ! mmol / m3 -> umol / kg 956 f_TALK = (zalk / f_rhosw) * 1000. ! meq / m3 -> ueq / kg 957 f_dcf = f_rhosw 958 # else 854 959 iters = 0 855 960 !! 856 961 !! carbon dioxide (CO2); Jerry Blackford code (ostensibly OCMIP-2, but not) 857 CALL trc_co2_medusa( ztmp, zsal, zdic, zalk, 0.0, f_ wind, f_pco2a, &! inputs858 f_ph, f_pco2w, f_h2co3, f_hco3, f_co3, f_omcal(ji,jj), &! outputs859 f_omarg(ji,jj), f_co2flux, f_TDIC, f_TALK, f_dcf, f_henry, iters ) 962 CALL trc_co2_medusa( ztmp, zsal, zdic, zalk, 0.0, f_kw660, f_pco2a, & ! inputs 963 f_ph, f_pco2w, f_h2co3, f_hco3, f_co3, f_omcal(ji,jj), & ! outputs 964 f_omarg(ji,jj), f_co2flux, f_TDIC, f_TALK, f_dcf, f_henry, iters ) ! outputs 860 965 !! 861 966 !! AXY (09/01/14): removed iteration and NaN checks; these have … … 869 974 iters, ' AT (', ji, ', ', jj, ', ', jk, ') AT ', kt 870 975 endif 976 # endif 871 977 # else 872 978 !! AXY (18/04/13): switch off carbonate chemistry calculations; provide … … 882 988 f_TDIC = zdic 883 989 f_TALK = zalk 884 f_dcf = 1. 990 f_dcf = 1.026 885 991 f_henry = 1. 992 !! AXY (23/06/15): add in some extra MOCSY diagnostics 993 f_fco2w = fpco2a 994 f_BetaD = 1. 995 f_rhosw = 1.026 996 f_opres = 0. 997 f_insitut = ztmp 998 f_pco2atm = fpco2a 999 f_fco2atm = fpco2a 1000 f_schmidtco2 = 660. 1001 f_kwco2 = 0. 1002 f_K0 = 0. 1003 f_co2starair = fpco2a 1004 f_dpco2 = 0. 886 1005 # endif 887 1006 !! 888 !! already in right units; correct for sea-ice; divide through by layer thickness889 f_co2flux = (1. - fr_i(ji,jj)) * f_co2flux / fthk1007 !! mmol/m2/s -> mmol/m3/d; correct for sea-ice; divide through by layer thickness 1008 f_co2flux = (1. - fr_i(ji,jj)) * f_co2flux * 86400. / fthk 890 1009 !! 891 1010 !! oxygen (O2); OCMIP-2 code 892 CALL trc_oxy_medusa( ztmp, zsal, f_uwind, f_vwind, f_pp0, zoxy / 1000., fthk, & ! inputs 893 f_kw660, f_o2flux, f_o2sat ) ! outputs 894 !! 895 !! mol/m3/s -> mmol/m3/d; correct for sea-ice 896 f_o2flux = (1. - fr_i(ji,jj)) * f_o2flux * 1000. * 60. * 60. * 24. 897 f_o2sat = f_o2sat * 1000. 1011 !! AXY (23/06/15): amend input list for oxygen to account for common gas 1012 !! transfer velocity 1013 !! CALL trc_oxy_medusa( ztmp, zsal, f_uwind, f_vwind, f_pp0, zoxy / 1000., fthk, & ! inputs 1014 !! f_kw660, f_o2flux, f_o2sat ) ! outputs 1015 CALL trc_oxy_medusa( ztmp, zsal, f_kw660, f_pp0, zoxy, & ! inputs 1016 f_kwo2, f_o2flux, f_o2sat ) ! outputs 1017 !! 1018 !! mmol/m2/s -> mol/m3/d; correct for sea-ice; divide through by layer thickness 1019 f_o2flux = (1. - fr_i(ji,jj)) * f_o2flux * 86400. / fthk 898 1020 !! 899 1021 !! Jpalm (08-2014) … … 911 1033 !! 912 1034 IF (jdms .eq. 1) THEN 913 !! 914 !! CALL trc_dms_medusa( zchn, zchd, hmld(ji,jj), & !! inputs 915 !! dms_surf ) !! outputs 916 CALL trc_dms_medusa( zchn, zchd, hmld(ji,jj), qsr(ji,jj), zdin, & !! inputs 917 dms_surf, dms_andr, dms_simo, dms_aran, dms_hall ) !! outputs 1035 !! 1036 !! feed in correct inputs 1037 if (jdms_input .eq. 0) then 1038 !! use instantaneous inputs 1039 CALL trc_dms_medusa( zchn, zchd, hmld(ji,jj), qsr(ji,jj), zdin, & ! inputs 1040 dms_andr, dms_simo, dms_aran, dms_hall ) ! outputs 1041 else 1042 !! use diel-average inputs 1043 CALL trc_dms_medusa( zn_dms_chn(ji,jj), zn_dms_chd(ji,jj), & ! inputs 1044 zn_dms_mld(ji,jj), zn_dms_qsr(ji,jj), zn_dms_din(ji,jj), & ! inputs 1045 dms_andr, dms_simo, dms_aran, dms_hall ) ! outputs 1046 endif 1047 !! 1048 !! assign correct output to variable passed to atmosphere 1049 if (jdms_model .eq. 1) then 1050 dms_surf = dms_andr 1051 elseif (jdms_model .eq. 2) then 1052 dms_surf = dms_simo 1053 elseif (jdms_model .eq. 3) then 1054 dms_surf = dms_aran 1055 elseif (jdms_model .eq. 4) then 1056 dms_surf = dms_hall 1057 endif 918 1058 ENDIF 919 1059 !! End DMS Loop … … 953 1093 !! alkalinity are derived from continent-scale DIC estimates (Huang et al., 954 1094 !! 2012) and some Arctic river alkalinity estimates (Katya?) 955 !! 1095 !! 956 1096 !! as of 19/07/12, riverine nutrients can now be spread vertically across 957 1097 !! several grid cells rather than just poured into the surface box; this … … 3317 3457 ENDDO 3318 3458 ENDDO 3459 !! AXY (07/07/15): temporary hijacking 3460 # if defined key_roam 3461 trc2d(:,:,126) = zn_dms_chn(:,:) 3462 trc2d(:,:,127) = zn_dms_chd(:,:) 3463 trc2d(:,:,128) = zn_dms_mld(:,:) 3464 trc2d(:,:,129) = zn_dms_qsr(:,:) 3465 trc2d(:,:,130) = zn_dms_din(:,:) 3466 # endif 3319 3467 3320 3468 if (ibenthic.eq.2) then … … 3322 3470 !! effectively commented out because it does not work as 3323 3471 !! anticipated; it can be deleted at a later date 3324 if (jorgben.eq.1) then 3325 za_sed_n(:,:) = ( f_sbenin_n(:,:) + f_fbenin_n(:,:) - f_benout_n(:,:) ) * rdt 3326 za_sed_fe(:,:) = ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - f_benout_fe(:,:) ) * rdt 3327 za_sed_c(:,:) = ( f_sbenin_c(:,:) + f_fbenin_c(:,:) - f_benout_c(:,:) ) * rdt 3328 endif 3329 if (jinorgben.eq.1) then 3330 za_sed_si(:,:) = ( f_fbenin_si(:,:) - f_benout_si(:,:) ) * rdt 3331 za_sed_ca(:,:) = ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) * rdt 3332 endif 3333 !! 3334 !! Leap-frog scheme - only in explicit case, otherwise the time stepping 3335 !! is already being done in trczdf 3336 !! IF( l_trczdf_exp .AND. (ln_trcadv_cen2 .OR. ln_trcadv_tvd) ) THEN 3337 !! zfact = 2. * rdttra(jk) * FLOAT( ndttrc ) 3338 !! IF( neuler == 0 .AND. kt == nittrc000 ) zfact = rdttra(jk) * FLOAT(ndttrc) 3339 !! if (jorgben.eq.1) then 3340 !! za_sed_n(:,:) = zb_sed_n(:,:) + ( zfact * za_sed_n(:,:) ) 3341 !! za_sed_fe(:,:) = zb_sed_fe(:,:) + ( zfact * za_sed_fe(:,:) ) 3342 !! za_sed_c(:,:) = zb_sed_c(:,:) + ( zfact * za_sed_c(:,:) ) 3343 !! endif 3344 !! if (jinorgben.eq.1) then 3345 !! za_sed_si(:,:) = zb_sed_si(:,:) + ( zfact * za_sed_si(:,:) ) 3346 !! za_sed_ca(:,:) = zb_sed_ca(:,:) + ( zfact * za_sed_ca(:,:) ) 3347 !! endif 3348 !! ENDIF 3349 !! 3350 !! Time filter and swap of arrays 3351 IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN ! centred or tvd scheme 3352 IF( neuler == 0 .AND. kt == nittrc000 ) THEN 3472 if (jorgben.eq.1) then 3473 za_sed_n(:,:) = ( f_sbenin_n(:,:) + f_fbenin_n(:,:) - f_benout_n(:,:) ) * rdt 3474 za_sed_fe(:,:) = ( f_sbenin_fe(:,:) + f_fbenin_fe(:,:) - f_benout_fe(:,:) ) * rdt 3475 za_sed_c(:,:) = ( f_sbenin_c(:,:) + f_fbenin_c(:,:) - f_benout_c(:,:) ) * rdt 3476 endif 3477 if (jinorgben.eq.1) then 3478 za_sed_si(:,:) = ( f_fbenin_si(:,:) - f_benout_si(:,:) ) * rdt 3479 za_sed_ca(:,:) = ( f_fbenin_ca(:,:) - f_benout_ca(:,:) ) * rdt 3480 endif 3481 !! 3482 !! Leap-frog scheme - only in explicit case, otherwise the time stepping 3483 !! is already being done in trczdf 3484 !! IF( l_trczdf_exp .AND. (ln_trcadv_cen2 .OR. ln_trcadv_tvd) ) THEN 3485 !! zfact = 2. * rdttra(jk) * FLOAT( ndttrc ) 3486 !! IF( neuler == 0 .AND. kt == nittrc000 ) zfact = rdttra(jk) * FLOAT(ndttrc) 3487 !! if (jorgben.eq.1) then 3488 !! za_sed_n(:,:) = zb_sed_n(:,:) + ( zfact * za_sed_n(:,:) ) 3489 !! za_sed_fe(:,:) = zb_sed_fe(:,:) + ( zfact * za_sed_fe(:,:) ) 3490 !! za_sed_c(:,:) = zb_sed_c(:,:) + ( zfact * za_sed_c(:,:) ) 3491 !! endif 3492 !! if (jinorgben.eq.1) then 3493 !! za_sed_si(:,:) = zb_sed_si(:,:) + ( zfact * za_sed_si(:,:) ) 3494 !! za_sed_ca(:,:) = zb_sed_ca(:,:) + ( zfact * za_sed_ca(:,:) ) 3495 !! endif 3496 !! ENDIF 3497 !! 3498 !! Time filter and swap of arrays 3499 IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN ! centred or tvd scheme 3500 IF( neuler == 0 .AND. kt == nittrc000 ) THEN 3501 if (jorgben.eq.1) then 3502 zb_sed_n(:,:) = zn_sed_n(:,:) 3503 zn_sed_n(:,:) = za_sed_n(:,:) 3504 za_sed_n(:,:) = 0.0 3505 !! 3506 zb_sed_fe(:,:) = zn_sed_fe(:,:) 3507 zn_sed_fe(:,:) = za_sed_fe(:,:) 3508 za_sed_fe(:,:) = 0.0 3509 !! 3510 zb_sed_c(:,:) = zn_sed_c(:,:) 3511 zn_sed_c(:,:) = za_sed_c(:,:) 3512 za_sed_c(:,:) = 0.0 3513 endif 3514 if (jinorgben.eq.1) then 3515 zb_sed_si(:,:) = zn_sed_si(:,:) 3516 zn_sed_si(:,:) = za_sed_si(:,:) 3517 za_sed_si(:,:) = 0.0 3518 !! 3519 zb_sed_ca(:,:) = zn_sed_ca(:,:) 3520 zn_sed_ca(:,:) = za_sed_ca(:,:) 3521 za_sed_ca(:,:) = 0.0 3522 endif 3523 ELSE 3524 if (jorgben.eq.1) then 3525 zb_sed_n(:,:) = (atfp * ( zb_sed_n(:,:) + za_sed_n(:,:) )) + (atfp1 * zn_sed_n(:,:) ) 3526 zn_sed_n(:,:) = za_sed_n(:,:) 3527 za_sed_n(:,:) = 0.0 3528 !! 3529 zb_sed_fe(:,:) = (atfp * ( zb_sed_fe(:,:) + za_sed_fe(:,:) )) + (atfp1 * zn_sed_fe(:,:)) 3530 zn_sed_fe(:,:) = za_sed_fe(:,:) 3531 za_sed_fe(:,:) = 0.0 3532 !! 3533 zb_sed_c(:,:) = (atfp * ( zb_sed_c(:,:) + za_sed_c(:,:) )) + (atfp1 * zn_sed_c(:,:) ) 3534 zn_sed_c(:,:) = za_sed_c(:,:) 3535 za_sed_c(:,:) = 0.0 3536 endif 3537 if (jinorgben.eq.1) then 3538 zb_sed_si(:,:) = (atfp * ( zb_sed_si(:,:) + za_sed_si(:,:) )) + (atfp1 * zn_sed_si(:,:)) 3539 zn_sed_si(:,:) = za_sed_si(:,:) 3540 za_sed_si(:,:) = 0.0 3541 !! 3542 zb_sed_ca(:,:) = (atfp * ( zb_sed_ca(:,:) + za_sed_ca(:,:) )) + (atfp1 * zn_sed_ca(:,:)) 3543 zn_sed_ca(:,:) = za_sed_ca(:,:) 3544 za_sed_ca(:,:) = 0.0 3545 endif 3546 ENDIF 3547 ELSE ! case of smolar scheme or muscl 3353 3548 if (jorgben.eq.1) then 3354 zb_sed_n(:,:) = z n_sed_n(:,:)3549 zb_sed_n(:,:) = za_sed_n(:,:) 3355 3550 zn_sed_n(:,:) = za_sed_n(:,:) 3356 3551 za_sed_n(:,:) = 0.0 3357 3552 !! 3358 zb_sed_fe(:,:) = z n_sed_fe(:,:)3553 zb_sed_fe(:,:) = za_sed_fe(:,:) 3359 3554 zn_sed_fe(:,:) = za_sed_fe(:,:) 3360 3555 za_sed_fe(:,:) = 0.0 3361 3556 !! 3362 zb_sed_c(:,:) = z n_sed_c(:,:)3557 zb_sed_c(:,:) = za_sed_c(:,:) 3363 3558 zn_sed_c(:,:) = za_sed_c(:,:) 3364 3559 za_sed_c(:,:) = 0.0 3365 3560 endif 3366 3561 if (jinorgben.eq.1) then 3367 zb_sed_si(:,:) = z n_sed_si(:,:)3562 zb_sed_si(:,:) = za_sed_si(:,:) 3368 3563 zn_sed_si(:,:) = za_sed_si(:,:) 3369 3564 za_sed_si(:,:) = 0.0 3370 3565 !! 3371 zb_sed_ca(:,:) = zn_sed_ca(:,:) 3372 zn_sed_ca(:,:) = za_sed_ca(:,:) 3373 za_sed_ca(:,:) = 0.0 3374 endif 3375 ELSE 3376 if (jorgben.eq.1) then 3377 zb_sed_n(:,:) = (atfp * ( zb_sed_n(:,:) + za_sed_n(:,:) )) + (atfp1 * zn_sed_n(:,:) ) 3378 zn_sed_n(:,:) = za_sed_n(:,:) 3379 za_sed_n(:,:) = 0.0 3380 !! 3381 zb_sed_fe(:,:) = (atfp * ( zb_sed_fe(:,:) + za_sed_fe(:,:) )) + (atfp1 * zn_sed_fe(:,:)) 3382 zn_sed_fe(:,:) = za_sed_fe(:,:) 3383 za_sed_fe(:,:) = 0.0 3384 !! 3385 zb_sed_c(:,:) = (atfp * ( zb_sed_c(:,:) + za_sed_c(:,:) )) + (atfp1 * zn_sed_c(:,:) ) 3386 zn_sed_c(:,:) = za_sed_c(:,:) 3387 za_sed_c(:,:) = 0.0 3388 endif 3389 if (jinorgben.eq.1) then 3390 zb_sed_si(:,:) = (atfp * ( zb_sed_si(:,:) + za_sed_si(:,:) )) + (atfp1 * zn_sed_si(:,:)) 3391 zn_sed_si(:,:) = za_sed_si(:,:) 3392 za_sed_si(:,:) = 0.0 3393 !! 3394 zb_sed_ca(:,:) = (atfp * ( zb_sed_ca(:,:) + za_sed_ca(:,:) )) + (atfp1 * zn_sed_ca(:,:)) 3566 zb_sed_ca(:,:) = za_sed_ca(:,:) 3395 3567 zn_sed_ca(:,:) = za_sed_ca(:,:) 3396 3568 za_sed_ca(:,:) = 0.0 3397 3569 endif 3398 3570 ENDIF 3399 ELSE ! case of smolar scheme or muscl3400 if (jorgben.eq.1) then3401 zb_sed_n(:,:) = za_sed_n(:,:)3402 zn_sed_n(:,:) = za_sed_n(:,:)3403 za_sed_n(:,:) = 0.03404 !!3405 zb_sed_fe(:,:) = za_sed_fe(:,:)3406 zn_sed_fe(:,:) = za_sed_fe(:,:)3407 za_sed_fe(:,:) = 0.03408 !!3409 zb_sed_c(:,:) = za_sed_c(:,:)3410 zn_sed_c(:,:) = za_sed_c(:,:)3411 za_sed_c(:,:) = 0.03412 endif3413 if (jinorgben.eq.1) then3414 zb_sed_si(:,:) = za_sed_si(:,:)3415 zn_sed_si(:,:) = za_sed_si(:,:)3416 za_sed_si(:,:) = 0.03417 !!3418 zb_sed_ca(:,:) = za_sed_ca(:,:)3419 zn_sed_ca(:,:) = za_sed_ca(:,:)3420 za_sed_ca(:,:) = 0.03421 endif3422 ENDIF3423 3571 endif 3424 3572 3425 3573 IF( ln_diatrc ) THEN 3426 3574 !!---------------------------------------------------------------------- -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcco2_medusa.F90
r5726 r5841 47 47 !======================================================================= 48 48 ! 49 SUBROUTINE trc_co2_medusa( Temp, Sal, DIC, ALK, Depth, Wnd, pCO2a, &49 SUBROUTINE trc_co2_medusa( Temp, Sal, DIC, ALK, Depth, xkw, pCO2a, & 50 50 pH, pCO2w, h2co3, hco3, co3, om_cal, om_arg, co2flux, TDIC, TALK, & 51 51 dcf, henry, iters ) … … 72 72 ! 17/02/2010. Update calculation of K1, K2, Kb to make consistant with the OCMIP protocols. 73 73 ! 29/07/2011. Merged into MEDUSA with a raft of changes to this subroutine; less elsewhere 74 ! 23/06/2015. Modified to take gas transfer velocity as an input (rather than wind speed); 75 ! alter CO2 flux to /s rather than /d for consistency with other schemes 74 76 ! 75 77 ! Changes for MEDUSA include: … … 85 87 REAL(wp), INTENT( in ) :: ALK ! meq / m3 86 88 REAL(wp), INTENT( in ) :: Depth ! m 87 REAL(wp), INTENT( in ) :: Wnd ! m / s 89 ! REAL(wp), INTENT( in ) :: Wnd ! m / s 90 REAL(wp), INTENT( in ) :: xkw ! m / s 88 91 REAL(wp), INTENT( in ) :: pCO2a ! uatm 89 92 !---------------------------------------------------------------------- … … 95 98 REAL(wp), INTENT( inout ) :: om_cal ! normalised 96 99 REAL(wp), INTENT( inout ) :: om_arg ! normalised 97 REAL(wp), INTENT( inout ) :: co2flux ! mmol / m2 / d100 REAL(wp), INTENT( inout ) :: co2flux ! mmol / m2 / s 98 101 REAL(wp), INTENT( inout ) :: TDIC ! umol / kg 99 102 REAL(wp), INTENT( inout ) :: TALK ! ueq / kg … … 129 132 ! (i.e. surface calculations being performed) 130 133 if (Depth .eq. 0.0) then 131 call Air_sea_exchange( Temp, Wnd, pCO2w, pCO2a, henry, dcf, & ! inputs134 call Air_sea_exchange( Temp, xkw, pCO2w, pCO2a, henry, dcf, & ! inputs 132 135 co2flux ) ! output 133 136 else … … 145 148 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: zdic =', DIC 146 149 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: zalk =', ALK 147 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: f_ wind =', Wnd150 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: f_kw660 =', xkw 148 151 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: f_ph =', ph 149 152 IF(lwp) WRITE(numout,*) ' trc_co2_medusa: f_pco2w =', pCO2w … … 162 165 & ' DIC', DIC, ' ALK', ALK 163 166 if (lwp) write (numout,'(a,a,f10.3,a,f10.3)') 'CO2FLUX-NAN', & 164 & ' WND', Wnd, ' PH ', ph167 & ' XKW', xkw, ' PH ', ph 165 168 if (lwp) write (numout,'(a,a,i6)') 'CO2FLUX-NAN', & 166 169 & ' ITERS', iters … … 196 199 ! WRITE(*,'(A27,F10.3)') " Omega calcite (~) = ", om_cal 197 200 ! WRITE(*,'(A27,F10.3)') " Omega aragonite (~) = ", om_arg 198 ! WRITE(*,'(A27,F10.3)') " air sea flux(mmol/m2/ d) = ", flux201 ! WRITE(*,'(A27,F10.3)') " air sea flux(mmol/m2/s) = ", flux 199 202 ! WRITE(*,*) " " 200 203 … … 287 290 !======================================================================= 288 291 ! 289 SUBROUTINE Air_sea_exchange( T, Wnd, pco2w, pco2a, henry, dcf, &292 SUBROUTINE Air_sea_exchange( T, xkw, pco2w, pco2a, henry, dcf, & 290 293 flux ) 291 294 ! … … 302 305 ! pCO2a partial pressure of CO2 in the atmosphere (usually external forcing). 303 306 ! T temperature (C) 304 ! Wnd wind speed, metres 307 ! Wnd wind speed, metres (DELETED) 308 ! xkw gas transfer velocity 305 309 ! Henry henry's constant 306 310 ! density the density of water for conversion between mmol/m3 and umol/kg … … 312 316 IMPLICIT NONE 313 317 314 REAL(wp), INTENT( in ) :: T, wnd, pco2w, pco2a, henry, dcf ! INPUT PARAMETERS:318 REAL(wp), INTENT( in ) :: T, xkw, pco2w, pco2a, henry, dcf ! INPUT PARAMETERS: 315 319 !----------------------------------------------------------------------- 316 320 REAL(wp), INTENT( inout ) :: flux ! OUTPUT Variables … … 320 324 ! calculate the Schmidt number and unit conversions 321 325 sc = 2073.1-125.62*T+3.6276*T**2.0-0.0432190*T**3.0 322 fwind = (0.222d0 * wnd**2d0 + 0.333d0 * wnd)*(sc/660.d0)**(-0.5) 326 ! fwind = (0.222d0 * wnd**2d0 + 0.333d0 * wnd)*(sc/660.d0)**(-0.5) 327 fwind = xkw * (sc/660.d0)**(-0.5) 323 328 fwind = fwind*24.d0/100.d0 ! convert to m/day 324 329 … … 326 331 ! here it is rescaled to mmol/m2/d 327 332 flux = fwind * henry * ( pco2a - pco2w ) * dcf 333 334 ! AXY (23/06/15): let's get it from /d to /s 335 flux = flux / ( 86400. ) 328 336 329 337 RETURN -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcdms_medusa.F90
r5726 r5841 39 39 ! 40 40 SUBROUTINE trc_dms_medusa( chn, chd, mld, xqsr, xdin, & !! inputs 41 & dms_ surf, dms_andr, dms_simo, dms_aran, dms_hall )!! outputs41 & dms_andr, dms_simo, dms_aran, dms_hall ) !! outputs 42 42 ! 43 43 !======================================================================= … … 70 70 !! published (and different from the above) 71 71 !! 72 !! AXY (08/07/15): amend to remove Julien's original calculation 73 !! as this is now superfluous; the four schemes 74 !! are calculated and one is chosen to be passed 75 !! to the atmosphere in trc_bio_medusa 76 !! 72 77 !======================================================================= 73 78 … … 79 84 REAL(wp), INTENT( in ) :: xqsr !! surface irradiance (W/m2) 80 85 REAL(wp), INTENT( in ) :: xdin !! surface DIN (mmol N/m3) 81 REAL(wp), INTENT( inout ) :: dms_surf !! DMS surface concentration (mol/m3)82 86 REAL(wp), INTENT( inout ) :: dms_andr !! DMS surface concentration (mol/m3) 83 87 REAL(wp), INTENT( inout ) :: dms_simo !! DMS surface concentration (mol/m3) … … 89 93 !! temporary variables 90 94 REAL(wp) :: fq1,fq2,fq3 91 !92 !! IJT (30/03/13): DMS calc needs this93 !! Julien : in Simo & Dachs, GBC, 2002, DMS is derived from94 !! CHL/MLD ratio in mg/m4 (i.e. CHL is in mg/m395 !! MLD in m).96 !! In MEDUSA, we already have CHL in mg/m3 for both97 !! Diatoms and non-diatoms (zchn,zchd); and mld from98 !! NEMO (hmld) in m.99 CHL = 0.0100 !!101 !! CHL = mask * TT(I,J,1,PHYTO_TRACER) &102 !! & * c2n_p * mw_carbon / CCHL_P(I,J,1,1)103 CHL = chn+chd !! mg/m3104 !!105 !! ------------------------------------------------106 !! Calculate the DMS concentration in nM (nanomol/litre)107 !! from Simo & Dachs, GBC, 2002, modified to be positive-definite108 !! for MLD>182.536m, using DMS=90./MLD (Aranami & Tsunogai, JGR, 2004)109 !! Multiply by 1.0E-6 to convert nM to (mol/m3)110 !! cmr = fm(i,1)*chl/mld(i)111 !! IF (cmr .lt. 0.02) THEN112 !! IF (mld(i) .le. 182.536) THEN113 !! csurf(i,astf_dms) = 1.0e-6*fm(i,1)*(-LOG(mld(i)) + 5.7)114 !! ELSE115 !! csurf(i,astf_dms) = 1.0e-6*fm(i,1)*(90./mld(i))116 !! ENDIF117 !! ELSE118 !! csurf(i,astf_dms) = 1.0e-6*fm(i,1)*(55.8*cmr + 0.6)119 !! ENDIF120 !!121 cmr = CHL / mld122 ! sw_dms = 0.5 + SIGN( 0.5, cmr - 0.02 )123 !! Jpalm (11-08-2014)124 !! Explanation about the SIGN function :125 !! not easy to read, but maybe "more elegant and efficient")126 !! here for example:127 !! sw_dms = 1 if cmr is greater than 0.02,128 !! 0 if cmr lower than 0.02129 !! then130 !! if cmr < 0.02131 !! dms_surf = 1.0e-6 * 90.0 / mld132 !! or = 1.0e-6 * 5.7 - LOG(mld)133 !! and if cmr > 0.02134 !! dms_surf = 1.0e-6 * ( 55.8 * cmr + 0.6 )135 !! what is equivalent to the IF loops formulations.136 !! difference is on the stresholds between mld = 182.536m137 !! (strange value...)138 !! and the Max function... that stay uncertain.139 !!140 ! dms_surf = 1.0e-6 * ( sw_dms * &141 ! & ( 55.8 * cmr + 0.6 ) + ( 1.0 - sw_dms ) * &142 ! & ( MAX( 90.0 / mld, 5.7 - LOG(mld) ) ) )143 !144 ! AXY (12/01/15): the DMS equation donated by the UKMO does not match145 ! that reported in Halloran et al. (2010); amend the146 ! equations appropriately147 !148 if (cmr .lt. 0.02) then149 dms_surf = (-1.0 * log(mld)) + 5.7150 else151 dms_surf = (55.8 * cmr) + 0.6152 endif153 !154 if (mld > 182.5) then155 dms_surf = (90.0 / mld)156 endif157 !158 dms_surf = 1.0e-6 * dms_surf159 160 95 ! 161 96 !======================================================================= … … 163 98 ! AXY (13/03/15): per remarks above, the following calculations estimate 164 99 ! DMS using all of the schemes examined for UKESM1 100 ! 101 CHL = 0.0 102 CHL = chn+chd !! mg/m3 103 cmr = CHL / mld 165 104 ! 166 105 ! AXY (13/03/15): Anderson et al. (2001) … … 180 119 ! 181 120 ! AXY (13/03/15): Simo & Dachs (2002) 182 cmr = CHL / mld183 121 fq1 = (-1 * log(mld)) + 5.7 184 122 fq2 = (55.8 * cmr) + 0.6 … … 191 129 ! 192 130 ! AXY (13/03/15): Aranami & Tsunogai (2004) 193 cmr = CHL / mld194 131 fq1 = 60.0 / mld 195 132 fq2 = (55.8 * cmr) + 0.6 … … 202 139 ! 203 140 ! AXY (13/03/15): Halloran et al. (2010) 204 cmr = CHL / mld205 141 fq1 = (-1 * log(mld)) + 5.7 206 142 fq2 = (55.8 * cmr) + 0.6 -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcini_medusa.F90
r5741 r5841 38 38 !! JPALM (14/09/15) 39 39 LOGICAL, PUBLIC :: & 40 ln_ccd = . FALSE.40 ln_ccd = .TRUE. 41 41 42 42 INTEGER :: & … … 236 236 237 237 !!---------------------------------------------------------------------- 238 !! Averaged properties for DMS calculations (various units) 239 !!---------------------------------------------------------------------- 240 !! 241 !! these store temporally averaged properties for DMS calculations (AXY, 07/07/15) 242 zb_dms_chn(:,:) = 0.0 !! CHN 243 zn_dms_chn(:,:) = 0.0 244 za_dms_chn(:,:) = 0.0 245 zb_dms_chd(:,:) = 0.0 !! CHD 246 zn_dms_chd(:,:) = 0.0 247 za_dms_chd(:,:) = 0.0 248 zb_dms_mld(:,:) = 0.0 !! MLD 249 zn_dms_mld(:,:) = 0.0 250 za_dms_mld(:,:) = 0.0 251 zb_dms_qsr(:,:) = 0.0 !! QSR 252 zn_dms_qsr(:,:) = 0.0 253 za_dms_qsr(:,:) = 0.0 254 zb_dms_din(:,:) = 0.0 !! DIN 255 zn_dms_din(:,:) = 0.0 256 za_dms_din(:,:) = 0.0 257 !! 258 IF(lwp) WRITE(numout,*) ' trc_ini_medusa: average fields for DMS initialised to zero' 259 IF(lwp) CALL flush(numout) 260 261 !!---------------------------------------------------------------------- 238 262 !! AXY (04/11/13): initialise fields previously done by trc_sed_medusa 239 263 !!---------------------------------------------------------------------- -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcnam_medusa.F90
r5736 r5841 83 83 & xsdiss, & 84 84 & vsed,xhr, & 85 & sedlam,sedlostpoc,jpkb,jdms 85 & sedlam,sedlostpoc,jpkb,jdms,jdms_input,jdms_model 86 86 #if defined key_roam 87 87 NAMELIST/natroam/ xthetaphy,xthetazoo,xthetanit, & … … 138 138 IF( ( .NOT.lk_iomput .AND. ln_diatrc ) .OR. ( ln_diatrc .AND. lk_medusa ) ) THEN 139 139 ! 140 ! Namelist nam pisdia140 ! Namelist nammeddia 141 141 ! ------------------- 142 REWIND( numnatp_ref ) ! Namelist nam pisdia in reference namelist : Piscesdiagnostics142 REWIND( numnatp_ref ) ! Namelist nammeddia in reference namelist : MEDUSA diagnostics 143 143 READ ( numnatp_ref, nammeddia, IOSTAT = ios, ERR = 901) 144 144 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammeddia in reference namelist', lwp ) 145 145 146 REWIND( numnatp_cfg ) ! Namelist nam pisdia in configuration namelist : Piscesdiagnostics146 REWIND( numnatp_cfg ) ! Namelist nammeddia in configuration namelist : MEDUSA diagnostics 147 147 READ ( numnatp_cfg, nammeddia, IOSTAT = ios, ERR = 902 ) 148 148 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nammeddia in configuration namelist', lwp ) … … 338 338 jpkb = 0. 339 339 jdms = 0 340 jdms_input = 0 341 jdms_input = 3 340 342 341 343 !REWIND(numnatm) … … 343 345 ! Namelist natbio 344 346 ! ------------------- 345 REWIND( numnatp_ref ) ! Namelist na mpisdia in reference namelist : Piscesdiagnostics347 REWIND( numnatp_ref ) ! Namelist natbio in reference namelist : MEDUSA diagnostics 346 348 READ ( numnatp_ref, natbio, IOSTAT = ios, ERR = 903) 347 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain reference namelist', lwp )348 349 REWIND( numnatp_cfg ) ! Namelist na mpisdia in configuration namelist : Piscesdiagnostics349 903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natbio in reference namelist', lwp ) 350 351 REWIND( numnatp_cfg ) ! Namelist natbio in configuration namelist : MEDUSA diagnostics 350 352 READ ( numnatp_cfg, natbio, IOSTAT = ios, ERR = 904 ) 351 904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain configuration namelist', lwp )353 904 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natbio in configuration namelist', lwp ) 352 354 IF(lwm) WRITE ( numonp, natbio ) 353 355 … … 488 490 !! UKESM1 - new diagnostics !! Jpalm 489 491 !! jdms : include dms diagnostics 490 !! 491 !! 492 !! 493 492 !! jdms_input : use instant (0) or diel-avg (1) inputs 493 !! jdms_model : choice of DMS model passed to atmosphere 494 !! 1 = ANDR, 2 = SIMO, 3 = ARAN, 4 = HALL 495 !! 494 496 IF(lwp) THEN 495 497 !! … … 510 512 & ' key_axy_carbchem = INACTIVE' 511 513 #endif 514 #if defined key_mocsy 515 WRITE(numout,*) & 516 & ' key_mocsy = ACTIVE' 517 #else 518 WRITE(numout,*) & 519 & ' key_mocsy = INACTIVE' 520 #endif 521 #if defined key_avgqsr_medusa 522 WRITE(numout,*) & 523 & ' key_avgqsr_medusa = ACTIVE' 524 #else 525 WRITE(numout,*) & 526 & ' key_avgqsr_medusa = INACTIVE' 527 #endif 512 528 #if defined key_bs_axy_zforce 513 529 WRITE(numout,*) & … … 544 560 WRITE(numout,*) & 545 561 & ' key_axy_pi_co2 = INACTIVE' 562 # endif 563 # if defined key_debug_medusa 564 WRITE(numout,*) & 565 & ' key_debug_medusa = ACTIVE' 566 #else 567 WRITE(numout,*) & 568 & ' key_debug_medusa = INACTIVE' 546 569 # endif 547 570 WRITE(numout,*) ' ' … … 971 994 & ' Vert layer for diagnostic of vertical flux, jpkp = ', jpkb 972 995 !! 973 !! UKESM1 - new diagnostics !! Jpalm 996 !! UKESM1 - new diagnostics !! Jpalm; AXY (08/07/15) 974 997 WRITE(numout,*) '=== UKESM1-related parameters' 975 998 WRITE(numout,*) & 976 999 & ' include DMS diagnostic?, jdms = ', jdms 1000 if (jdms_input .eq. 0) then 1001 WRITE(numout,*) & 1002 & ' use instant (0) or diel-avg (1) inputs, jdms_input = instantaneous' 1003 else 1004 WRITE(numout,*) & 1005 & ' use instant (0) or diel-avg (1) inputs, jdms_input = diel-average' 1006 endif 1007 if (jdms_model .eq. 1) then 1008 WRITE(numout,*) & 1009 & ' choice of DMS model passed to atmosphere, jdms_model = Anderson et al. (2001)' 1010 elseif (jdms_model .eq. 2) then 1011 WRITE(numout,*) & 1012 & ' choice of DMS model passed to atmosphere, jdms_model = Simo & Dachs (2002)' 1013 elseif (jdms_model .eq. 3) then 1014 WRITE(numout,*) & 1015 & ' choice of DMS model passed to atmosphere, jdms_model = Aranami & Tsunogai (2004)' 1016 elseif (jdms_model .eq. 4) then 1017 WRITE(numout,*) & 1018 & ' choice of DMS model passed to atmosphere, jdms_model = Halloran et al. (2010)' 1019 endif 977 1020 !! 978 1021 ENDIF … … 1032 1075 1033 1076 !READ(numnatm,natroam) 1034 ! Namelist nat bio1077 ! Namelist natroam 1035 1078 ! ------------------- 1036 REWIND( numnatp_ref ) ! Namelist na mpisdia in reference namelist : Piscesdiagnostics1037 READ ( numnatp_ref, nat bio, IOSTAT = ios, ERR = 905)1038 905 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain reference namelist', lwp )1039 1040 REWIND( numnatp_cfg ) ! Namelist na mpisdia in configuration namelist : Piscesdiagnostics1041 READ ( numnatp_cfg, nat bio, IOSTAT = ios, ERR = 906 )1042 906 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain configuration namelist', lwp )1043 IF(lwm) WRITE ( numonp, nat bio)1079 REWIND( numnatp_ref ) ! Namelist natroam in reference namelist : MEDUSA diagnostics 1080 READ ( numnatp_ref, natroam, IOSTAT = ios, ERR = 905) 1081 905 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natroam in reference namelist', lwp ) 1082 1083 REWIND( numnatp_cfg ) ! Namelist natroam in configuration namelist : MEDUSA diagnostics 1084 READ ( numnatp_cfg, natroam, IOSTAT = ios, ERR = 906 ) 1085 906 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natroam in configuration namelist', lwp ) 1086 IF(lwm) WRITE ( numonp, natroam ) 1044 1087 1045 1088 !! ROAM carbon, alkalinity and oxygen cycle parameters … … 1086 1129 ! Namelist natopt 1087 1130 ! ------------------- 1088 REWIND( numnatp_ref ) ! Namelist na mpisdia in reference namelist : Piscesdiagnostics1131 REWIND( numnatp_ref ) ! Namelist natopt in reference namelist : MEDUSA diagnostics 1089 1132 READ ( numnatp_ref, natopt, IOSTAT = ios, ERR = 907) 1090 907 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain reference namelist', lwp )1091 1092 REWIND( numnatp_cfg ) ! Namelist na mpisdia in configuration namelist : Piscesdiagnostics1133 907 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natopt in reference namelist', lwp ) 1134 1135 REWIND( numnatp_cfg ) ! Namelist natopt in configuration namelist : MEDUSA diagnostics 1093 1136 READ ( numnatp_cfg, natopt, IOSTAT = ios, ERR = 908 ) 1094 908 IF( ios /= 0 ) CALL ctl_nam ( ios , 'na mmeddiain configuration namelist', lwp )1137 908 IF( ios /= 0 ) CALL ctl_nam ( ios , 'natopt in configuration namelist', lwp ) 1095 1138 IF(lwm) WRITE ( numonp, natopt ) 1096 1139 -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcopt_medusa.F90
r5726 r5841 73 73 ! determination of surface irradiance 74 74 ! ----------------------------------- 75 zpar0m (:,:) = qsr (:,:) * 0.43 75 ! AXY (23/07/15): the inclusion of empirical DMS calculations requires 76 ! daily averages of a series of properties that are 77 ! used as inputs; these include surface irradiance; 78 ! here, this is taken advantage of to allow MEDUSA to 79 ! base its submarine light field on daily average 80 ! rather than "instantaneous" irradiance; largely 81 ! because MEDUSA was originally formulated to work 82 ! with diel average irradiance rather than a diel 83 ! cycle; using key_avgqsr_medusa activates this 84 ! functionality, while its absence gives default 85 ! MEDUSA (which is whatever is supplied by NEMO) 86 # if defined key_avgqsr_medusa 87 ! surface irradiance input is rolling average irradiance 88 zpar0m (:,:) = zn_dms_qsr(:,:) * 0.43 89 # else 90 ! surface irradiance input is instantaneous irradiance 91 zpar0m (:,:) = qsr(:,:) * 0.43 92 # endif 76 93 ! AXY (22/08/14): when zpar0m = 0, zpar100 is also zero and calculating 77 94 ! euphotic depth is not possible (cf. the Arctic Octopus); … … 92 109 zparr (:,:,1) = 0.5 * zpar0m(:,:) 93 110 zparg (:,:,1) = 0.5 * zpar0m(:,:) 94 95 111 96 112 ! determination of xpar … … 146 162 ENDDO 147 163 148 149 164 IF(ln_ctl) THEN ! print mean trends (used for debugging) 150 165 WRITE(charout, FMT="('opt')") -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcoxy_medusa.F90
r5726 r5841 34 34 ! The following is a map of the subroutines contained within this module 35 35 ! - trc_oxy_medusa 36 ! - CALLS gas_transfer37 36 ! - CALLS oxy_schmidt 38 37 ! - CALLS oxy_sato … … 42 41 !======================================================================= 43 42 ! 44 SUBROUTINE trc_oxy_medusa( pt, ps, uwind, vwind, pp0, o2, dz,& !! inputs45 kw 660, o2flux, o2sat )!! outputs43 SUBROUTINE trc_oxy_medusa( pt, ps, kw660, pp0, o2, & !! inputs 44 kwo2, o2flux, o2sat ) !! outputs 46 45 ! 47 46 !======================================================================= … … 57 56 !! number of oxygen) and oxy_sato.f (calculates oxygen saturation 58 57 !! concentration at 1 atm). 58 !! 59 !! AXY (23/06/15): revised to allow common gas transfer velocity 60 !! to be used for CO2 and O2; outputs of this 61 !! routine amended to mmol/m3 from mol/m3 59 62 !! 60 63 !! Function inputs are (in order) : 61 64 !! pt temperature (degrees C) 62 65 !! ps salinity (PSU) 63 !! uwind u-wind velocity (m/s) 64 !! vwind v-wind velocity (m/s) 66 !! kw660 gas transfer velocity (m/s) 65 67 !! pp0 surface pressure (divided by 1 atm) 66 !! o2 surface O2 concentration (mol/m3) 67 !! dz surface layer thickness (m) 68 !! (*) kw660 gas transfer velocity (m/s) 69 !! (*) o2flux exchange rate of oxygen (mol/m3/s) 70 !! (+) o2sat oxygen saturation concentration (mol/m3) 68 !! o2 surface O2 concentration (mmol/m3) 69 !! (+) kwo2 gas transfer velocity for O2 (m/s) 70 !! (*) o2flux exchange rate of oxygen (mmol/m2/s) 71 !! (+) o2sat oxygen saturation concentration (mmol/m3) 71 72 !! 72 73 !! Where (*) is the function output (note its units). … … 78 79 REAL(wp), INTENT( in ) :: pt 79 80 REAL(wp), INTENT( in ) :: ps 80 REAL(wp), INTENT( in ) :: uwind 81 REAL(wp), INTENT( in ) :: vwind 81 REAL(wp), INTENT( in ) :: kw660 82 82 REAL(wp), INTENT( in ) :: pp0 83 83 REAL(wp), INTENT( in ) :: o2 84 REAL(wp), INTENT( in ) :: dz 85 REAL(wp), INTENT( inout ) :: kw660, o2flux, o2sat 86 ! 87 REAL(wp) :: scl_wind, kwo2, o2schmidt, o2sato 88 ! 89 ! Calculate gas transfer 90 ! 91 call gas_transfer(uwind, vwind, scl_wind, kw660) 84 REAL(wp), INTENT( out ) :: kwo2, o2flux, o2sat 85 ! 86 REAL(wp) :: o2schmidt, o2sato, mol_o2 87 ! 88 ! Oxygen to mol / m3 89 ! 90 mol_o2 = o2 / 1000. 92 91 ! 93 92 ! Calculate oxygen Schmidt number … … 106 105 ! Calculate time rate of change of O2 due to gas exchange (mol/m3/s) 107 106 ! 108 o2flux = kwo2 * (o2sat - o2) / dz 107 o2flux = kwo2 * (o2sat - mol_o2) 108 ! 109 ! Oxygen flux and saturation to mmol / m3 110 ! 111 o2sat = o2sat * 1000. 112 o2flux = o2flux * 1000. 109 113 ! 110 114 END SUBROUTINE trc_oxy_medusa … … 130 134 !! are taken from Keeling et al. (1998, GBC, 12, 141-163). 131 135 !! 136 !! AXY (23/06/2015) 137 !! UPDATED: revised formulation from Wanninkhof (2014) for 138 !! consistency with MOCSY 139 !! 140 !! Winninkhof, R. (2014). Relationship between wind speed and gas 141 !! exchange over the ocean revisited. LIMNOLOGY AND OCEANOGRAPHY-METHODS 142 !! 12, 351-362, doi:10.4319/lom.2014.12.351 143 !! 132 144 !! Function inputs are (in order) : 133 145 !! t temperature (degrees C) … … 141 153 ! 142 154 REAL(wp) :: pt, o2_schmidt 143 REAL(wp) :: a0, a1, a2, a3 144 ! 145 data a0 / 1638.0 / 146 data a1 / -81.83 / 147 data a2 / 1.483 / 148 data a3 / -0.008004 / 149 ! 150 o2_schmidt = a0 + pt*(a1 + pt*(a2 + pt*a3)) 155 REAL(wp) :: a0, a1, a2, a3, a4 156 ! 157 ! AXY (23/06/15): OCMIP-2 coefficients 158 ! data a0 / 1638.0 / 159 ! data a1 / -81.83 / 160 ! data a2 / 1.483 / 161 ! data a3 / -0.008004 / 162 ! 163 ! AXY (23/06/15): Wanninkhof (2014) coefficients 164 data a0 / 1920.4 / 165 data a1 / -135.6 / 166 data a2 / 5.2121 / 167 data a3 / -0.10939 / 168 data a4 / 0.00093777 / 169 ! 170 ! o2_schmidt = a0 + pt*(a1 + pt*(a2 + pt*a3)) 171 o2_schmidt = a0 + pt*(a1 + pt*(a2 + pt*(a3 + pt*a4))) 151 172 ! 152 173 END SUBROUTINE oxy_schmidt … … 231 252 !======================================================================= 232 253 233 !=======================================================================234 !235 SUBROUTINE gas_transfer( uwind, vwind, & !! input236 scl_wind, k ) !! output237 !238 !=======================================================================239 !!240 !! Title : Calculates gas transfer velocity241 !! Author : Andrew Yool242 !! Date : 15/10/04 (revised 04/08/2011)243 !!244 !! This subroutine uses near-surface wind speed to calculate gas245 !! transfer velocity for use in CO2 and O2 exchange calculations.246 !!247 !! Note that the parameterisation of Wanninkhof quoted here is a248 !! truncation of the original equation. It excludes a chemical249 !! enhancement function (based on temperature), although such250 !! temperature dependence is reported negligible by Etcheto &251 !! Merlivat (1988).252 !!253 !! Note also that in calculating scalar wind, the variance of the254 !! wind over the period of a timestep is ignored. Some authors,255 !! for instance OCMIP-2, favour including some reference to the256 !! variability of wind. However, their wind fields are averaged257 !! over relatively long time periods, and so this issue may be258 !! safely (!) ignored here.259 !!260 !! Subroutine inputs are (in order) :261 !! uwind wind u velocity at 10 m (m/s)262 !! vwind wind v velocity at 10 m (m/s)263 !! (+) scl_wind scalar wind velocity at 10 m (m/s)264 !! (*) k gas transfer velocity (m/s)265 !! Where (*) is the function output and (+) is a diagnostic output.266 !!267 !=======================================================================268 269 implicit none270 !271 ! Input variables272 REAL(wp) :: uwind, vwind273 !274 ! Output variables275 REAL(wp) :: scl_wind, k, tmp_k276 !277 ! Choice of parameterisation278 INTEGER :: eqn279 !280 ! Coefficients for various parameterisations281 REAL(wp), DIMENSION(6) :: a282 REAL(wp), DIMENSION(6) :: b283 !284 ! Values of coefficients285 data a(1) / 0.166 / ! Liss & Merlivat (1986) [approximated]286 data a(2) / 0.3 / ! Wanninkhof (1992) [sans enhancement]287 data a(3) / 0.23 / ! Nightingale et al. (2000) [good]288 data a(4) / 0.23 / ! Nightingale et al. (2000) [better]289 data a(5) / 0.222 / ! Nightingale et al. (2000) [best]290 data a(6) / 0.337 / ! OCMIP-2 [sans variability]291 !292 data b(1) / 0.133 /293 data b(2) / 0.0 /294 data b(3) / 0.0 /295 data b(4) / 0.1 /296 data b(5) / 0.333 /297 data b(6) / 0.0 /298 !299 ! Which parameterisation is to be used?300 eqn = 2301 !302 ! Calculate scalar wind (m/s)303 scl_wind = (uwind**2 + vwind**2)**0.5304 !305 ! Calculate gas transfer velocity (cm/h)306 tmp_k = (a(eqn) * scl_wind**2) + (b(eqn) * scl_wind)307 !308 ! Convert tmp_k from cm/h to m/s309 k = tmp_k / (100. * 3600.)310 !311 END SUBROUTINE gas_transfer312 313 !=======================================================================314 !=======================================================================315 !=======================================================================316 317 254 #else 318 255 !!====================================================================== … … 322 259 CONTAINS 323 260 324 SUBROUTINE trc_oxy_medusa( pt, ps, uwind, vwind, pp0, o2, dz,& !! inputs325 kw660, o2flux, o2sat )!! outputs261 SUBROUTINE trc_oxy_medusa( pt, ps, kw660, pp0, o2, & !! inputs 262 o2flux, o2sat ) !! outputs 326 263 USE par_kind 327 264 328 265 REAL(wp), INTENT( in ) :: pt 329 266 REAL(wp), INTENT( in ) :: ps 267 REAL(wp), INTENT( in ) :: kw660 330 268 REAL(wp), INTENT( in ) :: pp0 331 269 REAL(wp), INTENT( in ) :: o2 332 REAL(wp), INTENT( in ) :: dz 333 REAL(wp), INTENT( inout ) :: kw660, o2flux, o2sat 270 REAL(wp), INTENT( inout ) :: o2flux, o2sat 334 271 335 272 WRITE(*,*) 'trc_oxy_medusa: You should not have seen this print! error?', kt -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcsed_medusa.F90
r5742 r5841 40 40 !! AXY (10/02/09) 41 41 LOGICAL, PUBLIC :: & 42 bdustfer = . FALSE.42 bdustfer = .TRUE. 43 43 REAL(wp), PUBLIC :: & 44 44 sedfeinput = 1.e-9_wp , & … … 322 322 CALL prihre( dustmo(:,:,1),jpi,jpj,1,jpi,20,1,jpj,10,1e9,numout ) 323 323 ENDIF 324 324 325 325 ENDIF 326 326 -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcsms_medusa.F90
r5726 r5841 21 21 USE trcopt_medusa 22 22 USE trcsed_medusa 23 USE trcavg_medusa 23 24 24 25 … … 46 47 INTEGER, INTENT(in) :: kt ! ocean time-step index 47 48 49 # if defined key_debug_medusa 50 IF(lwp) WRITE(numout,*) ' MEDUSA inside trc_sms_medusa' 51 CALL flush(numout) 52 # endif 53 48 54 IF( kt == nittrc000 ) THEN 49 55 IF(lwp) WRITE(numout,*) … … 52 58 ENDIF 53 59 60 CALL trc_avg_medusa( kt ) ! rolling average module 61 # if defined key_debug_medusa 62 IF(lwp) WRITE(numout,*) ' MEDUSA done trc_avg_medusa' 63 CALL flush(numout) 64 # endif 65 54 66 CALL trc_opt_medusa( kt ) ! optical model 67 # if defined key_debug_medusa 68 IF(lwp) WRITE(numout,*) ' MEDUSA done trc_opt_medusa' 69 CALL flush(numout) 70 # endif 55 71 56 72 # if defined key_kill_medusa … … 60 76 # else 61 77 CALL trc_bio_medusa( kt ) ! biological model 78 # if defined key_debug_medusa 79 IF(lwp) WRITE(numout,*) ' MEDUSA done trc_bio_medusa' 80 CALL flush(numout) 81 # endif 62 82 63 83 CALL trc_sed_medusa( kt ) ! sedimentation model 84 # if defined key_debug_medusa 85 IF(lwp) WRITE(numout,*) ' MEDUSA done trc_sed_medusa' 86 CALL flush(numout) 87 # endif 64 88 # endif 65 89 -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/par_trc.F90
r5726 r5841 8 8 !! 1.0 ! 2004-03 (C. Ethe) Free form and module 9 9 !! 2.0 ! 2007-12 (C. Ethe, G. Madec) revised architecture 10 !! - ! 2014-06 (A. Yool, J. Palmieri) adding MEDUSA-2 10 11 !!---------------------------------------------------------------------- 11 12 USE par_kind ! kind parameters -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/trcnam.F90
r5726 r5841 11 11 !! - ! 2001-01 (E Kestenare) suppress ndttrc=1 for CEN2 and TVD schemes 12 12 !! 1.0 ! 2005-03 (O. Aumont, A. El Moussaoui) F90 13 !! - ! 2014-06 (A. Yool, J. Palmieri) adding MEDUSA-2 13 14 !!---------------------------------------------------------------------- 14 15 #if defined key_top … … 56 57 !! ** Method : - read passive tracer namelist 57 58 !! - read namelist of each defined SMS model 58 !! ( (PISCES, CFC, MY_TRC )59 !! ( (PISCES, CFC, MY_TRC, MEDUSA, IDTRA ) 59 60 !!--------------------------------------------------------------------- 60 61 INTEGER :: jn, jk ! dummy loop indice … … 231 232 ENDIF 232 233 ! 233 # if defined key_debug_medusa 234 CALL flush(numout) 235 IF (lwp) write (numout,*) '------------------------------' 236 IF (lwp) write (numout,*) 'Jpalm - debug' 237 IF (lwp) write (numout,*) 'CALL trc_nam_pisces -- OK' 234 235 IF( lk_cfc ) THEN ; CALL trc_nam_cfc ! CFC tracers 236 ELSE ; IF(lwp) WRITE(numout,*) ' CFC not used' 237 ENDIF 238 239 IF( lk_c14b ) THEN ; CALL trc_nam_c14b ! C14 bomb tracers 240 ELSE ; IF(lwp) WRITE(numout,*) ' C14 not used' 241 ENDIF 242 243 IF( lk_my_trc ) THEN ; CALL trc_nam_my_trc ! MY_TRC tracers 244 ELSE ; IF(lwp) WRITE(numout,*) ' MY_TRC not used' 245 ENDIF 246 ! 247 # if defined key_debug_medusa 248 CALL flush(numout) 249 IF (lwp) write (numout,*) '------------------------------' 250 IF (lwp) write (numout,*) 'Jpalm - debug' 251 IF (lwp) write (numout,*) 'CALL trc_nam_pisces - CFC - C14 - my_trc -- OK' 238 252 IF (lwp) write (numout,*) 'in trc_nam - just before CALL trc_nam_medusa' 239 253 IF (lwp) write (numout,*) ' ' … … 262 276 IF (lwp) write (numout,*) 'Jpalm - debug' 263 277 IF (lwp) write (numout,*) 'CALL trc_nam_idtra -- OK' 264 IF (lwp) write (numout,*) 'in trc_nam - just before CALL trc_nam_cfc' 265 IF (lwp) write (numout,*) ' ' 266 # endif 267 ! 268 IF( lk_cfc ) THEN ; CALL trc_nam_cfc ! CFC tracers 269 ELSE ; IF(lwp) WRITE(numout,*) ' CFC not used' 270 ENDIF 271 272 IF( lk_c14b ) THEN ; CALL trc_nam_c14b ! C14 bomb tracers 273 ELSE ; IF(lwp) WRITE(numout,*) ' C14 not used' 274 ENDIF 275 276 IF( lk_my_trc ) THEN ; CALL trc_nam_my_trc ! MY_TRC tracers 277 ELSE ; IF(lwp) WRITE(numout,*) ' MY_TRC not used' 278 ENDIF 278 IF (lwp) write (numout,*) 'in trc_nam - CALL trc_nam OK' 279 IF (lwp) write (numout,*) ' ' 280 # endif 279 281 ! 280 282 IF(lwp) CALL flush(numout) -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/trcsms.F90
r5726 r5841 51 51 ! 52 52 IF( lk_pisces ) CALL trc_sms_pisces ( kt ) ! main program of PISCES 53 IF( lk_medusa ) CALL trc_sms_medusa ( kt ) ! MEDUSA tracers54 IF( lk_idtra ) CALL trc_sms_idtra ( kt ) ! radioactive decay of Id. tracer55 53 IF( lk_cfc ) CALL trc_sms_cfc ( kt ) ! surface fluxes of CFC 56 54 IF( lk_c14b ) CALL trc_sms_c14b ( kt ) ! surface fluxes of C14 57 55 IF( lk_my_trc ) CALL trc_sms_my_trc ( kt ) ! MY_TRC tracers 56 IF( lk_medusa ) CALL trc_sms_medusa ( kt ) ! MEDUSA tracers 57 IF( lk_idtra ) CALL trc_sms_idtra ( kt ) ! radioactive decay of Id. tracer 58 58 59 59 IF(ln_ctl) THEN ! print mean trends (used for debugging) -
branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/trcwri.F90
r5726 r5841 5 5 !!====================================================================== 6 6 !! History : 1.0 ! 2009-05 (C. Ethe) Original code 7 !! - ! 2014-06 (A. Yool, J. Palmieri) adding MEDUSA-2 7 8 !!---------------------------------------------------------------------- 8 9 #if defined key_top && defined key_iomput … … 58 59 ! --------------------------------------- 59 60 IF( lk_pisces ) CALL trc_wri_pisces ! PISCES 61 IF( lk_cfc ) CALL trc_wri_cfc ! surface fluxes of CFC 62 IF( lk_c14b ) CALL trc_wri_c14b ! surface fluxes of C14 63 IF( lk_my_trc ) CALL trc_wri_my_trc ! MY_TRC tracers 60 64 ! 61 65 # if defined key_debug_medusa … … 81 85 !!! JPALM 82 86 !!! don't forget to add idtra 83 IF( lk_cfc ) CALL trc_wri_cfc ! surface fluxes of CFC84 IF( lk_c14b ) CALL trc_wri_c14b ! surface fluxes of C1485 IF( lk_my_trc ) CALL trc_wri_my_trc ! MY_TRC tracers86 87 ! 87 88 IF( nn_timing == 1 ) CALL timing_stop('trc_wri')
Note: See TracChangeset
for help on using the changeset viewer.