[MITgcm-devel] bug: division by zero

Martin Losch Martin.Losch at awi.de
Wed May 30 12:40:44 EDT 2007


Hi there,

do you have an opinion on this? What should/can I do about it?

Martin

On 25 May 2007, at 09:16, Martin Losch wrote:

> Hi,
>
> this is probably a bug:
>
> on our new sx8 I get a division by zeros in mom_calc_visc.F for  
> curvilinear grids on the cubed sphere/llc and variable viscosity.  
> The culprits are recip_dxV and recip_dyU, because they are the only  
> ones that are not exchanged in ini_curvilinear_grid.F
> Why are they not exchanged? I gather from the comments in  
> ini_curvilinear_grid, that for cube sphere exchanges this has not  
> been sorted out and some sort of hack has been applied which is now  
> catching up with me. What can we do?
> 1. use the wrong exchanges so that dxV/dyU and their reciprocals  
> have values other than zero in the overlaps and keep our fingers  
> crossed that no-one will notice (I assume that the viscosity values  
> in the overlaps are not that important (o:)
> 2. modify mom_calc_visc, so that it does not compute anything in  
> the overlaps (will that work?)
>
> solution 1 changes the results of testreport for  
> global_ocean.cs32x15, if I uncomment these two lines in  
> ini_curvilinear_grid.F:
> c     CALL EXCH_Z_3D_RS( dxV, 1, myThid )
> c     CALL EXCH_Z_3D_RS( dyU, 1, myThid )
> I could do this:
> #ifdef TARGET_SX_NEC
>      CALL EXCH_Z_3D_RS( dxV, 1, myThid )
>      CALL EXCH_Z_3D_RS( dyU, 1, myThid )
> #endif
> which will fix the problem for me but as soon as you run into a  
> compiler that is as picky as the sxf90 you'll have the problem again.
>
> (BTW. sxf90 has options to initial all stack/heap with nan. If I do  
> that the model stops in mom_vi_hdissip.F, havent' figured out yet,  
> why.)
>
> What do you think?
>
> Martin
>
> PS. I hate going to new machines, all these problems that you  
> discover ...
>
> _______________________________________________
> MITgcm-devel mailing list
> MITgcm-devel at mitgcm.org
> http://mitgcm.org/mailman/listinfo/mitgcm-devel




More information about the MITgcm-devel mailing list