Initial population.

1 files changed, 820 insertions, 0 deletions

diff --git a/src/interp/fortcall.boot.pamphlet b/src/interp/fortcall.boot.pamphlet
new file mode 100644
index 00000000..9513e313
--- /dev/null
+++ b/src/interp/fortcall.boot.pamphlet
@@ -0,0 +1,820 @@
+\documentclass{article}
+\usepackage{axiom}
+\begin{document}
+\title{\$SPAD/src/interp fortcall.boot}
+\author{The Axiom Team}
+\maketitle
+\begin{abstract}
+\end{abstract}
+\eject
+\tableofcontents
+\eject
+\section{License}
+<<license>>=
+-- Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
+-- All rights reserved.
+--
+-- Redistribution and use in source and binary forms, with or without
+-- modification, are permitted provided that the following conditions are
+-- met:
+--
+--     - Redistributions of source code must retain the above copyright
+--       notice, this list of conditions and the following disclaimer.
+--
+--     - Redistributions in binary form must reproduce the above copyright
+--       notice, this list of conditions and the following disclaimer in
+--       the documentation and/or other materials provided with the
+--       distribution.
+--
+--     - Neither the name of The Numerical ALgorithms Group Ltd. nor the
+--       names of its contributors may be used to endorse or promote products
+--       derived from this software without specific prior written permission.
+--
+-- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+-- IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+-- TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+-- PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
+-- OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+-- EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+-- PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+-- PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+-- LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+-- NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+-- SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+@
+<<*>>=
+<<license>>
+
+makeFort(name,args,decls,results,returnType,aspInfo) ==
+  -- Create an executable Fortran file to call a given library function,
+  -- and a stub Axiom function to process its arguments.
+  -- the following is a list of objects for which values need not be
+  -- passed by the user.
+  dummies := [SECOND(u) for u in args | EQUAL(car u,0)]
+  args := [untangle2(u) for u in args] -- lose spad Union representation
+    where untangle2 u ==
+      atom (v := rest(u)) => v
+      first(v)
+  userArgs := [u for u in args | not member(u,dummies)]  -- Temporary
+  decls := [untangle(u) for u in decls] -- lose spad Union representation
+    where untangle u ==
+      [if atom(rest(v)) then rest(v) else _
+        [if atom(w) then w else rest(w) for w in rest(v)] for v in u]
+  makeFort1(name,args,userArgs,dummies,decls,results,returnType,aspInfo)
+
+makeFort1(name,args,userArgs,dummies,decls,results,returnType,aspInfo) ==
+  asps := [first(u) for u in aspInfo]
+  -- Now reorder the arguments so that all the scalars come first, so
+  -- that when we come to deal with arrays we know all the dimensions.
+  scalarArgs := [u for u in args | atom getFortranType(u,decls)]
+  arrayArgs := [u for u in args | not member(u,scalarArgs)]
+  orderedArgs := [:scalarArgs,:arrayArgs]
+  file := if $fortranDirectory then
+    STRCONC($fortranDirectory,"/",STRINGIMAGE name)
+  else
+    STRINGIMAGE name
+  makeFortranFun(name,orderedArgs,args,dummies,decls,results,file,
+                 $fortranDirectory,returnType,asps)
+  makeSpadFun(name,userArgs,orderedArgs,dummies,decls,results,returnType,asps,
+              aspInfo,file)
+  name
+
+makeFortranFun(name,args,fortranArgs,dummies,decls,results,file,dir,
+               returnType,asps) ==
+  -- Create a C file to call the library function, and compile it.
+  fp := MAKE_-OUTSTREAM(STRCONC(file,".c"))
+  writeCFile(name,args,fortranArgs,dummies,decls,results,returnType,asps,fp)
+  if null dir then dir := '"."
+  asps => SYSTEM STRCONC("cc -c ",file,".c ; mv ",file,".o ",dir)
+  SYSTEM STRCONC("cc ",file,".c -o ",file,".spadexe ",$fortranLibraries)
+
+writeCFile(name,args,fortranArgs,dummies,decls,results,returnType,asps,fp) ==
+  WRITE_-LINE('"#include <stdio.h>",fp)
+  WRITE_-LINE('"#include <sys/select.h>",fp)
+  WRITE_-LINE('"#include <rpc/rpc.h>",fp)
+  WRITE_-LINE('"#ifndef NULL",fp)
+  WRITE_-LINE('"#define NULL 0",fp)
+  WRITE_-LINE('"#endif  NULL",fp)
+  WRITE_-LINE('"#define MAX__ARRAY(x) (x ? x :  20000)",fp)
+  WRITE_-LINE('"#define CHECK(x) if (!x) {fprintf(stderr,_"xdr failed_"); exit(1);}",fp)
+  WRITE_-LINE('"void main()",fp)
+  WRITE_-LINE('"{",fp)
+  WRITE_-LINE('"  XDR xdrs;",fp)
+  WRITE_-LINE('"  {",fp)
+  if $addUnderscoreToFortranNames then
+    routineName := STRCONC(name,STRING(95))
+  else
+    routineName := name
+  -- If it is a function then give it somewhere to stick its result:
+  if returnType then
+    returnName := INTERN STRCONC(name,"__result")
+    wl(['"    ",getCType returnType,'" ",returnName,'",",routineName,'"();"],fp)
+  -- print out type declarations for the Fortran parameters, and build an
+  -- ordered list of pairs [<parameter> , <type>]
+  argList := nil
+  for a in args repeat
+    argList := [[a, getCType getFortranType(a,decls)], :argList]
+    printDec(SECOND first argList,a,asps,fp)
+  argList := nreverse argList;
+  -- read in the data
+  WRITE_-LINE('"    xdrstdio__create(&xdrs, stdin, XDR__DECODE);",fp)
+  for a in argList repeat
+    if LISTP SECOND a then writeMalloc(first a,first SECOND a,rest SECOND a,fp)
+    not MEMQ(first a,[:dummies,:asps]) => writeXDR(a,'"&xdrs",fp)
+  -- now call the Library routine.  FORTRAN names may have an underscore
+  -- appended.
+  if returnType then
+    wt(['"    ",returnName,'"="],fp)
+  else
+    wt(['"    "],fp)
+  wt([routineName,'"("],fp)
+  if first fortranArgs then
+    printCName(first fortranArgs,isPointer?(first fortranArgs,decls),asps,fp)
+  for a in rest fortranArgs repeat
+    PRINC('",",fp)
+    printCName(a,isPointer?(a,decls),asps,fp)
+  writeStringLengths(fortranArgs,decls,fp)
+  WRITE_-LINE('");",fp)
+  -- now export the results.
+  WRITE_-LINE('"    xdrstdio__create(&xdrs, stdout, XDR__ENCODE);",fp)
+  if returnType then
+    writeXDR([returnName,getCType returnType],'"&xdrs",fp)
+  for r in results repeat
+    writeXDR([r,getCType getFortranType(r,decls)],'"&xdrs",fp)
+  WRITE_-LINE('"    exit(0);",fp)
+  WRITE_-LINE('"  }",fp)
+  WRITE_-LINE('"}",fp)
+
+writeStringLengths(fortranArgs,decls,fp) ==
+  for a in fortranArgs repeat
+    if isString?(a,decls) then wt(['",&",a,'"__length"],fp)
+
+isString?(u,decls) ==
+  EQUAL(ty := getFortranType(u,decls),"character") or
+    LISTP(ty) and EQUAL(first ty,"character")
+
+isPointer?(u,decls) ==
+  ty := getFortranType(u,decls)
+  LISTP(ty) or ty in ["character","complex","double complex"]
+
+printCName(u,ispointer,asps,fp) ==
+  member(u,asps) =>
+    PRINC(u,fp)
+    if $addUnderscoreToFortranNames then PRINC(STRING(95),fp)
+  if not ispointer then PRINC('"&",fp)
+  PRINC(u,fp)
+
+getFortranType(u,decls) ==
+  -- find u in decls, return the given (Fortran) type.
+  result := nil
+  for d in decls repeat for dec in rest d repeat
+    atom(dec) and dec=u =>
+      return( result := first d )
+    LISTP(dec) and first(dec)=u =>
+      return( result := [first d,:rest dec] )
+  result => result
+  error ['"Undeclared Fortran parameter: ",u]
+
+getCType t ==
+  -- Return the equivalent C type.
+  LISTP(t) =>
+    --[if first(t)="character" then '"char" else getCType first t,:rest t]
+    first(t)="character" => ['"char",:rest t]
+    first(t)="complex" => ['"float",2,:rest t]
+    first(t)="double complex" =>  ['"double",2,:rest t]
+    [getCType first t,:rest t]
+  t="double" => '"double"
+  t="double precision" => '"double"
+  t="integer" => '"int"
+  t="real" => '"float"
+  t="logical" => '"int"
+  t="character" => ['"char",1]
+  t="complex" => ['"float",2] --'"Complex" -- we use our own typedef
+  t="double complex" =>  ['"double",2] --'"DComplex" -- we use our own typedef
+  error ['"Unrecognised Fortran type: ",t]
+
+XDRFun t ==
+  LISTP(ty := SECOND t) =>
+    if first(ty)='"char" then '"wrapstring" else '"array"
+  ty
+
+printDec(type,dec,asps,fp) ==
+  wt(['"    ",if LISTP(type) then first(type) else type,'" "],fp)
+  member(dec,asps) =>
+    if $addUnderscoreToFortranNames then
+      wl([dec,STRING(95),'"();"],fp)
+    else
+      wl([dec,'"();"],fp)
+  LISTP(type) =>
+    wl(['"*",dec,'" = NULL;"],fp)
+    wl(['"    u__int ",dec, '"__length = 0;"],fp)
+  type = '"char" =>
+    wl(['"*",dec,'" = NULL;"],fp)
+  wl([dec, '";"],fp)
+
+writeXDR(v,str,fp) ==
+  -- Generate the calls to the filters which will read from the temp
+  -- file.  The CHECK macro ensures that the translation worked.
+  underscore := STRING CHAR("__:",0) -- to avoid a compiler bug which won't
+                                     -- parse " ... __" properly.
+  wt(['"    CHECK(xdr",underscore, XDRFun(v), '"(", str, '",&", first(v)],fp)
+  if (LISTP (ty :=SECOND v)) and not EQUAL(first ty,'"char") then
+    wt(['",&",first(v),'"__length,MAX__ARRAY(",first(v),'"__length),"],fp)
+    wt(['"sizeof(",first(ty),'"),xdr",underscore,first ty],fp)
+  wl(['"));"],fp)
+
+prefix2Infix(l) ==
+  atom(l) => [l]
+  #l=2 => [first l,"(",:prefix2Infix SECOND l,")"]
+  #l=3 => ["(",:prefix2Infix SECOND l,first l,:prefix2Infix THIRD l,")"]
+  error '"Function in array dimensions with more than two arguments"
+
+writeMalloc(name,type,dims,fp) ==
+  -- Write out a malloc for array arguments
+  -- Need the size as well
+  wl(['"    ",name,'"__length=",prefix2Infix first dims,:[:["*",:prefix2Infix u]
+      for u in rest dims],'";"], fp)
+  type = '"char" =>
+    wl(['"    ",name,'"=(",type," *)malloc((1+",name,
+       '"__length)*sizeof(",type,'"));"],fp)
+  wl(['"    ",name,'"=(",type," *)malloc(",name,
+     '"__length*sizeof(",type,'"));"],fp)
+
+wl (l,fp) ==
+  for u in l repeat PRINC(u,fp)
+  TERPRI(fp)
+
+wt (l,fp) ==
+  for u in l repeat PRINC(u,fp)
+
+-- spadRecordType(v,decs) ==
+--   -- Build a lisp representation of the declaration of a spad record.
+--   -- This will be the returned type of the spad function which calls the
+--   -- Fortran code.
+--   ["Record",:[spadRecordType1(u,decs) for u in v]]
+-- 
+-- spadRecordType1(u,decls) ==
+--   -- Create a list of the form '( |:| u <spadTypeTTT u>)
+--   [":",u,spadTypeTTT getFortranType(u,decls)]
+
+spadTypeTTT u ==
+  -- Return the spad domain equivalent to the given Fortran type.
+  -- Changed by MCD 8/4/94 to reflect correct format for domains in 
+  -- current system.
+  LISTP u =>
+    first(u)="character" => ["String"]
+    first(u)="logical" and #u=2 => ["List",["Boolean"]]
+    first(u)="logical" => ["List",["List",["Boolean"]]]
+    #u=2 => ["Matrix",spadTypeTTT first u]
+    #u=3 => ["Matrix",spadTypeTTT first u]
+    #u=4 => ["ThreeDimensionalMatrix",spadTypeTTT first u]
+    error '"Can only handle one-, two- and three-dimensional matrices"
+  u = "double" => ["DoubleFloat"]
+  u = "double precision" => ["DoubleFloat"]
+  u = "real" => ["DoubleFloat"]
+  u = "integer" => ["Integer"]
+  u = "logical" => ["Boolean"]
+  u = "character" => ["String"]
+  u = "complex" => ["Complex",["DoubleFloat"]]
+  u = "double complex" => ["Complex",["DoubleFloat"]]
+  error ['"Unrecognised Fortran type: ",u]
+
+mkQuote l ==
+ [addQuote(u)for u in l] where
+    addQuote u ==
+      atom u => ['QUOTE,u]
+      ["construct",:[addQuote(v) for v in u]]
+
+makeLispList(l) ==
+  outputList := []
+  for u in l repeat
+    outputList := [:outputList, _
+                  if atom(u) then ['QUOTE,u] else [["$elt","Lisp","construct"],_
+                  :makeLispList(u)]]
+  outputList
+
+makeSpadFun(name,userArgs,args,dummies,decls,results,returnType,asps,aspInfo,
+            file) ==
+  -- Create an interpreter function for the user to call.
+
+  fType := ["List", ["Record" , [":","key","Symbol"], [":","entry","Any"]]]
+
+  -- To make sure the spad interpreter isn't confused:
+  if returnType then
+    returnName := INTERN STRCONC(name,"Result")
+    decls := [[returnType,returnName], :decls]
+    results := [returnName, :results]
+  argNames := [INTERN STRCONC(STRINGIMAGE(u),'"__arg") for u in userArgs]
+  aType := [axiomType(a,decls,asps,aspInfo) for a in userArgs]
+  aspTypes := [SECOND NTH(POSITION(u,userArgs),aType) for u in asps]
+  nilLst := MAKE_-LIST(#args+1)
+  decPar := [["$elt","Lisp","construct"],:makeLispList decls]
+  fargNames := [INTERN STRCONC(STRINGIMAGE(u),'"__arg") for u in args |
+                 not (MEMQ(u,dummies) or MEMQ(u,asps)) ]
+  for u in asps repeat
+    fargNames := delete(INTERN STRCONC(STRINGIMAGE(u),'"__arg"),fargNames)
+  resPar := ["construct",["@",["construct",:fargNames],_
+             ["List",["Any"]]]]
+  call := [["$elt","Lisp","invokeFortran"],STRCONC(file,".spadexe"),_
+           [["$elt","Lisp","construct"],:mkQuote args],_
+           [["$elt","Lisp","construct"],:mkQuote union(asps,dummies)], decPar,_
+           [["$elt","Lisp","construct"],:mkQuote results],resPar]
+  if asps then
+    -- Make a unique(ish) id for asp files
+    aspId := STRCONC(getEnv('"SPADNUM"), GENTEMP('"NAG"))
+    body := ["SEQ",:makeAspGenerators(asps,aspTypes,aspId),_
+             makeCompilation(asps,file,aspId),_
+             ["pretend",call,fType] ]
+  else
+    body := ["pretend",call,fType]
+  interpret ["DEF",[name,:argNames],["Result",:aType],nilLst,_
+             [["$elt","Result","construct"],body]]
+
+stripNil u ==
+  [CAR(u), ["construct",:CADR(u)], if CADDR(u) then "true" else "false"]
+
+makeUnion aspType ==
+  -- The argument is the type of the asp to be generated.  We would like to
+  -- allow the user to be able to provide a fileName as an alternative
+  -- argument, so this builds the Union of aspType and FileName.
+  ["Union",[":","fp",aspType],[":","fn","FileName"]]
+
+axiomType(a,decls,asps,aspInfo) ==
+  a in asps =>
+    entry := first [u for u in aspInfo | first(u) = a]
+    ftc := ["$elt","FortranType","construct"]
+    rc  := ["$elt", _
+             ["Record",[":","key","Symbol"],[":","entry","FortranType"]], _
+            "construct"]
+    makeUnion ["FortranProgram",_
+      a,_
+      CADR(entry),_
+      ["construct",:mkQuote CADDR entry], _
+      [ ["$elt", "SymbolTable","symbolTable"],_
+          ["construct",_
+            :[[rc,first(v),[ftc,:stripNil rest(v)]] for v in CADDDR entry]]_
+    ] ]
+  spadTypeTTT(getFortranType(a,decls))
+
+makeAspGenerators(asps,types,aspId) ==
+-- The code generated here will manipulate the Fortran output stack and write
+-- the asps out as Fortran.
+  [:makeAspGenerators1(u,v,aspId) for u in asps for v in types]
+
+makeAspGenerators1(asp,type,aspId) ==
+  [[["$elt","FOP","pushFortranOutputStack"] ,_
+    ["filename",'"",STRCONC(STRINGIMAGE asp,aspId),'"f"]] , _
+   makeOutputAsFortran INTERN STRCONC(STRINGIMAGE(asp),'"__arg"), _
+   [["$elt","FOP","popFortranOutputStack"]]   _
+  ]
+
+makeOutputAsFortran arg ==
+  ["IF",["case",arg,"fn"],["outputAsFortran",[arg,"fn"]],_
+                          ["outputAsFortran",[arg,"fp"]] ]
+
+makeCompilation(asps,file,aspId) ==
+  [["$elt","Lisp","compileAndLink"],_
+   ["construct",:[STRCONC(STRINGIMAGE a,aspId,'".f") for a in asps]], _
+   $fortranCompilerName,_
+   STRCONC(file,'".o"),_
+   STRCONC(file,'".spadexe"),_
+   $fortranLibraries]
+
+
+compileAndLink(fortFileList,fortCompiler,cFile,outFile,linkerArgs) ==
+  SYSTEM STRCONC (fortCompiler, addSpaces fortFileList,_
+                  cFile, " -o ",outFile," ",linkerArgs)
+
+addSpaces(stringList) ==
+  l := " "
+  for s in stringList repeat l := STRCONC(l,s," ")
+  l
+
+complexRows z ==
+-- Take a list of lists of complexes (i.e. pairs of floats) and
+-- make them look like a Fortran vector!
+  [:[:pair2list(u.i) for u in z] for i in 0..#(z.0)-1]
+
+pair2list u == [car u,cdr u]
+vec2Lists1 u == [ELT(u,i) for i in 0..#u-1]
+vec2Lists u == [vec2Lists1 ELT(u,i) for i in 0..#u-1]
+
+spad2lisp(u) ==
+  -- Turn complexes into arrays of floats
+  first first(u)="Complex" =>
+    makeVector([makeVector([CADR u,CDDR u],'DOUBLE_-FLOAT)],NIL)
+  -- Turn arrays of complexes into arrays of floats so that tarnsposing
+  -- them puts them in the correct fortran order
+  first first(u)="Matrix" and first SECOND first(u) = "Complex" =>
+    makeVector([makeVector(complexRows vec2Lists rest u,'DOUBLE_-FLOAT)],NIL)
+  rest(u)
+
+invokeFortran(objFile,args,dummies,decls,results,actual) ==
+  actual := [spad2lisp(u) for u in first actual]
+  returnedValues := spadify( _
+        fortCall(objFile,prepareData(args,dummies,actual,decls),_
+                 prepareResults(results,args,dummies,actual,decls)),_
+                            results,decls,inFirstNotSecond(args,dummies),actual)
+
+--  -- If there are one or two elements in returnedValues we must return a
+--  -- cons cell, otherwise a vector.  This is to match the internal
+--  -- representation of an Axiom Record.
+--  #returnedValues = 1 => returnedValues
+--  #returnedValues = 2 => CONS(first returnedValues,SECOND returnedValues)
+--  makeVector(returnedValues,nil)
+
+int2Bool u ==
+  -- Return something which looks like an axiom boolean
+  u=1 => "TRUE"
+  NIL
+
+makeResultRecord(name,type,value) ==
+  -- Take an object returned by the NAG routine and make it into an AXIOM
+  -- object of type Record(key:Symbol,entry:Any) for use by Result.
+  CONS(name,CONS(spadTypeTTT type,value))
+
+spadify(l,results,decls,names,actual) ==
+  -- The elements of list l are the output forms returned from the Fortran
+  -- code: integers, floats and vectors.  Return spad forms of these, of 
+  -- type Record(key:Symbol,entry:Any) (for use with the Result domain).
+  SETQ(RESULTS,l)
+  spadForms := nil
+  for i in 0..(#l -1) repeat
+    fort := NTH(i,l)
+    name := NTH(i,results)
+    ty := getFortranType(name,decls)
+    -- Result is a string
+    STRINGP fort =>
+      spadForms := [makeResultRecord(name,ty,fort), :spadForms]
+    -- Result is a Complex Scalar
+    ty in ["double complex" , "complex"] =>
+       spadForms := [makeResultRecord(name,ty, _
+                                     CONS(ELT(fort,0),ELT(fort,1)) ),:spadForms]
+    -- Result is a Complex vector or array
+    LISTP(ty) and first(ty) in ["double complex" , "complex"] =>
+      dims := [getVal(u,names,actual) for u in rest ty]
+      els := nil
+      if #dims=1 then
+        els := [makeVector([CONS(ELT(fort,2*i),ELT(fort,2*i+1)) _
+                for i in 0..(first(dims)-1)],nil)]
+      else if #dims=2 then
+        for r in 0..(first(dims) - 1) repeat
+          innerEls := nil
+          for c in 0..(SECOND(dims) - 1) repeat
+            offset := 2*(c*first(dims)+r)
+            innerEls := [CONS(ELT(fort,offset),ELT(fort,offset+1)),:innerEls]
+          els := [makeVector(NREVERSE innerEls,nil),:els]
+      else
+         error ['"Can't cope with complex output dimensions higher than 2"]
+      spadForms := [makeResultRecord(name,ty,makeVector(NREVERSE els,nil)),
+                    :spadForms]
+    -- Result is a Boolean vector or array
+    LISTP(ty) and first(ty)="logical" and #ty=2 =>
+      dim := getVal(first rest ty,names,actual)
+      spadForms := [makeResultRecord(name,ty,_
+                          [int2Bool ELT(fort,i) for i in 0..dim-1]), :spadForms]
+    LISTP(ty) and first(ty)="logical" =>
+      dims := [getVal(u,names,actual) for u in rest ty]
+      els := nil
+      if #dims=2 then
+        for r in 0..(first(dims) - 1) repeat
+          innerEls := nil
+          for c in 0..(SECOND(dims) - 1) repeat
+            innerEls := [int2Bool ELT(fort,c*first(dims)+r),:innerEls]
+          els := [NREVERSE innerEls,:els]
+      else
+         error ['"Can't cope with logical output dimensions higher than 2"]
+      spadForms := [makeResultRecord(name,ty,NREVERSE els), :spadForms]
+    -- Result is a vector or array
+    VECTORP fort =>
+      dims := [getVal(u,names,actual) for u in rest ty]
+      els := nil
+      -- Check to see whether we are dealing with a dummy (0-dimensional) array.
+      if MEMQ(0,dims) then
+        els := [[]]
+      else if #dims=1 then
+        els := [makeVector([ELT(fort,i) for i in 0..(first(dims)-1)],nil)]
+      else if #dims=2 then
+        for r in 0..(first(dims) - 1) repeat
+          innerEls := nil
+          for c in 0..(SECOND(dims) - 1) repeat
+            innerEls := [ELT(fort,c*first(dims)+r),:innerEls]
+          els := [makeVector(NREVERSE innerEls,nil),:els]
+      else if #dims=3 then
+        iDim := first(dims)
+        jDim := SECOND dims
+        kDim := THIRD dims
+        for r in 0..(iDim - 1) repeat
+          middleEls := nil
+          for c in 0..(jDim - 1) repeat
+            innerEls := nil
+            for p in 0..(kDim - 1) repeat
+              offset := p*jDim + c*kDim + r
+              innerEls := [ELT(fort,offset),:innerEls]
+            middleEls := [makeVector(NREVERSE innerEls,nil),:middleEls]
+          els := [makeVector(NREVERSE middleEls,nil),:els]
+      else
+         error ['"Can't cope with output dimensions higher than 3"]
+      if not MEMQ(0,dims) then els := makeVector(NREVERSE els,nil)
+      spadForms := [makeResultRecord(name,ty,els), :spadForms]
+    -- Result is a Boolean Scalar
+    atom fort and ty="logical" =>
+      spadForms := [makeResultRecord(name,ty,int2Bool fort), :spadForms]
+    -- Result is a Scalar
+    atom fort => 
+      spadForms := [makeResultRecord(name,ty,fort),:spadForms]
+    error ['"Unrecognised output format: ",fort]
+  NREVERSE spadForms
+
+lispType u ==
+  -- Return the lisp type equivalent to the given Fortran type.
+  LISTP u => lispType first u
+  u = "real" => "SHORT-FLOAT"
+  u = "double" => "DOUBLE-FLOAT"
+  u = "double precision" => "DOUBLE-FLOAT"
+  u = "integer" => "FIXNUM"
+  u = "logical" => "BOOLEAN"
+  u = "character" => "CHARACTER"
+  u = "complex" => "SHORT-FLOAT"
+  u = "double complex" => "DOUBLE-FLOAT"
+  error ['"Unrecognised Fortran type: ",u]
+
+getVal(u,names,values) ==
+  -- if u is the i'th element of names, return the i'th element of values,
+  -- otherwise if it is an arithmetic expression evaluate it.
+  NUMBERP(u) => u
+  LISTP(u) => eval [first(u), :[getVal(v,names,values) for v in rest u]]
+  (place := POSITION(u,names)) => NTH(place,values)
+  error ['"No value found for parameter: ",u]
+
+
+prepareData(args,dummies,values,decls) ==
+-- TTT: we don't 
+-- writeData handles all the mess
+   [args,dummies,values,decls]
+
+
+checkForBoolean u ==
+  u = "BOOLEAN" => "FIXNUM"
+  u
+
+prepareResults(results,args,dummies,values,decls) ==
+  -- Create the floating point zeros (boot doesn't like 0.0d0, 0.0D0 etc)
+  shortZero : fluid := COERCE(0.0,'SHORT_-FLOAT)
+  longZero : fluid := COERCE(0.0,'DOUBLE_-FLOAT)
+  data := nil
+  for u in results repeat
+    type := getFortranType(u,decls)
+    data := [defaultValue(type,inFirstNotSecond(args,dummies),values),:data]
+      where defaultValue(type,argNames,actual) ==
+        LISTP(type) and first(type)="character" => MAKE_-STRING(1)
+        LISTP(type) and first(type) in ["complex","double complex"] =>
+          makeVector(  makeList(
+            2*APPLY('_*,[getVal(tt,argNames,actual) for tt in rest(type)]),_
+            if first(type)="complex" then shortZero else longZero),_
+          if first(type)="complex" then "SHORT-FLOAT" else "DOUBLE-FLOAT" )
+        LISTP type => makeVector(_
+          makeList(
+            APPLY('_*,[getVal(tt,argNames,actual) for tt in rest(type)]),_
+            defaultValue(first type,argNames,actual)),_
+          checkForBoolean lispType first(type) )
+        type = "integer" => 0
+        type = "real" => shortZero
+        type = "double" => longZero
+        type = "double precision" => longZero
+        type = "logical" => 0
+        type = "character" => MAKE_-STRING(1)
+        type = "complex" => makeVector([shortZero,shortZero],'SHORT_-FLOAT)
+        type = "double complex" => makeVector([longZero,longZero],'LONG_-FLOAT)
+        error ['"Unrecognised Fortran type: ",type]
+  NREVERSE data
+
+-- TTT this is dead code now
+--	transposeVector(u,type) ==
+--	  -- Take a vector of vectors and return a single vector which is in column
+--	  -- order (i.e. swap from C to Fortran order).
+--	  els  := nil
+--	  rows := CAR ARRAY_-DIMENSIONS(u)-1
+--	  cols := CAR ARRAY_-DIMENSIONS(ELT(u,0))-1
+--	  -- Could be a 3D Matrix
+--	  if VECTORP ELT(ELT(u,0),0) then
+--	    planes := CAR ARRAY_-DIMENSIONS(ELT(ELT(u,0),0))-1
+--	    for k in 0..planes repeat for j in 0..cols repeat for i in 0..rows repeat
+--	      els := [ELT(ELT(ELT(u,i),j),k),:els]
+--	  else
+--	    for j in 0..cols repeat for i in 0..rows repeat
+--	      els := [ELT(ELT(u,i),j),:els]
+--	  makeVector(NREVERSE els,type)
+
+
+writeData(tmpFile,indata) ==
+  -- Write the elements of the list data to a temporary file.  Return the
+  -- name of that file.
+  --
+  str := MAKE_-OUTSTREAM(tmpFile)
+  xstr := xdrOpen(str,true)
+  [args,dummies,values,decls] := indata
+  for v in values repeat
+	-- the two Boolean values
+	v = "T" => 
+		xdrWrite(xstr,1)
+	NULL v =>   
+		xdrWrite(xstr,0)
+	-- characters  
+	STRINGP v => 
+		xdrWrite(xstr,v)
+	-- some array
+	VECTORP v =>  
+		rows := CAR ARRAY_-DIMENSIONS(v)
+		-- is it 2d or more (most likely) ?
+		VECTORP ELT(v,0) =>	
+			cols := CAR ARRAY_-DIMENSIONS(ELT(v,0))
+			-- is it 3d ?
+			VECTORP ELT(ELT(v,0),0) =>
+				planes := CAR ARRAY_-DIMENSIONS(ELT(ELT(v,0),0))
+				-- write 3d array
+				xdrWrite(xstr,rows*cols*planes)
+				for k in 0..planes-1 repeat 
+					for j in 0..cols-1 repeat 
+						for i in 0..rows-1 repeat 
+							xdrWrite(xstr,ELT(ELT(ELT(v,i),j),k))
+			-- write 2d array
+			xdrWrite(xstr,rows*cols)
+			for j in 0..cols-1 repeat
+				for i in 0..rows-1 repeat xdrWrite(xstr,ELT(ELT(v,i),j))
+		-- write 1d array
+		xdrWrite(xstr,rows)
+		for i in 0..rows-1 repeat xdrWrite(xstr,ELT(v,i))
+	-- this is used for lists of booleans apparently in f01
+     	LISTP v => 
+		xdrWrite(xstr,LENGTH v)
+	        for el in v repeat 
+			if el then xdrWrite(xstr,1) else xdrWrite(xstr,0) 
+	-- integers
+	INTEGERP v => 
+		xdrWrite(xstr,v)
+	-- floats
+	FLOATP v => 
+		xdrWrite(xstr,v)
+  SHUT(str)
+  tmpFile
+
+readData(tmpFile,results) ==
+  -- read in the results from tmpFile.  The list results is a list of
+  -- dummy objects of the correct type which will receive the data.
+  str := MAKE_-INSTREAM(tmpFile)
+  xstr := xdrOpen(str,false)
+  results := [xdrRead1(xstr,r) for r in results] where
+    xdrRead1(x,dummy) ==
+      VECTORP(dummy) and ZEROP(LENGTH dummy) => dummy
+      xdrRead(x,dummy)
+  SHUT(str)
+  results
+
+generateDataName()==STRCONC($fortranTmpDir,getEnv('"HOST"),
+    getEnv('"SPADNUM"), GENTEMP('"NAG"),'"data")
+generateResultsName()==STRCONC($fortranTmpDir,getEnv('"HOST"),
+    getEnv('"SPADNUM"), GENTEMP('"NAG"),'"results")
+
+
+fortCall(objFile,data,results) ==
+  tmpFile1 := writeData(generateDataName(),data)
+  tmpFile2 := generateResultsName()
+  SYSTEM STRCONC(objFile," < ",tmpFile1," > ",tmpFile2)
+  results := readData(tmpFile2,results)
+  -- SYSTEM STRCONC("rm -f ",tmpFile1," ",tmpFile2)
+  PROBE_-FILE(tmpFile1) and DELETE_-FILE(tmpFile1)
+  PROBE_-FILE(tmpFile2) and DELETE_-FILE(tmpFile2)
+  results
+
+invokeNagman(objFiles,nfile,args,dummies,decls,results,actual) ==
+  actual := [spad2lisp(u) for u in first actual]
+  result := spadify(protectedNagCall(objFiles,nfile, _
+                 prepareData(args,dummies,actual,decls),_
+                 prepareResults(results,args,dummies,actual,decls)),_
+                 results,decls,inFirstNotSecond(args,dummies),actual)
+  -- Tidy up asps
+  -- if objFiles then SYSTEM STRCONC("rm -f ",addSpaces objFiles)
+  for fn in objFiles repeat PROBE_-FILE(fn) and DELETE_-FILE(fn)
+  result
+
+
+nagCall(objFiles,nfile,data,results,tmpFiled,tmpFiler) ==
+  nagMessagesString :=
+     $nagMessages => '"on"
+     '"off"
+  writeData(tmpFiled,data)
+  toSend:=STRCONC($nagHost," ",nfile," ",tmpFiler," ",tmpFiled," ",_
+      STRINGIMAGE($fortPersistence)," ", nagMessagesString," ",addSpaces objFiles)
+  sockSendString(8,toSend)
+  if sockGetInt(8)=1 then
+    results := readData(tmpFiler,results)
+  else
+    error ['"An error was detected while reading data: ", _
+           '"perhaps an incorrect array index was given ?"]
+  results
+
+protectedNagCall(objFiles,nfile,data,results) ==
+ errors :=true
+ val:=NIL
+ td:=generateDataName()
+ tr:=generateResultsName()
+ UNWIND_-PROTECT( (val:=nagCall(objFiles,nfile,data,results,td,tr) ;errors :=NIL),
+	errors =>( resetStackLimits(); sendNagmanErrorSignal();cleanUpAfterNagman(td,tr,objFiles)))
+ val
+
+
+cleanUpAfterNagman(f1,f2,listf)==
+  PROBE_-FILE(f1) and DELETE_-FILE(f1)
+  PROBE_-FILE(f2) and DELETE_-FILE(f2)
+  for fn in listf repeat PROBE_-FILE(fn) and DELETE_-FILE(fn)
+
+sendNagmanErrorSignal()==
+-- excite nagman's signal handler!
+ sockSendSignal(8,15)
+
+
+-- Globals
+-- $fortranDirectory := nil
+-- $fortranLibraries := '"-L/usr/local/lib/f90 -lf90 -L/usr/local/lib -lnag -lm"
+-- $fortranTmpDir := '"/tmp/"
+-- $addUnderscoreToFortranNames := true
+-- $fortranCompilerName := '"f90"
+
+inFirstNotSecond(f,s)==
+ [i for i in f | not i in s]
+
+-- Code for use in the Windows version of the AXIOM/NAG interface.
+
+multiToUnivariate f ==
+  -- Take an AnonymousFunction, replace the bound variables by references to
+  -- elements of a vector, and compile it.
+  (first f) ^= "+->" => error "in multiToUnivariate: not an AnonymousFunction"
+  if PAIRP CADR f then
+    vars := CDADR f -- throw away 'Tuple at start of variable list
+  else
+    vars := [CADR f]
+  body := COPY_-TREE CADDR f
+  newVariable := GENSYM()
+  for index in 0..#vars-1 repeat
+    -- Remember that AXIOM lists, vectors etc are indexed from 1
+    body := NSUBST(["elt",newVariable,index+1],vars.(index),body)
+  -- We want a Vector DoubleFloat -> DoubleFloat
+  target := [["DoubleFloat"],["Vector",["DoubleFloat"]]]
+  rest interpret ["ADEF",[newVariable],target,[[],[]],body]
+
+functionAndJacobian f ==
+  -- Take a mapping into n functions of n variables, produce code which will
+  -- evaluate function and jacobian values.
+  (first f) ^= "+->" => error "in functionAndJacobian: not an AnonymousFunction"
+  if PAIRP CADR f then
+    vars := CDADR f -- throw away 'Tuple at start of variable list
+  else
+    vars := [CADR f]
+  #(vars) ^= #(CDADDR f) => 
+    error "number of variables should equal number of functions"
+  funBodies := COPY_-TREE CDADDR f
+  jacBodies := [:[DF(f,v) for v in vars] for f in funBodies] where
+    DF(fn,var) == 
+      ["@",["convert",["differentiate",fn,var]],"InputForm"]
+  jacBodies := CDDR interpret [["$elt",["List",["InputForm"]],"construct"],:jacBodies]
+  newVariable := GENSYM()
+  for index in 0..#vars-1 repeat
+    -- Remember that AXIOM lists, vectors etc are indexed from 1
+    funBodies := NSUBST(["elt",newVariable,index+1],vars.(index),funBodies)
+    jacBodies := NSUBST(["elt",newVariable,index+1],vars.(index),jacBodies)
+  target := [["Vector",["DoubleFloat"]],["Vector",["DoubleFloat"]],["Integer"]]
+  rest interpret
+    ["ADEF",[newVariable,"flag"],target,[[],[],[]],_
+            ["IF", ["=","flag",1],_
+                   ["vector",["construct",:funBodies]],_
+                   ["vector",["construct",:jacBodies]]]]
+
+
+vectorOfFunctions f ==
+  -- Take a mapping into n functions of m variables, produce code which will
+  -- evaluate function values.
+  (first f) ^= "+->" => error "in vectorOfFunctions: not an AnonymousFunction"
+  if PAIRP CADR f then
+    vars := CDADR f -- throw away 'Tuple at start of variable list
+  else
+    vars := [CADR f]
+  funBodies := COPY_-TREE CDADDR f
+  newVariable := GENSYM()
+  for index in 0..#vars-1 repeat
+    -- Remember that AXIOM lists, vectors etc are indexed from 1
+    funBodies := NSUBST(["elt",newVariable,index+1],vars.(index),funBodies)
+  target := [["Vector",["DoubleFloat"]],["Vector",["DoubleFloat"]]]
+  rest interpret ["ADEF",[newVariable],target,[[],[]],["vector",["construct",:funBodies]]]
+
+
+
+
+
+@
+\eject
+\begin{thebibliography}{99}
+\bibitem{1} nothing
+\end{thebibliography}
+\end{document}