Changeset 7824 for branches/2015
- Timestamp:
- 2017-03-23T14:10:34+01:00 (7 years ago)
- Location:
- branches/2015/dev_r5803_UKMO_AGRIF_Vert_interp/NEMOGCM/NEMO/NST_SRC
- Files:
-
- 2 edited
-
agrif_opa_interp.F90 (modified) (9 diffs)
-
agrif_opa_update.F90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5803_UKMO_AGRIF_Vert_interp/NEMOGCM/NEMO/NST_SRC/agrif_opa_interp.F90
r6930 r7824 714 714 ! VERTICAL REFINEMENT BEGIN 715 715 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child 716 REAL(wp), DIMENSION(k1:k2,n1:n2-1) :: tabin 716 717 REAL(wp) :: h_in(k1:k2) 717 718 REAL(wp) :: h_out(1:jpk) 718 719 INTEGER :: N_in, N_out 719 720 REAL(wp) :: h_diff 721 REAL(wp) :: zrhoxy 720 722 ! VERTICAL REFINEMENT END 721 723 724 zrhoxy = Agrif_rhox()*Agrif_rhoy() 722 725 IF (before) THEN 723 ptab(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 726 DO jn = n1,n2-1 727 DO jk=k1,k2 728 DO jj=j1,j2 729 DO ji=i1,i2 730 ptab(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * e1e2t(ji,jj) * e3t_n(ji,jj,jk) 731 END DO 732 END DO 733 END DO 734 END DO 735 DO jk=k1,k2 736 DO jj=j1,j2 737 DO ji=i1,i2 738 ptab(ji,jj,jk,n2) = tmask(ji,jj,jk) * e1e2t(ji,jj) * e3t_n(ji,jj,jk) 739 END DO 740 END DO 741 END DO 742 724 743 ELSE 725 744 ! VERTICAL REFINEMENT BEGIN … … 728 747 do jj=j1,j2 729 748 do ji=i1,i2 730 h_in(k1:k2) = interp_scales_t(ji,jj,k1:k2) 731 h_out(1:jpk) = e3t_n(ji,jj,1:jpk) 732 h_diff = sum(h_out(1:jpk-1))-sum(h_in(k1:k2-1)) 733 N_in = k2-1 749 N_in = 0 750 DO jk=k1,k2 !k2 = jpk of parent grid 751 IF (ptab(ji,jj,jk,n2) == 0) EXIT 752 N_in = N_in + 1 753 tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)/ptab(ji,jj,jk,n2) 754 h_in(N_in) = ptab(ji,jj,jk,n2)/(e1e2t(ji,jj)*zrhoxy) 755 END DO 734 756 N_out = jpk-1 735 if (h_diff > 0) then 736 h_in(N_in+1) = h_diff 737 N_in = N_in + 1 738 else 739 h_out(N_out+1) = -h_diff 740 N_out = N_out + 1 741 endif 742 ptab(ji,jj,k2,:) = ptab(ji,jj,k2-1,:) 757 DO jk=1,jpk ! jpk of child grid 758 ! IF (tmask(ji,jj,jk) == 0) EXIT ! TODO: Will not work with ISF. !This doesn't seem to work at the moment. Is it just a GYRE issue??????? 759 ! N_out = N_out + 1 760 h_out(jk) = e3t_n(ji,jj,jk) !Child grid scale factors. Could multiply by e1e2t here instead of division above 761 ENDDO 762 IF (N_in > 0) THEN 763 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 764 ! if (h_diff > 0) then 765 ! h_in(N_in+1) = h_diff 766 ! N_in = N_in + 1 767 ! else 768 ! h_out(N_out+1) = -h_diff 769 ! N_out = N_out + 1 770 ! endif 771 ptab(ji,jj,k2,:) = ptab(ji,jj,k2-1,:) !what is this line for????? 743 772 do jn=1,jpts 744 call reconstructandremap( ptab(ji,jj,1:N_in,jn),h_in,ptab_child(ji,jj,1:N_out,jn),h_out,N_in,N_out)773 call reconstructandremap(tabin(1:N_in,jn),h_in,ptab_child(ji,jj,1:N_out,jn),h_out,N_in,N_out) 745 774 enddo 746 775 ! if (abs(h_diff) > 1000.) then … … 751 780 ! enddo 752 781 ! endif 782 ENDIF 753 783 enddo 754 784 enddo 785 755 786 756 787 ! VERTICAL REFINEMENT END … … 950 981 DO jj=j1,j2 951 982 DO ji=i1,i2 952 ptab(ji,jj,jk,1) = e2u(ji,jj) * un(ji,jj,jk) * e3u_n(ji,jj,jk)/zrhoy953 ptab(ji,jj,jk,2) = e2u(ji,jj) * e3u_n(ji,jj,jk)/zrhoy983 ptab(ji,jj,jk,1) = e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk) 984 ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk) 954 985 END DO 955 986 END DO … … 972 1003 N_in = N_in + 1 973 1004 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 974 h_in(N_in) = ptab(ji,jj,jk,2)/ e2u(ji,jj)1005 h_in(N_in) = ptab(ji,jj,jk,2)/(e2u(ji,jj) * zrhoy) 975 1006 ENDDO 976 1007 … … 992 1023 ENDIF 993 1024 994 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 995 IF (abs(h_diff) > 0.) THEN 996 CALL ctl_stop 997 ENDIF 1025 ! IF (N_in * N_out > 0) THEN 1026 ! h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 1027 ! Should be able to remove the next IF/ELSEIF statement once scale factors are dealt with properly 1028 ! if (h_diff < 0.) then 1029 ! print *,'CHECK YOUR BATHY ...', h_diff, sum(h_out(1:N_out)), sum(h_in(1:N_in)) 1030 ! stop 1031 ! endif 1032 ! ENDIF 998 1033 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 999 1034 … … 1009 1044 ua(i1:i2,jj,jk) = ptab_child(i1:i2,jj,jk) 1010 1045 !/(zrhoy*e2u(i1:i2,jj))) 1011 write(numout,*) 'ua=',ua(i1:i2,jj,jk)1012 1046 END DO 1013 1047 END DO … … 1043 1077 !!--------------------------------------------- 1044 1078 ! 1079 zrhox = Agrif_rhox() 1045 1080 IF (before) THEN 1046 !interpv entre 1 et k2 et interpv2d en jpkp1 1047 DO jk=k1,jpk 1081 DO jk=k1,k2 1048 1082 DO jj=j1,j2 1049 1083 DO ji=i1,i2 1050 ptab(ji,jj,jk,1) = e1v(ji,jj) * vn(ji,jj,jk) 1051 ptab(ji,jj,jk,1) = ptab(ji,jj,jk,1) * e3v_n(ji,jj,jk) 1052 ptab(ji,jj,jk,2) = e3v_n(ji,jj,jk) 1084 ptab(ji,jj,jk,1) = e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk) 1085 ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v_n(ji,jj,jk) 1053 1086 END DO 1054 1087 END DO … … 1060 1093 northern_side = (nb == 2).AND.(ndir == 2) 1061 1094 do jj=j1,j2 1062 jref = jj1063 IF (southern_side) jref = 21064 IF (northern_side) jref = nlcj-21065 do ji=i1,i21066 1067 N_in = 01068 do jk=k1,k21069 if (interp_scales_v(ji,jj,jk) == 0) EXIT1070 N_in = N_in + 11071 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2)1072 h_in(N_in) = ptab(ji,jj,jk,2)1073 enddo1074 IF (N_in == 0) THEN1075 ptab_child(ji,jj,:) = 0.1076 CYCLE1077 ENDIF1095 jref = jj 1096 IF (southern_side) jref = 2 1097 IF (northern_side) jref = nlcj-2 1098 do ji=i1,i2 1099 1100 N_in = 0 1101 do jk=k1,k2 1102 if (ptab(ji,jj,jk,2) == 0) EXIT 1103 N_in = N_in + 1 1104 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 1105 h_in(N_in) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox) 1106 enddo 1107 IF (N_in == 0) THEN 1108 ptab_child(ji,jj,:) = 0. 1109 CYCLE 1110 ENDIF 1078 1111 1079 N_out = 0 1080 do jk=1,jpk 1081 if (vmask(ji,jref,jk) == 0) EXIT 1082 N_out = N_out + 1 1083 h_out(N_out) = e3v_n(ji,jj,jk) 1084 enddo 1085 IF (N_out == 0) THEN 1086 ptab_child(ji,jj,:) = 0. 1087 CYCLE 1088 ENDIF 1112 N_out = 0 1113 do jk=1,jpk 1114 if (vmask(ji,jref,jk) == 0) EXIT 1115 N_out = N_out + 1 1116 h_out(N_out) = e3v_n(ji,jj,jk) 1117 enddo 1118 IF (N_out == 0) THEN 1119 ptab_child(ji,jj,:) = 0. 1120 CYCLE 1121 ENDIF 1122 1123 ! IF (N_in * N_out > 0) THEN 1124 ! h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 1125 ! Should be able to remove the next IF/ELSEIF statement once scale factors are dealt with properly 1126 ! if (h_diff < 0.) then 1127 ! print *,'CHECK YOUR BATHY interpvn...', h_diff, sum(h_out(1:N_out)), sum(h_in(1:N_in)) 1128 ! stop 1129 ! endif 1130 ! ENDIF 1131 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 1089 1132 1090 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in)) 1091 IF (h_diff > 0.) THEN 1092 N_in = N_in + 1 1093 h_in(N_in) = h_diff 1094 tabin(N_in) = 0. 1095 ELSE 1096 h_out(N_out) = h_out(N_out) - h_diff 1097 ENDIF 1098 1099 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 1100 1101 ptab_child(ji,jj,N_out) = ptab_child(ji,jj,N_out) * h_out(N_out) / e3v_n(ji,jj,N_out) 1102 1103 enddo 1133 enddo 1104 1134 enddo 1105 1135 ! in the following 1106 ! remove division of va by fs e3v (already done)1136 ! remove division of va by fs e3v, zrhox and e1v (already done) 1107 1137 ! VERTICAL REFINEMENT END 1108 zrhox= Agrif_Rhox()1109 1138 DO jk=1,jpkm1 1110 1139 DO jj=j1,j2 1111 va(i1:i2,jj,jk) = (ptab_child(i1:i2,jj,jk)/(zrhox*e1v(i1:i2,jj)))1140 va(i1:i2,jj,jk) = ptab_child(i1:i2,jj,jk) 1112 1141 END DO 1113 1142 END DO -
branches/2015/dev_r5803_UKMO_AGRIF_Vert_interp/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90
r6929 r7824 1 # undefTWO_WAY /* TWO WAY NESTING */1 #define TWO_WAY /* TWO WAY NESTING */ 2 2 #undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/ 3 3 … … 121 121 ! 122 122 IF (MOD(nbcline,nbclineupdate) == 0) THEN 123 WRITE(numout,*) 'TG print 1'124 CALL FLUSH(numout)125 123 # if ! defined DECAL_FEEDBACK 126 124 CALL Agrif_Update_Variable(tsn_id, procname=updateTS) … … 128 126 CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS) 129 127 # endif 130 WRITE(numout,*) 'TG print 2' 131 CALL FLUSH(numout) 132 ELSE 133 WRITE(numout,*) 'TG print 3' 134 CALL FLUSH(numout) 128 ELSE 135 129 # if ! defined DECAL_FEEDBACK 136 130 CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS) … … 138 132 CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS) 139 133 # endif 140 WRITE(numout,*) 'TG print 4'141 CALL FLUSH(numout)142 134 ENDIF 143 135 ! … … 380 372 ENDIF 381 373 ! 382 WRITE(numout,*) 'I got to end of updateTS before=',before383 CALL FLUSH(numout)384 374 END SUBROUTINE updateTS 385 375 … … 396 386 REAL(wp) :: zrhoy 397 387 ! VERTICAL REFINEMENT BEGIN 398 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk ,2) :: ptab_child388 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: ptab_child 399 389 REAL(wp) :: h_in(k1:k2) 400 390 REAL(wp) :: h_out(1:jpk) … … 405 395 !!--------------------------------------------- 406 396 ! 407 WRITE(numout,*) 'TG print 5: Start of updateu before = ',before408 CALL FLUSH(numout)409 397 IF( before ) THEN 410 398 zrhoy = Agrif_Rhoy() … … 412 400 DO jj=j1,j2 413 401 DO ji=i1,i2 414 ptab(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u_n(ji,jj,jk) * u n(ji,jj,jk)402 ptab(ji,jj,jk,1) = zrhoy * e2u(ji,jj) * e3u_n(ji,jj,jk) * umask(ji,jj,jk) * un(ji,jj,jk) 415 403 ptab(ji,jj,jk,2) = umask(ji,jj,jk) * zrhoy * e2u(ji,jj) * e3u_n(ji,jj,jk) 416 404 END DO … … 419 407 ELSE 420 408 ! VERTICAL REFINEMENT BEGIN 421 ptab_child(:,:,: ,:) = 0.409 ptab_child(:,:,:) = 0. 422 410 423 411 DO jj=j1,j2 … … 427 415 IF (ptab(ji,jj,jk,2) == 0) EXIT 428 416 N_in = N_in + 1 429 tabin(jk) = ptab(ji,jj,jk )/ptab(ji,jj,jk,2)417 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 430 418 h_in(N_in) = ptab(ji,jj,jk,2)/e2u(ji,jj) 431 419 ENDDO … … 447 435 ! h_in(N_in) = h_diff 448 436 endif 449 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out ,1),h_out(1:N_out),N_in,N_out)437 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 450 438 ENDIF 451 439 ENDDO … … 466 454 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 467 455 ub(ji,jj,jk) = ub(ji,jj,jk) & 468 & + atfp * ( ptab_child(ji,jj,jk ,1) - un(ji,jj,jk) ) * umask(ji,jj,jk)456 & + atfp * ( ptab_child(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk) 469 457 ENDIF 470 458 ! 471 un(ji,jj,jk) = ptab_child(ji,jj,jk ,1) * umask(ji,jj,jk)459 un(ji,jj,jk) = ptab_child(ji,jj,jk) * umask(ji,jj,jk) 472 460 END DO 473 461 END DO … … 475 463 ENDIF 476 464 ! 477 WRITE(numout,*) 'TG print 6: End of updateu before = ',before478 CALL FLUSH(numout)479 465 END SUBROUTINE updateu 480 466 … … 491 477 REAL(wp) :: zrhox 492 478 ! VERTICAL REFINEMENT BEGIN 493 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk ,2) :: ptab_child479 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: ptab_child 494 480 REAL(wp) :: h_in(k1:k2) 495 481 REAL(wp) :: h_out(1:jpk) … … 500 486 !!--------------------------------------------- 501 487 ! 502 WRITE(numout,*) 'TG print 7: Start of updatev before = ',before503 CALL FLUSH(numout)504 488 IF (before) THEN 505 489 zrhox = Agrif_Rhox() … … 507 491 DO jj=j1,j2 508 492 DO ji=i1,i2 509 ptab(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * v n(ji,jj,jk)493 ptab(ji,jj,jk,1) = zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) * vmask(ji,jj,jk) * vn(ji,jj,jk) 510 494 ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * zrhox * e1v(ji,jj) * e3v_n(ji,jj,jk) 511 495 END DO … … 520 504 N_in = 0 521 505 DO jk=k1,k2 522 IF ( update_scales_v(ji,jj,jk) == 0) EXIT506 IF (ptab(ji,jj,jk,2) == 0) EXIT 523 507 N_in = N_in + 1 524 tabin(jk) = ptab(ji,jj, 1)/ptab(ji,jj,jk,2)508 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 525 509 h_in(N_in) = ptab(ji,jj,jk,2)/e1v(ji,jj) 526 510 ENDDO … … 542 526 ! h_in(N_in) = h_diff 543 527 endif 544 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out ,1),h_out(1:N_out),N_in,N_out)528 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out) 545 529 ENDIF 546 530 ENDDO … … 553 537 ! VERTICAL REFINEMENT END 554 538 555 DO jk=k1, k2539 DO jk=k1,jpk 556 540 DO jj=j1,j2 557 541 DO ji=i1,i2 … … 561 545 IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 562 546 vb(ji,jj,jk) = vb(ji,jj,jk) & 563 & + atfp * ( ptab_child(ji,jj,jk ,1) - vn(ji,jj,jk) ) * vmask(ji,jj,jk)547 & + atfp * ( ptab_child(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk) 564 548 ENDIF 565 549 ! 566 vn(ji,jj,jk) = ptab_child(ji,jj,jk,1) * vmask(ji,jj,jk) 567 END DO 568 END DO 569 END DO 570 ENDIF 571 WRITE(numout,*) 'TG print 8: End of updatev before = ',before 572 ! 550 vn(ji,jj,jk) = ptab_child(ji,jj,jk) * vmask(ji,jj,jk) 551 END DO 552 END DO 553 END DO 554 ENDIF 573 555 END SUBROUTINE updatev 574 556
Note: See TracChangeset
for help on using the changeset viewer.
