ecearth3/sources/lpjg/offline/ccycle_coupling.f90

program CC_COUP

! First implementation: May 2020 (Etienne Tourigny) - dummy tm5 coupling - co2 ppm + fluxes, need to set co2 ppm value in namelist to activate this (284 for PI), or set to 0 to deactivate
! update: June 2020 (Etienne Tourigny) - always send co2, option to couple to pisces (disabled)

  use mod_oasis
  use mpi

  implicit none

  integer,parameter :: TM5GRID = 10800 ! CTM3
  integer,parameter :: t_max_year = 365
  integer,parameter :: t_max_leap = 366

  integer :: num_years,start_year
  double precision :: co2_ppm
  logical :: coupled_to_lpjguess, coupled_to_pisces
  namelist /NAMCC_COUP/ num_years,start_year,co2_ppm
  character(len=*),parameter :: NamelistFileName  = "namelist.ccycle_coupling"
  integer,parameter :: NAMCC_COUPU = 999
  logical :: CC_COUPnotFree

  character(len=6),parameter :: model_name    = "CC_COUP"

! dummy tm5 coupling - co2 ppm + fluxes
  character(len=128),parameter :: f_lco2_name = "LCO2_TM5"
  integer                    :: f_lco2_id   = -77
  double precision           :: f_lco2_data(TM5GRID,1)

  character(len=128),parameter :: f_oco2_name = "OCO2_TM5"
  integer                    :: f_oco2_id   = -77
  double precision           :: f_oco2_data(TM5GRID,1)

!  character(len=128),parameter :: f_cflx_name = "TM5_LandCFLX"
!  integer                    :: f_cflx_id   = -77
!  double precision           :: f_cflx_data(TM5GRID,1)

  character(len=128),parameter :: f_cnat_name = "TM5_LandCNAT"
  integer                    :: f_cnat_id   = -77
  double precision           :: f_cnat_data(TM5GRID,1)

  character(len=128),parameter :: f_cant_name = "TM5_LandCANT"
  integer                    :: f_cant_id   = -77
  double precision           :: f_cant_data(TM5GRID,1)

  character(len=128),parameter :: f_cnpp_name = "TM5_LandCNPP"
  integer                    :: f_cnpp_id   = -77
  double precision           :: f_cnpp_data(TM5GRID,1)

! TODO add option for lpjg,pisces,and feedback
  character(len=128),parameter :: f_ocecflx_name = "TM5_OceCFLX"
  integer                    :: f_ocecflx_id   = -77
  double precision           :: f_ocecflx_data(TM5GRID,1)

  integer :: comp_id = -77
  integer :: part_id_tm5 = -77
  integer :: t,t_step,yr,lyr,t_max,t_step_full,cell,loopyr
  integer :: ierror

  integer::ix,iy

  logical :: isleapyear = .false.
  character yearstr*4

  integer::localcomm, cplcomm, icpl;

! *** START ***

  write (*,'(A)') "*II* ccycle_coupling: Hello"
  write (*,*) "*II* ccycle_coupling: TM5GRID=",TM5GRID

  ! Read RunLengthSec,TimeStepSec from namelist
  inquire(NAMCC_COUPU,opened=CC_COUPnotFree)
  if (CC_COUPnotFree) call ERROR('Namelist CCOUP not free',1)
  open(NAMCC_COUPU,file=NamelistFileName,status='OLD')
  read(NAMCC_COUPU,nml=NAMCC_COUP)
  close(NAMCC_COUPU)

  ! by default we are coupled to lpjg and not pisces
  coupled_to_lpjguess = .true.
  coupled_to_pisces = .false.


  write (*,'(A)') "*II* ccycle_coupling: Now initialising Ccycle_coupling using oasis_..."

  call oasis_init_comp(comp_id,model_name,ierror)
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_init_comp returned  ierror=",ierror

  icpl = 1
  call oasis_create_couplcomm(icpl,MPI_COMM_SELF, cplcomm, ierror)
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_create_couplcomm  ierror=",ierror
  write (*,'(A,I12)') "*II* ccycle_coupling: oasis_create_couplcomm returned cplcomm =",cplcomm

