!$Id: sulfint.F90 104 2008-12-23 10:28:51Z acosce $ !! ========================================================================= !! INCA - INteraction with Chemistry and Aerosols !! !! Copyright Laboratoire des Sciences du Climat et de l'Environnement (LSCE) !! Unite mixte CEA-CNRS-UVSQ !! !! Contributors to this INCA subroutine: !! !! Didier Hauglustaine, LSCE, hauglustaine@cea.fr !! !! Anne Cozic, LSCE, anne.cozic@cea.fr !! Yann Meurdesoif, LSCE, yann.meurdesoif@cea.fr !! !! This software is a computer program whose purpose is to simulate the !! atmospheric gas phase and aerosol composition. The model is designed to be !! used within a transport model or a general circulation model. This version !! of INCA was designed to be coupled to the LMDz GCM. LMDz-INCA accounts !! for emissions, transport (resolved and sub-grid scale), photochemical !! transformations, and scavenging (dry deposition and washout) of chemical !! species and aerosols interactively in the GCM. Several versions of the INCA !! model are currently used depending on the envisaged applications with the !! chemistry-climate model. !! !! This software is governed by the CeCILL license under French law and !! abiding by the rules of distribution of free software. You can use, !! modify and/ or redistribute the software under the terms of the CeCILL !! license as circulated by CEA, CNRS and INRIA at the following URL !! "http://www.cecill.info". !! !! As a counterpart to the access to the source code and rights to copy, !! modify and redistribute granted by the license, users are provided only !! with a limited warranty and the software's author, the holder of the !! economic rights, and the successive licensors have only limited !! liability. !! !! In this respect, the user's attention is drawn to the risks associated !! with loading, using, modifying and/or developing or reproducing the !! software by the user in light of its specific status of free software, !! that may mean that it is complicated to manipulate, and that also !! therefore means that it is reserved for developers and experienced !! professionals having in-depth computer knowledge. Users are therefore !! encouraged to load and test the software's suitability as regards their !! requirements in conditions enabling the security of their systems and/or !! data to be ensured and, more generally, to use and operate it in the !! same conditions as regards security. !! !! The fact that you are presently reading this means that you have had !! knowledge of the CeCILL license and that you accept its terms. !! ========================================================================= #include SUBROUTINE SULF_INTI (filename) !---------------------------------------------------------------------- ! ... SULFATE climatologies initialize ! Didier Hauglustaine, IPSL, 2000. !---------------------------------------------------------------------- USE AC_SULF USE INCA_DIM USE MOD_INCA_PARA IMPLICIT NONE INCLUDE 'netcdf.inc' !---------------------------------------------------------------------- ! ... Dummy args !---------------------------------------------------------------------- CHARACTER(len=*), INTENT(in) :: filename !---------------------------------------------------------------------- ! ... Local variables !---------------------------------------------------------------------- INTEGER :: iret INTEGER :: varid INTEGER :: error ! ... Open the file !$OMP MASTER IF (is_mpi_root) THEN iret = nf_open(filename, 0, ncid) CALL check_err(iret, 'SULF_INTI', 'problem to open file') ! ... Get lengths iret = nf_inq_dimid(ncid, 'time', varid) CALL check_err(iret, 'SULF_INTI', 'problem to check dimid time') iret = nf_inq_dimlen(ncid, varid, ntimes) CALL check_err(iret, 'SULF_INTI', 'problem to check dimlen time') iret = nf_inq_dimid(ncid, 'lev', varid) CALL check_err(iret, 'SULF_INTI', 'problem to check dimid lev') iret = nf_inq_dimlen(ncid, varid, nlevs) CALL check_err(iret, 'SULF_INTI', 'problem to check dimlen lev') iret = nf_inq_dimid(ncid, 'vector', varid) CALL check_err(iret, 'SULF_INTI', 'problem to check dimid vector') iret = nf_inq_dimlen(ncid, varid, klons) CALL check_err(iret, 'SULF_INTI', 'problem to check dimlen vector') ENDIF !$OMP END MASTER CALL bcast(ntimes) CALL bcast(nlevs) CALL bcast(klons) IF (klons/=nbp_glo) THEN write(lunout,*) 'sulf_inti : [klons] [nbp_glo]', klons, nbp_glo call print_err(3, 'SULF_INTI','there a problem of vector size', 'check klons and nbp_glo', '') ELSE klons=PLON ENDIF ! ... Allocate variables ALLOCATE(times(ntimes), STAT=error) IF (error /= 0) call print_err(3, 'SULF_INTI', 'Space requested not possible for variable times ', '', '') ALLOCATE(levs(nlevs), STAT=error) IF (error /= 0) CALL print_err(3, 'SULF_INTI', 'Space requested not possible for variable levs', '', '') ALLOCATE(so4bd(klons,nlevs,2), STAT=error) IF (error /= 0) CALL print_err(3, 'SULF_INTI', 'Space requested not possible for variable so4bd ', '', '') ! ... Get the coordinates !$OMP MASTER IF (is_mpi_root) THEN iret = nf_inq_varid(ncid, 'time', varid) CALL check_err(iret, 'SULF_INTI', 'problem to check varid time') iret = nf_get_var_double(ncid, varid, times) CALL check_err(iret, 'SULF_INTI', 'problem to read time') iret = nf_inq_varid(ncid, 'lev', varid) CALL check_err(iret, 'SULF_INTI', 'problem to check varid lev') iret = nf_get_var_double(ncid, varid, levs) CALL check_err(iret, 'SULF_INTI', 'problem to read lev') ! ... Variable ids iret = nf_inq_varid(ncid, 'so4', so4_id) CALL check_err(iret, 'SULF_INTI', 'problem to check varid so4') ENDIF !$OMP END MASTER CALL bcast(times) CALL bcast(so4_id) call bcast(levs) END SUBROUTINE SULF_INTI SUBROUTINE SULF_READ(calday) !---------------------------------------------------------------------- ! ... Read SO4 distributions ! Didier Hauglustaine, IPSL, 2000. !---------------------------------------------------------------------- USE AC_SULF USE INCA_DIM USE MOD_INCA_PARA USE PRINT_INCA IMPLICIT NONE INCLUDE 'netcdf.inc' !---------------------------------------------------------------------- ! ... Dummy args !---------------------------------------------------------------------- REAL, INTENT(IN) :: calday !---------------------------------------------------------------------- ! ... Local variables !---------------------------------------------------------------------- INTEGER :: iret INTEGER :: varid INTEGER :: oldlotime, oldhitime INTEGER, DIMENSION(3) :: start3, count3 INTEGER, DIMENSION(2) :: start2, count2 REAL :: so4bd_glo(nbp_glo,nlevs,2) ! ... Check to see if model time still bounded by data set oldlotime = lotime oldhitime = hitime CALL FINDPLB (times, ntimes, calday, lotime) IF ( cyclical ) THEN hitime = MOD(lotime,ntimes) + 1 ELSE hitime = lotime + 1 END IF ! ... Read new data if needed IF ( hitime /= oldhitime ) THEN !$OMP MASTER IF (is_mpi_root) THEN loind = hiind hiind = MOD( loind, 2) + 1 count2(1) = nbp_glo count2(2) = 1 start2(1) = 1 count3(1) = nbp_glo count3(2) = nlevs count3(3) = 1 start3(1) = 1 start3(2) = 1 WRITE(lunout, *) 'SULF_READ : read new data for time ', times(hitime) start3(3) = hitime start2(2) = hitime iret = nf_get_vara_double(ncid, so4_id, start3, count3, so4bd_glo(1,1,hiind)) CALL check_err(iret, 'SULF_READ', 'problem to read so4bd_glo hiind') ENDIF !$OMP END MASTER CALL scatter(so4bd_glo(:,:,hiind),so4bd(:,:,hiind)) IF ( lotime /= oldlotime ) THEN !$OMP MASTER IF (is_mpi_root) THEN WRITE(lunout, *) 'SULF_READ : read new data for time ', times(lotime) start3(3) = lotime start2(2) = lotime iret = nf_get_vara_double(ncid, so4_id, start3, count3, so4bd_glo(1,1,loind)) CALL check_err(iret, 'SULF_READ', 'problem to read so4bd_glo loind') ENDIF !$OMP END MASTER CALL scatter(so4bd_glo(:,:,loind),so4bd(:,:,loind)) END IF END IF END SUBROUTINE SULF_READ SUBROUTINE SULF_INTERP(calday, pmid, zmid) !---------------------------------------------------------------------- ! ... Interpolate SO4 on model time ! Didier Hauglustaine, IPSL, 2000. !---------------------------------------------------------------------- USE AC_SULF USE INCA_DIM IMPLICIT NONE !---------------------------------------------------------------------- ! ... Dummy args !---------------------------------------------------------------------- REAL, INTENT(IN) :: calday REAL, INTENT(IN) :: pmid(PLON,PLEV) REAL, INTENT(IN) :: zmid(PLON,PLEV) !---------------------------------------------------------------------- ! ... Local variables !---------------------------------------------------------------------- REAL, PARAMETER :: dzs = 1. !specific to Considine's input file INTEGER :: k, i, m, j, l INTEGER :: lonlev REAL :: dt, dt1, tint REAL :: wght REAL, DIMENSION(PLON,PLEV) :: so4wrk lonlev = PLON * PLEV !---------------------------------------------------------------------- ! Note: 365 day year inconsistent with LMDz (360 days) !!! !---------------------------------------------------------------------- ! ... interpolate linearly on current model time IF ( times(hitime) < times(lotime) ) THEN dt = 365. - times(lotime) + times(hitime) IF (calday <= times(hitime) ) THEN dt1 = 365. - times(lotime) + calday ELSE dt1 = calday - times(lotime) END IF ELSE dt = times(hitime) - times(lotime) dt1 = calday - times(lotime) END IF tint = dt1/dt ! ... Interpolate to local time CALL linintp (& lonlev, & 0., & 1., & tint, & so4bd(1,1,loind), & so4bd(1,1,hiind), & so4wrk) !-------------------------------------------------------- ! ... Change units (from kg/m3 to g/cm3) !-------------------------------------------------------- DO i = 1, PLON DO l = 1, PLEV so4(i,l) = so4wrk(i,l) *1.e-3 END DO END DO END SUBROUTINE SULF_INTERP