AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
AbstractLinAlgPack_MA28_CppDecl.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
5 // Copyright (2003) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 //
42 // C declarations for MA28 functions. These declarations should not have to change
43 // for different platforms. As long as the fortran object code uses capitalized
44 // names for its identifers then the declarations in Teuchos_F77_wrappers.h should be
45 // sufficent for portability.
46 
47 #ifndef MA28_CPPDECL_H
48 #define MA28_CPPDECL_H
49 
50 #include "Teuchos_F77_wrappers.h"
51 
52 namespace MA28_CppDecl {
53 
54 // Declarations that will link to the fortran object file.
55 // These may change for different platforms
56 
57 using FortranTypes::f_int; // INTEGER
58 using FortranTypes::f_real; // REAL
59 using FortranTypes::f_dbl_prec; // DOUBLE PRECISION
60 using FortranTypes::f_logical; // LOGICAL
61 
62 namespace Fortran {
63 extern "C" {
64 
65 // analyze and factorize a matrix
66 FORTRAN_FUNC_DECL_UL(void,MA28AD,ma28ad) (const f_int& n, const f_int& nz, f_dbl_prec a[], const f_int& licn
67  , f_int irn[], const f_int& lirn, f_int icn[], const f_dbl_prec& u, f_int ikeep[], f_int iw[]
68  , f_dbl_prec w[], f_int* iflag);
69 
70 // factor using previous analyze
71 FORTRAN_FUNC_DECL_UL(void,MA28BD,ma28bd) (const f_int& n, const f_int& nz, f_dbl_prec a[], const f_int& licn
72  , const f_int ivect[], const f_int jvect[], const f_int icn[], const f_int ikeep[], f_int iw[]
73  , f_dbl_prec w[], f_int* iflag);
74 
75 // solve for rhs using internally stored factorized matrix
76 FORTRAN_FUNC_DECL_UL(void,MA28CD,ma28cd) (const f_int& n, const f_dbl_prec a[], const f_int& licn, const f_int icn[]
77  , const f_int ikeep[], f_dbl_prec rhs[], f_dbl_prec w[], const f_int& mtype);
78 
79 // /////////////////////////////////////////////////////////////////////////////////////////
80 // Declare structs that represent the MA28 common blocks.
81 // These are the common block variables that are ment to be accessed by the user
82 // Some are used to set the options of MA28 and others return information
83 // about the attempts to solve the system.
84 // I want to provide the access functions that allow all of those common block
85 // variables that are ment to be accessed by the user to be accessable.
86 // For each of the common data items there will be a get operation that
87 // returns the variable value. For those items that are ment to be
88 // set by the user there will also be set operations.
89 
90 // COMMON /MA28ED/ LP, MP, LBLOCK, GROW
91 // INTEGER LP, MP
92 // LOGICAL LBLOCK, GROW
93 struct MA28ED_struct {
94  f_int lp;
95  f_int mp;
96  f_logical lblock;
97  f_logical grow;
98 };
99 extern MA28ED_struct FORTRAN_NAME_UL(MA28ED,ma28ed); // link to fortan common block
100 
101 // COMMON /MA28FD/ EPS, RMIN, RESID, IRNCP, ICNCP, MINIRN, MINICN,
102 // * IRANK, ABORT1, ABORT2
103 // INTEGER IRNCP, ICNCP, MINIRN, MINICN, IRANK
104 // LOGICAL ABORT1, ABORT2
105 // REAL EPS, RMIN, RESID
106 struct MA28FD_struct {
107  f_dbl_prec eps;
108  f_dbl_prec rmin;
109  f_dbl_prec resid;
110  f_int irncp;
111  f_int icncp;
112  f_int minirn;
113  f_int minicn;
114  f_int irank;
115  f_logical abort1;
116  f_logical abort2;
117 };
118 extern MA28FD_struct FORTRAN_NAME_UL(MA28FD,ma28fd); // link to fortan common block
119 
120 
121 // COMMON /MA28GD/ IDISP
122 // INTEGER IDISP
123 struct MA28GD_struct {
124  f_int idisp[2];
125 };
126 extern MA28GD_struct FORTRAN_NAME_UL(MA28GD,ma28gd); // link to fortan common block
127 
128 // COMMON /MA28HD/ TOL, THEMAX, BIG, DXMAX, ERRMAX, DRES, CGCE,
129 // * NDROP, MAXIT, NOITER, NSRCH, ISTART, LBIG
130 // INTEGER NDROP, MAXIT, NOITER, NSRCH, ISTART
131 // LOGICAL LBIG
132 // REAL TOL, THEMAX, BIG, DXMAX, ERRMAX, DRES, CGCE
133 struct MA28HD_struct {
134  f_dbl_prec tol;
135  f_dbl_prec themax;
136  f_dbl_prec big;
137  f_dbl_prec dxmax;
138  f_dbl_prec errmax;
139  f_dbl_prec dres;
140  f_dbl_prec cgce;
141  f_int ndrop;
142  f_int maxit;
143  f_int noiter;
144  f_int nsrch;
145  f_int istart;
146  f_logical lbig;
147 };
148 extern MA28HD_struct FORTRAN_NAME_UL(MA28HD,ma28hd); // link to fortan common block
149 
150 // COMMON /MA30ED/ LP, ABORT1, ABORT2, ABORT3
151 // INTEGER LP
152 // LOGICAL ABORT1, ABORT2, ABORT3
153 struct MA30ED_struct {
154  f_int lp;
155  f_logical abort1;
156  f_logical abort2;
157  f_logical abort3;
158 };
159 extern MA30ED_struct FORTRAN_NAME_UL(MA30ED,ma30ed); // link to fortan common block
160 
161 // COMMON /MA30FD/ IRNCP, ICNCP, IRANK, IRN, ICN
162 // INTEGER IRNCP, ICNCP, IRANK, IRN, ICN
163 struct MA30FD_struct {
164  f_int irncp;
165  f_int icncp;
166  f_int irank;
167  f_int minirn;
168  f_int minicn;
169 };
170 extern MA30FD_struct FORTRAN_NAME_UL(MA30FD,ma30fd); // link to fortan common block
171 
172 // COMMON /MA30GD/ EPS, RMIN
173 // DOUBLE PRECISION EPS, RMIN
174 struct MA30GD_struct {
175  f_dbl_prec eps;
176  f_dbl_prec rmin;
177 };
178 extern MA30GD_struct FORTRAN_NAME_UL(MA30GD,ma30gd); // link to fortan common block
179 
180 // COMMON /MA30HD/ RESID
181 // DOUBLE PRECISION RESID
182 struct MA30HD_struct {
183  f_dbl_prec resid;
184 };
185 extern MA30HD_struct FORTRAN_NAME_UL(MA30HD,ma30hd); // link to fortan common block
186 
187 // COMMON /MA30ID/ TOL, BIG, NDROP, NSRCH, LBIG
188 // INTEGER NDROP, NSRCH
189 // LOGICAL LBIG
190 // DOUBLE PRECISION TOL, BIG
191 struct MA30ID_struct {
192  f_dbl_prec tol;
193  f_dbl_prec big;
194  f_int ndrop;
195  f_int nsrch;
196  f_logical lbig;
197 };
198 extern MA30ID_struct FORTRAN_NAME_UL(MA30ID,ma30id); // link to fortan common block
199 
200 // COMMON /MC23BD/ LP,NUMNZ,NUM,LARGE,ABORT
201 // INTEGER LP, NUMNZ, NUM, LARGE
202 // LOGICAL ABORT
203 struct MC23BD_struct {
204  f_int lp;
205  f_int numnz;
206  f_int num;
207  f_int large;
208  f_logical abort;
209 };
210 extern MC23BD_struct FORTRAN_NAME_UL(MC23BD,mc23bd); // link to fortan common block
211 
212 
213 } // end extern "C"
214 } // end namespace Fortran
215 
216 
217 /* * @name {\bf MA28 C++ Declarations}.
218  *
219  * These the C++ declarations for MA28 functions and common block data.
220  * These declarations will not change for different platforms.
221  * All of these functions are in the C++ namespace #MA28_CppDecl#.
222  *
223  * These functions perform the three phases that are normally associated
224  * with solving sparse systems of linear equations; analyze (and factorize)
225  * , factorize, and solve. The MA28 interface uses a coordinate format
226  * (aij, i, j) for the sparse matrix.
227  *
228  * There are three interface routienes that perform these steps:
229  * \begin{description}
230  * \item[#ma28ad#] Analyzes and factorizes a sparse matrix stored in #a#, #irn#
231  * , #icn# and returns the factorized matrix data structures in #a#, #icn#
232  * , and #ikeep#.
233  * \item[#ma28bd#] Factorizes a matrix with the same sparsity structure that was
234  * previously factorized by #ma28ad#. Information about the row and column
235  * permutations, fill-in elements etc. from the previous analyze and factorization
236  * is passed in the arguments #icn#, and #ikeep#. The matrix to be
237  * factorized is passed in #a#, #ivect#, and #jvect# and the non-zero
238  * elements of the factorization are returned in #a#.
239  * \item[#ma28cd#] Solves for a dense right hand side (rhs) given a matrix factorized
240  * by #ma28ad# or #ma28bd#. The rhs is passed in #rhs# and the solution
241  * is returned in #rhs#. The factorized matrix is passed in by #a#, #icn#
242  * and #ikeep#. The transposed or the non-transposed system can be solved
243  * for by passing in #mtype != 1# and #mtype == 1# respectively.
244  */
245 
246 // @{
247 // begin MA28 C++ Declarations
248 
249 /* * @name {\bf MA28 / MA30 Common Block Access}.
250  *
251  * These are references to structures that allow C++ users to set and retrive
252  * values of the MA28xD and MA30xD common blocks. Some of the common block
253  * items listed below for MA28xD are also present in MA30xD. The control
254  * parameters (abort1, eps, etc.) for MA28xD are transfered to the equivalent
255  * common block variables in the #ma28ad# function but not in any of the other
256  * functions.
257  *
258  * The internal states for MA28, MA30, and MC23 are determined by the
259  * values in these common block variables as there are no #SAVE# variables
260  * in any of the functions. So to use MA28 with more than
261  * one sparse matrix at a time you just have to keep copies of these
262  * common block variable for each system and then set them when every
263  * you want to work with that system agian. This is very good news.
264  *
265  * These common block variables are:
266  * \begin{description}
267  * \item[lp, mp]
268  * Integer: Used by the subroutine as the unit numbers for its warning
269  * and diagnostic messages. Default value for both is 6 (for line
270  * printer output). the user can either reset them to a different
271  * stream number or suppress the output by setting them to zero.
272  * While #lp# directs the output of error diagnostics from the
273  * principal subroutines and internally called subroutines, #mp#
274  * controls only the output of a message which warns the user that he
275  * has input two or more non-zeros a(i), . . ,a(k) with the same row
276  * and column indices. The action taken in this case is to proceed
277  * using a numerical value of a(i)+...+a(k). in the absence of other
278  * errors, #iflag# will equal -14 on exit.
279  * \item[lblock]
280  * Logical: Controls an option of first
281  * preordering the matrix to block lower triangular form (using
282  * harwell subroutine mc23a). The preordering is performed if #lblock#
283  * is equal to its default value of #true# if #lblock# is set to
284  * #false# , the option is not invoked and the space allocated to
285  * #ikeep# can be reduced to 4*n+1.
286  * \item[grow]
287  * Logical: If it is left at its default value of
288  * #true# , then on return from ma28a/ad or ma28b/bd, w(1) will give
289  * an estimate (an upper bound) of the increase in size of elements
290  * encountered during the decomposition. If the matrix is well
291  * scaled, then a high value for w(1), relative to the largest entry
292  * in the input matrix, indicates that the LU decomposition may be
293  * inaccurate and the user should be wary of his results and perhaps
294  * increase u for subsequent runs. We would like to emphasise that
295  * this value only relates to the accuracy of our LU decomposition
296  * and gives no indication as to the singularity of the matrix or the
297  * accuracy of the solution. This upper bound can be a significant
298  * overestimate particularly if the matrix is badly scaled. If an
299  * accurate value for the growth is required, #lbig# (q.v.) should be
300  * set to #true#
301  * \item[eps, rmin]
302  * Double Precision: If on entry to ma28b/bd, #eps# is less
303  * than one, then #rmin# will give the smallest ratio of the pivot to
304  * the largest element in the corresponding row of the upper
305  * triangular factor thus monitoring the stability of successive
306  * factorizations. if rmin becomes very large and w(1) from
307  * ma28b/bd is also very large, it may be advisable to perform a
308  * new decomposition using ma28a/ad.
309  * \item[resid]
310  * Double Precision: On exit from ma28c/cd gives the value
311  * of the maximum residual over all the equations unsatisfied because
312  * of dependency (zero pivots).
313  * \item[irncp,icncp]
314  * Integer: Monitors the adequacy of "elbow
315  * room" in #irn# and #a#/#icn# respectively. If either is quite large (say
316  * greater than n/10), it will probably pay to increase the size of
317  * the corresponding array for subsequent runs. if either is very low
318  * or zero then one can perhaps save storage by reducing the size of
319  * the corresponding array.
320  * \item[minirn, minicn]
321  * Integer: In the event of a
322  * successful return (#iflag# >= 0 or #iflag# = -14) give the minimum size
323  * of #irn# and #a#/#icn# respectively which would enable a successful run
324  * on an identical matrix. On an exit with #iflag# equal to -5, #minicn#
325  * gives the minimum value of #icn# for success on subsequent runs on
326  * an identical matrix. in the event of failure with #iflag# = -6, -4,
327  * -3, -2, or -1, then #minicn# and #minirn# give the minimum value of
328  * #licn# and #lirn# respectively which would be required for a
329  * successful decomposition up to the point at which the failure
330  * occurred.
331  * \item[irank]
332  * Integer: Gives an upper bound on the rank of the matrix.
333  * \item[abort1]
334  * Logical: Default value #true#. If #abort1# is
335  * set to #false# then ma28a/ad will decompose structurally singular
336  * matrices (including rectangular ones).
337  * \item[abort2]
338  * Logical: Default value #true#. If #abort2# is
339  * set to #false# then ma28a/ad will decompose numerically singular
340  * matrices.
341  * \item[idisp]
342  * Integer[2]: On output from ma28a/ad, the
343  * indices of the diagonal blocks of the factors lie in positions
344  * idisp(1) to idisp(2) of #a#/#icn#. This array must be preserved
345  * between a call to ma28a/ad and subsequent calls to ma28b/bd,
346  * ma28c/cd or ma28i/id.
347  * \item[tol]
348  * Double Precision: If it is set to a positive value, then any
349  * non-zero whose modulus is less than #tol# will be dropped from the
350  * factorization. The factorization will then require less storage
351  * but will be inaccurate. After a run of ma28a/ad with #tol# positive
352  * it is not possible to use ma28b/bd and the user is recommended to
353  * use ma28i/id to obtain the solution. The default value for #tol# is
354  * 0.0.
355  * \item[themax]
356  * Double Precision: On exit from ma28a/ad, it will hold the
357  * largest entry of the original matrix.
358  * \item[big]
359  * Double Precision: If #lbig# has been set to #true#, #big# will hold
360  * the largest entry encountered during the factorization by ma28a/ad
361  * or ma28b/bd.
362  * \item[dxmax]
363  * Double Precision: On exit from ma28i/id, #dxmax# will be set to
364  * the largest component of the solution.
365  * \item[errmax]
366  * Double Precision: On exit from ma28i/id, If #maxit# is
367  * positive, #errmax# will be set to the largest component in the
368  * estimate of the error.
369  * \item[dres]
370  * Double Precision: On exit from ma28i/id, if #maxit# is positive,
371  * #dres# will be set to the largest component of the residual.
372  * \item[cgce]
373  * Double Precision: It is used by ma28i/id to check the
374  * convergence rate. if the ratio of successive corrections is
375  * not less than #cgce# then we terminate since the convergence
376  * rate is adjudged too slow.
377  * \item[ndrop]
378  * Integer: If #tol# has been set positive, on exit
379  * from ma28a/ad, #ndrop# will hold the number of entries dropped from
380  * the data structure.
381  * \item[maxit]
382  * Integer: It is the maximum number of iterations
383  * performed by ma28i/id. It has a default value of 16.
384  * \item[noiter]
385  * Integer: It is set by ma28i/id to the number of
386  * iterative refinement iterations actually used.
387  * \item[nsrch]
388  * Integer: If #nsrch# is set to a value less than #n#,
389  * then a different pivot option will be employed by ma28a/ad. This
390  * may result in different fill-in and execution time for ma28a/ad.
391  * If #nsrch# is less than or equal to #n#, the workspace array #iw# can be
392  * reduced in length. The default value for nsrch is 32768.
393  * \item[istart]
394  * Integer: If #istart# is set to a value other than
395  * zero, then the user must supply an estimate of the solution to
396  * ma28i/id. The default value for istart is zero.
397  * \item[lbig]
398  * Logical: If #lbig# is set to #true#, the value of the
399  * largest element encountered in the factorization by ma28a/ad or
400  * ma28b/bd is returned in #big#. setting #lbig# to #true# will
401  * increase the time for ma28a/ad marginally and that for ma28b/bd
402  * by about 20%. The default value for #lbig# is #false#.
403  */
404 
405 // @{
406 // begin MA28 Common Block Access
407 
408 // / Common block with members: #lp#, #mp#, #lblock#, #grow#
409 static MA28ED_struct &ma28ed_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA28ED,ma28ed);
410 // / Common block with members: #eps#, #rmin#, #resid#, #irncp#, #icncp#, #minirc#, #minicn#, #irank#, #abort1#, #abort2#
411 static MA28FD_struct &ma28fd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA28FD,ma28fd);
412 // / Common block with members: #idisp#
413 static MA28GD_struct &ma28gd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA28GD,ma28gd);
414 // / Common block with members: #tol#, #themax#, #big#, #bxmax#, #errmax#, #dres#, #cgce#, #ndrop#, #maxit#, #noiter#, #nsrch#, #istart#, #lbig#
415 static MA28HD_struct &ma28hd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA28HD,ma28hd);
416 // / Common block with members: #lp#, #abort1#, #abort2#, #abort3#
417 static MA30ED_struct &ma30ed_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA30ED,ma30ed);
418 // / Common block with members: #irncp#, #icncp#, #irank#, #irn#, #icn#
419 static MA30FD_struct &ma30fd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA30FD,ma30fd);
420 // / Common block with members: #eps#, #rmin#
421 static MA30GD_struct &ma30gd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA30GD,ma30gd);
422 // / Common block with members: #resid#
423 static MA30HD_struct &ma30hd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA30HD,ma30hd);
424 // / Common block with members: #tol#, #big#, #ndrop#, #nsrch#, #lbig#
425 static MA30ID_struct &ma30id_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MA30ID,ma30id);
426 // / Common block with members: #lp#, #numnz#, #num#, #large#, #abort#
427 static MC23BD_struct &mc23bd_cb = FORTRAN_COMMMON_BLOCK_NAME_UL(MC23BD,mc23bd);
428 
429  // The reason that these are declared static is because I need to
430  // make sure that these references are initialized before they are
431  // used in other global defintions. This means that every translation
432  // unit will have their own copy of this data. To reduce this code
433  // blot you could declair them as pointers then set then using the
434  // trick of an initialization class (Myers).
435 
436 // end MA28 Common Block Access
437 // @}
438 
439 // /
440 /* * Analyze and factor a sparse matrix.
441  *
442  * This function analyzes (determines row and column pivots to minimize
443  * fill-in and result in a better conditioned factorization) then factors
444  * (calculates upper and lower triangular factors for the determined
445  * row and column pivots) the permuted system. The function takes a sparse
446  * matrix stored in coordinate form ([Aij, i, j] => [#a#, #irn#, #icn#]) and
447  * factorizes it. On entry, the first #nz# elements of #a#, #irn#, and
448  * #icn# hold the matrix elements. The remaining entires of #a#, #irn#
449  * , and #icn# hold the fill-in entries after the matrix factorization is
450  * complete.
451  *
452  * The amount of fill-in is influenced by #u#. A value of #u = 1.0# gives partial
453  * pivoting where fill-in is sacrificed for the sake of numerical stability and
454  * #u = 0.0# gives pivoting that strictly minimizes fill-in.
455  *
456  * The parameters #ikeep#, #iw#, and #w# are used as workspace but
457  * #ikeep# contains important information about the factoriation
458  * on return.
459  *
460  * The parameter #iflag# is used return error information about the attempt
461  * to factorized the system.
462  *
463  * @param n [input] Order of the system being solved.
464  * @param nz [input] Number of non-zeros of the input matrix (#nz >= n#).
465  * The ratio #nz/(n*n)# equals the sparsity fraction for the matrix.
466  * @param a [input/output] Length = #licn#. The first #nz# entries hold the
467  * non-zero entries of the input matrix on input and
468  * the non-zero entries of the factorized matrix on exit.
469  * @param licn [input] length of arrays #a# and #icn#. This
470  * is the total amount of storage advalable for the
471  * non-zero entries of the factorization of #a#.
472  * therefore #licn# must be greater than #nz#. How
473  * much greater depends on the amount of fill-in.
474  * @param irn [input/modifed] Length = #ircn#. The first #nz# entries hold
475  * the row indices of the matrix #a# on input.
476  * @param ircn [input] Lenght of irn.
477  * @param icn [input/output] Length = #licn#. Holds column indices of #a# on
478  * entry and the column indices of the reordered
479  * #a# on exit.
480  * @param u [input] Controls partial pivoting.
481  * \begin{description}
482  * \item[#* u >= 1.0#]
483  * Uses partial pivoting for maximum numerical stability
484  * at the expense of some extra fill-in.
485  * \item[#* 1.0 < u < 0.0#]
486  * Balances numerical stability and fill-in with #u# near 1.0
487  * favoring stability and #u# near 0.0 favoring less fill-in.
488  * \item[#* u <= 0.0#]
489  * Determines row and column pivots to minimize fill-in
490  * irrespective of the numerical stability of the
491  * resulting factorization.
492  * \end{description}
493  * @param ikeep [output] Length = 5 * #n#. On exist contains information about
494  * the factorization.
495  * \begin{description}
496  * \item[#* #ikeep(:,1)]
497  * Holds the total length of the part of row i in
498  * the diagonal block.
499  * \item[#* #ikeep(:,2)]
500  * Holds the row pivots. Row #ikeep(i,2)# of the
501  * input matrix is the ith row of the pivoted matrix
502  * which is factorized.
503  * \item[#* #ikeep(:,3)]
504  * Holds the column pivots. Column #ikeep(i,3)# of the
505  * input matrix is the ith column of the pivoted matrix
506  * which is factorized.
507  * \item[#* #ikeep(:,4)]
508  * Holds the length of the part of row i in the L
509  * part of the L/U decomposition.
510  * \item[#* #ikeep(:,5)]
511  * Holds the length of the part of row i in the
512  * off-diagonal blocks. If there is only one
513  * diagonal block, #ikeep(i,5)# is set to -1.
514  * \end{description}
515  * @param iw [] Length = #8*n#. Integer workspace.
516  * @param w [] Length = #n#. Real workspace.
517  * @param iflag [output] Used to return error condtions.
518  * \begin{description}
519  * \item[#* >= 0#] Success.
520  * \item[#* < 0#] Some error has occured.
521  * \end{description}
522  */
523 inline void ma28ad(const f_int& n, const f_int& nz, f_dbl_prec a[], const f_int& licn
524  , f_int irn[], const f_int& lirn, f_int icn[], const f_dbl_prec& u, f_int ikeep[], f_int iw[]
525  , f_dbl_prec w[], f_int* iflag)
526 { Fortran::FORTRAN_FUNC_CALL_UL(MA28AD,ma28ad) (n,nz,a,licn,irn,lirn,icn,u,ikeep,iw,w,iflag); }
527 
528 // /
529 /* * Factor a sparse matrix using previous analyze pivots.
530  *
531  * This function uses the pivots determined from a previous factorization
532  * to factorize the matrix #a# again. The assumption is that the
533  * sparsity structure of #a# has not changed but only its numerical
534  * values. It is therefore possible that the refactorization may
535  * become unstable.
536  *
537  * The matrix to be refactorized on order #n# with #nz# non-zero elements
538  * is input in coordinate format in #a#, #ivect# and #jvect#.
539  *
540  * Information about the factorization is contained in the #icn# and
541  * #ikeep# arrays returned from \Ref{ma28ad}.
542  *
543  * @param n [input] Order of the system being solved.
544  * @param nz [input] Number of non-zeros of the input matrix
545  * @param a [input/output] Length = #licn#. The first #nz# entries hold the
546  * non-zero entries of the input matrix on input and
547  * the non-zero entries of the factorized matrix on exit.
548  * @param licn [input] length of arrays #a# and #icn#. This
549  * is the total amount of storage avalable for the
550  * non-zero entries of the factorization of #a#.
551  * therefore #licn# must be greater than #nz#. How
552  * much greater depends on the amount of fill-in.
553  * @param icn [input] Length = #licn#. Same array output from #ma28ad#.
554  * It contains information about the analyze step.
555  * @param ikeep [input] Length = 5 * #n#. Same array output form #ma28ad#.
556  * It contains information about the analyze step.
557  * @param iw [] Length = #8*n#. Integer workspace.
558  * @param w [] Length = #n#. Real workspace.
559  * @param iflag [output] Used to return error condtions.
560  * \begin{description}
561  * \item[#* >= 0#] Success
562  * \item[#* < 0#] Some error has occured.
563  * \end{description}
564  */
565 inline void ma28bd(const f_int& n, const f_int& nz, f_dbl_prec a[], const f_int& licn
566  , const f_int ivect[], const f_int jvect[], const f_int icn[], const f_int ikeep[], f_int iw[]
567  , f_dbl_prec w[], f_int* iflag)
568 { Fortran::FORTRAN_FUNC_CALL_UL(MA28BD,ma28bd) (n,nz,a,licn,ivect,jvect,icn,ikeep,iw,w,iflag); }
569 
570 // /
571 /* * Solve for a rhs using a factorized matrix.
572  *
573  * This function solves for a rhs given a matrix factorized by
574  * #ma28ad# or #ma28bd#. The right hand side (rhs) is passed
575  * in in #rhs# and the solution is return in #rhs#. The
576  * factorized matrix is passed in in #a#, #icn#, and #ikeep#
577  * which were set by #ma28ad# and/or #ma28bd#. The
578  *
579  * The matrix or its transpose can be solved for by selecting
580  * #mtype == 1# or #mtype != 1# respectively.
581  *
582  * @param n [input] Order of the system being solved.
583  * @param a [input] Length = #licn#. Contains the non-zero
584  * elements of the factorized matrix.
585  * @param licn [input] length of arrays #a# and #icn#. This
586  * is the total amount of storage avalable for the
587  * non-zero entries of the factorization of #a#.
588  * therefore #licn# must be greater than #nz#. How
589  * much greater depends on the amount of fill-in.
590  * @param icn [input] Length = #licn#. Same array output from #ma28ad#.
591  * It contains information about the analyze step.
592  * @param ikeep [input] Length = 5 * #n#. Same array output form #ma28ad#.
593  * It contains information about the analyze step.
594  * @param w [] Length = #n#. Real workspace.
595  * @param mtype [input] Instructs to solve using the matrix or its transpoze.
596  * \begin{description}
597  * \item[#* mtype == 1#] Solve using the non-transposed matrix.
598  * \item[#* mtype != 1#] Solve using the transposed matrix.
599  * \end{description}
600  */
601 inline void ma28cd(const f_int& n, const f_dbl_prec a[], const f_int& licn, const f_int icn[]
602  , const f_int ikeep[], f_dbl_prec rhs[], f_dbl_prec w[], const f_int& mtype)
603 { Fortran::FORTRAN_FUNC_CALL_UL(MA28CD,ma28cd) (n,a,licn,icn,ikeep,rhs,w,mtype); }
604 
605 // end MA28 C++ Declarations
606 // @}
607 
608 } // end namespace MA28_CDecl
609 
610 #endif // MA28_CPPDECL_H