- Timestamp:
- 2017-03-20T17:21:42+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/nemo_v3_6_STABLE/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
r6963 r7814 27 27 !! - et_i(jpi,jpj) !total ice thermal content 28 28 !! - smt_i(jpi,jpj) !mean ice salinity 29 !! - ot_i(jpi,jpj) !average ice age29 !! - tm_i (jpi,jpj) !mean ice temperature 30 30 !!====================================================================== 31 31 !! History : - ! 2006-01 (M. Vancoppenolle) Original code … … 400 400 ! 401 401 END SUBROUTINE lim_var_salprof 402 402 403 403 404 SUBROUTINE lim_var_bv … … 654 655 INTEGER :: ji, jk, jl ! dummy loop indices 655 656 INTEGER :: ijpij, i_fill, jl0 656 REAL(wp) :: zarg, zV, zconv, zdh 657 REAL(wp) :: zarg, zV, zconv, zdh, zdv 657 658 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zai ! input ice/snow variables 658 659 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zht_i, zht_s, za_i ! output ice/snow variables … … 675 676 IF( zhti(ji) > 0._wp ) THEN 676 677 677 ! initialisation of tests 678 itest(:) = 0 678 ! find which category (jl0) the input ice thickness falls into 679 jl0 = jpl 680 DO jl = 1, jpl 681 IF ( ( zhti(ji) >= hi_max(jl-1) ) .AND. ( zhti(ji) < hi_max(jl) ) ) THEN 682 jl0 = jl 683 CYCLE 684 ENDIF 685 END DO 686 687 ! initialisation of tests 688 itest(:) = 0 679 689 680 i_fill = jpl + 1 !==================================== 681 DO WHILE ( ( SUM( itest(:) ) /= 4 ) .AND. ( i_fill >= 2 ) ) ! iterative loop on i_fill categories 682 ! iteration !==================================== 683 i_fill = i_fill - 1 690 i_fill = jpl + 1 !==================================== 691 DO WHILE ( ( SUM( itest(:) ) /= 4 ) .AND. ( i_fill >= 2 ) ) ! iterative loop on i_fill categories 692 ! iteration !==================================== 693 i_fill = i_fill - 1 694 695 ! initialisation of ice variables for each try 696 zht_i(ji,1:jpl) = 0._wp 697 za_i (ji,1:jpl) = 0._wp 698 itest(:) = 0 699 700 ! *** case very thin ice: fill only category 1 701 IF ( i_fill == 1 ) THEN 702 zht_i(ji,1) = zhti(ji) 703 za_i (ji,1) = zai (ji) 704 705 ! *** case ice is thicker: fill categories >1 706 ELSE 707 708 ! Fill ice thicknesses in the (i_fill-1) cat by hmean 709 DO jl = 1, i_fill - 1 710 zht_i(ji,jl) = hi_mean(jl) 711 END DO 712 713 ! Concentrations in the (i_fill-1) categories 714 za_i(ji,jl0) = zai(ji) / SQRT(REAL(jpl)) 715 DO jl = 1, i_fill - 1 716 IF ( jl /= jl0 ) THEN 717 zarg = ( zht_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp ) 718 za_i(ji,jl) = za_i (ji,jl0) * EXP(-zarg**2) 719 ENDIF 720 END DO 721 722 ! Concentration in the last (i_fill) category 723 za_i(ji,i_fill) = zai(ji) - SUM( za_i(ji,1:i_fill-1) ) 724 725 ! Ice thickness in the last (i_fill) category 726 zV = SUM( za_i(ji,1:i_fill-1) * zht_i(ji,1:i_fill-1) ) 727 zht_i(ji,i_fill) = ( zhti(ji) * zai(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 ) 728 729 ! clem: correction if concentration of upper cat is greater than lower cat 730 ! (it should be a gaussian around jl0 but sometimes it is not) 731 IF ( jl0 /= jpl ) THEN 732 DO jl = jpl, jl0+1, -1 733 IF ( za_i(ji,jl) > za_i(ji,jl-1) ) THEN 734 zdv = zht_i(ji,jl) * za_i(ji,jl) 735 zht_i(ji,jl ) = 0._wp 736 za_i (ji,jl ) = 0._wp 737 za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 ) 738 END IF 739 ENDDO 740 ENDIF 741 742 ENDIF ! case ice is thick or thin 684 743 685 ! initialisation of ice variables for each try 686 zht_i(ji,1:jpl) = 0._wp 687 za_i (ji,1:jpl) = 0._wp 688 itest(:) = 0 689 690 ! *** case very thin ice: fill only category 1 691 IF ( i_fill == 1 ) THEN 692 zht_i(ji,1) = zhti(ji) 693 za_i (ji,1) = zai (ji) 694 695 ! *** case ice is thicker: fill categories >1 696 ELSE 697 698 ! Fill ice thicknesses except the last one (i_fill) by hmean 699 DO jl = 1, i_fill - 1 700 zht_i(ji,jl) = hi_mean(jl) 701 END DO 744 !--------------------- 745 ! Compatibility tests 746 !--------------------- 747 ! Test 1: area conservation 748 zconv = ABS( zai(ji) - SUM( za_i(ji,1:jpl) ) ) 749 IF ( zconv < epsi06 ) itest(1) = 1 750 751 ! Test 2: volume conservation 752 zconv = ABS( zhti(ji)*zai(ji) - SUM( za_i(ji,1:jpl)*zht_i(ji,1:jpl) ) ) 753 IF ( zconv < epsi06 ) itest(2) = 1 702 754 703 ! find which category (jl0) the input ice thickness falls into 704 jl0 = i_fill 755 ! Test 3: thickness of the last category is in-bounds ? 756 IF ( zht_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1 757 758 ! Test 4: positivity of ice concentrations 759 itest(4) = 1 705 760 DO jl = 1, i_fill 706 IF ( ( zhti(ji) >= hi_max(jl-1) ) .AND. ( zhti(ji) < hi_max(jl) ) ) THEN 707 jl0 = jl 708 CYCLE 709 ENDIF 710 END DO 711 712 ! Concentrations in the (i_fill-1) categories 713 za_i(ji,jl0) = zai(ji) / SQRT(REAL(jpl)) 714 DO jl = 1, i_fill - 1 715 IF ( jl == jl0 ) CYCLE 716 zarg = ( zht_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp ) 717 za_i(ji,jl) = za_i (ji,jl0) * EXP(-zarg**2) 718 END DO 719 720 ! Concentration in the last (i_fill) category 721 za_i(ji,i_fill) = zai(ji) - SUM( za_i(ji,1:i_fill-1) ) 722 723 ! Ice thickness in the last (i_fill) category 724 zV = SUM( za_i(ji,1:i_fill-1) * zht_i(ji,1:i_fill-1) ) 725 zht_i(ji,i_fill) = ( zhti(ji) * zai(ji) - zV ) / za_i(ji,i_fill) 726 727 ENDIF ! case ice is thick or thin 728 729 !--------------------- 730 ! Compatibility tests 731 !--------------------- 732 ! Test 1: area conservation 733 zconv = ABS( zai(ji) - SUM( za_i(ji,1:jpl) ) ) 734 IF ( zconv < epsi06 ) itest(1) = 1 735 736 ! Test 2: volume conservation 737 zconv = ABS( zhti(ji)*zai(ji) - SUM( za_i(ji,1:jpl)*zht_i(ji,1:jpl) ) ) 738 IF ( zconv < epsi06 ) itest(2) = 1 739 740 ! Test 3: thickness of the last category is in-bounds ? 741 IF ( zht_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1 742 743 ! Test 4: positivity of ice concentrations 744 itest(4) = 1 745 DO jl = 1, i_fill 746 IF ( za_i(ji,jl) < 0._wp ) itest(4) = 0 747 END DO 748 !============================ 749 END DO ! end iteration on categories 750 !============================ 761 IF ( za_i(ji,jl) < 0._wp ) itest(4) = 0 762 END DO 763 ! !============================ 764 END DO ! end iteration on categories 765 ! !============================ 751 766 ENDIF ! if zhti > 0 752 767 END DO ! i loop
Note: See TracChangeset
for help on using the changeset viewer.