! dummy tm5 coupling
  call oasis_def_partition(part_id_tm5,(/ 0,0,TM5GRID /),ierror)
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_def_partition returned part_id_co2 =",part_id_tm5
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_def_partition returned  ierror =",ierror
  if ( coupled_to_lpjguess ) then
     call oasis_def_var( f_lco2_id, f_lco2_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_Out, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
     call oasis_def_var( f_cnat_id, f_cnat_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_In, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
     call oasis_def_var( f_cant_id, f_cant_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_In, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
     call oasis_def_var( f_cnpp_id, f_cnpp_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_In, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
  endif
  if ( coupled_to_pisces ) then
     call oasis_def_var( f_oco2_id, f_oco2_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_Out, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
     call oasis_def_var( f_ocecflx_id, f_ocecflx_name,   &
          part_id_tm5, (/ 2,1 /), PRISM_In, (/ 1,TM5GRID,1,1 /), PRISM_Real, ierror )
  endif

  write (*,'(A)') "*II* ccycle_coupling: before call oasis_enddef(ierror)"
  call oasis_enddef(ierror)
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_enddef returned ierror =",ierror

  write (*,'(A)') "*II* ccycle_coupling: Beginning time loop"

  t_step_full = 0

! Repetitions of forcing
!  do lyr = 1, num_loops

! Year loop
     do loopyr = start_year, start_year+num_years-1

        write (*,*) "*III* ccycle_coupling: Beginning year",loopyr

        yr = loopyr

! Leap year?
        isleapyear = .false.
        if (mod(yr,4) .eq. 0)   isleapyear = .true.
        if (mod(yr,100) .eq. 0) isleapyear = .false.
        if (mod(yr,400) .eq. 0) isleapyear = .true.

        if (isleapyear .eq. .true.) then
           t_max = t_max_leap
        else
           t_max = t_max_year
        endif

        ! Day/6-hr loop for each year

        write (*,'(A,I5)') "*II* ccycle_coupling: t_max =",t_max

        do t_step=0,t_max-1
           write (*,*) "*III* ccycle_coupling: Beginning dayloop ",t_step

           t = t_step_full*86400

           ! -----------------------------------------------------------------------------------
           ! *** PUT variables
           ! -----------------------------------------------------------------------------------

           write (*,'(A,I6,A,I12)') "*II* ccycle_couplingPUT : Time step t = ",t_step," - time t=",t

           write (*,*) "*II* ccycle_coupling: co2_ppm = ",co2_ppm

           if ( coupled_to_lpjguess ) then
              f_lco2_data(:,:) = co2_ppm
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_put with f_lco2_id =",f_lco2_id
              call oasis_put(f_lco2_id,t,f_lco2_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_put returned ierror =",ierror
           endif

           if ( coupled_to_pisces ) then
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_put with f_oco2_id =",f_oco2_id
              f_oco2_data(:,:) = co2_ppm
              call oasis_put(f_oco2_id,t,f_oco2_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_put returned ierror =",ierror
           endif

           ! -----------------------------------------------------------------------------------
           ! *** GET variables
           ! -----------------------------------------------------------------------------------

           write (*,'(A,I6,A,I12)') "*II* ccycle_couplingGET : Time step t = ",t_step," - time t=",t

           if ( coupled_to_lpjguess ) then
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_get with f_cnat_id =",f_cnat_id
              call oasis_get(f_cnat_id,t,f_cnat_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_get returned ierror =",ierror
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_get with f_cant_id =",f_cant_id
              call oasis_get(f_cant_id,t,f_cant_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_get returned ierror =",ierror
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_get with f_cnpp_id =",f_cnpp_id
              call oasis_get(f_cnpp_id,t,f_cnpp_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_get returned ierror =",ierror
              write(*,*)"*II* ccycle_coupling: CNAT obtained from OASIS for test gridcell (g=8000) is ",f_cnat_data(8000,1)
              write(*,*)"*II* ccycle_coupling: CANT obtained from OASIS for test gridcell (g=8000) is ",f_cant_data(8000,1)
              write(*,*)"*II* ccycle_coupling: CNPP obtained from OASIS for test gridcell (g=8000) is ",f_cnpp_data(8000,1)
           endif
           if ( coupled_to_pisces ) then
              write (*,'(A,I3)') "*II* ccycle_coupling: calling oasis_get with f_ocecflx_id =",f_ocecflx_id
              call oasis_get(f_ocecflx_id,t,f_ocecflx_data,ierror)
              write (*,'(A,I3)') "*II* ccycle_coupling: oasis_get returned ierror =",ierror
              write(*,*)"*II* ccycle_coupling: OCECFLX obtained from OASIS for test gridcell (g=8000) is ",f_ocecflx_data(8000,1)
           endif

           write (*,*) "*III* ccycle_coupling: Finished dayloop!",t_step, t_max

           ! Increase the full simulation counter
           t_step_full = t_step_full + 1

           ! End of dayloop
        enddo

        write (*,*) "*II* ccycle_coupling: End of year!",loopyr

        ! End of year loop
     enddo

     write (*,'(A)') "*II* ccycle_coupling: Finished time loop!"

     ! End of repetition loop
!  enddo

  write (*,'(A)') "*II* ccycle_coupling: oasis_terminate..."
  call oasis_terminate(ierror)
  write (*,'(A,I3)') "*II* ccycle_coupling: oasis_terminate returned ierror =",ierror

end program CC_COUP


subroutine ERROR(msg,status)
    character(len=*),intent(in) :: msg
    integer,intent(in) :: status

    write (*,'("*EE*",A)') msg
    call EXIT(status)
end subroutine ERROR