Ignore:
Timestamp:
03/18/13 15:44:08 (11 years ago)
Author:
ymipsl
Message:

Various optimisations

  • dissipation is not called every timestep (similar way than LMDZ)
  • transfert size of halos has been reduced : need to synchronise redondant data between tiles at itau_sync timestep

YM

File:
1 edited

Legend:

Unmodified
Added
Removed
  • codes/icosagcm/trunk/src/advect.f90

    r141 r148  
    1111  !========================================================================== 
    1212 
    13   SUBROUTINE init_advect(normal,tangent) 
     13  SUBROUTINE init_advect(normal,tangent,one_over_sqrt_leng) 
    1414    IMPLICIT NONE 
    1515    REAL(rstd),INTENT(OUT) :: normal(3*iim*jjm,3) 
    1616    REAL(rstd),INTENT(OUT) :: tangent(3*iim*jjm,3) 
     17    REAL(rstd),INTENT(OUT) :: one_over_sqrt_leng(iim*jjm) 
    1718 
    1819    INTEGER :: i,j,n 
     
    3940          tangent(n+u_ldown,:)=xyz_v(n+z_down,:)-xyz_v(n+z_ldown,:) 
    4041          tangent(n+u_ldown,:)=tangent(n+u_ldown,:)/sqrt(sum(tangent(n+u_ldown,:)**2)+1e-50) 
    41        END DO 
     42 
     43          one_over_sqrt_leng(n) = 1./sqrt(max(sum((xyz_v(n+z_up,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_down,:) - xyz_i(n,:))**2), & 
     44                                          sum((xyz_v(n+z_rup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_rdown,:) - xyz_i(n,:))**2),  & 
     45                                          sum((xyz_v(n+z_lup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_ldown,:) - xyz_i(n,:))**2)) ) 
     46       ENDDO 
    4247    ENDDO 
    4348 
     
    4651  !======================================================================================= 
    4752 
    48   SUBROUTINE compute_gradq3d(qi,gradq3d) 
     53  SUBROUTINE compute_gradq3d(qi,one_over_sqrt_leng,gradq3d) 
     54    USE trace 
    4955    IMPLICIT NONE 
    5056    REAL(rstd),INTENT(IN)  :: qi(iim*jjm,llm) 
     57    REAL(rstd),INTENT(IN)  :: one_over_sqrt_leng(iim*jjm) 
    5158    REAL(rstd),INTENT(OUT) :: gradq3d(iim*jjm,llm,3)  
    5259    REAL(rstd) :: maxq,minq,minq_c,maxq_c  
    53     REAL(rstd) :: alphamx,alphami,alpha ,maggrd,leng 
     60    REAL(rstd) :: alphamx,alphami,alpha ,maggrd 
    5461    REAL(rstd) :: leng1,leng2 
    5562    REAL(rstd) :: arr(2*iim*jjm) 
     63    REAL(rstd) :: ar(iim*jjm) 
    5664    REAL(rstd) :: gradtri(2*iim*jjm,llm,3)  
    57     REAL(rstd) :: ar 
    5865    INTEGER :: i,j,n,k,ind,l 
     66 
     67    CALL trace_start("compute_gradq3d") 
    5968 
    6069    ! TODO : precompute ar, drop arr as output argument of gradq ? 
     
    6776          DO i=ii_begin-1,ii_end+1 
    6877             n=(j-1)*iim+i    
    69              CALL gradq(n,n+t_rup,n+t_lup,n+z_up,qi(:,l),gradtri(n+z_up,l,:),arr(n+z_up)) 
    70              CALL gradq(n,n+t_ldown,n+t_rdown,n+z_down,qi(:,l),gradtri(n+z_down,l,:),arr(n+z_down)) 
     78             CALL gradq(n,l,n+t_rup,n+t_lup,n+z_up,qi,gradtri(n+z_up,l,:),arr(n+z_up)) 
     79             CALL gradq(n,l,n+t_ldown,n+t_rdown,n+z_down,qi,gradtri(n+z_down,l,:),arr(n+z_down)) 
    7180          END DO 
    7281       END DO 
    7382    END DO 
    7483 
    75     !      Do l =1,llm 
    76     DO j=jj_begin,jj_end 
    77        DO i=ii_begin,ii_end 
    78           n=(j-1)*iim+i 
    79           gradq3d(n,:,:) = gradtri(n+z_up,:,:) + gradtri(n+z_down,:,:)   +   & 
    80                gradtri(n+z_rup,:,:) +  gradtri(n+z_ldown,:,:)   +   &  
    81                gradtri(n+z_lup,:,:)+    gradtri(n+z_rdown,:,:)   
    82           ar = arr(n+z_up)+arr(n+z_lup)+arr(n+z_ldown)+arr(n+z_down)+arr(n+z_rdown)+arr(n+z_rup) 
    83           gradq3d(n,:,:) = gradq3d(n,:,:)/(ar+1.e-50)  
    84        END DO 
    85     END DO 
    86     !    END DO 
     84    DO j=jj_begin-1,jj_end+1 
     85      DO i=ii_begin-1,ii_end+1 
     86         n=(j-1)*iim+i    
     87         ar(n) = arr(n+z_up)+arr(n+z_lup)+arr(n+z_ldown)+arr(n+z_down)+arr(n+z_rdown)+arr(n+z_rup)+1.e-50 
     88      ENDDO 
     89    ENDDO 
     90       
     91    DO k=1,3 
     92      DO l =1,llm 
     93        DO j=jj_begin,jj_end 
     94           DO i=ii_begin,ii_end 
     95              n=(j-1)*iim+i 
     96              gradq3d(n,l,k) = ( gradtri(n+z_up,l,k) + gradtri(n+z_down,l,k)   +   & 
     97                                 gradtri(n+z_rup,l,k) +  gradtri(n+z_ldown,l,k)   +   &  
     98                                 gradtri(n+z_lup,l,k)+    gradtri(n+z_rdown,l,k) ) / ar(n)   
     99           END DO 
     100        END DO 
     101      END DO 
     102    ENDDO 
    87103 
    88104    !============================================================================================= LIMITING  
     
    94110             maggrd =  dot_product(gradq3d(n,l,:),gradq3d(n,l,:)) 
    95111             maggrd = sqrt(maggrd)  
    96              leng = max(sum((xyz_v(n+z_up,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_down,:) - xyz_i(n,:))**2), & 
    97                   sum((xyz_v(n+z_rup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_rdown,:) - xyz_i(n,:))**2),  & 
    98                   sum((xyz_v(n+z_lup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_ldown,:) - xyz_i(n,:))**2)) 
    99              maxq_c = qi(n,l) + maggrd*sqrt(leng) 
    100              minq_c = qi(n,l) - maggrd*sqrt(leng) 
     112!             leng = max(sum((xyz_v(n+z_up,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_down,:) - xyz_i(n,:))**2), & 
     113!                  sum((xyz_v(n+z_rup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_rdown,:) - xyz_i(n,:))**2),  & 
     114!                  sum((xyz_v(n+z_lup,:) - xyz_i(n,:))**2),sum((xyz_v(n+z_ldown,:) - xyz_i(n,:))**2)) 
     115!             maxq_c = qi(n,l) + maggrd*sqrt(leng(n)) 
     116!             minq_c = qi(n,l) - maggrd*sqrt(leng(n)) 
     117             maxq_c = qi(n,l) + maggrd*one_over_sqrt_leng(n) 
     118             minq_c = qi(n,l) - maggrd*one_over_sqrt_leng(n) 
    101119             maxq = max(qi(n,l),qi(n+t_right,l),qi(n+t_lup,l),qi(n+t_rup,l),qi(n+t_left,l), & 
    102120                  qi(n+t_rdown,l),qi(n+t_ldown,l)) 
     
    112130       END DO 
    113131    END DO 
     132 
     133  CALL trace_end("compute_gradq3d") 
    114134  END SUBROUTINE compute_gradq3d 
    115135 
    116136  ! Backward trajectories, for use with Miura approach 
    117137  SUBROUTINE compute_backward_traj(normal,tangent,ue,tau, cc) 
     138    USE trace 
    118139    IMPLICIT NONE 
    119140    REAL(rstd),INTENT(IN)    :: normal(3*iim*jjm,3) 
     
    125146    REAL(rstd) :: v_e(3), up_e, qe, ed(3)     
    126147    INTEGER :: i,j,n,l 
     148 
     149    CALL trace_start("compute_backward_traj") 
    127150 
    128151    ! TODO : compute normal displacement ue*tau as hfluxt / mass(upwind) then reconstruct tangential displacement 
     
    181204       ENDDO 
    182205    END DO 
     206 
     207    CALL trace_end("compute_backward_traj") 
     208 
    183209  END SUBROUTINE compute_backward_traj 
    184210 
     
    186212  ! Slope-limited van Leer approach with hexagons 
    187213  SUBROUTINE compute_advect_horiz(update_mass,hfluxt,cc,gradq3d, mass,qi)  
     214    USE trace 
    188215    IMPLICIT NONE 
    189216    LOGICAL, INTENT(IN)       :: update_mass 
     
    197224    REAL(rstd) :: qflux(3*iim*jjm,llm) 
    198225    INTEGER :: i,j,n,k,l 
     226 
     227    CALL trace_start("compute_advect_horiz") 
    199228 
    200229    ! evaluate tracer field at point cc using piecewise linear reconstruction 
     
    262291    END DO 
    263292 
     293    CALL trace_end("compute_advect_horiz") 
    264294  CONTAINS 
    265295 
     
    275305 
    276306  !========================================================================== 
    277   PURE SUBROUTINE gradq(n0,n1,n2,n3,q,dq,det) 
    278     IMPLICIT NONE 
    279     INTEGER, INTENT(IN) :: n0,n1,n2,n3 
    280     REAL(rstd), INTENT(IN)     :: q(iim*jjm) 
     307  PURE SUBROUTINE gradq(n0,l,n1,n2,n3,q,dq,det) 
     308    IMPLICIT NONE 
     309    INTEGER, INTENT(IN) :: n0,l,n1,n2,n3 
     310    REAL(rstd), INTENT(IN)     :: q(iim*jjm,llm) 
    281311    REAL(rstd), INTENT(OUT)    :: dq(3), det 
    282312     
     
    293323    A(3,1)=xyz_v(n3,1);              A(3,2)= xyz_v(n3,2);            A(3,3)=xyz_v(n3,3) 
    294324 
    295     dq(1) = q(n1)-q(n0) 
    296     dq(2) = q(n2)-q(n0) 
     325    dq(1) = q(n1,l)-q(n0,l) 
     326    dq(2) = q(n2,l)-q(n0,l) 
    297327    dq(3) = 0.0 
    298328    !    CALL DGESV(3,1,A,3,IPIV,dq(:),3,info) 
Note: See TracChangeset for help on using the changeset viewer.