! $Id: leapfrog_loc.F90 5116 2024-07-24 12:54:37Z abarral $ SUBROUTINE leapfrog_loc(ucov0, vcov0, teta0, ps0, & masse0, phis0, q0, time_0) USE misc_mod USE parallel_lmdz USE times USE mod_hallo USE Bands USE strings_mod, ONLY: int2str USE Write_Field_p USE vampir USE lmdz_timer_filtre, ONLY: print_filtre_timer USE infotrac USE guide_loc_mod, ONLY: guide_main USE getparam USE control_mod USE lmdz_filtreg_p USE write_field_loc USE allocate_field_mod USE call_dissip_mod, ONLY: call_dissip USE lmdz_call_calfis, ONLY: call_calfis USE leapfrog_mod, ONLY: ucov, vcov, teta, ps, masse, phis, q, dq & , ucovm1, vcovm1, tetam1, massem1, psm1, p, pks, pk, pkf, flxw & , pbaru, pbarv, du, dv, dteta, phi, dp, w & , leapfrog_allocate, leapfrog_switch_caldyn, leapfrog_switch_dissip use exner_hyb_loc_m, ONLY: exner_hyb_loc use exner_milieu_loc_m, ONLY: exner_milieu_loc USE comconst_mod, ONLY: cpp, dtvr, ihf USE comvert_mod, ONLY: ap, bp, pressure_exner USE logic_mod, ONLY: iflag_phys, ok_guide, forward, leapf, apphys, & statcl, conser, apdiss, purmats, ok_strato USE temps_mod, ONLY: itaufin, jD_ref, jH_ref, day_ini, & day_ref, start_time, dt USE mod_xios_dyn3dmem, ONLY: dyn3d_ctx_handle USE lmdz_xios, ONLY: xios_update_calendar, & xios_set_current_context, & using_xios USE lmdz_cppkeys_wrapper, ONLY: CPPKEY_DEBUGIO USE lmdz_description, ONLY: descript IMPLICIT NONE ! ...... Version du 10/01/98 .......... ! avec coordonnees verticales hybrides ! avec nouveaux operat. dissipation * ( gradiv2,divgrad2,nxgraro2 ) !======================================================================= ! ! Auteur: P. Le Van /L. Fairhead/F.Hourdin ! ------- ! ! Objet: ! ------ ! ! GCM LMD nouvelle grille ! !======================================================================= ! ! ... Dans inigeom , nouveaux calculs pour les elongations cu , cv ! et possibilite d'appeler une fonction f(y) a derivee tangente ! hyperbolique a la place de la fonction a derivee sinusoidale. ! ... Possibilite de choisir le shema pour l'advection de ! q , en modifiant iadv dans traceur.def (10/02) . ! ! Pour Van-Leer + Vapeur d'eau saturee, iadv(1)=4. (F.Codron,10/99) ! Pour Van-Leer iadv=10 ! !----------------------------------------------------------------------- ! Declarations: ! ------------- include "dimensions.h" include "paramet.h" include "comdissnew.h" include "comgeom.h" include "iniprint.h" include "academic.h" REAL, INTENT(IN) :: time_0 ! not used ! dynamical variables: REAL, INTENT(IN) :: ucov0(ijb_u:ije_u, llm) ! zonal covariant wind REAL, INTENT(IN) :: vcov0(ijb_v:ije_v, llm) ! meridional covariant wind REAL, INTENT(IN) :: teta0(ijb_u:ije_u, llm) ! potential temperature REAL, INTENT(IN) :: q0(ijb_u:ije_u, llm, nqtot) ! advected tracers REAL, INTENT(IN) :: ps0(ijb_u:ije_u) ! surface pressure (Pa) REAL, INTENT(IN) :: masse0(ijb_u:ije_u, llm) ! air mass REAL, INTENT(IN) :: phis0(ijb_u:ije_u) ! geopotentiat at the surface REAL :: zqmin, zqmax ! REAL,SAVE,ALLOCATABLE :: p (:,: ) ! pression aux interfac.des couches ! REAL,SAVE,ALLOCATABLE :: pks(:) ! exner au sol ! REAL,SAVE,ALLOCATABLE :: pk(:,:) ! exner au milieu des couches ! REAL,SAVE,ALLOCATABLE :: pkf(:,:) ! exner filt.au milieu des couches ! REAL,SAVE,ALLOCATABLE :: phi(:,:) ! geopotentiel ! REAL,SAVE,ALLOCATABLE :: w(:,:) ! vitesse verticale ! variables dynamiques intermediaire pour le transport ! REAL,SAVE,ALLOCATABLE :: pbaru(:,:),pbarv(:,:) !flux de masse ! variables dynamiques au pas -1 ! REAL,SAVE,ALLOCATABLE :: vcovm1(:,:),ucovm1(:,:) ! REAL,SAVE,ALLOCATABLE :: tetam1(:,:),psm1(:) ! REAL,SAVE,ALLOCATABLE :: massem1(:,:) ! tendances dynamiques ! REAL,SAVE,ALLOCATABLE :: dv(:,:),du(:,:) ! REAL,SAVE,ALLOCATABLE :: dteta(:,:),dp(:) ! REAL,DIMENSION(:,:,:), ALLOCATABLE, SAVE :: dq ! tendances de la dissipation ! REAL,SAVE,ALLOCATABLE :: dvdis(:,:),dudis(:,:) ! REAL,SAVE,ALLOCATABLE :: dtetadis(:,:) ! tendances physiques REAL, SAVE, ALLOCATABLE :: dvfi(:, :), dufi(:, :) REAL, SAVE, ALLOCATABLE :: dtetafi(:, :) REAL, SAVE, ALLOCATABLE :: dpfi(:) REAL, DIMENSION(:, :, :), ALLOCATABLE, SAVE :: dqfi ! variables pour le fichier histoire REAL :: dtav ! intervalle de temps elementaire REAL :: tppn(iim), tpps(iim), tpn, tps ! INTEGER :: itau, itaufinp1, iav ! INTEGER iday ! jour julien REAL :: time REAL :: SSUM ! REAL,SAVE,ALLOCATABLE :: finvmaold(:,:) !ym LOGICAL lafin LOGICAL :: lafin INTEGER :: ij, iq, l INTEGER :: ik REAL :: time_step, t_wrt, t_ops ! jD_cur: jour julien courant ! jH_cur: heure julienne courante REAL :: jD_cur, jH_cur INTEGER :: an, mois, jour REAL :: secondes logical :: physic LOGICAL :: first, callinigrads data callinigrads/.TRUE./ CHARACTER(LEN = 10) :: string10 ! REAL,SAVE,ALLOCATABLE :: flxw(:,:) ! flux de masse verticale !+jld variables test conservation energie ! REAL,SAVE,ALLOCATABLE :: ecin(:,:),ecin0(:,:) ! Tendance de la temp. potentiel d (theta)/ d t due a la ! tansformation d'energie cinetique en energie thermique ! cree par la dissipation ! REAL,SAVE,ALLOCATABLE :: dtetaecdt(:,:) ! REAL,SAVE,ALLOCATABLE :: vcont(:,:),ucont(:,:) ! REAL,SAVE,ALLOCATABLE :: vnat(:,:),unat(:,:) REAL :: d_h_vcol, d_qt, d_qw, d_ql, d_ec CHARACTER(len = 15) :: ztit !! INTEGER ip_ebil_dyn ! PRINT level for energy conserv. diag. ! SAVE ip_ebil_dyn ! DATA ip_ebil_dyn/0/ !-jld CHARACTER(LEN = 80) :: dynhist_file, dynhistave_file CHARACTER(LEN = *), parameter :: modname = "leapfrog_loc" CHARACTER(LEN = 80) :: abort_message logical, PARAMETER :: dissip_conservative = .TRUE. INTEGER :: testita PARAMETER (testita = 9) logical, parameter :: flag_verif = .FALSE. ! declaration liees au parallelisme INTEGER :: ierr LOGICAL :: FirstCaldyn LOGICAL :: FirstPhysic INTEGER :: ijb, ije, j, i type(Request) :: TestRequest type(Request) :: Request_Dissip type(Request) :: Request_physic INTEGER :: true_itau INTEGER :: iapptrac INTEGER :: AdjustCount ! INTEGER :: var_time LOGICAL :: ok_start_timer = .FALSE. LOGICAL, SAVE :: firstcall = .TRUE. TYPE(distrib), SAVE :: new_dist CALL check_isotopes(q0, ijb_u, ije_u, 'leapfrog204: debut') !$OMP MASTER ItCount = 0 !$OMP END MASTER true_itau = 0 FirstCaldyn = .TRUE. FirstPhysic = .TRUE. iapptrac = 0 AdjustCount = 0 lafin = .FALSE. if (nday>=0) THEN itaufin = nday * day_step else itaufin = -nday endif itaufinp1 = itaufin + 1 CALL check_isotopes(q0, ijb_u, ije_u, 'leapfrog 226') itau = 0 physic = .TRUE. if (iflag_phys==0.or.iflag_phys==2) physic = .FALSE. CALL init_nan CALL leapfrog_allocate ucov = ucov0 vcov = vcov0 teta = teta0 ps = ps0 masse = masse0 phis = phis0 q = q0 CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 239') ! iday = day_ini+itau/day_step ! time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 ! IF(time.GT.1.) THEN ! time = time-1. ! iday = iday+1 ! ENDIF ! Allocate variables depending on dynamic variable nqtot !$OMP MASTER if (firstcall) THEN ! ALLOCATE(p(ijb_u:ije_u,llmp1)) ! ALLOCATE(pks(ijb_u:ije_u)) ! ALLOCATE(pk(ijb_u:ije_u,llm)) ! ALLOCATE(pkf(ijb_u:ije_u,llm)) ! ALLOCATE(phi(ijb_u:ije_u,llm)) ! ALLOCATE(w(ijb_u:ije_u,llm)) ! ALLOCATE(pbaru(ip1jmp1,llm),pbarv(ip1jm,llm)) ! ALLOCATE(vcovm1(ijb_v:ije_v,llm),ucovm1(ijb_u:ije_u,llm)) ! ALLOCATE(tetam1(ijb_u:ije_u,llm),psm1(ijb_u:ije_u)) ! ALLOCATE(massem1(ijb_u:ije_u,llm)) ! ALLOCATE(dv(ijb_v:ije_v,llm),du(ijb_u:ije_u,llm)) ! ALLOCATE(dteta(ijb_u:ije_u,llm),dp(ijb_u:ije_u)) ! ALLOCATE(dvdis(ijb_v:ije_v,llm),dudis(ijb_u:ije_u,llm)) ! ALLOCATE(dtetadis(ijb_u:ije_u,llm)) ALLOCATE(dvfi(ijb_v:ije_v, llm), dufi(ijb_u:ije_u, llm)) ALLOCATE(dtetafi(ijb_u:ije_u, llm)) ALLOCATE(dpfi(ijb_u:ije_u)) ! ALLOCATE(dq(ijb_u:ije_u,llm,nqtot)) ALLOCATE(dqfi(ijb_u:ije_u, llm, nqtot)) ! ALLOCATE(dqfi_tmp(iip1,llm,nqtot)) ! ALLOCATE(finvmaold(ijb_u:ije_u,llm)) ! ALLOCATE(flxw(ijb_u:ije_u,llm)) ! ALLOCATE(ecin(ijb_u:ije_u,llm),ecin0(ijb_u:ije_u,llm)) ! ALLOCATE(dtetaecdt(ijb_u:ije_u,llm)) ! ALLOCATE(vcont(ijb_v:ije_v,llm),ucont(ijb_u:ije_u,llm)) ! ALLOCATE(vnat(ijb_v:ije_v,llm),unat(ijb_u:ije_u,llm)) endif !$OMP END MASTER !$OMP BARRIER ! CALL dynredem1_loc("restart.nc",0.0, ! & vcov,ucov,teta,q,masse,ps) !----------------------------------------------------------------------- ! On initialise la pression et la fonction d'Exner : ! -------------------------------------------------- !$OMP MASTER dq(:, :, :) = 0. CALL pression (ijnb_u, ap, bp, ps, p) !$OMP END MASTER if (pressure_exner) THEN CALL exner_hyb_loc(ijnb_u, ps, p, pks, pk, pkf) else CALL exner_milieu_loc(ijnb_u, ps, p, pks, pk, pkf) endif !----------------------------------------------------------------------- ! Debut de l'integration temporelle: ! ---------------------------------- ! et du parallelisme !! 1 CONTINUE ! Matsuno Forward step begins here ! date: (NB: date remains unchanged for Backward step) ! ----- jD_cur = jD_ref + day_ini - day_ref + & (itau + 1) / day_step jH_cur = jH_ref + start_time + & mod(itau + 1, day_step) / float(day_step) if (jH_cur > 1.0) THEN jD_cur = jD_cur + 1. jH_cur = jH_cur - 1. endif CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 321') if (ok_guide) THEN CALL guide_main(itau,ucov,vcov,teta,q,masse,ps) !$OMP BARRIER endif ! ! IF( MOD( itau, 10* day_step ).EQ.0 ) THEN ! CALL test_period ( ucov,vcov,teta,q,p,phis ) ! PRINT *,' ---- Test_period apres continue OK ! -----', itau ! ENDIF ! !ym CALL SCOPY( ijmllm ,vcov , 1, vcovm1 , 1 ) !ym CALL SCOPY( ijp1llm,ucov , 1, ucovm1 , 1 ) !ym CALL SCOPY( ijp1llm,teta , 1, tetam1 , 1 ) !ym CALL SCOPY( ijp1llm,masse, 1, massem1, 1 ) !ym CALL SCOPY( ip1jmp1, ps , 1, psm1 , 1 ) if (FirstCaldyn) THEN !$OMP MASTER ucovm1 = ucov vcovm1 = vcov tetam1 = teta massem1 = masse psm1 = ps ! Ehouarn: finvmaold is actually not used ! finvmaold = masse !$OMP END MASTER !$OMP BARRIER ! CALL filtreg_p ( finvmaold ,jjb_u,jje_u,jjb_u,jje_u,jjp1, llm, ! & -2,2, .TRUE., 1 ) else ! Save fields obtained at previous time step as '...m1' ijb = ij_begin ije = ij_end !$OMP MASTER psm1 (ijb:ije) = ps (ijb:ije) !$OMP END MASTER !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) DO l = 1, llm ije = ij_end ucovm1 (ijb:ije, l) = ucov (ijb:ije, l) tetam1 (ijb:ije, l) = teta (ijb:ije, l) massem1 (ijb:ije, l) = masse (ijb:ije, l) ! finvmaold(ijb:ije,l)=masse(ijb:ije,l) if (pole_sud) ije = ij_end - iip1 vcovm1(ijb:ije, l) = vcov (ijb:ije, l) ENDDO !$OMP ENDDO ! Ehouarn: finvmaold not used ! CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1, ! . llm, -2,2, .TRUE., 1 ) endif ! of if (FirstCaldyn) forward = .TRUE. leapf = .FALSE. dt = dtvr ! ... P.Le Van .26/04/94 .... !ym CALL SCOPY ( ijp1llm, masse, 1, finvmaold, 1 ) !ym CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 ) !ym ne sert a rien !ym CALL minmax(ijp1llm,q(:,:,3),zqmin,zqmax) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 400') 2 CONTINUE ! Matsuno backward or leapfrog step begins here CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 402') !$OMP MASTER ItCount = ItCount + 1 if (MOD(ItCount, 1)==1) THEN debug = .TRUE. else debug = .FALSE. endif !$OMP END MASTER !----------------------------------------------------------------------- ! date: (NB: only leapfrog step requires recomputing date) ! ----- IF (leapf) THEN jD_cur = jD_ref + day_ini - day_ref + & (itau + 1) / day_step jH_cur = jH_ref + start_time + & mod(itau + 1, day_step) / float(day_step) if (jH_cur > 1.0) THEN jD_cur = jD_cur + 1. jH_cur = jH_cur - 1. endif ENDIF ! gestion des appels de la physique et des dissipations: ! ------------------------------------------------------ ! ! ... P.Le Van ( 6/02/95 ) .... apphys = .FALSE. statcl = .FALSE. conser = .FALSE. apdiss = .FALSE. IF(purmats) THEN ! Purely Matsuno time stepping IF(MOD(itau, iconser) ==0.AND. forward) conser = .TRUE. IF(MOD(itau, dissip_period)==0.AND..NOT.forward) & apdiss = .TRUE. IF(MOD(itau, iphysiq)==0.AND..NOT.forward & .and. physic) apphys = .TRUE. ELSE ! Leapfrog/Matsuno time stepping IF(MOD(itau, iconser) == 0) conser = .TRUE. IF(MOD(itau + 1, dissip_period)==0 .AND. .NOT. forward) & apdiss = .TRUE. IF(MOD(itau + 1, iphysiq)==0.AND.physic) apphys = .TRUE. END IF ! Ehouarn: for Shallow Water case (ie: 1 vertical layer), ! supress dissipation step if (llm==1) THEN apdiss = .FALSE. endif !ym ---> Pour le moment !ym apphys = .FALSE. statcl = .FALSE. ! conser = .FALSE. ! ie: no output of control variables to stdout in // if (firstCaldyn) THEN !$OMP MASTER CALL Set_Distrib(distrib_caldyn) !$OMP END MASTER !$OMP BARRIER firstCaldyn = .FALSE. !ym CALL InitTime !$OMP MASTER CALL Init_timer !$OMP END MASTER endif !$OMP MASTER IF (ok_start_timer) THEN CALL InitTime ok_start_timer = .FALSE. ENDIF !$OMP END MASTER CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 471') !ym PAS D'AJUSTEMENT POUR LE MOMENT if (Adjust) THEN AdjustCount = AdjustCount + 1 ! if (iapptrac==iapp_tracvl .and. (forward .OR. leapf) ! & .and. itau/iphysiq>2 .and. Adjustcount>30) THEN if (Adjustcount>1) THEN AdjustCount = 0 !$OMP MASTER CALL allgather_timer_average if (prt_level > 9) THEN print *, '*********************************' print *, '****** TIMER CALDYN ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_caldyn(i), & ' : temps moyen :', & timer_average(jj_nb_caldyn(i), timer_caldyn, i), & '+-', timer_delta(jj_nb_caldyn(i), timer_caldyn, i) enddo print *, '*********************************' print *, '****** TIMER VANLEER ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_vanleer(i), & ' : temps moyen :', & timer_average(jj_nb_vanleer(i), timer_vanleer, i), & '+-', timer_delta(jj_nb_vanleer(i), timer_vanleer, i) enddo print *, '*********************************' print *, '****** TIMER DISSIP ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_dissip(i), & ' : temps moyen :', & timer_average(jj_nb_dissip(i), timer_dissip, i), & '+-', timer_delta(jj_nb_dissip(i), timer_dissip, i) enddo ! if (mpi_rank==0) CALL WriteBands endif CALL AdjustBands_caldyn(new_dist) !$OMP END MASTER !$OMP BARRIER CALL leapfrog_switch_caldyn(new_dist) !$OMP BARRIER !$OMP MASTER distrib_caldyn = new_dist CALL set_distrib(distrib_caldyn) !$OMP END MASTER !$OMP BARRIER ! CALL Register_SwapFieldHallo(ucov,ucov,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(ucovm1,ucovm1,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(vcov,vcov,ip1jm,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(vcovm1,vcovm1,ip1jm,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(teta,teta,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(tetam1,tetam1,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(masse,masse,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(massem1,massem1,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(ps,ps,ip1jmp1,1, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(psm1,psm1,ip1jmp1,1, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(pkf,pkf,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(pk,pk,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(pks,pks,ip1jmp1,1, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(phis,phis,ip1jmp1,1, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(phi,phi,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! CALL Register_SwapFieldHallo(finvmaold,finvmaold,ip1jmp1,llm, ! & jj_Nb_caldyn,0,0,TestRequest) ! do j=1,nqtot ! CALL Register_SwapFieldHallo(q(:,:,j),q(:,:,j),ip1jmp1,llm, ! & jj_nb_caldyn,0,0,TestRequest) ! enddo ! CALL Set_Distrib(distrib_caldyn) ! CALL SendRequest(TestRequest) ! CALL WaitRequest(TestRequest) !$OMP MASTER CALL AdjustBands_dissip(new_dist) !$OMP END MASTER !$OMP BARRIER CALL leapfrog_switch_dissip(new_dist) !$OMP BARRIER !$OMP MASTER distrib_dissip = new_dist !$OMP END MASTER !$OMP BARRIER ! CALL AdjustBands_physic !$OMP MASTER if (mpi_rank==0) CALL WriteBands !$OMP END MASTER endif endif CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 589') !----------------------------------------------------------------------- ! calcul des tendances dynamiques: ! -------------------------------- !$OMP BARRIER !$OMP MASTER CALL VTb(VThallo) !$OMP END MASTER CALL Register_Hallo_u(ucov, llm, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_v(vcov, llm, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_u(teta, llm, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_u(ps, 1, 1, 2, 2, 1, TestRequest) CALL Register_Hallo_u(pkf, llm, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_u(pk, llm, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_u(pks, 1, 1, 1, 1, 1, TestRequest) CALL Register_Hallo_u(p, llmp1, 1, 1, 1, 1, TestRequest) ! do j=1,nqtot ! CALL Register_Hallo(q(1,1,j),ip1jmp1,llm,1,1,1,1, ! * TestRequest) ! enddo CALL SendRequest(TestRequest) !$OMP BARRIER CALL WaitRequest(TestRequest) !$OMP MASTER CALL VTe(VThallo) !$OMP END MASTER !$OMP BARRIER if (debug) THEN CALL WriteField_u('ucov', ucov) CALL WriteField_v('vcov', vcov) CALL WriteField_u('teta', teta) CALL WriteField_u('ps', ps) CALL WriteField_u('masse', masse) CALL WriteField_u('pk', pk) CALL WriteField_u('pks', pks) CALL WriteField_u('pkf', pkf) CALL WriteField_u('phis', phis) do iq = 1, nqtot CALL WriteField_u('q' // trim(int2str(iq)), & q(:, :, iq)) enddo endif True_itau = True_itau + 1 !$OMP MASTER IF (prt_level>9) THEN WRITE(lunout, *)"leapfrog_p: Iteration No", True_itau ENDIF CALL start_timer(timer_caldyn) ! compute geopotential phi() CALL geopot_loc (ip1jmp1, teta, pk, pks, phis, phi) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 651') CALL VTb(VTcaldyn) !$OMP END MASTER ! var_time=time+iday-day_ini !$OMP BARRIER ! CALL FTRACE_REGION_BEGIN("caldyn") time = jD_cur + jH_cur CALL caldyn_loc & (itau, ucov, vcov, teta, ps, masse, pk, pkf, phis, & phi, conser, du, dv, dteta, dp, w, pbaru, pbarv, time) ! CALL FTRACE_REGION_END("caldyn") !$OMP MASTER if (mpi_rank==0.AND.conser) THEN WRITE(lunout, *) 'leapfrog_loc, Time step: ', itau, ' Day:', time ENDIF CALL VTe(VTcaldyn) !$OMP END MASTER IF (CPPKEY_DEBUGIO) THEN CALL WriteField_u('du', du) CALL WriteField_v('dv', dv) CALL WriteField_u('dteta', dteta) CALL WriteField_u('dp', dp) CALL WriteField_u('w', w) CALL WriteField_u('pbaru', pbaru) CALL WriteField_v('pbarv', pbarv) CALL WriteField_u('p', p) CALL WriteField_u('masse', masse) CALL WriteField_u('pk', pk) END IF !----------------------------------------------------------------------- ! calcul des tendances advection des traceurs (dont l'humidite) ! ------------------------------------------------------------- CALL check_isotopes(q, ijb_u, ije_u, & 'leapfrog 686: avant caladvtrac') IF(forward .OR. leapf) THEN ! Ehouarn: NB: fields sent to advtrac are those at the beginning of the time step !WRITE(*,*) 'leapfrog 679: avant CALL caladvtrac_loc' CALL caladvtrac_loc(q, pbaru, pbarv, & p, masse, dq, teta, & flxw, pk, iapptrac) ! CALL creation of mass flux IF (offline .AND. .NOT. adjust) THEN CALL fluxstokenc_p(pbaru, pbarv, masse, teta, phi) ENDIF !WRITE(*,*) 'leapfrog 719' CALL check_isotopes(q, ijb_u, ije_u, & 'leapfrog 698: apres caladvtrac') ! do j=1,nqtot ! CALL WriteField_u('qadv'//trim(int2str(j)),q(:,:,j)) ! enddo ! Ehouarn: Storage of mass flux for off-line tracers... not implemented... ENDIF ! of IF( forward .OR. leapf ) !----------------------------------------------------------------------- ! integrations dynamique et traceurs: ! ---------------------------------- !$OMP MASTER CALL VTb(VTintegre) !$OMP END MASTER IF (CPPKEY_DEBUGIO) THEN if (true_itau>20) THEN CALL WriteField_u('ucovm1', ucovm1) CALL WriteField_v('vcovm1', vcovm1) CALL WriteField_u('tetam1', tetam1) CALL WriteField_u('psm1', psm1) CALL WriteField_u('ucov_int', ucov) CALL WriteField_v('vcov_int', vcov) CALL WriteField_u('teta_int', teta) CALL WriteField_u('ps_int', ps) endif END IF !$OMP BARRIER ! CALL FTRACE_REGION_BEGIN("integrd") !WRITE(*,*) 'leapfrog 720' CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 756') ! CRisi: pourquoi aller jusqu'à 2 et non pas jusqu'à nqtot?? CALL integrd_loc (nqtot, vcovm1, ucovm1, tetam1, psm1, massem1, & dv, du, dteta, dq, dp, vcov, ucov, teta, q, ps, masse, phis) ! $ finvmaold ) ! !WRITE(*,*) 'leapfrog 724' CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 762') ! CALL FTRACE_REGION_END("integrd") !$OMP BARRIER IF (CPPKEY_DEBUGIO) THEN CALL WriteField_u('ucovm1', ucovm1) CALL WriteField_v('vcovm1', vcovm1) CALL WriteField_u('tetam1', tetam1) CALL WriteField_u('psm1', psm1) CALL WriteField_u('ucov_int', ucov) CALL WriteField_v('vcov_int', vcov) CALL WriteField_u('teta_int', teta) CALL WriteField_u('ps_int', ps) END IF CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 775') ! do j=1,nqtot ! CALL WriteField_p('q'//trim(int2str(j)), ! . reshape(q(:,:,j),(/iip1,jmp1,llm/))) ! CALL WriteField_p('dq'//trim(int2str(j)), ! . reshape(dq(:,:,j),(/iip1,jmp1,llm/))) ! enddo !$OMP MASTER CALL VTe(VTintegre) !$OMP END MASTER ! .P.Le Van (26/04/94 ajout de finvpold dans l'appel d'integrd) ! !----------------------------------------------------------------------- ! calcul des tendances physiques: ! ------------------------------- ! ######## P.Le Van ( Modif le 6/02/95 ) ########### ! IF(purmats) THEN IF(itau==itaufin.AND..NOT.forward) lafin = .TRUE. ELSE IF(itau + 1 == itaufin) lafin = .TRUE. ENDIF !c$OMP END PARALLEL ! ! IF(apphys) THEN CALL call_calfis(itau, lafin, ucov, vcov, teta, masse, ps, phis, q, flxw) ! c-jld !$OMP MASTER if (FirstPhysic) THEN ok_start_timer = .TRUE. FirstPhysic = .FALSE. endif !$OMP END MASTER ENDIF ! of IF( apphys ) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1132') !WRITE(*,*) 'leapfrog 1134: iflag_phys=',iflag_phys IF(iflag_phys==2) THEN ! "Newtonian" case !$OMP MASTER if (FirstPhysic) THEN ok_start_timer = .TRUE. FirstPhysic = .FALSE. endif !$OMP END MASTER ! Calcul academique de la physique = Rappel Newtonien + fritcion ! -------------------------------------------------------------- !ym teta(:,:)=teta(:,:) !ym s -iphysiq*dtvr*(teta(:,:)-tetarappel(:,:))/taurappel ijb = ij_begin ije = ij_end !LF teta(ijb:ije,:)=teta(ijb:ije,:) !LF s -iphysiq*dtvr*(teta(ijb:ije,:)-tetarappel(ijb:ije,:))/taurappel !$OMP DO SCHEDULE(STATIC,OMP_CHUNK) do l = 1, llm teta(ijb:ije, l) = teta(ijb:ije, l) - dtvr * & (teta(ijb:ije, l) - tetarappel(ijb:ije, l)) * & (knewt_g + knewt_t(l) * clat4(ijb:ije)) enddo !$OMP END DO !$OMP MASTER if (planet_type=="giant") THEN ! add an intrinsic heat flux at the base of the atmosphere teta(ijb:ije, 1) = teta(ijb:ije, 1) & + dtvr * aire(ijb:ije) * ihf / cpp / masse(ijb:ije, 1) endif !$OMP END MASTER !$OMP BARRIER CALL Register_Hallo_u(ucov, llm, 0, 1, 1, 0, Request_Physic) CALL Register_Hallo_v(vcov, llm, 1, 1, 1, 1, Request_Physic) CALL SendRequest(Request_Physic) !$OMP BARRIER CALL WaitRequest(Request_Physic) !$OMP BARRIER CALL friction_loc(ucov, vcov, dtvr) !$OMP BARRIER ! Sponge layer (if any) IF (ok_strato) THEN CALL top_bound_loc(vcov, ucov, teta, masse, dtvr) !$OMP BARRIER ENDIF ! of IF (ok_strato) ENDIF ! of IF(iflag_phys.EQ.2) CALL pression_loc (ip1jmp1, ap, bp, ps, p) !$OMP BARRIER if (pressure_exner) THEN CALL exner_hyb_loc(ijnb_u, ps, p, pks, pk, pkf) else CALL exner_milieu_loc(ijnb_u, ps, p, pks, pk, pkf) endif !$OMP BARRIER CALL massdair_loc(p, masse) !$OMP BARRIER !c$OMP END PARALLEL CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1196') !----------------------------------------------------------------------- ! dissipation horizontale et verticale des petites echelles: ! ---------------------------------------------------------- IF(apdiss) THEN CALL call_dissip(ucov, vcov, teta, p, pk, ps) END IF CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1430') ! ******************************************************************** ! .... fin de l'integration dynamique et physique pour le pas itau .. ! ******************************************************************** ! preparation du pas d'integration suivant ...... !$OMP MASTER CALL stop_timer(timer_caldyn) !$OMP END MASTER IF (itau==itaumax) THEN !$OMP MASTER CALL allgather_timer_average CALL barrier if (mpi_rank==0) THEN print *, '*********************************' print *, '****** TIMER CALDYN ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_caldyn(i), & ' : temps moyen :', & timer_average(jj_nb_caldyn(i), timer_caldyn, i) enddo print *, '*********************************' print *, '****** TIMER VANLEER ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_vanleer(i), & ' : temps moyen :', & timer_average(jj_nb_vanleer(i), timer_vanleer, i) enddo print *, '*********************************' print *, '****** TIMER DISSIP ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_dissip(i), & ' : temps moyen :', & timer_average(jj_nb_dissip(i), timer_dissip, i) enddo print *, '*********************************' print *, '****** TIMER PHYSIC ******' do i = 0, mpi_size - 1 print *, 'proc', i, ' : Nb Bandes :', jj_nb_physic(i), & ' : temps moyen :', & timer_average(jj_nb_physic(i), timer_physic, i) enddo endif CALL barrier print *, 'Taille du Buffer MPI (REAL*8)', MaxBufferSize print *, 'Taille du Buffer MPI utilise (REAL*8)', MaxBufferSize_Used print *, 'Temps total ecoule sur la parallelisation :', DiffTime() print *, 'Temps CPU ecoule sur la parallelisation :', DiffCpuTime() CALL print_filtre_timer !$OMP END MASTER CALL dynredem1_loc("restart.nc", 0.0, & vcov, ucov, teta, q, masse, ps) !$OMP MASTER CALL fin_getparam !$OMP END MASTER if (ok_guide) THEN ! set ok_guide to false to avoid extra output ! in following forward step ok_guide = .FALSE. endif IF (ANY(type_trac == ['inca', 'inco'])) THEN CALL finalize_inca ! switching back to LMDZDYN context !$OMP MASTER IF (ok_dyn_xios) THEN CALL xios_set_current_context(dyn3d_ctx_handle) ENDIF !$OMP END MASTER ENDIF #ifdef REPROBUS if (type_trac == 'repr') CALL finalize_reprobus #endif !$OMP MASTER CALL finalize_parallel !$OMP END MASTER !$OMP BARRIER RETURN ENDIF CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1509') IF (.NOT.purmats) THEN ! ........................................................ ! .............. schema matsuno + leapfrog .............. ! ........................................................ IF(forward .OR. leapf) THEN itau = itau + 1 ! iday= day_ini+itau/day_step ! time= REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 ! IF(time.GT.1.) THEN ! time = time-1. ! iday = iday+1 ! ENDIF ENDIF IF(itau == itaufinp1) THEN if (flag_verif) THEN WRITE(79, *) 'ucov', ucov WRITE(80, *) 'vcov', vcov WRITE(81, *) 'teta', teta WRITE(82, *) 'ps', ps WRITE(83, *) 'q', q WRITE(85, *) 'q1 = ', q(:, :, 1) WRITE(86, *) 'q3 = ', q(:, :, 3) endif !$OMP MASTER CALL fin_getparam !$OMP END MASTER IF (ANY(type_trac == ['inca', 'inco'])) THEN CALL finalize_inca ! switching back to LMDZDYN context !$OMP MASTER IF (ok_dyn_xios) THEN CALL xios_set_current_context(dyn3d_ctx_handle) ENDIF !$OMP END MASTER ENDIF #ifdef REPROBUS if (type_trac == 'repr') CALL finalize_reprobus #endif !$OMP MASTER CALL finalize_parallel !$OMP END MASTER abort_message = 'Simulation finished' CALL abort_gcm(modname, abort_message, 0) RETURN ENDIF !----------------------------------------------------------------------- ! ecriture du fichier histoire moyenne: ! ------------------------------------- IF(MOD(itau, iperiod)==0 .OR. itau==itaufin) THEN !$OMP BARRIER IF(itau==itaufin) THEN iav = 1 ELSE iav = 0 ENDIF ! Ehouarn: re-compute geopotential for outputs !$OMP BARRIER !$OMP MASTER CALL geopot_loc(ip1jmp1, teta, pk, pks, phis, phi) !$OMP END MASTER !$OMP BARRIER IF (ok_dynzon) THEN CALL bilan_dyn_loc(2,dtvr*iperiod,dtvr*day_step*periodav, & ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q) ENDIF !ok_dynzon IF (ok_dyn_ave) THEN CALL writedynav_loc(itau,vcov, & ucov,teta,pk,phi,q,masse,ps,phis) ENDIF ENDIF CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1584') !----------------------------------------------------------------------- ! ecriture de la bande histoire: ! ------------------------------ IF(MOD(itau, iecri)==0) THEN ! Ehouarn: output only during LF or Backward Matsuno if (leapf.or.(.not.leapf.and.(.not.forward))) THEN !$OMP BARRIER !$OMP MASTER CALL geopot_loc(ip1jmp1, teta, pk, pks, phis, phi) !$OMP END MASTER !$OMP BARRIER if (ok_dyn_ins) THEN CALL writehist_loc(itau,vcov,ucov,teta,pk,phi,q, & masse,ps,phis) endif IF (ok_dyn_xios) THEN !$OMP MASTER CALL xios_update_calendar(itau) !$OMP END MASTER !$OMP BARRIER CALL writedyn_xios(vcov, & ucov, teta, pk, phi, q, masse, ps, phis) ENDIF endif ! of if (leapf.or.(.not.leapf.and.(.not.forward))) ENDIF ! of IF(MOD(itau,iecri).EQ.0) IF(itau==itaufin) THEN !$OMP BARRIER ! if (planet_type.eq."earth") THEN ! Write an Earth-format restart file CALL dynredem1_loc("restart.nc", 0.0, & vcov, ucov, teta, q, masse, ps) ! endif ! of if (planet_type.eq."earth") if (ok_guide) THEN ! set ok_guide to false to avoid extra output ! in following forward step ok_guide = .FALSE. endif ! CLOSE(99) ENDIF ! of IF (itau.EQ.itaufin) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1624') !----------------------------------------------------------------------- ! gestion de l'integration temporelle: ! ------------------------------------ IF(MOD(itau, iperiod)==0) THEN GO TO 1 ELSE IF (MOD(itau - 1, iperiod) == 0) THEN IF(forward) THEN ! fin du pas forward et debut du pas backward forward = .FALSE. leapf = .FALSE. GO TO 2 ELSE ! fin du pas backward et debut du premier pas leapfrog leapf = .TRUE. dt = 2. * dtvr GO TO 2 END IF ELSE ! ...... pas leapfrog ..... leapf = .TRUE. dt = 2. * dtvr GO TO 2 END IF ! of IF (MOD(itau,iperiod).EQ.0) ! ELSEIF (MOD(itau-1,iperiod).EQ.0) ELSE ! of IF (.not.purmats) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1664') ! ........................................................ ! .............. schema matsuno ............... ! ........................................................ IF(forward) THEN itau = itau + 1 ! iday = day_ini+itau/day_step ! time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0 ! IF(time.GT.1.) THEN ! time = time-1. ! iday = iday+1 ! ENDIF forward = .FALSE. IF(itau == itaufinp1) THEN !$OMP MASTER CALL fin_getparam !$OMP END MASTER IF (ANY(type_trac == ['inca', 'inco'])) THEN CALL finalize_inca ! switching back to LMDZDYN context !$OMP MASTER IF (ok_dyn_xios) THEN CALL xios_set_current_context(dyn3d_ctx_handle) ENDIF !$OMP END MASTER ENDIF #ifdef REPROBUS if (type_trac == 'repr') CALL finalize_reprobus #endif !$OMP MASTER CALL finalize_parallel !$OMP END MASTER abort_message = 'Simulation finished' CALL abort_gcm(modname, abort_message, 0) RETURN ENDIF GO TO 2 ELSE ! of IF(forward) i.e. backward step CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1698') IF(MOD(itau, iperiod)==0 .OR. itau==itaufin) THEN IF(itau==itaufin) THEN iav = 1 ELSE iav = 0 ENDIF ! Ehouarn: re-compute geopotential for outputs !$OMP BARRIER !$OMP MASTER CALL geopot_loc(ip1jmp1,teta,pk,pks,phis,phi) !$OMP END MASTER !$OMP BARRIER IF (ok_dynzon) THEN CALL bilan_dyn_loc(2,dtvr*iperiod,dtvr*day_step*periodav, & ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q) ENDIF IF (ok_dyn_ave) THEN CALL writedynav_loc(itau,vcov, & ucov,teta,pk,phi,q,masse,ps,phis) ENDIF ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) IF(MOD(itau, iecri)==0) THEN !$OMP BARRIER !$OMP MASTER CALL geopot_loc(ip1jmp1, teta, pk, pks, phis, phi) !$OMP END MASTER !$OMP BARRIER if (ok_dyn_ins) THEN CALL writehist_loc(itau,vcov,ucov,teta,pk,phi,q, & masse,ps,phis) endif ! of if (ok_dyn_ins) IF (ok_dyn_xios) THEN !$OMP MASTER CALL xios_update_calendar(itau) !$OMP END MASTER !$OMP BARRIER CALL writedyn_xios(vcov, & ucov, teta, pk, phi, q, masse, ps, phis) ENDIF ENDIF ! of IF(MOD(itau,iecri).EQ.0) IF(itau==itaufin) THEN ! if (planet_type.eq."earth") THEN CALL dynredem1_loc("restart.nc", 0.0, & vcov, ucov, teta, q, masse, ps) ! endif ! of if (planet_type.eq."earth") if (ok_guide) THEN ! set ok_guide to false to avoid extra output ! in following forward step ok_guide = .FALSE. endif ENDIF ! of IF(itau.EQ.itaufin) forward = .TRUE. GO TO 1 ENDIF ! of IF (forward) CALL check_isotopes(q, ijb_u, ije_u, 'leapfrog 1750') END IF ! of IF(.not.purmats) !$OMP MASTER CALL fin_getparam !$OMP END MASTER IF (ANY(type_trac == ['inca', 'inco'])) THEN CALL finalize_inca ! switching back to LMDZDYN context !$OMP MASTER IF (ok_dyn_xios) THEN CALL xios_set_current_context(dyn3d_ctx_handle) ENDIF !$OMP END MASTER ENDIF #ifdef REPROBUS if (type_trac == 'repr') CALL finalize_reprobus #endif !$OMP MASTER CALL finalize_parallel !$OMP END MASTER abort_message = 'Simulation finished' CALL abort_gcm(modname, abort_message, 0) END SUBROUTINE leapfrog_loc