[ut.na] NAgMAg v89 #17

krj@na.toronto.edu (Ken Jackson) (07/10/89)

From:	nagmag%ukc.ac.uk@NSFNET-RELAY.AC.UK
Date:	Thu, 6 Jul 89 17:59:59 EDT
Subject: NAgMAg v89 #17

NAgMAg		Thursday,  July 6 1989	Volume 89   Issue 17

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%   The official electronic digest of the NAG Users Association   %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Today's Topics
                          Queries about D03EEF
                          Minimization Problem
                     Multiple Precision Arithmetic
                   Bug in NAG Graphics routine J06HGF

---------------------------

Date: Fri, 30 Jun 1989 12:29:06 PDT
From: Proskurowski <proskuro@castor.usc.edu>
To: nagmag@uk.ac.ukc
Subject: Queries about D03EEF

This spring I gave the following problem as a part of the take home
quiz to my graduate class in numerical PDEs:
	The NAG software package contains a multigrid routine (D03EEF)
	for solving general 2nd order elliptic PDEs on rectangular
regions.  It has two options for approximating first derivatives by
using a) central, or b) forward differences.  The following test
example (from p.8 of the Mark 13 Release NAG Manual) was run on the SUN
computer in double precision:
		-Du+100(u_x+u_y)=f  in [0,1]^2          with
Dirichlet boundary conditions corresponding to the exact solution
u=x^2+y^2.  Here D denotes the Laplacian, and u_x partial wrt x. The
obtained results were:

# of levels	# of iterations		  ||erorr||_2
		a)	   b)		  a)	     b)
   3		33	   7		5e-11	   4.0e-2
   4		14	   8		3e-12	   1.6e-2
   5		9	   8		7e-13	   3.0e-3

Explain the behavior of the rate of convergence and the error as a
function of N=2^(#of levels) in both cases.  Why such extreme
differences in the results occur (note that the method is implemented
correctly and there are no bugs in the program)?


Comments and questions to NAG.

1. Only the result for 3 levels is given in the manual.  Moreover,
given there are actually squares of the 2-norm (additionally, without
normalization), so the numbers read: a) 1.7e-19 and b) 1.28e-1.

2. Who cursorily looking at these rusults would want to use the
routine: one option is super accurate but extremely expensive, the
other fast but gives practically no useful information (compare with
the norm of the solution!)

3. I hope you want NAG to be used not only by numerical analysts who
have time to invest (or students to do the work) to find out that the
routine works well only the test example is ill chosen, especially in
the complete absence of proper explanations.


---------------------------


From: AlBert DeKnuydt <prlb2!kulcs!kulesat!deknuydt@uunet.UU.NET>
Date: Thu, 22 Jun 89 18:11:27 GMT
Subject: Mininization Problem

This MAIL contains a description of a problem on numerical
minimisation and numerical differentiation.  The domain of
the problem is image processing.


1) Function description :

The rather complicated function to minimize depends on 5
variables.  These 5 variables define an image of which a
histogram is calculated. The function result is the following
operation on this histogram. (In fact the entropy)

N(i)   =   number of occurrences of value i
N      =   total number of occurrences

function_result =  -  SUM    [N(i)/N] * blog [(N(i)/N)]
                   all i with
                   N(i) <>  0

One single function evaluation takes about 3 minutes CPU on a
VAX 8530.

2) What we tried up to now :

We tried the following NAG routines.

E04JAF  :  Quasi Newton algorithm using function values only,
           easy-to-use version.

           Problem : accurracy of result nor function
           evaluation not controllable.  Routine doesn't use
           large enough steps towards the minimum.

E04JBF  :  Quasi Newton algorithm using function values only,
           comprehensive (= not "easy-to-use" ?) version.

           Problem : steps too small.  Routine E04HBF, used
           to determine initial step length for making
           difference approximations to the partial
           derivatives of the target function, expects target
           function to be of machine accuracy. Because our
           function is essentially discrete, this is
           questionable.

E04VCF  :  Sequential QP (Quadratic Programming) method,
           using first derivatives.

           Problem is now with the routine to calculate these
           derivatives.  It almost always returns with  a
           diagnostic complaining about accuracy (forward and
           central difference estimates don't agree to half a
           decimal place).  When we check these derivatives
           with E04ZCF, it says it doubts about the
           correctness of the derivatives.  Probably, the
           cause is the slightly discrete character of  the
           target function.

(We use MARK 12 release of the NAG library, if this is of any
 importance)

3)Questions

Anybody has an idea :

           What the reason of the failure of the NAG routines
           might be ? And how to solve this ?

           Of another way to minimize this kind of function ?


eMAIL deknuydt@kulesat.uucp           UUCP
      deknuydt%kulesat.uucp@blekul60  BITNET

