![(please configure the [header_logo] section in trac.ini)](/dynamico/chrome/site/your_project_logo.png){kind=logo}
Changeset 148
- Timestamp:
- 03/18/13 15:44:08 (11 years ago)
- Location:
- codes/icosagcm/trunk/src
- Files:
-
- 7 edited
-
advect.f90 (modified) (13 diffs)
-
advect_tracer.f90 (modified) (19 diffs)
-
caldyn_gcm_opt.f90 (modified) (1 diff)
-
dissip_gcm.f90 (modified) (18 diffs)
-
timeloop_gcm.f90 (modified) (12 diffs)
-
transfert.F90 (modified) (1 diff)
-
transfert_mpi.f90 (modified) (8 diffs)
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/src/advect.f90
r141 r148 11 11 !========================================================================== 12 12 13 SUBROUTINE init_advect(normal,tangent )13 SUBROUTINE init_advect(normal,tangent,one_over_sqrt_leng) 14 14 IMPLICIT NONE 15 15 REAL(rstd),INTENT(OUT) :: normal(3*iim*jjm,3) 16 16 REAL(rstd),INTENT(OUT) :: tangent(3*iim*jjm,3) 17 REAL(rstd),INTENT(OUT) :: one_over_sqrt_leng(iim*jjm) 17 18 18 19 INTEGER :: i,j,n … … 39 40 tangent(n+u_ldown,:)=xyz_v(n+z_down,:)-xyz_v(n+z_ldown,:) 40 41 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 42 47 ENDDO 43 48 … … 46 51 !======================================================================================= 47 52 48 SUBROUTINE compute_gradq3d(qi,gradq3d) 53 SUBROUTINE compute_gradq3d(qi,one_over_sqrt_leng,gradq3d) 54 USE trace 49 55 IMPLICIT NONE 50 56 REAL(rstd),INTENT(IN) :: qi(iim*jjm,llm) 57 REAL(rstd),INTENT(IN) :: one_over_sqrt_leng(iim*jjm) 51 58 REAL(rstd),INTENT(OUT) :: gradq3d(iim*jjm,llm,3) 52 59 REAL(rstd) :: maxq,minq,minq_c,maxq_c 53 REAL(rstd) :: alphamx,alphami,alpha ,maggrd ,leng60 REAL(rstd) :: alphamx,alphami,alpha ,maggrd 54 61 REAL(rstd) :: leng1,leng2 55 62 REAL(rstd) :: arr(2*iim*jjm) 63 REAL(rstd) :: ar(iim*jjm) 56 64 REAL(rstd) :: gradtri(2*iim*jjm,llm,3) 57 REAL(rstd) :: ar58 65 INTEGER :: i,j,n,k,ind,l 66 67 CALL trace_start("compute_gradq3d") 59 68 60 69 ! TODO : precompute ar, drop arr as output argument of gradq ? … … 67 76 DO i=ii_begin-1,ii_end+1 68 77 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)) 71 80 END DO 72 81 END DO 73 82 END DO 74 83 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 87 103 88 104 !============================================================================================= LIMITING … … 94 110 maggrd = dot_product(gradq3d(n,l,:),gradq3d(n,l,:)) 95 111 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) 101 119 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), & 102 120 qi(n+t_rdown,l),qi(n+t_ldown,l)) … … 112 130 END DO 113 131 END DO 132 133 CALL trace_end("compute_gradq3d") 114 134 END SUBROUTINE compute_gradq3d 115 135 116 136 ! Backward trajectories, for use with Miura approach 117 137 SUBROUTINE compute_backward_traj(normal,tangent,ue,tau, cc) 138 USE trace 118 139 IMPLICIT NONE 119 140 REAL(rstd),INTENT(IN) :: normal(3*iim*jjm,3) … … 125 146 REAL(rstd) :: v_e(3), up_e, qe, ed(3) 126 147 INTEGER :: i,j,n,l 148 149 CALL trace_start("compute_backward_traj") 127 150 128 151 ! TODO : compute normal displacement ue*tau as hfluxt / mass(upwind) then reconstruct tangential displacement … … 181 204 ENDDO 182 205 END DO 206 207 CALL trace_end("compute_backward_traj") 208 183 209 END SUBROUTINE compute_backward_traj 184 210 … … 186 212 ! Slope-limited van Leer approach with hexagons 187 213 SUBROUTINE compute_advect_horiz(update_mass,hfluxt,cc,gradq3d, mass,qi) 214 USE trace 188 215 IMPLICIT NONE 189 216 LOGICAL, INTENT(IN) :: update_mass … … 197 224 REAL(rstd) :: qflux(3*iim*jjm,llm) 198 225 INTEGER :: i,j,n,k,l 226 227 CALL trace_start("compute_advect_horiz") 199 228 200 229 ! evaluate tracer field at point cc using piecewise linear reconstruction … … 262 291 END DO 263 292 293 CALL trace_end("compute_advect_horiz") 264 294 CONTAINS 265 295 … … 275 305 276 306 !========================================================================== 277 PURE SUBROUTINE gradq(n0, n1,n2,n3,q,dq,det)278 IMPLICIT NONE 279 INTEGER, INTENT(IN) :: n0, n1,n2,n3280 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) 281 311 REAL(rstd), INTENT(OUT) :: dq(3), det 282 312 … … 293 323 A(3,1)=xyz_v(n3,1); A(3,2)= xyz_v(n3,2); A(3,3)=xyz_v(n3,3) 294 324 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) 297 327 dq(3) = 0.0 298 328 ! CALL DGESV(3,1,A,3,IPIV,dq(:),3,info) -
codes/icosagcm/trunk/src/advect_tracer.f90
r146 r148 8 8 TYPE(t_field),POINTER :: f_gradq3d(:) 9 9 TYPE(t_field),POINTER :: f_cc(:) ! starting point of backward-trajectory (Miura approach) 10 10 TYPE(t_field),POINTER :: f_one_over_sqrt_leng(:) 11 11 12 REAL(rstd), PARAMETER :: pente_max=2.0 ! for vlz 12 13 … … 19 20 REAL(rstd),POINTER :: tangent(:,:) 20 21 REAL(rstd),POINTER :: normal(:,:) 22 REAL(rstd),POINTER :: one_over_sqrt_leng(:) 21 23 INTEGER :: ind 22 24 … … 25 27 CALL allocate_field(f_gradq3d,field_t,type_real,llm,3, name='gradq3d') 26 28 CALL allocate_field(f_cc,field_u,type_real,llm,3, name='cc') 29 CALL allocate_field(f_one_over_sqrt_leng,field_t,type_real, name='one_over_sqrt_leng') 27 30 28 31 DO ind=1,ndomain … … 31 34 normal=f_normal(ind) 32 35 tangent=f_tangent(ind) 33 CALL init_advect(normal,tangent) 36 one_over_sqrt_leng=f_one_over_sqrt_leng(ind) 37 CALL init_advect(normal,tangent,one_over_sqrt_leng) 34 38 END DO 35 39 … … 49 53 TYPE(t_field),POINTER :: f_rhodz(:) ! mass field at beginning of macro time step 50 54 51 REAL(rstd),POINTER :: q(:,:,:), normal(:,:), tangent(:,:), gradq3d(:,:,:), cc(:,:,:)55 REAL(rstd),POINTER :: q(:,:,:), normal(:,:), tangent(:,:), one_over_sqrt_leng(:), gradq3d(:,:,:), cc(:,:,:) 52 56 REAL(rstd),POINTER :: hfluxt(:,:), wfluxt(:,:) 53 57 REAL(rstd),POINTER :: rhodz(:,:), u(:,:) … … 106 110 rhodz = f_rhodz(ind) 107 111 wfluxt = f_wfluxt(ind) 112 108 113 DO k = 1, nqtot 109 CALL vlz(k==nqtot,0.5, wfluxt,rhodz,q(:,:,k)) 110 END DO 114 CALL vlz(k==nqtot,0.5, wfluxt,rhodz,q(:,:,k),1) 115 END DO 116 111 117 CALL compute_backward_traj(tangent,normal,u,0.5*dt*itau_adv, cc) 112 118 END DO 113 119 114 CALL transfert_request(f_q,req_i1) ! necessary ?115 CALL transfert_request(f_rhodz,req_i1) ! necessary ?120 ! CALL transfert_request(f_q,req_i1) ! necessary ? 121 ! CALL transfert_request(f_rhodz,req_i1) ! necessary ? 116 122 117 123 ! horizontal transport - split in two to place transfer of gradq3d … … 123 129 q = f_q(ind) 124 130 gradq3d = f_gradq3d(ind) 125 CALL compute_gradq3d(q(:,:,k),gradq3d) 131 one_over_sqrt_leng=f_one_over_sqrt_leng(ind) 132 CALL compute_gradq3d(q(:,:,k),one_over_sqrt_leng,gradq3d) 126 133 END DO 127 134 128 135 CALL transfert_request(f_gradq3d,req_i1) 136 137 129 138 130 139 DO ind=1,ndomain … … 140 149 END DO 141 150 142 CALL transfert_request(f_q,req_i1) ! necessary ?143 CALL transfert_request(f_rhodz,req_i1) ! necessary ?151 ! CALL transfert_request(f_q,req_i1) ! necessary ? 152 ! CALL transfert_request(f_rhodz,req_i1) ! necessary ? 144 153 145 154 ! 1/2 vertical transport … … 151 160 wfluxt = f_wfluxt(ind) 152 161 DO k = 1,nqtot 153 CALL vlz(k==nqtot, 0.5,wfluxt,rhodz, q(:,:,k) )162 CALL vlz(k==nqtot, 0.5,wfluxt,rhodz, q(:,:,k),0) 154 163 END DO 155 164 END DO … … 159 168 END SUBROUTINE advect_tracer 160 169 161 SUBROUTINE vlz(update_mass, fac,wfluxt,mass, q )170 SUBROUTINE vlz(update_mass, fac,wfluxt,mass, q, halo) 162 171 ! 163 172 ! Auteurs: P.Le Van, F.Hourdin, F.Forget, T. Dubos … … 168 177 ! wfluxt >0 for upward transport 169 178 ! ******************************************************************** 179 USE trace 170 180 IMPLICIT NONE 171 181 LOGICAL, INTENT(IN) :: update_mass … … 173 183 REAL(rstd), INTENT(INOUT) :: mass(iim*jjm,llm) 174 184 REAL(rstd), INTENT(INOUT) :: q(iim*jjm,llm) 185 INTEGER, INTENT(IN) :: halo 175 186 176 187 REAL(rstd) :: dq(iim*jjm,llm), & ! increase of q … … 182 193 INTEGER :: i,ij,l,j 183 194 195 CALL trace_start("vlz") 196 184 197 ! finite difference of q 185 198 DO l=2,llm 186 DO j=jj_begin- 1,jj_end+1187 DO i=ii_begin- 1,ii_end+1199 DO j=jj_begin-halo,jj_end+halo 200 DO i=ii_begin-halo,ii_end+halo 188 201 ij=(j-1)*iim+i 189 202 dzqw(ij,l)=q(ij,l)-q(ij,l-1) … … 196 209 ! dzq = slope*dz, i.e. the reconstructed q varies by dzq inside level l 197 210 DO l=2,llm-1 198 DO j=jj_begin- 1,jj_end+1199 DO i=ii_begin- 1,ii_end+1211 DO j=jj_begin-halo,jj_end+halo 212 DO i=ii_begin-halo,ii_end+halo 200 213 ij=(j-1)*iim+i 201 214 IF(dzqw(ij,l)*dzqw(ij,l+1).gt.0.) THEN … … 211 224 212 225 ! 0 slope in top and bottom layers 213 DO j=jj_begin- 1,jj_end+1214 DO i=ii_begin- 1,ii_end+1226 DO j=jj_begin-halo,jj_end+halo 227 DO i=ii_begin-halo,ii_end+halo 215 228 ij=(j-1)*iim+i 216 229 dzq(ij,1)=0. … … 222 235 ! then amount of q leaving level l/l+1 = wq = w * qq 223 236 DO l = 1,llm-1 224 DO j=jj_begin- 1,jj_end+1225 DO i=ii_begin- 1,ii_end+1237 DO j=jj_begin-halo,jj_end+halo 238 DO i=ii_begin-halo,ii_end+halo 226 239 ij=(j-1)*iim+i 227 240 w = fac*wfluxt(ij,l+1) … … 238 251 END DO 239 252 ! wq = 0 at top and bottom 240 DO j=jj_begin- 1,jj_end+1241 DO i=ii_begin- 1,ii_end+1253 DO j=jj_begin-halo,jj_end+halo 254 DO i=ii_begin-halo,ii_end+halo 242 255 ij=(j-1)*iim+i 243 256 wq(ij,llm+1)=0. … … 248 261 ! update q, mass is updated only after all q's have been updated 249 262 DO l=1,llm 250 DO j=jj_begin- 1,jj_end+1251 DO i=ii_begin- 1,ii_end+1263 DO j=jj_begin-halo,jj_end+halo 264 DO i=ii_begin-halo,ii_end+halo 252 265 ij=(j-1)*iim+i 253 266 newmass = mass(ij,l) + fac*(wfluxt(ij,l)-wfluxt(ij,l+1)) … … 258 271 END DO 259 272 273 CALL trace_end("vlz") 274 260 275 END SUBROUTINE vlz 261 276 -
codes/icosagcm/trunk/src/caldyn_gcm_opt.f90
r147 r148 111 111 112 112 CALL trace_start("caldyn") 113 CALL transfert_request(f_phis,req_i1)114 113 CALL transfert_request(f_ps,req_i1) 115 114 CALL transfert_request(f_theta_rhodz,req_i1) -
codes/icosagcm/trunk/src/dissip_gcm.f90
r146 r148 27 27 REAL, save :: rayleigh_tau 28 28 29 INTEGER,SAVE :: idissip29 ! INTEGER,SAVE :: itau_dissip 30 30 REAL,SAVE :: dtdissip 31 31 … … 57 57 USE mpi_mod 58 58 USE mpipara 59 USE time_mod 59 60 60 61 IMPLICIT NONE … … 75 76 INTEGER :: l 76 77 CHARACTER(len=255) :: rayleigh_friction_key 78 REAL(rstd) :: mintau 77 79 78 80 … … 415 417 tau_divgrad(l) = tau_divgrad(l)/zvert 416 418 ENDDO 419 420 mintau=tau_graddiv(1) 421 DO l=1,llm 422 mintau=MIN(mintau,tau_graddiv(l)) 423 mintau=MIN(mintau,tau_gradrot(l)) 424 mintau=MIN(mintau,tau_divgrad(l)) 425 ENDDO 417 426 418 ! idissip=INT(MIN(tetagdiv,tetagrot)/(2.*dt))419 ! idissip=MAX(1,idissip)420 ! dtdissip=idissip*dt421 ! PRINT *,"idissip",idissip," dtdissip ",dtdissip427 itau_dissip=INT(mintau/dt) 428 itau_dissip=MAX(1,itau_dissip) 429 dtdissip=itau_dissip*dt 430 IF (is_mpi_root) PRINT *,"mintau ",mintau,"itau_dissip",itau_dissip," dtdissip ",dtdissip 422 431 423 432 END SUBROUTINE init_dissip … … 431 440 USE geopotential_mod 432 441 USE trace 442 USE time_mod 433 443 IMPLICIT NONE 434 444 TYPE(t_field),POINTER :: f_ue(:) … … 475 485 n=(j-1)*iim+i 476 486 477 due(n+u_right,l) = -0.5*( due_diss1(n+u_right,l)/tau_graddiv(l) + due_diss2(n+u_right,l)/tau_gradrot(l) )478 due(n+u_lup,l) = -0.5*( due_diss1(n+u_lup,l) /tau_graddiv(l) + due_diss2(n+u_lup,l) /tau_gradrot(l) )479 due(n+u_ldown,l) = -0.5*( due_diss1(n+u_ldown,l)/tau_graddiv(l) + due_diss2(n+u_ldown,l)/tau_gradrot(l) )480 481 dtheta_rhodz(n,l) = -0.5*dtheta_rhodz_diss(n,l)/tau_divgrad(l) 487 due(n+u_right,l) = -0.5*( due_diss1(n+u_right,l)/tau_graddiv(l) + due_diss2(n+u_right,l)/tau_gradrot(l))*itau_dissip 488 due(n+u_lup,l) = -0.5*( due_diss1(n+u_lup,l) /tau_graddiv(l) + due_diss2(n+u_lup,l) /tau_gradrot(l))*itau_dissip 489 due(n+u_ldown,l) = -0.5*( due_diss1(n+u_ldown,l)/tau_graddiv(l) + due_diss2(n+u_ldown,l)/tau_gradrot(l))*itau_dissip 490 491 dtheta_rhodz(n,l) = -0.5*dtheta_rhodz_diss(n,l)/tau_divgrad(l)*itau_dissip 482 492 ENDDO 483 493 ENDDO … … 736 746 REAL(rstd),INTENT(IN) :: ue(iim*3*jjm,ll) 737 747 REAL(rstd),INTENT(OUT) :: gradivu_e(iim*3*jjm,ll) 738 REAL(rstd) ,SAVE,ALLOCATABLE :: divu_i(:,:)748 REAL(rstd) :: divu_i(iim*jjm,ll) 739 749 740 750 INTEGER :: i,j,n,l 741 751 742 !$OMP BARRIER743 !$OMP MASTER744 ALLOCATE(divu_i(iim*jjm,ll))745 !$OMP END MASTER746 !$OMP BARRIER747 748 752 DO l=1,ll 749 !$OMP DO750 753 DO j=jj_begin,jj_end 751 754 DO i=ii_begin,ii_end … … 762 765 763 766 DO l=1,ll 764 !$OMP DO765 767 DO j=jj_begin,jj_end 766 768 DO i=ii_begin,ii_end … … 779 781 780 782 DO l=1,ll 781 !$OMP DO782 783 DO j=jj_begin,jj_end 783 784 DO i=ii_begin,ii_end … … 790 791 ENDDO 791 792 792 !$OMP BARRIER 793 !$OMP MASTER 794 DEALLOCATE(divu_i) 795 !$OMP END MASTER 796 !$OMP BARRIER 797 798 793 799 794 END SUBROUTINE compute_gradiv 800 795 … … 805 800 REAL(rstd),INTENT(IN) :: theta(iim*jjm,ll) 806 801 REAL(rstd),INTENT(OUT) :: divgrad_i(iim*jjm,ll) 807 REAL(rstd) ,SAVE,ALLOCATABLE :: grad_e(:,:)802 REAL(rstd) :: grad_e(3*iim*jjm,ll) 808 803 809 804 INTEGER :: i,j,n,l 810 805 811 !$OMP BARRIER 812 !$OMP MASTER 813 ALLOCATE(grad_e(3*iim*jjm,ll)) 814 !$OMP END MASTER 815 !$OMP BARRIER 816 806 817 807 DO l=1,ll 818 !$OMP DO819 808 DO j=jj_begin-1,jj_end+1 820 809 DO i=ii_begin-1,ii_end+1 … … 834 823 835 824 DO l=1,ll 836 !$OMP DO837 825 DO j=jj_begin,jj_end 838 826 DO i=ii_begin,ii_end … … 849 837 850 838 DO l=1,ll 851 !$OMP DO852 839 DO j=jj_begin,jj_end 853 840 DO i=ii_begin,ii_end … … 857 844 ENDDO 858 845 ENDDO 859 860 !$OMP BARRIER861 !$OMP MASTER862 DEALLOCATE(grad_e)863 !$OMP END MASTER864 !$OMP BARRIER865 846 866 847 END SUBROUTINE compute_divgrad … … 873 854 REAL(rstd),INTENT(IN) :: ue(iim*3*jjm,ll) 874 855 REAL(rstd),INTENT(OUT) :: gradrot_e(iim*3*jjm,ll) 875 REAL(rstd) ,SAVE,ALLOCATABLE :: rot_v(:,:)856 REAL(rstd) :: rot_v(2*iim*jjm,ll) 876 857 877 858 INTEGER :: i,j,n,l 878 879 !$OMP BARRIER880 !$OMP MASTER881 ALLOCATE(rot_v(2*iim*jjm,ll))882 !$OMP END MASTER883 !$OMP BARRIER884 859 885 860 DO l=1,ll 886 !$OMP DO887 861 DO j=jj_begin-1,jj_end+1 888 862 DO i=ii_begin-1,ii_end+1 … … 902 876 903 877 DO l=1,ll 904 !$OMP DO905 878 DO j=jj_begin,jj_end 906 879 DO i=ii_begin,ii_end … … 918 891 919 892 DO l=1,ll 920 !$OMP DO921 893 DO j=jj_begin,jj_end 922 894 DO i=ii_begin,ii_end … … 930 902 ENDDO 931 903 ENDDO 932 933 !$OMP BARRIER934 !$OMP MASTER935 DEALLOCATE(rot_v)936 !$OMP END MASTER937 !$OMP BARRIER938 904 939 905 END SUBROUTINE compute_gradrot -
codes/icosagcm/trunk/src/timeloop_gcm.f90
r146 r148 6 6 7 7 INTEGER, PARAMETER :: euler=1, rk4=2, mlf=3 8 INTEGER :: itau_sync=10 8 9 9 10 CONTAINS … … 19 20 USE mpipara 20 21 USE trace 22 USE transfert_mod 21 23 IMPLICIT NONE 22 24 TYPE(t_field),POINTER :: f_phis(:) … … 47 49 CHARACTER(len=255) :: scheme_name 48 50 LOGICAL :: fluxt_zero(ndomain) ! set to .TRUE. to start accumulating fluxes in time 51 LOGICAL, PARAMETER :: check=.FALSE. 49 52 ! INTEGER :: itaumax 50 53 ! REAL(rstd) ::write_period … … 132 135 CALL init_advect_tracer 133 136 CALL init_physics 134 137 138 135 139 CALL etat0(f_ps,f_phis,f_theta_rhodz,f_u, f_q) 136 140 CALL writefield("phis",f_phis,once=.TRUE.) … … 155 159 END DO 156 160 161 CALL transfert_request(f_phis,req_i0) 162 CALL transfert_request(f_phis,req_i1) 163 157 164 DO it=0,itaumax 158 165 IF (MOD(it,itau_sync)==0) THEN 166 CALL transfert_request(f_ps,req_i0) 167 CALL transfert_request(f_theta_rhodz,req_i0) 168 CALL transfert_request(f_u,req_e0_vect) 169 CALL transfert_request(f_q,req_i0) 170 ENDIF 171 159 172 IF (is_mpi_root) PRINT *,"It No :",It," t :",dt*It 160 173 IF (mod(it,itau_out)==0 ) THEN … … 188 201 END DO 189 202 190 CALL dissip(f_u,f_du,f_ps,f_phis, f_theta_rhodz,f_dtheta_rhodz) 191 CALL euler_scheme(.FALSE.) 203 IF (MOD(it+1,itau_dissip)==0) THEN 204 CALL dissip(f_u,f_du,f_ps,f_phis, f_theta_rhodz,f_dtheta_rhodz) 205 CALL euler_scheme(.FALSE.) 206 ENDIF 192 207 193 208 IF(MOD(it+1,itau_adv)==0) THEN 194 CALL transfert_request(f_wfluxt,req_i1) ! FIXME209 ! CALL transfert_request(f_wfluxt,req_i1) ! FIXME 195 210 ! CALL transfert_request(f_hfluxt,req_e1) ! FIXME 196 211 … … 199 214 200 215 ! FIXME : check that rhodz is consistent with ps 201 CALL transfert_request(f_rhodz,req_i1) 202 CALL transfert_request(f_ps,req_i1) 203 CALL transfert_request(f_dps,req_i1) ! FIXME 204 CALL transfert_request(f_wflux,req_i1) ! FIXME 205 DO ind=1,ndomain 206 CALL swap_dimensions(ind) 207 CALL swap_geometry(ind) 208 rhodz=f_rhodz(ind); ps=f_ps(ind); dps=f_dps(ind); 209 wflux=f_wflux(ind); wfluxt=f_wfluxt(ind) 210 CALL compute_rhodz(.FALSE., ps, rhodz) 211 END DO 216 IF (check) THEN 217 CALL transfert_request(f_rhodz,req_i1) 218 CALL transfert_request(f_ps,req_i1) 219 CALL transfert_request(f_dps,req_i1) ! FIXME 220 CALL transfert_request(f_wflux,req_i1) ! FIXME 221 DO ind=1,ndomain 222 CALL swap_dimensions(ind) 223 CALL swap_geometry(ind) 224 rhodz=f_rhodz(ind); ps=f_ps(ind); dps=f_dps(ind); 225 wflux=f_wflux(ind); wfluxt=f_wfluxt(ind) 226 CALL compute_rhodz(.FALSE., ps, rhodz) 227 END DO 228 ENDIF 212 229 213 230 END IF … … 222 239 LOGICAL :: with_dps 223 240 INTEGER :: ind 224 241 INTEGER :: i,j,ij,l 225 242 CALL trace_start("Euler_scheme") 226 243 … … 229 246 CALL swap_geometry(ind) 230 247 IF(with_dps) THEN 231 ps=f_ps(ind) ; dps=f_dps(ind) ; 232 ps(:)=ps(:)+dt*dps(:) 233 hflux=f_hflux(ind); hfluxt=f_hfluxt(ind) 234 wflux=f_wflux(ind); wfluxt=f_wfluxt(ind) 235 CALL accumulate_fluxes(hflux,wflux,hfluxt,wfluxt,dt,fluxt_zero(ind)) 248 ps=f_ps(ind) ; dps=f_dps(ind) ; 249 250 DO j=jj_begin,jj_end 251 DO i=ii_begin,ii_end 252 ij=(j-1)*iim+i 253 ps(ij)=ps(ij)+dt*dps(ij) 254 ENDDO 255 ENDDO 256 hflux=f_hflux(ind); hfluxt=f_hfluxt(ind) 257 wflux=f_wflux(ind); wfluxt=f_wfluxt(ind) 258 CALL accumulate_fluxes(hflux,wflux,hfluxt,wfluxt,dt,fluxt_zero(ind)) 236 259 END IF 260 237 261 u=f_u(ind) ; theta_rhodz=f_theta_rhodz(ind) 238 262 du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind) 239 u(:,:)=u(:,:)+dt*du(:,:) 240 theta_rhodz(:,:)=theta_rhodz(:,:)+dt*dtheta_rhodz(:,:) 263 264 DO l=1,llm 265 DO j=jj_begin,jj_end 266 DO i=ii_begin,ii_end 267 ij=(j-1)*iim+i 268 u(ij+u_right,l)=u(ij+u_right,l)+dt*du(ij+u_right,l) 269 u(ij+u_lup,l)=u(ij+u_lup,l)+dt*du(ij+u_lup,l) 270 u(ij+u_ldown,l)=u(ij+u_ldown,l)+dt*du(ij+u_ldown,l) 271 theta_rhodz(ij,l)=theta_rhodz(ij,l)+dt*dtheta_rhodz(ij,l) 272 ENDDO 273 ENDDO 274 ENDDO 241 275 ENDDO 242 276 … … 250 284 REAL(rstd), DIMENSION(4), PARAMETER :: coef = (/ .25, 1./3., .5, 1. /) 251 285 REAL(rstd) :: tau 286 INTEGER :: i,j,ij,l 252 287 253 288 CALL trace_start("RK_scheme") … … 263 298 264 299 IF (stage==1) THEN ! first stage : save model state in XXm1 265 psm1(:)=ps(:) ; um1(:,:)=u(:,:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) 300 301 DO j=jj_begin,jj_end 302 DO i=ii_begin,ii_end 303 ij=(j-1)*iim+i 304 psm1(ij)=ps(ij) 305 ENDDO 306 ENDDO 307 308 DO l=1,llm 309 DO j=jj_begin,jj_end 310 DO i=ii_begin,ii_end 311 ij=(j-1)*iim+i 312 um1(ij+u_right,l)=u(ij+u_right,l) 313 um1(ij+u_lup,l)=u(ij+u_lup,l) 314 um1(ij+u_ldown,l)=u(ij+u_ldown,l) 315 theta_rhodzm1(ij,l)=theta_rhodz(ij,l) 316 ENDDO 317 ENDDO 318 ENDDO 319 266 320 END IF 267 321 ! updates are of the form x1 := x0 + tau*f(x1) 268 ps(:)=psm1(:)+tau*dps(:) 269 u(:,:)=um1(:,:)+tau*du(:,:) 270 theta_rhodz(:,:)=theta_rhodzm1(:,:)+tau*dtheta_rhodz(:,:) 322 DO j=jj_begin,jj_end 323 DO i=ii_begin,ii_end 324 ij=(j-1)*iim+i 325 ps(ij)=psm1(ij)+tau*dps(ij) 326 ENDDO 327 ENDDO 328 329 DO l=1,llm 330 DO j=jj_begin,jj_end 331 DO i=ii_begin,ii_end 332 ij=(j-1)*iim+i 333 u(ij+u_right,l)=um1(ij+u_right,l)+tau*du(ij+u_right,l) 334 u(ij+u_lup,l)=um1(ij+u_lup,l)+tau*du(ij+u_lup,l) 335 u(ij+u_ldown,l)=um1(ij+u_ldown,l)+tau*du(ij+u_ldown,l) 336 theta_rhodz(ij,l)=theta_rhodzm1(ij,l)+tau*dtheta_rhodz(ij,l) 337 ENDDO 338 ENDDO 339 ENDDO 340 271 341 IF(stage==nb_stage) THEN ! accumulate mass fluxes at last stage 272 342 hflux=f_hflux(ind); hfluxt=f_hfluxt(ind) … … 450 520 REAL(rstd), INTENT(IN) :: tau 451 521 LOGICAL, INTENT(INOUT) :: fluxt_zero 522 INTEGER :: l,i,j,ij 523 452 524 IF(fluxt_zero) THEN 453 525 ! PRINT *, 'Accumulating fluxes (first)' 454 526 fluxt_zero=.FALSE. 455 hfluxt = tau*hflux 456 wfluxt = tau*wflux 527 DO l=1,llm 528 DO j=jj_begin-1,jj_end+1 529 DO i=ii_begin-1,ii_end+1 530 ij=(j-1)*iim+i 531 hfluxt(ij+u_right,l) = tau*hflux(ij+u_right,l) 532 hfluxt(ij+u_lup,l) = tau*hflux(ij+u_lup,l) 533 hfluxt(ij+u_ldown,l) = tau*hflux(ij+u_ldown,l) 534 ENDDO 535 ENDDO 536 ENDDO 537 538 DO l=1,llm+1 539 DO j=jj_begin,jj_end 540 DO i=ii_begin,ii_end 541 ij=(j-1)*iim+i 542 wfluxt(ij,l) = tau*wflux(ij,l) 543 ENDDO 544 ENDDO 545 ENDDO 457 546 ELSE 458 547 ! PRINT *, 'Accumulating fluxes (next)' 459 hfluxt = hfluxt + tau*hflux 460 wfluxt = wfluxt + tau*wflux 548 DO l=1,llm 549 DO j=jj_begin-1,jj_end+1 550 DO i=ii_begin-1,ii_end+1 551 ij=(j-1)*iim+i 552 hfluxt(ij+u_right,l) = hfluxt(ij+u_right,l)+tau*hflux(ij+u_right,l) 553 hfluxt(ij+u_lup,l) = hfluxt(ij+u_lup,l)+tau*hflux(ij+u_lup,l) 554 hfluxt(ij+u_ldown,l) = hfluxt(ij+u_ldown,l)+tau*hflux(ij+u_ldown,l) 555 ENDDO 556 ENDDO 557 ENDDO 558 559 DO l=1,llm+1 560 DO j=jj_begin,jj_end 561 DO i=ii_begin,ii_end 562 ij=(j-1)*iim+i 563 wfluxt(ij,l) = wfluxt(ij,l)+tau*wflux(ij,l) 564 ENDDO 565 ENDDO 566 ENDDO 567 461 568 END IF 462 569 END SUBROUTINE accumulate_fluxes -
codes/icosagcm/trunk/src/transfert.F90
r146 r148 3 3 #ifdef CPP_USING_MPI 4 4 USE transfert_mpi_mod, ONLY : init_transfert, transfert_request=>transfert_request_mpi, req_i1,req_e1_vect, & 5 req_e1_scal, request_add_point, create_request, gather_field5 req_e1_scal, req_i0, req_e0_vect, req_e0_scal, request_add_point, create_request, gather_field 6 6 #else 7 7 USE transfert_mpi_mod, ONLY : init_transfert, transfert_request, req_i1,req_e1_vect, & -
codes/icosagcm/trunk/src/transfert_mpi.f90
r146 r148 37 37 TYPE(t_request),POINTER :: req_e1_vect(:) 38 38 39 TYPE(t_request),POINTER :: req_i0(:) 40 TYPE(t_request),POINTER :: req_e0_scal(:) 41 TYPE(t_request),POINTER :: req_e0_vect(:) 42 39 43 40 44 CONTAINS … … 69 73 70 74 DO i=ii_begin,ii_end 71 CALL request_add_point(ind,i,jj_begin,req_i1)75 ! CALL request_add_point(ind,i,jj_begin,req_i1) 72 76 ENDDO 73 77 74 78 DO j=jj_begin,jj_end 75 CALL request_add_point(ind,ii_end,j,req_i1)79 ! CALL request_add_point(ind,ii_end,j,req_i1) 76 80 ENDDO 77 81 78 82 DO i=ii_begin,ii_end 79 CALL request_add_point(ind,i,jj_end,req_i1)83 ! CALL request_add_point(ind,i,jj_end,req_i1) 80 84 ENDDO 81 85 82 86 DO j=jj_begin,jj_end 83 CALL request_add_point(ind,ii_begin,j,req_i1)87 ! CALL request_add_point(ind,ii_begin,j,req_i1) 84 88 ENDDO 85 89 … … 87 91 88 92 CALL finalize_request(req_i1) 89 93 94 95 CALL create_request(field_t,req_i0) 96 97 DO ind=1,ndomain 98 CALL swap_dimensions(ind) 99 100 DO i=ii_begin,ii_end 101 CALL request_add_point(ind,i,jj_begin,req_i0) 102 ENDDO 103 104 DO j=jj_begin,jj_end 105 CALL request_add_point(ind,ii_end,j,req_i0) 106 ENDDO 107 108 DO i=ii_begin,ii_end 109 CALL request_add_point(ind,i,jj_end,req_i0) 110 ENDDO 111 112 DO j=jj_begin,jj_end 113 CALL request_add_point(ind,ii_begin,j,req_i0) 114 ENDDO 115 116 ENDDO 117 118 CALL finalize_request(req_i0) 119 120 90 121 CALL create_request(field_u,req_e1_scal) 91 122 DO ind=1,ndomain … … 112 143 113 144 DO i=ii_begin+1,ii_end-1 114 CALL request_add_point(ind,i,jj_begin,req_e1_scal,right)115 CALL request_add_point(ind,i,jj_end,req_e1_scal,right)145 ! CALL request_add_point(ind,i,jj_begin,req_e1_scal,right) 146 ! CALL request_add_point(ind,i,jj_end,req_e1_scal,right) 116 147 ENDDO 117 148 118 149 DO j=jj_begin+1,jj_end-1 119 CALL request_add_point(ind,ii_begin,j,req_e1_scal,rup)120 CALL request_add_point(ind,ii_end,j,req_e1_scal,rup)150 ! CALL request_add_point(ind,ii_begin,j,req_e1_scal,rup) 151 ! CALL request_add_point(ind,ii_end,j,req_e1_scal,rup) 121 152 ENDDO 122 153 123 CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_scal,left)124 CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_scal,ldown)125 CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_scal,left)126 CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_scal,ldown)154 ! CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_scal,left) 155 ! CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_scal,ldown) 156 ! CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_scal,left) 157 ! CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_scal,ldown) 127 158 128 159 ENDDO 129 160 130 161 CALL finalize_request(req_e1_scal) 162 163 164 CALL create_request(field_u,req_e0_scal) 165 DO ind=1,ndomain 166 CALL swap_dimensions(ind) 167 168 169 DO i=ii_begin+1,ii_end-1 170 CALL request_add_point(ind,i,jj_begin,req_e0_scal,right) 171 CALL request_add_point(ind,i,jj_end,req_e0_scal,right) 172 ENDDO 173 174 DO j=jj_begin+1,jj_end-1 175 CALL request_add_point(ind,ii_begin,j,req_e0_scal,rup) 176 CALL request_add_point(ind,ii_end,j,req_e0_scal,rup) 177 ENDDO 178 179 CALL request_add_point(ind,ii_begin+1,jj_begin,req_e0_scal,left) 180 CALL request_add_point(ind,ii_begin,jj_begin+1,req_e0_scal,ldown) 181 CALL request_add_point(ind,ii_begin+1,jj_end,req_e0_scal,left) 182 CALL request_add_point(ind,ii_end,jj_begin+1,req_e0_scal,ldown) 183 184 ENDDO 185 186 CALL finalize_request(req_e0_scal) 187 188 131 189 132 190 CALL create_request(field_u,req_e1_vect,.TRUE.) … … 154 212 155 213 DO i=ii_begin+1,ii_end-1 156 CALL request_add_point(ind,i,jj_begin,req_e1_vect,right)157 CALL request_add_point(ind,i,jj_end,req_e1_vect,right)214 ! CALL request_add_point(ind,i,jj_begin,req_e1_vect,right) 215 ! CALL request_add_point(ind,i,jj_end,req_e1_vect,right) 158 216 ENDDO 159 217 160 218 DO j=jj_begin+1,jj_end-1 161 CALL request_add_point(ind,ii_begin,j,req_e1_vect,rup)162 CALL request_add_point(ind,ii_end,j,req_e1_vect,rup)219 ! CALL request_add_point(ind,ii_begin,j,req_e1_vect,rup) 220 ! CALL request_add_point(ind,ii_end,j,req_e1_vect,rup) 163 221 ENDDO 164 222 165 CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_vect,left) 166 CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_vect,ldown) 167 CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_vect,left) 168 CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_vect,ldown) 169 223 ! CALL request_add_point(ind,ii_begin+1,jj_begin,req_e1_vect,left) 224 ! CALL request_add_point(ind,ii_begin,jj_begin+1,req_e1_vect,ldown) 225 ! CALL request_add_point(ind,ii_begin+1,jj_end,req_e1_vect,left) 226 ! CALL request_add_point(ind,ii_end,jj_begin+1,req_e1_vect,ldown) 170 227 171 228 ENDDO 172 229 173 230 CALL finalize_request(req_e1_vect) 231 232 233 CALL create_request(field_u,req_e0_vect,.TRUE.) 234 DO ind=1,ndomain 235 CALL swap_dimensions(ind) 236 237 DO i=ii_begin+1,ii_end-1 238 CALL request_add_point(ind,i,jj_begin,req_e0_vect,right) 239 CALL request_add_point(ind,i,jj_end,req_e0_vect,right) 240 ENDDO 241 242 DO j=jj_begin+1,jj_end-1 243 CALL request_add_point(ind,ii_begin,j,req_e0_vect,rup) 244 CALL request_add_point(ind,ii_end,j,req_e0_vect,rup) 245 ENDDO 246 247 CALL request_add_point(ind,ii_begin+1,jj_begin,req_e0_vect,left) 248 CALL request_add_point(ind,ii_begin,jj_begin+1,req_e0_vect,ldown) 249 CALL request_add_point(ind,ii_begin+1,jj_end,req_e0_vect,left) 250 CALL request_add_point(ind,ii_end,jj_begin+1,req_e0_vect,ldown) 251 252 ENDDO 253 254 CALL finalize_request(req_e0_vect) 255 174 256 175 257 END SUBROUTINE init_transfert … … 533 615 USE mpi_mod 534 616 USE mpipara 617 USE trace 535 618 IMPLICIT NONE 536 619 TYPE(t_field),POINTER :: field(:) … … 553 636 INTEGER :: dim3,dim4 554 637 638 CALL trace_start("transfert_mpi") 639 555 640 IF (field(1)%data_type==type_real) THEN 556 641 IF (field(1)%ndim==2) THEN … … 747 832 748 833 ENDIF 834 835 CALL trace_end("transfert_mpi") 749 836 750 837 END SUBROUTINE transfert_request_mpi
Note: See TracChangeset
for help on using the changeset viewer.
