From b69b6d2dea8536b85c68f762dbe524e30c0243bb Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Mon, 10 Sep 2018 09:37:38 -0600 Subject: [PATCH 1/6] ozphys_2015.f --- physics/ozphys_2015.f | 108 +++++++++++++++++++++++++++++++++++++ physics/ozphys_2015.f.ORIG | 108 +++++++++++++++++++++++++++++++++++++ 2 files changed, 216 insertions(+) create mode 100755 physics/ozphys_2015.f create mode 100755 physics/ozphys_2015.f.ORIG diff --git a/physics/ozphys_2015.f b/physics/ozphys_2015.f new file mode 100755 index 000000000..1d7cad57c --- /dev/null +++ b/physics/ozphys_2015.f @@ -0,0 +1,108 @@ + subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, + & prsl, prdout, pl_coeff, delp, ldiag3d, + & ozp,me) +! +! this code assumes that both prsl and po3 are from bottom to top +! as are all other variables +! This code is specifically for NRL parameterization and +! climatological T and O3 are in location 5 and 6 of prdout array +! June 2015 - Shrinivas Moorthi +! + use machine , only : kind_phys + use physcons, only : grav => con_g + implicit none +! + real, parameter :: gravi=1.0/grav + integer im, ix, levs, ko3, pl_coeff,me + real(kind=kind_phys) ozi(ix,levs), ozo(ix,levs), po3(ko3), + & prsl(ix,levs), tin(ix,levs), delp(ix,levs), + & prdout(ix,ko3,pl_coeff), + & ozp(ix,levs,4), dt +! + integer k,kmax,kmin,l,i,j + logical ldiag3d, flg(im) + real(kind=kind_phys) pmax, pmin, tem, temp + real(kind=kind_phys) wk1(im), wk2(im), wk3(im), prod(im,pl_coeff), + & ozib(im), colo3(im,levs+1), coloz(im,levs+1) +! + colo3(:,levs+1) = 0.0 + coloz(:,levs+1) = 0.0 +! + do l=levs,1,-1 + pmin = 1.0e10 + pmax = -1.0e10 +! + do i=1,im + wk1(i) = log(prsl(i,l)) + pmin = min(wk1(i), pmin) + pmax = max(wk1(i), pmax) + prod(i,:) = 0.0 + enddo + kmax = 1 + kmin = 1 + do k=1,ko3-1 + if (pmin < po3(k)) kmax = k + if (pmax < po3(k)) kmin = k + enddo +! + do k=kmin,kmax + temp = 1.0 / (po3(k) - po3(k+1)) + do i=1,im + flg(i) = .false. + if (wk1(i) < po3(k) .and. wk1(i) >= po3(k+1)) then + flg(i) = .true. + wk2(i) = (wk1(i) - po3(k+1)) * temp + wk3(i) = 1.0 - wk2(i) + endif + enddo + do j=1,pl_coeff + do i=1,im + if (flg(i)) then + prod(i,j) = wk2(i) * prdout(i,k,j) + & + wk3(i) * prdout(i,k+1,j) + endif + enddo + enddo + enddo +! + do j=1,pl_coeff + do i=1,im + if (wk1(i) < po3(ko3)) then + prod(i,j) = prdout(i,ko3,j) + endif + if (wk1(i) >= po3(1)) then + prod(i,j) = prdout(i,1,j) + endif + enddo + enddo + do i=1,im + colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l)*gravi + coloz(i,l) = coloz(i,l+1) + prod(i,6) * delp(i,l)*gravi + prod(i,2) = min(prod(i,2), 0.0) + enddo +! write(1000+me,*) ' colo3=',colo3(1,l),' coloz=',coloz(1,l) +! &,' l=',l + do i=1,im + ozib(i) = ozi(i,l) ! no filling + tem = prod(i,1) - prod(i,2) * prod(i,6) + & + prod(i,3) * (tin(i,l) - prod(i,5)) + & + prod(i,4) * (colo3(i,l)-coloz(i,l)) + +! if (me .eq. 0) print *,'ozphys_2015 tem=',tem,' prod=',prod(i,:) +! &,' ozib=',ozib(i),' l=',l,' tin=',tin(i,l),'colo3=',colo3(i,l+1) + + ozo(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) + enddo + if (ldiag3d) then ! ozone change diagnostics + do i=1,im + ozp(i,l,1) = ozp(i,l,1) + (prod(i,1)-prod(i,2)*prod(i,6))*dt + ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) + ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*(tin(i,l)-prod(i,5))*dt + ozp(i,l,4) = ozp(i,l,4) + prod(i,4) + & * (colo3(i,l)-coloz(i,l))*dt + enddo + endif + enddo ! vertical loop +! + return + end diff --git a/physics/ozphys_2015.f.ORIG b/physics/ozphys_2015.f.ORIG new file mode 100755 index 000000000..1d7cad57c --- /dev/null +++ b/physics/ozphys_2015.f.ORIG @@ -0,0 +1,108 @@ + subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, + & prsl, prdout, pl_coeff, delp, ldiag3d, + & ozp,me) +! +! this code assumes that both prsl and po3 are from bottom to top +! as are all other variables +! This code is specifically for NRL parameterization and +! climatological T and O3 are in location 5 and 6 of prdout array +! June 2015 - Shrinivas Moorthi +! + use machine , only : kind_phys + use physcons, only : grav => con_g + implicit none +! + real, parameter :: gravi=1.0/grav + integer im, ix, levs, ko3, pl_coeff,me + real(kind=kind_phys) ozi(ix,levs), ozo(ix,levs), po3(ko3), + & prsl(ix,levs), tin(ix,levs), delp(ix,levs), + & prdout(ix,ko3,pl_coeff), + & ozp(ix,levs,4), dt +! + integer k,kmax,kmin,l,i,j + logical ldiag3d, flg(im) + real(kind=kind_phys) pmax, pmin, tem, temp + real(kind=kind_phys) wk1(im), wk2(im), wk3(im), prod(im,pl_coeff), + & ozib(im), colo3(im,levs+1), coloz(im,levs+1) +! + colo3(:,levs+1) = 0.0 + coloz(:,levs+1) = 0.0 +! + do l=levs,1,-1 + pmin = 1.0e10 + pmax = -1.0e10 +! + do i=1,im + wk1(i) = log(prsl(i,l)) + pmin = min(wk1(i), pmin) + pmax = max(wk1(i), pmax) + prod(i,:) = 0.0 + enddo + kmax = 1 + kmin = 1 + do k=1,ko3-1 + if (pmin < po3(k)) kmax = k + if (pmax < po3(k)) kmin = k + enddo +! + do k=kmin,kmax + temp = 1.0 / (po3(k) - po3(k+1)) + do i=1,im + flg(i) = .false. + if (wk1(i) < po3(k) .and. wk1(i) >= po3(k+1)) then + flg(i) = .true. + wk2(i) = (wk1(i) - po3(k+1)) * temp + wk3(i) = 1.0 - wk2(i) + endif + enddo + do j=1,pl_coeff + do i=1,im + if (flg(i)) then + prod(i,j) = wk2(i) * prdout(i,k,j) + & + wk3(i) * prdout(i,k+1,j) + endif + enddo + enddo + enddo +! + do j=1,pl_coeff + do i=1,im + if (wk1(i) < po3(ko3)) then + prod(i,j) = prdout(i,ko3,j) + endif + if (wk1(i) >= po3(1)) then + prod(i,j) = prdout(i,1,j) + endif + enddo + enddo + do i=1,im + colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l)*gravi + coloz(i,l) = coloz(i,l+1) + prod(i,6) * delp(i,l)*gravi + prod(i,2) = min(prod(i,2), 0.0) + enddo +! write(1000+me,*) ' colo3=',colo3(1,l),' coloz=',coloz(1,l) +! &,' l=',l + do i=1,im + ozib(i) = ozi(i,l) ! no filling + tem = prod(i,1) - prod(i,2) * prod(i,6) + & + prod(i,3) * (tin(i,l) - prod(i,5)) + & + prod(i,4) * (colo3(i,l)-coloz(i,l)) + +! if (me .eq. 0) print *,'ozphys_2015 tem=',tem,' prod=',prod(i,:) +! &,' ozib=',ozib(i),' l=',l,' tin=',tin(i,l),'colo3=',colo3(i,l+1) + + ozo(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) + enddo + if (ldiag3d) then ! ozone change diagnostics + do i=1,im + ozp(i,l,1) = ozp(i,l,1) + (prod(i,1)-prod(i,2)*prod(i,6))*dt + ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) + ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*(tin(i,l)-prod(i,5))*dt + ozp(i,l,4) = ozp(i,l,4) + prod(i,4) + & * (colo3(i,l)-coloz(i,l))*dt + enddo + endif + enddo ! vertical loop +! + return + end From 49fb209f3e7a177aa89f60a3c7ae56e754e12c78 Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Thu, 13 Sep 2018 16:15:43 -0600 Subject: [PATCH 2/6] add ccpp-compliant ozphys_2015.f --- physics/ozphys_2015.f | 103 ++++++++++++++++++++++++++++++++++++------ 1 file changed, 89 insertions(+), 14 deletions(-) diff --git a/physics/ozphys_2015.f b/physics/ozphys_2015.f index 1d7cad57c..f606810cd 100755 --- a/physics/ozphys_2015.f +++ b/physics/ozphys_2015.f @@ -1,6 +1,60 @@ - subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, - & prsl, prdout, pl_coeff, delp, ldiag3d, - & ozp,me) +!> \file ozphys_2015.f +!! This file is ozone sources and sinks. + + +!> This module contains the CCPP-compliant Ozone 2015 photochemistry scheme. + module ozphys_2015 + + contains + +!> \section arg_table_ozphys_2015_init Argument Table +!! + subroutine ozphys_2015_init() + end subroutine ozphys_2015_init + +! \brief Brief description of the subroutine +! +!> \section arg_table_ozphys_2015_finalize Argument Table +!! + subroutine ozphys_2015_finalize() + end subroutine ozphys_2015_finalize + + +!>\defgroup GFS_ozphys_2015 GFS ozphys_2015 Main +!! @{ +!! \brief The operational GFS currently parameterizes ozone production +!and +!! destruction based on monthly mean coefficients ( +!! \c ozprdlos_2015_new_sbuvO3_tclm15_nuchem.f77) provided by Naval +!! Research Laboratory through CHEM2D chemistry model +!! (McCormack et al. (2006) \cite mccormack_et_al_2006). +!! \section arg_table_ozphys_2015_run Argument Table +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |----------------|---------------------------------------------------|---------------------------------------------------|---------|------|-----------|-----------|--------|----------| +!! | ix | horizontal_dimension | horizontal dimension | count | 0 | integer | | in | F | +!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | +!! | levs | vertical_dimension | number of vertical layers | count | 0 | integer | | in | F | +!! | ko3 | vertical_dimension_of_ozone_forcing_data | number of vertical layers in ozone forcing data | count | 0 | integer | | in | F | +!! | dt | time_step_for_physics | physics time step | s | 0 | real | kind_phys | in | F | +!! | oz | ozone_concentration_updated_by_physics | ozone concentration updated by physics | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | tin | air_temperature_updated_by_physics | updated air temperature | K | 2 | real | kind_phys | in | F | +!! | po3 | natural_log_of_ozone_forcing_data_pressure_levels | natural log of ozone forcing data pressure levels | log(Pa) | 1 | real | kind_phys | in | F | +!! | prsl | air_pressure | mid-layer pressure | Pa | 2 | real | kind_phys | in | F | +!! | prdout | ozone_forcing | ozone forcing data | various | 3 | real | kind_phys | in | F | +!! | pl_coeff | number_of_coefficients_in_ozone_forcing_data | number of coefficients in ozone forcing data | index | 0 | integer | | in | F | +!! | delp | air_pressure_difference_between_midlayers | difference between mid-layer pressures | Pa | 2 | real | kind_phys | in | F | +!! | ldiag3d | flag_diagnostics_3D | flag for calculating 3-D diagnostic fields | flag | 0 | logical | | in | F | +!! | ozp | change_in_ozone_concentration | change in ozone concentration | kg kg-1 | 3 | real | kind_phys | inout | F | +!! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! +!> \section genal_ozphys_2015 GFS ozphys_2015_run General Algorithm +!! @{ + subroutine ozphys_2015_run ( & + & ix, im, levs, ko3, dt, oz, tin, po3, & + & prsl, prdout, pl_coeff, delp, ldiag3d, & + & ozp,me, errmsg, errflg) ! ! this code assumes that both prsl and po3 are from bottom to top ! as are all other variables @@ -13,18 +67,32 @@ subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, implicit none ! real, parameter :: gravi=1.0/grav - integer im, ix, levs, ko3, pl_coeff,me - real(kind=kind_phys) ozi(ix,levs), ozo(ix,levs), po3(ko3), - & prsl(ix,levs), tin(ix,levs), delp(ix,levs), - & prdout(ix,ko3,pl_coeff), - & ozp(ix,levs,4), dt -! + integer, intent(in) :: im, ix, levs, ko3, pl_coeff,me + real(kind=kind_phys), intent(in) :: po3(ko3), & + & prsl(ix,levs), tin(ix,levs), & + & delp(ix,levs), & + & prdout(ix,ko3,pl_coeff), dt + real(kind=kind_phys), intent(inout) :: ozp(ix,levs,4) + real(kind=kind_phys), intent(inout) :: oz(ix,levs) + + + character(len=*), intent(out) :: errmsg + integer, intent(out) :: errflg + integer k,kmax,kmin,l,i,j logical ldiag3d, flg(im) real(kind=kind_phys) pmax, pmin, tem, temp - real(kind=kind_phys) wk1(im), wk2(im), wk3(im), prod(im,pl_coeff), - & ozib(im), colo3(im,levs+1), coloz(im,levs+1) + real(kind=kind_phys) wk1(im), wk2(im), wk3(im),prod(im,pl_coeff), & + & ozib(im), colo3(im,levs+1), coloz(im,levs+1),& + & ozi(ix,levs) ! + ! Initialize CCPP error handling variables + errmsg = '' + errflg = 0 + +!ccpp: save input oz in ozi + ozi = oz + colo3(:,levs+1) = 0.0 coloz(:,levs+1) = 0.0 ! @@ -91,12 +159,14 @@ subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, ! if (me .eq. 0) print *,'ozphys_2015 tem=',tem,' prod=',prod(i,:) ! &,' ozib=',ozib(i),' l=',l,' tin=',tin(i,l),'colo3=',colo3(i,l+1) - ozo(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) +!ccpp ozo(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) + oz(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) enddo if (ldiag3d) then ! ozone change diagnostics do i=1,im ozp(i,l,1) = ozp(i,l,1) + (prod(i,1)-prod(i,2)*prod(i,6))*dt - ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) +!ccpp ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) + ozp(i,l,2) = ozp(i,l,2) + (oz(i,l) - ozib(i)) ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*(tin(i,l)-prod(i,5))*dt ozp(i,l,4) = ozp(i,l,4) + prod(i,4) & * (colo3(i,l)-coloz(i,l))*dt @@ -105,4 +175,9 @@ subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, enddo ! vertical loop ! return - end + end subroutine ozphys_2015_run + +!! @} +!! @} + + end module ozphys_2015 From 102a248f734042a3b9286b3b1f2d411044824b75 Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Wed, 19 Sep 2018 13:48:26 -0600 Subject: [PATCH 3/6] seperate diag%dq3dt(:,:6:9) --- physics/ozphys.f | 83 +++++++++++++++++++++++++------------------ physics/ozphys_2015.f | 54 +++++++++++++++------------- 2 files changed, 77 insertions(+), 60 deletions(-) diff --git a/physics/ozphys.f b/physics/ozphys.f index ccfe39527..6cae531aa 100644 --- a/physics/ozphys.f +++ b/physics/ozphys.f @@ -29,53 +29,61 @@ end subroutine ozphys_finalize !! Research Laboratory through CHEM2D chemistry model !! (McCormack et al. (2006) \cite mccormack_et_al_2006). !! \section arg_table_ozphys_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |----------------|---------------------------------------------------|---------------------------------------------------|---------|------|-----------|-----------|--------|----------| -!! | ix | horizontal_dimension | horizontal dimension | count | 0 | integer | | in | F | -!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | -!! | levs | vertical_dimension | number of vertical layers | count | 0 | integer | | in | F | -!! | ko3 | vertical_dimension_of_ozone_forcing_data | number of vertical layers in ozone forcing data | count | 0 | integer | | in | F | -!! | dt | time_step_for_physics | physics time step | s | 0 | real | kind_phys | in | F | -!! | oz | ozone_concentration_updated_by_physics | ozone concentration updated by physics | kg kg-1 | 2 | real | kind_phys | inout | F | -!! | tin | air_temperature_updated_by_physics | updated air temperature | K | 2 | real | kind_phys | in | F | -!! | po3 | natural_log_of_ozone_forcing_data_pressure_levels | natural log of ozone forcing data pressure levels | log(Pa) | 1 | real | kind_phys | in | F | -!! | prsl | air_pressure | mid-layer pressure | Pa | 2 | real | kind_phys | in | F | -!! | prdout | ozone_forcing | ozone forcing data | various | 3 | real | kind_phys | in | F | -!! | oz_coeff | number_of_coefficients_in_ozone_forcing_data | number of coefficients in ozone forcing data | index | 0 | integer | | in | F | -!! | delp | air_pressure_difference_between_midlayers | difference between mid-layer pressures | Pa | 2 | real | kind_phys | in | F | -!! | ldiag3d | flag_diagnostics_3D | flag for calculating 3-D diagnostic fields | flag | 0 | logical | | in | F | -!! | ozp | change_in_ozone_concentration | change in ozone concentration | kg kg-1 | 3 | real | kind_phys | inout | F | -!! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |----------------|---------------------------------------------------------------------------|----------------------------------------------------------------------------|---------|------|-----------|-----------|--------|----------| +!! | ix | horizontal_dimension | horizontal dimension | count | 0 | integer | | in | F | +!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | +!! | levs | vertical_dimension | number of vertical layers | count | 0 | integer | | in | F | +!! | ko3 | vertical_dimension_of_ozone_forcing_data | number of vertical layers in ozone forcing data | count | 0 | integer | | in | F | +!! | dt | time_step_for_physics | physics time step | s | 0 | real | kind_phys | in | F | +!! | oz | ozone_concentration_updated_by_physics | ozone concentration updated by physics | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | tin | air_temperature_updated_by_physics | updated air temperature | K | 2 | real | kind_phys | in | F | +!! | po3 | natural_log_of_ozone_forcing_data_pressure_levels | natural log of ozone forcing data pressure levels | log(Pa) | 1 | real | kind_phys | in | F | +!! | prsl | air_pressure | mid-layer pressure | Pa | 2 | real | kind_phys | in | F | +!! | prdout | ozone_forcing | ozone forcing coefficients | various | 3 | real | kind_phys | in | F | +!! | oz_coeff | number_of_coefficients_in_ozone_forcing_data | number of coefficients in ozone forcing data | index | 0 | integer | | in | F | +!! | delp | air_pressure_difference_between_midlayers | difference between mid-layer pressures | Pa | 2 | real | kind_phys | in | F | +!! | ldiag3d | flag_diagnostics_3D | flag for calculating 3-D diagnostic fields | flag | 0 | logical | | in | F | +!! | ozp1 | cumulative_change_in_ozone_concentration_due_to_production_and_loss_rate | cumulative change in ozone concentration due to production and loss rate | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp2 | cumulative_change_in_ozone_concentration_due_to_ozone_mixing_ratio | cumulative change in ozone concentration due to ozone mixing ratio | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp3 | cumulative_change_in_ozone_concentration_due_to_temperature | cumulative change in ozone concentration due to temperature | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp4 | cumulative_change_in_ozone_concentration_due_to_overhead_ozone_column | cumulative change in ozone concentration due to overhead ozone column | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | con_g | gravitational_acceleration | gravitational acceleration | m s-2 | 0 | real | kind_phys | in | F | +!! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | !! !> \section genal_ozphys GFS ozphys_run General Algorithm !! @{ subroutine ozphys_run ( & & ix, im, levs, ko3, dt, oz, tin, po3, & & prsl, prdout, oz_coeff, delp, ldiag3d, & - & ozp, me, errmsg, errflg) + & ozp1, ozp2, ozp3, ozp4, con_g, me, errmsg, errflg) ! ! this code assumes that both prsl and po3 are from bottom to top ! as are all other variables ! use machine , only : kind_phys - use physcons, only : grav => con_g +! use physcons, only : grav => con_g implicit none ! ! Interface variables integer, intent(in) :: im, ix, levs, ko3, oz_coeff, me real(kind=kind_phys), intent(inout) :: & - & oz(ix,levs), ozp(ix,levs,oz_coeff) + & oz(ix,levs), & + & ozp1(ix,levs), ozp2(ix,levs), ozp3(ix,levs), & + & ozp4(ix,levs) real(kind=kind_phys), intent(in) :: & & dt, po3(ko3), prdout(ix,ko3,oz_coeff), & - & prsl(ix,levs), tin(ix,levs), delp(ix,levs) + & prsl(ix,levs), tin(ix,levs), delp(ix,levs), & + & con_g logical, intent(in) :: ldiag3d + character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg ! ! Local variables - real, parameter :: gravi=1.0/grav +! real, parameter :: gravi=1.0/grav integer k,kmax,kmin,l,i,j logical flg(im) real(kind=kind_phys) pmax, pmin, tem, temp @@ -88,13 +96,18 @@ subroutine ozphys_run ( & ! ! save input oz in ozi ozi = oz + + if (me .eq. 0) write(6,*) 'CCPP DUBUG in ozphys_run: oz_coeff, ko3& + & =', oz_coeff, ko3 + if (me .eq. 0) write(6,*) 'CCPP DUBUG in ozphys_run: po3(ko3) = ',& + & po3(:) ! !> - Calculate vertical integrated column ozone values. if (oz_coeff > 2) then colo3(:,levs+1) = 0.0 do l=levs,1,-1 do i=1,im - colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l) * gravi + colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l) * 1.0/con_g enddo enddo endif @@ -156,16 +169,16 @@ subroutine ozphys_run ( & ! if (ldiag3d) then ! ozone change diagnostics do i=1,im - ozp(i,l,1) = ozp(i,l,1) + prod(i,1)*dt - ozp(i,l,2) = ozp(i,l,2) + (oz(i,l) - ozib(i)) + ozp1(i,l) = ozp1(i,l) + prod(i,1)*dt + ozp2(i,l) = ozp2(i,l) + (oz(i,l) - ozib(i)) enddo endif endif !> - Calculate the 4 terms of prognostic ozone change during time \a dt: -!! - ozp(:,:,1) - Ozone production at model layers -!! - ozp(:,:,2) - Ozone tendency at model layers -!! - ozp(:,:,3) - Ozone production from temperature term at model layers -!! - ozp(:,:,4) - Ozone production from column ozone term at model layers +!! - ozp1(:,:) - Ozone production from production/loss ratio +!! - ozp2(:,:) - Ozone production from ozone mixing ratio +!! - ozp3(:,:) - Ozone production from temperature term at model layers +!! - ozp4(:,:) - Ozone production from column ozone term at model layers if (oz_coeff == 4) then do i=1,im ozib(i) = ozi(i,l) ! no filling @@ -177,10 +190,10 @@ subroutine ozphys_run ( & enddo if (ldiag3d) then ! ozone change diagnostics do i=1,im - ozp(i,l,1) = ozp(i,l,1) + prod(i,1)*dt - ozp(i,l,2) = ozp(i,l,2) + (oz(i,l) - ozib(i)) - ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*tin(i,l)*dt - ozp(i,l,4) = ozp(i,l,4) + prod(i,4)*colo3(i,l+1)*dt + ozp1(i,l) = ozp1(i,l) + prod(i,1)*dt + ozp2(i,l) = ozp2(i,l) + (oz(i,l) - ozib(i)) + ozp3(i,l) = ozp3(i,l) + prod(i,3)*tin(i,l)*dt + ozp4(i,l) = ozp4(i,l) + prod(i,4)*colo3(i,l+1)*dt enddo endif endif diff --git a/physics/ozphys_2015.f b/physics/ozphys_2015.f index f606810cd..da0c7f931 100755 --- a/physics/ozphys_2015.f +++ b/physics/ozphys_2015.f @@ -29,32 +29,35 @@ end subroutine ozphys_2015_finalize !! Research Laboratory through CHEM2D chemistry model !! (McCormack et al. (2006) \cite mccormack_et_al_2006). !! \section arg_table_ozphys_2015_run Argument Table -!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | -!! |----------------|---------------------------------------------------|---------------------------------------------------|---------|------|-----------|-----------|--------|----------| -!! | ix | horizontal_dimension | horizontal dimension | count | 0 | integer | | in | F | -!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | -!! | levs | vertical_dimension | number of vertical layers | count | 0 | integer | | in | F | -!! | ko3 | vertical_dimension_of_ozone_forcing_data | number of vertical layers in ozone forcing data | count | 0 | integer | | in | F | -!! | dt | time_step_for_physics | physics time step | s | 0 | real | kind_phys | in | F | -!! | oz | ozone_concentration_updated_by_physics | ozone concentration updated by physics | kg kg-1 | 2 | real | kind_phys | inout | F | -!! | tin | air_temperature_updated_by_physics | updated air temperature | K | 2 | real | kind_phys | in | F | -!! | po3 | natural_log_of_ozone_forcing_data_pressure_levels | natural log of ozone forcing data pressure levels | log(Pa) | 1 | real | kind_phys | in | F | -!! | prsl | air_pressure | mid-layer pressure | Pa | 2 | real | kind_phys | in | F | -!! | prdout | ozone_forcing | ozone forcing data | various | 3 | real | kind_phys | in | F | -!! | pl_coeff | number_of_coefficients_in_ozone_forcing_data | number of coefficients in ozone forcing data | index | 0 | integer | | in | F | -!! | delp | air_pressure_difference_between_midlayers | difference between mid-layer pressures | Pa | 2 | real | kind_phys | in | F | -!! | ldiag3d | flag_diagnostics_3D | flag for calculating 3-D diagnostic fields | flag | 0 | logical | | in | F | -!! | ozp | change_in_ozone_concentration | change in ozone concentration | kg kg-1 | 3 | real | kind_phys | inout | F | -!! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | -!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | -!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | +!! | local_name | standard_name | long_name | units | rank | type | kind | intent | optional | +!! |----------------|--------------------------------------------------------------------------|----------------------------------------------------------------------------|---------|------|-----------|-----------|--------|----------| +!! | ix | horizontal_dimension | horizontal dimension | count | 0 | integer | | in | F | +!! | im | horizontal_loop_extent | horizontal loop extent | count | 0 | integer | | in | F | +!! | levs | vertical_dimension | number of vertical layers | count | 0 | integer | | in | F | +!! | ko3 | vertical_dimension_of_ozone_forcing_data | number of vertical layers in ozone forcing data | count | 0 | integer | | in | F | +!! | dt | time_step_for_physics | physics time step | s | 0 | real | kind_phys | in | F | +!! | oz | ozone_concentration_updated_by_physics | ozone concentration updated by physics | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | tin | air_temperature_updated_by_physics | updated air temperature | K | 2 | real | kind_phys | in | F | +!! | po3 | natural_log_of_ozone_forcing_data_pressure_levels | natural log of ozone forcing data pressure levels | log(Pa) | 1 | real | kind_phys | in | F | +!! | prsl | air_pressure | mid-layer pressure | Pa | 2 | real | kind_phys | in | F | +!! | prdout | ozone_forcing | ozone forcing data | various | 3 | real | kind_phys | in | F | +!! | pl_coeff | number_of_coefficients_in_ozone_forcing_data | number of coefficients in ozone forcing data | index | 0 | integer | | in | F | +!! | delp | air_pressure_difference_between_midlayers | difference between mid-layer pressures | Pa | 2 | real | kind_phys | in | F | +!! | ldiag3d | flag_diagnostics_3D | flag for calculating 3-D diagnostic fields | flag | 0 | logical | | in | F | +!! | ozp1 | cumulative_change_in_ozone_concentration_due_to_production_and_loss_rate | cumulative change in ozone concentration due to production and loss rate | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp2 | cumulative_change_in_ozone_concentration_due_to_ozone_mixing_ratio | cumulative change in ozone concentration due to ozone mixing ratio | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp3 | cumulative_change_in_ozone_concentration_due_to_temperature | cumulative change in ozone concentration due to temperature | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | ozp4 | cumulative_change_in_ozone_concentration_due_to_overhead_ozone_column | cumulative change in ozone concentration due to overhead ozone column | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | +!! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | +!! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | !! !> \section genal_ozphys_2015 GFS ozphys_2015_run General Algorithm !! @{ subroutine ozphys_2015_run ( & & ix, im, levs, ko3, dt, oz, tin, po3, & & prsl, prdout, pl_coeff, delp, ldiag3d, & - & ozp,me, errmsg, errflg) + & ozp1,ozp2,ozp3,ozp4,me, errmsg, errflg) ! ! this code assumes that both prsl and po3 are from bottom to top ! as are all other variables @@ -72,7 +75,8 @@ subroutine ozphys_2015_run ( & & prsl(ix,levs), tin(ix,levs), & & delp(ix,levs), & & prdout(ix,ko3,pl_coeff), dt - real(kind=kind_phys), intent(inout) :: ozp(ix,levs,4) + real(kind=kind_phys), intent(inout) :: ozp1(ix,levs), & + & ozp2(ix,levs), ozp3(ix,levs),ozp4(ix,levs) real(kind=kind_phys), intent(inout) :: oz(ix,levs) @@ -164,11 +168,11 @@ subroutine ozphys_2015_run ( & enddo if (ldiag3d) then ! ozone change diagnostics do i=1,im - ozp(i,l,1) = ozp(i,l,1) + (prod(i,1)-prod(i,2)*prod(i,6))*dt + ozp1(i,l) = ozp1(i,l) + (prod(i,1)-prod(i,2)*prod(i,6))*dt !ccpp ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) - ozp(i,l,2) = ozp(i,l,2) + (oz(i,l) - ozib(i)) - ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*(tin(i,l)-prod(i,5))*dt - ozp(i,l,4) = ozp(i,l,4) + prod(i,4) + ozp2(i,l) = ozp2(i,l) + (oz(i,l) - ozib(i)) + ozp3(i,l) = ozp3(i,l) + prod(i,3)*(tin(i,l)-prod(i,5))*dt + ozp4(i,l) = ozp4(i,l) + prod(i,4) & * (colo3(i,l)-coloz(i,l))*dt enddo endif From ada85a40bc081da20d4d591e9a587fb1c22f85ff Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Wed, 19 Sep 2018 14:41:36 -0600 Subject: [PATCH 4/6] clean up codes --- physics/ozinterp.f90 | 7 ++++++- physics/ozphys.f | 5 ----- physics/ozphys_2015.f | 11 ++++++++--- 3 files changed, 14 insertions(+), 9 deletions(-) diff --git a/physics/ozinterp.f90 b/physics/ozinterp.f90 index e7704924e..dad0ef9fa 100644 --- a/physics/ozinterp.f90 +++ b/physics/ozinterp.f90 @@ -67,6 +67,11 @@ SUBROUTINE read_o3data (ntoz, me, master) oz_pres(:) = log(100.0*oz_pres(:)) oz_lat(:) = oz_lat4(:) oz_time(:) = oz_time4(:) + if (me == master) then + write(*,*) ' oz_pres(:) = ', oz_pres(:) + write(*,*) ' oz_lat(:) = ', oz_lat(:) + write(*,*) ' oz_time(:) = ', oz_time(:) + endif deallocate (oz_lat4, oz_pres4, oz_time4) !--- read in ozplin which is in order of (lattitudes, ozone levels, coeff number, time) @@ -203,4 +208,4 @@ SUBROUTINE ozinterpol(me,npts,IDATE,FHOUR,jindx1,jindx2,ozplout,ddy) RETURN END -end module ozinterp \ No newline at end of file +end module ozinterp diff --git a/physics/ozphys.f b/physics/ozphys.f index 6cae531aa..0d69c364a 100644 --- a/physics/ozphys.f +++ b/physics/ozphys.f @@ -96,11 +96,6 @@ subroutine ozphys_run ( & ! ! save input oz in ozi ozi = oz - - if (me .eq. 0) write(6,*) 'CCPP DUBUG in ozphys_run: oz_coeff, ko3& - & =', oz_coeff, ko3 - if (me .eq. 0) write(6,*) 'CCPP DUBUG in ozphys_run: po3(ko3) = ',& - & po3(:) ! !> - Calculate vertical integrated column ozone values. if (oz_coeff > 2) then diff --git a/physics/ozphys_2015.f b/physics/ozphys_2015.f index da0c7f931..d1d683e90 100755 --- a/physics/ozphys_2015.f +++ b/physics/ozphys_2015.f @@ -48,6 +48,7 @@ end subroutine ozphys_2015_finalize !! | ozp2 | cumulative_change_in_ozone_concentration_due_to_ozone_mixing_ratio | cumulative change in ozone concentration due to ozone mixing ratio | kg kg-1 | 2 | real | kind_phys | inout | F | !! | ozp3 | cumulative_change_in_ozone_concentration_due_to_temperature | cumulative change in ozone concentration due to temperature | kg kg-1 | 2 | real | kind_phys | inout | F | !! | ozp4 | cumulative_change_in_ozone_concentration_due_to_overhead_ozone_column | cumulative change in ozone concentration due to overhead ozone column | kg kg-1 | 2 | real | kind_phys | inout | F | +!! | con_g | gravitational_acceleration | gravitational acceleration | m s-2 | 0 | real | kind_phys | in | F | !! | me | mpi_rank | rank of the current MPI task | index | 0 | integer | | in | F | !! | errmsg | ccpp_error_message | error message for error handling in CCPP | none | 0 | character | len=* | out | F | !! | errflg | ccpp_error_flag | error flag for error handling in CCPP | flag | 0 | integer | | out | F | @@ -57,7 +58,8 @@ end subroutine ozphys_2015_finalize subroutine ozphys_2015_run ( & & ix, im, levs, ko3, dt, oz, tin, po3, & & prsl, prdout, pl_coeff, delp, ldiag3d, & - & ozp1,ozp2,ozp3,ozp4,me, errmsg, errflg) + & ozp1,ozp2,ozp3,ozp4,con_g, & + & me, errmsg, errflg) ! ! this code assumes that both prsl and po3 are from bottom to top ! as are all other variables @@ -66,10 +68,12 @@ subroutine ozphys_2015_run ( & ! June 2015 - Shrinivas Moorthi ! use machine , only : kind_phys - use physcons, only : grav => con_g + !use physcons, only : grav => con_g implicit none ! - real, parameter :: gravi=1.0/grav + !real, parameter :: gravi=1.0/grav + real(kind=kind_phys),intent(in) :: con_g + real :: gravi integer, intent(in) :: im, ix, levs, ko3, pl_coeff,me real(kind=kind_phys), intent(in) :: po3(ko3), & & prsl(ix,levs), tin(ix,levs), & @@ -96,6 +100,7 @@ subroutine ozphys_2015_run ( & !ccpp: save input oz in ozi ozi = oz + gravi=1.0/con_g colo3(:,levs+1) = 0.0 coloz(:,levs+1) = 0.0 From 23294fb7f733d5239cc5f24f010d51cd191ef4bf Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Thu, 20 Sep 2018 15:38:54 -0600 Subject: [PATCH 5/6] codes cleanup --- physics/ozinterp.f90 | 5 -- physics/ozphys_2015.f.ORIG | 108 ------------------------------------- 2 files changed, 113 deletions(-) delete mode 100755 physics/ozphys_2015.f.ORIG diff --git a/physics/ozinterp.f90 b/physics/ozinterp.f90 index dad0ef9fa..d90bd61fa 100644 --- a/physics/ozinterp.f90 +++ b/physics/ozinterp.f90 @@ -67,11 +67,6 @@ SUBROUTINE read_o3data (ntoz, me, master) oz_pres(:) = log(100.0*oz_pres(:)) oz_lat(:) = oz_lat4(:) oz_time(:) = oz_time4(:) - if (me == master) then - write(*,*) ' oz_pres(:) = ', oz_pres(:) - write(*,*) ' oz_lat(:) = ', oz_lat(:) - write(*,*) ' oz_time(:) = ', oz_time(:) - endif deallocate (oz_lat4, oz_pres4, oz_time4) !--- read in ozplin which is in order of (lattitudes, ozone levels, coeff number, time) diff --git a/physics/ozphys_2015.f.ORIG b/physics/ozphys_2015.f.ORIG deleted file mode 100755 index 1d7cad57c..000000000 --- a/physics/ozphys_2015.f.ORIG +++ /dev/null @@ -1,108 +0,0 @@ - subroutine ozphys_2015 (ix, im, levs, ko3, dt, ozi, ozo, tin, po3, - & prsl, prdout, pl_coeff, delp, ldiag3d, - & ozp,me) -! -! this code assumes that both prsl and po3 are from bottom to top -! as are all other variables -! This code is specifically for NRL parameterization and -! climatological T and O3 are in location 5 and 6 of prdout array -! June 2015 - Shrinivas Moorthi -! - use machine , only : kind_phys - use physcons, only : grav => con_g - implicit none -! - real, parameter :: gravi=1.0/grav - integer im, ix, levs, ko3, pl_coeff,me - real(kind=kind_phys) ozi(ix,levs), ozo(ix,levs), po3(ko3), - & prsl(ix,levs), tin(ix,levs), delp(ix,levs), - & prdout(ix,ko3,pl_coeff), - & ozp(ix,levs,4), dt -! - integer k,kmax,kmin,l,i,j - logical ldiag3d, flg(im) - real(kind=kind_phys) pmax, pmin, tem, temp - real(kind=kind_phys) wk1(im), wk2(im), wk3(im), prod(im,pl_coeff), - & ozib(im), colo3(im,levs+1), coloz(im,levs+1) -! - colo3(:,levs+1) = 0.0 - coloz(:,levs+1) = 0.0 -! - do l=levs,1,-1 - pmin = 1.0e10 - pmax = -1.0e10 -! - do i=1,im - wk1(i) = log(prsl(i,l)) - pmin = min(wk1(i), pmin) - pmax = max(wk1(i), pmax) - prod(i,:) = 0.0 - enddo - kmax = 1 - kmin = 1 - do k=1,ko3-1 - if (pmin < po3(k)) kmax = k - if (pmax < po3(k)) kmin = k - enddo -! - do k=kmin,kmax - temp = 1.0 / (po3(k) - po3(k+1)) - do i=1,im - flg(i) = .false. - if (wk1(i) < po3(k) .and. wk1(i) >= po3(k+1)) then - flg(i) = .true. - wk2(i) = (wk1(i) - po3(k+1)) * temp - wk3(i) = 1.0 - wk2(i) - endif - enddo - do j=1,pl_coeff - do i=1,im - if (flg(i)) then - prod(i,j) = wk2(i) * prdout(i,k,j) - & + wk3(i) * prdout(i,k+1,j) - endif - enddo - enddo - enddo -! - do j=1,pl_coeff - do i=1,im - if (wk1(i) < po3(ko3)) then - prod(i,j) = prdout(i,ko3,j) - endif - if (wk1(i) >= po3(1)) then - prod(i,j) = prdout(i,1,j) - endif - enddo - enddo - do i=1,im - colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l)*gravi - coloz(i,l) = coloz(i,l+1) + prod(i,6) * delp(i,l)*gravi - prod(i,2) = min(prod(i,2), 0.0) - enddo -! write(1000+me,*) ' colo3=',colo3(1,l),' coloz=',coloz(1,l) -! &,' l=',l - do i=1,im - ozib(i) = ozi(i,l) ! no filling - tem = prod(i,1) - prod(i,2) * prod(i,6) - & + prod(i,3) * (tin(i,l) - prod(i,5)) - & + prod(i,4) * (colo3(i,l)-coloz(i,l)) - -! if (me .eq. 0) print *,'ozphys_2015 tem=',tem,' prod=',prod(i,:) -! &,' ozib=',ozib(i),' l=',l,' tin=',tin(i,l),'colo3=',colo3(i,l+1) - - ozo(i,l) = (ozib(i) + tem*dt) / (1.0 - prod(i,2)*dt) - enddo - if (ldiag3d) then ! ozone change diagnostics - do i=1,im - ozp(i,l,1) = ozp(i,l,1) + (prod(i,1)-prod(i,2)*prod(i,6))*dt - ozp(i,l,2) = ozp(i,l,2) + (ozo(i,l) - ozib(i)) - ozp(i,l,3) = ozp(i,l,3) + prod(i,3)*(tin(i,l)-prod(i,5))*dt - ozp(i,l,4) = ozp(i,l,4) + prod(i,4) - & * (colo3(i,l)-coloz(i,l))*dt - enddo - endif - enddo ! vertical loop -! - return - end From 1935234188592f3fd30d000062d9ad8c33895233 Mon Sep 17 00:00:00 2001 From: "Man.Zhang" Date: Mon, 24 Sep 2018 16:07:38 -0600 Subject: [PATCH 6/6] calculate gravi out of loop per Grant's suggestion and clean up code. --- physics/ozphys.f | 6 +++--- physics/ozphys_2015.f | 2 -- 2 files changed, 3 insertions(+), 5 deletions(-) diff --git a/physics/ozphys.f b/physics/ozphys.f index 0d69c364a..29ad2d0cf 100644 --- a/physics/ozphys.f +++ b/physics/ozphys.f @@ -64,7 +64,6 @@ subroutine ozphys_run ( & ! as are all other variables ! use machine , only : kind_phys -! use physcons, only : grav => con_g implicit none ! ! Interface variables @@ -77,13 +76,13 @@ subroutine ozphys_run ( & & dt, po3(ko3), prdout(ix,ko3,oz_coeff), & & prsl(ix,levs), tin(ix,levs), delp(ix,levs), & & con_g + real :: gravi logical, intent(in) :: ldiag3d character(len=*), intent(out) :: errmsg integer, intent(out) :: errflg ! ! Local variables -! real, parameter :: gravi=1.0/grav integer k,kmax,kmin,l,i,j logical flg(im) real(kind=kind_phys) pmax, pmin, tem, temp @@ -96,13 +95,14 @@ subroutine ozphys_run ( & ! ! save input oz in ozi ozi = oz + gravi=1.0/con_g ! !> - Calculate vertical integrated column ozone values. if (oz_coeff > 2) then colo3(:,levs+1) = 0.0 do l=levs,1,-1 do i=1,im - colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l) * 1.0/con_g + colo3(i,l) = colo3(i,l+1) + ozi(i,l) * delp(i,l) * gravi enddo enddo endif diff --git a/physics/ozphys_2015.f b/physics/ozphys_2015.f index d1d683e90..c0f63f794 100755 --- a/physics/ozphys_2015.f +++ b/physics/ozphys_2015.f @@ -68,10 +68,8 @@ subroutine ozphys_2015_run ( & ! June 2015 - Shrinivas Moorthi ! use machine , only : kind_phys - !use physcons, only : grav => con_g implicit none ! - !real, parameter :: gravi=1.0/grav real(kind=kind_phys),intent(in) :: con_g real :: gravi integer, intent(in) :: im, ix, levs, ko3, pl_coeff,me