Skip to content

Commit

Permalink
Merge pull request #31 from mvertens/feature/UpdateNoresmBranchRebase
Browse files Browse the repository at this point in the history 
Updates  to oslo-aero including performance

- updates oslo-aero (and in particular oslo_aero/src_cam) to be consistent with the changes that will be brought into CAM with an updated branch (feature/UpdateNorESMBranch)
- implements significant performance updates to oslo-aero
- Testing: regression testing was carried out using an updated test suite - BFB
  • Loading branch information
@gold2718
gold2718 authored Aug 28, 2024
2 parents a93b80d + 773d82e commit ba624de
Show file tree
Hide file tree
Showing 29 changed files with 5,552 additions and 7,863 deletions.
267 changes: 134 additions & 133 deletions src/aero_model.F90
View file

Large diffs are not rendered by default.

1,548 changes: 779 additions & 769 deletions src/oslo_aero_aerocom.F90
View file

Large diffs are not rendered by default.

974 changes: 487 additions & 487 deletions src/oslo_aero_aerocom_tables.F90
View file

Large diffs are not rendered by default.

604 changes: 302 additions & 302 deletions src/oslo_aero_aerodry_tables.F90
View file

Large diffs are not rendered by default.

237 changes: 127 additions & 110 deletions src/oslo_aero_coag.F90
View file

Large diffs are not rendered by default.

593 changes: 292 additions & 301 deletions src/oslo_aero_conc.F90
View file

Large diffs are not rendered by default.

271 changes: 143 additions & 128 deletions src/oslo_aero_condtend.F90
View file

Large diffs are not rendered by default.

750 changes: 380 additions & 370 deletions src/oslo_aero_depos.F90
View file

Large diffs are not rendered by default.

