Skip to content

Fix Tsnice calculation #542

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
eclare108213 merged 4 commits into from
Nov 3, 2025
Merged

Fix Tsnice calculation #542

Changes from 1 commit
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
e9d626f
Fix Tsnice calculation
dabail10 Sep 19, 2025
1b47cf1
Fix some bug with Tsnice
dabail10 Sep 24, 2025
453743d
Add zero initialization for Tsnice in driver
dabail10 Sep 24, 2025
ba7a2d3
Fix the computation of Tsnice to be a straight line
dabail10 Oct 30, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions columnphysics/icepack_therm_vertical.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -378,9 +378,9 @@ subroutine thermo_vertical (dt, aicen, &
Tsnice = c0
Copy link
Contributor

@anton-seaice anton-seaice Sep 24, 2025 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this Tsnice = c0 is in the wrong location.

thermo_vertical is called for once for each thickness category, which means Tsnice=0 gets reset every time and the result would be just the temperature for the ncat = 5

Copy link
Contributor

@eclare108213 eclare108213 Sep 24, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

agreed

Copy link
Contributor Author

@dabail10 dabail10 Sep 24, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Good catch. Fixed.

Copy link
Contributor Author

@dabail10 dabail10 Sep 24, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Tsnice is now initialized in icedrv_flux.F90 and ice_flux.F90.

if ((hslyr+hilyr) > puny) then
if (hslyr > puny) then
Tsnice = (hslyr*zTsn(nslyr) + hilyr*zTin(1)) / (hslyr+hilyr)
Tsnice = Tsnice + aicen*((hslyr*zTsn(nslyr) + hilyr*zTin(1)) / (hslyr+hilyr))
Copy link
Contributor

@NickSzapiro-NOAA NickSzapiro-NOAA Sep 22, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can I check that I understand? Does this aicen factor get "canceled" in ice_history in:

! Only average for timesteps when ice present
if (avail_hist_fields(n)%avg_ice_present) then
   do j = jlo, jhi
   do i = ilo, ihi
      if (tmask(i,j,iblk)) then
            a2D(i,j,n,iblk) = &
            a2D(i,j,n,iblk)*avgct(ns)*ravgip(i,j)
      endif

Copy link
Contributor Author

@dabail10 dabail10 Sep 22, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No, this is aggregating the category Tsnice into the grid cell Tsnice. The time averaging in ice_history.F90 will come in a CICE PR!

Copy link
Contributor

@NickSzapiro-NOAA NickSzapiro-NOAA Sep 23, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks! Then isn't a Tsnice = Tsnice/sum(aicen) factor missing? For example, for single category aicen almost puny

Copy link
Contributor Author

@dabail10 dabail10 Sep 23, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry this is so confusing. I should really call this:

Tsnice = sum(aicen(n)*Tsnicen(n)) for n=1,5

this is done for every point in space within icepack. The history averaging in CICE is then:

sum(aice(t)*Tsnice(t)) / sum(aice(t)), for t over the averaging period. The CICE changes are coming in a future PR!

Copy link
Contributor

@anton-seaice anton-seaice Sep 23, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm also confused :)

Tsnice = sum(aicen(n)*Tsnicen(n)) for n=1,5

sum(aicen(n)) is not 1, so Tsnice ends up weighted by aice here right ?

... The history averaging in CICE is then:

sum(aice(t)*Tsnice(t)) / sum(aice(t)), for t over the averaging period. The CICE changes are coming in a future PR!

Doesn't this mean its weighted by aice twice ?

Copy link
Contributor

@anton-seaice anton-seaice Sep 24, 2025 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh - maybe this line would be correct now:

worka(i,j) = aice(i,j,iblk)*(Tsnice(i,j,iblk)/aice_init(i,j,iblk)+Tffresh)

https://github.com/CICE-Consortium/CICE/blob/b01032506d71d86ca51cef8c3c28c1ab87532833/cicecore/cicedyn/analysis/ice_history.F90#L2768

the divide by aice_init corrects it for the aice weighting applied during thermo vertical

Copy link
Contributor

@eclare108213 eclare108213 Sep 24, 2025 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

No, this is aggregating the category Tsnice into the grid cell Tsnice.

This is true only if Tsnice=0 for the open water area. Assuming that, then the grid cell average Tavg = sum(aicen(n)*Tsnice(n))/sum(aicen(n)) =sum(aicen(n)*Tsnice(n))/1 where the sums are for n=0,ncat mathematically but only need to be taken for n=1,ncat in the code. Confusing, yes, but possibly okay here. (Explicit code comments to this effect would be helpful.) A better assumption for Tsnice for open water might be the ocean temperature. On the other hand, if the data request is for Tsnice only over the ice-covered area, then it had better be
sum(aicen(n)*Tsnice(n))/sum(aicen(n)) for n=1,ncat, where aicen might actually be aicen_init at this point in the code.

I'm generally suspicious of any 'mixed' ice area scaling using both aice and aice_init. Those will need to be considered on a case-by-case basis, in my opinion.

else
Tsnice = Tsf
Tsnice = Tsnice + aicen*Tsf
endif
endif

Expand Down