Henry, the combined form below fails when b is 0. Kip
r
3 : 0
'a b c' =. y
((a %~ ]) , (c % ])) _0.5 * (%: (*:b) - 4*a*c) ((* *) + ]) b
)
r 1 3 2
_2 _1
r 1 2 1
_1 _1
r 1 0 1
0 _
r 1 0 _1
0 __
roots 1 0 1
0j1 0j_1
roots 1 0 _1
1 _1
Sent from my iPad
On Dec 12, 2012, at 6:25 AM, Henry Rich <[email protected]> wrote:
> I have to disagree. That's way cool.
>
> I assumed you used Newton's method to shine up the final values of the roots,
> but I never knew you tried rational results too.
>
> It is possible to produce accurate results for quadratics without such heroic
> measures. You just can't use the quadratic formula like you learned it in
> highschool.
>
> There's the usual quadratic formula:
>
> (2*a) %~ (-b) (+,-) %: (b^2) - 4*a*c
>
> and there's a less-well-known alternative form:
>
> (2*c) % (-b) (-,+) %: (b^2) - 4*a*c
>
> Each version produces two roots. The roots created are more accurate when
> +-(-b) is positive, so you can get one accurate root from each formula.
>
> To calculate, use the combined form
>
> ((a %~ ]) , (c % ])) _0.5 * (%: (*:b) - 4*a*c) ((* *) + ]) b
>
> Accuracy has a coolness all its own.
>
> Henry Rich
>
> On 12/12/2012 3:13 AM, Roger Hui wrote:
>>> ... produce more accurate results in some difficult cases.
>>
>> c=: p. <20$1.5
>> c
>> 3325.26 _44336.8 280799 _1.1232e6 3.18239e6 _6.78911e6 1.13152e7 _1.50869e7
>> 1.63441e7 _1.45281e7 1.0654e7 _6.45695e6 3.22847e6 _1.3245e6 441501 _117734
>> 24527.8 _3847.5 427.5 _30 1
>> p. c
>> ┌─┬───────────────────────────────────────────────────────────────────────────────┐
>> │1│3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2 3r2
>> 3r2 3r2│
>> └─┴───────────────────────────────────────────────────────────────────────────────┘
>>
>>
>>
>> On Tue, Dec 11, 2012 at 11:29 PM, Roger Hui <[email protected]>wrote:
>>
>>> You may or may not know that p. employs some extraordinary measures which
>>> produce more accurate results in some difficult cases. But those
>>> extraordinary measures are not "cool". For example:
>>>
>>> w=: p. <1+i.20 NB. Wilkinson's
>>> polynomial<http://en.wikipedia.org/wiki/Wilkinson_polynomial>
>>> w
>>> 2432902008176640000 _8752948036761600000 13803759753640704000
>>> _12870931245150988800 8037811822645051776 _3599979517947607200
>>> 1206647803780373360 _311333643161390640 63030812099294896
>>> _10142299865511450 1307535010540395 _135585182899530 11310276995381 _7561...
>>>
>>> p. w
>>> ┌─┬──────────────────────────────────────────────────┐
>>> │1│20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1│
>>> └─┴──────────────────────────────────────────────────┘
>>>
>>>
>>>
>>>
>>>
>>> On Tue, Dec 11, 2012 at 3:40 PM, Henry Rich <[email protected]> wrote:
>>>
>>>> None of these cute ways is very accurate in tough cases:
>>>>
>>>> 0j15 ": (2*a) %~ (-b) (+,-) %: (b^2) - 4*a*c [ 'a b c' =. 1e_6 1e6 1e_6
>>>> 0.000000000000000 _1000000000000.000000000000000
>>>>
>>>> But p. does better:
>>>>
>>>> 0j15 ": 1 {:: p. c,b,a
>>>> _1000000000000.000000000000000 _0.000000000001000
>>>>
>>>> Henry Rich
>>>>
>>>>
>>>>
>>>> On 12/11/2012 6:29 PM, Roger Hui wrote:
>>>>
>>>>> There are some cheeky (or is it cheesy?) versions:
>>>>>
>>>>> (2*a) %~ (-b) (+,-) %: (b^2) - 4*a*c NB. Kip Murray
>>>>> (2*a) %~ - b (+,-) %: (b^2) - 4*a*c
>>>>> (+:a) %~ - b (+,-) %: (*:b) - 4*a*c
>>>>> -: a %~ - b (+,-) %: (*:b) - 4*a*c
>>>>>
>>>>>
>>>>>
>>>>> On Tue, Dec 11, 2012 at 11:37 AM, km <[email protected]> wrote:
>>>>>
>>>>> It appears this could be translated into J as the rather cool
>>>>>>
>>>>>> (2*a) %~ (-b) (+,-) %: (b^2) - 4*a*c
>>>>>>
>>>>>> Sent from my iPad
>>>>>>
>>>>>>
>>>>>> On Dec 11, 2012, at 12:59 PM, Roger Hui <[email protected]>
>>>>>> wrote:
>>>>>>
>>>>>> Example from the Iverson and McDonnell *Phrasal
>>>>>>> Forms*<http://www.jsoftware.**com/papers/fork.htm<http://www.jsoftware.com/papers/fork.htm>>paper
>>>>>>> (which
>>>>>>> introduced fork):
>>>>>>>
>>>>>>> (-b)(+,-)√((b*2)-4×a×c)÷2×a
>>>>>>>
>>>>>>> √ is a postulated APL primitive, spelled %: in J.
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On Tue, Dec 11, 2012 at 10:49 AM, km <[email protected]> wrote:
>>>>>>>
>>>>>>> What is the coolest way of programming the quadratic formula in J? We
>>>>>>> are
>>>>>>
>>>>>>> finding the roots of polynomial c + x*(b + x*a) without using p. . I
>>>>>>> offer
>>>>>>
>>>>>>>
>>>>>>>> roots
>>>>>>>> 3 : 0
>>>>>>>> 'a b c' =. y
>>>>>>>> q =. %: (b^2) - 4*a*c
>>>>>>>> (2*a) %~ (-b) + q,-q
>>>>>>>> )
>>>>>>>> roots 1 3 2
>>>>>>>> _1 _2
>>>>>>>> roots 1 0 1
>>>>>>>> 0j1 0j_1
>>>>>>>> roots 1 _2 1
>>>>>>>> 1 1
>>>>>>>>
>>>>>>>> partly as problem definition. I am looking for cool roots verbs!
>>>>>>>>
>>>>>>>> Kip Murray
>>>>>>>>
>>>>>>>> Sent from my iPad
>>>>>>>>
>>>>>>>> ------------------------------**------------------------------**
>>>>>>>> ----------
>>>>>>>> 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<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
----------------------------------------------------------------------
For information about J forums see http://www.jsoftware.com/forums.htm
