Changeset 3963 for branches/2013/dev_r3406_CNRS_LIM3
- Timestamp:
- 2013-07-09T17:41:20+02:00 (11 years ago)
- Location:
- branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO
- Files:
-
- 17 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limdiahsb.F90
r3938 r3963 127 127 CALL iom_put( 'ibgfrcemp',frc_vol * 1.e-9 ) !vol - forcing : km3 (equivalent liquid) 128 128 CALL iom_put( 'ibgfrcemps',frc_sal * 1.e-9 ) !sal - forcing : psu*km3 129 CALL iom_put( 'ibgemps',zbg_emps *86400 )130 CALL iom_put( 'ibgemp',zbg_emp *86400 )131 CALL iom_put( 'ibgfsbri',zbg_fsbri *86400 )132 CALL iom_put( 'ibgfseqv',zbg_fseqv *86400 )133 CALL iom_put( 'ibgfsaltres',zbg_fsalt_res *86400 )134 CALL iom_put( 'ibgfsaltrpo',zbg_fsalt_rpo *86400 )129 CALL iom_put( 'ibgemps',zbg_emps *86400.) 130 CALL iom_put( 'ibgemp',zbg_emp *86400.) 131 CALL iom_put( 'ibgfsbri',zbg_fsbri *86400.) 132 CALL iom_put( 'ibgfseqv',zbg_fseqv *86400.) 133 CALL iom_put( 'ibgfsaltres',zbg_fsalt_res *86400.) 134 CALL iom_put( 'ibgfsaltrpo',zbg_fsalt_rpo *86400.) 135 135 CALL iom_put( 'ibggrpme',bg_grpme * rhoic / rau0 * 1.e-9 ) ! km3 (equivalent liquid) 136 136 ! -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90
r3938 r3963 66 66 REAL(wp), POINTER, DIMENSION(:,:) :: zu_io, zv_io ! ice-ocean velocity 67 67 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 68 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 68 69 !!--------------------------------------------------------------------- 69 70 … … 234 235 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 235 236 237 zchk_vmin = glob_min(v_i) 238 zchk_amax = glob_max(SUM(a_i,dim=3)) 239 zchk_amin = glob_min(a_i) 240 236 241 IF(lwp) THEN 237 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limdyn) = ',(zchk_v_i * 86400.) 238 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limdyn) = ',(zchk_smv * 86400.) 239 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limdyn) = ',(MINVAL(v_i) * 1.e-3) 240 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+1.e-10 ) WRITE(numout,*) 'violation a_i>amax (limdyn) = ',MAXVAL(SUM(a_i,dim=3)) 242 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limdyn) = ',(zchk_v_i * 86400.) 243 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limdyn) = ',(zchk_smv * 86400.) 244 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limdyn) = ',(zchk_vmin * 1.e-3) 245 !IF ( zchk_amax > amax+1.e-10 ) WRITE(numout,*) 'violation a_i>amax (limdyn) = ',zchk_amax 246 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limdyn) = ',zchk_amin 241 247 ENDIF 242 248 !CALL ctl_opn( numhsb , cl_name , 'UNKNOWN' , 'FORMATTED' , 'SEQUENTIAL' , 1 , numout , lwp , 1 ) -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r3938 r3963 142 142 REAL(wp), POINTER, DIMENSION(:,:) :: vt_i_init, vt_i_final ! ice volume summed over categories 143 143 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 144 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 144 145 ! mass and salt flux (clem) 145 146 REAL(wp), POINTER, DIMENSION(:,:,:) :: zviold, zvsold, zsmvold ! old ice volume... … … 473 474 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 474 475 476 zchk_vmin = glob_min(v_i) 477 zchk_amax = glob_max(SUM(a_i,dim=3)) 478 zchk_amin = glob_min(a_i) 479 475 480 IF(lwp) THEN 476 IF ( 477 IF ( 478 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_me) = ',(MINVAL(v_i)* 1.e-3)479 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_me) = ',MAXVAL(SUM(a_i,dim=3))480 IF ( MINVAL( a_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation a_i<0 (limitd_me) = ',MINVAL(a_i)481 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limitd_me) = ',(zchk_v_i * 86400.) 482 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limitd_me) = ',(zchk_smv * 86400.) 483 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_me) = ',(zchk_vmin * 1.e-3) 484 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_me) = ',zchk_amax 485 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limitd_me) = ',zchk_amin 481 486 ENDIF 482 487 ENDIF -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90
r3938 r3963 68 68 INTEGER :: jl, ja, jm, jbnd1, jbnd2 ! ice types dummy loop index 69 69 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 70 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 70 71 !!------------------------------------------------------------------ 71 72 ! ------------------------------- … … 165 166 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 166 167 168 zchk_vmin = glob_min(v_i) 169 zchk_amax = glob_max(SUM(a_i,dim=3)) 170 zchk_amin = glob_min(a_i) 171 167 172 IF(lwp) THEN 168 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limitd_th) = ',(zchk_v_i * 86400.) 169 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limitd_th) = ',(zchk_smv * 86400.) 170 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_th) = ',(MINVAL(v_i) * 1.e-3) 171 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_th) = ',MAXVAL(SUM(a_i,dim=3)) 173 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limitd_th) = ',(zchk_v_i * 86400.) 174 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limitd_th) = ',(zchk_smv * 86400.) 175 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limitd_th) = ',(zchk_vmin * 1.e-3) 176 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limitd_th) = ',zchk_amax 177 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limitd_th) = ',zchk_amin 172 178 ENDIF 173 179 ENDIF … … 206 212 INTEGER :: zji, zjj, nd ! local integer 207 213 REAL(wp) :: zx1, zwk1, zdh0, zetamin, zdamax ! local scalars 208 REAL(wp) :: zx2, zwk2, zda0, zetamax , zhimin! - -214 REAL(wp) :: zx2, zwk2, zda0, zetamax ! - - 209 215 REAL(wp) :: zx3, zareamin, zindb ! - - 210 216 CHARACTER (len = 15) :: fieldid … … 242 248 CALL wrk_alloc( jpi,jpj, zhb0,zhb1,vt_i_init,vt_i_final,vt_s_init,vt_s_final,et_i_init,et_i_final,et_s_init,et_s_final ) 243 249 244 zhimin = 0.1 !minimum ice thickness tolerated by the model245 250 zareamin = epsi10 !minimum area in thickness categories tolerated by the conceptors of the model 246 251 … … 527 532 zji = nind_i(ji) 528 533 zjj = nind_j(ji) 529 IF ( ( zhimin .GT. 0.0 ) .AND. & 530 ( ( a_i(zji,zjj,1) .GT. epsi10 ) .AND. ( ht_i(zji,zjj,1) .LT. zhimin ) ) & 531 ) THEN 532 a_i(zji,zjj,1) = a_i(zji,zjj,1) * ht_i(zji,zjj,1) / zhimin 533 ht_i(zji,zjj,1) = zhimin 534 v_i(zji,zjj,1) = a_i(zji,zjj,1)*ht_i(zji,zjj,1) !clem@useless 534 IF ( ( a_i(zji,zjj,1) .GT. epsi10 ) .AND. ( ht_i(zji,zjj,1) .LT. hiclim ) ) THEN 535 a_i(zji,zjj,1) = a_i(zji,zjj,1) * ht_i(zji,zjj,1) / hiclim 536 ht_i(zji,zjj,1) = hiclim 537 v_i(zji,zjj,1) = a_i(zji,zjj,1) * ht_i(zji,zjj,1) !clem@useless 535 538 ENDIF 536 539 END DO !ji -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r3938 r3963 94 94 REAL(wp), POINTER, DIMENSION(:,:) :: zqlbsbq ! link with lead energy budget qldif 95 95 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 96 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 96 97 !!------------------------------------------------------------------- 97 98 … … 476 477 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 477 478 479 zchk_vmin = glob_min(v_i) 480 zchk_amax = glob_max(SUM(a_i,dim=3)) 481 zchk_amin = glob_min(a_i) 482 478 483 IF(lwp) THEN 479 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limthd) = ',(zchk_v_i * 86400.) 480 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limthd) = ',(zchk_smv * 86400.) 481 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limthd) = ',(MINVAL(v_i) * 1.e-3) 482 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limthd) = ',MAXVAL(SUM(a_i,dim=3)) 484 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limthd) = ',(zchk_v_i * 86400.) 485 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limthd) = ',(zchk_smv * 86400.) 486 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limthd) = ',(zchk_vmin * 1.e-3) 487 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limthd) = ',zchk_amax 488 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limthd) = ',zchk_amin 483 489 ENDIF 484 490 ENDIF -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r3938 r3963 121 121 122 122 ! mass and salt flux (clem) 123 REAL(wp) :: zdvres 123 REAL(wp) :: zdvres, zdvsur, zdvbot 124 124 REAL(wp), POINTER, DIMENSION(:) :: zviold, zvsold ! old ice volume... 125 125 … … 255 255 zhn = 1.0 - MAX( zzero , SIGN( zone , - zhsnew ) ) 256 256 ht_s_b(ji) = MAX( zzero , zhsnew ) 257 ! we recompute dh_s_tot (clem) 258 dh_s_tot (ji) = ht_s_b(ji) - zhsold(ji) 257 259 ! Volume and mass variations of snow 258 260 ! dvsbq_1d (ji) = a_i_b(ji) * ( ht_s_b(ji) - zhsold(ji) - zdh_s_mel(ji) ) … … 283 285 zdq_i (ji) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) / rdt_ice 284 286 ! 285 ! IOVINO 286 !zfsalt_melt(ji) = zfsalt_melt(ji) - sm_i_b(ji) * a_i_b(ji) & 287 ! & * MIN( zdeltah(ji,jk) , 0.e0 ) * rhoic / rdt_ice 287 ! clem 288 288 fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) & 289 289 & * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic / rdt_ice 290 ! fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdeltah(ji,jk) * rhoic / rdt_ice 290 291 END DO 291 292 END DO … … 449 450 ! Salinity update 450 451 ! entrapment during bottom growth 451 dsm_i_se_1d(ji) = ( s_i_new(ji) * dh_i_bott(ji) + sm_i_b(ji) * ht_i_b(ji) ) &452 & / MAX( ht_i_b(ji) + dh_i_bott(ji) ,epsi13 ) - sm_i_b(ji)452 !clem:movedown dsm_i_se_1d(ji) = ( s_i_new(ji) * dh_i_bott(ji) + sm_i_b(ji) * ht_i_b(ji) ) & 453 !clem:movedown & / MAX( ht_i_b(ji) + dh_i_bott(ji) ,epsi13 ) - sm_i_b(ji) 453 454 fseqv_1d(ji) = fseqv_1d(ji) - s_i_new(ji) * a_i_b(ji) * dh_i_bott(ji) * rhoic / rdt_ice 454 455 ENDIF ! heat budget … … 489 490 zdq_i(ji) = zdq_i(ji) + zdeltah(ji,jk) * q_i_b(ji,jk) / rdt_ice 490 491 ENDIF 491 ! IOVINO contribution to salt flux 492 !zfsalt_melt(ji) = zfsalt_melt(ji) - sm_i_b(ji) * a_i_b(ji) & 493 ! & * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic / rdt_ice 492 ! clem 494 493 fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) & 495 494 & * MIN( zdeltah(ji,jk) , 0._wp ) * rhoic / rdt_ice 495 ! fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdeltah(ji,jk) * rhoic / rdt_ice 496 496 ENDIF 497 497 END DO ! ji … … 549 549 ! ! excessive energy is sent to lateral ablation 550 550 fsup (ji) = rhoic * lfus * at_i_b(ji) / MAX( 1.0 - at_i_b(ji) , epsi13 ) * zdvres / rdt_ice 551 ! !since ice volume is only used for outputs, we keep it global for all categories552 dvbbq_1d (ji) = a_i_b(ji) * dh_i_bott(ji)553 ! !new ice thickness554 zhgnew (ji) = ht_i_b(ji) + dh_i_surf(ji) + dh_i_bott(ji)555 551 END DO 556 552 … … 564 560 ! Adapt the remaining energy if too much ice melts 565 561 !-------------------------------------------------- 566 zdvres = MAX( 0._wp, - zhgnew(ji) ) 567 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) !1 if ice 568 dh_i_bott (ji) = dh_i_bott(ji) + zdvres ! clem@bug 569 zhgnew (ji) = zihgnew * zhgnew(ji) ! ice thickness is put to 0 570 ! remaining heat 562 zdvres = MAX( 0._wp, - ht_i_b(ji) - dh_i_surf(ji) - dh_i_bott(ji) ) 563 zdvsur = MIN( 0._wp, dh_i_surf(ji) + zdvres ) - dh_i_surf(ji) ! fill the surface first 564 zdvbot = MAX( 0._wp, zdvres - zdvsur ) ! then the bottom 565 dh_i_surf (ji) = dh_i_surf(ji) + zdvsur ! clem 566 dh_i_bott (ji) = dh_i_bott(ji) + zdvbot ! clem 567 568 ! new ice thickness (clem) 569 zhgnew(ji) = ht_i_b(ji) + dh_i_surf(ji) + dh_i_bott(ji) 570 zihgnew = 1.0 - MAX( zzero , SIGN( zone , - zhgnew(ji) ) ) !1 if ice 571 zhgnew(ji) = zihgnew * zhgnew(ji) ! ice thickness is put to 0 572 573 ! !since ice volume is only used for outputs, we keep it global for all categories 574 dvbbq_1d (ji) = a_i_b(ji) * dh_i_bott(ji) 575 576 ! remaining heat 571 577 zfdt_final(ji) = ( 1.0 - zihgnew ) * ( zqfont_su(ji) + zqfont_bo(ji) ) 572 578 … … 582 588 ht_s_b(ji) = MAX( zzero , zhnfi ) 583 589 zqt_s(ji) = zqt_s(ji) * ht_s_b(ji) 590 ! we recompute dh_s_tot (clem) 591 dh_s_tot (ji) = ht_s_b(ji) - zhsold(ji) 584 592 585 593 ! Mass variations of ice and snow … … 598 606 !--------------------------------- 599 607 focea(ji) = - zfdt_final(ji) / rdt_ice ! focea is in W.m-2 * dt 600 ! salt flux 601 fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdvres * rhoic / rdt_ice 608 609 ! residual salt flux (clem) 610 !-------------------------- 611 ! surface 612 fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdvsur * rhoic / rdt_ice 613 ! bottom 614 IF ( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) >= 0._wp ) THEN ! melting 615 fseqv_1d(ji) = fseqv_1d(ji) - sm_i_b(ji) * a_i_b(ji) * zdvbot * rhoic / rdt_ice 616 ELSE ! growth 617 fseqv_1d(ji) = fseqv_1d(ji) - s_i_new(ji) * a_i_b(ji) * zdvbot * rhoic / rdt_ice 618 ENDIF 619 ! 602 620 ! ! diagnostic ( bottom ice growth ) 603 621 zji = MOD( npb(ji) - 1, jpi ) + 1 … … 700 718 !ENDIF 701 719 ! entrapment during snow ice formation 702 i_ice_switch = 1.0 - MAX( 0.e0 , SIGN( 1.0 , - ht_i_b(ji) + epsi13 ) )703 isnowic = 1.0 - MAX( 0.e0 , SIGN( 1.0 , - dh_snowice(ji) ) ) * i_ice_switch704 720 705 721 !clem IF( num_sal == 2 .OR. num_sal == 4 ) & … … 707 723 !clem & + sm_i_b(ji) * ht_i_b(ji) / MAX( ht_i_b(ji) + dh_snowice(ji), epsi13) & 708 724 !clem & - sm_i_b(ji) ) * isnowic 709 IF ( num_sal == 2 .OR. num_sal == 4 ) & 710 & dsm_i_si_1d(ji) = ( ( zsm_snowice * dh_snowice(ji) + sm_i_b(ji) * ht_i_b(ji) ) & 711 & / MAX( ht_i_b(ji) + dh_snowice(ji), epsi13 ) - sm_i_b(ji) ) * isnowic 725 726 ! clem: new salinity difference stored (to be used in limthd_ent.F90) 727 IF ( num_sal == 2 .OR. num_sal == 4 ) THEN 728 i_ice_switch = 1.0 - MAX( 0.e0 , SIGN( 1.0 , - zhgnew(ji) + epsi13 ) ) 729 ! salinity dif due to snow-ice formation 730 dsm_i_si_1d(ji) = ( zsm_snowice - sm_i_b(ji) ) * dh_snowice(ji) / MAX( zhgnew(ji), epsi13 ) * i_ice_switch 731 ! salinity dif due to bottom growth 732 IF ( ( fc_bo_i(ji) + fbif_1d(ji) + qlbbq_1d(ji) ) < 0._wp ) THEN 733 dsm_i_se_1d(ji) = ( s_i_new(ji) - sm_i_b(ji) ) * dh_i_bott(ji) / MAX( zhgnew(ji), epsi13 ) * i_ice_switch 734 ENDIF 735 ENDIF 712 736 713 737 ! Actualize new snow and ice thickness. -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90
r3938 r3963 81 81 INTEGER :: layer, nbpac ! local integers 82 82 INTEGER :: zji, zjj, iter ! - - 83 REAL(wp) :: ztmelts, zdv, zqold, zfrazb, zweight, zalphai, zindb, z de ! local scalars83 REAL(wp) :: ztmelts, zdv, zqold, zfrazb, zweight, zalphai, zindb, zinda, zde ! local scalars 84 84 REAL(wp) :: zgamafr, zvfrx, zvgx, ztaux, ztwogp, zf , zhicol_new ! - - 85 85 REAL(wp) :: ztenagm, zvfry, zvgy, ztauy, zvrel2, zfp, zsqcd , zhicrit ! - - … … 454 454 zji = MOD( npac(ji) - 1, jpi ) + 1 455 455 zjj = ( npac(ji) - 1 ) / jpi + 1 456 diag_lat_gr(zji,zjj) = zv_newice(ji) / rdt_ice456 diag_lat_gr(zji,zjj) = diag_lat_gr(zji,zjj) + zv_newice(ji) / rdt_ice ! clem 457 457 END DO !ji 458 458 … … 547 547 DO ji = 1, nbpac 548 548 zindb = MAX( 0._wp, SIGN( 1._wp , zdv_res(ji) ) ) 549 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + zindb * zdv_res(ji) * za_i_ac(ji,jl) / MAX( zat_i_lev(ji) , epsi06 ) 549 zinda = MAX( 0._wp, SIGN( 1._wp , zat_i_lev(ji) - epsi06 ) ) ! clem 550 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + zindb * zinda * zdv_res(ji) * za_i_ac(ji,jl) / MAX( zat_i_lev(ji) , epsi06 ) 550 551 END DO 551 552 END DO … … 623 624 ! Update salinity 624 625 !----------------- 625 IF( num_sal == 2 .OR. num_sal == 4 ) THEN626 !clem IF( num_sal == 2 .OR. num_sal == 4 ) THEN 626 627 DO jl = 1, jpl 627 628 DO ji = 1, nbpac 628 629 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zv_i_ac(ji,jl) + epsi10 ) ) ! 0 if no ice and 1 if yes 629 630 zdv = zv_i_ac(ji,jl) - zv_old(ji,jl) 630 zsmv_i_ac(ji,jl) = ( zsmv_i_ac(ji,jl) + zdv * zs_newice(ji) ) * zindb631 zsmv_i_ac(ji,jl) = zsmv_i_ac(ji,jl) + zdv * zs_newice(ji) * zindb ! clem modif 631 632 END DO 632 633 END DO 633 ENDIF634 !clem ENDIF 634 635 635 636 !-------------------------------- … … 640 641 zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - zv_i_ac(ji,jl) + epsi10 ) ) ! 0 if no ice and 1 if yes 641 642 zdv = zv_i_ac(ji,jl) - zv_old(ji,jl) 642 rdmicif_1d(ji) = rdmicif_1d(ji) + zdv * rhoic !* zindb643 rdmicif_1d(ji) = rdmicif_1d(ji) + zdv * rhoic * zindb 643 644 fseqv_1d(ji) = fseqv_1d(ji) - zdv * rhoic * zs_newice(ji) / rdt_ice * zindb 644 645 END DO … … 652 653 CALL tab_1d_2d( nbpac, v_i (:,:,jl), npac(1:nbpac), zv_i_ac (1:nbpac,jl), jpi, jpj ) 653 654 CALL tab_1d_2d( nbpac, oa_i(:,:,jl), npac(1:nbpac), zoa_i_ac(1:nbpac,jl), jpi, jpj ) 654 IF ( num_sal == 2 .OR. num_sal == 4 ) &655 !clem IF ( num_sal == 2 .OR. num_sal == 4 ) & 655 656 CALL tab_1d_2d( nbpac, smv_i (:,:,jl), npac(1:nbpac), zsmv_i_ac(1:nbpac,jl) , jpi, jpj ) 656 657 DO jk = 1, nlay_i -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90
r3938 r3963 127 127 i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 128 128 sm_i_b(ji) = i_ice_switch * sm_i_b(ji) + s_i_min * ( 1._wp - i_ice_switch ) 129 END DO ! ji 129 130 !---------------------------- 131 ! Heat flux - brine drainage 132 !---------------------------- 133 fhbri_1d(ji) = 0._wp 134 135 !---------------------------- 136 ! Salt flux - brine drainage 137 !---------------------------- 138 fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) * ( sm_i_b(ji) - zsiold(ji) ) / rdt_ice 139 IF( num_sal == 4 ) fsbri_1d(ji) = 0._wp 140 141 END DO ! ji 130 142 131 143 ! Salinity profile … … 136 148 !---------------------------- 137 149 138 DO ji = kideb, kiut150 !clem:move DO ji = kideb, kiut 139 151 !!gm useless 140 152 ! iflush : 1 if summer 141 iflush = MAX( 0._wp , SIGN ( 1._wp , t_su_b(ji) - rtt ) )153 !clem iflush = MAX( 0._wp , SIGN ( 1._wp , t_su_b(ji) - rtt ) ) 142 154 ! igravdr : 1 if t_su lt t_bo 143 igravdr = MAX( 0._wp , SIGN ( 1._wp , t_bo_b(ji) - t_su_b(ji) ) )155 !clem igravdr = MAX( 0._wp , SIGN ( 1._wp , t_bo_b(ji) - t_su_b(ji) ) ) 144 156 ! iaccrbo : 1 if bottom accretion 145 iaccrbo = MAX( 0._wp , SIGN ( 1._wp , dh_i_bott(ji) ) )157 !clem iaccrbo = MAX( 0._wp , SIGN ( 1._wp , dh_i_bott(ji) ) ) 146 158 !!gm end useless 147 159 ! 148 fhbri_1d(ji) = 0._wp149 END DO ! ji160 !clem:move fhbri_1d(ji) = 0._wp 161 !clem:move END DO ! ji 150 162 151 163 !---------------------------- 152 164 ! Salt flux - brine drainage 153 165 !---------------------------- 154 DO ji = kideb, kiut155 i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) )156 fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) * ( sm_i_b(ji) - zsiold(ji) ) / rdt_ice166 !clem:move DO ji = kideb, kiut 167 !clem:move i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 168 !clem:move fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) * ( sm_i_b(ji) - zsiold(ji) ) / rdt_ice 157 169 !i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - zhiold(ji) ) ) 158 170 !fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * zhiold(ji) * ( sm_i_b(ji) - zsiold(ji) ) / rdt_ice 159 171 !clem fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) & 160 172 !clem & * ( MAX(dsm_i_gd_1d(ji) + dsm_i_fl_1d(ji), sm_i_b(ji) - zsiold(ji) ) ) / rdt_ice 161 IF( num_sal == 4 ) fsbri_1d(ji) = 0._wp162 END DO ! ji173 !clem:move IF( num_sal == 4 ) fsbri_1d(ji) = 0._wp 174 !clem:move END DO ! ji 163 175 164 176 ! Only necessary for conservation check since salinity is modified -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90
r3938 r3963 81 81 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zs0e 82 82 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 83 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 83 84 ! mass and salt flux (clem) 84 85 REAL(wp), POINTER, DIMENSION(:,:,:) :: zviold ! old ice volume... … … 170 171 zusnit = 1.0 / REAL( initad ) 171 172 IF( zcfl > 0.5 .AND. lwp ) & 172 WRITE(numout,*) 'lim_trp _2: CFL violation at day ', nday, ', cfl = ', zcfl, &173 WRITE(numout,*) 'lim_trp : CFL violation at day ', nday, ', cfl = ', zcfl, & 173 174 & ': the ice time stepping is split in two' 174 175 … … 430 431 DO ji = 1, jpi 431 432 432 IF ( v_i(ji,jj,jl) > 0. ) THEN433 IF ( v_i(ji,jj,jl) > 0._wp ) THEN 433 434 zvi = v_i(ji,jj,jl) 434 435 zvs = v_s(ji,jj,jl) … … 437 438 438 439 IF ( ( zdv > 0.0 .AND. ht_i(ji,jj,jl) > zhimax(ji,jj,jl) .AND. SUM( zaiold(ji,jj,1:jpl) ) < 0.80 ) .OR. & 439 440 & ( zdv < 0.0 .AND. ht_i(ji,jj,jl) > zhimax(ji,jj,jl) ) ) THEN 440 441 ht_i(ji,jj,jl) = MIN( zhimax(ji,jj,jl), hi_max(jl) ) 441 442 zindh = MAX( rzero, SIGN( rone, ht_i(ji,jj,jl) - epsi10 ) ) … … 447 448 ENDIF 448 449 449 ! zindh = MAX( rzero, SIGN( rone, ht_i(ji,jj,jl) - epsi10 ) )450 v_i(ji,jj,jl) = a_i(ji,jj,jl) * ht_i(ji,jj,jl)451 v_s(ji,jj,jl) = a_i(ji,jj,jl) * ht_s(ji,jj,jl)452 453 ! Update mass fluxes (clem)450 ! small correction due to *zindh for a_i 451 v_i(ji,jj,jl) = zindh * v_i(ji,jj,jl) 452 v_s(ji,jj,jl) = zindh * v_s(ji,jj,jl) 453 454 ! Update mass fluxes 454 455 rdmicif(ji,jj) = rdmicif(ji,jj) + ( v_i(ji,jj,jl) - zvi ) * rhoic 455 456 rdmsnif(ji,jj) = rdmsnif(ji,jj) + ( v_s(ji,jj,jl) - zvs ) * rhosn … … 587 588 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 588 589 590 zchk_vmin = glob_min(v_i) 591 zchk_amax = glob_max(SUM(a_i,dim=3)) 592 zchk_amin = glob_min(a_i) 593 589 594 IF(lwp) THEN 590 IF ( 591 IF ( 592 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limtrp) = ',(MINVAL(v_i)* 1.e-3)593 IF ( MINVAL( a_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation a_i<0 (limtrp) = ',MINVAL(a_i)595 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limtrp) = ',(zchk_v_i * 86400.) 596 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limtrp) = ',(zchk_smv * 86400.) 597 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limtrp) = ',(zchk_vmin * 1.e-3) 598 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limtrp) = ',zchk_amin 594 599 ENDIF 595 600 ENDIF -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90
r3938 r3963 86 86 INTEGER :: i_ice_switch 87 87 INTEGER :: ind_im, layer ! indices for internal melt 88 REAL(wp) :: zweight, zesum, z_da_i 89 REAL(wp) :: zind b, zindsn, zindic88 REAL(wp) :: zweight, zesum, z_da_i, zhimax 89 REAL(wp) :: zinda, zindb, zindsn, zindic 90 90 REAL(wp) :: zindg, zh, zdvres, zviold2 91 91 REAL(wp) :: zbigvalue, zvsold2, z_da_ex … … 94 94 REAL(wp), POINTER, DIMENSION(:) :: zthick0, zqm0 ! thickness of the layers and heat contents for 95 95 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 96 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 96 97 ! mass and salt flux (clem) 97 98 REAL(wp), POINTER, DIMENSION(:,:,:) :: zviold, zvsold, zsmvold ! old ice volume... … … 182 183 ! 2. Review of all pathological cases 183 184 !-------------------------------------- 185 186 !------------------------------------------- 187 ! 2.1) Advection of ice in an ice-free cell 188 !------------------------------------------- 189 ! should be removed since it is treated after dynamics now 190 191 ! !IF( ln_nicep ) THEN 192 ! WRITE(numout,*) ' limupdate1 - before h correction ' 193 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 194 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 195 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 196 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 197 ! !ENDIF 198 ! 199 zhimax = 1._wp 200 ! first category 201 DO jj = 1, jpj 202 DO ji = 1, jpi 203 !--- the thickness of such an ice is often out of bounds 204 !--- thus we recompute a new area while conserving ice volume 205 zat_i_old = SUM( old_a_i(ji,jj,:) ) 206 zindb = MAX( 0._wp, SIGN( 1._wp, ABS( d_a_i_trp(ji,jj,1)) - epsi10 ) ) 207 IF ( ( ABS(d_v_i_trp(ji,jj,1))/MAX(ABS(d_a_i_trp(ji,jj,1)),epsi10)*zindb .GT. zhimax) & 208 .AND.( ( v_i(ji,jj,1)/MAX(a_i(ji,jj,1),epsi10)*zindb) .GT. zhimax ) & 209 .AND.( zat_i_old .LT. 1.e-6 ) ) THEN ! new line 210 ht_i(ji,jj,1) = hi_max(1) / 2.0 211 a_i (ji,jj,1) = v_i(ji,jj,1) / ht_i(ji,jj,1) 212 ENDIF 213 END DO 214 END DO 215 216 217 ! !IF( ln_nicep ) THEN 218 ! at_i(:,:) = 0._wp 219 ! DO jl = 1, jpl 220 ! at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 221 ! END DO 222 ! WRITE(numout,*) ' limupdate1 - after h correction 1 ' 223 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 224 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 225 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 226 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 227 ! !ENDIF 228 229 zhimax = 10._wp 230 ! other categories 231 DO jl = 2, jpl 232 jm = ice_types(jl) 233 DO jj = 1, jpj 234 DO ji = 1, jpi 235 zindb = MAX( rzero, SIGN( rone, ABS(d_a_i_trp(ji,jj,jl)) - epsi10 ) ) 236 ! this correction is very tricky... sometimes, advection gets wrong i don't know why 237 ! it makes problems when the advected volume and concentration do not seem to be 238 ! related with each other 239 ! the new thickness is sometimes very big! 240 ! and sometimes d_a_i_trp and d_v_i_trp have different sign 241 ! which of course is plausible 242 ! but fuck! it fucks everything up :) 243 IF ( (ABS(d_v_i_trp(ji,jj,jl))/MAX(ABS(d_a_i_trp(ji,jj,jl)),epsi10)*zindb .GT. zhimax) & 244 .AND.(v_i(ji,jj,jl)/MAX(a_i(ji,jj,jl),epsi10)*zindb) .GT. zhimax ) THEN 245 ht_i(ji,jj,jl) = ( hi_max_typ(jl-ice_cat_bounds(jm,1),jm) + hi_max_typ(jl-ice_cat_bounds(jm,1)+1,jm) ) / 2.0 246 a_i (ji,jj,jl) = v_i(ji,jj,jl) / ht_i(ji,jj,jl) 247 ENDIF 248 END DO ! ji 249 END DO !jj 250 END DO !jl 184 251 185 252 at_i(:,:) = 0._wp … … 187 254 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 188 255 END DO 256 257 ! !IF( ln_nicep ) THEN 258 ! WRITE(numout,*) ' limupdate1 - after h correction 2 ' 259 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 260 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 261 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 262 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 263 ! !ENDIF 189 264 190 265 !---------------------------------------------------- … … 247 322 !---------------------------------------------------------------------------- 248 323 zindg = tms(ji,jj) * MAX( 0._wp, SIGN( 1._wp, -v_i(ji,jj,jl) ) ) 249 v_i(ji,jj,jl) = v_i(ji,jj,jl) + zindg * rhosn * v_s(ji,jj,jl) / rau0 250 v_s(ji,jj,jl) = ( 1._wp - zindg ) * v_s(ji,jj,jl) + zindg * v_i(ji,jj,jl) * ( rau0 - rhoic ) / rhosn 324 zdvres = zindg * rhosn * v_s(ji,jj,jl) / rau0 325 v_i(ji,jj,jl) = v_i(ji,jj,jl) + zdvres 326 327 zdvres = zindsn*zindb * ( - zindg * v_s(ji,jj,jl) + zindg * v_i(ji,jj,jl) * ( rau0 - rhoic ) / rhosn ) 328 v_s(ji,jj,jl) = v_s(ji,jj,jl) + zdvres 251 329 252 330 !--- 2.7 Correction to ice concentrations … … 286 364 287 365 !--- 2.13 ice concentration should not exceed amax 366 ! (it should not be the case) 288 367 !----------------------------------------------------- 289 368 DO jj = 1, jpj 290 369 DO ji = 1, jpi 291 370 z_da_ex = MAX( at_i(ji,jj) - amax , 0.0 ) 371 zindb = MAX( rzero, SIGN( rone, at_i(ji,jj) - epsi06 ) ) 292 372 DO jl = 1, jpl 293 zindb = MAX( rzero, SIGN( rone, v_i(ji,jj,jl) - epsi04 ) )294 373 z_da_i = a_i(ji,jj,jl) * z_da_ex / MAX( at_i(ji,jj), epsi06 ) * zindb 295 a_i(ji,jj,jl) = a_i(ji,jj,jl) - z_da_i 296 END DO 297 END DO 298 END DO 374 a_i(ji,jj,jl) = MAX( 0._wp, a_i(ji,jj,jl) - z_da_i ) 375 ! 376 zinda = MAX( rzero, SIGN( rone, a_i(ji,jj,jl) - epsi06 ) ) 377 ht_i(ji,jj,jl) = v_i(ji,jj,jl) / MAX( a_i(ji,jj,jl), epsi06 ) * zinda 378 !v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl) ! makes ice shrinken but should not be used 379 END DO 380 END DO 381 END DO 299 382 at_i(:,:) = a_i(:,:,1) 300 383 DO jl = 2, jpl 301 384 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 302 385 END DO 386 303 387 304 388 ! Final thickness distribution rebinning … … 311 395 ENDIF 312 396 END DO 397 313 398 314 399 !--------------------- … … 356 441 !- check conservation (C Rousset) 357 442 IF (ln_limdiahsb) THEN 443 358 444 zchk_fs = glob_sum( ( fsbri(:,:) + fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:) ) * area(:,:) * tms(:,:) ) - zchk_fs_b 359 445 zchk_fw = glob_sum( rdmicif(:,:) * area(:,:) * tms(:,:) ) - zchk_fw_b … … 362 448 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 363 449 450 zchk_vmin = glob_min(v_i) 451 zchk_amax = glob_max(SUM(a_i,dim=3)) 452 zchk_amin = glob_min(a_i) 453 364 454 IF(lwp) THEN 365 IF ( 366 IF ( 367 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate1) = ',(MINVAL(v_i)* 1.e-3)368 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate1) = ',MAXVAL(SUM(a_i,dim=3))369 IF ( MINVAL( a_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate1) = ',MINVAL(a_i)455 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limupdate1) = ',(zchk_v_i * 86400.) 456 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limupdate1) = ',(zchk_smv * 86400.) 457 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate1) = ',(zchk_vmin * 1.e-3) 458 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate1) = ',zchk_amax 459 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate1) = ',zchk_amin 370 460 ENDIF 371 461 ENDIF -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90
r3938 r3963 83 83 INTEGER :: ind_im, layer ! indices for internal melt 84 84 REAL(wp) :: zweight, zesum, zhimax, z_da_i 85 REAL(wp) :: zind b, zindsn, zindic85 REAL(wp) :: zinda, zindb, zindsn, zindic 86 86 REAL(wp) :: zindg, zh, zdvres, zviold2 87 87 REAL(wp) :: zbigvalue, zvsold2, z_da_ex … … 91 91 REAL(wp), POINTER, DIMENSION(:) :: zthick0, zqm0 ! thickness of the layers and heat contents for 92 92 REAL(wp) :: zchk_v_i, zchk_smv, zchk_fs, zchk_fw, zchk_v_i_b, zchk_smv_b, zchk_fs_b, zchk_fw_b ! Check conservation (C Rousset) 93 93 REAL(wp) :: zchk_vmin, zchk_amin, zchk_amax ! Check errors (C Rousset) 94 94 ! mass and salt flux (clem) 95 95 REAL(wp), POINTER, DIMENSION(:,:,:) :: zviold, zvsold, zsmvold ! old ice volume... … … 177 177 END DO 178 178 179 180 179 !-------------------------------------- 181 180 ! 2. Review of all pathological cases 182 181 !-------------------------------------- 183 182 !------------------------------------------- 183 ! 2.1) Advection of ice in an ice-free cell 184 !------------------------------------------- 185 ! should be removed since it is treated after dynamics now 186 187 ! !IF( ln_nicep ) THEN 188 ! WRITE(numout,*) ' limupdate2 - before h correction ' 189 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 190 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 191 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 192 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 193 ! !ENDIF 194 195 zhimax = 1._wp 196 ! first category 197 DO jj = 1, jpj 198 DO ji = 1, jpi 199 !--- the thickness of such an ice is often out of bounds 200 !--- thus we recompute a new area while conserving ice volume 201 zat_i_old = SUM( old_a_i(ji,jj,:) ) 202 zindb = MAX( 0._wp, SIGN( 1._wp, ABS( d_a_i_thd(ji,jj,1)) - epsi10 ) ) 203 IF ( ( ABS(d_v_i_thd(ji,jj,1))/MAX(ABS(d_a_i_thd(ji,jj,1)),epsi10)*zindb .GT. zhimax) & 204 .AND.( ( v_i(ji,jj,1)/MAX(a_i(ji,jj,1),epsi10)*zindb) .GT. zhimax ) & 205 .AND.( zat_i_old .LT. 1.e-6 ) ) THEN ! new line 206 ht_i(ji,jj,1) = hi_max(1) / 2.0 207 a_i (ji,jj,1) = v_i(ji,jj,1) / ht_i(ji,jj,1) 208 ENDIF 209 END DO 210 END DO 211 212 ! !IF( ln_nicep ) THEN 213 ! at_i(:,:) = 0._wp 214 ! DO jl = 1, jpl 215 ! at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 216 ! END DO 217 ! WRITE(numout,*) ' limupdate2 - after h correction 1 ' 218 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 219 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 220 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 221 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 222 ! !ENDIF 223 ! 224 zhimax = 10._wp 225 ! other categories 226 DO jl = 2, jpl 227 jm = ice_types(jl) 228 DO jj = 1, jpj 229 DO ji = 1, jpi 230 zindb = MAX( rzero, SIGN( rone, ABS(d_a_i_thd(ji,jj,jl)) - epsi10 ) ) 231 ! this correction is very tricky... sometimes, advection gets wrong i don't know why 232 ! it makes problems when the advected volume and concentration do not seem to be 233 ! related with each other 234 ! the new thickness is sometimes very big! 235 ! and sometimes d_a_i_trp and d_v_i_trp have different sign 236 ! which of course is plausible 237 ! but fuck! it fucks everything up :) 238 IF ( (ABS(d_v_i_thd(ji,jj,jl))/MAX(ABS(d_a_i_thd(ji,jj,jl)),epsi10)*zindb .GT. zhimax) & 239 .AND.(v_i(ji,jj,jl)/MAX(a_i(ji,jj,jl),epsi10)*zindb) .GT. zhimax ) THEN 240 ht_i(ji,jj,jl) = ( hi_max_typ(jl-ice_cat_bounds(jm,1),jm) + hi_max_typ(jl-ice_cat_bounds(jm,1)+1,jm) ) / 2.0 241 a_i (ji,jj,jl) = v_i(ji,jj,jl) / ht_i(ji,jj,jl) 242 ENDIF 243 END DO ! ji 244 END DO !jj 245 END DO !jl 246 247 at_i(:,:) = 0._wp 248 DO jl = 1, jpl 249 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 250 END DO 251 252 ! !IF( ln_nicep ) THEN 253 ! WRITE(numout,*) ' limupdate2 - after h correction 2 ' 254 ! WRITE(numout,*) ' a_i : ', a_i(12, 44, 1:jpl) 255 ! WRITE(numout,*) ' at_i : ', at_i(12,44) 256 ! WRITE(numout,*) ' v_i : ', v_i(12, 44, 1:jpl) 257 ! WRITE(numout,*) ' ht_i : ', ht_i(12, 44, 1:jpl) 258 ! !ENDIF 184 259 185 260 !---------------------------------------------------- … … 363 438 v_i(ji,jj,jl) = v_i(ji,jj,jl) + zdvres 364 439 365 zdvres = - zindg * v_s(ji,jj,jl) + zindg * v_i(ji,jj,jl) * ( rau0 - rhoic ) / rhosn440 zdvres = zindsn*zindb * ( - zindg * v_s(ji,jj,jl) + zindg * v_i(ji,jj,jl) * ( rau0 - rhoic ) / rhosn ) 366 441 v_s(ji,jj,jl) = v_s(ji,jj,jl) + zdvres 367 442 … … 421 496 422 497 !--- 2.13 ice concentration should not exceed amax 498 ! (it should not be the case) 423 499 !----------------------------------------------------- 424 500 DO jj = 1, jpj 425 501 DO ji = 1, jpi 426 502 z_da_ex = MAX( at_i(ji,jj) - amax , 0.0 ) 503 zindb = MAX( rzero, SIGN( rone, at_i(ji,jj) - epsi06 ) ) 427 504 DO jl = 1, jpl 428 zindb = MAX( rzero, SIGN( rone, v_i(ji,jj,jl) - epsi04 ) )429 505 z_da_i = a_i(ji,jj,jl) * z_da_ex / MAX( at_i(ji,jj), epsi06 ) * zindb 430 a_i(ji,jj,jl) = a_i(ji,jj,jl) - z_da_i 506 a_i(ji,jj,jl) = MAX( 0._wp, a_i(ji,jj,jl) - z_da_i ) 507 ! 508 zinda = MAX( rzero, SIGN( rone, a_i(ji,jj,jl) - epsi06 ) ) 509 ht_i(ji,jj,jl) = v_i(ji,jj,jl) / MAX( a_i(ji,jj,jl), epsi06 ) * zinda 510 !v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl) ! makes ice shrinken but should not be used 431 511 END DO 432 512 END DO … … 517 597 zchk_smv = ( glob_sum( SUM( smv_i(:,:,:), dim=3 ) * area(:,:) * tms(:,:) ) - zchk_smv_b ) / rdt_ice + ( zchk_fs / rhoic ) 518 598 599 zchk_vmin = glob_min(v_i) 600 zchk_amax = glob_max(SUM(a_i,dim=3)) 601 zchk_amin = glob_min(a_i) 602 519 603 IF(lwp) THEN 520 IF ( 521 IF ( 522 IF ( MINVAL( v_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate2) = ',(MINVAL(v_i)* 1.e-3)523 IF ( MAXVAL( SUM(a_i(:,:,:),dim=3) ) > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate2) = ',MAXVAL(SUM(a_i,dim=3))524 IF ( MINVAL( a_i(:,:,:) ) < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate2) = ',MINVAL(a_i)604 IF ( ABS( zchk_v_i ) > 1.e-5 ) WRITE(numout,*) 'violation volume [m3/day] (limupdate2) = ',(zchk_v_i * 86400.) 605 IF ( ABS( zchk_smv ) > 1.e-4 ) WRITE(numout,*) 'violation saline [psu*m3/day] (limupdate2) = ',(zchk_smv * 86400.) 606 IF ( zchk_vmin < 0. ) WRITE(numout,*) 'violation v_i<0 [mm] (limupdate2) = ',(zchk_vmin * 1.e-3) 607 IF ( zchk_amax > amax+epsi10 ) WRITE(numout,*) 'violation a_i>amax (limupdate2) = ',zchk_amax 608 IF ( zchk_amin < 0. ) WRITE(numout,*) 'violation a_i<0 (limupdate2) = ',zchk_amin 525 609 ENDIF 526 610 ENDIF -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
r3938 r3963 180 180 ht_s(ji,jj,jl) = v_s (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi10 ) * zindb 181 181 o_i(ji,jj,jl) = oa_i(ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi10 ) * zindb 182 a_i(ji,jj,jl) = a_i (ji,jj,jl) * zindb ! clem correction 182 183 END DO 183 184 END DO -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90
r3938 r3963 122 122 ! Normal file 123 123 !------------- 124 125 zsto = rdt_ice126 IF( ln_mskland ) THEN ; clop = "ave(only(x))" ! put 1.e+20 on land (very expensive!!)127 ELSE ; clop = "ave(x)" ! no use of the mask value (require less cpu time)128 ENDIF129 zout = nwrite * rdt_ice / nn_fsbc130 124 niter = ( nit000 - 1 ) / nn_fsbc 131 zdept(1) = 0.132 133 125 CALL ymds2ju ( nyear, nmonth, nday, rdt, zjulian ) 134 126 zjulian = zjulian - adatrj ! set calendar origin to the beginning of the experiment 135 CALL dia_nam ( clhstnam, nwrite, 'icemod_old' ) 136 CALL histbeg ( clhstnam, jpi, glamt, jpj, gphit, 1, jpi, 1, jpj, niter, zjulian, rdt_ice, & 137 & nhorid, nice, domain_id=nidom, snc4chunks=snc4set ) 138 CALL histvert( nice, "deptht", "Vertical T levels", "m", 1, zdept, ndepid, "down") 139 CALL wheneq ( jpij , tmask(:,:,1), 1, 1., ndex51, ndim) 140 141 DO jf = 1 , noumef 142 IF(lwp) WRITE(numout,*) 'jf', jf 143 IF ( nc(jf) == 1 ) THEN 144 CALL histdef( nice, nam(jf), titn(jf), uni(jf), jpi, jpj & 145 , nhorid, 1, 1, 1, -99, 32, clop, zsto, zout ) 146 IF(lwp) WRITE(numout,*) 'nice, nam(jf), titn(jf), uni(jf), nhorid, clop, zsto, zout' 147 IF(lwp) WRITE(numout,*) nice, nam(jf), titn(jf), uni(jf), nhorid, clop, zsto, zout 148 ENDIF 149 END DO 150 151 CALL histend(nice, snc4set) 152 127 !clem 128 ! zsto = rdt_ice 129 ! IF( ln_mskland ) THEN ; clop = "ave(only(x))" ! put 1.e+20 on land (very expensive!!) 130 ! ELSE ; clop = "ave(x)" ! no use of the mask value (require less cpu time) 131 ! ENDIF 132 ! zout = nwrite * rdt_ice / nn_fsbc 133 ! zdept(1) = 0. 134 ! 135 ! CALL dia_nam ( clhstnam, nwrite, 'icemod_old' ) 136 ! CALL histbeg ( clhstnam, jpi, glamt, jpj, gphit, 1, jpi, 1, jpj, niter, zjulian, rdt_ice, & 137 ! & nhorid, nice, domain_id=nidom, snc4chunks=snc4set ) 138 ! CALL histvert( nice, "deptht", "Vertical T levels", "m", 1, zdept, ndepid, "down") 139 ! CALL wheneq ( jpij , tmask(:,:,1), 1, 1., ndex51, ndim) 140 ! 141 ! DO jf = 1 , noumef 142 ! IF(lwp) WRITE(numout,*) 'jf', jf 143 ! IF ( nc(jf) == 1 ) THEN 144 ! CALL histdef( nice, nam(jf), titn(jf), uni(jf), jpi, jpj & 145 ! , nhorid, 1, 1, 1, -99, 32, clop, zsto, zout ) 146 ! IF(lwp) WRITE(numout,*) 'nice, nam(jf), titn(jf), uni(jf), nhorid, clop, zsto, zout' 147 ! IF(lwp) WRITE(numout,*) nice, nam(jf), titn(jf), uni(jf), nhorid, clop, zsto, zout 148 ! ENDIF 149 ! END DO 150 ! 151 ! CALL histend(nice, snc4set) 152 !clem 153 ! 153 154 !----------------- 154 155 ! ITD file output … … 280 281 ! ecriture d'un fichier netcdf 281 282 ! 283 282 284 niter = niter + 1 283 DO jf = 1 , noumef 284 ! 285 zfield(:,:) = zcmo(:,:,jf) * cmulti(jf) + cadd(jf) 286 ! 287 IF( jf == 7 .OR. jf == 8 .OR. jf == 15 .OR. jf == 16 ) THEN ; CALL lbc_lnk( zfield, 'T', -1. ) 288 ELSE ; CALL lbc_lnk( zfield, 'T', 1. ) 289 ENDIF 290 ! 291 IF( ln_nicep ) THEN 292 WRITE(numout,*) 293 WRITE(numout,*) 'nc(jf), nice, nam(jf), niter, ndim' 294 WRITE(numout,*) nc(jf), nice, nam(jf), niter, ndim 295 ENDIF 296 IF( nc(jf) == 1 ) CALL histwrite( nice, nam(jf), niter, zfield, ndim, ndex51 ) 297 ! 298 END DO 299 300 IF( ( nn_fsbc * niter ) >= nitend .OR. kindic < 0 ) THEN 301 IF( lwp) WRITE(numout,*) ' Closing the icemod file ' 302 CALL histclo( nice ) 303 ENDIF 304 285 !clem 286 ! DO jf = 1 , noumef 287 ! ! 288 ! zfield(:,:) = zcmo(:,:,jf) * cmulti(jf) + cadd(jf) 289 ! ! 290 ! IF( jf == 7 .OR. jf == 8 .OR. jf == 15 .OR. jf == 16 ) THEN ; CALL lbc_lnk( zfield, 'T', -1. ) 291 ! ELSE ; CALL lbc_lnk( zfield, 'T', 1. ) 292 ! ENDIF 293 ! ! 294 ! IF( ln_nicep ) THEN 295 ! WRITE(numout,*) 296 ! WRITE(numout,*) 'nc(jf), nice, nam(jf), niter, ndim' 297 ! WRITE(numout,*) nc(jf), nice, nam(jf), niter, ndim 298 ! ENDIF 299 ! IF( nc(jf) == 1 ) CALL histwrite( nice, nam(jf), niter, zfield, ndim, ndex51 ) 300 ! ! 301 ! END DO 302 ! 303 ! IF( ( nn_fsbc * niter ) >= nitend .OR. kindic < 0 ) THEN 304 ! IF( lwp) WRITE(numout,*) ' Closing the icemod file ' 305 ! CALL histclo( nice ) 306 ! ENDIF 307 !clem 308 ! 305 309 CALL iom_put ('iceconc', zcmo(:,:,1) ) ! field1: ice concentration 306 310 CALL iom_put ('icethic_cea', zcmo(:,:,2) ) ! field2: ice thickness (i.e. icethi(:,:)) … … 328 332 CALL iom_put ('isnowhc', zcmo(:,:,26) ) ! field 26: snow total heat content 329 333 CALL iom_put ('icest', zcmo(:,:,27) ) ! field 27: ice surface temperature 330 CALL iom_put ('fsbri', zcmo(:,:,28) * 86400 ) ! field 28: brine salt flux331 CALL iom_put ('fseqv', zcmo(:,:,29) * 86400 ) ! field 29: equivalent FW salt flux334 CALL iom_put ('fsbri', zcmo(:,:,28) * 86400. ) ! field 28: brine salt flux 335 CALL iom_put ('fseqv', zcmo(:,:,29) * 86400. ) ! field 29: equivalent FW salt flux 332 336 CALL iom_put ('ibrinv', zcmo(:,:,30) *100 ) ! field 30: brine volume 333 337 CALL iom_put ('icecolf', zcmo(:,:,31) ) ! field 31: frazil ice collection thickness … … 341 345 CALL iom_put ('icevolu', zcmo(:,:,39) ) ! field 39: ice volume 342 346 CALL iom_put ('snowvol', zcmo(:,:,40) ) ! field 40: snow volume 343 CALL iom_put ('fsrpo', zcmo(:,:,41) * 86400 ) ! field 41: salt flux from ridging rafting344 CALL iom_put ('fsres', zcmo(:,:,42) * 86400 ) ! field 42: salt flux from limupdate (resultant)347 CALL iom_put ('fsrpo', zcmo(:,:,41) * 86400. ) ! field 41: salt flux from ridging rafting 348 CALL iom_put ('fsres', zcmo(:,:,42) * 86400. ) ! field 42: salt flux from limupdate (resultant) 345 349 CALL iom_put ('icetrp', zcmo(:,:,43) ) ! field 43: ice volume transport 346 350 -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limwri_dimg.h90
r3938 r3963 14 14 !! modif : 03/06/98 15 15 !!------------------------------------------------------------------- 16 USE lib_fortran ! to use key_nosignedzero 16 17 USE diawri, ONLY : dia_wri_dimg 17 18 REAL(wp),DIMENSION(1) :: zdept -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r3938 r3963 41 41 USE sbc_ice ! Surface boundary condition: ice fields 42 42 #endif 43 USE lib_fortran ! to use key_nosignedzero 43 44 44 45 IMPLICIT NONE … … 530 531 DO ji = fs_2, fs_jpim1 ! vect. opt. 531 532 p_taui(ji,jj) = zcoef_wnorm2 * ( z_wnds_t(ji+1,jj ) + z_wnds_t(ji,jj) ) & 532 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - pui(ji,jj) )533 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * pui(ji,jj) ) 533 534 p_tauj(ji,jj) = zcoef_wnorm2 * ( z_wnds_t(ji,jj+1 ) + z_wnds_t(ji,jj) ) & 534 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - pvi(ji,jj) )535 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * pvi(ji,jj) ) 535 536 END DO 536 537 END DO … … 560 561 ! iovino 561 562 IF( ff(ji,jj) .GT. 0._wp ) THEN 562 z_qlw(ji,jj,jl) = ( 0.95 * sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) ! MV test 2 garde juste cette ligne ci563 z_qlw(ji,jj,jl) = ( 0.95 * sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) 563 564 ELSE 564 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) ! MV test 1 garde juste cette ligne ci565 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) 565 566 ENDIF 566 567 ! lw sensitivity … … 610 611 !CDIR COLLAPSE 611 612 p_spr(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! solid precipitation [kg/m2/s] 612 CALL iom_put( 'snowpre', p_spr * 86400 ) ! Snow precipitation613 CALL iom_put( 'precip', p_tpr * 86400 ) ! Total precipitation613 CALL iom_put( 'snowpre', p_spr * 86400. ) ! Snow precipitation 614 CALL iom_put( 'precip', p_tpr * 86400. ) ! Total precipitation 614 615 ! 615 616 IF(ln_ctl) THEN -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r3938 r3963 97 97 98 98 ! ! overwrite namelist parameter using CPP key information 99 !IF( Agrif_Root() ) THEN ! AGRIF zoom100 !IF( lk_lim2 ) nn_ice = 2101 !IF( lk_lim3 ) nn_ice = 3102 !IF( lk_cice ) nn_ice = 4103 !ENDIF99 IF( Agrif_Root() ) THEN ! AGRIF zoom 100 IF( lk_lim2 ) nn_ice = 2 101 IF( lk_lim3 ) nn_ice = 3 102 IF( lk_cice ) nn_ice = 4 103 ENDIF 104 104 IF( cp_cfg == 'gyre' ) THEN ! GYRE configuration 105 105 ln_ana = .TRUE. -
branches/2013/dev_r3406_CNRS_LIM3/NEMOGCM/NEMO/OPA_SRC/lib_fortran.F90
r3294 r3963 5 5 !!====================================================================== 6 6 !! History : 3.2 ! 2010-05 (M. Dunphy, R. Benshila) Original code 7 !! 3.4 ! 2013-06 (C. Rousset) add glob_min, glob_max 8 !! + 3d dim. of input is fexible (jpk, jpl...) 7 9 !!---------------------------------------------------------------------- 8 10 … … 23 25 24 26 PUBLIC glob_sum 27 PUBLIC glob_min, glob_max 25 28 #if defined key_nosignedzero 26 29 PUBLIC SIGN … … 29 32 INTERFACE glob_sum 30 33 MODULE PROCEDURE glob_sum_2d, glob_sum_3d,glob_sum_2d_a, glob_sum_3d_a 34 END INTERFACE 35 INTERFACE glob_min 36 MODULE PROCEDURE glob_min_2d, glob_min_3d,glob_min_2d_a, glob_min_3d_a 37 END INTERFACE 38 INTERFACE glob_max 39 MODULE PROCEDURE glob_max_2d, glob_max_3d,glob_max_2d_a, glob_max_3d_a 31 40 END INTERFACE 32 41 … … 47 56 48 57 #if ! defined key_mpp_rep 58 59 ! --- SUM --- 49 60 FUNCTION glob_sum_2d( ptab ) 50 61 !!----------------------------------------------------------------------- … … 61 72 ! 62 73 END FUNCTION glob_sum_2d 63 64 74 65 75 FUNCTION glob_sum_3d( ptab ) 66 76 !!----------------------------------------------------------------------- … … 73 83 !! 74 84 INTEGER :: jk 75 !!----------------------------------------------------------------------- 85 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 86 !!----------------------------------------------------------------------- 87 ! 88 zjpk = SIZE(ptab,3) 76 89 ! 77 90 glob_sum_3d = 0.e0 78 DO jk = 1, jpk91 DO jk = 1, zjpk 79 92 glob_sum_3d = glob_sum_3d + SUM( ptab(:,:,jk)*tmask_i(:,:) ) 80 93 END DO … … 111 124 !! 112 125 INTEGER :: jk 113 !!----------------------------------------------------------------------- 126 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 127 !!----------------------------------------------------------------------- 128 ! 129 zjpk = SIZE(ptab1,3) 114 130 ! 115 131 glob_sum_3d_a(:) = 0.e0 116 DO jk = 1, jpk132 DO jk = 1, zjpk 117 133 glob_sum_3d_a(1) = glob_sum_3d_a(1) + SUM( ptab1(:,:,jk)*tmask_i(:,:) ) 118 134 glob_sum_3d_a(2) = glob_sum_3d_a(2) + SUM( ptab2(:,:,jk)*tmask_i(:,:) ) … … 121 137 ! 122 138 END FUNCTION glob_sum_3d_a 139 140 141 ! --- MIN --- 142 FUNCTION glob_min_2d( ptab ) 143 !!----------------------------------------------------------------------- 144 !! *** FUNCTION glob_min_2D *** 145 !! 146 !! ** Purpose : perform a masked min on the inner global domain of a 2D array 147 !!----------------------------------------------------------------------- 148 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab ! input 2D array 149 REAL(wp) :: glob_min_2d ! global masked min 150 !!----------------------------------------------------------------------- 151 ! 152 glob_min_2d = MINVAL( ptab(:,:)*tmask_i(:,:) ) 153 IF( lk_mpp ) CALL mpp_min( glob_min_2d ) 154 ! 155 END FUNCTION glob_min_2d 156 157 FUNCTION glob_min_3d( ptab ) 158 !!----------------------------------------------------------------------- 159 !! *** FUNCTION glob_min_3D *** 160 !! 161 !! ** Purpose : perform a masked min on the inner global domain of a 3D array 162 !!----------------------------------------------------------------------- 163 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab ! input 3D array 164 REAL(wp) :: glob_min_3d ! global masked min 165 !! 166 INTEGER :: jk 167 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 168 !!----------------------------------------------------------------------- 169 ! 170 zjpk = SIZE(ptab,3) 171 ! 172 glob_min_3d = 0.e0 173 DO jk = 1, zjpk 174 glob_min_3d = glob_min_3d + MINVAL( ptab(:,:,jk)*tmask_i(:,:) ) 175 END DO 176 IF( lk_mpp ) CALL mpp_min( glob_min_3d ) 177 ! 178 END FUNCTION glob_min_3d 179 180 181 FUNCTION glob_min_2d_a( ptab1, ptab2 ) 182 !!----------------------------------------------------------------------- 183 !! *** FUNCTION glob_min_2D _a *** 184 !! 185 !! ** Purpose : perform a masked min on the inner global domain of two 2D array 186 !!----------------------------------------------------------------------- 187 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab1, ptab2 ! input 2D array 188 REAL(wp) , DIMENSION(2) :: glob_min_2d_a ! global masked min 189 !!----------------------------------------------------------------------- 190 ! 191 glob_min_2d_a(1) = MINVAL( ptab1(:,:)*tmask_i(:,:) ) 192 glob_min_2d_a(2) = MINVAL( ptab2(:,:)*tmask_i(:,:) ) 193 IF( lk_mpp ) CALL mpp_min( glob_min_2d_a, 2 ) 194 ! 195 END FUNCTION glob_min_2d_a 196 197 198 FUNCTION glob_min_3d_a( ptab1, ptab2 ) 199 !!----------------------------------------------------------------------- 200 !! *** FUNCTION glob_min_3D_a *** 201 !! 202 !! ** Purpose : perform a masked min on the inner global domain of two 3D array 203 !!----------------------------------------------------------------------- 204 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab1, ptab2 ! input 3D array 205 REAL(wp) , DIMENSION(2) :: glob_min_3d_a ! global masked min 206 !! 207 INTEGER :: jk 208 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 209 !!----------------------------------------------------------------------- 210 ! 211 zjpk = SIZE(ptab1,3) 212 ! 213 glob_min_3d_a(:) = 0.e0 214 DO jk = 1, zjpk 215 glob_min_3d_a(1) = glob_min_3d_a(1) + MINVAL( ptab1(:,:,jk)*tmask_i(:,:) ) 216 glob_min_3d_a(2) = glob_min_3d_a(2) + MINVAL( ptab2(:,:,jk)*tmask_i(:,:) ) 217 END DO 218 IF( lk_mpp ) CALL mpp_min( glob_min_3d_a, 2 ) 219 ! 220 END FUNCTION glob_min_3d_a 221 222 ! --- MAX --- 223 FUNCTION glob_max_2d( ptab ) 224 !!----------------------------------------------------------------------- 225 !! *** FUNCTION glob_max_2D *** 226 !! 227 !! ** Purpose : perform a masked max on the inner global domain of a 2D array 228 !!----------------------------------------------------------------------- 229 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab ! input 2D array 230 REAL(wp) :: glob_max_2d ! global masked max 231 !!----------------------------------------------------------------------- 232 ! 233 glob_max_2d = MAXVAL( ptab(:,:)*tmask_i(:,:) ) 234 IF( lk_mpp ) CALL mpp_max( glob_max_2d ) 235 ! 236 END FUNCTION glob_max_2d 237 238 FUNCTION glob_max_3d( ptab ) 239 !!----------------------------------------------------------------------- 240 !! *** FUNCTION glob_max_3D *** 241 !! 242 !! ** Purpose : perform a masked max on the inner global domain of a 3D array 243 !!----------------------------------------------------------------------- 244 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab ! input 3D array 245 REAL(wp) :: glob_max_3d ! global masked max 246 !! 247 INTEGER :: jk 248 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 249 !!----------------------------------------------------------------------- 250 ! 251 zjpk = SIZE(ptab,3) 252 ! 253 glob_max_3d = 0.e0 254 DO jk = 1, zjpk 255 glob_max_3d = glob_max_3d + MAXVAL( ptab(:,:,jk)*tmask_i(:,:) ) 256 END DO 257 IF( lk_mpp ) CALL mpp_max( glob_max_3d ) 258 ! 259 END FUNCTION glob_max_3d 260 261 262 FUNCTION glob_max_2d_a( ptab1, ptab2 ) 263 !!----------------------------------------------------------------------- 264 !! *** FUNCTION glob_max_2D _a *** 265 !! 266 !! ** Purpose : perform a masked max on the inner global domain of two 2D array 267 !!----------------------------------------------------------------------- 268 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab1, ptab2 ! input 2D array 269 REAL(wp) , DIMENSION(2) :: glob_max_2d_a ! global masked max 270 !!----------------------------------------------------------------------- 271 ! 272 glob_max_2d_a(1) = MAXVAL( ptab1(:,:)*tmask_i(:,:) ) 273 glob_max_2d_a(2) = MAXVAL( ptab2(:,:)*tmask_i(:,:) ) 274 IF( lk_mpp ) CALL mpp_max( glob_max_2d_a, 2 ) 275 ! 276 END FUNCTION glob_max_2d_a 277 278 279 FUNCTION glob_max_3d_a( ptab1, ptab2 ) 280 !!----------------------------------------------------------------------- 281 !! *** FUNCTION glob_max_3D_a *** 282 !! 283 !! ** Purpose : perform a masked max on the inner global domain of two 3D array 284 !!----------------------------------------------------------------------- 285 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab1, ptab2 ! input 3D array 286 REAL(wp) , DIMENSION(2) :: glob_max_3d_a ! global masked max 287 !! 288 INTEGER :: jk 289 INTEGER :: zjpk ! local variable: size of the 3d dimension of ptab 290 !!----------------------------------------------------------------------- 291 ! 292 zjpk = SIZE(ptab1,3) 293 ! 294 glob_max_3d_a(:) = 0.e0 295 DO jk = 1, zjpk 296 glob_max_3d_a(1) = glob_max_3d_a(1) + MAXVAL( ptab1(:,:,jk)*tmask_i(:,:) ) 297 glob_max_3d_a(2) = glob_max_3d_a(2) + MAXVAL( ptab2(:,:,jk)*tmask_i(:,:) ) 298 END DO 299 IF( lk_mpp ) CALL mpp_max( glob_max_3d_a, 2 ) 300 ! 301 END FUNCTION glob_max_3d_a 302 123 303 124 304 #else … … 127 307 !!---------------------------------------------------------------------- 128 308 309 ! --- SUM --- 129 310 FUNCTION glob_sum_2d( ptab ) 130 311 !!---------------------------------------------------------------------- … … 133 314 !! ** Purpose : perform a sum in calling DDPDD routine 134 315 !!---------------------------------------------------------------------- 135 REAL(wp), INTENT(in), DIMENSION( jpi,jpj) :: ptab136 REAL(wp) 316 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab 317 REAL(wp) :: glob_sum_2d ! global masked sum 137 318 !! 138 319 COMPLEX(wp):: ctmp 139 320 REAL(wp) :: ztmp 140 321 INTEGER :: ji, jj ! dummy loop indices 141 !!----------------------------------------------------------------------- 142 ! 143 ztmp = 0.e0 144 ctmp = CMPLX( 0.e0, 0.e0, wp ) 145 DO jj = 1, jpj 146 DO ji =1, jpi 322 INTEGER :: zjpi, zjpj ! local variables: size of ptab 323 !!----------------------------------------------------------------------- 324 zjpi = SIZE(ptab,1) 325 zjpj = SIZE(ptab,2) 326 ! 327 ztmp = 0.e0 328 ctmp = CMPLX( 0.e0, 0.e0, wp ) 329 DO jj = 1, zjpj 330 DO ji =1, zjpi 147 331 ztmp = ptab(ji,jj) * tmask_i(ji,jj) 148 332 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) … … 161 345 !! ** Purpose : perform a sum on a 3D array in calling DDPDD routine 162 346 !!---------------------------------------------------------------------- 163 REAL(wp), INTENT(in), DIMENSION( jpi,jpj,jpk) :: ptab164 REAL(wp) 347 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab 348 REAL(wp) :: glob_sum_3d ! global masked sum 165 349 !! 166 350 COMPLEX(wp):: ctmp 167 351 REAL(wp) :: ztmp 168 352 INTEGER :: ji, jj, jk ! dummy loop indices 169 !!----------------------------------------------------------------------- 170 ! 171 ztmp = 0.e0 172 ctmp = CMPLX( 0.e0, 0.e0, wp ) 173 DO jk = 1, jpk 174 DO jj = 1, jpj 175 DO ji =1, jpi 353 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 354 !!----------------------------------------------------------------------- 355 ! 356 zjpi = SIZE(ptab,1) 357 zjpj = SIZE(ptab,2) 358 zjpk = SIZE(ptab,3) 359 ! 360 ztmp = 0.e0 361 ctmp = CMPLX( 0.e0, 0.e0, wp ) 362 DO jk = 1, zjpk 363 DO jj = 1, zjpj 364 DO ji =1, zjpi 176 365 ztmp = ptab(ji,jj,jk) * tmask_i(ji,jj) 177 366 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) … … 191 380 !! ** Purpose : perform a sum on two 2D arrays in calling DDPDD routine 192 381 !!---------------------------------------------------------------------- 193 REAL(wp), INTENT(in), DIMENSION( jpi,jpj) :: ptab1, ptab2194 REAL(wp) 382 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab1, ptab2 383 REAL(wp) :: glob_sum_2d_a ! global masked sum 195 384 !! 196 385 COMPLEX(wp):: ctmp 197 386 REAL(wp) :: ztmp 198 387 INTEGER :: ji, jj ! dummy loop indices 199 !!----------------------------------------------------------------------- 200 ! 201 ztmp = 0.e0 202 ctmp = CMPLX( 0.e0, 0.e0, wp ) 203 DO jj = 1, jpj 204 DO ji =1, jpi 388 INTEGER :: zjpi, zjpj ! local variables: size of ptab 389 !!----------------------------------------------------------------------- 390 ! 391 zjpi = SIZE(ptab1,1) 392 zjpj = SIZE(ptab1,2) 393 ! 394 ztmp = 0.e0 395 ctmp = CMPLX( 0.e0, 0.e0, wp ) 396 DO jj = 1, zjpj 397 DO ji =1, zjpi 205 398 ztmp = ptab1(ji,jj) * tmask_i(ji,jj) 206 399 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) … … 221 414 !! ** Purpose : perform a sum on two 3D array in calling DDPDD routine 222 415 !!---------------------------------------------------------------------- 223 REAL(wp), INTENT(in), DIMENSION( jpi,jpj,jpk) :: ptab1, ptab2224 REAL(wp) 416 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab1, ptab2 417 REAL(wp) :: glob_sum_3d_a ! global masked sum 225 418 !! 226 419 COMPLEX(wp):: ctmp 227 420 REAL(wp) :: ztmp 228 421 INTEGER :: ji, jj, jk ! dummy loop indices 229 !!----------------------------------------------------------------------- 230 ! 231 ztmp = 0.e0 232 ctmp = CMPLX( 0.e0, 0.e0, wp ) 233 DO jk = 1, jpk 234 DO jj = 1, jpj 235 DO ji =1, jpi 422 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 423 !!----------------------------------------------------------------------- 424 ! 425 zjpi = SIZE(ptab1,1) 426 zjpj = SIZE(ptab1,2) 427 zjpk = SIZE(ptab1,3) 428 ! 429 ztmp = 0.e0 430 ctmp = CMPLX( 0.e0, 0.e0, wp ) 431 DO jk = 1, zjpk 432 DO jj = 1, zjpj 433 DO ji =1, zjpi 236 434 ztmp = ptab1(ji,jj,jk) * tmask_i(ji,jj) 237 435 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) … … 247 445 248 446 447 ! --- MIN --- 448 FUNCTION glob_min_2d( ptab ) 449 !!---------------------------------------------------------------------- 450 !! *** FUNCTION glob_min_2d *** 451 !! 452 !! ** Purpose : perform a min in calling DDPDD routine 453 !!---------------------------------------------------------------------- 454 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab 455 REAL(wp) :: glob_min_2d ! global masked min 456 !! 457 COMPLEX(wp):: ctmp 458 REAL(wp) :: ztmp 459 INTEGER :: ji, jj ! dummy loop indices 460 INTEGER :: zjpi, zjpj ! local variables: size of ptab 461 !!----------------------------------------------------------------------- 462 zjpi = SIZE(ptab,1) 463 zjpj = SIZE(ptab,2) 464 ! 465 ztmp = 0.e0 466 ctmp = CMPLX( 0.e0, 0.e0, wp ) 467 DO jj = 1, zjpj 468 DO ji =1, zjpi 469 ztmp = ptab(ji,jj) * tmask_i(ji,jj) 470 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 471 END DO 472 END DO 473 IF( lk_mpp ) CALL mpp_min( ctmp ) ! min over the global domain 474 glob_min_2d = REAL(ctmp,wp) 475 ! 476 END FUNCTION glob_min_2d 477 478 479 FUNCTION glob_min_3d( ptab ) 480 !!---------------------------------------------------------------------- 481 !! *** FUNCTION glob_min_3d *** 482 !! 483 !! ** Purpose : perform a min on a 3D array in calling DDPDD routine 484 !!---------------------------------------------------------------------- 485 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab 486 REAL(wp) :: glob_min_3d ! global masked min 487 !! 488 COMPLEX(wp):: ctmp 489 REAL(wp) :: ztmp 490 INTEGER :: ji, jj, jk ! dummy loop indices 491 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 492 !!----------------------------------------------------------------------- 493 ! 494 zjpi = SIZE(ptab,1) 495 zjpj = SIZE(ptab,2) 496 zjpk = SIZE(ptab,3) 497 ! 498 ztmp = 0.e0 499 ctmp = CMPLX( 0.e0, 0.e0, wp ) 500 DO jk = 1, zjpk 501 DO jj = 1, zjpj 502 DO ji =1, zjpi 503 ztmp = ptab(ji,jj,jk) * tmask_i(ji,jj) 504 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 505 END DO 506 END DO 507 END DO 508 IF( lk_mpp ) CALL mpp_min( ctmp ) ! min over the global domain 509 glob_min_3d = REAL(ctmp,wp) 510 ! 511 END FUNCTION glob_min_3d 512 513 514 FUNCTION glob_min_2d_a( ptab1, ptab2 ) 515 !!---------------------------------------------------------------------- 516 !! *** FUNCTION glob_min_2d_a *** 517 !! 518 !! ** Purpose : perform a min on two 2D arrays in calling DDPDD routine 519 !!---------------------------------------------------------------------- 520 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab1, ptab2 521 REAL(wp) :: glob_min_2d_a ! global masked min 522 !! 523 COMPLEX(wp):: ctmp 524 REAL(wp) :: ztmp 525 INTEGER :: ji, jj ! dummy loop indices 526 INTEGER :: zjpi, zjpj ! local variables: size of ptab 527 !!----------------------------------------------------------------------- 528 ! 529 zjpi = SIZE(ptab1,1) 530 zjpj = SIZE(ptab1,2) 531 ! 532 ztmp = 0.e0 533 ctmp = CMPLX( 0.e0, 0.e0, wp ) 534 DO jj = 1, zjpj 535 DO ji =1, zjpi 536 ztmp = ptab1(ji,jj) * tmask_i(ji,jj) 537 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 538 ztmp = ptab2(ji,jj) * tmask_i(ji,jj) 539 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 540 END DO 541 END DO 542 IF( lk_mpp ) CALL mpp_min( ctmp ) ! min over the global domain 543 glob_min_2d_a = REAL(ctmp,wp) 544 ! 545 END FUNCTION glob_min_2d_a 546 547 548 FUNCTION glob_min_3d_a( ptab1, ptab2 ) 549 !!---------------------------------------------------------------------- 550 !! *** FUNCTION glob_min_3d_a *** 551 !! 552 !! ** Purpose : perform a min on two 3D array in calling DDPDD routine 553 !!---------------------------------------------------------------------- 554 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab1, ptab2 555 REAL(wp) :: glob_min_3d_a ! global masked min 556 !! 557 COMPLEX(wp):: ctmp 558 REAL(wp) :: ztmp 559 INTEGER :: ji, jj, jk ! dummy loop indices 560 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 561 !!----------------------------------------------------------------------- 562 ! 563 zjpi = SIZE(ptab1,1) 564 zjpj = SIZE(ptab1,2) 565 zjpk = SIZE(ptab1,3) 566 ! 567 ztmp = 0.e0 568 ctmp = CMPLX( 0.e0, 0.e0, wp ) 569 DO jk = 1, zjpk 570 DO jj = 1, zjpj 571 DO ji =1, zjpi 572 ztmp = ptab1(ji,jj,jk) * tmask_i(ji,jj) 573 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 574 ztmp = ptab2(ji,jj,jk) * tmask_i(ji,jj) 575 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 576 END DO 577 END DO 578 END DO 579 IF( lk_mpp ) CALL mpp_min( ctmp ) ! min over the global domain 580 glob_min_3d_a = REAL(ctmp,wp) 581 ! 582 END FUNCTION glob_min_3d_a 583 584 585 ! --- MAX --- 586 FUNCTION glob_max_2d( ptab ) 587 !!---------------------------------------------------------------------- 588 !! *** FUNCTION glob_max_2d *** 589 !! 590 !! ** Purpose : perform a max in calling DDPDD routine 591 !!---------------------------------------------------------------------- 592 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab 593 REAL(wp) :: glob_max_2d ! global masked max 594 !! 595 COMPLEX(wp):: ctmp 596 REAL(wp) :: ztmp 597 INTEGER :: ji, jj ! dummy loop indices 598 INTEGER :: zjpi, zjpj ! local variables: size of ptab 599 !!----------------------------------------------------------------------- 600 zjpi = SIZE(ptab,1) 601 zjpj = SIZE(ptab,2) 602 ! 603 ztmp = 0.e0 604 ctmp = CMPLX( 0.e0, 0.e0, wp ) 605 DO jj = 1, zjpj 606 DO ji =1, zjpi 607 ztmp = ptab(ji,jj) * tmask_i(ji,jj) 608 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 609 END DO 610 END DO 611 IF( lk_mpp ) CALL mpp_max( ctmp ) ! max over the global domain 612 glob_max_2d = REAL(ctmp,wp) 613 ! 614 END FUNCTION glob_max_2d 615 616 617 FUNCTION glob_max_3d( ptab ) 618 !!---------------------------------------------------------------------- 619 !! *** FUNCTION glob_max_3d *** 620 !! 621 !! ** Purpose : perform a max on a 3D array in calling DDPDD routine 622 !!---------------------------------------------------------------------- 623 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab 624 REAL(wp) :: glob_max_3d ! global masked max 625 !! 626 COMPLEX(wp):: ctmp 627 REAL(wp) :: ztmp 628 INTEGER :: ji, jj, jk ! dummy loop indices 629 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 630 !!----------------------------------------------------------------------- 631 ! 632 zjpi = SIZE(ptab,1) 633 zjpj = SIZE(ptab,2) 634 zjpk = SIZE(ptab,3) 635 ! 636 ztmp = 0.e0 637 ctmp = CMPLX( 0.e0, 0.e0, wp ) 638 DO jk = 1, zjpk 639 DO jj = 1, zjpj 640 DO ji =1, zjpi 641 ztmp = ptab(ji,jj,jk) * tmask_i(ji,jj) 642 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 643 END DO 644 END DO 645 END DO 646 IF( lk_mpp ) CALL mpp_max( ctmp ) ! max over the global domain 647 glob_max_3d = REAL(ctmp,wp) 648 ! 649 END FUNCTION glob_max_3d 650 651 652 FUNCTION glob_max_2d_a( ptab1, ptab2 ) 653 !!---------------------------------------------------------------------- 654 !! *** FUNCTION glob_max_2d_a *** 655 !! 656 !! ** Purpose : perform a max on two 2D arrays in calling DDPDD routine 657 !!---------------------------------------------------------------------- 658 REAL(wp), INTENT(in), DIMENSION(:,:) :: ptab1, ptab2 659 REAL(wp) :: glob_max_2d_a ! global masked max 660 !! 661 COMPLEX(wp):: ctmp 662 REAL(wp) :: ztmp 663 INTEGER :: ji, jj ! dummy loop indices 664 INTEGER :: zjpi, zjpj ! local variables: size of ptab 665 !!----------------------------------------------------------------------- 666 ! 667 zjpi = SIZE(ptab1,1) 668 zjpj = SIZE(ptab1,2) 669 ! 670 ztmp = 0.e0 671 ctmp = CMPLX( 0.e0, 0.e0, wp ) 672 DO jj = 1, zjpj 673 DO ji =1, zjpi 674 ztmp = ptab1(ji,jj) * tmask_i(ji,jj) 675 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 676 ztmp = ptab2(ji,jj) * tmask_i(ji,jj) 677 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 678 END DO 679 END DO 680 IF( lk_mpp ) CALL mpp_max( ctmp ) ! max over the global domain 681 glob_max_2d_a = REAL(ctmp,wp) 682 ! 683 END FUNCTION glob_max_2d_a 684 685 686 FUNCTION glob_max_3d_a( ptab1, ptab2 ) 687 !!---------------------------------------------------------------------- 688 !! *** FUNCTION glob_max_3d_a *** 689 !! 690 !! ** Purpose : perform a max on two 3D array in calling DDPDD routine 691 !!---------------------------------------------------------------------- 692 REAL(wp), INTENT(in), DIMENSION(:,:,:) :: ptab1, ptab2 693 REAL(wp) :: glob_max_3d_a ! global masked max 694 !! 695 COMPLEX(wp):: ctmp 696 REAL(wp) :: ztmp 697 INTEGER :: ji, jj, jk ! dummy loop indices 698 INTEGER :: zjpi, zjpj, zjpk ! local variables: size of ptab 699 !!----------------------------------------------------------------------- 700 ! 701 zjpi = SIZE(ptab1,1) 702 zjpj = SIZE(ptab1,2) 703 zjpk = SIZE(ptab1,3) 704 ! 705 ztmp = 0.e0 706 ctmp = CMPLX( 0.e0, 0.e0, wp ) 707 DO jk = 1, zjpk 708 DO jj = 1, zjpj 709 DO ji =1, zjpi 710 ztmp = ptab1(ji,jj,jk) * tmask_i(ji,jj) 711 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 712 ztmp = ptab2(ji,jj,jk) * tmask_i(ji,jj) 713 CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 714 END DO 715 END DO 716 END DO 717 IF( lk_mpp ) CALL mpp_max( ctmp ) ! max over the global domain 718 glob_max_3d_a = REAL(ctmp,wp) 719 ! 720 END FUNCTION glob_max_3d_a 721 722 249 723 SUBROUTINE DDPDD( ydda, yddb ) 250 724 !!----------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.