Changeset 8435 for branches/UKMO
- Timestamp:
- 2017-08-14T12:19:09+02:00 (7 years ago)
- Location:
- branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO
- Files:
-
- 4 edited
-
LIM_SRC_2/dom_ice_2.F90 (modified) (1 diff)
-
OPA_SRC/ASM/bias.F90 (modified) (16 diffs)
-
OPA_SRC/DYN/dynhpg.F90 (modified) (2 diffs)
-
OPA_SRC/step.F90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/LIM_SRC_2/dom_ice_2.F90
r7065 r8435 1 1 MODULE dom_ice_2 2 !test3 2 !!====================================================================== 4 3 !! *** MODULE dom_ice *** -
branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/ASM/bias.F90
r7667 r8435 211 211 !! * Local declarations 212 212 213 CHARACTER(len=100) :: cn_dir ! dir for location ofline bias213 CHARACTER(len=100) :: cn_dir ! dir for location ofline bias 214 214 INTEGER :: ierror 215 INTEGER :: ios ! Local integer output status for namelist read216 REAL(wp) :: eft_rlx, & ! efolding time (bias memory) in days217 & eft_asm, & ! " "215 INTEGER :: ios ! Local integer output status for namelist read 216 REAL(wp) :: eft_rlx, & ! efolding time (bias memory) in days 217 & eft_asm, & ! " " 218 218 & log2, & 219 & lenscl_bias, & ! lengthscale of the pressure bias decay between minlat and maxlat.220 & minlat_bias, & !used in ipc221 & maxlat_bias !used in ipc219 & lenscl_bias, & ! lengthscale pressure bias decay between minlat and maxlat. 220 & minlat_bias, & ! used in ipc 221 & maxlat_bias ! used in ipc 222 222 223 223 NAMELIST/nambias/ ln_bias, ln_bias_asm, ln_bias_rlx, ln_bias_ofl, & … … 396 396 tbias_i(:,:,:) = 0.0_wp 397 397 sbias_i(:,:,:) = 0.0_wp 398 IF(lwp) WRITE(numout,*) 'gru_pc and grv_pc are set to zero '399 398 gru_pc(:,:) = 0.0_wp 400 399 grv_pc(:,:) = 0.0_wp 400 401 401 IF ( ln_bias_rlx ) THEN 402 402 tbias_rlx(:,:,:) = 0.0_wp 403 403 sbias_rlx(:,:,:) = 0.0_wp 404 404 ENDIF 405 405 406 IF ( ln_bias_asm ) THEN !now rlx and asm bias in same file 406 tbias_asm(:,:,:) = 0.0_wp407 sbias_asm(:,:,:) = 0.0_wp408 tbias_asm_out(:,:,:) = 0.0_wp409 sbias_asm_out(:,:,:) = 0.0_wp407 tbias_asm(:,:,:) = 0.0_wp 408 sbias_asm(:,:,:) = 0.0_wp 409 tbias_asm_out(:,:,:) = 0.0_wp 410 sbias_asm_out(:,:,:) = 0.0_wp 410 411 ENDIF 411 412 412 413 IF ( ln_incpc ) THEN !incr pressure correction 413 tbias_asm_stscale(:,:,:) = 0.0_wp414 sbias_asm_stscale(:,:,:) = 0.0_wp415 tbias_asm_stscale_out(:,:,:) = 0.0_wp416 sbias_asm_stscale_out(:,:,:) = 0.0_wp414 tbias_asm_stscale(:,:,:) = 0.0_wp 415 sbias_asm_stscale(:,:,:) = 0.0_wp 416 tbias_asm_stscale_out(:,:,:) = 0.0_wp 417 sbias_asm_stscale_out(:,:,:) = 0.0_wp 417 418 ENDIF 418 419 … … 458 459 ! Get the T and S bias data 459 460 CALL iom_get( numbias_tot, jpdom_autoglo, 'tbias_asm_stscale', tbias_asm_stscale ) 460 461 CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_asm_stscale', sbias_asm_stscale ) 461 462 ELSE 462 463 CALL ctl_stop( 'Short time scale bias assim variables not found in ',cn_bias_tot ) 463 464 ENDIF 464 465 ENDIF … … 550 551 551 552 IF ( ln_incpc) THEN 552 ! not sure if this should be a special case of nn_lat_ramp == 1553 553 minlat_bias = 3.0_wp 554 554 maxlat_bias = 8.0_wp … … 616 616 & fbcoef(:,:) * fb_t_asm_max ) 617 617 zcof2(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_p_asm 618 619 618 620 619 IF ( ln_itdecay ) THEN !decay the pressure correction at each time step 621 620 622 ztsday = rday / real(rdt)623 624 zdecay = (1-ztscale)**(1/ real(ztsday)) ! used in ipc625 zfrac1 = zdecay** real(kt) ! used in ipc621 ztsday = rday / REAL(rdt,wp) 622 623 zdecay = (1-ztscale)**(1/REAL(ztsday,wp)) ! used in ipc 624 zfrac1 = zdecay**REAL(kt,wp) ! used in ipc 626 625 IF ( zfrac1 <= 0.0 ) zfrac1 = 0.0_wp 627 626 628 629 627 IF( ln_asmiau .and. ln_trainc ) THEN ! nudge in increments and decay historical contribution 630 631 628 632 629 IF ( kt <= nitiaufin ) THEN ! During IAU calculate the fraction of increments to apply at each time step 633 630 634 ztfrac = real(kt) / real(nitiaufin) ! nudging factor during the IAU 635 631 ztfrac = REAL(kt,wp) / REAL(nitiaufin,wp) ! nudging factor during the IAU 636 632 637 633 IF (lwp) THEN … … 639 635 WRITE(numout,*) '~~~~~~~~~~~~' 640 636 WRITE(numout,* ) "proportion of increment applied in pcbias ",ztfrac 641 WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**( real(kt)/ztsday)637 WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(REAL(kt,wp)/ztsday) 642 638 ENDIF 643 639 644 645 DO jk = 1, jpkm1 646 640 DO jk = 1, jpkm1 647 641 tbias(:,:,jk) = tbias(:,:,jk) + & 648 & ( t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) + &642 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) + & 649 643 & ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:) 650 644 sbias(:,:,jk) = sbias(:,:,jk) + & 651 & ( t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) + &645 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) + & 652 646 & ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:) 653 647 654 tbias_p(:,:,jk) = tbias_p(:,:,jk) + &655 & ( t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) + &648 tbias_p(:,:,jk) = tbias_p(:,:,jk) + & 649 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) + & 656 650 & ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:) 657 651 sbias_p(:,:,jk) = sbias_p(:,:,jk) + & 658 & ( t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) + &652 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) + & 659 653 & ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:) 660 654 ENDDO 661 655 662 663 656 IF (ln_incpc) THEN 664 665 657 666 658 IF (lwp) THEN … … 687 679 688 680 DO jk = 1, jpkm1 689 tbias_asm_out(:,:,jk) = t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) + &681 tbias_asm_out(:,:,jk) = t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) + & 690 682 & ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) 691 sbias_asm_out(:,:,jk) = t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) + &683 sbias_asm_out(:,:,jk) = t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) + & 692 684 & ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) 693 685 END DO … … 707 699 IF ( kt == nitiaufin ) THEN 708 700 DO jk = 1, jpkm1 709 tbias_asm(:,:,jk) = t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) + &701 tbias_asm(:,:,jk) = t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) + & 710 702 & ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) 711 sbias_asm(:,:,jk) = t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) + &703 sbias_asm(:,:,jk) = t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) + & 712 704 & ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) 713 705 END DO … … 721 713 722 714 ENDIF 723 724 715 725 716 ELSE ! decay pressure correction from combined historical component and increments after IAU 726 717 727 ztfrac=t_asm_mem**(real(kt-nitiaufin)/real(nitiaufin)) ! decay from end of IAU 728 718 ztfrac=t_asm_mem**(REAL(kt-nitiaufin,wp)/REAL(nitiaufin,wp)) ! decay from end of IAU 729 719 730 720 DO jk = 1, jpkm1 … … 737 727 sbias_p(:,:,jk) = sbias_p(:,:,jk) + & 738 728 & ( ztfrac * sbias_asm(:,:,jk) ) * zcof2(:,:) 739 740 729 ENDDO 741 730 742 731 IF (ln_incpc) THEN 743 732 744 zfrac = zdecay** real(kt-nitiaufin)733 zfrac = zdecay**REAL(kt-nitiaufin,wp) 745 734 IF ( zfrac <= 0.0 ) zfrac = 0.0_wp 746 735 747 748 736 DO jk = 1, jpkm1 749 737 tbias_i(:,:,jk) = tbias_asm_stscale(:,:,jk) * zfrac * (1.0 - fbcoef_stscale(:,:)) … … 752 740 753 741 IF (lwp) THEN 754 WRITE(numout,*) 'tra_bias : bias weights'755 WRITE(numout,*) '~~~~~~~~~~~~'756 WRITE(numout,* ) "IPC - proportion of increments + historical pcbias applied ",zfrac742 WRITE(numout,*) 'tra_bias : bias weights' 743 WRITE(numout,*) '~~~~~~~~~~~~' 744 WRITE(numout,* ) "IPC - proportion of increments + historical pcbias applied ",zfrac 757 745 ENDIF 758 746 759 747 ENDIF 760 748 761 IF (lwp) THEN 762 WRITE(numout,*) 'tra_bias : bias weights' 763 WRITE(numout,*) '~~~~~~~~~~~~' 764 WRITE(numout,* ) "proportion of increments + historical pcbias applied ",ztfrac 765 ENDIF 766 767 IF ( ln_trainc .and. .not.ln_bsyncro ) THEN 768 IF ( kt == nn_bias_itwrt ) THEN 769 DO jk = 1, jpkm1 770 tbias_asm_out(:,:,jk) = ztfrac * tbias_asm(:,:,jk) 771 sbias_asm_out(:,:,jk) = ztfrac * sbias_asm(:,:,jk) 772 END DO 773 774 IF ( ln_incpc ) THEN 775 IF (lwp) WRITE(numout,*) 'after end of IAU - IPC - saving s/tbias_asm_stscale' 776 DO jk = 1, jpkm1 777 tbias_asm_stscale_out(:,:,jk) = tbias_asm_stscale(:,:,jk) * zfrac 778 sbias_asm_stscale_out(:,:,jk) = sbias_asm_stscale(:,:,jk) * zfrac 779 ENDDO 780 ENDIF 781 782 ENDIF 783 784 ENDIF 785 786 ENDIF 787 749 IF (lwp) THEN 750 WRITE(numout,*) 'tra_bias : bias weights' 751 WRITE(numout,*) '~~~~~~~~~~~~' 752 WRITE(numout,* ) "proportion of increments + historical pcbias applied ",ztfrac 753 ENDIF 754 755 IF ( ln_trainc .and. .not.ln_bsyncro ) THEN 756 IF ( kt == nn_bias_itwrt ) THEN 757 DO jk = 1, jpkm1 758 tbias_asm_out(:,:,jk) = ztfrac * tbias_asm(:,:,jk) 759 sbias_asm_out(:,:,jk) = ztfrac * sbias_asm(:,:,jk) 760 END DO 761 762 IF ( ln_incpc ) THEN 763 IF (lwp) WRITE(numout,*) 'after end of IAU - IPC - saving s/tbias_asm_stscale' 764 DO jk = 1, jpkm1 765 tbias_asm_stscale_out(:,:,jk) = tbias_asm_stscale(:,:,jk) * zfrac 766 sbias_asm_stscale_out(:,:,jk) = sbias_asm_stscale(:,:,jk) * zfrac 767 ENDDO 768 ENDIF 769 ENDIF 770 ENDIF 771 ENDIF 788 772 789 773 ELSE ! no assimilation increments, simply decay pressure correction (e.g for forecasts) … … 792 776 DO jk = 1, jpkm1 793 777 tbias(:,:,jk) = tbias(:,:,jk) + & 794 & ( t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) ) * zcof1(:,:)778 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) ) * zcof1(:,:) 795 779 sbias(:,:,jk) = sbias(:,:,jk) + & 796 & ( t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) ) * zcof1(:,:)780 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) ) * zcof1(:,:) 797 781 tbias_p(:,:,jk) = tbias_p(:,:,jk) + & 798 & ( t_asm_mem**( real(kt)/ztsday) * tbias_asm(:,:,jk) ) * zcof2(:,:)782 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) ) * zcof2(:,:) 799 783 sbias_p(:,:,jk) = sbias_p(:,:,jk) + & 800 & ( t_asm_mem**( real(kt)/ztsday) * sbias_asm(:,:,jk) ) * zcof2(:,:)784 & ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) ) * zcof2(:,:) 801 785 ENDDO 802 786 … … 804 788 WRITE(numout,*) 'tra_bias : bias weights' 805 789 WRITE(numout,*) '~~~~~~~~~~~~' 806 WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**( real(kt)/ztsday)790 WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(REAL(kt,wp)/ztsday) 807 791 ENDIF 808 809 810 792 811 793 IF (ln_incpc) THEN … … 957 939 tsw(:,:,:,jp_sal) = tsa(:,:,:,jp_sal) - sbias_p(:,:,:) - sbias_i(:,:,:) 958 940 IF ( ln_incpc_only ) THEN 959 IF(lwp) WRITE(numout,*) 'SUM(tbias_i) is', SUM(tbias_i), 'before to be added to tsw'960 941 IF(lwp) WRITE(numout,*) 'ln_incpc_only =', ln_incpc_only, 'tsw updated with IPC only and ln_dynhpg_imp = ',ln_dynhpg_imp 961 942 tsw(:,:,:,jp_tem) = tsa(:,:,:,jp_tem) - tbias_i(:,:,:) … … 970 951 tsw(:,:,:,jp_sal) = tsb(:,:,:,jp_sal) - sbias_p(:,:,:) - sbias_i(:,:,:) 971 952 IF ( ln_incpc_only ) THEN 972 IF(lwp) WRITE(numout,*) 'SUM(tbias_i) is', SUM(tbias_i), 'before to be added to tsw'973 953 IF(lwp) WRITE(numout,*) 'ln_incpc_only =', ln_incpc_only, 'tsw updated with IPC only and ln_dynhpg_imp = ',ln_dynhpg_imp 974 954 tsw(:,:,:,jp_tem) = tsb(:,:,:,jp_tem) - tbias_i(:,:,:) … … 978 958 ENDIF 979 959 980 981 982 960 IF(lwp) WRITE(numout,*) 'dyn_bias( kt ) calculating rhd_pc, kt =', kt 983 !! GO5: CALL eos( tsw, rhd_pc, rhop ) 984 CALL eos( tsw, rhd_pc, rhop, fsdept_n(:,:,:) ) 961 CALL eos( tsw, rhd_pc, rhop, fsdept_n(:,:,:) ) 985 962 986 ! is this needed?987 IF(lwp) WRITE(numout,*) 'CALL lbc_lnk( rhd_pc, T, 1.0_wp )'988 963 CALL lbc_lnk( rhd_pc, 'T', 1.0_wp ) 989 964 990 965 ! Partial steps: now horizontal gradient of t,s,rd 991 966 ! at the bottom ocean level 992 IF(lwp) WRITE(numout,*) 'horizontal gradient of t,s,rd ' 967 993 968 IF( ln_zps ) THEN 994 969 CALL zps_hde( kt, jpts, tsw, gtsu, gtsv, & -
branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90
r7658 r8435 99 99 ! 100 100 IF ( ln_bias .AND. ln_bias_pc_app ) THEN 101 IF(lwp) THEN102 WRITE(numout,*) " ! store orig density, use pressure corrected density using gru_pc and grv_pc"103 ENDIF104 101 z_rhd_st(:,:,:) = rhd(:,:,:) ! store orig density 105 102 rhd(:,:,:) = rhd_pc(:,:,:) ! use pressure corrected density … … 110 107 ENDIF 111 108 112 IF(lwp) WRITE(numout,*) " computing pressure gradient with nhpg = ",nhpg113 109 SELECT CASE ( nhpg ) ! Hydrostatic pressure gradient computation 114 110 CASE ( 0 ) ; CALL hpg_zco ( kt ) ! z-coordinate -
branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/step.F90
r7665 r8435 273 273 CALL tra_nxt( kstp ) ! tracer fields at next time step 274 274 IF( ln_sto_eos ) CALL sto_pts( tsn ) ! Random T/S fluctuations 275 276 275 CALL eos ( tsa, rhd, rhop, fsdept_n(:,:,:) ) ! Time-filtered in situ density for hpg computation 276 CALL lbc_lnk( rhd, 'T', 1.0_wp ) 277 277 IF( ln_zps .AND. .NOT. ln_isfcav) & 278 278 & CALL zps_hde ( kstp, jpts, tsa, gtsu, gtsv, & ! Partial steps: before horizontal gradient … … 283 283 & gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi ) ! of t, s, rd at the last ocean level 284 284 IF( ln_bias ) CALL dyn_bias( kstp ) 285 ELSE ! centered hpg (eos then time stepping)286 IF ( .NOT. lk_dynspg_ts ) THEN ! eos already called in time-split case287 IF( ln_sto_eos ) CALL sto_pts( tsn ) ! Random T/S fluctuations288 289 285 ELSE ! centered hpg (eos then time stepping) 286 IF ( .NOT. lk_dynspg_ts ) THEN ! eos already called in time-split case 287 IF( ln_sto_eos ) CALL sto_pts( tsn ) ! Random T/S fluctuations 288 CALL eos ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 289 CALL lbc_lnk( rhd, 'T', 1.0_wp ) 290 290 IF( ln_zps .AND. .NOT. ln_isfcav) & 291 291 & CALL zps_hde ( kstp, jpts, tsn, gtsu, gtsv, & ! Partial steps: before horizontal gradient
Note: See TracChangeset
for help on using the changeset viewer.
