Re: Strange behavior in algebric computations

Mateus,

I've tracked down the source of the problem, and I think it's a bug
on our end, which I've reported.

The issue is that when dh=2114, and the two v's equal -10.7 and -4.8,
then this part of the calcuation is equal to exactly -32767, which is
the missing value for "dh":

     dh * ( v(c,t,k-1,j,i) + v(c,t,k,j,i) )

For some reason this makes the rest of the calculation missing.

When you changed it so that it was multiplying by 0.5 first, then 0.5
got multiplied to dh before the rest of the calculation was done, so
you didn't end up with -32767.

You've already found one work-around. Another one is to change the
missing value of "dh" to something other than -32767, like 1e20.

I'll let you know if I hear about a fix, or whether this is a not a bug.

--Mary

On Oct 1, 2009, at 10:18 AM, Mateus Teixeira wrote:

> Hi Mary,
>
> Here are the results of your codes:
>
> ncl 0> su = 0.
> ncl 1> dh = 1271
> ncl 2> v1 = 0.699997
> ncl 3> v2 = -4.8
> ncl 4> su1 = su + dh * (v1 + v2) / 2
> ncl 5> su2 = su + 0.5 * dh * (v1 + v2)
> ncl 6> print("su1 = " + su1)
> (0) su1 = -2605.55
> ncl 7> print("su2 = " + su2)
> (0) su2 = -2605.55
>
> ncl 10> su = -2605.55
> ncl 11> dh = 2114
> ncl 12> v1 = -4.8
> ncl 13> v2 = -10.7
> ncl 14> su1 = su + dh * (v1 + v2) / 2
> ncl 15> su2 = su + 0.5 * dh * (v1 + v2)
> ncl 16> print("su1 = " + su1)
> (0) su1 = -18989.1
> ncl 17> print("su2 = " + su2)
> (0) su2 = -18989.1
>
> As I mentioned in the previous email, with both formulae the result is
> the same, when computed by hand, but in the loop it goes to a missing
> value.
>
> And here are the outputs from printVarSummary for all variables used
> in the formula:
>
>
> Variable: dh
> Type: float
> Total Size: 4 bytes
> 1 values
> Number of Dimensions: 1
> Dimensions and sizes: [1]
> Coordinates:
> Number Of Attributes: 1
> _FillValue : -32767
>
> Variable: sdh
> Type: float
> Total Size: 4 bytes
> 1 values
> Number of Dimensions: 1
> Dimensions and sizes: [1]
> Coordinates:
> Number Of Attributes: 1
> _FillValue : -32767
>
>
> Variable: sv
> Type: float
> Total Size: 4 bytes
> 1 values
> Number of Dimensions: 1
> Dimensions and sizes: [1]
> Coordinates:
> Number Of Attributes: 1
> _FillValue : -32767
>
>
> Variable: v (subsection)
> Type: float
> Total Size: 4 bytes
> 1 values
> Number of Dimensions: 1
> Dimensions and sizes: [1]
> Coordinates:
> Number Of Attributes: 8
> lon : 302.5
> lat : -35
> nivel : 250
> tempo : 1
> padrao : 1
> _FillValue : -32767
> long_name : vento meridional
> units : m/s
>
>
> Variable: hgt (subsection)
> Type: float
> Total Size: 4 bytes
> 1 values
> Number of Dimensions: 1
> Dimensions and sizes: [1]
> Coordinates:
> Number Of Attributes: 8
> lon : 302.5
> lat : -35
> nivel : 300
> tempo : 1
> padrao : 1
> _FillValue : -32767
> long_name : altura geopotencial
> units : m
>
>
> And here is the loop
>
> sdh = 0.
> su = 0.
> sv = 0.
> dh = 0.
> do k = k10000, k3000
> dh = hgt(c,t,k-1,j,i) - hgt(c,t,k,j,i)
> sdh = sdh + dh
> su = su + 0.5 * dh * ( u(c,t,k-1,j,i) + u(c,t,k,j,i) )
> sv = sv + 0.5 * dh * ( v(c,t,k-1,j,i) + v(c,t,k,j,i) )
> end do
>
> And all these lines are inside other loops for c, t, i, and j. I
> searched for missing values in the data (hgt, u, and v), founding
> nothing. It is very strange. I thought it was due a error
> accumulations, but the variables involved are equaled to zero before k
> loop.
>
> Best regards,
>
> Mateus
>
>
>
> 2009/10/1 Mary Haley <haley_at_ucar.edu>:
>> Hi Mateus,
>>
>> I was unable to reproduce your problem. Can you try this code and
>> let me
>> know what you get:
>>
>> su = 0.
>> dh = 1271
>> v1 = 0.699997
>> v2 = -4.8
>> su1 = su + dh * (v1 + v2) / 2
>> su2 = su + 0.5 * dh * (v1 + v2)
>> print("su1 = " + su1)
>> print("su2 = " + su2)
>>
>> su = -2605.55
>> dh = 2114
>> v1 = -4.8
>> v2 = -10.7
>> su1 = su + dh * (v1 + v2) / 2
>> su2 = su + 0.5 * dh * (v1 + v2)
>> print("su1 = " + su1)
>> print("su2 = " + su2)
>>
>> My first thought is that this has something to do with the types
>> of your
>> variables, and the order of which the arithmetic is being done
>> on them.
>>
>> Can you do a "printVarSummary" on all your variables and send me
>> the output
>> so I have the correct types?
>>
>> --Mary
>>
>>
>> On Oct 1, 2009, at 8:13 AM, Mateus Teixeira wrote:
>>
>>> Dear NCL users,
>>>
>>> I'm doing some algebric computations that are having strange
>>> behavior
>>> depending on way the formula is written. When written this way
>>>
>>> su = su + dh * ( u(c,t,k-1,j,i) + u(c,t,k,j,i) ) / 2
>>>
>>> I'm getting a missing value (the formula is inside a do loop), as
>>> shown
>>> below:
>>>
>>> (0) sv(before)=0; dh=1271; v(k-1)=0.699997; v(k)=-4.8
>>> (0) sv(after)=-2605.55
>>> (0) sv(before)=-2605.55; dh=2114; v(k-1)=-4.8; v(k)=-10.7
>>> (0) sv(after)=-32767 <= MISSING VALUE
>>>
>>> But when the formula is written as
>>>
>>> su = su + 0.5 * dh * ( u(c,t,k-1,j,i) + u(c,t,k,j,i) )
>>>
>>> I haven't a missingl value, as shown below for the same stop
>>> point showed
>>> above
>>>
>>> (0) sv(before)=0; dh=1271; v(k-1)=0.699997; v(k)=-4.8
>>> (0) sv(after)=-2605.55
>>> (0) sv(before)=-2605.55; dh=2114; v(k-1)=-4.8; v(k)=-10.7
>>> (0) sv(after)=-18989.1 <=CORRECT VALUE
>>>
>>> I made these same computations by hand in NCL with both formulae and
>>> get the correct value given by the second way of the formula.
>>>
>>> Can anyone imagine what is causing this? I thought that the float
>>> representation of the 1/2 could be the source of the error, but the
>>> difference between missing value and the correct value is not
>>> little.
>>>
>>> Best regards,
>>>
>>> --
>>> Mateus da Silva Teixeira
>>> Meteorologista
>>> Instituto de Pesquisas Meteorológicas - IPMet
>>> Universidade Estadual Paulista - UNESP
>>> Av. Luis Edmundo Carrijo Coube, 14-01 - Cx. Postal 281 - CEP
>>> 17033-360
>>> Bauru - SP - Brasil
>>> Fone: +55 14 3103-6030 / fax: 3203-3649
>>>
>>> Registered Linux User #466740 (http://counter.li.org/)
>>> _______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>
>
>
> --
> Mateus da Silva Teixeira
> Meteorologista
> Instituto de Pesquisas Meteorológicas - IPMet
> Universidade Estadual Paulista - UNESP
> Av. Luis Edmundo Carrijo Coube, 14-01 - Cx. Postal 281 - CEP 17033-360
> Bauru - SP - Brasil
> Fone: +55 14 3103-6030 / fax: 3203-3649
>
> Registered Linux User #466740 (http://counter.li.org/)
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Oct 01 2009 - 15:17:42 MDT