- Timestamp:
- 02/09/15 20:18:34 (9 years ago)
- Location:
- codes/icosagcm/trunk
- Files:
-
- 16 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/trunk/make_icosa
r279 r327 16 16 arch="" 17 17 parallel="none" 18 physics="none" 18 19 CPP_KEY="CPP_NONE" 19 20 ICOSA_LIB="" … … 45 46 "-parallel") 46 47 parallel=$2 ; parallel_defined="TRUE"; shift ; shift ;; 48 49 "-physics") 50 physics=$2 ; physics_defined="TRUE"; shift ; shift ;; 47 51 48 52 "-job") … … 119 123 fi 120 124 125 if [[ "$physics" == "lmdz_generic" ]] 126 then 127 CPP_KEY="$CPP_KEY CPP_PHYSICS_LMDZ_GENERIC" 128 COMPIL_FFLAGS="$COMPIL_FFLAGS $PHYSICS_INCDIR" 129 ICOSA_LIB="$ICOSA_LIB $PHYSICS_LIBDIR $PHYSICS_LIB" 130 fi 131 121 132 if [[ "$with_xios_defined" == "TRUE" ]] 122 133 then … … 127 138 128 139 ICOSA_LIB="$ICOSA_LIB $NETCDF_LIBDIR $NETCDF_LIB $HDF5_LIBDIR $HDF5_LIB" 140 129 141 130 142 rm -f config.fcm -
codes/icosagcm/trunk/src/advect.f90
r252 r327 272 272 CALL trace_end("compute_gradq3d3") 273 273 274 gradq3d_out=gradq3d 274 275 DO k=1,3 276 DO l = ll_begin,ll_end 277 DO ij=ij_begin,ij_end 278 gradq3d_out(ij,l,k)=gradq3d(ij,l,k) 279 ENDDO 280 ENDDO 281 ENDDO 275 282 276 283 CONTAINS -
codes/icosagcm/trunk/src/advect_tracer.f90
r295 r327 95 95 96 96 CALL send_message(f_u,req_u) 97 CALL send_message(f_wfluxt,req_wfluxt) 98 CALL send_message(f_q,req_q) 99 CALL send_message(f_rhodz,req_rhodz) 100 97 101 CALL wait_message(req_u) 98 CALL send_message(f_wfluxt,req_wfluxt)99 102 CALL wait_message(req_wfluxt) 100 CALL send_message(f_q,req_q)101 103 CALL wait_message(req_q) 102 CALL send_message(f_rhodz,req_rhodz)103 104 CALL wait_message(req_rhodz) 104 105 ! CALL wait_message(req_u)106 ! CALL wait_message(req_wfluxt)107 ! CALL wait_message(req_q)108 ! CALL wait_message(req_rhodz)109 105 110 106 ! 1/2 vertical transport + back-trajectories … … 134 130 135 131 CALL send_message(f_cc,req_cc) 136 CALL wait_message(req_cc)137 132 138 133 139 134 ! horizontal transport - split in two to place transfer of gradq3d 140 !!$OMP BARRIER141 135 DO k = 1, nqtot 142 136 DO ind=1,ndomain … … 148 142 sqrt_leng=f_sqrt_leng(ind) 149 143 CALL compute_gradq3d(q(:,:,k),sqrt_leng,gradq3d,xyz_i,xyz_v) 144 150 145 END DO 151 146 152 147 CALL send_message(f_gradq3d,req_gradq3d) 153 !CALL wait_message(req_cc)148 CALL wait_message(req_cc) 154 149 CALL wait_message(req_gradq3d) 155 150 -
codes/icosagcm/trunk/src/caldyn_gcm.f90
r295 r327 201 201 IF(caldyn_eta==eta_mass) THEN 202 202 CALL send_message(f_ps,req_ps) 203 CALL wait_message(req_ps)204 203 ELSE 205 204 CALL send_message(f_mass,req_mass) 206 CALL wait_message(req_mass)207 205 END IF 208 206 207 CALL send_message(f_theta_rhodz,req_theta_rhodz) 209 208 CALL send_message(f_u,req_u) 210 CALL wait_message(req_u)211 CALL send_message(f_theta_rhodz,req_theta_rhodz)212 CALL wait_message(req_theta_rhodz)213 214 ! CALL wait_message(req_u)215 ! CALL wait_message(req_theta_rhodz)216 209 217 210 SELECT CASE(caldyn_conserv) … … 232 225 233 226 CALL send_message(f_qu,req_qu) 234 CALL wait_message(req_qu)227 ! CALL wait_message(req_qu) 235 228 236 229 DO ind=1,ndomain … … 364 357 365 358 IF(caldyn_eta==eta_mass) THEN 366 !CALL wait_message(req_ps)359 CALL wait_message(req_ps) 367 360 ELSE 368 !CALL wait_message(req_mass)361 CALL wait_message(req_mass) 369 362 END IF 370 !CALL wait_message(req_theta_rhodz)363 CALL wait_message(req_theta_rhodz) 371 364 372 365 IF(caldyn_eta==eta_mass) THEN ! Compute mass & theta 373 366 DO l = ll_begin,ll_end 374 !CALL test_message(req_u)367 CALL test_message(req_u) 375 368 !DIR$ SIMD 376 369 DO ij=ij_begin_ext,ij_end_ext … … 382 375 ELSE ! Compute only theta 383 376 DO l = ll_begin,ll_end 384 !CALL test_message(req_u)377 CALL test_message(req_u) 385 378 !DIR$ SIMD 386 379 DO ij=ij_begin_ext,ij_end_ext … … 390 383 END IF 391 384 392 !CALL wait_message(req_u)385 CALL wait_message(req_u) 393 386 394 387 !!! Compute shallow-water potential vorticity … … 446 439 INTEGER :: i,j,ij,l 447 440 REAL(rstd) :: p_ik, exner_ik 441 INTEGER,SAVE ::ij_omp_begin_ext, ij_omp_end_ext 442 !$OMP THREADPRIVATE(ij_omp_begin_ext, ij_omp_end_ext) 443 LOGICAL,SAVE :: first=.TRUE. 444 !$OMP THREADPRIVATE(first) 445 448 446 449 447 CALL trace_start("compute_geopot") 448 449 IF (first) THEN 450 first=.FALSE. 451 CALL distrib_level(ij_end_ext-ij_begin_ext+1,ij_omp_begin_ext,ij_omp_end_ext) 452 ij_omp_begin_ext=ij_omp_begin_ext+ij_begin_ext-1 453 ij_omp_end_ext=ij_omp_end_ext+ij_begin_ext-1 454 ENDIF 450 455 451 456 IF(caldyn_eta==eta_mass) THEN 452 457 453 458 !!! Compute exner function and geopotential 454 IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency455 459 DO l = 1,llm 456 ! !$OMP DO SCHEDULE(STATIC)457 460 !DIR$ SIMD 458 DO ij=ij_ begin_ext,ij_end_ext461 DO ij=ij_omp_begin_ext,ij_omp_end_ext 459 462 p_ik = ptop + mass_ak(l) + mass_bk(l)*ps(ij) ! FIXME : leave ps for the moment ; change ps to Ms later 460 463 ! p_ik = ptop + g*(mass_ak(l)+ mass_bk(l)*ps(i,j)) … … 465 468 ENDDO 466 469 ENDDO 467 ENDIF470 ! ENDIF 468 471 ELSE 469 472 ! We are using a Lagrangian vertical coordinate … … 474 477 IF(boussinesq) THEN ! compute only geopotential : pressure pk will be computed in compute_caldyn_horiz 475 478 ! specific volume 1 = dphi/g/rhodz 476 IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency 477 DO l = 1,llm 478 ! !$OMP DO SCHEDULE(STATIC) 479 !DIR$ SIMD 480 DO ij=ij_begin_ext,ij_end_ext 481 geopot(ij,l+1) = geopot(ij,l) + g*rhodz(ij,l) 482 ENDDO 479 ! IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency 480 DO l = 1,llm 481 !DIR$ SIMD 482 DO ij=ij_omp_begin_ext,ij_omp_end_ext 483 geopot(ij,l+1) = geopot(ij,l) + g*rhodz(ij,l) 483 484 ENDDO 484 END IF485 ENDDO 485 486 ELSE ! non-Boussinesq, compute geopotential and Exner pressure 486 487 ! uppermost layer 487 IF (is_omp_level_master) THEN ! no openMP on vertical due to dependency 488 488 489 !DIR$ SIMD 490 DO ij=ij_omp_begin_ext,ij_omp_end_ext 491 pk(ij,llm) = ptop + (.5*g)*rhodz(ij,llm) 492 END DO 493 ! other layers 494 DO l = llm-1, 1, -1 495 !DIR$ SIMD 496 DO ij=ij_omp_begin_ext,ij_omp_end_ext 497 pk(ij,l) = pk(ij,l+1) + (.5*g)*(rhodz(ij,l)+rhodz(ij,l+1)) 498 END DO 499 END DO 500 ! surface pressure (for diagnostics) 501 DO ij=ij_omp_begin_ext,ij_omp_end_ext 502 ps(ij) = pk(ij,1) + (.5*g)*rhodz(ij,1) 503 END DO 504 505 ! specific volume v = kappa*theta*pi/p = dphi/g/rhodz 506 DO l = 1,llm 489 507 !DIR$ SIMD 490 DO ij=ij_begin_ext,ij_end_ext 491 pk(ij,llm) = ptop + (.5*g)*rhodz(ij,llm) 492 END DO 493 ! other layers 494 DO l = llm-1, 1, -1 495 496 ! !$OMP DO SCHEDULE(STATIC) 497 !DIR$ SIMD 498 DO ij=ij_begin_ext,ij_end_ext 499 pk(ij,l) = pk(ij,l+1) + (.5*g)*(rhodz(ij,l)+rhodz(ij,l+1)) 500 END DO 501 END DO 502 ! surface pressure (for diagnostics) 503 DO ij=ij_begin_ext,ij_end_ext 504 ps(ij) = pk(ij,1) + (.5*g)*rhodz(ij,1) 505 END DO 506 507 ! specific volume v = kappa*theta*pi/p = dphi/g/rhodz 508 DO l = 1,llm 509 510 ! !$OMP DO SCHEDULE(STATIC) 511 !DIR$ SIMD 512 DO ij=ij_begin_ext,ij_end_ext 513 p_ik = pk(ij,l) 514 exner_ik = cpp * (p_ik/preff) ** kappa 515 geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik 516 pk(ij,l) = exner_ik 517 ENDDO 508 DO ij=ij_omp_begin_ext,ij_omp_end_ext 509 p_ik = pk(ij,l) 510 exner_ik = cpp * (p_ik/preff) ** kappa 511 geopot(ij,l+1) = geopot(ij,l) + (g*kappa)*rhodz(ij,l)*theta(ij,l)*exner_ik/p_ik 512 pk(ij,l) = exner_ik 518 513 ENDDO 519 END IF514 ENDDO 520 515 END IF 521 516 … … 562 557 DO l = ll_begin, ll_end 563 558 !!! Compute mass and theta fluxes 564 !IF (caldyn_conserv==energy) CALL test_message(req_qu)559 IF (caldyn_conserv==energy) CALL test_message(req_qu) 565 560 !DIR$ SIMD 566 561 DO ij=ij_begin_ext,ij_end_ext … … 602 597 CASE(energy) ! energy-conserving TRiSK 603 598 604 !CALL wait_message(req_qu)599 CALL wait_message(req_qu) 605 600 606 601 DO l=ll_begin,ll_end … … 796 791 REAL(rstd),INTENT(INOUT) :: wflux(iim*jjm,llm+1) ! vertical mass flux (kg/m2/s) 797 792 REAL(rstd),INTENT(INOUT) :: wwuu(iim*3*jjm,llm+1) 798 REAL(rstd),INTENT( OUT) :: du(iim*3*jjm,llm)799 REAL(rstd),INTENT( OUT) :: dtheta_rhodz(iim*jjm,llm)793 REAL(rstd),INTENT(INOUT) :: du(iim*3*jjm,llm) 794 REAL(rstd),INTENT(INOUT) :: dtheta_rhodz(iim*jjm,llm) 800 795 REAL(rstd),INTENT(OUT) :: dps(iim*jjm) 801 796 … … 803 798 INTEGER :: i,j,ij,l 804 799 REAL(rstd) :: p_ik, exner_ik 800 INTEGER,SAVE ::ij_omp_begin, ij_omp_end 801 !$OMP THREADPRIVATE(ij_omp_begin, ij_omp_end) 802 LOGICAL,SAVE :: first=.TRUE. 803 !$OMP THREADPRIVATE(first) 804 805 806 CALL trace_start("compute_geopot") 807 808 IF (first) THEN 809 first=.FALSE. 810 CALL distrib_level(ij_end-ij_begin+1,ij_omp_begin,ij_omp_end) 811 ij_omp_begin=ij_omp_begin+ij_begin-1 812 ij_omp_end=ij_omp_end+ij_begin-1 813 ENDIF 805 814 806 815 ! REAL(rstd) :: wwuu(iim*3*jjm,llm+1) ! tmp var, don't know why but gain 30% on the whole code in opemp … … 812 821 !$OMP BARRIER 813 822 !!! cumulate mass flux convergence from top to bottom 814 IF (is_omp_level_master) THEN823 ! IF (is_omp_level_master) THEN 815 824 DO l = llm-1, 1, -1 816 825 ! IF (caldyn_conserv==energy) CALL test_message(req_qu) … … 818 827 !!$OMP DO SCHEDULE(STATIC) 819 828 !DIR$ SIMD 820 DO ij=ij_ begin,ij_end829 DO ij=ij_omp_begin,ij_omp_end 821 830 convm(ij,l) = convm(ij,l) + convm(ij,l+1) 822 831 ENDDO 823 832 ENDDO 824 ENDIF833 ! ENDIF 825 834 826 835 !$OMP BARRIER -
codes/icosagcm/trunk/src/checksum.f90
r295 r327 14 14 TYPE(t_field), POINTER :: field(:) 15 15 INTEGER :: intval(2) 16 INTEGER :: ind,i,j,ij,l 16 INTEGER :: ind,i,j,ij,l,k 17 17 INTEGER :: tot_sum 18 18 … … 49 49 ENDDO 50 50 ENDDO 51 52 ELSE IF (field(ind)%ndim==4) THEN 53 54 DO k=1,size(field(ind)%rval4d,3) 55 DO l=1,size(field(ind)%rval4d,2) 56 DO j=jj_begin,jj_end 57 DO i=ii_begin,ii_end 58 ij=(j-1)*iim+i 59 IF (domain(ind)%own(i,j)) THEN 60 intval=transfer(field(ind)%rval4d(ij,l,k),intval,2) 61 tot_sum=tot_sum+intval(1)+intval(2) 62 ENDIF 63 ENDDO 64 ENDDO 65 ENDDO 66 ENDDO 51 67 52 68 ENDIF … … 56 72 57 73 !$OMP MASTER 58 PRINT*,"CheckSum Field ",field(1)%name,tot_sum74 PRINT*,"CheckSum Field : ",field(1)%name,tot_sum 59 75 !$OMP END MASTER 60 76 -
codes/icosagcm/trunk/src/domain.f90
r295 r327 402 402 SUBROUTINE assign_domain 403 403 USE mpipara 404 USE grid_param 404 405 IMPLICIT NONE 405 406 INTEGER :: nb_domain(0:mpi_size-1) 406 407 INTEGER :: rank, ind,ind_glo 408 INTEGER :: block_j,jb,i,j,nd_glo,n,nf 409 LOGICAL :: exit 407 410 408 411 DO rank=0,mpi_size-1 … … 415 418 ALLOCATE(domloc_glo_ind(ndomain)) 416 419 420 421 block_j=sqrt(nsplit_i*nsplit_j*nb_face*1./mpi_size) 422 exit=.FALSE. 423 jb=1 424 i=1 425 j=1 426 ind=1 427 nd_glo=0 417 428 rank=0 418 ind=0 419 DO ind_glo=1,ndomain_glo 420 ind=ind+1 429 DO WHILE (.NOT. exit) 430 431 IF (j==MIN(jb+block_j,nsplit_j*nb_face+1)) THEN 432 j=jb 433 i=i+1 434 ENDIF 435 436 IF (i>nsplit_i) THEN 437 i=1 438 jb=jb+block_j 439 j=jb 440 ENDIF 441 442 IF (ind>nb_domain(rank)) THEN 443 rank=rank+1 444 ind=1 445 ENDIF 446 ind_glo=(j-1)*nsplit_i+i 447 448 nd_glo=nd_glo+1 449 IF (nd_glo==ndomain_glo) THEN 450 451 exit=.TRUE. 452 IF (.NOT. (rank==mpi_size-1 .AND. ind==nb_domain(rank) )) THEN 453 PRINT *, "Distribution problem in assign_domain" 454 STOP 455 ENDIF 456 457 ENDIF 458 421 459 domglo_rank(ind_glo)=rank 422 460 domglo_loc_ind(ind_glo)=ind … … 426 464 ENDIF 427 465 428 IF (ind==nb_domain(rank)) THEN 429 rank=rank+1 430 ind=0 431 ENDIF 466 j=j+1 467 ind=ind+1 468 432 469 ENDDO 470 471 IF (is_mpi_master) THEN 472 473 ind_glo=0 474 WRITE(*,'') 475 PRINT*, ' MPI PROCESS DISTRIBUTION' 476 WRITE(*,'') 477 478 WRITE(*,"(' ')", ADVANCE='NO') 479 DO n=1,nsplit_i*7-1 480 WRITE(*,"('=')", ADVANCE='NO') 481 ENDDO 482 WRITE(*,'') 483 484 DO nf=1,nb_face 485 DO j=1,nsplit_j 486 IF (j>1) THEN 487 WRITE(*,"(' ')", ADVANCE='NO') 488 DO n=1,nsplit_i*7-1 489 WRITE(*,"('-')", ADVANCE='NO') 490 ENDDO 491 WRITE(*,'') 492 ENDIF 493 494 WRITE(*,"('|')", ADVANCE='NO') 495 DO i=1,nsplit_i 496 WRITE(*,"(' ',' ',' |')",ADVANCE='NO') 497 ENDDO 498 WRITE(*,'') 499 500 WRITE(*,"('|')", ADVANCE='NO') 501 DO i=1,nsplit_i 502 ind_glo=ind_glo+1 503 WRITE(*,"(' ',i4.4 ,' |')",ADVANCE='NO'),domglo_rank(ind_glo) 504 END DO 505 WRITE(*,'') 506 507 WRITE(*,"('|')", ADVANCE='NO') 508 DO i=1,nsplit_i 509 WRITE(*,"(' ',' ',' |')",ADVANCE='NO') 510 ENDDO 511 WRITE(*,'') 512 513 ENDDO 514 515 WRITE(*,"(' ')", ADVANCE='NO') 516 DO n=1,nsplit_i*7-1 517 WRITE(*,"('=')", ADVANCE='NO') 518 ENDDO 519 WRITE(*,'') 520 ENDDO 521 ENDIF 522 523 ! rank=0 524 ! ind=0 525 ! DO ind_glo=1,ndomain_glo 526 ! ind=ind+1 527 ! domglo_rank(ind_glo)=rank 528 ! domglo_loc_ind(ind_glo)=ind 529 ! IF (rank==mpi_rank) THEN 530 ! CALL copy_domain(domain_glo(ind_glo),domain(ind)) 531 ! domloc_glo_ind(ind)=ind_glo 532 ! ENDIF 533 ! 534 ! IF (ind==nb_domain(rank)) THEN 535 ! rank=rank+1 536 ! ind=0 537 ! ENDIF 538 ! ENDDO 433 539 434 540 !$OMP PARALLEL -
codes/icosagcm/trunk/src/etat0.f90
r325 r327 18 18 USE etat0_dcmip5_mod, ONLY : getin_etat0_dcmip5=>getin_etat0 19 19 USE etat0_williamson_mod, ONLY : getin_etat0_williamson=>getin_etat0 20 USE etat0_temperature_mod, ONLY: getin_etat0_temperature=>getin_etat0 20 21 ! Old interface 21 22 USE etat0_academic_mod, ONLY : etat0_academic=>etat0 … … 55 56 CASE ('isothermal') 56 57 CALL getin_etat0_isothermal 58 CALL etat0_collocated(f_phis,f_ps,f_mass,f_theta_rhodz,f_u, f_q) 59 CASE ('temperature_profile') 60 CALL getin_etat0_temperature 57 61 CALL etat0_collocated(f_phis,f_ps,f_mass,f_theta_rhodz,f_u, f_q) 58 62 CASE ('jablonowsky06') … … 161 165 USE etat0_dcmip5_mod, ONLY : compute_dcmip5 => compute_etat0 162 166 USE etat0_williamson_mod, ONLY : compute_w91_6 => compute_etat0 167 USE etat0_temperature_mod, ONLY: compute_etat0_temperature => compute_etat0 163 168 IMPLICIT NONE 164 169 REAL(rstd),INTENT(INOUT) :: ps(iim*jjm) … … 186 191 CALL compute_etat0_isothermal(iim*jjm, phis, ps, temp_i, ulon_i, ulat_i, q) 187 192 CALL compute_etat0_isothermal(3*iim*jjm, phis_e, ps_e, temp_e, ulon_e, ulat_e, q_e) 193 CASE ('temperature_profile') 194 CALL compute_etat0_temperature(iim*jjm, phis, ps, temp_i, ulon_i, ulat_i, q) 195 CALL compute_etat0_temperature(3*iim*jjm, phis_e, ps_e, temp_e, ulon_e, ulat_e, q_e) 188 196 CASE('jablonowsky06') 189 197 CALL compute_jablonowsky06(iim*jjm,lon_i,lat_i, phis, ps, temp_i, ulon_i, ulat_i) -
codes/icosagcm/trunk/src/icosa_gcm.f90
r295 r327 12 12 USE xios_mod 13 13 USE write_field 14 USE physics_mod 15 14 16 ! USE getin_mod 15 17 IMPLICIT NONE … … 42 44 CALL output_field_init 43 45 CALL init_timeloop 44 46 !$OMP END PARALLEL 47 48 CALL init_physics 49 50 !$OMP PARALLEL 45 51 CALL timeloop 46 52 CALL switch_omp_no_distrib_level -
codes/icosagcm/trunk/src/omp_para.F90
r324 r327 68 68 IMPLICIT NONE 69 69 LOGICAL, INTENT(IN) :: is_mpi_master 70 INTEGER :: ll_nb,i 70 INTEGER :: ll_nb,i,llb,lle 71 71 72 72 #ifdef CPP_USING_OMP … … 77 77 78 78 IF (using_openmp) THEN 79 !$OMP PARALLEL PRIVATE(ll_nb,i )79 !$OMP PARALLEL PRIVATE(ll_nb,i,llb,lle) 80 80 81 81 !$OMP MASTER … … 123 123 IF (omp_level_rank==omp_level_size-1) is_omp_last_level=.TRUE. 124 124 125 ll_end=0 125 lle=0 126 126 127 DO i=0,omp_level_rank 127 ll _begin=ll_end+1128 llb=lle+1 128 129 ll_nb=llm/omp_level_size 129 130 IF (MOD(llm,omp_level_size)>i) ll_nb=ll_nb+1 130 ll _end=ll_begin+ll_nb-1131 lle=llb+ll_nb-1 131 132 ENDDO 133 ll_begin=llb 134 ll_end=lle 132 135 133 136 ll_beginp1=ll_begin … … 160 163 161 164 ELSE 162 is_master=is_mpi_master163 165 omp_size=1 164 166 omp_level_size=1 … … 198 200 END SUBROUTINE init_omp_para 199 201 202 SUBROUTINE distrib_level(size,lbegin,lend) 203 IMPLICIT NONE 204 INTEGER,INTENT(IN) :: size 205 INTEGER,INTENT(OUT) :: lbegin 206 INTEGER,INTENT(OUT) :: lend 207 INTEGER :: div,rest 208 209 div=size/omp_level_size 210 rest=MOD(size,omp_level_size) 211 IF (omp_level_rank<rest) THEN 212 lbegin=(div+1)*omp_level_rank+1 213 lend=lbegin+div 214 ELSE 215 lbegin=(div+1)*rest + (omp_level_rank-rest)*div+1 216 lend=lbegin+div-1 217 ENDIF 218 END SUBROUTINE distrib_level 219 200 220 201 221 SUBROUTINE switch_omp_distrib_level -
codes/icosagcm/trunk/src/physics.f90
r325 r327 5 5 PRIVATE 6 6 7 INTEGER, PARAMETER :: phys_none=0, phys_HS94=1, phys_DCMIP=2, phys_ LB2012=37 INTEGER, PARAMETER :: phys_none=0, phys_HS94=1, phys_DCMIP=2, phys_lmdz_generic=3, phys_LB2012=4 8 8 9 9 INTEGER :: phys_type … … 26 26 USE physics_dcmip_mod, ONLY : init_physics_dcmip=>init_physics 27 27 USE etat0_venus_mod, ONLY : init_phys_venus=>init_physics 28 USE physics_lmdz_generic_mod, ONLY : init_physics_lmdz_generic=>init_physics 28 29 IMPLICIT NONE 29 30 … … 40 41 phys_type = phys_LB2012 41 42 CALL init_phys_venus 43 44 CASE ('phys_lmdz_generic') 45 CALL init_physics_lmdz_generic 46 phys_type=phys_lmdz_generic 42 47 CASE ('dcmip') 43 48 CALL allocate_field(f_dulon,field_t,type_real,llm, name='dulon') … … 50 55 CASE DEFAULT 51 56 IF(is_mpi_root) PRINT*, 'init_physics : Bad selector for variable physics <',& 52 TRIM(physics_type), '> options are <none>, <held_suarez>, <Lebonnois2012>, <dcmip> '57 TRIM(physics_type), '> options are <none>, <held_suarez>, <Lebonnois2012>, <dcmip>, <phys_lmdz_generic>' 53 58 STOP 54 59 END SELECT … … 57 62 END SUBROUTINE init_physics 58 63 59 SUBROUTINE physics(it,f_phis, f_ps, f_theta_rhodz, f_ue, f_q) 60 USE icosa 61 USE physics_interface_mod 64 SUBROUTINE physics(it,f_phis, f_ps, f_theta_rhodz, f_ue, f_wflux, f_q) 65 USE icosa 66 USE physics_interface_mod 67 USE physics_lmdz_generic_mod, ONLY : physics_lmdz_generic => physics 62 68 USE physics_dcmip_mod, ONLY : write_physics_dcmip => write_physics 63 69 USE etat0_heldsz_mod … … 69 75 TYPE(t_field),POINTER :: f_theta_rhodz(:) 70 76 TYPE(t_field),POINTER :: f_ue(:) 77 TYPE(t_field),POINTER :: f_wflux(:) 71 78 TYPE(t_field),POINTER :: f_q(:) 72 79 REAL(rstd),POINTER :: phis(:) … … 87 94 CASE(phys_HS94) 88 95 CALL held_suarez(f_ps,f_theta_rhodz,f_ue) 96 CASE (phys_lmdz_generic) 97 CALL physics_lmdz_generic(it ,f_phis, f_ps, f_theta_rhodz, f_ue, f_wflux, f_q) 89 98 CASE(phys_LB2012) 90 99 CALL phys_venus(f_ps,f_theta_rhodz,f_ue) -
codes/icosagcm/trunk/src/restart.f90
r266 r327 21 21 USE spherical_geom_mod 22 22 USE transfert_mod 23 USE disvert_mod 23 24 24 25 IMPLICIT NONE … … 38 39 CHARACTER(LEN=255) :: restart_file_name 39 40 INTEGER,PARAMETER :: nvert=6 40 INTEGER :: ncid, cellId, levId, edgeId, vertid, lonId, latId, bounds_lonId, bounds_latId, nqId 41 INTEGER :: ncid, cellId, levId, edgeId, vertid, lonId, latId, bounds_lonId, bounds_latId, nqId, levAxisId 41 42 INTEGER :: ind,ind_glo,i,j,k,nf 42 43 INTEGER :: status … … 91 92 status = NF90_DEF_VAR(ncid,'bounds_lon',NF90_DOUBLE,(/ vertId,cellId /),bounds_lonId) 92 93 status = NF90_DEF_VAR(ncid,'bounds_lat',NF90_DOUBLE,(/ vertId,cellId /),bounds_latId) 94 status = NF90_DEF_VAR(ncid,'lev',NF90_DOUBLE,(/ levId /),levAxisId) 95 status = NF90_PUT_ATT(ncid,levAxisId,"axis","Z") 96 status = NF90_PUT_ATT(ncid,levAxisId,"units","Pa") 97 status = NF90_PUT_ATT(ncid,levAxisId,"positive","down") 93 98 94 99 DO nf=1,nfield … … 97 102 IF (field(1)%ndim==2) THEN 98 103 status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(field(1)%name)),NF90_DOUBLE,(/ cellId /),fieldId(nf)) 104 status = NF90_PUT_ATT(ncid,FieldId(nf),"coordinates","lon lat") 99 105 ELSE IF (field(1)%ndim==3) THEN 100 106 status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(field(1)%name)),NF90_DOUBLE,(/ cellId, levId /),fieldId(nf)) 107 status = NF90_PUT_ATT(ncid,FieldId(nf),"coordinates","lev lon lat") 101 108 ELSE IF (field(1)%ndim==4) THEN 102 109 status = NF90_DEF_VAR(ncid,TRIM(ADJUSTL(field(1)%name)),NF90_DOUBLE,(/ cellId, levId,nqId /),fieldId(nf)) 110 status = NF90_PUT_ATT(ncid,FieldId(nf),"coordinates","nq lev lon lat") 103 111 ENDIF 104 status = NF90_PUT_ATT(ncid,FieldId(nf),"coordinates","lon lat")105 112 ELSE IF (field(1)%field_type==field_U) THEN 106 113 IF (field(1)%ndim==2) THEN … … 139 146 status=NF90_PUT_VAR(ncid,bounds_lonId,REAL(bounds_lon,r8),start=(/ 1,1 /),count=(/ nvert,ncell_glo /)) 140 147 status=NF90_PUT_VAR(ncid,bounds_latId,REAL(bounds_lat,r8),start=(/ 1,1 /),count=(/ nvert,ncell_glo /)) 148 status=NF90_PUT_VAR(ncid,levAxisId,REAL(presnivs,r8),start=(/ 1 /),count=(/ llm /)) 141 149 ENDIF 142 150 -
codes/icosagcm/trunk/src/time.f90
r295 r327 92 92 itau_physics=1 93 93 CALL getin('itau_physics',itau_physics) 94 if (itau_physics<=0) itau_physics = HUGE(itau_physics) 94 95 95 96 itau_check_conserv=HUGE(itau_check_conserv) -
codes/icosagcm/trunk/src/timeloop_gcm.f90
r326 r327 40 40 USE write_field 41 41 USE theta2theta_rhodz_mod 42 USE sponge_mod 42 43 IMPLICIT NONE 43 44 … … 111 112 f_psm2 => f_phis 112 113 END IF 114 CASE ('none') 115 nb_stage=0 113 116 114 117 CASE default … … 120 123 CALL init_theta2theta_rhodz 121 124 CALL init_dissip 125 CALL init_sponge 122 126 CALL init_caldyn 123 127 CALL init_guided 124 128 CALL init_advect_tracer 125 129 CALL init_check_conserve 126 CALL init_physics127 128 130 129 131 CALL etat0(f_ps,f_mass,f_phis,f_theta_rhodz,f_u, f_q) … … 144 146 USE icosa 145 147 USE dissip_gcm_mod 148 USE sponge_mod 146 149 USE disvert_mod 147 150 USE caldyn_mod … … 176 179 INTEGER :: rate_clock 177 180 INTEGER :: l 181 LOGICAL,SAVE :: first_physic=.TRUE. 182 !$OMP THREADPRIVATE(first_physic) 178 183 179 180 ! CALL write_etat0(f_ps, f_phis,f_theta_rhodz,f_u,f_q) 181 ! CALL read_start(f_ps,f_mass,f_phis,f_theta_rhodz,f_u,f_q) 182 ! CALL write_restart(f_ps,f_mass,f_phis,f_theta_rhodz,f_u,f_q) 183 184 184 185 CALL switch_omp_distrib_level 185 186 CALL caldyn_BC(f_phis, f_wflux) ! set constant values in first/last interfaces … … 228 229 CALL wait_message(req_q0) 229 230 230 ! CALL wait_message(req_ps0)231 ! CALL wait_message(req_mass0)232 ! CALL wait_message(req_theta_rhodz0)233 ! CALL wait_message(req_u0)234 ! CALL wait_message(req_q0)235 231 ENDIF 236 232 … … 264 260 265 261 IF (MOD(it,itau_dissip)==0) THEN 266 ! CALL send_message(f_ps,req_ps)267 ! CALL wait_message(req_ps)268 262 269 263 IF(caldyn_eta==eta_mass) THEN … … 278 272 END DO 279 273 ENDIF 280 ! CALL send_message(f_mass,req_mass) 281 ! CALL wait_message(req_mass) 274 282 275 CALL dissip(f_u,f_du,f_mass,f_phis, f_theta_rhodz,f_dtheta_rhodz) 283 276 284 ! CALL send_message(f_mass,req_mass)285 ! CALL wait_message(req_mass)286 277 CALL euler_scheme(.FALSE.) ! update only u, theta 278 IF (iflag_sponge > 0) THEN 279 CALL sponge(f_u,f_du,f_theta_rhodz,f_dtheta_rhodz) 280 CALL euler_scheme(.FALSE.) ! update only u, theta 281 ENDIF 287 282 END IF 288 283 … … 308 303 END IF 309 304 310 CALL physics(it,f_phis, f_ps, f_theta_rhodz, f_u, f_q)311 305 312 306 IF (MOD(it,itau_check_conserv)==0) THEN 313 CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it) 307 CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it) 308 ENDIF 309 310 IF (MOD(it,itau_physics)==0) THEN 311 CALL physics(it,f_phis, f_ps, f_theta_rhodz, f_u, f_wflux, f_q) 312 313 !$OMP MASTER 314 IF (first_physic) CALL SYSTEM_CLOCK(start_clock) 315 !$OMP END MASTER 316 first_physic=.FALSE. 314 317 ENDIF 315 318 … … 354 357 ELSE ! update mass 355 358 mass=f_mass(ind) ; dmass=f_dmass(ind) ; 356 DO l= 1,llm359 DO l=ll_begin,ll_end 357 360 !$SIMD 358 361 DO ij=ij_begin,ij_end … … 606 609 END SUBROUTINE accumulate_fluxes 607 610 608 ! FUNCTION maxval_i(p)609 ! USE icosa610 ! IMPLICIT NONE611 ! REAL(rstd), DIMENSION(iim*jjm) :: p612 ! REAL(rstd) :: maxval_i613 ! INTEGER :: j, ij614 !615 ! maxval_i=p((jj_begin-1)*iim+ii_begin)616 !617 ! DO j=jj_begin-1,jj_end+1618 ! ij=(j-1)*iim619 ! maxval_i = MAX(maxval_i, MAXVAL(p(ij+ii_begin:ij+ii_end)))620 ! END DO621 ! END FUNCTION maxval_i622 623 ! FUNCTION maxval_ik(p)624 ! USE icosa625 ! IMPLICIT NONE626 ! REAL(rstd) :: p(iim*jjm, llm)627 ! REAL(rstd) :: maxval_ik(llm)628 ! INTEGER :: l,j, ij629 !630 ! DO l=1,llm631 ! maxval_ik(l)=p((jj_begin-1)*iim+ii_begin,l)632 ! DO j=jj_begin-1,jj_end+1633 ! ij=(j-1)*iim634 ! maxval_ik(l) = MAX(maxval_ik(l), MAXVAL(p(ij+ii_begin:ij+ii_end,l)))635 ! END DO636 ! END DO637 ! END FUNCTION maxval_ik638 639 611 END MODULE timeloop_gcm_mod -
codes/icosagcm/trunk/src/trace.F90
r186 r327 43 43 44 44 END SUBROUTINE trace_off 45 45 46 46 47 SUBROUTINE trace_start(name) … … 59 60 END SUBROUTINE trace_start 60 61 62 61 63 SUBROUTINE trace_end(name) 62 64 IMPLICIT NONE … … 74 76 75 77 END SUBROUTINE trace_end 78 79 76 80 77 81 SUBROUTINE trace_start2(name) -
codes/icosagcm/trunk/src/transfert_mpi.f90
r295 r327 159 159 DO j=jj_begin,jj_end 160 160 CALL request_add_point(ind,ii_end+1,j,req_e1_scal,left) 161 ENDDO 162 DO j=jj_begin,jj_end 161 163 CALL request_add_point(ind,ii_end+1,j-1,req_e1_scal,lup) 162 164 ENDDO … … 169 171 DO j=jj_begin,jj_end 170 172 CALL request_add_point(ind,ii_begin-1,j,req_e1_scal,right) 173 ENDDO 174 DO j=jj_begin,jj_end 171 175 CALL request_add_point(ind,ii_begin-1,j+1,req_e1_scal,rdown) 172 176 ENDDO … … 213 217 DO j=jj_begin,jj_end 214 218 CALL request_add_point(ind,ii_end+1,j,req_e1_vect,left) 219 ENDDO 220 DO j=jj_begin,jj_end 215 221 CALL request_add_point(ind,ii_end+1,j-1,req_e1_vect,lup) 216 222 ENDDO … … 223 229 DO j=jj_begin,jj_end 224 230 CALL request_add_point(ind,ii_begin-1,j,req_e1_vect,right) 231 ENDDO 232 DO j=jj_begin,jj_end 225 233 CALL request_add_point(ind,ii_begin-1,j+1,req_e1_vect,rdown) 226 234 ENDDO … … 1096 1104 INTEGER,POINTER :: sign(:) 1097 1105 INTEGER :: offset,msize,rank 1098 1099 CALL trace_start("transfert_mpi") 1106 INTEGER :: lbegin, lend 1107 1108 ! CALL trace_start("send_message_mpi") 1100 1109 1101 1110 !$OMP BARRIER … … 1192 1201 1193 1202 DO ind=1,ndomain 1194 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1203 IF (.NOT. assigned_domain(ind) ) CYCLE 1195 1204 1196 1205 dim3=size(field(ind)%rval3d,2) 1206 CALL distrib_level(dim3,lbegin,lend) 1207 1197 1208 rval3d=>field(ind)%rval3d 1198 1209 req=>message%request(ind) … … 1205 1216 ireq=send%ireq 1206 1217 buffer_r=>message%buffers(ireq)%r 1207 offset=send%offset*dim3 1208 1209 DO d3=1,dim3 1218 1219 offset=send%offset*dim3 + (lbegin-1)*send%size 1220 1221 CALL trace_start("copy_to_buffer") 1222 1223 DO d3=lbegin,lend 1210 1224 DO n=1,send%size 1211 1225 buffer_r(n+offset)=rval3d(value(n),d3) … … 1213 1227 offset=offset+send%size 1214 1228 ENDDO 1215 1216 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1217 !$OMP CRITICAL 1218 CALL MPI_ISEND(buffer_r,send%size*dim3,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1219 !$OMP END CRITICAL 1220 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1221 CALL MPI_ISEND(buffer_r,send%size*dim3,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1229 CALL trace_end("copy_to_buffer") 1230 1231 IF (is_omp_level_master) THEN 1232 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1233 !$OMP CRITICAL 1234 CALL MPI_ISEND(buffer_r,send%size*dim3,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1235 !$OMP END CRITICAL 1236 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1237 CALL MPI_ISEND(buffer_r,send%size*dim3,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1238 ENDIF 1222 1239 ENDIF 1223 1240 ENDIF … … 1226 1243 1227 1244 DO ind=1,ndomain 1228 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1245 IF (.NOT. assigned_domain(ind) ) CYCLE 1229 1246 dim3=size(field(ind)%rval3d,2) 1247 CALL distrib_level(dim3,lbegin,lend) 1230 1248 rval3d=>field(ind)%rval3d 1231 1249 req=>message%request(ind) … … 1240 1258 sgn=>recv%sign 1241 1259 1242 DO n=1,recv%size 1243 rval3d(value(n),:)=src_rval3d(src_value(n),:)*sgn(n) 1260 CALL trace_start("copy_data") 1261 1262 DO d3=lbegin,lend 1263 DO n=1,recv%size 1264 rval3d(value(n),d3)=src_rval3d(src_value(n),d3)*sgn(n) 1265 ENDDO 1244 1266 ENDDO 1267 CALL trace_end("copy_data") 1245 1268 1246 1269 ELSE … … 1248 1271 buffer_r=>message%buffers(ireq)%r 1249 1272 1250 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1251 !$OMP CRITICAL 1252 CALL MPI_IRECV(buffer_r,recv%size*dim3,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1253 !$OMP END CRITICAL 1254 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1255 CALL MPI_IRECV(buffer_r,recv%size*dim3,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1273 IF (is_omp_level_master) THEN 1274 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1275 !$OMP CRITICAL 1276 CALL MPI_IRECV(buffer_r,recv%size*dim3,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1277 !$OMP END CRITICAL 1278 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1279 CALL MPI_IRECV(buffer_r,recv%size*dim3,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1280 ENDIF 1256 1281 ENDIF 1257 ENDIF 1282 ENDIF 1258 1283 ENDDO 1259 1284 … … 1263 1288 1264 1289 DO ind=1,ndomain 1265 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1290 IF (.NOT. assigned_domain(ind) ) CYCLE 1266 1291 1267 1292 dim3=size(field(ind)%rval4d,2) 1293 CALL distrib_level(dim3,lbegin,lend) 1268 1294 dim4=size(field(ind)%rval4d,3) 1269 1295 rval4d=>field(ind)%rval4d … … 1278 1304 ireq=send%ireq 1279 1305 buffer_r=>message%buffers(ireq)%r 1280 offset=send%offset*dim3*dim4 1281 1306 1307 CALL trace_start("copy_to_buffer") 1282 1308 DO d4=1,dim4 1283 DO d3=1,dim3 1309 offset=send%offset*dim3*dim4 + dim3*send%size*(d4-1) + (lbegin-1)*send%size 1310 1311 DO d3=lbegin,lend 1284 1312 DO n=1,send%size 1285 1313 buffer_r(n+offset)=rval4d(value(n),d3,d4) … … 1288 1316 ENDDO 1289 1317 ENDDO 1290 1291 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1292 !$OMP CRITICAL 1293 CALL MPI_ISEND(buffer_r,send%size*dim3*dim4,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1294 !$OMP END CRITICAL 1295 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1296 CALL MPI_ISEND(buffer_r,send%size*dim3*dim4,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1318 CALL trace_end("copy_to_buffer") 1319 1320 IF (is_omp_level_master) THEN 1321 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1322 !$OMP CRITICAL 1323 CALL MPI_ISEND(buffer_r,send%size*dim3*dim4,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1324 !$OMP END CRITICAL 1325 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1326 CALL MPI_ISEND(buffer_r,send%size*dim3*dim4,MPI_REAL8,send%rank,send%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1327 ENDIF 1297 1328 ENDIF 1298 1329 … … 1302 1333 1303 1334 DO ind=1,ndomain 1304 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1335 IF (.NOT. assigned_domain(ind) ) CYCLE 1305 1336 1306 1337 dim3=size(field(ind)%rval4d,2) 1338 CALL distrib_level(dim3,lbegin,lend) 1307 1339 dim4=size(field(ind)%rval4d,3) 1308 1340 rval4d=>field(ind)%rval4d … … 1317 1349 sgn=>recv%sign 1318 1350 1319 DO n=1,recv%size 1320 rval4d(value(n),:,:)=src_rval4d(src_value(n),:,:)*sgn(n) 1351 CALL trace_start("copy_data") 1352 DO d4=1,dim4 1353 DO d3=lbegin,lend 1354 DO n=1,recv%size 1355 rval4d(value(n),d3,d4)=src_rval4d(src_value(n),d3,d4)*sgn(n) 1356 ENDDO 1357 ENDDO 1321 1358 ENDDO 1359 1360 CALL trace_end("copy_data") 1322 1361 1323 1362 ELSE … … 1325 1364 ireq=recv%ireq 1326 1365 buffer_r=>message%buffers(ireq)%r 1327 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1328 !$OMP CRITICAL 1329 CALL MPI_IRECV(buffer_r,recv%size*dim3*dim4,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1330 !$OMP END CRITICAL 1331 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1332 CALL MPI_IRECV(buffer_r,recv%size*dim3*dim4,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1366 IF (is_omp_level_master) THEN 1367 IF (mpi_threading_mode==MPI_THREAD_SERIALIZED) THEN 1368 !$OMP CRITICAL 1369 CALL MPI_IRECV(buffer_r,recv%size*dim3*dim4,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1370 !$OMP END CRITICAL 1371 ELSE IF (mpi_threading_mode==MPI_THREAD_MULTIPLE) THEN 1372 CALL MPI_IRECV(buffer_r,recv%size*dim3*dim4,MPI_REAL8,recv%rank,recv%tag+1000*message%number,comm_icosa, message%mpi_req(ireq),ierr) 1373 ENDIF 1333 1374 ENDIF 1334 1335 1375 ENDIF 1336 1376 ENDDO … … 1362 1402 1363 1403 !$OMP BARRIER 1364 CALL trace_end("transfert_mpi")1404 ! CALL trace_end("send_message_mpi") 1365 1405 1366 1406 END SUBROUTINE send_message_mpi … … 1402 1442 INTEGER :: ireq,nreq 1403 1443 INTEGER :: ind,n,l,m,i 1404 INTEGER :: dim3,dim4,d3,d4 1444 INTEGER :: dim3,dim4,d3,d4,lbegin,lend 1405 1445 INTEGER :: offset 1406 1446 1407 1447 IF (.NOT. message%pending) RETURN 1408 1448 1409 CALL trace_start("transfert_mpi")1449 ! CALL trace_start("wait_message_mpi") 1410 1450 1411 1451 field=>message%field … … 1452 1492 1453 1493 DO ind=1,ndomain 1454 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1494 IF (.NOT. assigned_domain(ind) ) CYCLE 1455 1495 1456 1496 rval3d=>field(ind)%rval3d … … 1465 1505 1466 1506 dim3=size(rval3d,2) 1467 offset=recv%offset*dim3 1468 DO d3=1,dim3 1469 DO n=1,recv%size 1470 rval3d(value(n),d3)=buffer_r(n+offset)*sgn(n) 1507 1508 CALL distrib_level(dim3,lbegin,lend) 1509 offset=recv%offset*dim3 + (lbegin-1)*recv%size 1510 CALL trace_start("copy_from_buffer") 1511 1512 IF (req%vector) THEN 1513 DO d3=lbegin,lend 1514 DO n=1,recv%size 1515 rval3d(value(n),d3)=buffer_r(n+offset)*sgn(n) 1516 ENDDO 1517 offset=offset+recv%size 1471 1518 ENDDO 1472 offset=offset+recv%size 1473 ENDDO 1519 ELSE 1520 DO d3=lbegin,lend 1521 DO n=1,recv%size 1522 rval3d(value(n),d3)=buffer_r(n+offset) 1523 ENDDO 1524 offset=offset+recv%size 1525 ENDDO 1526 ENDIF 1527 1528 CALL trace_end("copy_from_buffer") 1474 1529 ENDIF 1475 1530 ENDDO … … 1485 1540 1486 1541 DO ind=1,ndomain 1487 IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE1542 IF (.NOT. assigned_domain(ind) ) CYCLE 1488 1543 1489 1544 rval4d=>field(ind)%rval4d … … 1498 1553 1499 1554 dim3=size(rval4d,2) 1555 CALL distrib_level(dim3,lbegin,lend) 1500 1556 dim4=size(rval4d,3) 1501 offset=recv%offset*dim3*dim41557 CALL trace_start("copy_from_buffer") 1502 1558 DO d4=1,dim4 1503 DO d3=1,dim3 1559 offset=recv%offset*dim3*dim4 + dim3*recv%size*(d4-1) + (lbegin-1)*recv%size 1560 DO d3=lbegin,lend 1504 1561 DO n=1,recv%size 1505 1562 rval4d(value(n),d3,d4)=buffer_r(n+offset)*sgn(n) … … 1508 1565 ENDDO 1509 1566 ENDDO 1567 CALL trace_end("copy_from_buffer") 1510 1568 ENDIF 1511 1569 ENDDO … … 1521 1579 !$OMP END MASTER 1522 1580 1523 CALL trace_end("transfert_mpi")1581 ! CALL trace_end("wait_message_mpi") 1524 1582 !$OMP BARRIER 1525 1583 … … 1728 1786 1729 1787 END SUBROUTINE scatter_field 1730 1731 1732 1788 1733 1789 SUBROUTINE trace_in 1734 1790 USE trace … … 1752 1808 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 1753 1809 1754 !! -- Les chaine de charact ère -- !!1810 !! -- Les chaine de charactï¿œre -- !! 1755 1811 1756 1812 SUBROUTINE bcast_mpi_c(var1) … … 1949 2005 IF (.NOT. using_mpi) RETURN 1950 2006 1951 CALL MPI_BCAST(Var,nb,MPI_REAL ,mpi_master,comm_icosa,ierr)2007 CALL MPI_BCAST(Var,nb,MPI_REAL8,mpi_master,comm_icosa,ierr) 1952 2008 1953 2009 END SUBROUTINE bcast_mpi_rgen -
codes/icosagcm/trunk/src/transfert_omp.f90
r295 r327 102 102 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 103 103 104 !! -- Les chaine de charact ère -- !!104 !! -- Les chaine de charactï¿œre -- !! 105 105 106 106 SUBROUTINE bcast_omp_c(var)
Note: See TracChangeset
for help on using the changeset viewer.