[I wrote this earlier, and sent it, but when I came back it turned out
not to have been sent. It's a bit stale now, but perhaps still useful
for someone.]
When reading this thread, my instinctual reaction to the requirement
was best matched by chris burke's solution. I approve of and recommend
that solution.
As for the speculations about optimizations, keep in mind that J is
currently built on C, which uses whatever the hardware has implemented
for IEEE-754 floating point arithmetic for division (or perhaps
integer arithmetic in some cases involving special code).
The C "function" for this implementation is divDD which is defined on
line 43 of ve.c. That line looks like this:
APFX( divDB, D,D,B, DIVI ) APFX( divDI, D,D,I, DIVI) ANAN(
divDD, D,D,D, DIV)
After macro expansion using gcc -E (I'm not attempting to investigate
if there's any relevant platform specific code), the divDD definition
looks like this:
void divDD(J jt,B b,I m, I n,D*z,D*x,D*y){D u;D v; (_clearfp());
if(1==n) {I i=0,_n=(m); for(;i<_n;i++){*z++=((*x)||(*y)?(*x)/(*y):0);
x++; y++;}} else if(b){I i=0,_n=(m); for(;i<_n;i++){u=*x++; {I
i=0,_n=(n); for(;i<_n;i++){*z++=((u)||(*y)?(u)/(*y):0); y++;}}}} else
{I i=0,_n=(m); for(;i<_n;i++){v=*y++; {I i=0,_n=(n);
for(;i<_n;i++){*z++=((*x)||(v)?(*x)/(v):0); x++;}}}};
{if(0x0001&_clearfp()){jtjsignal(jt,(33)); return ;}}; }
Restructuring that into k&r style, using indent (I borrowed a recipe
from stack overflow and then manually replaced tabs with spaces
because composing in gmail seemed to eat the tabs), I get:
void divDD(J jt, B b, I m, I n, D * z, D * x, D * y) {
D u;
D v;
(_clearfp());
if (1 == n) {
I i = 0, _n = (m);
for (; i < _n; i++) {
*z++ = ((*x) || (*y) ? (*x) / (*y) : 0);
x++;
y++;
}
} else if (b) {
I i = 0, _n = (m);
for (; i < _n; i++) {
u = *x++; {
I i = 0, _n = (n);
for (; i < _n; i++) {
*z++ = ((u) || (*y) ? (u) / (*y) : 0);
y++;
}
}
}
} else {
I i = 0, _n = (m);
for (; i < _n; i++) {
v = *y++; {
I i = 0, _n = (n);
for (; i < _n; i++) {
*z++ = ((*x) || (v) ? (*x) / (v) : 0);
x++;
}
}
}
}; {
if (0x0001 & _clearfp()) {
jtjsignal(jt, (33));
return;
}
};
}
I'm not sure I see anything here which I would classify as an
"optimization", though of course it's traditional to have gcc compile
with optimizations enabled.
That said, taking a step back, on many modern cpus (for example,
intel) an unpredicted branch is an order of magnitude more expensive
than a machine instruction (such as floating point division). At the
hardware level, the machine is optimistically deciding what memory to
fetch before the calculations that control that decision have
completed.
Adding conditional code can also create other problems. Generally you
want to make sure to only use conditional code where it's important
for correctness.
Did I mention that I liked chris burke's solution?
Thanks,
--
Raul
On Sat, May 18, 2013 at 2:13 AM, Robert Knight
<[email protected]> wrote:
> Thanks, Henry. I also defer to Devon's more expert solution for a similar
> approach (which I had not seen while composing my reply to PT).
>
> Your advice about premature optimization strategy is cogent: I can't help
> wonder whether J's low-level arithmetic routines have "common sense"
> optimizations already built in.
>
> For example: Does J's "addition" procedure *automatically* return an
> infinity result when either of the addends is infinity -- rather than
> wasting cycles on the logically needless floating point addition?
>
>
> *-Robert Knight
> *
> *
> *
> *
> *
> On Fri, May 17, 2013 at 10:28 PM, Henry Rich <[email protected]> wrote:
>
>
>> Devon's & Robert's ideas are good. Don't worry about performance until
>> you get something working. Then use the performance monitor to tune.
>> Almost certainly you will be fixing the wrong thing if you try to optimize
>> first.
>>
>> Henry Rich
>>
>>
>> On 5/17/2013 10:13 PM, P T wrote:
>>
>>> Thanks Devon and Robert.
>>>
>>> Sure, we can replace zeros with infinity and vice versa before and after
>>> division. But, can it be avoided? What I am looking for is conditional
>>> division without writing an explicit loop. My intention is to represent an
>>> electrical transmission network with matrices and they can potentially be
>>> large. I am trying to minimizing the number of operations.
>>>
>>> I am reading on sparse matrices in "Learning J" by Roger Stokes and looks
>>> like it can do what I am looking for. But, I am not sure if it actually
>>> does not perform the unnecessary operations (i.e. division by zeros) or
>>> just a display issue.
>>>
>>> x =: 1 $. 6 6 NB. an empty 6x6 matrix
>>> x =: 4 5 6 7 ( 0 0 ; 1 1; 2 2; 3 3) } x NB. insert some data
>>> x
>>> 0 0 │ 4
>>> 1 1 │ 5
>>> 2 2 │ 6
>>> 3 3 │ 7
>>>
>>> y =: 1 $. 6 6
>>> y =: 1 2 3 4 ( 0 0 ; 1 1; 2 2; 3 3) } y
>>> y
>>> 0 0 │ 1
>>> 1 1 │ 2
>>> 2 2 │ 3
>>> 3 3 │ 4
>>>
>>> x%y
>>> 0 0 │ 4
>>> 1 1 │ 2.5
>>> 2 2 │ 2
>>> 3 3 │ 1.75
>>>
>>> Thanks,
>>> PT
>>>
>>>
>>>
>>> On Fri, May 17, 2013 at 6:45 PM, Robert Knight <[email protected]
>>> >**wrote:
>>>
>>> *PT-*
>>>>
>>>> To modify the (infinite) "zero-division" result from *infinity* to
>>>> *zero*...
>>>>
>>>> How about adding infinity to the divisor's zero-elements?
>>>>
>>>> *z =: x%y+_*y=0*
>>>>
>>>>
>>>> ]x =. 2 2 $ 2
>>>>
>>>> 2 2
>>>> 2 2
>>>>
>>>> ]y =. 2 2 $ i.4
>>>>
>>>> 0 1
>>>> 2 3
>>>>
>>>> ]z =. x%y+_*y=0
>>>>
>>>> 0 2
>>>> 1 0.666667
>>>>
>>>> *-Robert Knight*
>>>> (Also a J-newbie)
>>>>
>>>>
>>>> On Fri, May 17, 2013 at 6:44 PM, P T <[email protected]> wrote:
>>>>
>>>>
>>>> I am learning J (J602) and dividing one table with another. When
>>>>> division
>>>>> by zero occurs, I want the the value to be zero instead of infinity.
>>>>> For
>>>>> example, in the results below, I want the first element to be 0.0
>>>>> instead
>>>>> of _
>>>>>
>>>>> ]x =. 2 2 $ 2
>>>>> 2 2
>>>>> 2 2
>>>>>
>>>>> ]y=. 2 2 $ i.4
>>>>> 0 1
>>>>> 2 3
>>>>>
>>>>> x%y
>>>>> _ 2
>>>>> 1 0.666667
>>>>>
>>>>>
>>>>> May be I can replace all occurrences of _ with 0.0. But, can I avoid
>>>>> this
>>>>> additional step?
>>>>>
>>>>> Thanks,
>>>>> PT
>>>>> ------------------------------**------------------------------**
>>>>> ----------
>>>>> For information about J forums see http://www.jsoftware.com/**
>>>>> forums.htm <http://www.jsoftware.com/forums.htm>
>>>>>
>>>>> ------------------------------**------------------------------**
>>>> ----------
>>>> For information about J forums see
>>>> http://www.jsoftware.com/**forums.htm<http://www.jsoftware.com/forums.htm>
>>>>
>>>> ------------------------------**------------------------------**
>>> ----------
>>> For information about J forums see
>>> http://www.jsoftware.com/**forums.htm<http://www.jsoftware.com/forums.htm>
>>>
>>> ------------------------------**------------------------------**
>> ----------
>> For information about J forums see
>> http://www.jsoftware.com/**forums.htm<http://www.jsoftware.com/forums.htm>
>>
> ----------------------------------------------------------------------
> For information about J forums see http://www.jsoftware.com/forums.htm
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