B. DeKnuydt & J. Smolders
K.U.Leuven ESAT/MI2
Kardinaal Mercierlaan 94
B-3030 Heverlee-Leuven
B E L G I U M

%% This message was forwarded from the news network. Please reply
%% to nagmag and I'll forward responses to the original senders.
%% Tim

---------------------------

Date:		Wed,  5 JUL 89 10:17:40 GMT
From:		CAROLINE@vax.nag.co.uk
Subject: Multiple Precision Arithmetic

Replies to issue 16
-------------------

re: Multiple Precision Arithmetic
    -----------------------------

The question from Richard Steele about the modulus for large integers
can be easily solved in terms of the add, subtract and multiply he has
already written, together with a "number-of-digits" (nod) function
which should be easy to write.

To calculate A mod B (both assumed positive): 

if A < B 
then
  return A
else 
  NODD := nod(A) - nod(B)
  while NODD >= 0
    C := B * 10**NODD
    while A > 0 
      A := A - C
    A := A + C
    NODD := NODD - 1
  return A

B * 10**NODD, of course, can be achieved by multiplying by 10, NODD times

I'm sure it can be done much more efficiently than this (check and see
how it is done in Computer Algebra packages) but, if a quick and dirty
solution is adequate, try coding the above algorithm in "C".

Mike Richardson

%% I have mailed somebody regarding the availability of a multiple
%% precision package called BIGNUM apparently written in C. I'll
%% summarize when I get a reply. Tim

---------------------------

Date: Thu, 06 Jul 89 12:51:20 +0100
From: drm@ukc.ac.uk
Subject: Bug in NAG Graphics routine J06HGF


    There is a bug in NAG Graphics routine J06HGF, the routine which
draws a perspective view of a three-dimensional histogram.  The bug
appears when you try to label some, but not all of the rows (or columns)
of the histogram, by putting spaces in the array for those rows which
shouldn't be labelled. 

    The easiest way to demonstrate the bug is with the example program
from the manual.  The program does not need to be altered, but the data
file does by changing one of the labels to spaces, as follows.  Change
line 12 of the data file from: 1978197919801981198219831984
                           to: 1978    19801981198219831984
In other words, don't put a label on the second row of the histogram
(labels being four characters long), but do put a label on all the other
rows.  The example program crashes with the following error:

%FOR-F, adjustable array dimension error
%TRACE-F, symbolic stack dump follows
module name     routine name                     line       rel PC    abs PC

                                                           0028A345  0028A345
J06HEY          J06HEY                            216      00000616  0000AB02
J06HEZ          J06HEZ                            381      00000DAB  00009167
J06HGY          J06HGY                            159      00000519  00008309
J06HGZ          J06HGZ                            122      0000044A  00007D9A
J06HGF          J06HGF                             63      00000142  0000773A
BUG             BUG                                75      00000285  00007485
(I had to call the example program something, and BUG seemed appropriate)

    What I was actually trying to do was to put labels on the first and
last rows and columns of the histogram, but to leave the intervening
labels blank.  The obvious solution is to set up the arrays for labels
on each row and column, with only the first and last elements given
labels, the others were filled in with spaces.  But this causes the
routine to crash... 

    As I can't do what I want to, and this routine puts in the text at
an angle anyway, for the time being I'll put in my axis labels another
way.  Anyone know any good ways of putting text in vertically, in the
same way as it appears on an axis label?  Its easy using calls directly
to GKS, but the only way I can see of doing it in NAG Graphics is to use
the axis labelling routine and to keep moving the viewport around until
the text label ends up in the right place.  This technique works, but it
takes a lot of trial and error getting the viewport in the right place!!

David Morse,
Computing Lab.,
University of Kent at Canterbury.

P.S. You'd never believe the trouble I had sending this missive to NAGMAG,
the machine refused to acknowledge its existence!!

%% Why didn't you ask a local guru :-) Tim


---------------------------

%%   For further information about the NAG Users Association please contact:
%% Janet Bentley, Administrator NAGUA,
%% Shore Lane Farm, Blackstone Edge Old Road,
%% LITTLEBOROUGH, Lancashire, OL15 0LQ, UK.

%%
%%   Replies or submissions to          nagmag@uk.ac.ukc
%%   Distribution changes to    nagmag-request@uk.ac.ukc
%%
%%   END OF ISSUE 

-----------------------------

Reposted by

-- 
Kenneth R. Jackson,            krj@na.toronto.edu   (on Internet, CSNet, 
Computer Science Dept.,                              ARPAnet, BITNET)
University of Toronto,         krj@na.utoronto.ca   (CDNnet and other 
Toronto, Canada  M5S 1A4                             X.400 nets (Europe))
(Phone: 416-978-7075)          ...!{uunet,pyramid,watmath,ubc-cs}!utai!krj