Skip to content
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

Jump to bottom

Undefined behaviour in calculate_Tin_from_qin #482

Closed
kieranricardo opened this issue Jan 9, 2024 · 7 comments · Fixed by #502
Closed

Undefined behaviour in calculate_Tin_from_qin #482

kieranricardo opened this issue Jan 9, 2024 · 7 comments · Fixed by #502
Assignees
@apcraig
Labels
Priority: High Type: Bug

Comments

@kieranricardo
Copy link

kieranricardo commented Jan 9, 2024 •
edited
Loading

Hi there,

I've come across some undefined behaviour in calculate_Tin_from_qin. Depending on the compiler and compilation options, when the sqrt argument is negative the following will either yield a NaN or Tmltk. Replacing min with ieee_min_num fixes this for my use case, and will consistently return Tmltk when the sqrt argument is negative, independent of the compiler.

https://github.com/CICE-Consortium/Icepack/blob/f6ff8f7c4d4cb6feabe3651b13204cf43fc948e3/columnphysics/icepack_therm_shared.F90#L80C1-L81

Tin =  min((-bb1 - sqrt(bb1*bb1 - c4*aa1*cc1)) /  &
                         (c2*aa1),Tmltk)

Happy to submit a PR from my fork https://github.com/kieranricardo/Icepack/tree/bugfix-calculate_Tin_from_qin-undefined-behaviour :)
@eclare108213
Copy link
Contributor

The quantity under the square root probably should not ever be negative. What are the values of aa1, bb1, cc1 when this fails? Is it a roundoff problem or a bigger issue that needs to be looked into?

@dabail10
Copy link
Contributor

dabail10 commented Jul 5, 2024

Any updates on this one?

@anton-seaice
Copy link

Kieran is on leave this week, we resolved why the quantity was negative.

We could just ignore this, or if we think returning Tmltk is better than returning Nan when the first argument is Nan, then make the change?

@dabail10
Copy link
Contributor

dabail10 commented Jul 9, 2024 •
edited
Loading

I am happy to put in a fix for this. Maybe we should just put in a check for the argument in the square root and set it to zero when negative. Then if the result is greater than Tmltk, then you just get that. I suspect this is a roundoff issue as highlighted by @eclare108213 .

@anton-seaice
Copy link

@apcraig - Do you have any experience with ieee_min_num and if it is better than min? I am not sure about compiler support for it?

@dabail10
Copy link
Contributor

I did some looking into this and ieee_min_num is a fortran 2018 feature, and so we need to be careful here. I think the typical approach here is to do the following instead:

Tin =  min((-bb1 - sqrt(max(bb1*bb1 - c4*aa1*cc1),puny)) /  &
                         (c2*aa1),Tmltk)

@anton-seaice
Copy link

@kieranricardo resolved why the quantity was negative some time ago, we think it was an issue with how the coupling was getting set-up.

My view is we shouldn't change anything, because we don't think its valid for the quantity under the square root to be negative, and giving an error helps with debugging.

I guess a more graceful handling would be to check if bb1*bb1 - c4*aa1*cc1 is negative, and abort with a meaningful error if it is.

@apcraig apcraig self-assigned this Oct 2, 2024
apcraig added a commit to apcraig/Icepack that referenced this issue Oct 2, 2024
@apcraig
Update Derecho Cray, rename subroutines, minor cleanup
74ef9cd 
- Update Derecho Cray, add -Rp to -O2 standard optimization to address
SAME_TBS compiler bug that is triggered in limited cases.  Reported
issue to NCAR.  -O2 -Rp changes answers relative to -O2.

- Rename two public Icepack subroutines, but maintain backward compatibility,
      use icepack_therm_shared  , only: icepack_init_thermo => icepack_init_salinity
      use icepack_therm_shared  , only: icepack_init_trcr => icepack_init_enthalpy
Update the calls in Icepack driver, but also confirmed it works fine with the old
names as well as in CICE with the old names.
Closes CICE-Consortium#255

- Add a check and abort for negative values in the sqrt in computation of Tin in
function calculate_Tin_from_qin.
Closes CICE-Consortium#482

- Refactor calls to icepack_aggregate to make them consistent.  This was part of the
testing for the Derecho Cray bug, and decided to keep the implementation.

- Update comments associated with floeshape constant attribution.  Change from
Steele to Rothrock 1984.
Closes CICE-Consortium#479

- Clean up some of the variable declarations in subroutine set_state_var and module
icedrv_state, merge multiple lines to one line and shift to assumed shape arrays
where appropriate.

- Clean up implementation error in icedrv_restart.F90, subroutine restartfile.  This subroutine
was using a parameter, ntrcr, directly from icepack_tracers.  Switched that to a call to
icepack_query_tracer_sizes.
@apcraig apcraig mentioned this issue Oct 2, 2024
Merged
16 tasks
dabail10 pushed a commit to ESCOMP/Icepack that referenced this issue Oct 25, 2024
@apcraig @dabail10
Update Derecho Cray, rename subroutines, minor cleanup (CICE-Consorti…
04ba407 
…um#502)

* Update Derecho Cray, rename subroutines, minor cleanup

Update Derecho Cray, add -Rp to -O2 standard optimization to address SAME_TBS compiler bug that is triggered in limited cases. Reported issue to NCAR. -O2 -Rp changes answers relative to -O2.

Rename two public Icepack subroutines, but maintain backward compatibility. Update the calls in Icepack driver, but also confirmed it works fine with the old names as well as in CICE with the old names.
Closes CICE-Consortium#255 

    use icepack_therm_shared  , only: icepack_init_thermo => icepack_init_salinity 
    use icepack_therm_shared  , only: icepack_init_trcr => icepack_init_enthalpy 

Add a check and abort for negative values in the sqrt in computation of Tin in function calculate_Tin_from_qin.
Closes CICE-Consortium#482

Refactor calls to icepack_aggregate to make them consistent. This was part of the testing for the Derecho Cray bug, and decided to keep the implementation.

Update comments associated with floeshape constant attribution. Change from Steele to Rothrock 1984.
Closes CICE-Consortium#479

Clean up some of the variable declarations in subroutine set_state_var and module icedrv_state, merge multiple lines to one line and shift to assumed shape arrays where appropriate.

Clean up implementation error in icedrv_restart.F90, subroutine restartfile. This subroutine was using a parameter, ntrcr, directly from icepack_tracers. Switched that to a call to icepack_query_tracer_sizes.

Update documentation of kice noting it's use with BL99 and MU71. See CICE-Consortium#447.

Generate updated interface documentation (./icepack.setup --docintfc)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees

@apcraig apcraig

Labels
Priority: High Type: Bug
Projects
None yet
Milestone
No milestone
Development

Successfully merging a pull request may close this issue.

Update Derecho Cray, rename subroutines, minor cleanup apcraig/Icepack
5 participants
@eclare108213 @apcraig @dabail10 @kieranricardo @anton-seaice