MODULE diagflux_mod USE icosa USE omp_para IMPLICIT NONE SAVE PRIVATE TYPE(t_field), POINTER, PUBLIC :: & f_masst(:), f_qmasst(:), & ! time-averaged mass, tracer mass, f_massfluxt(:), f_qfluxt(:), & ! time-integrated mass flux and tracer flux f_qfluxt_lon(:), f_qfluxt_lat(:), & ! scalar flux reconstructed at cell centers f_epot(:), f_ekin(:), f_enthalpy(:), & ! time-averaged potential E, kinetic E and enthalpy f_epotfluxt(:), f_ekinfluxt(:), f_enthalpyfluxt(:) ! time averaged 'fluxes' of epot, ekin and enthalpy LOGICAL :: diagflux_on !$OMP THREADPRIVATE(diagflux_on) PUBLIC :: diagflux_on, init_diagflux, zero_qfluxt, qflux_centered_lonlat, flux_centered_lonlat, diagflux_energy CONTAINS SUBROUTINE init_diagflux USE getin_mod INTEGER :: ll diagflux_on = .FALSE. CALL getin("diagflux", diagflux_on) ll = MERGE(llm,1,diagflux_on) CALL allocate_field(f_masst, field_t,type_real,ll, name="masst") CALL allocate_field(f_epot, field_t,type_real,ll, name="epot") CALL allocate_field(f_ekin, field_t,type_real,ll, name="ekin") CALL allocate_field(f_enthalpy, field_t,type_real,ll, name="enthalpy") CALL allocate_field(f_qmasst, field_t,type_real,ll,nqtot, name="qmasst") CALL allocate_field(f_massfluxt, field_u,type_real,ll, name="massfluxt") CALL allocate_field(f_epotfluxt, field_u,type_real,ll, name="epotfluxt") CALL allocate_field(f_ekinfluxt, field_u,type_real,ll, name="ekinfluxt") CALL allocate_field(f_enthalpyfluxt, field_u,type_real,ll, name="enthalpyfluxt") CALL allocate_field(f_qfluxt, field_u,type_real,ll,nqtot, name="qfluxt") CALL allocate_field(f_qfluxt_lon, field_t,type_real,ll,nqtot, name="qfluxt_lon") CALL allocate_field(f_qfluxt_lat, field_t,type_real,ll,nqtot, name="qfluxt_lat") IF(diagflux_on) CALL zero_qfluxt END SUBROUTINE init_diagflux #define ZERO2(field) buf2=field(ind) ; buf2(:,ll_begin:ll_end)=0. #define ZERO3(field) buf3=field(ind) ; buf3(:,ll_begin:ll_end,:)=0. SUBROUTINE zero_qfluxt INTEGER :: ind REAL(rstd), POINTER :: buf2(:,:),buf3(:,:,:) DO ind=1,ndomain IF (.NOT. assigned_domain(ind)) CYCLE CALL swap_dimensions(ind) ZERO2(f_masst) ZERO2(f_epot) ZERO2(f_ekin) ZERO2(f_enthalpy) ZERO3(f_qmasst) ZERO2(f_massfluxt) ZERO2(f_epotfluxt) ZERO2(f_ekinfluxt) ZERO2(f_enthalpyfluxt) ZERO3(f_qfluxt) END DO END SUBROUTINE zero_qfluxt !------------------------------------ Reconstruct fluxes at cell centers --------------------------------------- SUBROUTINE qflux_centered_lonlat(scale, f_flux, f_flux_lon, f_flux_lat) REAL(rstd), INTENT(IN) :: scale TYPE(t_field),POINTER :: f_flux(:), f_flux_lon(:), f_flux_lat(:) REAL(rstd), POINTER :: flux(:,:,:), flux_lon(:,:,:), flux_lat(:,:,:) INTEGER :: ind, itrac DO ind=1,ndomain IF (.NOT. assigned_domain(ind)) CYCLE CALL swap_dimensions(ind) CALL swap_geometry(ind) flux=f_flux(ind) flux_lon=f_flux_lon(ind) flux_lat=f_flux_lat(ind) DO itrac=1,nqtot CALL compute_flux_centered_lonlat(scale, flux(:,:,itrac), flux_lon(:,:,itrac), flux_lat(:,:,itrac)) END DO END DO END SUBROUTINE qflux_centered_lonlat SUBROUTINE flux_centered_lonlat(scale, f_flux, f_flux_lon, f_flux_lat) REAL(rstd), INTENT(IN) :: scale TYPE(t_field),POINTER :: f_flux(:), f_flux_lon(:), f_flux_lat(:) REAL(rstd), POINTER :: flux(:,:), flux_lon(:,:), flux_lat(:,:) INTEGER :: ind DO ind=1,ndomain IF (.NOT. assigned_domain(ind)) CYCLE CALL swap_dimensions(ind) CALL swap_geometry(ind) flux=f_flux(ind) flux_lon=f_flux_lon(ind) flux_lat=f_flux_lat(ind) CALL compute_flux_centered_lonlat(scale, flux, flux_lon, flux_lat) END DO END SUBROUTINE flux_centered_lonlat SUBROUTINE compute_flux_centered_lonlat(scale, flux, flux_lon, flux_lat) USE wind_mod REAL(rstd), INTENT(IN) :: scale REAL(rstd), INTENT(IN) :: flux(3*iim*jjm,llm) REAL(rstd), INTENT(OUT) :: flux_lon(iim*jjm,llm), flux_lat(iim*jjm,llm) REAL(rstd) :: flux_3d(iim*jjm,llm,3) CALL compute_flux_centered(scale, flux, flux_3d) CALL compute_wind_centered_lonlat_compound(flux_3d, flux_lon, flux_lat) END SUBROUTINE compute_flux_centered_lonlat !------------------------------------ Compute energy fluxes --------------------------------------- SUBROUTINE diagflux_energy(frac, f_phis,f_rhodz,f_theta_rhodz,f_u, f_geopot,f_theta, f_hfluxt) REAL(rstd), INTENT(IN) :: frac TYPE(t_field),POINTER :: f_phis(:),f_rhodz(:),f_theta_rhodz(:),f_u(:), f_geopot(:), f_theta(:), f_hfluxt(:) REAL(rstd), POINTER :: phis(:), rhodz(:,:), theta_rhodz(:,:,:), u(:,:), & geopot(:,:), pk(:,:,:), hfluxt(:,:), & epot(:,:), ekin(:,:), enthalpy(:,:), & epotflux(:,:), ekinflux(:,:), enthalpyflux(:,:) INTEGER :: ind DO ind=1,ndomain IF (.NOT. assigned_domain(ind)) CYCLE CALL swap_dimensions(ind) CALL swap_geometry(ind) hfluxt = f_hfluxt(ind) phis = f_phis(ind) rhodz = f_rhodz(ind) theta_rhodz = f_theta_rhodz(ind) u = f_u(ind) geopot = f_geopot(ind) pk = f_theta(ind) ! buffer epot = f_epot(ind) ekin = f_ekin(ind) enthalpy = f_enthalpy(ind) epotflux = f_epotfluxt(ind) ekinflux = f_ekinfluxt(ind) enthalpyflux = f_enthalpyfluxt(ind) CALL compute_diagflux_energy(frac,hfluxt, phis,rhodz,theta_rhodz,u, geopot,pk, epot,ekin,enthalpy, epotflux, ekinflux, enthalpyflux) END DO END SUBROUTINE diagflux_energy SUBROUTINE compute_diagflux_energy(frac, massflux, phis,rhodz,theta_rhodz,u, geopot,pk, epot,ekin,enthalpy, epot_flux, ekin_flux, enthalpy_flux) USE disvert_mod, ONLY : ptop REAL(rstd), INTENT(IN) :: frac REAL(rstd), INTENT(IN) :: massflux(3*iim*jjm,llm), u(3*iim*jjm,llm),& phis(iim*jjm), rhodz(iim*jjm,llm), theta_rhodz(iim*jjm,llm,nqtot) REAL(rstd), INTENT(INOUT) :: geopot(iim*jjm,llm+1), pk(iim*jjm,llm) ! pk = buffer REAL(rstd), INTENT(INOUT), DIMENSION(iim*jjm, llm) :: epot, ekin, enthalpy REAL(rstd), INTENT(INOUT), DIMENSION(3*iim*jjm, llm) :: epot_flux, ekin_flux, enthalpy_flux REAL(rstd) :: energy, p_ik, theta_ik, temp_ik, gv, Rd INTEGER :: ij, l, ij_omp_begin_ext, ij_omp_end_ext Rd = kappa*cpp CALL distrib_level(ij_end_ext-ij_begin_ext+1,ij_omp_begin_ext,ij_omp_end_ext) ij_omp_begin_ext = ij_omp_begin_ext+ij_begin_ext-1 ij_omp_end_ext = ij_omp_end_ext+ij_begin_ext-1 #include "../kernels/energy_fluxes.k90" END SUBROUTINE compute_diagflux_energy END MODULE diagflux_mod