Changeset 4594
- Timestamp:
- 2014-03-26T11:47:55+01:00 (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/TRA/zpshde_tam.F90
r3660 r4594 35 35 PUBLIC zps_hde_tan ! routine called by step_tam.F90 36 36 PUBLIC zps_hde_adj ! routine called by step_tam.F90 37 !PUBLIC zps_hde_adj_tst ! routine called by tst.F9037 PUBLIC zps_hde_adj_tst ! routine called by tst.F90 38 38 39 39 !! * Substitutions … … 104 104 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT( in ) :: & 105 105 pta, pta_tl ! 3D T, S and rd direct fields 106 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in out), OPTIONAL :: &106 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ), OPTIONAL :: & 107 107 prd_tl 108 108 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out ), OPTIONAL :: & … … 261 261 END SUBROUTINE zps_hde_tan 262 262 263 SUBROUTINE zps_hde_adj ( kt, kjpt, pta, pgtu, pgtv,&263 SUBROUTINE zps_hde_adj ( kt, kjpt, pta, & 264 264 & pta_ad, prd_ad, & 265 265 & pgtu_ad, pgru_ad, & … … 325 325 prd_ad ! 3D T, S and rd tangent fields 326 326 REAL(wp), DIMENSION(jpi,jpj,kjpt), INTENT( inout ), OPTIONAL :: & 327 pgtu , pgtv, pgtu_ad, pgtv_ad ! 3D T, S and rd tangent fields327 pgtu_ad, pgtv_ad ! 3D T, S and rd tangent fields 328 328 REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ), OPTIONAL :: & 329 329 pgru_ad, & ! horizontal grad. of T, S and rd at u- … … 368 368 ! interpolated values of tracers 369 369 zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 370 ! gradient of tracers371 pgtu(ji,jj,jn) = umask(ji,jj,1) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) )372 370 ELSE ! case 2 373 371 zmaxu = -ze3wu / fse3w(ji,jj,iku) 374 372 ! interpolated values of tracers 375 373 zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 376 ! gradient of tracers377 pgtu(ji,jj,jn) = umask(ji,jj,1) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) )378 374 ENDIF 379 375 ! … … 383 379 ! interpolated values of tracers 384 380 ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 385 ! gradient of tracers386 pgtv(ji,jj,jn) = vmask(ji,jj,1) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) )387 381 ELSE ! case 2 388 382 zmaxv = -ze3wv / fse3w(ji,jj,ikv) 389 383 ! interpolated values of tracers 390 384 ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 391 ! gradient of tracers392 pgtv(ji,jj,jn) = vmask(ji,jj,1) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) )393 385 ENDIF 394 386 # if ! defined key_vectopt_loop … … 396 388 # endif 397 389 END DO 398 CALL lbc_lnk( pgtu(:,:,jn), 'U', -1. ) ; CALL lbc_lnk( pgtv(:,:,jn), 'V', -1. ) ! Lateral boundary cond.399 390 ! 400 391 END DO … … 566 557 ! 567 558 END SUBROUTINE zps_hde_adj 568 !SUBROUTINE zps_hde_adj_tst( kumadt ) 569 !!!----------------------------------------------------------------------- 570 !!! 571 !!! *** ROUTINE zps_hde_adj_tst *** 572 !!! 573 !!! ** Purpose : Test the adjoint routine. 574 !!! 575 !!! ** Method : Verify the scalar product 576 !!! 577 !!! ( L dx )^T W dy = dx^T L^T W dy 578 !!! 579 !!! where L = tangent routine 580 !!! L^T = adjoint routine 581 !!! W = diagonal matrix of scale factors 582 !!! dx = input perturbation (random field) 583 !!! dy = L dx 584 !!! 585 !!! 586 !!! History : 587 !!! ! 08-08 (A. Vidard) 588 !!!----------------------------------------------------------------------- 589 !!! * Modules used 590 591 !!! * Arguments 592 !INTEGER, INTENT(IN) :: & 593 !& kumadt ! Output unit 594 595 !INTEGER :: & 596 !& ji, & ! dummy loop indices 597 !& jj, & 598 !& jk, & 599 !& kt, & 600 !& jn 601 602 !!! * Local declarations 603 !REAL(KIND=wp), DIMENSION(:,:,:,:), ALLOCATABLE :: & 604 !& zts, & ! Direct field : temperature 605 !& zts_tlin, & ! Tangent input: temperature 606 !& zts_adout, & ! Adjoint output: temperature 607 !& zats ! 3D random field for t 608 609 !REAL(KIND=wp), DIMENSION(:,:,:), ALLOCATABLE :: & 610 !& zgtu_tlout, & ! Tangent output: horizontal gradient 611 !& zgtv_tlout, & ! Tangent output: horizontal gradient 612 !& zrd_adout, & ! Adjoint output: 613 !& zar, & ! 3D random field for rd 614 !& zrd_tlin, & ! Tangent input: 615 !& zgtu_adin, & ! Adjoint input : horizontal gradient 616 !& zgtv_adin ! Adjoint input : horizontal gradient 617 618 !REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: & 619 !& zgru_tlout, & ! Tangent output: horizontal gradient 620 !& zgrv_tlout, & ! Tangent output: horizontal gradient 621 !& zgru_adin, & ! Adjoint input : horizontal gradient 622 !& zgrv_adin ! Adjoint input : horizontal gradient 623 624 !REAL(KIND=wp) :: & 625 !! random field standard deviation for: 626 !& zsp1, & ! scalar product involving the tangent routine 627 !& zsp1_1, & ! scalar product components 628 !& zsp1_2, & 629 !& zsp1_3, & ! scalar product components 630 !& zsp1_4, & 631 !& zsp1_5, & ! scalar product components 632 !& zsp1_6, & 633 !& zsp2, & ! scalar product involving the adjoint routine 634 !& zsp2_1, & ! scalar product components 635 !& zsp2_2, & 636 !& zsp2_3 637 !CHARACTER (LEN=14) :: & 638 !& cl_name 639 640 !kt = nit000 641 !! Allocate memory 642 !ALLOCATE( & 643 !& zts(jpi,jpj,jpk,jpts), & 644 !& zts_tlin(jpi,jpj,jpk,jpts), & 645 !& zrd_tlin(jpi,jpj,jpk), & 646 !& zts_adout(jpi,jpj,jpk,jpts), & 647 !& zrd_adout(jpi,jpj,jpk), & 648 !& zar(jpi,jpj,jpk), & 649 !& zats(jpi,jpj,jpk,jpts), & 650 !& zgtu_tlout(jpi,jpj,jpts), & 651 !& zgtv_tlout(jpi,jpj,jpts), & 652 !& zgru_tlout(jpi,jpj), & 653 !& zgrv_tlout(jpi,jpj), & 654 !& zgtu_adin(jpi,jpj,jpts), & 655 !& zgtv_adin(jpi,jpj,jpts), & 656 !& zgru_adin(jpi,jpj), & 657 !& zgrv_adin(jpi,jpj) & 658 !& ) 659 !! Initialize random field standard deviationsthe reference state 660 !zts = tsn(:,:,:,:) 661 662 !!============================================================= 663 !! 1) dx = ( T ) and dy = ( T ) 664 !!============================================================= 665 666 !!-------------------------------------------------------------------- 667 !! Reset the tangent and adjoint variables 668 !!-------------------------------------------------------------------- 669 !zts_tlin(:,:,:,:) = 0.0_wp 670 !zrd_tlin(:,:,:) = 0.0_wp 671 !zts_adout(:,:,:,:) = 0.0_wp 672 !zrd_adout(:,:,:) = 0.0_wp 673 !zgtu_tlout(:,:,:) = 0.0_wp 674 !zgtv_tlout(:,:,:) = 0.0_wp 675 !zgru_tlout(:,:) = 0.0_wp 676 !zgrv_tlout(:,:) = 0.0_wp 677 !zgtu_adin(:,:,:) = 0.0_wp 678 !zgtv_adin(:,:,:) = 0.0_wp 679 !zgru_adin(:,:) = 0.0_wp 680 !zgrv_adin(:,:) = 0.0_wp 681 682 !!-------------------------------------------------------------------- 683 !! Initialize the tangent input with random noise: dx 684 !!-------------------------------------------------------------------- 685 !DO jn = 1, jpts 686 !CALL grid_random( zats(:,:,:,jn), 'T', 0.0_wp, stdt ) 687 !END DO 688 !CALL grid_random( zar, 'T', 0.0_wp, stdr ) 689 690 !zts_tlin(:,:,:,:) = zats(:,:,:,:) 691 !zrd_tlin(:,:,:) = zar(:,:,:) 692 !CALL zps_hde_tan ( nit000, jpts, zts, & 693 !& zts_tlin , zrd_tlin , & 694 !& zgtu_tlout, zgru_tlout, & 695 !& zgtv_tlout, zgrv_tlout ) 696 !DO jn = 1, jpts 697 !DO jj = nldj, nlej 698 !DO ji = nldi, nlei 699 !jk = mbku(ji,jj) 700 !zgtu_adin(ji,jj,jn) = zgtu_tlout(ji,jj,jn) & 701 !& * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) & 702 !& * umask(ji,jj,jk) 703 !jk = mbkv(ji,jj) 704 !zgtv_adin(ji,jj,jn) = zgtv_tlout(ji,jj,jn) & 705 !& * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) & 706 !& * vmask(ji,jj,jk) 707 !END DO 708 !END DO 709 !END DO 710 !DO jj = nldj, nlej 711 !DO ji = nldi, nlei 712 !jk = mbku(ji,jj) 713 !zgru_adin(ji,jj) = zgru_tlout(ji,jj) & 714 !& * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) & 715 !& * umask(ji,jj,jk) 716 !jk = mbkv(ji,jj) 717 !zgrv_adin(ji,jj) = zgrv_tlout(ji,jj) & 718 !& * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) & 719 !& * vmask(ji,jj,jk) 720 !END DO 721 !END DO 559 SUBROUTINE zps_hde_adj_tst( kumadt ) 560 !!----------------------------------------------------------------------- 561 !! 562 !! *** ROUTINE zps_hde_adj_tst *** 563 !! 564 !! ** Purpose : Test the adjoint routine. 565 !! 566 !! ** Method : Verify the scalar product 567 !! 568 !! ( L dx )^T W dy = dx^T L^T W dy 569 !! 570 !! where L = tangent routine 571 !! L^T = adjoint routine 572 !! W = diagonal matrix of scale factors 573 !! dx = input perturbation (random field) 574 !! dy = L dx 575 !! 576 !! 577 !! History : 578 !! ! 08-08 (A. Vidard) 579 !!----------------------------------------------------------------------- 580 !! * Modules used 581 !! * Arguments 582 INTEGER, INTENT(IN) :: & 583 & kumadt ! Output unit 584 585 INTEGER :: & 586 & ji, & ! dummy loop indices 587 & jj, & 588 & jk, & 589 & kt, & 590 & jn 591 592 !! * Local declarations 593 REAL(KIND=wp), DIMENSION(:,:,:,:), POINTER :: & 594 & zts, & ! Direct field : temperature 595 & zts_tlin, & ! Tangent input: temperature 596 & zts_adout, & ! Adjoint output: temperature 597 & zats ! 3D random field for t 598 599 REAL(KIND=wp), DIMENSION(:,:,:), POINTER :: & 600 & zgtu_tlout, & ! Tangent output: horizontal gradient 601 & zgtv_tlout, & ! Tangent output: horizontal gradient 602 & zrd_adout, & ! Adjoint output: 603 & zar, & ! 3D random field for rd 604 & zrd_tlin, & ! Tangent input: 605 & zgtu_adin, & ! Adjoint input : horizontal gradient 606 & zgtv_adin ! Adjoint input : horizontal gradient 607 608 REAL(KIND=wp), DIMENSION(:,:), POINTER :: & 609 & zgru_tlout, & ! Tangent output: horizontal gradient 610 & zgrv_tlout, & ! Tangent output: horizontal gradient 611 & zgru_adin, & ! Adjoint input : horizontal gradient 612 & zgrv_adin ! Adjoint input : horizontal gradient 613 614 REAL(KIND=wp) :: & 615 ! random field standard deviation for: 616 & zsp1, & ! scalar product involving the tangent routine 617 & zsp1_1, & ! scalar product components 618 & zsp1_2, & 619 & zsp1_3, & ! scalar product components 620 & zsp1_4, & 621 & zsp1_5, & ! scalar product components 622 & zsp1_6, & 623 & zsp2, & ! scalar product involving the adjoint routine 624 & zsp2_1, & ! scalar product components 625 & zsp2_2, & 626 & zsp2_3 627 CHARACTER (LEN=14) :: & 628 & cl_name 629 630 kt = nit000 631 ! Allocate memory 632 CALL wrk_alloc(jpi, jpj, jpk, jpts, zts, zts_tlin,zts_adout, zats ) 633 634 CALL wrk_alloc(jpi, jpj, jpts, zgtu_tlout, zgtv_tlout, zgtu_adin, zgtv_adin ) 635 636 CALL wrk_alloc(jpi, jpj, jpk , zrd_adout, zrd_tlin, zar ) 637 638 CALL wrk_alloc(jpi, jpj, zgru_tlout, zgrv_tlout, zgru_adin, zgrv_adin ) 639 640 ! Initialize random field standard deviationsthe reference state 641 zts = tsn(:,:,:,:) 642 643 !============================================================= 644 ! 1) dx = ( T ) and dy = ( T ) 645 !============================================================= 646 647 !-------------------------------------------------------------------- 648 ! Reset the tangent and adjoint variables 649 !-------------------------------------------------------------------- 650 zts_tlin(:,:,:,:) = 0.0_wp 651 zrd_tlin(:,:,:) = 0.0_wp 652 zts_adout(:,:,:,:) = 0.0_wp 653 zrd_adout(:,:,:) = 0.0_wp 654 zgtu_tlout(:,:,:) = 0.0_wp 655 zgtv_tlout(:,:,:) = 0.0_wp 656 zgru_tlout(:,:) = 0.0_wp 657 zgrv_tlout(:,:) = 0.0_wp 658 zgtu_adin(:,:,:) = 0.0_wp 659 zgtv_adin(:,:,:) = 0.0_wp 660 zgru_adin(:,:) = 0.0_wp 661 zgrv_adin(:,:) = 0.0_wp 662 663 !-------------------------------------------------------------------- 664 ! Initialize the tangent input with random noise: dx 665 !-------------------------------------------------------------------- 666 DO jn = 1, jpts 667 CALL grid_random( zats(:,:,:,jn), 'T', 0.0_wp, stdt ) 668 END DO 669 CALL grid_random( zar, 'T', 0.0_wp, stdr ) 670 671 DO jj = nldj, nlej 672 DO ji = nldi, nlei 673 zts_tlin(ji,jj,:,:) = zats(ji,jj,:,:) 674 zrd_tlin(ji,jj,:) = zar (ji,jj,:) 675 END DO 676 END DO 677 678 CALL zps_hde_tan ( nit000, jpts, zts, & 679 & zts_tlin , zrd_tlin , & 680 & zgtu_tlout, zgru_tlout, & 681 & zgtv_tlout, zgrv_tlout ) 682 683 DO jn = 1, jpts 684 DO jj = nldj, nlej 685 DO ji = nldi, nlei 686 jk = mbku(ji,jj) 687 zgtu_adin(ji,jj,jn) = zgtu_tlout(ji,jj,jn) & 688 & * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) & 689 & * umask(ji,jj,jk) 690 jk = mbkv(ji,jj) 691 zgtv_adin(ji,jj,jn) = zgtv_tlout(ji,jj,jn) & 692 & * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) & 693 & * vmask(ji,jj,jk) 694 END DO 695 END DO 696 END DO 697 DO jj = nldj, nlej 698 DO ji = nldi, nlei 699 jk = mbku(ji,jj) 700 zgru_adin(ji,jj) = zgru_tlout(ji,jj) & 701 & * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) & 702 & * umask(ji,jj,jk) 703 jk = mbkv(ji,jj) 704 zgrv_adin(ji,jj) = zgrv_tlout(ji,jj) & 705 & * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) & 706 & * vmask(ji,jj,jk) 707 END DO 708 END DO 722 709 !!-------------------------------------------------------------------- 723 710 !! Compute the scalar product: ( L dx )^T W dy 724 711 !!-------------------------------------------------------------------- 725 712 726 !zsp1_1 = DOT_PRODUCT( zgtu_adin(:,:,jp_tem), zgtu_tlout(:,:,jp_tem) ) 727 !zsp1_2 = DOT_PRODUCT( zgtu_adin(:,:,jp_sal), zgtu_tlout(:,:,jp_sal) ) 728 !zsp1_3 = DOT_PRODUCT( zgru_adin, zgru_tlout ) 729 !zsp1_4 = DOT_PRODUCT( zgtv_adin(:,:,jp_tem), zgtv_tlout(:,:,jp_tem) ) 730 !zsp1_5 = DOT_PRODUCT( zgtv_adin(:,:,jp_sal), zgtv_tlout(:,:,jp_sal) ) 731 !zsp1_6 = DOT_PRODUCT( zgrv_adin, zgrv_tlout ) 732 !zsp1 = zsp1_1 + zsp1_2 + zsp1_3 + zsp1_4 + zsp1_5 + zsp1_6 733 734 735 !!-------------------------------------------------------------------- 736 !! Call the adjoint routine: dx^* = L^T dy^* 737 !!-------------------------------------------------------------------- 738 !CALL zps_hde_adj ( kt, jpts, zts, gtsu, gtsv , & 739 !& tsa_ad, , & 740 !& zgtu_adin , zgru_adin , & 741 !& zgtv_adin , zgrv_adin ) 742 743 !zsp2_1 = DOT_PRODUCT( zts_tlin(:,:,:,jp_tem), zts_adout(:,:,:,jp_tem) ) 744 !zsp2_2 = DOT_PRODUCT( zts_tlin(:,:,:,jp_sal), zts_adout(:,:,:,jp_sal) ) 745 !zsp2_3 = DOT_PRODUCT( zrd_tlin , zrd_adout ) 746 !zsp2 = zsp2_1 + zsp2_2 + zsp2_3 747 748 !! Compare the scalar products 749 750 !cl_name = 'zps_hde_adj ' 751 !CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 ) 752 753 !! Deallocate memory 754 !DEALLOCATE( & 755 !& zts, & 756 !& zts_tlin, & 757 !& zrd_tlin, & 758 !& zts_adout, & 759 !& zrd_adout, & 760 !& zar, & 761 !& zats, & 762 !& zgtu_tlout, & 763 !& zgtv_tlout, & 764 !& zgru_tlout, & 765 !& zgrv_tlout, & 766 !& zgtu_adin, & 767 !& zgtv_adin, & 768 !& zgru_adin, & 769 !& zgrv_adin & 770 !& ) 771 772 !END SUBROUTINE zps_hde_adj_tst 713 zsp1_1 = DOT_PRODUCT( zgtu_adin(:,:,jp_tem), zgtu_tlout(:,:,jp_tem) ) 714 zsp1_2 = DOT_PRODUCT( zgtu_adin(:,:,jp_sal), zgtu_tlout(:,:,jp_sal) ) 715 zsp1_3 = DOT_PRODUCT( zgru_adin, zgru_tlout ) 716 zsp1_4 = DOT_PRODUCT( zgtv_adin(:,:,jp_tem), zgtv_tlout(:,:,jp_tem) ) 717 zsp1_5 = DOT_PRODUCT( zgtv_adin(:,:,jp_sal), zgtv_tlout(:,:,jp_sal) ) 718 zsp1_6 = DOT_PRODUCT( zgrv_adin, zgrv_tlout ) 719 zsp1 = zsp1_1 + zsp1_2 + zsp1_3 + zsp1_4 + zsp1_5 + zsp1_6 720 721 722 !-------------------------------------------------------------------- 723 ! Call the adjoint routine: dx^* = L^T dy^* 724 !-------------------------------------------------------------------- 725 CALL zps_hde_adj ( kt, jpts, zts , & 726 & zts_adout , zrd_adout , & 727 & zgtu_adin , zgru_adin , & 728 & zgtv_adin , zgrv_adin ) 729 730 zsp2_1 = DOT_PRODUCT( zts_tlin(:,:,:,jp_tem), zts_adout(:,:,:,jp_tem) ) 731 zsp2_2 = DOT_PRODUCT( zts_tlin(:,:,:,jp_sal), zts_adout(:,:,:,jp_sal) ) 732 zsp2_3 = DOT_PRODUCT( zrd_tlin , zrd_adout ) 733 zsp2 = zsp2_1 + zsp2_2 + zsp2_3 734 735 ! Compare the scalar products 736 737 cl_name = 'zps_hde_adj ' 738 CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 ) 739 740 ! Deallocate memory 741 CALL wrk_dealloc(jpi, jpj, jpk, jpts, zts, zts_tlin,zts_adout, zats ) 742 743 CALL wrk_dealloc(jpi, jpj, jpts, zgtu_tlout, zgtv_tlout, zgtu_adin, zgtv_adin ) 744 745 CALL wrk_dealloc(jpi, jpj, jpk , zrd_adout, zrd_tlin, zar ) 746 747 CALL wrk_dealloc(jpi, jpj, zgru_tlout, zgrv_tlout, zgru_adin, zgrv_adin ) 748 749 750 END SUBROUTINE zps_hde_adj_tst 773 751 #endif 774 752 END MODULE zpshde_tam
Note: See TracChangeset
for help on using the changeset viewer.