94 changes: 47 additions & 47 deletions src/oslo_aero_diurnal_var.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,8 +29,8 @@ subroutine set_diurnal_invariants(invariants,dtc,ncol,lchnk,inv_oh,inv_ho2,id_oh
integer , intent(in) :: id_oh, id_ho2, id_no3
real(r8) , intent(inout) :: invariants(ncol,pver,nfs)

integer :: i ! column index
integer :: k ! height index
integer :: icol ! column index
integer :: ilev ! height index
integer :: iriseset ! sunrise/set flag
integer :: day, mon, yr, jyr ! date stuff
integer :: j ! working var
Expand Down Expand Up @@ -58,14 +58,14 @@ subroutine set_diurnal_invariants(invariants,dtc,ncol,lchnk,inv_oh,inv_ho2,id_oh
! at any rate only valid between 1950 and 2050, so important years e.g. 1850
! is out of boundary

do i=1,ncol
do icol=1,ncol

fdiurn_oxid=1._r8
fdiurn_no3oxid=1._r8

deglat = rlats(i)*180._r8/pi
deglat = rlats(icol)*180._r8/pi
deglat = max( -89.9999_r8, min( +89.9999_r8, deglat ) )
deglon = rlons(i)*180._r8/pi
deglon = rlons(icol)*180._r8/pi

! get sunrise and sunset times in UTC hours
call sunrisesetxx( deglon, deglat, jyr, mon, day, iriseset, trise, tset, solardec )
Expand All @@ -83,7 +83,7 @@ subroutine set_diurnal_invariants(invariants,dtc,ncol,lchnk,inv_oh,inv_ho2,id_oh
end if
else
trise = 0._r8
if (abs(deglat+solardec) .ge. 90._r8) then
if (abs(deglat+solardec) >= 90._r8) then
tset = 1._r8
else
tset = 0._r8
Expand All @@ -93,7 +93,7 @@ subroutine set_diurnal_invariants(invariants,dtc,ncol,lchnk,inv_oh,inv_ho2,id_oh

! if all day or all night (or very close to it), set fdiurn = 1.0
! Also in periods with all night, we put the mean value for all night steps
if ((tlight .ge. 0.99_r8) .or. (tlight .le. 0.01_r8)) then
if ((tlight >= 0.99_r8) .or. (tlight <= 0.01_r8)) then
fdiurn_oxid = 1._r8
fdiurn_no3oxid = 1._r8 !++IH
! otherwise determine overlap between current timestep and daylight times
Expand Down Expand Up @@ -129,22 +129,22 @@ subroutine set_diurnal_invariants(invariants,dtc,ncol,lchnk,inv_oh,inv_ho2,id_oh
end if

if (inv_oh) then
do k=1,pver
invariants(i,k,id_oh)=invariants(i,k,id_oh)*fdiurn_oxid
do ilev=1,pver
invariants(icol,ilev,id_oh)=invariants(icol,ilev,id_oh)*fdiurn_oxid
end do
end if
if (inv_ho2) then
do k=1,pver
invariants(i,k,id_ho2)=invariants(i,k,id_ho2)*fdiurn_oxid
do ilev=1,pver
invariants(icol,ilev,id_ho2)=invariants(icol,ilev,id_ho2)*fdiurn_oxid
end do
end if
if (inv_no3) then
do k=1,pver
invariants(i,k,id_no3)=invariants(i,k,id_no3)*fdiurn_no3oxid
do ilev=1,pver
invariants(icol,ilev,id_no3)=invariants(icol,ilev,id_no3)*fdiurn_no3oxid
end do
end if

end do ! i= 1,ncol
end do ! icol= 1,ncol
end subroutine set_diurnal_invariants


Expand Down Expand Up @@ -184,23 +184,23 @@ subroutine sunrisesetxx( xlong, ylat, iyear, imonth, iday, &
! local
real(r8) :: sunrise, sunset, ap_dec
real(r8) :: xlongb
integer :: iriseset,i
integer :: iriseset

! need xlong between -180 and +180
xlongb = xlong

if (xlongb .lt. -180.) then
if (xlongb < -180.) then
xlongb = xlongb + 360._r8
else if (xlongb .gt. 180._r8) then
else if (xlongb > 180._r8) then
xlongb = xlongb - 360._r8
end if

call srisesetxx( iyear, imonth, iday, ylat, xlongb, iriseset,sunrise, sunset, ap_dec)

iflag = iriseset
if (iflag .eq. 0) then
if (iflag == 0) then
iflag = 1
if (abs(sunrise+100_r8) .le. 0.01_r8) iflag = 0
if (abs(sunrise+100_r8) <= 0.01_r8) iflag = 0
end if
trise = sunrise
tset = sunset
Expand Down Expand Up @@ -358,49 +358,49 @@ subroutine srisesetxx(iyear, month, iday, rlat, rlong, iriseset,sunrise, sunset,
iriseset = 0

! check latitude, longitude, dates for proper range before calculating dates.
if (((rlat .lt. -90._r8) .or. (rlat .gt. 90._r8)) .or. &
((rlong .lt. -180._r8) .or. (rlong .gt. 180._r8))) then
if (((rlat < -90._r8) .or. (rlat > 90._r8)) .or. &
((rlong < -180._r8) .or. (rlong > 180._r8))) then
iriseset = -1
return
end if

! Year assumed to be betweeen 1950 and 2049. As the model is outside these
! boundary in many cases. year 2000 is assumed for this version of the
! model
! if (iyear .lt. 1950 .or. iyear .gt. 2049) then
! if (iyear < 1950 .or. iyear > 2049) then
! iriseset = -1
! return
! end if
! if (((month .lt. 1) .or. (month .gt. 12)) .or. &
! ((iday .lt. 0) .or. (iday .gt. 32))) then
! if (((month < 1) .or. (month > 12)) .or. &
! ((iday < 0) .or. (iday > 32))) then
! iriseset = -1
! return
! end if
! determine julian day number

! there is no year 0 in the Gregorian calendar and the leap year cycle
! changes for earlier years.
! if (iyear .lt. 1) then
! if (iyear < 1) then
! iriseset = -1
! return
! end if
! leap years are divisible by 4, except for centurial years not divisible by 400.

! year = real (iyear)
! if ((amod(year,4.) .eq. 0.0) .and. (amod(year,100.) .ne. 0.0)) &
! if ((amod(year,4.) == 0.0) .and. (amod(year,100.) /= 0.0)) &
! leapyr = 1
! if(amod(year,400.) .eq. 0.0) leapyr = 1
! if(amod(year,400.) == 0.0) leapyr = 1

jday = iimonth(month) + iday
! if ((leapyr .eq. 1) .and. (month .gt. 2)) jday = jday + 1
! if ((leapyr == 1) .and. (month > 2)) jday = jday + 1

! construct Julian centuries since J2000 at 0 hours UT of date,
! days.fraction since J2000, and UT hours.
delta_years = iyear - 2000._r8

! delta_days is days from 2000/01/00 (1900's are negative).
delta_days = delta_years * 365._r8 + delta_years / 4._r8 + jday
if (iyear .gt. 2000) delta_days = delta_days + 1._r8
if (iyear > 2000) delta_days = delta_days + 1._r8

! J2000 is 2000/01/01.5
days_j2000 = delta_days - 1.5_r8
Expand All @@ -419,45 +419,45 @@ subroutine srisesetxx(iyear, month, iday, rlat, rlong, iriseset,sunrise, sunset,

! tangent of ecliptic_long separated into sine and cosine parts for ap_ra.
f_ap_ra = atan2(cos(mean_obliquity) * sin(ecliptic_long), cos(ecliptic_long))

! change range of ap_ra from -pi -> pi to 0 -> 2 pi.
if (f_ap_ra .lt. 0.0) f_ap_ra = f_ap_ra + twopi
if (f_ap_ra < 0.0) f_ap_ra = f_ap_ra + twopi

! put ap_ra in the range 0 -> 24 hours.
ap_ra = (f_ap_ra / twopi - int(f_ap_ra /twopi)) * 24._r8
ap_dec = asin(sin(mean_obliquity) * sin(ecliptic_long))

! calculate local mean sidereal time.
! A. A. 1990, B6-B7.
! horner's method of polynomial exponent expansion used for gmst0h.
f_gmst0h = 24110.54841_r8 + cent_j2000 * (8640184.812866_r8 &
+cent_j2000 * (0.093104_r8 - cent_j2000 * 6.2e-6_r8))

! convert gmst0h from seconds to hours and put in the range 0 -> 24.
! 24 hours = 86400 seconds
gmst0h = (f_gmst0h / 86400._r8 - int(f_gmst0h / 86400._r8)) * 24._r8
if (gmst0h .lt. 0._r8) gmst0h = gmst0h + 24._r8
if (gmst0h < 0._r8) gmst0h = gmst0h + 24._r8

! convert latitude to radians.
rlat_r = rlat * deg_rad

! avoid tangent overflow at +-90 degrees.
! 1.57079615 radians is equal to 89.99999 degrees.
if (abs(rlat_r) .lt. 1.57079615_r8) then
if (abs(rlat_r) < 1.57079615_r8) then
tan_lat = tan(rlat_r)
else
tan_lat = 6.0e6_r8
end if
if (abs(ap_dec) .lt. 1.57079615_r8) then
if (abs(ap_dec) < 1.57079615_r8) then
tan_dec = tan(ap_dec)
else
tan_dec = 6.0e6_r8
end if

! compute UTs of sunrise and sunset.
! A. A. 1990, A12.
tangterm = tan_lat * tan_dec
if (abs(tangterm) .gt. 1.0_r8) then
if (abs(tangterm) > 1.0_r8) then
sunrise = -100._r8
sunset = -100._r8
else
Expand All @@ -470,19 +470,19 @@ subroutine srisesetxx(iyear, month, iday, rlat, rlong, iriseset,sunrise, sunset,
! put sunrise and sunset in the range 0 to 24 hours.
!ec inserted following statement since in some latitudes timeterm
!ec minus tangterm is less than -25
if (sunrise .le. -24._r8) sunrise = sunrise + 48._r8
if (sunrise .lt. 0._r8) sunrise = sunrise + 24._r8
if (sunrise .ge. 24._r8) sunrise = sunrise - 24._r8
if (sunset .lt. 0._r8) sunset = sunset + 24._r8
if (sunset .ge. 24._r8) sunset = sunset - 24._r8
if (sunrise <= -24._r8) sunrise = sunrise + 48._r8
if (sunrise < 0._r8 ) sunrise = sunrise + 24._r8
if (sunrise >= 24._r8 ) sunrise = sunrise - 24._r8
if (sunset < 0._r8 ) sunset = sunset + 24._r8
if (sunset >= 24._r8 ) sunset = sunset - 24._r8

! mean sidereal day is 0.99727 mean solar days.
sunrise = sunrise * 0.99727_r8
sunset = sunset * 0.99727_r8
end if
! convert ap_dec to degrees.
ap_dec = ap_dec * rad_deg
return

end subroutine srisesetxx

end module oslo_aero_diurnal_var
28 changes: 14 additions & 14 deletions src/oslo_aero_dust.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,7 @@ end subroutine oslo_aero_dust_readnl
subroutine oslo_aero_dust_init()

! local variables
integer :: i
integer :: imode

call soil_erod_init( dust_emis_fact, soil_erod_file )

Expand All @@ -97,16 +97,16 @@ subroutine oslo_aero_dust_init()
tracerMap(2) = l_dst_a3

dust_names(:)=" "
do i=1,numberOfDustModes
dust_names(i) = cnst_name(tracerMap(i))
do imode=1,numberOfDustModes
dust_names(imode) = cnst_name(tracerMap(imode))
end do

end subroutine oslo_aero_dust_init

!===============================================================================
subroutine oslo_aero_dust_emis(lchnk, ncol, dstflx, cflx)

!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
! Purpose: Interface to emission of all dusts.
! Notice that the mobilization is calculated in the land model and
! the soil erodibility factor is applied here.
Expand All @@ -115,37 +115,37 @@ subroutine oslo_aero_dust_emis(lchnk, ncol, dstflx, cflx)
! Arguments:
integer , intent(in) :: lchnk
integer , intent(in) :: ncol
real(r8) , intent(in) :: dstflx(pcols,4)
real(r8) , intent(in) :: dstflx(pcols,4)
real(r8) , intent(inout) :: cflx(pcols,pcnst) ! Surface fluxes

! Local variables
integer :: i,n
integer :: icol,imode
real(r8) :: soil_erod_tmp(pcols)
real(r8) :: totalEmissionFlux(pcols)

! Filter away unreasonable values for soil erodibility
! (using low values e.g. gives emissions in greenland..)
where(soil_erodibility(:,lchnk) .lt. 0.1_r8)
where(soil_erodibility(:,lchnk) < 0.1_r8)
soil_erod_tmp(:)=0.0_r8
elsewhere
soil_erod_tmp(:)=soil_erodibility(:,lchnk)
end where

totalEmissionFlux(:) = 0.0_r8
do i=1,ncol
totalEmissionFlux(i) = totalEmissionFlux(i) + sum(dstflx(i,:))
do icol=1,ncol
totalEmissionFlux(icol) = totalEmissionFlux(icol) + sum(dstflx(icol,:))
end do

! Note that following CESM use of "dust_emis_fact", the emissions are
! Note that following CESM use of "dust_emis_fact", the emissions are
! scaled by the INVERSE of the factor!!
! There is another random scale factor of 1.15 there. Adapting the exact
! same formulation as MAM now and tune later
! As of NE-380: Oslo dust emissions are 2/3 of CAM emissions
! gives better AOD close to dust sources

do n = 1,numberOfDustModes
cflx(:ncol, tracerMap(n)) = -1.0_r8*emis_fraction_in_mode(n) &
*totalEmissionFlux(:ncol)*soil_erod_tmp(:ncol)/(dust_emis_fact)*1.15_r8
do imode = 1,numberOfDustModes
cflx(:ncol, tracerMap(imode)) = -1.0_r8*emis_fraction_in_mode(imode) &
*totalEmissionFlux(:ncol)*soil_erod_tmp(:ncol)/(dust_emis_fact)*1.15_r8
end do

end subroutine oslo_aero_dust_emis
Expand Down Expand Up @@ -180,7 +180,7 @@ subroutine soil_erod_init( dust_emis_fact, soil_erod_file )

! read in soil erodibility factors, similar to Zender's boundary conditions

! Get file name.
! Get file name.
call getfil(soil_erod_file, infile, 0)
call cam_pio_openfile (ncid, trim(infile), PIO_NOWRITE)

Expand Down
Loading

0 comments on commit ba624de

Please sign in to comment.